-15 -10 -5 0 5 10 15 -15-10-5 0 5 10 15 0
0.2 0.4 0.6 0.8 1
-15 -10 -5 0 5 10 15 -15-10-5 0 5 10 15 0
0.2 0.4 0.6 0.8 1
-15 -10 -5 0 5 10 15 -15-10-5 0 5 10 15 0
0.2 0.4 0.6 0.8 1
(a)falloff = 10−4 (b)falloff = 10−3 (c) falloff = 10−2 Figure 8.1: Various weightings of the image area surrounding the search block. A larger falloff factors of an exponential distribution allows selecting between uniform con-stant weighting where every pixel are equally important and point weighting where only the center pixel contributes to the evaluation of local image dissimilarity.
To enforce smoothness, hence, regularity of the deformation, a weighted regular-ization term is added which includes the square of the Laplacians. This regulariza-tion scheme was already proposed in the seminal paper of Horn and Schunk.
The regularized cost function can be expressed by Q=X
(f −m)2+λaTLTLa, (8.10) whereLstands for a discretized Laplacian operator andλis a user-defined parameter driving the strength of the regularization term.
In the implementation of ShIRT, the registration problem is iteratively solved by a conjugate gradient (CG) method. The current estimate for the parameters a at iterationt≥1 is noted byat. The accuracy of the estimate is iteratively updated asat+1 =at+ ∆at with
∆at=
TTT +λLTL−1
TT(f −mt)−λLTLat
. (8.11)
For the present work, the block is round-shaped and the coordinates of the pixels Bi,∀i≤ |B|belonging to a compact circular support is noted by
B ≡n
(x, y)i : k(x, y)ik ≥ (x, y)j
, ∀i > jo
, (8.12)
where|B|is the number of pixels included in the setB. In a similar way, the search region in the target image is limited to a round neighborhood around the center of the block. This set of displacement vectors Nj,∀j ≤ |N| is also ordered according to magnitudes such that
N ≡n
(u, v)i: k(u, v)ik ≥ (u, v)j
, ∀i > jo
, (8.13)
where|N|is the number of pixels inside the search regionN.
If displacements are small, it can be assumed that the larger the magnitude of the displacement is, the less probable a block of better similarity will be found further away from the best match found so far. Therefore, a weighting is introduced within the image dissimilarity metric to control the influence of pixels far away from the block center.
The weights are associated to each pixel of the blockB and are precomputed in a look-up table (LUT) using the formula
wi = exp h
−falloff · kBik2i
, ∀i≤ |B|. (8.14) Furthermore, weights are normalized such that the integral over every pixel of the block equals to one. Thus, normalized weights are finally assigned with the following normalization procedure:
wi ← wi P
k≤|B|wk
, ∀i≤ |B|. (8.15)
The weights follow an exponential function dependent on the squared Euclidean distance from the search center and a constant falloff factor. The spatial distribution of function for different values offalloff in a square of 33×33 pixels can be visualized in figure 8.1. It can be seen that the number of pixels with a significant weight decreases with increasing falloff factors, while the shape of the weighting function is almost flat for a small constant.
For computing the dissimilarity between a block extracted from position (x, y) inside the target imagef and another block from the source imageg, the algorithm tests first if the search block covers image elements which are outside the image borders. Then, the dissimilarity integrated over pixels that belong to the image support is normalized according to the sum of the weightings of the pixels inside the image.
A simple computation of the image distance is used in this work, while multi-modal registration requires appropriate image similarity metrics such as the pop-ular mutual information (Pluim et al., 2000, 2003) or the cross correlation ( An-dronacheet al., 2008). The sum of absolute differences (SAD) is used for
measur-ing the dissimilarity such as Q= 1
W X
i≤|B|
wi|f((x, y) +Bi)−g(Nj+ (x, y) +Bi)|, (8.16) where every difference has been weighted. The normalization term is given by
W = X
i≤|B|
wi. (8.17)
Note that if the block is completely inside the image, thenW = 1 and the evaluation of the normalization term can be ignored.
Additionally, in order to favor small displacements, a shift-variant penalty term dependent on the squared length of displacement vectors is added to the dissimilarity metric. The penalization is weighted by the factorλ≥0 such that
Q= 1 W
X
i≤|B|
wi|f((x, y) +Bi)−g(Nj+ (x, y) +Bi)|+λkNjk2. (8.18)
The displacement vector Nj ∈N that minimizes the cost function Qis selected by an exhaustive search strategy for every control point located in (x, y). The rationale is that, since the number of candidate displacement vectors |N| is small, an exhaustive search can outperform more complex optimization strategies.
In comparison to the clasical and popular sum of squared differences (SSD) metric, the SAD metric emphasizes much less small relative differences between bright pixels. Many alternative criteria have been proposed for computing image dissimilarities, such as the relative entropy (Kullback-Liebler divergence) that is based on elements of information theory. Future work could evaluate the robustness of the proposed block matching scheme when using different metrics.
In the current implementation, the displacement vectors are estimated on a coarser grid. Using a smaller displacement field decreases the number of compu-tations that is needed and therefore speeds up the registration (Plishker et al., 2007). Additionally, the displacement field is scaled up to the full size using cubic interpolation and a smooth displacement vector field is obtained. Since the motion is computed at a coarser resolution, small local deformations can not be captured (Robinson and Milanfar, 2004).
Computing resources are especially needed for computing the dissimilarity be-tween blocks for each possible displacement vector candidate. The number of blocks to be processed can be reduced further by previously thresholding the image and determining which pixels are close to object borders (Malschet al., 2006). Another possibility to increase the speed of the algorithm is to alter the search strategy (Zhu and Ma, 2000).
A full search algorithm is used in standard block matching, that means that in a predefined search area every possible displacement candidate is evaluated. In addition, multiresolution image representation (Adelson et al., 1984) can both speed-up the registration process and regularize the results (Bajcsyand Kovacic, 1989; Kostelecet al., 1998).
(a) Acquired projection (b) Selected reference (c) Approx. reference Figure 8.2: Two series of pairs of projections are used for experiments but only the frontal projection is shown as illustration. The sharp acquired projection (a) is registered on either a sharp projection (b) from a selected reference motion state or a more blurry approximate reference projection (c), projected from a scout reconstruction. The approximate projection contains motion blur artifacts that are well visible at the border of the diaphragm.
(a) ShIRT (b) Block matching λ= 0
(c) Block matching λ= 100
Figure 8.3: Effect of the penalization parameterλon the regularity of the grid produced by block matching. In comparison to the grid produced by ShIRT (a), the non-penalized block matching scheme produces irregular deformations (b). However, irregularities are largely prevented by using a slight penalization (c) within the image dissimilarity metric.