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

affine_registration_fft.hxx VIGRA

1 /************************************************************************/
2 /* */
3 /* Copyright 2007-2014 by Benjamin Seppke */
4 /* */
5 /* This file is part of the VIGRA computer vision library. */
6 /* The VIGRA Website is */
7 /* http://hci.iwr.uni-heidelberg.de/vigra/ */
8 /* Please direct questions, bug reports, and contributions to */
9 /* ullrich.koethe@iwr.uni-heidelberg.de or */
10 /* vigra@informatik.uni-hamburg.de */
11 /* */
12 /* Permission is hereby granted, free of charge, to any person */
13 /* obtaining a copy of this software and associated documentation */
14 /* files (the "Software"), to deal in the Software without */
15 /* restriction, including without limitation the rights to use, */
16 /* copy, modify, merge, publish, distribute, sublicense, and/or */
17 /* sell copies of the Software, and to permit persons to whom the */
18 /* Software is furnished to do so, subject to the following */
19 /* conditions: */
20 /* */
21 /* The above copyright notice and this permission notice shall be */
22 /* included in all copies or substantial portions of the */
23 /* Software. */
24 /* */
25 /* THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND */
26 /* EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES */
27 /* OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND */
28 /* NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT */
29 /* HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, */
30 /* WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING */
31 /* FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR */
32 /* OTHER DEALINGS IN THE SOFTWARE. */
33 /* */
34 /************************************************************************/
35 
36 #ifndef VIGRA_AFFINE_REGISTRATION_FFT_HXX
37 #define VIGRA_AFFINE_REGISTRATION_FFT_HXX
38 
39 #include "mathutil.hxx"
40 #include "matrix.hxx"
41 #include "linear_solve.hxx"
42 #include "tinyvector.hxx"
43 #include "splineimageview.hxx"
44 #include "imagecontainer.hxx"
45 #include "multi_shape.hxx"
46 #include "affinegeometry.hxx"
47 #include "correlation.hxx"
48 
49 #include <cmath>
50 
51 namespace vigra {
52 
53 /** \addtogroup Registration Image Registration
54 */
55 //@{
56 
57 /********************************************************/
58 /* */
59 /* transformToPolarCoordinates */
60 /* */
61 /********************************************************/
62 
63 /** \brief Transforms a given image to its (image-centered) polar coordinates representation.
64 
65  This algorithm transforms a given image (by means of an spline image view) to its
66  image-centered polar coordinates reprensentation. The sampling of the polar coordinate system
67  is determined by the shape of the dest. image.
68 
69  <b> Declarations:</b>
70 
71  <b>\#include</b> <vigra/affine_registration_fft.hxx><br>
72  Namespace: vigra
73 
74  pass 2D array views:
75  \code
76  namespace vigra {
77  template <class SplineImage,
78  class T1, class S1>
79  void
80  transformToPolarCoordinates(SplineImage const & src,
81  MultiArrayView<2, T1, S1> dest);
82  }
83  \endcode
84 
85  \deprecatedAPI{estimateTranslation}
86  pass \ref ImageIterators and \ref DataAccessors :
87  \code
88  namespace vigra {
89  template <class SplineImage,
90  class DestIterator, class DestAccessor>
91  void
92  transformToPolarCoordinates(SplineImage const & src,
93  DestIterator dul, DestIterator dlr, DestAccessor dest)
94  }
95  \endcode
96  use argument objects in conjunction with \ref ArgumentObjectFactories :
97  \code
98  namespace vigra {
99  template <class SplineImage,
100  class DestIterator, class DestAccessor>
101  void
102  transformToPolarCoordinates(SplineImage const & src,
103  triple<DestIterator, DestIterator, DestAccessor> dest)
104  }
105  \endcode
106  \deprecatedEnd
107 */
108 
109 template <class SplineImage,
110  class DestIterator, class DestAccessor>
111 void
112 transformToPolarCoordinates(SplineImage const & src,
113  DestIterator d_ul, DestIterator d_lr, DestAccessor d_acc)
114 {
115  typename DestIterator::difference_type d_shape = (d_lr - d_ul);
116 
117  int s_w = src.width(),
118  s_h = src.height();
119 
120  int s_size = min(s_w, s_h);
121 
122  int d_w = d_shape.x,
123  d_h = d_shape.y;
124 
125  double r_max = s_size / 2.0;
126 
127  DestIterator yd = d_ul;
128  DestIterator xd = yd;
129 
130  for (int t_step = 0; t_step < d_h; t_step++, yd.y++)
131  {
132  xd = yd;
133  for (int r_step = 0; r_step < d_w; r_step++, xd.x++)
134  {
135  double theta = 2.0 * M_PI * double(t_step) / double(d_h);
136  double r = r_max * double(r_step) / double(d_w);
137  double u = r * cos(theta) + r_max;
138  double v = r * -sin(theta) + r_max;
139 
140  if ( u >= 0 && u < s_size
141  && v >= 0 && v < s_size)
142  {
143  d_acc.set(src(u, v), xd);
144  }
145  }
146  }
147 }
148 
149 template <class SplineImage,
150  class DestIterator, class DestAccessor>
151 inline void
152 transformToPolarCoordinates(SplineImage const & src,
153  vigra::triple<DestIterator, DestIterator, DestAccessor> dest)
154 {
155  transformToPolarCoordinates(src, dest.first, dest.second, dest.third);
156 }
157 
158 template <class SplineImage,
159  class T1, class S1>
160 void
161 transformToPolarCoordinates(SplineImage const & src,
162  MultiArrayView<2, T1, S1> dest)
163 {
164  transformToPolarCoordinates(src, srcImageRange(dest));
165 }
166 
167 
168 
169 
170 namespace detail
171 {
172  template <class SrcIterator, class SrcAccessor,
173  class DestIterator, class DestAccessor>
174  void
175  maximumFastNCC(SrcIterator s_ul, SrcIterator s_lr, SrcAccessor s_acc,
176  DestIterator d_ul, DestIterator d_lr, DestAccessor d_acc,
177  TinyVector<int,2> & position,
178  double & correlation_coefficent)
179  {
180  typename DestIterator::difference_type s_shape = s_lr - s_ul;
181  typename DestIterator::difference_type d_shape = d_lr - d_ul;
182 
183  MultiArray<2, float> src(s_shape.x, s_shape.y), dest(d_shape.x, d_shape.y), ncc(d_shape.x, d_shape.y);
184  BasicImageView<float> src_view = makeBasicImageView(src);
185  BasicImageView<float> dest_view = makeBasicImageView(dest);
186 
187  copyImage(srcIterRange(s_ul, s_lr, s_acc), destImage(src_view));
188  copyImage(srcIterRange(d_ul, d_lr, d_acc), destImage(dest_view));
189 
190  fastNormalizedCrossCorrelation(dest, src, ncc);
191 
192  int max_x = 0;
193  int max_y = 0;
194  float val = 0.0;
195  float max_val = -1.0;
196 
197  for (int y = 0; y < ncc.height()-s_shape.y; y++)
198  {
199  for (int x = 0; x < ncc.width()-s_shape.x; x++)
200  {
201  val = ncc(x+s_shape.x/2, y+s_shape.y/2);
202 
203  if (val > max_val)
204  {
205  max_x = x;
206  max_y = y;
207  max_val = val;
208  }
209  }
210  }
211 
212  position[0] = max_x;
213  position[1] = max_y;
214 
215  correlation_coefficent = max_val;
216  }
217 
218  template <class SrcIterator, class SrcAccessor,
219  class DestIterator, class DestAccessor>
220  inline void
221  maximumFastNCC(triple<SrcIterator, SrcIterator, SrcAccessor> src,
222  triple<DestIterator, DestIterator, DestAccessor> dest,
223  TinyVector<int,2> & position,
224  double & correlation_coefficent)
225  {
226  maximumFastNCC(src.first, src.second, src.third,
227  dest.first, dest.second, dest.third,
228  position,
229  correlation_coefficent);
230  }
231 
232  template <class SrcIterator, class SrcAccessor,
233  class DestIterator, class DestAccessor>
234  void fourierLogAbsSpectrumInPolarCoordinates(SrcIterator s_ul, SrcIterator s_lr, SrcAccessor s_acc,
235  DestIterator d_ul, DestIterator d_lr, DestAccessor d_acc)
236  {
237  using namespace vigra;
238 
239  typename SrcIterator::difference_type shape = s_lr - s_ul;
240 
241  FFTWComplexImage fourier(shape);
242  FImage fft_mag(shape);
243 
244  fourierTransform(srcIterRange(s_ul, s_lr, s_acc), destImage(fourier));
245  moveDCToCenter(srcImageRange(fourier, FFTWMagnitudeAccessor<>()), destImage(fft_mag));
246 
247  vigra::SplineImageView<2, double> spl(srcImageRange(fft_mag));
248 
250  destIterRange(d_ul, d_lr, d_acc));
251 
252  transformImage(srcIterRange(d_ul,d_lr,d_acc), destIter(d_ul,d_acc), log(abs(functor::Arg1())));
253  }
254 
255  template <class SrcIterator, class SrcAccessor,
256  class DestIterator, class DestAccessor>
257  void fourierLogAbsSpectrumInPolarCoordinates(triple<SrcIterator, SrcIterator, SrcAccessor> src,
258  triple<DestIterator, DestIterator, DestAccessor> dest)
259  {
260  fourierLogAbsSpectrumInPolarCoordinates(src.first, src.second, src.third, dest.first, dest.second, dest.third);
261  }
262 } //namespace detail
263 
264 
265 
266 
267 /********************************************************/
268 /* */
269 /* estimateGlobalRotation */
270 /* */
271 /********************************************************/
272 
273 /** \brief Estimate the rotation between two images by means of a normalized cross correlation matching of the FFT spectra.
274 
275  This algorithm uses the fast normalized cross correlation to determine a global rotation
276  between two images (from image2 to image1). To derive the rotation, the algorithm performs the following steps:
277  <ol>
278  <li>Transforming both images to the frequency domain using FFTW</li>
279  <li>Create LogAbs PolarCoordinate representations for each spectrum.</li>
280  <li>Determining the final Rotation by performing a fast normalized cross correlation
281  based on the polar representations.</li>
282  </ol>
283  The images are cropped to the corresponding images center-squared before the estimation
284  takes place.
285 
286  <b> Declarations:</b>
287 
288  <b>\#include</b> <vigra/affine_registration_fft.hxx><br>
289  Namespace: vigra
290 
291  pass 2D array views:
292  \code
293  namespace vigra {
294  template <class T1, class S1,
295  class T2, class S2>
296  void
297  estimateGlobalRotation(MultiArrayView<2, T1, S1> const & src,
298  MultiArrayView<2, T2, S2> dest,
299  Matrix<double> & affineMatrix,
300  double & correlation_coefficent);
301  }
302  \endcode
303 
304  \deprecatedAPI{estimateGlobalRotation}
305  pass \ref ImageIterators and \ref DataAccessors :
306  \code
307  namespace vigra {
308  template <class SrcIterator, class SrcAccessor,
309  class DestIterator, class DestAccessor>
310  void
311  estimateGlobalRotation(SrcIterator sul, SrcIterator slr, SrcAccessor src,
312  DestIterator dul, DestIterator dlr, DestAccessor dest,
313  Matrix<double> & affineMatrix,
314  double & correlation_coefficent)
315  }
316  \endcode
317  use argument objects in conjunction with \ref ArgumentObjectFactories :
318  \code
319  namespace vigra {
320  template <class SrcIterator, class SrcAccessor,
321  class DestIterator, class DestAccessor>
322  void
323  estimateGlobalRotation(triple<SrcIterator, SrcIterator, SrcAccessor> src,
324  triple<DestIterator, DestIterator, DestAccessor> dest,
325  Matrix<double> & affineMatrix,
326  double & correlation_coefficent)
327  }
328  \endcode
329  \deprecatedEnd
330 */
332 
333 template <class SrcIterator, class SrcAccessor,
334  class DestIterator, class DestAccessor>
335 void
336 estimateGlobalRotation(SrcIterator s_ul, SrcIterator s_lr, SrcAccessor s_acc,
337  DestIterator d_ul, DestIterator d_lr, DestAccessor d_acc,
338  Matrix<double> & affineMatrix,
339  double & correlation_coefficient)
340 {
341  //determine squared centers of both images without consuming any additional memory!
342  typename SrcIterator::difference_type s_shape = s_lr - s_ul;
343  Diff2D s2_shape(min(s_shape.x, s_shape.y),min(s_shape.x, s_shape.y));
344  Diff2D s2_offset = (s_shape-s2_shape)/2;
345 
346  typename DestIterator::difference_type d_shape = d_lr - d_ul;
347  Diff2D d2_shape(min(d_shape.x, d_shape.y),min(d_shape.x, d_shape.y));
348  Diff2D d2_offset = (d_shape-d2_shape)/2;
349 
350  //Determine Shape for united polar coordinate representation
351  Diff2D mean_shape = (s_shape + d_shape)/2;
352 
353  int size = min(mean_shape.x, mean_shape.y);
354  if(size %2 == 0)
355  size++;
356 
357  const int pc_w = size,
358  pc_h = size*6+1;
359 
360 
361  DImage s_polCoords(pc_w, pc_h/2),
362  d_polCoords(pc_w, pc_h),
363  ncc(pc_w, pc_h);
364 
365  detail::fourierLogAbsSpectrumInPolarCoordinates(srcIterRange(s_ul+s2_offset, s_ul+s2_offset+s2_shape, s_acc),
366  destImageRange(d_polCoords));
367  copyImage(srcIterRange(d_polCoords.upperLeft(), d_polCoords.upperLeft() + vigra::Diff2D(pc_w, pc_h/2), d_polCoords.accessor()),
368  destImage(s_polCoords));
369 
370  detail::fourierLogAbsSpectrumInPolarCoordinates(srcIterRange(d_ul+d2_offset, d_ul+d2_offset+d2_shape, d_acc),
371  destImageRange(d_polCoords));
372 
373  //Basic Cross correlation is assumed to be faster here, as there are only pc_h comparisons possible...
374  normalizedCrossCorrelation(srcImageRange(d_polCoords), srcImageRange(s_polCoords), destImage(ncc));
375 
376  int max_idx = 0;
377  double max_val = -1;
378 
379  const int x=pc_w/2;
380  double val;
381 
382  //Only look at a stripe for the maximum angle of rotation
383  //at the image center, at find the best fitting angle...
384  for (int y=0; y<pc_h/2; y++)
385  {
386  val = ncc(x,y+pc_h/4);
387 
388  if (val > max_val)
389  {
390  max_idx = y;
391  max_val = val;
392  }
393  }
394 
395  double theta = double(max_idx) / pc_h * 360.0;
396 
397  affineMatrix = rotationMatrix2DDegrees(theta, vigra::TinyVector<double,2>(s_shape.x/2.0, s_shape.y/2.0));
398  correlation_coefficient = max_val;
399 }
400 
401 template <class SrcIterator, class SrcAccessor,
402  class DestIterator, class DestAccessor>
403 inline void
404 estimateGlobalRotation(triple<SrcIterator, SrcIterator, SrcAccessor> src,
405  triple<DestIterator, DestIterator, DestAccessor> dest,
406  Matrix<double> & affineMatrix,
407  double & correlation_coefficient)
408 {
409  estimateGlobalRotation(src.first, src.second, src.third,
410  dest.first, dest.second, dest.third,
411  affineMatrix,
412  correlation_coefficient);
413 }
414 
415 template <class T1, class S1,
416  class T2, class S2>
417 inline void
420  Matrix<double> & affineMatrix,
421  double & correlation_coefficent)
422 {
423  estimateGlobalRotation(srcImageRange(src),
424  destImageRange(dest),
425  affineMatrix,
426  correlation_coefficent);
427 }
428 
429 /********************************************************/
430 /* */
431 /* estimateGlobalTranslation */
432 /* */
433 /********************************************************/
434 
435 /** \brief Estimate the translation between two images by means of a normalized cross correlation matching.
436 
437  This algorithm uses the fast normalized cross correlation to determine a global translation
438  between two images (from image2 to image1). To derive the translation, the algorithm consists of differents steps:
439  <ol>
440  <li>Separation of the second image<br/>
441  The second image (the one, for which the translation shall be determined) is cut into five
442  subregions: UpperLeft, UpperRight, Center, LowerLeft and LowerRight, each of 1/4 the size of
443  the original image. Using a border > 0 results in (all) overlapping regions.</li>
444  <li>Cross-Correlation of the subimages to the first image<br/>
445  The subimages are normalized cross-correlated to the (complete) first image.
446  The resulting maximum-likelihood translations and the correlation coefficients are stored for the next step.</li>
447  <li>Determining the final Translation by voting<br/>
448  Each correlation vector gets one vote at the beginning. For each equality of derived motion vectors,
449  the voting to these vectors is incremented. If the maximum number of votes is larger than 1, the vector with the
450  most votes is chosen. If the maximum number of votes is 1, the vector with the maximum likelihood is choosen.</li>
451  </ol>
452  <b> Declarations:</b>
453 
454  <b>\#include</b> <vigra/affine_registration_fft.hxx><br>
455  Namespace: vigra
456 
457  pass 2D array views:
458  \code
459  namespace vigra {
460  template <class T1, class S1,
461  class T2, class S2>
462  void
463  estimateGlobalTranslation(MultiArrayView<2, T1, S1> const & src,
464  MultiArrayView<2, T2, S2> dest,
465  Matrix<double> & affineMatrix,
466  double & correlation_coefficent,
467  Diff2D border = Diff2D(0,0));
468  }
469  \endcode
470 
471  \deprecatedAPI{estimateGlobalTranslation}
472  pass \ref ImageIterators and \ref DataAccessors :
473  \code
474  namespace vigra {
475  template <class SrcIterator, class SrcAccessor,
476  class DestIterator, class DestAccessor>
477  void
478  estimateGlobalTranslation(SrcIterator sul, SrcIterator slr, SrcAccessor src,
479  DestIterator dul, DestIterator dlr, DestAccessor dest,
480  Matrix<double> & affineMatrix,
481  double & correlation_coefficent,
482  Diff2D border = Diff2D(0,0))
483  }
484  \endcode
485  use argument objects in conjunction with \ref ArgumentObjectFactories :
486  \code
487  namespace vigra {
488  template <class SrcIterator, class SrcAccessor,
489  class DestIterator, class DestAccessor>
490  void
491  estimateGlobalTranslation(triple<SrcIterator, SrcIterator, SrcAccessor> src,
492  triple<DestIterator, DestIterator, DestAccessor> dest,
493  Matrix<double> & affineMatrix,
494  double & correlation_coefficent,
495  Diff2D border = Diff2D(0,0))
496  }
497  \endcode
498  \deprecatedEnd
499 */
501 
502 template <class SrcIterator, class SrcAccessor,
503  class DestIterator, class DestAccessor>
504 void estimateGlobalTranslation(SrcIterator s_ul, SrcIterator s_lr, SrcAccessor s_acc,
505  DestIterator d_ul, DestIterator d_lr, DestAccessor d_acc,
506  Matrix<double> & affine_matrix,
507  double & correlation_coefficent,
508  Diff2D border = Diff2D(0,0))
509 {
510  ignore_argument(d_lr);
511  typename SrcIterator::difference_type s_shape = s_lr - s_ul;
512 
513  //determine matrix by using 5 quater-matches and a maximum likelihood decision:
514  Diff2D q_shape = (s_shape - border - border)/2;
515  if (q_shape.x % 2 == 0) q_shape.x--;
516  if (q_shape.y % 2 == 0) q_shape.y--;
517 
518  Diff2D q_offsets[5];
519  q_offsets[0] = border;
520  q_offsets[1] = Diff2D(s_shape.x, 0)/2 + border;
521  q_offsets[2] = s_shape/4;
522  q_offsets[3] = Diff2D(0, s_shape.y)/2 + border;
523  q_offsets[4] = s_shape/2 + border;
524 
525  TinyVector<int,2> translation_vectors[5];
526  double translation_correlations[5] = {0.0,0.0,0.0,0.0,0.0};
527  int translation_votes[5] = {1,1,1,1,1};
528 
529  int max_corr_idx=0;
530 
531  for (int q=0; q!=5; q++)
532  {
533  Diff2D offset = q_offsets[q];
534 
535  //we are searching a transformation from img2 -> transformed image1, thus we switch dest and tmp
536  detail::maximumFastNCC(srcIterRange(d_ul+offset, d_ul+offset+q_shape, d_acc),
537  srcIterRange(s_ul, s_lr, s_acc),
538  translation_vectors[q],
539  translation_correlations[q]);
540 
541  translation_vectors[q] = translation_vectors[q] - TinyVector<int,2>(offset);
542 
543  if(translation_correlations[q] > translation_correlations[max_corr_idx])
544  {
545  max_corr_idx=q;
546  }
547 
548  for (int q_old=0; q_old!=q; q_old++)
549  {
550  if (translation_vectors[q] == translation_vectors[q_old])
551  {
552  translation_votes[q_old]++;
553  }
554  }
555  }
556 
557  int max_votes_idx=0;
558 
559  for (int q=0; q!=5; q++)
560  {
561  if(translation_votes[q] > translation_votes[max_votes_idx])
562  {
563  max_votes_idx=q;
564  }
565  }
566 
567  int best_idx=0;
568  if(translation_votes[max_votes_idx] > 1)
569  {
570  best_idx = max_votes_idx;
571  }
572  else
573  {
574  best_idx = max_corr_idx;
575  }
576 
577  affine_matrix = translationMatrix2D(translation_vectors[best_idx]);
578  correlation_coefficent = translation_correlations[best_idx];
579 }
580 
581 template <class SrcIterator, class SrcAccessor,
582  class DestIterator, class DestAccessor>
583 inline void
584 estimateGlobalTranslation(triple<SrcIterator, SrcIterator, SrcAccessor> src,
585  triple<DestIterator, DestIterator, DestAccessor> dest,
586  Matrix<double> & affineMatrix,
587  double & correlation_coefficent,
588  Diff2D border = Diff2D(0,0))
589 {
590  estimateGlobalTranslation(src.first, src.second, src.third,
591  dest.first, dest.second, dest.third,
592  affineMatrix,
593  correlation_coefficent,
594  border);
595 }
596 
597 template <class T1, class S1,
598  class T2, class S2>
599 inline void
602  Matrix<double> & affineMatrix,
603  double & correlation_coefficent,
604  Diff2D border = Diff2D(0,0))
605 {
606  estimateGlobalTranslation(srcImageRange(src),
607  destImageRange(dest),
608  affineMatrix,
609  correlation_coefficent,
610  border);
611 }
612 
613 /********************************************************/
614 /* */
615 /* estimateGlobalRotationTranslation */
616 /* */
617 /********************************************************/
618 
619 /** \brief Estimate the (global) rotation and translation between two images by means a normalized cross correlation matching.
620 
621  This algorithm use the functions \ref estimateGlobalRotation() and
622  \ref estimateGlobalTranslation() to estimate a matrix which describes
623  the global rotation and translation from the second to the first image.
624 
625  <b> Declarations:</b>
626 
627  <b>\#include</b> <vigra/affine_registration_fft.hxx><br>
628  Namespace: vigra
629 
630  pass 2D array views:
631  \code
632  namespace vigra {
633  template <class T1, class S1,
634  class T2, class S2>
635  void
636  estimateGlobalRotationTranslation(MultiArrayView<2, T1, S1> const & src,
637  MultiArrayView<2, T2, S2> dest,
638  Matrix<double> & affineMatrix,
639  double & rotation_correlation,
640  double & translation_correlation,
641  Diff2D border = Diff2D(0,0));
642  }
643  \endcode
644 
645  \deprecatedAPI{estimateGlobalRotationTranslation}
646  pass \ref ImageIterators and \ref DataAccessors :
647  \code
648  namespace vigra {
649  template <class SrcIterator, class SrcAccessor,
650  class DestIterator, class DestAccessor>
651  void
652  estimateGlobalRotationTranslation(SrcIterator sul, SrcIterator slr, SrcAccessor src,
653  DestIterator dul, DestIterator dlr, DestAccessor dest,
654  Matrix<double> & affineMatrix,
655  double & rotation_correlation,
656  double & translation_correlation,
657  Diff2D border = Diff2D(0,0))
658  }
659  \endcode
660  use argument objects in conjunction with \ref ArgumentObjectFactories :
661  \code
662  namespace vigra {
663  template <class SrcIterator, class SrcAccessor,
664  class DestIterator, class DestAccessor>
665  void
666  estimateGlobalRotationTranslation(triple<SrcIterator, SrcIterator, SrcAccessor> src,
667  triple<DestIterator, DestIterator, DestAccessor> dest,
668  Matrix<double> & affineMatrix,
669  double & rotation_correlation,
670  double & translation_correlation,
671  Diff2D border = Diff2D(0,0))
672  }
673  \endcode
674  \deprecatedEnd
675 */
677 template <class SrcIterator, class SrcAccessor,
678  class DestIterator, class DestAccessor>
679 void
680 estimateGlobalRotationTranslation(SrcIterator s_ul, SrcIterator s_lr, SrcAccessor s_acc,
681  DestIterator d_ul, DestIterator d_lr, DestAccessor d_acc,
682  Matrix<double> & affineMatrix,
683  double & rotation_correlation,
684  double & translation_correlation,
685  Diff2D border = Diff2D(0,0))
686 {
687  typename DestIterator::difference_type d_shape = d_lr - d_ul;
688 
689  //First step: Estimate rotation from img2 -> img1.
690  Matrix<double> rotation_matrix;
691  estimateGlobalRotation(srcIterRange(s_ul+border, s_lr-border, s_acc),
692  srcIterRange(d_ul+border, d_lr-border, d_acc),
693  rotation_matrix,
694  rotation_correlation);
695 
696  //Second step: correct image according to the estimated rotation:
697  FImage tmp(d_shape);
698  SplineImageView<3, double> spl(srcIterRange(s_ul, s_lr, s_acc));
699  affineWarpImage(spl, destImageRange(tmp), rotation_matrix);
700 
701  //Third step: find rotation between temp image (of step 2) and dest:
702  Matrix<double> translation_matrix;
703  estimateGlobalTranslation(srcImageRange(tmp),
704  srcIterRange(d_ul, d_lr, d_acc),
705  translation_matrix,
706  translation_correlation,
707  border);
708 
709  affineMatrix = rotation_matrix * translation_matrix;
710 }
711 
712 template <class SrcIterator, class SrcAccessor,
713  class DestIterator, class DestAccessor>
714 inline void
715 estimateGlobalRotationTranslation(triple<SrcIterator, SrcIterator, SrcAccessor> src,
716  triple<DestIterator, DestIterator, DestAccessor> dest,
717  Matrix<double> & affineMatrix,
718  double & rotation_correlation,
719  double & translation_correlation,
720  Diff2D border = Diff2D(0,0))
721 {
722  estimateGlobalRotationTranslation(src.first, src.second, src.third,
723  dest.first, dest.second, dest.third,
724  affineMatrix,
725  rotation_correlation,
726  translation_correlation,
727  border);
728 }
729 
730 template <class T1, class S1,
731  class T2, class S2>
732 inline void
735  Matrix<double> & affineMatrix,
736  double & rotation_correlation,
737  double & translation_correlation,
738  Diff2D border = Diff2D(0,0))
739 {
740  estimateGlobalRotationTranslation(srcImageRange(src),
741  destImageRange(dest),
742  affineMatrix,
743  rotation_correlation,
744  translation_correlation,
745  border);
746 }
747 
748 //@}
749 
750 } // namespace vigra
751 
752 
753 #endif /* VIGRA_AFFINE_REGISTRATION_FFT_HXX */
void fastNormalizedCrossCorrelation(...)
This function performes a fast normalized cross-correlation.
void moveDCToCenter(...)
Rearrange the quadrants of a Fourier image so that the origin is in the image center.
int y
Definition: diff2d.hxx:392
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.
Definition: affine_registration_fft.hxx:112
linalg::TemporaryMatrix< T > sin(MultiArrayView< 2, T, C > const &v)
void transformImage(...)
Apply unary point transformation to each pixel.
int x
Definition: diff2d.hxx:385
linalg::TemporaryMatrix< double > rotationMatrix2DDegrees(double angle)
Create homogeneous matrix representing a 2D rotation about the coordinate origin. ...
Definition: affinegeometry.hxx:133
Two dimensional difference vector.
Definition: diff2d.hxx:185
void estimateGlobalRotationTranslation(...)
Estimate the (global) rotation and translation between two images by means a normalized cross correla...
linalg::TemporaryMatrix< double > translationMatrix2D(TinyVector< double, 2 > const &shift)
Create homogeneous matrix representing a 2D translation.
Definition: affinegeometry.hxx:64
void estimateGlobalRotation(...)
Estimate the rotation between two images by means of a normalized cross correlation matching of the F...
BasicImageView< T > makeBasicImageView(MultiArrayView< 2, T, Stride > const &array)
Definition: multi_array.hxx:3476
doxygen_overloaded_function(template<...> void separableConvolveBlockwise) template< unsigned int N
Separated convolution on ChunkedArrays.
void estimateGlobalTranslation(...)
Estimate the translation between two images by means of a normalized cross correlation matching...
void copyImage(...)
Copy source image into destination image.
Definition: fftw3.hxx:1377
linalg::TemporaryMatrix< T > log(MultiArrayView< 2, T, C > const &v)
Fundamental class template for images.
Definition: basicimage.hxx:475
FFTWComplex< R >::NormType abs(const FFTWComplex< R > &a)
absolute value (= magnitude)
Definition: fftw3.hxx:1002
Base class for, and view to, vigra::MultiArray.
Definition: multi_array.hxx:704
void affineWarpImage(...)
Warp an image according to an affine transformation.
Create a continuous view onto a discrete image using splines.
Definition: splineimageview.hxx:99
linalg::TemporaryMatrix< T > cos(MultiArrayView< 2, T, C > const &v)
void normalizedCrossCorrelation(...)
This function performes a (slow) normalized cross-correlation.
void fourierTransform(...)
Compute forward and inverse Fourier transforms.

© 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)