FFTW 结果与 MATLAB 中的 FFT 不同答案

【问题标题】：FFTW results differ from FFT in MATLABFFTW 结果与 MATLAB 中的 FFT 不同
【发布时间】：2019-06-24 14:40:45
【问题描述】：

我使用 FFTW 和 MATLAB fft 比较前向 FFT。输入信号是高斯信号。代码：

FFTW 使用 C:

float *signal; /*input signal*/
int nt; /*length of the signal*/
fftw_complex *in, *out;
fftw_plan plan1;

in = fftw_malloc(nt*sizeof(fftw_complex)); 
out = fftw_malloc(nt*sizeof(fftw_complex));
for (j=0;j<nt;j++){
        in[j][0]=(double)signal[j];
        in[j][1]=0.0;
}    
plan1 = fftw_plan_dft_1d(nt, in, out, -1, FFTW_ESTIMATE);
fftw_execute(plan1);        
fftw_destroy_plan(plan1);

for (j=0;j<nt;j++){
        real[j]=(float)out[j][0];
        imag[j]=(float)out[j][1];
}

fft MATLAB 中的函数：

fft(signal);

我绘制了两个结果的实部和虚部：

实部的值几乎相同，而虚部的值完全不同。如何解决这个问题？

【问题讨论】：

Imag 的差异在 10 到 -15 次方的范围内。所以并没有那么不同。
另一个事实是，因为我的信号是实数，但在 FFTW 中我创建了一个复数 DFT，这意味着该信号是具有零虚值的复数。对于真实信号，out[i] 是 out[n-i] 的共轭，如 matlab fft 结果所示。但是 FFTW 的图像不会这样。
一个很好的测试是在这两种情况下应用逆变换并比较结果，这将显示虚部的差异是否相关（我也认为它们不相关）。跨度>
“这意味着信号是具有零虚值的复数” DFT 的虚值为零，因为输入是对称的，而不是因为输入是实值。此外，知道虚值应该为零，您应该已经注意到它们非常接近于零。您只是在这里查看浮点舍入错误。舍入误差不需要是对称的。
您可以通过将所有小于 epsilon（例如 10^-14）的虚数值归零来解决“非精确共轭”部分，因为您知道它无论如何都是浮点噪声。

标签： c matlab fft fftw

【解决方案1】：

您应该在“Imag”图的左侧查看图的比例因子。它说 10^-15。这与实际信号幅度（至少大于 10^1 的较大部分）相比非常小，因此结果非常相似。

浮点算法通常不会提供完全相同的结果，只要它们没有以完全相同的方式实现。（即便如此，它们也会因舍入的不同选项而有所不同。

此 QA 可能会提供一些见解： Floating point inaccuracy examples

【讨论】：

在这种情况下，“相当小”只能相对于其他东西来判断。绝对量级是不够的。然而，在这种情况下，实信号的几个特征（最大值、平均值、绝对范围）中的任何一个都为比较提供了适当的基础，以得出结论，是的，在这种情况下虚信号可以忽略不计。
@JohnBollinger 完全正确，我只是根据频谱的真实部分推断出这一点。有机会让这一点更清楚。
一般来说，浮点算法往往不会提供完全相同的结果，只要它们没有以完全相同的方式实现。不仅如此，FFTW 还可以产生不同的结果结果取决于选择的确切解决方案（FFTW 的“计划”）方法。例如，如果将 40 元素数据数组拆分为 4x10 或 5x8 子元素，计算 DFT 的结果会略有不同。在不同的计算机上运行它或者只是碰巧得到不同的内存布局，你可以得到不同的结果。要从 FFTW 获得一致的*结果，需要保存计划并重复使用它们。
* - 与“完全相同的输入总是产生完全相同的输出直到最后一位”一致。
@AndrewHenle 很高兴知道我从未与 FFTW 合作过，所以我无法触及这个问题，但这通常属于我所说的关于不同算法实现的内容，因为 FFTW 似乎有一些选择from.This有人可能不会想到这种行为，所以谢谢你的补充。

【解决方案2】：

四舍五入到最接近的满量程（实数）的 0.001%，注意虚数全为零。

【讨论】：