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

details Helper Functions for FFT VIGRA


Functions

template<...> void moveDCToCenter (SrcImageIterator src_upperleft, SrcImageIterator src_lowerright, SrcAccessor sa, DestImageIterator dest_upperleft, DestAccessor da)
 Rearrange the quadrants of a Fourier image so that the origin is in the image center.

template<class SrcImageIterator, class SrcAccessor, class DestImageIterator, class DestAccessor> void moveDCToUpperLeft (SrcImageIterator src_upperleft, SrcImageIterator src_lowerright, SrcAccessor sa, DestImageIterator dest_upperleft, DestAccessor da)
 Rearrange the quadrants of a Fourier image so that the origin is in the image's upper left.

template<...> void applyFourierFilter (triple< SrcImageIterator, SrcImageIterator, SrcAccessor > src, pair< FilterImageIterator, FilterAccessor > filter, pair< DestImageIterator, DestAccessor > dest)
 Apply a filter (defined in the frequency domain) to an image.

template<...> void applyFourierFilterFamily (triple< SrcImageIterator, SrcImageIterator, SrcAccessor > src, const ImageArray< FilterType > &filters, ImageArray< DestImage > &results)
 Apply an array of filters (defined in the frequency domain) to an image.



Detailed Description


Rearrange quadrants in a Fourier image or apply filters on images using the FFTW.


Function Documentation


  void applyFourierFilter (...)
 
 

Apply a filter (defined in the frequency domain) to an image.

After transferring the image into the frequency domain, it is multiplied pixel-wise with the filter and transformed back. The result is always put into the given FFTWComplexImage destImg which must have the right size. Finally, the result will be normalized to compensate for the two FFTs.

The input and filter images can be scalar images (more precisely, SrcAccessor::value_type must be scalar) or FFTWComplexImages (in this case, one conversion is saved).

The DC entry of the filter must be in the upper left, which is the position where FFTW expects it (see moveDCToUpperLeft()).

You can optionally pass two fftwnd_plans as last parameters, the forward plan and the in-place backward plan. Both must have been created for the right image size (see sample code).

Declarations:

pass arguments explicitly:

    namespace vigra {
        template <class SrcImageIterator, class SrcAccessor,
            class FilterImageIterator, class FilterAccessor>
        void applyFourierFilter(SrcImageIterator srcUpperLeft,
                                SrcImageIterator srcLowerRight, SrcAccessor sa,
                                FilterImageIterator filterUpperLeft, FilterAccessor fa,
                                FFTWComplexImage & destImg)

        template <class SrcImageIterator, class SrcAccessor,
            class FilterImageIterator, class FilterAccessor>
        void applyFourierFilter(SrcImageIterator srcUpperLeft,
                                SrcImageIterator srcLowerRight, SrcAccessor sa,
                                FilterImageIterator filterUpperLeft, FilterAccessor fa,
                                FFTWComplexImage & destImg,
                                const fftwnd_plan & forwardPlan, const fftwnd_plan & backwardPlan)
    }

use argument objects in conjuction with Argument Object Factories:

    namespace vigra {
        template <class SrcImageIterator, class SrcAccessor,
            class FilterImageIterator, class FilterAccessor>
        inline
        void applyFourierFilter(triple<SrcImageIterator, SrcImageIterator, SrcAccessor> src,
                                pair<FilterImageIterator, FilterAccessor> filter,
                                FFTWComplexImage &destImg)

        template <class SrcImageIterator, class SrcAccessor,
            class FilterImageIterator, class FilterAccessor>
        inline
        void applyFourierFilter(triple<SrcImageIterator, SrcImageIterator, SrcAccessor> src,
                                pair<FilterImageIterator, FilterAccessor> filter,
                                FFTWComplexImage &destImg,
                                const fftwnd_plan &forwardPlan, const fftwnd_plan &backwardPlan)
    }

Usage:

#include "vigra/fftw.hxx"
Namespace: vigra

    // create a Gaussian filter in Fourier space
    vigra::FImage gaussFilter(w, h), filter(w, h);
    for(int y=0; y<h; ++y)
        for(int x=0; x<w; ++x)
        {
            xx = float(x - w / 2) / w;
            yy = float(y - h / 2) / h;

            gaussFilter(x,y) = std::exp(-(xx*xx + yy*yy) / 2.0 * scale);
        }

    // applyFourierFilter() expects the filter's DC in the upper left
    moveDCToUpperLeft(srcImageRange(gaussFilter), destImage(filter));

    vigra::FFTWComplexImage result(w, h);

    vigra::applyFourierFilter(srcImageRange(image), srcImage(filter), result);

For inspection of the result, FFTWMagnitudeAccessor might be useful.

As mentioned, you can pass along two FFTW-plans. This is useful to speed up the application of many filters of the same size. Look into the FFTW documentation for details. You can create optimized plans like this:

    // note that the height comes before the width here!
    fftwnd_plan forwardPlan= fftw2d_create_plan(h, w, FFTW_FORWARD, FFTW_MEASURE );
    fftwnd_plan backwardPlan= fftw2d_create_plan(h, w, FFTW_BACKWARD, FFTW_MEASURE | FFTW_IN_PLACE);

    applyFourierFilter(srcImageRange(image), srcImage(filter), result,
                       forwardPlan, backwardPlan);

Note: The creation of the plans in this way can take several seconds - the FFTW library will measure different possible implementations to decide which one is the fastest for this specific environment. The result of those measurements is called "wisdom" in the FFTW slang and there are functions to save it to disk.


  void applyFourierFilterFamily (...)
 
 

Apply an array of filters (defined in the frequency domain) to an image.

This provides the same functionality as applyFourierFilter(), but applying several filters at once allows to avoid repeated Fourier transforms of the source image.

Filters and result images must be stored in vigra::ImageArray data structures. In contrast to applyFourierFilter(), this function adjusts the size of the result images and the the length of the array.

Declarations:

pass arguments explicitly:

    namespace vigra {
        template <class SrcImageIterator, class SrcAccessor, class FilterType>
        void applyFourierFilterFamily(SrcImageIterator srcUpperLeft,
                                      SrcImageIterator srcLowerRight, SrcAccessor sa,
                                      const ImageArray<FilterType> &filters,
                                      ImageArray<FFTWComplexImage> &results)

        template <class SrcImageIterator, class SrcAccessor, class FilterType>
        void applyFourierFilterFamily(SrcImageIterator srcUpperLeft,
                                      SrcImageIterator srcLowerRight, SrcAccessor sa,
                                      const ImageArray<FilterType> &filters,
                                      ImageArray<FFTWComplexImage> &results,
                                      const fftwnd_plan &forwardPlan, const fftwnd_plan &backwardPlan)
    }

use argument objects in conjuction with Argument Object Factories:

    namespace vigra {
        template <class SrcImageIterator, class SrcAccessor, class FilterType>
        inline
        void applyFourierFilterFamily(triple<SrcImageIterator, SrcImageIterator, SrcAccessor> src,
                                      const ImageArray<FilterType> &filters,
                                      ImageArray<FFTWComplexImage> &results)

        template <class SrcImageIterator, class SrcAccessor, class FilterType>
        inline
        void applyFourierFilterFamily(triple<SrcImageIterator, SrcImageIterator, SrcAccessor> src,
                                      const ImageArray<FilterType> &filters,
                                      ImageArray<FFTWComplexImage> &results,
                                      const fftwnd_plan &forwardPlan, const fftwnd_plan &backwardPlan)
    }

Usage:

#include "vigra/fftw.hxx"
Namespace: vigra

    // assuming the presence of a real-valued image named "image" to
    // be filtered in this example

    vigra::ImageArray<vigra::FImage> filters(16, image.size());

    for (int i=0; i<filters.size(); i++)
         // create some meaningful filters here
         createMyFilterOfScale(i, destImage(filters[i]));

    vigra::ImageArray<vigra::FFTWComplexImage> results();

    vigra::applyFourierFilterFamily(srcImageRange(image), filters, results);

For details about the FFTW-plans, see the applyFourierFilter() example.


  void moveDCToCenter (...)
 
 

Rearrange the quadrants of a Fourier image so that the origin is in the image center.

FFTW produces fourier images where the DC component (origin of fourier space) is located in the upper left corner of the image. The quadrants are placed like this (using a 4x4 image for example):

            DC 4 3 3
             4 4 3 3
             1 1 2 2
             1 1 2 2

After applying the function, the quadrants are at their usual places:

            2 2  1 1
            2 2  1 1
            3 3 DC 4
            3 3  4 4

This transformation can be reversed by moveDCToUpperLeft(). Note that the transformation must not be executed in place - input and output images must be different.

Declarations:

pass arguments explicitly:

    namespace vigra {
        template <class SrcImageIterator, class SrcAccessor,
          class DestImageIterator, class DestAccessor>
        void moveDCToCenter(SrcImageIterator src_upperleft,
                               SrcImageIterator src_lowerright, SrcAccessor sa,
                               DestImageIterator dest_upperleft, DestAccessor da)
    }

use argument objects in conjuction with Argument Object Factories:

    namespace vigra {
        template <class SrcImageIterator, class SrcAccessor,
                  class DestImageIterator, class DestAccessor>
        inline void moveDCToCenter(
            triple<SrcImageIterator, SrcImageIterator, SrcAccessor> src,
            pair<DestImageIterator, DestAccessor> dest)
    }

Usage:

#include "vigra/fftw.hxx"
Namespace: vigra

    vigra::FFTWComplexImage spatial(width,height), fourier(width,height);
    ... // fill image with data

    // create a plan for optimal performance
    fftwnd_plan forwardPlan=
        fftw2d_create_plan(height, width, FFTW_FORWARD, FFTW_ESTIMATE );

    // calculate FFT
    fftwnd_one(forwardPlan, spatial.begin(), fourier.begin());

    vigra::FFTWComplexImage rearrangedFourier(width, height);
    moveDCToCenter(srcImageRange(fourier), destImage(rearrangedFourier));

    //delete the plan
    fftwnd_destroy_plan(forwardPlan);


void moveDCToUpperLeft SrcImageIterator    src_upperleft,
SrcImageIterator    src_lowerright,
SrcAccessor    sa,
DestImageIterator    dest_upperleft,
DestAccessor    da

 

Rearrange the quadrants of a Fourier image so that the origin is in the image's upper left.

This function is the inversion of moveDCToCenter(). See there for declarations and a usage example.

© Ullrich Köthe (koethe@informatik.uni-hamburg.de)
Cognitive Systems Group, University of Hamburg, Germany

html generated using doxygen and Python
VIGRA 1.2.0 (7 Aug 2003)