【问题标题】:Evaluating Fourier transform of a Gaussian with FFTW3 and C++使用 FFTW3 和 C++ 评估高斯的傅里叶变换
【发布时间】:2013-11-27 19:11:16
【问题描述】:

我尝试在 C++ 中使用 FFTW3 计算高斯函数的傅立叶变换。这是我的代码的主要部分

main(int argc, char** argv)
{
   fftw_plan p;
   complex<double> *in,*out;
   long N=8;

   //allocation of in and the fftw plan called 
   in=(complex<double>*) calloc(N,sizeof(complex<double>));
   p=fftw_plan_dft_1d(N,(fftw_complex*)in,(fftw_complex*)in,FFTW_BACKWARD,FFTW_ESTIMATE);

   //initialize the array in with the values of a Gaussian function
   gaussianf(in,N);
   //Fourier transform in
   fftw_execute(p);  
   //write the result into a file
   writeft(in,N);
   fftw_destroy_plan(p);
}

由于数组已经用高斯的值初始化,我希望输出也是高斯的,但实际上只有包络具有高斯形状。正如我在下面的数据中所示,可以看到出现了一些负值。

#input values
#x       real part     imag part

-10     3.72008e-44     0
-7.5    3.72336e-25     0
-5      1.38879e-11     0
-2.5    0.00193045      0
0       1       0
2.5     0.00193045      0
5       1.38879e-11     0
7.5     3.72336e-25     0

#output
#x       real part     imag part
-10     1.00386 0
-7.5    -1.00273        0
-5      1       0
-2.5    -0.99727        0
0       0.996139        0
2.5     -0.99727        0
5       1       0
7.5     -1.00273        0

谁能告诉我我做错了什么?我将衷心感谢您的帮助。非常感谢。

【问题讨论】:

  • 您是否刚刚使用了如图所示的开箱即用命令?我使用的是同一个库,但在转换实数和偶数数据时得到非零虚部。

标签: c++ fftw


【解决方案1】:

就 C 编程或 FFTW 调用而言,您没有做错任何事情:这些都是正确的值。高斯曲线的 FFT 的实部确实在零附近振荡。如果您绘制 absolute 值,这可能看起来更像您的预期。

【讨论】:

  • 非常感谢您的回答。也许我误解了一些东西,但是如果通过通常的积分计算 FT,则可以得到另一个没有负值的高斯。为什么不一样?
  • 如果你乘以cos(n*theta) 和中间的theta=0(即高斯峰值所在的位置),积分总是正数。但这不是 DFT 假设的:它假设 theta=0 在信号的开头。所以现在,根据n,余弦波和高斯峰值之间会有不同的相位关系。
  • 谢谢,这是有道理的,现在很清楚了。现在的问题是如何使用 FFTW 在数值上计算“真实”FT。
  • 在组成输入信号时采用相同的约定:第一个样本是 t=0,第二个是 t=1;最后一个是 t=-1,倒数第二个是 t=-2。我相信 fftw 具有某种循环移位功能,可让您在此和 t=0-in-the-middle 约定之间来回转换信号(matlab 将其称为 FFTSHIFT)。即使没有,手动编程也是一个非常简单的例程。
【解决方案2】:

这些波动是意料之中的。在实践中,需要使用窗口函数来减少这些振荡。 http://en.wikipedia.org/wiki/Window_function

【讨论】:

    猜你喜欢
    • 2019-04-25
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    相关资源
    最近更新 更多