时间戳

首页 > 解决方案 > dft（离散傅立叶变换）与 C++ 代码

c++ - dft（离散傅立叶变换）与 C++ 代码

问题描述

我在 Visual Studio C++ 中编写了 C。

在 Visual Studio C++ 中，我使用 C 代码语法而不是 C++ 操作复数<complex>。

a[]是 4001 数组，所以使用b[4001]存储操作值，最后返回a[]。 NXm从 main 定义为 4001。

当我与 matlab 的 fft 的结果进行比较时，差异发生在第 169 到第 4000 个值之间。

你看看有没有错误？或者是什么原因？

感谢您阅读问题。

  1. 我尝试从 NXm 到 0 反向推进“for”。
  2. 我试过double ak = (double)k * (double)n * (2.0 * M_PI / (double)NXm); 改成double ak = k * n * (2.0 * M_PI / (double)NXm);
  3. 当我尝试操作短字段时，功能运行良好。
    这是代码。
void fft(complex<double> a[], int NXm)
{

    complex<double> sum = 0.0 + 0.0*I; ;
    complex<double> c = 1.0*I;
    complex<double> b[4001] = { 0 };
    for (int k = 0; k < NXm; k++)
    {   
        sum = 0.0 + 0.0*I;
        for (int n = 0; n < NXm; n++)
        {
            double ak = (double)k * (double)n * (2.0 * M_PI / (double)NXm);
            sum = sum + a[n] * exp(-c * ak);
        }
        b[k] = sum;

    }
    for (int i = 0; i < NXm; i++)
    {
        a[i] = b[i];
    }
}

我希望得到与 matlab 的 fft 几乎相同的结果。epsilon 级别的轻微错误是可以的

标签： c++fftdft

解决方案

复指数是：

exp(x + I*y) = exp(x) ( cos(y) + I*sin(y)) = exp(x) * cis(y)

问题可能是由处理器的 FPU 部分的缺陷引起的，如果y相当小，则内置一个内在函数以sin(y)获得比|y - sin(y)|. 欧拉公式的幼稚实现可能会因小参数而出错，草率的可能会假设sin(y)为零。

最好sin()用余割内在函数代替，它是倒数正弦 (1/sin(x))。

另一个不太常见的问题是当三角函数由大于 的值提供时 2*PI*n，其中n自然大于 1。精度随着 和 的增加而n下降sin (x) != sin(2*PI*n + x)

第三，根据编译器的不同，可能会有优化标志来规范编译器如何处理数学函数和浮点数学的使用。他们可能默认进行不精确但快速的计算。这可能会在 C++ 代码中使用 C 标头取代。

推荐阅读