递归 IFFT 实现中的缩放不正确

Question

我正在用 C++ 开发一个简单的递归 FFT/IFFT 实现。代码如下：

#include <complex>
#include <numbers>
#include <ranges>
#include <valarray>
#include <vector>

#include "assert.hpp"

template<typename T = float>
std::valarray<std::complex<T>> fft(std::valarray<std::complex<T>> const& data)
{
    const auto N = data.size();
    ASSERT(!(N & (N - 1)), "FFT length has to be power of two");

    if (N == 1)
        return data;

    auto even = fft(data[std::slice(0, N / 2, 2)]);
    auto odd = fft(data[std::slice(1, N / 2, 2)]);

    auto result = std::valarray<std::complex<T>>(N);
    for (size_t k = 0; k < N / 2; ++k)
    {
        auto t = std::exp(std::complex<T>(0, -2 * std::numbers::pi * T(k) / N));
        result[k] = even[k] + odd[k] * t;
        result[k + N / 2] = even[k] - odd[k] * t;
    }

    return result;
}

template<typename T = float>
std::valarray<std::complex<T>> ifft(std::valarray<std::complex<T>> const& data)
{
    const auto N = data.size();
    ASSERT(!(N & (N - 1)), "FFT length has to be a power of two");

    if (N == 1)
        return data;

    auto even = ifft(data[std::slice(0, N / 2, 2)]);
    auto odd = ifft(data[std::slice(1, N / 2, 2)]);

    auto result = std::valarray<std::complex<T>>(N);
    for (size_t k = 0; k < N / 2; ++k)
    {
        auto t = std::exp(std::complex<T>(0, 2 * std::numbers::pi * T(k) / N));
        result[k] = (even[k] + odd[k] * t) / T(N);
        result[k + N / 2] = (even[k] - odd[k] * t) / T(N);
    }

    return result;
}

虽然正向 FFT 工作正常，但反向 FFT 会错误地缩放其输出（对于大小为 8 的 FFT，缩放比例为 1/8）。 我正在使用此代码来测试我的实现：

std::valarray<std::complex<float>> data = { 3, 5, 8, 9, 1, 4, 7, 2 };

std::cout << "Data: ";
for (auto const& e : data)
    std::cout << e << " ";
std::cout << "\n";

auto spectrum = fft(data);

std::cout << "Forward FFT: ";
for (auto const& e : spectrum)
    std::cout << e << " ";
std::cout << "\n";

std::cout << "Inverse FFT: ";
auto inverse = ifft(spectrum);
for (auto const& e : inverse)
    std::cout << e << " ";
std::cout << "\n";

产生以下输出：

Data: (3,0) (5,0) (8,0) (9,0) (1,0) (4,0) (7,0) (2,0) 
Forward FFT: (39,0) (-2.24264,-6.65685) (-11,2) (6.24264,-4.65685) (-1,0) (6.24264,4.65685) (-11,-2) (-2.24264,6.65685) 
Inverse FFT: (0.375,0) (0.625,6.58126e-09) (1,2.73196e-09) (1.125,-1.21692e-08) (0.125,0) (0.5,-1.20452e-08) (0.875,-2.73196e-09) (0.25,1.76331e-08)

根据numpy，前向FFT是正确的，而逆向FFT当然应该是原始数据，但事实并非如此。

是什么原因导致这种行为？我该如何解决它？

Answer 1

我发现我做错了什么：

1/N

的缩放系数不正确。使用

1/sqrt(N)

时效果非常好。

IFFT 函数必须更改为：

template<typename T = float>
std::valarray<std::complex<T>> ifft(std::valarray<std::complex<T>> const& data)
{
    const auto N = data.size();
    ASSERT(!(N & (N - 1)), "FFT length has to be a power of two");

    if (N == 1)
        return data;

    auto even = ifft(data[std::slice(0, N / 2, 2)]);
    auto odd = ifft(data[std::slice(1, N / 2, 2)]);

    auto result = std::valarray<std::complex<T>>(N);
    for (size_t k = 0; k < N / 2; ++k)
    {
        auto t = std::exp(std::complex<T>(0, 2 * std::numbers::pi * T(k) / N));
        result[k] = (even[k] + odd[k] * t) / T(sqrt(N));
        result[k + N / 2] = (even[k] - odd[k] * t) / T(sqrt(N));
    }

    return result;
}