[ VIGRA Homepage | Function Index | Class Index | Namespaces | File List | Main Page ]

details Image Registration VIGRA

Classes

class  AffineMotionEstimationOptions< SPLINEORDER >
 Option object for affine registration functions. More...
 
struct  DistancePowerFunctor< N >
 
struct  ThinPlateSplineFunctor
 

Functions

template<class SrcIterator , class DestIterator >
linalg::TemporaryMatrix< double > affineMatrix2DFromCorrespondingPoints (SrcIterator s, SrcIterator send, DestIterator d)
 Create homogeneous matrix that maps corresponding points onto each other. More...
 
template<... >
void estimateAffineTransform (...)
 Estimate the optical flow between two images according to an affine transform model (e.g. translation, rotation, non-uniform scaling, and shearing). More...
 
template<... >
void estimateGlobalRotation (...)
 Estimate the rotation between two images by means of a normalized cross correlation matching of the FFT spectra. More...
 
template<... >
void estimateGlobalRotationTranslation (...)
 Estimate the (global) rotation and translation between two images by means a normalized cross correlation matching. More...
 
template<... >
void estimateGlobalTranslation (...)
 Estimate the translation between two images by means of a normalized cross correlation matching. More...
 
template<... >
void estimateSimilarityTransform (...)
 Estimate the optical flow between two images according to a similarity transform model (e.g. translation, rotation, and uniform scaling). More...
 
template<... >
void estimateTranslation (...)
 Estimate the optical flow between two images according to a translation model. More...
 
template<int PolynomOrder, class SrcPointIterator , class DestPointIterator >
linalg::TemporaryMatrix< double > polynomialMatrix2DFromCorrespondingPoints (SrcPointIterator s, SrcPointIterator s_end, DestPointIterator d)
 Create polynomial matrix of a certain degree that maps corresponding points onto each other. More...
 
template<... >
void polynomialWarpImage (...)
 Warp an image according to an polynomial transformation. More...
 
std::vector< double > polynomialWarpWeights (double x, double y, unsigned int polynom_order)
 
template<class SrcPointIterator , class DestPointIterator >
linalg::TemporaryMatrix< double > projectiveMatrix2DFromCorrespondingPoints (SrcPointIterator s, SrcPointIterator send, DestPointIterator d)
 Create homogeneous matrix that maps corresponding points onto each other. More...
 
template<... >
void projectiveWarpImage (...)
 Warp an image according to an projective transformation. More...
 
template<class RadialBasisFunctor , class SrcPointIterator , class DestPointIterator >
linalg::TemporaryMatrix< double > rbfMatrix2DFromCorrespondingPoints (SrcPointIterator s, SrcPointIterator s_end, DestPointIterator d, RadialBasisFunctor const &rbf)
 Create a matrix that maps corresponding points onto each other using a given RBF. More...
 
template<int ORDER, class T , class DestIterator , class DestAccessor , class DestPointIterator , class C , class RadialBasisFunctor >
void rbfWarpImage (SplineImageView< ORDER, T > const &src, DestIterator dul, DestIterator dlr, DestAccessor dest, DestPointIterator d, DestPointIterator d_end, MultiArrayView< 2, double, C > const &W, RadialBasisFunctor rbf)
 Warp an image according to an radial basis function based transformation. More...
 
template<class SplineImage , class DestIterator , class DestAccessor >
void transformToPolarCoordinates (SplineImage const &src, DestIterator d_ul, DestIterator d_lr, DestAccessor d_acc)
 Transforms a given image to its (image-centered) polar coordinates representation. More...
 

Detailed Description

Transform different image into a common coordinate system.

Function Documentation

linalg::TemporaryMatrix<double> vigra::affineMatrix2DFromCorrespondingPoints ( SrcIterator  s,
SrcIterator  send,
DestIterator  d 
)

Create homogeneous matrix that maps corresponding points onto each other.

For use with affineWarpImage(). When only two corresponding points are given, the matrix will only represent a similarity transform (translation, rotation, and uniform scaling). When only one point pair is given, the result will be a pure translation.

void vigra::estimateTranslation (   ...)

Estimate the optical flow between two images according to a translation model.

This function applies the same algorithm as estimateAffineTransform() with the additional constraint that the motion model must be a translation rather than affine.

Declarations:

#include <vigra/affine_registration.hxx>
Namespace: vigra

pass 2D array views:

namespace vigra {
template <class T1, class S1,
class T2, class S2,
int SPLINEORDER>
void
estimateTranslation(MultiArrayView<2, T1, S1> const & src,
MultiArrayView<2, T2, S2> dest,
Matrix<double> & affineMatrix,
AffineMotionEstimationOptions<SPLINEORDER> const & options = AffineMotionEstimationOptions<SPLINEORDER>());
}

show deprecated declarations

void vigra::estimateSimilarityTransform (   ...)

Estimate the optical flow between two images according to a similarity transform model (e.g. translation, rotation, and uniform scaling).

This function applies the same algorithm as estimateAffineTransform() with the additional constraint that the motion model must be a similarity transform rather than affine.

Declarations:

#include <vigra/affine_registration.hxx>
Namespace: vigra

pass 2D array views:

namespace vigra {
template <class T1, class S1,
class T2, class S2,
int SPLINEORDER>
void
estimateSimilarityTransform(MultiArrayView<2, T1, S1> const & src,
MultiArrayView<2, T2, S2> dest,
Matrix<double> & affineMatrix,
AffineMotionEstimationOptions<SPLINEORDER> const & options =
AffineMotionEstimationOptions<SPLINEORDER>());
}

show deprecated declarations

void vigra::estimateAffineTransform (   ...)

Estimate the optical flow between two images according to an affine transform model (e.g. translation, rotation, non-uniform scaling, and shearing).

This function implements the algorithm described in

J.R. Bergen, P. Anandan, K.J. Hanna, R. Hingorani: "Hierarchical model-based motion estimation", ECCV 1992

Specifically, it minimizes the squared loss between the images I at two consecutive time points t-1 and t:

\[ \min_{\theta} \sum_{\mathbf{x}} \left(I(\mathbf{x}, t) - I(\mathbf{x} - \mathbf{u}_{\theta}(\mathbf{x}), t-1)\right)^2 \]

where $\mathbf{x}$ are the pixel coordinates and $\mathbf{u}_{\theta}(\mathbf{x})$ is an affine motion model parameterized by $\theta$. Since the objective is non-linear, it is linearized by first-order Taylor expansion w.r.t. $\theta$, and a local optimum is determined iteratively by the Gauss-Newton method. To handle larger displacements, the algorithm employs a coarse-to-fine strategy, where the motion is first estimated on downsampled versions of the images and then refined at consecutively higher resolutions.

The algorithm's parameters can be controlled by the option object vigra::AffineMotionEstimationOptions. In particular, one can determine if

  • I in the objective refers to the original images (default) or their second derivatives (options.useLaplacianPyramid() – makes motion estimation invariant against additive intensity offsets);
  • the highest pyramid level to be used (options.highestPyramidLevel(h) – images are downsampled to 2-h times their original size, default: h=4);
  • the number of Gauss-Newton iterations per resolution level (options.iterationsPerLevel(i), default: i=4);
  • the interpolation order to compute subpixel intensities $I(\mathbf{x} - \mathbf{u}_{\theta}(\mathbf{x}), t-1)$ (default: 2)

The resulting affine model is stored in parameter affineMatrix, which can be used by affineWarpImage() to apply the transformation to time frame t-1. See documentation there for the precise meaning of the matrix elements.

Declarations:

#include <vigra/affine_registration.hxx>
Namespace: vigra

pass 2D array views:

namespace vigra {
template <class T1, class S1,
class T2, class S2,
int SPLINEORDER>
void
estimateAffineTransform(MultiArrayView<2, T1, S1> const & src,
MultiArrayView<2, T2, S2> dest,
Matrix<double> & affineMatrix,
AffineMotionEstimationOptions<SPLINEORDER> const & options =
AffineMotionEstimationOptions<SPLINEORDER>());
}

show deprecated declarations

void vigra::transformToPolarCoordinates ( SplineImage const &  src,
DestIterator  d_ul,
DestIterator  d_lr,
DestAccessor  d_acc 
)

Transforms a given image to its (image-centered) polar coordinates representation.

This algorithm transforms a given image (by means of an spline image view) to its image-centered polar coordinates reprensentation. The sampling of the polar coordinate system is determined by the shape of the dest. image.

Declarations:

#include <vigra/affine_registration_fft.hxx>
Namespace: vigra

pass 2D array views:

namespace vigra {
template <class SplineImage,
class T1, class S1>
void
transformToPolarCoordinates(SplineImage const & src,
MultiArrayView<2, T1, S1> dest);
}

show deprecated declarations

void vigra::estimateGlobalRotation (   ...)

Estimate the rotation between two images by means of a normalized cross correlation matching of the FFT spectra.

This algorithm uses the fast normalized cross correlation to determine a global rotation between two images (from image2 to image1). To derive the rotation, the algorithm performs the following steps:

  1. Transforming both images to the frequency domain using FFTW
  2. Create LogAbs PolarCoordinate representations for each spectrum.
  3. Determining the final Rotation by performing a fast normalized cross correlation based on the polar representations.

The images are cropped to the corresponding images center-squared before the estimation takes place.

Declarations:

#include <vigra/affine_registration_fft.hxx>
Namespace: vigra

pass 2D array views:

namespace vigra {
template <class T1, class S1,
class T2, class S2>
void
estimateGlobalRotation(MultiArrayView<2, T1, S1> const & src,
MultiArrayView<2, T2, S2> dest,
Matrix<double> & affineMatrix,
double & correlation_coefficent);
}

show deprecated declarations

void vigra::estimateGlobalTranslation (   ...)

Estimate the translation between two images by means of a normalized cross correlation matching.

This algorithm uses the fast normalized cross correlation to determine a global translation between two images (from image2 to image1). To derive the translation, the algorithm consists of differents steps:

  1. Separation of the second image
    The second image (the one, for which the translation shall be determined) is cut into five subregions: UpperLeft, UpperRight, Center, LowerLeft and LowerRight, each of 1/4 the size of the original image. Using a border > 0 results in (all) overlapping regions.
  2. Cross-Correlation of the subimages to the first image
    The subimages are normalized cross-correlated to the (complete) first image. The resulting maximum-likelihood translations and the correlation coefficients are stored for the next step.
  3. Determining the final Translation by voting
    Each correlation vector gets one vote at the beginning. For each equality of derived motion vectors, the voting to these vectors is incremented. If the maximum number of votes is larger than 1, the vector with the most votes is chosen. If the maximum number of votes is 1, the vector with the maximum likelihood is choosen.

Declarations:

#include <vigra/affine_registration_fft.hxx>
Namespace: vigra

pass 2D array views:

namespace vigra {
template <class T1, class S1,
class T2, class S2>
void
estimateGlobalTranslation(MultiArrayView<2, T1, S1> const & src,
MultiArrayView<2, T2, S2> dest,
Matrix<double> & affineMatrix,
double & correlation_coefficent,
Diff2D border = Diff2D(0,0));
}

show deprecated declarations

void vigra::estimateGlobalRotationTranslation (   ...)

Estimate the (global) rotation and translation between two images by means a normalized cross correlation matching.

This algorithm use the functions estimateGlobalRotation() and estimateGlobalTranslation() to estimate a matrix which describes the global rotation and translation from the second to the first image.

Declarations:

#include <vigra/affine_registration_fft.hxx>
Namespace: vigra

pass 2D array views:

namespace vigra {
template <class T1, class S1,
class T2, class S2>
void
estimateGlobalRotationTranslation(MultiArrayView<2, T1, S1> const & src,
MultiArrayView<2, T2, S2> dest,
Matrix<double> & affineMatrix,
double & rotation_correlation,
double & translation_correlation,
Diff2D border = Diff2D(0,0));
}

show deprecated declarations

std::vector<double> vigra::polynomialWarpWeights ( double  x,
double  y,
unsigned int  polynom_order 
)

Iterative function for determinination of the polynom weights:

Example: order=2, x, y --—> [1, x, y, x^2, x*y, y^2]

This function is needed, because the polynomial transformation Matrix has the the same number of rows. the target position is then determined by multiplying each x- and y-transformation result value with the corresponding weight for the current x- and y-coordinate, given by this function.

linalg::TemporaryMatrix<double> vigra::polynomialMatrix2DFromCorrespondingPoints ( SrcPointIterator  s,
SrcPointIterator  s_end,
DestPointIterator  d 
)

Create polynomial matrix of a certain degree that maps corresponding points onto each other.

For use with polynomialWarpImage() of same degree.

Since polynoms are usually non-linear functions, a special semantics is embedded to define a matrix here. Each matrix consist of two rows, containing x- and y-factors of the polynom.

The meaning of the matrix is explained at the example of a polynom of 2nd order:

First Row = [a_x b_x c_x d_x e_x f_x] Second Row = [a_y b_y c_y d_y e_y f_y]

The transformed coordinate p'=[x' y'] of a position p=[x y] is then:

x' = a_x + b_x*x + c_x*y + d_x*x^2 + e_x*x*y + f_x*y^2 y' = a_y + b_y*x + c_y*y + d_y*x^2 + e_y*x*y + f_y*y^2

Note that the order of the polynom's factors is directly influenced by the polynomialWarpWeights() function and follows the intuitive scheme.

void vigra::polynomialWarpImage (   ...)

Warp an image according to an polynomial transformation.

To get more information about the structure of the matrix, see polynomialMatrix2DFromCorrespondingPoints().

#include <vigra/polynomial_registration.hxx>
Namespace: vigra

pass 2D array views:

namespace vigra {
template <int ORDER, class T,
class T2, class S2,
class C>
void
polynomialWarpImage(SplineImageView<ORDER, T> const & src,
MultiArrayView<2, T2, S2> dest,
MultiArrayView<2, double, C> const & polynomialMatrix);
}

show deprecated declarations

linalg::TemporaryMatrix<double> vigra::projectiveMatrix2DFromCorrespondingPoints ( SrcPointIterator  s,
SrcPointIterator  send,
DestPointIterator  d 
)

Create homogeneous matrix that maps corresponding points onto each other.

For use with projectiveWarpImage(). Since four corresponding points are needed to be given, the matrix will compute a full projective transform.

void vigra::projectiveWarpImage (   ...)

Warp an image according to an projective transformation.

The matrix can be computed from a set of correspondung points using projectiveMatrix2DFromCorrespondingPoints().

Declarations:

#include <vigra/projective_registration.hxx>
Namespace: vigra

pass 2D array views:

namespace vigra {
template <int ORDER, class T,
class T2, class S2,
class C>
void
projectiveWarpImage(SplineImageView<ORDER, T> const & src,
MultiArrayView<2, T2, S2> dest,
MultiArrayView<2, double, C> const & projectiveMatrix);
}

show deprecated declarations

The algorithm applies the given projectiveMatrix to the destination coordinates and copies the image value from the resulting source coordinates, using the given SplineImageView src for interpolation. If the resulting coordinate is outside the source image, nothing will be written at that destination point.

for all dest pixels:
currentSrcCoordinate = projectiveMatrix * currentDestCoordinate;
if src.isInside(currentSrcCoordinate):
dest[currentDestCoordinate] = src[currentSrcCoordinate]; // copy an interpolated value

The matrix represent a 2-dimensional projective transform by means of homogeneous coordinates, i.e. it must be a 3x3 matrix whose last row is (p1,p2,1).

Required Interface:

DestImageIterator dest_upperleft;
double x = ..., y = ...;
if (spline.isInside(x,y))
dest_accessor.set(spline(x, y), dest_upperleft);

See also: Functions to specify projective transformation: translationMatrix2D(), scalingMatrix2D(), shearMatrix2D(), rotationMatrix2DRadians(), rotationMatrix2DDegrees() and projectiveMatrix2DFromCorrespondingPoints()

linalg::TemporaryMatrix<double> vigra::rbfMatrix2DFromCorrespondingPoints ( SrcPointIterator  s,
SrcPointIterator  s_end,
DestPointIterator  d,
RadialBasisFunctor const &  rbf 
)

Create a matrix that maps corresponding points onto each other using a given RBF.

For use with rbfWarpImage(). For n given (corresponding) points, the matrix will be of size (n+3,2). Note that the representation of this matrix is exactly the same as the "W" matrix of Bookstein. More information can be found in the following article:

Fred L. Bookstein. Principal Warps: Thin-Plate Splines and the Decomposition of Deformations. IEEE PAMI, Vol 11, No 8. 1989

void vigra::rbfWarpImage ( SplineImageView< ORDER, T > const &  src,
DestIterator  dul,
DestIterator  dlr,
DestAccessor  dest,
DestPointIterator  d,
DestPointIterator  d_end,
MultiArrayView< 2, double, C > const &  W,
RadialBasisFunctor  rbf 
)

Warp an image according to an radial basis function based transformation.

To get more information about the structure of the matrix, see rbfMatrix2DFromCorrespondingPoints()

#include <vigra/rbf_registration.hxx>
Namespace: vigra

pass 2D array views:

namespace vigra {
template <int ORDER, class T,
class T2, class S2,
class DestPointIterator,
class C,
class RadialBasisFunctor>
void
rbfWarpImage(SplineImageView<ORDER, T> const & src,
MultiArrayView<2, T2, S2> dest,
DestPointIterator d, DestPointIterator d_end,
MultiArrayView<2, double, C> const & W,
RadialBasisFunctor rbf);
}

show deprecated declarations

© Ullrich Köthe (ullrich.koethe@iwr.uni-heidelberg.de)
Heidelberg Collaboratory for Image Processing, University of Heidelberg, Germany

html generated using doxygen and Python
vigra 1.11.1 (Fri May 19 2017)