A short Outline of the Method

In the following, an image is given as a mapping $I:\Omega \to V$ from its coordinate domain $\Omega \in {\mathbb{R}}^3$ to its intensity range $V\in {\mathbb{R}}$. Given a coordinate $\vec{x}\in \Omega$, and the intensity of the image $I$ at this coordinate $I(\vec{x})$, the ordered pair $(\vec{x}, I(\vec{x}))$ is referred to as a voxel (volume element). Using a transformation $T:\Omega \to \Omega$, an image can be changed according to $I_T := I(T(\vec{x}))$. The set of all these transformations is called the transformation space $\Gamma$.

In this paper, the transformations correspond to spatial displacements of voxels and are described in the so-called Eulerian reference frame. Here the voxels are tracked by their position: A voxel originates at time $t_0=0$ at coordinate $\vec{x}\in \Omega$. As it moves through $\Omega$, the displacement of a voxel $(\vec{x}, I(\vec{x}))$ at time t is given as a vector $\mathbf{u}(\vec{x},t)$. The set of the displacements of all voxels of an image is called a displacement field over domain $\Omega$, and its value at time $t$ is denoted as $\mathbf{u}(t)$. The corresponding transformation T can be given coordinate-wise:

$$T_{t}(\vec{x}):=\vec{x}- \mathbf{u}(\vec{x},t)\: \forall \: \vec{x}\in \Omega .$$ (1)

The concatenation of transformations is then given as

$$T_{1}\circ T_{2} := \vec{x}-\mathbf{u}_{1}(\vec{x}-\mathbf{u}_{2}(\vec{x}))-\mathbf{u}_{2}(\vec{x}),$$ (2)

The focus of the registration of one (study) image $S:\Omega \to V$ to another (reference) image $R:\Omega \to V$ is to find a transformation $T_{min} \in \Gamma$ that minimizes a given cost function $F(R,S_T)$ describing the similarity between transformed study image $S$ and reference image $R$ in conjunction with an energy normalization (smoothness) term $E(T)$ that enforces topology preservation:

$$T_{\min}:= \arg \min_{T \in \Gamma } \left( F(R,S_T) + \kappa E(T) \right).$$ (3)

$\kappa$ is a Lagrangian multiplier to balance between registration accuracy and transformation smoothness. Minimizing (3) can be done in terms of its first order derivative:

$$\kappa \frac{\partial}{\partial T}E(T) = - \frac{\partial}{\partial T} F(T,S,R).$$ (4)

In the non-rigid registration software I use the sum of squared differences as a cost function:

$$F_{c}(T,S,R):=\frac{1}{2}\int _{\Omega}\left[ R(\vec{x})-S(T(\vec{x}))\right] ^{2}d\vec{x},$$ (5)

and fluid dynamics as smoothness measure.

Then the first order derivative of the cost function (5) can be used to estimate a deforming force:

$$\mathbf{f}(\vec{x},t) := - [S(T(\vec{x},t))-R(\vec{x})]\left. \nabla S\right\vert _{T(\vec{x},t)},$$ (6)

and with $\kappa=1.0$, this force (6), and fluid dynamics energy regularisation, (4) can be written as

$$\left( \mu \nabla^2 + ( \mu + \lambda ) \nabla ( \nabla \cdot ) \right) \mathbf{v}(\vec{x}, t) = - \mathbf{f}(\vec{x},\mathbf{u}(\vec{x},t)).$$ (7)

In order to solve the registration problem, (7) is solved for constant time, and the deformation field $\mathbf{u}(t)$ is updated from the estimated velocity field using a time integration step with step-width $\Delta t$:

$$\mathbf{u}(\vec{x},t + \Delta t) := \mathbf{u}(\vec{x},t) + \Delt... ...mathbf{v}(\vec{x},t)- \nabla \mathbf{u}(\vec{x},t)\mathbf{v}(\vec{x},t)\right].$$ (8)

The solution of the registration problem is summarized in algorithm 1