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.