Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

FFTW - computing real 2D FFT, special requirements

Tags:

c++

fft

fftw

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"?

like image 428
tir38 Avatar asked Oct 05 '22 05:10

tir38


2 Answers

  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.

like image 144
Carl Norum Avatar answered Oct 07 '22 19:10

Carl Norum


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

like image 42
Chiel Avatar answered Oct 07 '22 19:10

Chiel