fftw中虚部为零的实数fft和复数fft之间的区别?

Difference between real fft and complex fft with imaginary part of zero in fftw?

我有一个真正的二维矩阵。我正在使用 fftw 获取它的 fft。但是使用实数到复数 fft 的结果不同于复数(虚部等于零)到复数 fft。

实矩阵

 0     1     2
 3     4     5
 6     7     8

实数到复数 fft 的结果

36 -4.5+2.59808i  -13.5+7.79423i 
0  -13.5-7.79423i 0 
0  0              0 

代码:

int r = 3, c = 3;
int sz = r * c;
double *in = (double*) malloc(sizeof(double) * sz);
fftw_complex *out = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * sz);
fftw_plan p = fftw_plan_dft_r2c_2d(r, c, in, out, FFTW_MEASURE);
for ( int i=0; i<r; ++i ){
    for ( int j=0; j<c; ++j ){
        in[i*c+j] = i*c + j;
    }
}
fftw_execute(p);

使用虚部为零的复数矩阵

复杂矩阵

 0+0i     1+0i     2+0i
 3+0i     4+0i     5+0i
 6+0i     7+0i     8+0i

复数到复数 fft 的结果

36               -4.5 + 2.59808i  -4.5 - 2.59808i 
-13.5 + 7.79423i 0               0 
-13.5 - 7.79423i 0               0  

代码:

int r = 3, c = 3;
int sz = r * c;
fftw_complex *out = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * sz);
fftw_complex *inc = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * sz);
p = fftw_plan_dft_2d( r,c, inc, out, FFTW_FORWARD,FFTW_MEASURE);
for ( int i=0; i<r; ++i ){
    for ( int j=0; j<c; ++j ){
        inc[i*c+j][0] = i*c+j;
        inc[i*c+j][1] = 0;
    }
}
fftw_execute(p);

我正在寻找复杂到复杂的 fft 的结果。但是从真实到复杂的 fft 速度要快得多,而且我的数据是真实的。我是否犯了编程错误或结果应该不同?

FFTW documentation

所示

Then, after an r2c transform, the output is an n0 × n1 × n2 × … × (nd-1/2 + 1) array of fftw_complex values in row-major order

换句话说,样本实矩阵的实数到复数变换的输出实际上是:

 36            -4.5+2.59808i
-13.5+7.79423i  0
-13.5-7.79423i  0 

您可能会注意到这两列与复杂到复杂转换的前两列完全匹配。缺失的列从实数到复数的变换中被省略,因为它由于对称性而冗余。因此,包含缺失列的完整 3x3 矩阵可以使用以下方法构建:

fftw_complex *outfull = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * sz);
int outc = (c/2+1);
for ( int i=0; i<r; ++i ){
    // copy existing columns
    for ( int j=0; j<outc; ++j ){
        outfull[i*c+j][0] = out[i*outc+j][0];
        outfull[i*c+j][1] = out[i*outc+j][1];
    }
    // generate missing column(s) from symmetry
    for ( int j=outc; j<c; ++j){
        int row = (r-i)%r;
        int col = c-j;
        outfull[i*c+j][0] =  out[row*outc+col][0];
        outfull[i*c+j][1] = -out[row*outc+col][1];
    }
}