使用 numpy 进行 FFT 归一化答案

【问题标题】：FFT normalization with numpy使用 numpy 进行 FFT 归一化
【发布时间】：2016-09-01 11:02:30
【问题描述】：

刚开始使用 numpy 包，并从计算输入信号的 FFT 的简单任务开始。代码如下：

import numpy as np
import matplotlib.pyplot as plt

#Some constants
L = 128
p = 2
X = 20
x = np.arange(-X/2,X/2,X/L)
fft_x = np.linspace(0,128,128, True)

fwhl = 1

fwhl_y = (2/fwhl) \
*(np.log([2])/np.pi)**0.5*np.e**(-(4*np.log([2]) \
*x**2)/fwhl**2)

fft_fwhl = np.fft.fft(fwhl_y, norm='ortho')

ampl_fft_fwhl = np.abs(fft_fwhl)

plt.bar(fft_x, ampl_fft_fwhl, width=.7, color='b')

plt.show()

由于我使用的是一个指数函数，其前面有一些常数除以 pi，因此我希望在傅里叶空间中得到指数函数，其中 FFT 的常数部分始终等于 1（零频率）。但是我使用 numpy 获得的那个组件的值更大（大约是 1,13）。在这里，我有一个幅度谱，它由 1/(number_of_counts)**0.5 归一化（这就是我在 numpy 文档中读到的）。我不明白出了什么问题...有人可以帮帮我吗？

谢谢！

[EDITED] 看来问题已经解决了，你需要得到与傅里叶积分和 FFT 相同的结果就是将 FFT 乘以步长（在我的例子中是 X/L）。至于作为 numpy.fft.fft(..., norm='ortho') 选项的归一化，它仅用于保存变换的比例，否则您需要将逆 FFT 的结果除以样本数。感谢大家的帮助！

【问题讨论】：

小心除法 - 如果您在 python 中将两个整数相除，它们不会返回浮点数。如果你将 15 除以 10，你将得到 1，而不是 1.5。同样，20/128= 0，而不是 0.15
@flyingmeatball 在 python 2.x 上确实如此。在 python 3.x 上，它们被强制转换为明确地浮动；）
零频率对应于输入的平均值：fft_fwhl[0] # (0.56568542494923801+0j) 匹配 np.mean(fwhl_y) * np.sqrt(len(fwhl_y)) # 0.56568542494923812（它是 sqrt，因为您使用的是 norm="ortho"）。
我发现从 DFT 定义的角度考虑 FFT 比连续傅里叶变换更有帮助：docs.scipy.org/doc/numpy/reference/… 通常很难让 DFT/FFT 给出相同的答案作为连续傅里叶变换。
对于希望更全面地介绍 FFT 缩放（和窗口化！）的每个人，我可以推荐 Heinzel 等人的这篇论文。这对我理解 FFT 缩放有很大帮助。它比大多数介绍性文本更面向应用程序，并且仍然有足够的理论来超越简单的配方：Spectrum and spectral density estimation by the Discrete Fouriertransform (DFT), including a comprehensive list of window functions and some new flat-top windows, G. Heinzel et al., 2002

标签： python numpy fft

【解决方案1】：

我终于解决了我的问题。将 FFT 与傅立叶积分结合所需的只是将变换 (FFT) 的结果乘以步长（在我的情况下为 X/L，FFTX/L），它通常可以工作。就我而言，它有点复杂，因为我有一个额外的规则来转换函数。我要确保曲线下面积等于1，因为它是δ函数的模型，所以由于步长是不可改变的，所以我要满足步长sum(fwhl_y)=1的条件，即X/L=1/总和（fwhl_y）。所以为了得到正确的结果，我必须做以下事情：

计算 FFT fft_fwhl = np.fft.fft(fwhl_y)
去除fwhl_y函数对称性带来的相位分量，即定义在[-T/2,T/2]区间的函数，其中T是周期，np .fft.fft 操作认为我的函数定义在 [0,T] 区间。因此，仅获取幅度谱（这就是我需要的）我只需使用 np.abs(FFT)
要获得我期望的值，我应该将上一步得到的结果乘以 X/L，即 np.abs(FFT)*X/L
我对曲线下的面积有一个额外的条件，所以它是 X/L*sum(fwhl_y)=1，我终于来到 np.abs(FFT)*X /L = np.abs(FFT)/sum(fwhl_y)

希望它至少可以帮助任何人。

【讨论】：

【解决方案2】：

以下是解决您问题的可能方法：

import numpy as np
import matplotlib.pyplot as plt
from scipy import fft
from numpy import log, pi, e

# Signal setup
Fs = 150
Ts = 1.0 / Fs
t = np.arange(0, 1, Ts)
ff = 50
fwhl = 1
y = (2 / fwhl) * (log([2]) / pi)**0.5 * e**(-(4 * log([2]) * t**2) / fwhl**2)

# Plot original signal
plt.subplot(2, 1, 1)
plt.plot(t, y, 'k-')
plt.xlabel('time')
plt.ylabel('amplitude')

# Normalized FFT
plt.subplot(2, 1, 2)
n = len(y)
k = np.arange(n)
T = n / Fs
frq = k / T
freq = frq[range(n / 2)]

Y = np.fft.fft(y) / n
Y = Y[range(n / 2)]

plt.plot(freq, abs(Y), 'r-')
plt.xlabel('freq (Hz)')
plt.ylabel('|Y(freq)|')

plt.show()

使用 fwhl=1：

使用 fwhl=0.1：

您可以在上图中看到当 fwhl 接近 0 时指数和 FFT 图如何变化

【讨论】：

但是 FFT 的零分量仍然不同，尽管它应该始终保持不变。我的意思是，如果我求解积分以找到this function 的 FFT，我得到this result。所以无论参数的值是多少，我总是在零点有 1。在你的情况下，FFT 的零分量总是会有不同的值。
@GuntherOnFire 哦，是的，DFT/FFT 是在周期性的离散序列上定义的。一直以来的傅立叶积分与离散化函数的 FFT 不易相关。
@Ahmed Fasih，似乎我终于找到了如何将傅里叶积分与 FFT/DFT 连接起来的解决方案，您只需将 FFT 的结果乘以您选择的步长。就我而言，它是 X/L。
@GuntherOnFire 不错！你能把它写成答案并接受吗？
@Ahmed Fasih，终于成功了。希望对你有用。