4

我正在使用 numpy 和 matplotlib 从我的模拟中分析数据输出。有一个（明显的）不一致之处，我找不到其根源。如下：

我有一个具有给定能量 a^2~1 的信号。当我使用 rfft 进行 FFT 并计算傅里叶空间中的能量时，它会明显更大。为了避免提供我的数据等的详细信息，这里有一个简单的正弦波示例：

from pylab import *
xx=np.linspace(0.,2*pi,128)
a=np.zeros(128)
for i in range(0,128):
    a[i]=sin(xx[i])
aft=rfft(a)
print mean(abs(aft)**2),mean(a**2)

原则上，这两个数字应该是相同的（至少在数字意义上），但这就是我从这段代码中得到的：

62.523081632 0.49609375

我试图通过 numpy.fft 文档但找不到任何东西。这里的搜索给出了以下内容，但我无法理解那里的解释：

现有（合成）信号与滤波信号之间的 FFT 幅度差异很大

我错过了什么/误解了什么？在这方面的任何帮助/指针将不胜感激。

谢谢！

4

2 回答 2

9

Henry 在非规范化部分是正确的，但还有更多内容，因为您使用的是rfft，而不是fft。以下与他的回答一致：

>>> x = np.linspace(0, 2 * np.pi, 128)
>>> y = 1 - np.sin(x)
>>> fft = np.fft.fft(y)
>>> np.mean((fft * fft.conj()).real)
191.49999999999991
>>> np.mean(y**2)
1.4960937500000004
>>> fft = fft / np.sqrt(len(fft))
>>> np.mean((fft * fft.conj()).real)
1.4960937499999991

但是，如果您现在尝试对 进行相同rfft操作，事情就不太顺利了：

>>> rfft = np.fft.rfft(y)
>>> np.mean((rfft * rfft.conj()).real)
314.58462009358772
>>> rfft /= np.sqrt(len(rfft))
>>> np.mean((rfft * rfft.conj()).real)
4.8397633860551954
65
>>> np.mean((rfft * rfft.conj()).real) / len(rfft)
4.8397633860551954

但是，以下确实可以正常工作：

>>> (rfft[0] * rfft[0].conj() +
...  2 * np.sum(rfft[1:] * rfft[1:].conj())).real / len(y)
1.4960937873636722

当你使用rfft你得到的东西时，你的数据的 DFT 并不正确，而只是它的正半部分，因为负值与它是对称的。要计算平均值，您需要将除 DC 分量之外的每个值都考虑两次，这是最后一行代码所做的。

于 2013-03-01T00:17:19.383 回答
3

在大多数 FFT 库中，各种 DFT 风格不是正交的。该numpy.fft库仅在逆变换期间应用必要的归一化。

考虑DFT 的 Wikipedia 描述；逆 DFT 具有 DFT 没有的 1/N 项（其中 N 是变换的长度）。要制作 DFT 的正交版本，您需要将未归一化的 DFT 的结果缩放 1/sqrt(N)。在这种情况下，变换是正交的（也就是说，如果我们将正交 DFT 定义为 F，那么逆 DFT 就是 F 的共轭或厄米转置）。

aft在您的情况下，您可以通过简单地缩放来获得正确的答案1.0/sqrt(len(a))（注意 N 是从变换的长度中找到的；真正的 FFT 只会抛出大约一半的值，所以长度a很重要）。

我怀疑将归一化保留到最后的原因是，在大多数情况下，这并不重要，因此您可以节省进行两次归一化的计算成本。事实上，非常快速的FFTW 库不会在任何一个方向上进行任何规范化，而是完全由用户来处理。

编辑：为了清楚起见，上面的解释并不完全正确。使用这种简单的缩放不会得出正确的答案，因为在您的情况下，DC 分量将被添加两次，尽管1.0/sqrt(len(a))仍然是产生单一变换的正确缩放。

于 2013-02-28T23:42:34.317 回答