0

I'm trying to figure out, how can I calculate fft row by row of the following matrix (stored as an array) using fftw3:

double signal[9] = {
        1, 2, 3, 
        4, 5, 6, 
        7, 8, 9  
}

to obtain the matrix with rows (pseudo-code)

    fft([signal[0], signal[1], signal[2]])
    fft([signal[3], signal[4], signal[5]])
    fft([signal[6], signal[7], signal[8]])

I can check the result using Python just applying np.fft.fft([[1, 2, 3], [4, 5, 6], [7, 8, 9]], axis=1) and I want to get the same result with fftw3. My code:


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


int main() 
{
    // it is needed to compute fft row-by-row of the following matrix
    double signal[9] { 
        1, 2, 3, 
        4, 5, 6, 
        7, 8, 9  
    }; // 3x3 matrix
    // so that the output will have the form (pseudo-code)
    // fft(signal[0], signal[1], signal[2])
    // fft(signal[3], signal[4], signal[5])
    // fft(signal[6], signal[7], signal[8])
    fftw_complex result[9]{  };

    // prepare parameters
    int rank = 1;
    int n[] = { 3 };
    int howmany = 3;
    int idist = 3;
    int odist = 3;
    int istride = 1;
    int ostride = 1;

    auto plan_ = fftw_plan_many_dft_r2c(
        rank,       // 1D problem
        n,          // 1D transforms of length n[0] = 3
        howmany,    // 3 1D transforms
        signal,     // input matrix (array)
        NULL,       // will be assigned to n
        istride,    // distance between two elements in the same row
        idist,      // distance between the first elements of neighboring rows
        result,     // output
        NULL,       // will be assigned to n
        ostride,    // distance between two elements in the same row
        odist,      // distance between the firsn elements of neighboring rows
        FFTW_ESTIMATE
    );

    fftw_execute(plan_);
    for (int i = 0; i < 3; ++i)
    {
        for (int j = 0; j < 3; ++j)
        {
            std::cout << result[i * 3 + j][0] << " + 1j*" << result[i * 3 + j][1] << '\t';
        }
        std::cout << '\n';
    }
    fftw_destroy_plan(plan_);

    return 0;
}

The corresponding explanations are in the code. Python gives the following answer:

array([[ 6. +0.j       , -1.5+0.8660254j, -1.5-0.8660254j],
       [15. +0.j       , -1.5+0.8660254j, -1.5-0.8660254j],
       [24. +0.j       , -1.5+0.8660254j, -1.5-0.8660254j]])

My code gives the following

6 + 1j*0        -1.5 + 1j*0.866025      0 + 1j*0
15 + 1j*0       -1.5 + 1j*0.866025      0 + 1j*0
24 + 1j*0       -1.5 + 1j*0.866025      0 + 1j*0

I see that the output almost the same as with Python, but third elements of each row are zeros. Can someone help me?

Thank you.

  • Hmmm, `int one_dimensional_index = (row * MAXIMUM_COLUMNS) + column;` I'm sure you could use this in your indexing (or derive the row and column based on a 1D index). – Thomas Matthews Nov 11 '20 at 01:21
  • @ThomasMatthews, corrected. Many thanx! It does not affect the result, but it's good for cout. –  Nov 11 '20 at 01:28
  • 1
    the posted code is C++ not C. They are two different languages. Please remove the `c` tag – user3629249 Nov 12 '20 at 14:01

0 Answers0