36#ifndef VIGRA_CORRELATION_HXX
37#define VIGRA_CORRELATION_HXX
39#include "stdimage.hxx"
40#include "inspectimage.hxx"
41#include "functorexpression.hxx"
42#include "multi_math.hxx"
43#include "multi_fft.hxx"
44#include "integral_image.hxx"
47#include "applywindowfunction.hxx"
116template <
class T1,
class S1,
125 vigra_precondition(
in.shape() ==
out.shape(),
126 "vigra::fastCrossCorrelation(): shape mismatch between input and output.");
139template<
class T1,
class S1>
140inline double integralMultiArrayWindowMean(MultiArrayView<1, T1, S1>
const & in,
144 return in[right]-in[left];
147template<
class T1,
class S1>
148inline double integralMultiArrayWindowMean(MultiArrayView<2, T1, S1>
const & in,
152 return in[lr] - in(lr[0],ul[1]) - in(ul[0],lr[1]) + in[ul];
155template<
class T1,
class S1>
156inline double integralMultiArrayWindowMean(MultiArrayView<3, T1, S1>
const & in,
160 return (in[lr] - in(lr[0],ul[1],lr[2]) - in(ul[0],lr[1],lr[2]) + in(ul[0],ul[1],lr[2]))
161 - (in(lr[0],lr[1],ul[2]) - in(lr[0],ul[1],ul[2]) - in(ul[0],lr[1],ul[2]) + in[ul]);
223template <
class T1,
class S1,
235 vigra_precondition(
in.shape() ==
out.shape(),
236 "vigra::fastNormalizedCrossCorrelation(): shape mismatch between input and output.");
238 vigra_precondition(N>0 && N<=3,
239 "vigra::fastNormalizedCrossCorrelation(): only implemented for arrays of 1, 2 or 3 dimensions.");
241 for(
unsigned int dim=0; dim<N; dim++)
243 vigra_precondition(mask.shape()[dim] % 2 == 1,
"vigra::fastNormalizedCrossCorrelation(): Mask width has to be of odd size!");
244 vigra_precondition(
in.shape()[dim] >= mask.shape()[dim] ,
"vigra::fastNormalizedCrossCorrelation(): Mask is larger than image!");
308template<
class MaskIterator,
class MaskAccessor>
309class CorrelationFunctor
320 : m_mask_ul(
m.first),
326 template <
class SrcIterator,
class SrcAccessor,
class DestIterator,
class DestAccessor>
329 using namespace vigra;
356 Diff2D windowShape()
const
358 return m_mask_lr - m_mask_ul;
453 BorderTreatmentMode
border = BORDER_TREATMENT_AVOID)
465 BorderTreatmentMode
border = BORDER_TREATMENT_AVOID)
468 mask.first, mask.second, mask.third,
473template <
class T1,
class S1,
480 vigra_precondition(
in.shape() ==
out.shape(),
481 "vigra::crossCorrelation(): shape mismatch between input and output.");
487template<
class MaskIterator,
class MaskAccessor>
488class NormalizedCorrelationFunctor
503 : m_mask_ul(mask.first),
504 m_mask_lr(mask.second),
505 m_mask_acc(mask.third),
512 template <
class SrcIterator,
class SrcAccessor,
class DestIterator,
class DestAccessor>
515 using namespace vigra;
536 double s1=0,s2=0,
s12=0,
s22=0;
544 s12 += (s1-m_avg1)*(s2-average());
545 s22 += (s2-average())*(s2-average());
559 Diff2D windowShape()
const
561 return m_mask_lr - m_mask_ul;
569 inspectImage(srcIterRange(m_mask_ul, m_mask_lr, m_mask_acc), average);
576 for( ;
ym.y != m_mask_lr.y;
ym.y++)
578 for(
xm =
ym;
xm.x != m_mask_lr.x;
xm.x++)
580 m_s11 += (m_mask_acc(
xm)-m_avg1)*(m_mask_acc(
xm)-m_avg1);
680 BorderTreatmentMode
border = BORDER_TREATMENT_AVOID)
692 BorderTreatmentMode
border = BORDER_TREATMENT_AVOID)
695 mask.first, mask.second, mask.third,
700template <
class T1,
class S1,
707 vigra_precondition(
in.shape() ==
out.shape(),
708 "vigra::normalizedCrossCorrelation(): shape mismatch between input and output.");
Two dimensional difference vector.
Definition diff2d.hxx:186
MultiArrayShape< actual_dimension >::type difference_type
Definition multi_array.hxx:739
Class for a single RGB value.
Definition rgbvalue.hxx:128
TinyVectorView< VALUETYPE, TO-FROM > subarray() const
Definition tinyvector.hxx:887
Class for fixed size vectors.
Definition tinyvector.hxx:1008
int floor(FixedPoint< IntBits, FracBits > v)
rounding down.
Definition fixedpoint.hxx:667
void fastNormalizedCrossCorrelation(...)
This function performes a fast normalized cross-correlation.
void crossCorrelation(...)
This function performes a (slow) cross-correlation.
void correlateFFT(...)
Correlate an array with a kernel by means of the Fourier transform.
void inspectImage(...)
Apply read-only functor to every pixel in the image.
void fastCrossCorrelation(...)
This function performes a fast cross-correlation.
void applyWindowFunction(...)
Apply a window function to each pixels of a given image.
NumericTraits< V >::Promote prod(TinyVectorBase< V, SIZE, D1, D2 > const &l)
product of the vector's elements
Definition tinyvector.hxx:2097
void normalizedCrossCorrelation(...)
This function performes a (slow) normalized cross-correlation.
void initMultiArrayBorder(...)
Write values to the specified border values in the array.