Image Registration

Image registration is the task of bringing two images into the same coordinate frame: finding the transformation that maps each point of one image onto the point of the other that shows the same part of the scene or the same anatomy. The images may have been taken at different times, from different viewpoints, or by different sensors (Zitová and Flusser, 2003). Registration comes before many comparisons and combinations of images: stitching photographs into a panorama, stabilizing video, detecting change between satellite passes, fusing a CT scan with an MRI scan, and tracking a template through a sequence.

A registration method makes three choices: a transform model that says how the images may differ geometrically, a similarity measure that says how well they agree under a candidate transform, and a search that finds the best transform. Feature-based methods match distinctive points and fit the transform to the matches; intensity-based (direct) methods compare pixel values and optimize the transform directly.

Problem Definition

One image is held fixed, the reference or fixed image $I_R$; the other, the moving, sensed or floating image $I_M$, is transformed onto it. Registration seeks a transformation $T$ from a family $\mathcal{T}$ such that $I_M(T(\mathbf{x}))$ matches $I_R(\mathbf{x})$ for every point $\mathbf{x}$ of the overlap. Intensity-based methods write this as an optimization,

$$\hat{T} = \arg\min_{T \in \mathcal{T}} \; D\big(I_R,\; I_M \circ T\big) + \lambda\, \mathcal{R}(T),$$

where $D$ is a dissimilarity measure (or the negative of a similarity), $\mathcal{R}$ penalizes implausible transforms, and $\lambda$ balances the two. Sotiras, Davatzikos and Paragios (2013) use this form for deformable registration, with a deformation model, an objective function and an optimization method as the three components of every algorithm. For a few-parameter transform the regularizer is often dropped; for a deformable one it is essential, because a displacement vector per pixel has more unknowns than the image has intensities to constrain them.

Zitová and Flusser divide most registration methods into four steps: detect distinctive features, match them between the images, estimate the transform from the matches, and resample the moving image with it (image warping). Intensity-based (area-based) methods detect no features, so the first step is omitted and matching is done by the similarity measure.

Registration is ill-posed in general. Sotiras et al. give a simple example: two images of a uniform disk are aligned by any rotation about its center, so even the rigid transform is not unique. Where the images carry no structure, nothing can be estimated, and where the scene itself changed (a tumor grew, a building was demolished), registration must not remove the very difference the user is looking for.

Input and Output

• Input: two images (2D or 3D), a transform model, and often an initial guess, a region of interest or a mask of pixels to ignore; for feature-based methods, detected keypoints and descriptors.
• Output: the transform parameters (a matrix for parametric models, a displacement field for deformable ones), usually with the warped moving image, and ideally an estimate of the registration error.

Variants

Transform models

Planar parametric transforms form a nested family, written as $3 \times 3$ matrices in homogeneous coordinates in the image warping article:

Model	Degrees of freedom (2D)	Preserves	Typical cause
Translation	2	everything but position	camera panning slightly, scanner offset
Rigid (Euclidean)	3	lengths and angles	rotated document, same-subject medical scans
Similarity	4	angles, length ratios	zoom plus rotation
Affine	6	parallel lines, area ratios	distant scene, oblique view of a small region
Homography (projective)	8	straight lines	plane seen from two viewpoints, rotating camera
Deformable (non-rigid)	2 per pixel, or per control point	smoothness, if regularized	breathing, soft tissue, lens and terrain effects

In 3D, a rigid transform has 6 degrees of freedom and an affine one 12; Maes et al. (1997) registered brain volumes with the 6-parameter rigid model, and Viola and Wells (1997) aligned 2D MR images with a 6-parameter affine one. The model should be the least flexible one that explains the geometry: extra freedom fits noise and hides real differences. Zitová and Flusser also distinguish global models, with one set of parameters for the whole image, from local models fitted patch by patch, for example piecewise over a triangulation; radial basis functions such as thin-plate splines are global models that can still handle locally varying distortions.

Other variants

• Monomodal and multimodal registration. In monomodal registration both images come from the same kind of sensor, and corresponding pixels have similar intensities. In multimodal registration (CT and MRI, visible and thermal, optical and radar) the same structure can be bright in one image and dark in the other, so the similarity measure must not assume equal intensities.
• Pairwise, groupwise and atlas registration. Many images can be registered to one reference, to a population average (an atlas), or to each other jointly.
• Image-to-model registration. Viola and Wells (1997) aligned a 3D surface model of an object to a video image of it, and model-based alignment of this kind is a registration problem whose transform is a 3D pose.
• Point-set registration. When the data are points rather than images, as in LiDAR scans, alignment is done by Iterative Closest Point and related point-cloud registration methods.

Medical image registration builds on the same methods and adds scanner-specific geometry, anatomical priors and clinical validation.

Challenges

• Large motions. Intensity-based objectives are non-convex, so a local optimizer started too far from the solution converges to the wrong alignment.
• Different appearance. Illumination, exposure, season, and especially a different sensor change the intensities that a measure compares.
• Partial overlap and occlusion. Only part of each image has a counterpart, and moving objects or parallax violate a single global transform.
• Weak or repetitive structure. Flat regions constrain nothing; repeated patterns produce several plausible alignments.
• Real change versus misalignment. It is hard to tell a registration error from a genuine change in the scene (Zitová and Flusser, 2003).
• Cost. Deformable 3D registration has millions of unknowns; VoxelMorph’s authors cite tens of minutes to hours per pair for accurate CPU methods (Balakrishnan et al., 2019).

Approaches

Feature-based alignment

Feature-based methods detect keypoints in both images, describe them with descriptors such as SIFT, match descriptors (feature matching), and fit the transform to the matches. Because some matches are always wrong, the fit is made with a robust estimator, most often RANSAC, followed by a least-squares refinement on the inliers. For a homography, four correspondences give a minimal solution and the direct linear transform handles more.

Features tolerate large displacements, rotations and scale changes because matching is global: a keypoint can be found anywhere in the other image. Szeliski (2007) recounts preferring direct methods earlier, then finding that scale-space features with orientation-invariant descriptors matched images well enough for automatic panorama stitching where direct methods failed too often. Features fail where distinctive structure is missing or where the two modalities share no repeatable keypoints.

Intensity-based (direct) alignment

Direct methods use every pixel: they warp the moving image with the current parameters, measure how well it agrees with the reference, and update the parameters. The Lucas–Kanade algorithm, introduced by Lucas and Kanade (1981) as an iterative image registration technique, linearizes the residual and solves for a Gauss–Newton update of a translation or of a general warp such as an affine transform or a homography. Direct methods use all the information in the overlap and reach sub-pixel accuracy, but they converge only from a nearby start. Szeliski (2007) therefore suggests a combination: align with features, then refine with a direct method in a common frame.

Similarity measures

The measure encodes what “aligned” means. Writing $a_i = I_R(\mathbf{x}_i)$ and $b_i = I_M(T(\mathbf{x}_i))$ for the $n$ pixels of the overlap:

• Sum of squared differences (SSD), $\sum_i (a_i - b_i)^2$, is optimal for equal intensities plus Gaussian noise. Robust penalties down-weight outliers, and dividing by the overlap area removes a bias toward solutions with little overlap (Szeliski, 2007).
• Normalized cross-correlation (NCC), $$\mathrm{NCC} = \frac{\sum_i (a_i - \bar{a})(b_i - \bar{b})}{\sqrt{\sum_i (a_i - \bar{a})^2 \sum_i (b_i - \bar{b})^2}},$$ lies in $[-1, 1]$ and is unchanged by a positive gain and an offset applied to either image. It is undefined for a patch with zero variance and degrades in noisy, low-contrast regions (Szeliski, 2007).
• Mutual information (MI) measures statistical dependence instead of equality. With $p_{AB}$ the joint histogram of intensity pairs $(a_i, b_i)$, normalized to sum to 1, and $p_A$, $p_B$ its marginals, $$\mathrm{MI}(A, B) = H(A) + H(B) - H(A, B) = \sum_{a,b} p_{AB}(a, b)\, \log \frac{p_{AB}(a, b)}{p_A(a)\, p_B(b)},$$ where $H$ is entropy. When the images are aligned, a given intensity in one image predicts the intensity in the other, whatever the mapping between them is, and the joint histogram concentrates into a few sharp clusters; misalignment spreads it out.

MI came to registration in the mid-1990s. Viola and Wells presented alignment by maximization of MI at ICCV in June 1995, with a journal version in 1997, aligning MR images to each other and 3D object models to video; they estimated the entropies with Parzen windows from small random samples and maximized MI by stochastic gradient ascent (Viola and Wells, 1997). Zitová and Flusser (2003) name their paper as one of the first to propose the technique. The Leuven group’s 1997 paper states that it expands on ideas first presented by Collignon et al. at the IPMI conference, which its reference list dates to June 1995; it computed MI from the joint histogram of CT, MR and PET volumes, maximized it with Powell’s method, and reported subvoxel accuracy against a stereotactic reference (Maes et al., 1997). In Viola and Wells’s experiments MI still aligned two MR images after one was passed through a non-monotonic (quadratic) intensity mapping, a case that, they note, would defeat traditional correlation. Zitová and Flusser (2003) describe MI methods as the leading technique in multimodal registration.

MI is computed over the overlap, so its value changes with how much of each image the overlap contains. Normalized mutual information, $(H(A) + H(B))/H(A, B)$, was proposed by Studholme, Hill and Hawkes (1999), whose title calls it an overlap invariant entropy measure; Zitová and Flusser (2003) summarize the paper as comparing joint entropy, MI and normalized MI with histogram estimates. The histogram’s bin count, the interpolation used to resample the moving image, and smoothing all shape the MI surface. Maes et al. (1997) proposed partial-volume interpolation, which spreads each sample’s contribution over the histogram entries of its eight grid neighbors with trilinear weights, so that the histogram changes smoothly with the parameters.

Phase correlation

When the images differ by a translation, the Fourier transform turns the shift into a phase ramp, and the inverse transform of the normalized cross-power spectrum is ideally an impulse at the displacement. Phase correlation, for which Szeliski (2007) and the scikit-image documentation cite Kuglin and Hines (1975), finds a translation anywhere in the image’s range in $O(n \log n)$ time for $n$ pixels, with no initial guess, and the peak can be refined to sub-pixel precision; scikit-image’s implementation upsamples the cross-correlation around the peak by a matrix-multiply DFT (Guizar-Sicairos et al., 2008, as cited in its documentation). Whitening the spectrum helps when noise is confined to a few frequencies but can hurt when the signal-to-noise ratio is low at many frequencies, as in blurry images with high-frequency noise (Szeliski, 2007).

Rotation and scale break the single phase ramp. Resampling the Fourier magnitudes, which do not change under translation, onto a log-polar grid turns rotation and scale into shifts that phase correlation can find, a method Szeliski attributes to De Castro and Morandi (1987) and that he notes works only when the images overlap largely.

Optimization and coarse-to-fine search

The search depends on the transform and the measure. Exhaustive search over shifts is used for small translations, as in block matching. Gradient-based methods (nonlinear least squares solvers such as Gauss–Newton and Levenberg–Marquardt for SSD-type measures, gradient ascent for MI) need derivatives of the measure; derivative-free methods such as Powell’s direction-set search, used by Maes et al. (1997), need only function values, which (our reasoning) suits histogram-based MI, whose value jumps slightly as samples cross bin boundaries.

All of them benefit from coarse-to-fine estimation: registering smoothed, subsampled versions of the images from an image pyramid first, where large motions become small and the objective is smoother, then refining at finer levels. Szeliski (2007) notes that in practice only two or three pyramid levels can be used before important detail blurs away, which limits the capture range of direct methods. A good initialization, from features, phase correlation, sensor metadata or the previous frame of a video, is often what decides success.

Deformable registration

When the geometry differs locally, the transform is a displacement field $\mathbf{u}$, with $T(\mathbf{x}) = \mathbf{x} + \mathbf{u}(\mathbf{x})$.

• Free-form deformations (FFDs) place a coarse grid of control points over the image and interpolate their displacements with cubic B-splines, giving a smooth field with few parameters; each point is influenced only by nearby control points. Sotiras et al. (2013) describe FFDs as one of the most common deformation models in medical image registration, adopted widely once coupled with cubic B-splines, as in Rueckert et al.’s (1999) registration of breast MR images. The resulting deformations are smooth but need not preserve topology unless the displacements are constrained.
• Demons. Thirion (1998) treated matching as a diffusion process. Every pixel is a “demon” that pushes the moving image along the reference image’s gradient $\nabla s$ by an optical-flow-like displacement, $$\mathbf{v} = \frac{(m - s)\,\nabla s}{\lVert \nabla s \rVert^2 + (m - s)^2},$$ where $m$ and $s$ are the moving and reference intensities; the second term in the denominator keeps the step bounded where the gradient is small. The accumulated field is smoothed with a Gaussian after each iteration, which acts as the regularizer. Because the push is driven by the difference $m - s$, demons assume that corresponding points have equal intensities, as SSD does.
• Learned registration. VoxelMorph (Balakrishnan et al., 2019) trains a convolutional network that takes the two volumes and outputs a displacement field, warps the moving volume with a differentiable spatial transformer, and minimizes an image dissimilarity (mean squared error or local cross-correlation) plus a smoothness penalty, without ground-truth deformations. For 3D brain MRI registered to an atlas, the authors reported an average Dice overlap of anatomical structures on their test set of 0.753 for the cross-correlation version, against 0.749 for ANTs SyN and 0.755 for NiftyReg, and 0.584 for affine alignment alone. VoxelMorph took 0.45 s per pair on a GPU and 57 s on a CPU, against 9,059 s and 2,347 s on a CPU for ANTs and NiftyReg (their Table I). The network replaces a per-pair optimization with a single forward pass.

Deep networks are also used inside the feature-based pipeline, as learned keypoints and matchers; the feature matching article describes them.

Example

The code below simulates a multimodal pair with a known affine transform. A synthetic scene is cropped to a $128 \times 128$ reference; the second image sees the same scene shifted by $(6.4, -3.7)$ pixels, rotated by 8°, scaled by 1.06 and 0.96 along its axes and sheared by 0.05, with noise. Its intensities are then inverted and squared, so bright structures in the reference are dark in the second “modality”. The code recovers a pure translation by phase correlation, compares SSD, NCC and MI at the true transform and at the identity, and registers the multimodal pair by maximizing MI, and for comparison by minimizing SSD, with Powell’s method on Gaussian-smoothed images and then on the unsmoothed ones. It reports the recovered parameters and the mean target registration error (TRE) at five points spread over the image.

import numpy as np
from scipy import ndimage, optimize
from skimage.filters import window
from skimage.registration import phase_cross_correlation

rng = np.random.default_rng(0)
M, N, o = 192, 128, 32                       # scene size, image size, crop offset
Y, X = np.mgrid[0:M, 0:M].astype(float)

# A synthetic scene: blobs, bars, a ramp and fine texture (std of the smoothed noise ~0.18)
scene = 0.3 * X / M + 0.25 * ndimage.gaussian_filter(rng.standard_normal((M, M)), 1.5) / 0.18
for x0, y0, s, a in [(72, 67, 9, 1.0), (122, 82, 14, 0.7), (92, 127, 6, 0.9), (132, 132, 10, 0.5)]:
    scene += a * np.exp(-((X - x0) ** 2 + (Y - y0) ** 2) / (2 * s**2))
scene[52:62, 102:147] += 0.6
scene[87:142, 50:58] += 0.8
scene = ndimage.gaussian_filter(scene, 1.0)
scene = (scene - scene.min()) / (scene.max() - scene.min())

yy, xx = np.mgrid[0:N, 0:N].astype(float)
grid = np.stack([xx.ravel(), yy.ravel()])
c = np.array([(N - 1) / 2, (N - 1) / 2])     # transforms act about the image center
ref = scene[o:o + N, o:o + N]

def affine(p):
    """p = (tx, ty, rotation in degrees, sx, sy, shear) -> A, t"""
    tx, ty, th, sx, sy, k = p
    r = np.radians(th)
    R = np.array([[np.cos(r), -np.sin(r)], [np.sin(r), np.cos(r)]])
    return R @ np.array([[sx, k], [0, sy]]), np.array([tx, ty])

def mapped(p, pts):                          # T(x) = A (x - c) + c + t
    A, t = affine(p)
    return A @ (pts - c[:, None]) + (c + t)[:, None]

def sample(img, q, offset=0):                # cubic interpolation, NaN outside
    v = ndimage.map_coordinates(img, [q[1] + offset, q[0] + offset], order=3,
                                mode="constant", cval=np.nan)
    return v.reshape(N, N)

def capture(p, noise):
    """Second image of the scene: moving(T(x)) = ref(x), plus sensor noise."""
    A, t = affine(p)
    q = np.linalg.solve(A, grid - (c + t)[:, None]) + c[:, None]
    return sample(scene, q, offset=o) + noise * rng.standard_normal((N, N))

truth = np.array([6.4, -3.7, 8.0, 1.06, 0.96, 0.05])
shifted = capture([6.4, -3.7, 0, 1, 1, 0], 0.01)
moving = capture(truth, 0.01)
other = (1 - np.clip(moving, 0, 1)) ** 2 + 0.02 * rng.standard_normal((N, N))  # "other modality"

# 1. Phase correlation (Hann window against border effects)
W = window("hann", (N, N))
def phase_corr(img):
    d = phase_cross_correlation(ref * W, img * W, upsample_factor=100)[0]
    return -d[::-1]                          # (tx, ty) of the content
print("translation only, truth (6.40, -3.70):  estimate (%.2f, %.2f)" % tuple(phase_corr(shifted)))
print("affine, same modality:                  estimate (%.2f, %.2f)" % tuple(phase_corr(moving)))
print("affine, other modality:                 estimate (%.2f, %.2f)" % tuple(phase_corr(other)))

# 2. Similarity measures at the identity and at the true transform
def overlap(img, p, r=ref):
    w = sample(img, mapped(p, grid))
    m = ~np.isnan(w)
    return r[m], w[m]

def ssd(a, b):
    return np.mean((a - b) ** 2)

def ncc(a, b):
    a, b = a - a.mean(), b - b.mean()
    return (a * b).sum() / np.sqrt((a * a).sum() * (b * b).sum())

def mi(a, b, bins=32):
    p = np.histogram2d(a, b, bins=bins)[0]
    p /= p.sum()
    pa, pb = p.sum(1, keepdims=True), p.sum(0, keepdims=True)
    nz = p > 0
    return (p[nz] * np.log(p[nz] / (pa @ pb)[nz])).sum()

identity = np.array([0, 0, 0, 1, 1, 0.0])
print("\n                           SSD     NCC   MI (nats)")
for name, img in [("same modality", moving), ("other modality", other)]:
    for where, p in [("identity", identity), ("truth", truth)]:
        a, b = overlap(img, p)
        print(f"{name:14s} at {where:8s}  {ssd(a, b):.4f}  {ncc(a, b):+.3f}  {mi(a, b):.3f}")

# 3. Register the other-modality pair: Powell search, smoothed level then full detail
units = np.array([1, 1, 1, 0.01, 0.01, 0.01])
def register(img, cost, start):
    p = start.copy()
    for sigma in (2.0, 0.0):
        r_s = ndimage.gaussian_filter(ref, sigma) if sigma else ref
        i_s = ndimage.gaussian_filter(img, sigma) if sigma else img
        def f(z):
            a, b = overlap(i_s, z * units, r_s)
            return cost(a, b) if a.size > N * N / 2 else 1e3
        p = optimize.minimize(f, p / units, method="Powell",
                              options={"xtol": 1e-3, "ftol": 1e-6}).x * units
    return p

start = np.r_[phase_corr(other), 0, 1, 1, 0]
est_mi = register(other, lambda a, b: -mi(a, b), start)
est_ssd = register(other, ssd, start)

targets = np.array([[16, 16], [112, 16], [16, 112], [112, 112], [64, 64]], float).T
def tre(p):                                  # mean error at five target points (pixels)
    return np.linalg.norm(mapped(p, targets) - mapped(truth, targets), axis=0).mean()

print("\n          tx      ty    rot     sx     sy    shear   TRE (px)")
for name, p in [("truth", truth), ("start", start), ("MI", est_mi), ("SSD", est_ssd)]:
    print(f"{name:5s}  {p[0]:6.2f}  {p[1]:6.2f}  {p[2]:5.2f}  {p[3]:.3f}  {p[4]:.3f}  {p[5]:+.3f}"
          f"  {tre(p):7.2f}")

Output:

translation only, truth (6.40, -3.70):  estimate (6.38, -3.57)
affine, same modality:                  estimate (-3.24, 5.03)
affine, other modality:                 estimate (4.28, -6.73)

                           SSD     NCC   MI (nats)
same modality  at identity  0.0251  +0.381  0.120
same modality  at truth     0.0001  +0.998  2.279
other modality at identity  0.0863  -0.361  0.115
other modality at truth     0.1167  -0.971  1.810

          tx      ty    rot     sx     sy    shear   TRE (px)
truth    6.40   -3.70   8.00  1.060  0.960  +0.050     0.00
start    4.28   -6.73   0.00  1.000  1.000  +0.000     7.68
MI       6.41   -3.71   7.99  1.060  0.959  +0.050     0.02
SSD    -55.85    4.43  178.40  0.585  1.003  +0.014   118.82

Both images are multiplied by a Hann window before phase correlation, because their borders are not periodic. For a pure translation, phase correlation recovers the shift to within 0.13 pixels without any starting guess. Once the second image is also rotated and scaled, no single shift relates the two images, and the estimates are off by several pixels and differ between the two modalities. Our reading: the rotation spreads the correlation peak, so which part of the scene dominates it depends on the intensities. scikit-image (in the 0.25.2 source) locates the peak by its magnitude, so the global contrast inversion alone does not defeat it.

The measure table shows why SSD cannot register the multimodal pair: its value at the true transform is higher than at the identity, so it prefers misalignment. NCC is close to $-1$ at the truth, because the intensity mapping is decreasing; maximizing $|\mathrm{NCC}|$ would work for this monotonic mapping, but not for a non-monotonic one such as Viola and Wells’s quadratic example. MI rises from about 0.1 nats at the identity to 2.3 nats (same modality) and 1.8 nats (other modality) at the truth.

Starting from the phase-correlation translation, whose TRE is 7.7 pixels, MI maximization recovers all six parameters, with a TRE of 0.02 pixels. SSD from the same start runs to a transform rotated by about 180°, more than 100 pixels from the truth: with inverted intensities, its minimum is not at the true alignment. The TRE here is exact because the transform is synthetic; real evaluations have to measure it with markers or landmarks.

Datasets and Benchmarks

Registration is often evaluated on synthetic transforms like the one above, because they give exact ground truth. For real images, ground truth comes from external markers, landmarks annotated by experts, or labeled structures:

• Medical rigid registration. Zitová and Flusser (2003) call the Vanderbilt University project, a comparison of retrospective intermodality brain registration techniques published by West et al. in 1997, probably the most comprehensive case study in the validation of rigid-body registration. Maes et al. (1997) validated MI against a stereotactic reference solution.
• Deformable brain MRI. VoxelMorph was evaluated on 3,731 T1-weighted brain scans from eight public datasets, split into 3,231 training, 250 validation and 250 test volumes, and on the Buckner40 set with manual segmentations (Balakrishnan et al., 2019).
• Homographies between photographs. HPatches provides image sequences with ground-truth homographies, used to evaluate homography estimation from matched features (see feature matching).

Evaluation Metrics

• Target registration error (TRE) is the distance between where the estimated and the true transforms send a point of interest, $\lVert \hat{T}(\mathbf{x}) - T(\mathbf{x}) \rVert$, at points not used to compute the registration. Zitová and Flusser (2003) credit Fitzpatrick and colleagues with using anatomical points inside the region of interest as such “target points”.
• Landmark error. With corresponding landmarks $\mathbf{p}_k$ and $\mathbf{q}_k$ marked in the two images, the mean or maximum of $\lVert \hat{T}(\mathbf{p}_k) - \mathbf{q}_k \rVert$ estimates the TRE but includes the error of locating the landmarks.
• Control-point error versus test-point error. The residual at the points used to fit the transform measures only how well the model fits them; a flexible enough model drives it to zero (overfitting), so Zitová and Flusser recommend the error at held-out test points instead.
• Overlap of warped labels. In deformable registration, the Dice coefficient between a segmentation warped from the moving image and the reference segmentation measures anatomical agreement, as in VoxelMorph’s evaluation.
• Plausibility of the deformation. The fraction of voxels where the Jacobian determinant of the transform is not positive counts folds, where the deformation is not invertible; Balakrishnan et al. (2019) report it next to Dice.
• Corner error for homographies. The reprojection error of the image corners under the estimated homography, as in SuperGlue’s homography evaluation (see feature matching).
• Consistency checks and visual inspection. Without ground truth, the result can be compared with that of a gold-standard or otherwise different method, two images registered to the same reference can also be registered to each other and the two routes compared, and a domain expert can inspect the result visually (Zitová and Flusser, 2003).

A low intensity residual or a high similarity score after registration shows that the optimizer converged, not that the alignment is correct.

Applications

• Panoramas and image stitching. Photos from a rotating camera are registered by homographies and blended (Szeliski, 2007).
• Video stabilization. Consecutive frames are registered by a global parametric model, a form of motion estimation, to estimate and remove camera shake (Szeliski, 2007).
• Multi-frame enhancement. Image fusion and multi-frame super-resolution start by registering the frames (Zitová and Flusser, 2003).
• Remote sensing and cartography. Multispectral and multitemporal images are registered for classification, environmental monitoring, change detection and map updating (Zitová and Flusser, 2003).
• Medicine. CT and MRI of the same patient are fused, follow-up scans are compared to monitor treatment, and scans are mapped to anatomical atlases (Zitová and Flusser, 2003).
• Template tracking. A planar target is tracked frame to frame with a homography estimated by features or by direct alignment, for example to anchor content in augmented reality.
• Industrial inspection. Zitová and Flusser (2003) list automatic quality control among the computer vision uses: a part is aligned to a reference before it is compared with it.

Open Problems

Zitová and Flusser (2003) named registration with complex local distortions, multimodal registration and the cost of registering large 3D and higher-dimensional data as the hardest problems, and noted that validation of non-rigid methods was only beginning. They called the distinction between a local deformation and a real change of the scene a conceptual question. Learned methods are fast once trained, but a network learns from one kind of data, and (our reasoning) nothing guarantees that it transfers to another scanner, body region or camera without retraining. Zitová and Flusser also observed that MI methods, being area-based, fail or become very slow when multimodal images differ by large rotations or scale changes, and pointed to modality-insensitive features as the way forward.