FFT2D.hh

/* 
This file is part of the BIAS library (Basic ImageAlgorithmS).

Copyright (C) 2003-2009    (see file CONTACT for details)
Multimediale Systeme der Informationsverarbeitung
Institut fuer Informatik
Christian-Albrechts-Universitaet Kiel


BIAS is free software; you can redistribute it and/or modify
it under the terms of the GNU Lesser General Public License as published by
the Free Software Foundation; either version 2.1 of the License, or
(at your option) any later version.

BIAS is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
GNU Lesser General Public License for more details.

You should have received a copy of the GNU Lesser General Public License
along with BIAS; if not, write to the Free Software
Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
*/


#ifndef __2DFFT_hh__
#define __2DFFT_hh__

#include <bias_config.h>
#ifndef BIAS_HAVE_FFTW3
#  error You need BIAS_HAVE_FFTW3 to use FFT2D. Please enable USE_FFTW3 and recompile BIAS.
#endif

#include "Filter/FilterNTo2N.hh"

#include <complex>
#include <iostream>
#include <fftw3.h>

namespace BIAS {

  /** \class FFT2D
   * \brief Wrapper to the fftw3 library adapted for 2D image filtering
   * \note Use fast as follows:
   * - call Init() with appropriate size
   * - call Filter()
   * \author evers and koeser, jw, oniemann
   * \date 2004/01, 2006 */
  template <class InputStorageType, class OutputStorageType>
  class BIASImage_EXPORT  FFT2D   /*: public FilterNTo2N<InputStorageType, OutputStorageType>*/
  {
  public:
    FFT2D();
    virtual ~FFT2D();

    /** @brief initializes for forward transformation on complete image. 
    This may take some time 
    (i.e. several seconds) since some measurements are done. */
    void Init(int width, int height);

    /** call this before a second call to Init() */
    void Release();    

    /** dst.GetChannelCount()==2*src.GetCHannelCount() 
    The result is not normalized! */
    virtual int Filter(const Image<InputStorageType>& src, 
      Image<OutputStorageType>& dst);

    /** dstX.GetChannelCount()==src.GetCHannelCount() 
    The result is not normalized! */
    virtual int Filter(const Image<InputStorageType>& src, 
      Image<OutputStorageType>& dst1, 
      Image<OutputStorageType>& dst2);

    /** Transform forward and get absolute value from complex result
    The result is not normalized! */
    virtual int TransformAbs(const Image<InputStorageType>& src, 
      Image<OutputStorageType>& dst);
    virtual int TransformLogAbs(const Image<InputStorageType>& src, 
      Image<OutputStorageType>& dst);
    virtual int TransformPhase(const Image<InputStorageType>& src, 
      Image<OutputStorageType>& dst);

    /** Transform reverse, src _must_ be of _SizeX,_OutSizeY,2
    */
    virtual int TransformReverse(const Image<OutputStorageType>& src,
      Image<InputStorageType>& dst);

    /** Normalize output after transformation */
    virtual int Normalize(Image<OutputStorageType>& dst);

    virtual int CrossPowerSpectrum(const Image<InputStorageType>& src1,
      const Image<InputStorageType>& src2,
      Image<OutputStorageType>& dst1);


  protected:

    virtual int Forward_(const Image<InputStorageType>& src);

    // assumes src.GetChannelCount() == 2
    virtual int Reverse_(const Image<OutputStorageType>& src);

    virtual void GetBordersValid_(int& border_x, int& border_y) const;

    fftw_complex *_out;
    double *_in;
    fftw_plan _p_forward, _p_reverse;
    // size of input _Size = SizeX * _SizeY
    // Warning: _SizeX is height and _SizeY is width of image
    int _SizeX, _SizeY, _Size;
    int _OutSizeY;
    int _OutNum; // length of output = _SizeX * (_SizeY/2 + 1)
  };

} // namespace

#endif // __FFT2D_hh__