5

I'm using FFTW3 to compute 2D real FFT in c++. I've read the manual but have some questions. From the manual: http://www.fftw.org/fftw3_doc/One_002dDimensional-DFTs-of-Real-Data.html#One_002dDimensional-DFTs-of-Real-Data

In exchange for these speed and space advantages, the user sacrifices some of the simplicity of FFTW's complex transforms. First of all, the input and output arrays are of different sizes and types: the input is n real numbers, while the output is n/2+1 complex numbers (the non-redundant outputs); this also requires slight “padding” of the input array for in-place transforms. Second, the inverse transform (complex to real) has the side-effect of overwriting its input array, by default. Neither of these inconveniences should pose a serious problem for users, but it is important to be aware of them.

  1. I understand that I need to convert my input 2D matrix in into row-order 1D vector. But what does the output look like? What do the n/2 + 1 numbers mean? In other words, how do I reorder the output to get 2D matrix?

  2. What specifically do I have to do to create this "padding"?

tir38
  • 9,810
  • 10
  • 64
  • 107
  • Hm. In some other FFT libraries I've used, the Nyquist value is just packed into the imaginary component of `X_0` (the DC offset) since both values are real. Then the padding is not necessary. Evidently FFTW doesn't do this. – nneonneo Jan 28 '13 at 23:32
  • Just for reference, `such packing [..] does not generalize well to multi-dimensional transforms`. See the second last paragraph of [this FFTW doc](http://www.fftw.org/doc/One_002dDimensional-DFTs-of-Real-Data.html#One_002dDimensional-DFTs-of-Real-Data) for details. – Sveltely Mar 24 '14 at 08:09

2 Answers2

2

  1. If your input is already in a normal C++ 2D array, all you should need to do is typecast it:

    double twoDarray[10][10];
double *oneDarrayPointer = (double *)twoDarray;

    If your input is 100 (like it is in the above example), your output array is going to be 51 complex numbers. The format of those numbers should be described by your library, but probably is an array of 102 doubles - 51 entries times 2 (real/imaginary parts).

    Edit: Confirmed - fftw_complex is defined as:

    typedef double fftw_complex[2];

    So they're just consecutive pairs of doubles representing the real and imaginary parts of a complex number.

  2. If you don't want to do it in place, you don't have to pad anything - just allocate an appropriately sized output array. If you do need to do it in place, your input buffer has to have space for the 2 extra doubles vs the input size. Assuming the declarations above, you'd want something like:

    double *inPlaceFFTPointer = malloc(sizeof twoDarray + 2*sizeof(double));
memcpy(inPlaceFFTPointer, oneDarrayPointer, sizeof twoDarray);

    I'm not sure if you'd need to make sure to have 0.0 in the last two entries or not, but that's easy enough to add.

Carl Norum
  • 219,201
  • 40
  • 422
  • 469
  • 1
    you can use C++ std::complex as well as the C style or the above typedef with FFTW – pyCthon Jan 29 '13 at 03:48
  • @pyCthon, yup. I think they're pretty much all the same in-memory. – Carl Norum Jan 29 '13 at 04:27
  • "So they're just consecutive pairs of doubles representing the real and imaginary parts of a complex number." Right, but why only 51? Why not 100? Doesn't 2D FFT output the same size matrix as input? – tir38 Jan 29 '13 at 18:29
  • 2
    From the link in your post: "In many practical applications, the input data in[i] are purely real numbers, in which case the DFT output satisfies the “Hermitian” redundancy: `out[i]` is the conjugate of `out[n-i]`." So you only need to have half of the output in order to generate all of the output should you so wish. – Carl Norum Jan 29 '13 at 18:30
  • I kind of understand. This helped me understand the visual symmetry: "You may begin to notice there is a lot of symmetry. For all REAL (as opposed to IMAGINARY or COMPLEX) images, the FT is symmetrical about the origin so the 1st and 3rd quadrants are the same and the 2nd and 4th quadrants are the same. If the image is symmetrical about the x-axis (as the cosine images are) 4-fold symmetry results." http://www.cs.unm.edu/~brayer/vision/fourier.html – tir38 Feb 27 '13 at 21:54
2

You could have a look to the real-to-real transforms in FFTW3, that do exactly what you were asking for. These do not require padding and take into account that both the wavenumber 0 and the one that represents the Nyquist frequency have only a real component. Have a look here:

FFTW3 Real-to-Real Transforms

and for the layout in memory:

FFTW3 Real-to-Real Transform Kinds

Chiel
  • 6,006
  • 2
  • 32
  • 57