问题描述
我正在通过在16x16数据集中设置某些模式的实数分量来手动初始化2d频域中的状态。然后,我执行2D IDFT来获取真实域数据。这一切都按预期工作。
然后,我对实域数据执行DFT,以获取(应该是)与我手动初始化的频率模式相同的频率模式。但是,它们以振幅减半返回,并且垂直频率“镜像”。为了说明:
输入模式:
k [1,0]:32 + 0i
k [2,0]:16 + 0i
k [3,0]: 8 + 0i
k [4,0]: 4 + 0i
IDFT之后的输出模式-> DFT:
k [1,0]:16 + 0i
k [2,0]: 8 + 0i
k [3,0]: 4 + 0i
k [4,0]: 2 + 0i
k [12,0]: 2 + 0i
k [13,0]: 4 + 0i
k [14,0]: 8 + 0i
k [15,0]:16 + 0i
我的问题是,为什么DFT输出中的模式与IDFT的初始输入不同?
在某些情况下,我遇到的问题是我正在使用这些数据来“求解”热方程,并且高频信号会迅速缩小。因此,经过数个时间步长后,k [12,0]至k [15,0]模式实际上贡献不大。
重现问题的代码:
int N = 16; // Dimensions of the data
int logical_width = (N / 2) + 1; // Logical width of the frequency domain
double* real = new double[N * N];
fftw_complex* complex = (fftw_complex*)fftw_malloc(sizeof(fftw_complex) * N * logical_width);
fftw_plan plan = fftw_plan_dft_r2c_2d(N,N,real,complex,FFTW_ESTIMATE);
fftw_plan iplan = fftw_plan_dft_c2r_2d(N,FFTW_ESTIMATE);
// Initialize all real data to 0
for (int i = 0; i < N * N; i++) {
real[i] = 0.0;
}
// Initialize all complex data to 0
for (int i = 0; i < N * logical_width; i++) {
complex[i][REAL] = 0.0;
complex[i][IMAG] = 0.0;
}
// Set first 4 vertical modes
complex[1 * logical_width][REAL] = 32;
complex[2 * logical_width][REAL] = 16;
complex[3 * logical_width][REAL] = 8;
complex[4 * logical_width][REAL] = 4;
// Print before IDFT -> DFT
printComplex(complex,N);
// IDFT
fftw_execute(iplan);
// DFT back
fftw_execute(plan);
// Print after IDFT -> DFT
printComplex(complex,true); // Pass true to divide amplitudes by N*N
// Clean up
fftw_destroy_plan(plan);
fftw_destroy_plan(iplan);
delete[] real;
fftw_free(complex);
两个printComplex(...)
调用的输出可以在上面的问题中看到。
解决方法
您需要阅读离散傅立叶变换。
对于实值时域信号,DFT具有共轭对称性:
F(k)= conj(F(N-k)),
样本数为N。通过对非对称频域信号进行逆变换,可以得到复数值的时域信号,但是由于使用的是复数到实数的转换,因此实际上只计算了该结果的实数部分。您在这里丢掉了一半的数据。然后,前向变换返回该变换信号的DFT。因为您的时域信号现在是实值,所以您的频域结果具有共轭对称性。