Kissfft 缩放答案

【问题标题】：kissfft scalingKissfft 缩放
【发布时间】：2011-08-03 10:31:33
【问题描述】：

我希望使用 FFT 和 Kissfft 库计算快速相关性，并且缩放需要精确。什么缩放是必要的（向前和向后）以及我使用什么值来缩放我的数据？

【问题讨论】：

是的，我应该更具体一些。我希望有人能帮我省去自己发现缩放的麻烦（特定于 Kissfft），但是我自己做了，可以告诉你我必须将逆缩放 1/N。当然，最初出于效率目的而忽略了缩放，但在计算相关性时它很重要。自相关并不重要，因为结果是成比例的，但重要的是当你计算像 ASDF 这样的东西时，它有两个精确的项要添加到末尾。
更详细的信息，在 Kiss_fftr 上进行测试（实际 FFT 的实现）：正向 fttr 按nfft/2 缩放幅度，而反向 fftr 通过2 缩放幅度，与 nfft 无关。要获得原始幅度，您必须分别按2/nfft 和1/2 进行缩放。（这解释了为什么在串联进行正向和反向 FFT 时必须按 1/nfft 进行缩放。）

标签： signal-processing fft correlation kissfft

【解决方案1】：

将每个频率响应乘以 1/sqrt(N)，以获得 1/N 的整体缩放比例

在伪代码中：

ifft( fft(x)*conj( fft(y) )/N ) == circular_correlation(x,y)

至少对于浮点类型的 kisfft 是这样。

以下 c++ 示例代码的输出应该类似于

[1, 3i, 0 0 ....] 与自身的循环相关 = (10,0),(1.19796e-10,3),(-4.91499e-08,1.11519e- 15),(1.77301e-08,-1.19588e-08) ...

#include <complex>
#include <iostream>
#include "kiss_fft.h"
using namespace std;

int main()
{
    const int nfft=256;
    kiss_fft_cfg fwd = kiss_fft_alloc(nfft,0,NULL,NULL);
    kiss_fft_cfg inv = kiss_fft_alloc(nfft,1,NULL,NULL);

    std::complex<float> x[nfft];
    std::complex<float> fx[nfft];
    memset(x,0,sizeof(x));
    x[0] = 1;
    x[1] = std::complex<float>(0,3);

    kiss_fft(fwd,(kiss_fft_cpx*)x,(kiss_fft_cpx*)fx);
    for (int k=0;k<nfft;++k) {
        fx[k] = fx[k] * conj(fx[k]);
        fx[k] *= 1./nfft;
    }
    kiss_fft(inv,(kiss_fft_cpx*)fx,(kiss_fft_cpx*)x);
    cout << "the circular correlation of [1, 3i, 0 0 ....] with itself = ";
    cout
        << x[0] << ","
        << x[1] << ","
        << x[2] << ","
        << x[3] << " ... " << endl;
    kiss_fft_free(fwd);
    kiss_fft_free(inv);
    return 0;
}

【讨论】：

我使用的代码是：ifft( fft(x)*conj( fft(y) ) )/N == corr With Kissfft v1_2_9
是的。很抱歉最初误导您。我已经更新了描述以匹配我的代码实际在做什么。

【解决方案2】：

3 个最常见的 FFT 缩放因子是：

1.0 正向 FFT，1.0/N 反向 FFT
1.0/N 正向 FFT，1.0 反向 FFT
1.0/sqrt(N) 双向，FFT 和 IFFT

考虑到文档中可能存在的任何歧义，以及用户期望对其目的“正确”的任何缩放，最好只提供一个已知（1.0 浮点或 255 整数）幅度的纯正弦波，并且在FFT 长度到有问题的 FFT（和/或 IFFT），并查看缩放是否与上述之一匹配，可能与上述之一不同 2X 或 sqrt(2)，或者所需的缩放完全不同。

例如在您的环境中为您的数据类型编写一个 Kissfft 单元测试。

【讨论】：

这就是我最终做的事情，但我希望有人知道可以帮我省去麻烦。总的来说，这是一个很好的建议。
如果您不知道有这些不同的表示形式，您不会相信 FFT 文献是多么令人困惑 :)