自定义 DFT 实现返回奇数索引答案的逆序

问题描述 投票:0回答:1

这是我的 DFT 实现,使用 {0, 1, 2, 3} 进行测试。

#include <vector>
#include <complex>
#include <numbers>
#include <cmath>
#include <iostream>

std::vector<std::complex<double>> DFT(std::vector<std::complex<double>>&& P)
{
    int n = P.size();
    if (n == 1) {
        return P;
    }

    std::vector<std::complex<double>> Pe(n / 2);
    std::vector<std::complex<double>> Po(n / 2);
    for (int i = 0; i < n / 2; ++i) {
        Pe[i] = P[2 * i];
        Po[i] = P[2 * i + 1];
    }

    auto ye = DFT(std::move(Pe)), yo = DFT(std::move(Po));

    // in place algorithm, use input P to store output
    auto wi = std::complex<double>(1.0, 0.0);
    auto wn = std::complex<double>(std::cos(2.0 * std::numbers::pi_v<double> / (double)n),
                                   std::sin(2.0 * std::numbers::pi_v<double> / (double)n));
    for (int i = 0; i < n / 2; ++i) {
        P[i]         = ye[i] + wi * yo[i];
        P[i + n / 2] = ye[i] - wi * yo[i];
        wi           = wi * wn;
    }

    return P;
}

int main()
{
    int                               N = 4;
    std::vector<std::complex<double>> p_input(N);
    for (int i = 0; i < N; ++i)
    {
        p_input[i] = {(double)i, 0.0};
    }

    auto p_output = DFT(std::move(p_input));

    for (int i = 0; i < p_output.size(); ++i)
    {
        std::cout << p_output[i] << std::endl;
    }
}

测试结果是

(6,0)
(-2,-2)
(-2,0)
(-2,2)

但是 MATLAB 的答案是

>> fft(0:1:3)

ans =

   6.0000 + 0.0000i  -2.0000 + 2.0000i  -2.0000 + 0.0000i  -2.0000 - 2.0000i

我测试了更长的输入长度,并且我的结果总是在奇数索引位置具有相反的顺序。

不知道哪一部分是错误的?我的参考公式是:

P(x): [p_0, p_1, ..., p_{n-1}]

w: [w^0, w^1, ..., w^{n-1}]

Pe(x^2): [p_0, p_2, ..., p_{n-2}]

Po(x^2): [p_1, p_3, ..., p_{n-1}]

ye = [Pe(w^0), Pe(w^2), ..., Pe(w^{n-2})]

yo = [Po(w^0), Po(w^2), ..., Po(w^{n-2})]

P(w^j) = ye[j] + w^j yo[j]

P(w^{j+n/2}) = ye[j] - w^j yo[j]

y = [P(w^0), P(w^1), ..., P(w^{n-1})]

c++ fft
1个回答
0
投票

您只是使用了与其实现不同的 n 次单位根,即

cos(2*pi / n) + sin(2*pi / n) * i

而不是

cos(2*pi / n) - sin(2*pi / n) * i

就FFT而言,你使用哪一种并不重要,只要逆FFT与其一致即可。我对 FFT 不太熟悉,不知道按照惯例,其中一种方法是否优于另一种方法。

© www.soinside.com 2019 - 2024. All rights reserved.