在 Python 中绘制快速傅里叶变换

Question

我可以访问 NumPy 和 SciPy，并希望创建一个数据集的简单 FFT。我有两个列表，一个是 y 值，另一个是那些 y 值的时间戳。

将这些列表输入 SciPy 或 NumPy 方法并绘制结果 FFT 的最简单方法是什么？

我已经查找了示例，但它们都依赖于创建一组具有一定数量的数据点和频率等的假数据，并没有真正展示如何仅使用一组数据和相应的时间戳来做到这一点.

我尝试了以下示例：

 from scipy.fftpack import fft

# Number of samplepoints
N = 600

# Sample spacing
T = 1.0 / 800.0
x = np.linspace(0.0, N*T, N)
y = np.sin(50.0 * 2.0*np.pi*x) + 0.5*np.sin(80.0 * 2.0*np.pi*x)
yf = fft(y)
xf = np.linspace(0.0, 1.0/(2.0*T), N/2)
import matplotlib.pyplot as plt
plt.plot(xf, 2.0/N * np.abs(yf[0:N/2]))
plt.grid()
plt.show()

但是当我将 fft 的参数更改为我的数据集并绘制它时，我得到了非常奇怪的结果，并且频率的缩放似乎可能关闭了。我不确定。

这是我尝试进行 FFT 的数据的粘贴箱

http://pastebin.com/0WhjjMkb http://pastebin.com/ksM4FvZS

当我使用 fft() 时，它只是在零处有一个巨大的峰值，没有别的。

这是我的代码：

 ## Perform FFT with SciPy
signalFFT = fft(yInterp)

## Get power spectral density
signalPSD = np.abs(signalFFT) ** 2

## Get frequencies corresponding to signal PSD
fftFreq = fftfreq(len(signalPSD), spacing)

## Get positive half of frequencies
i = fftfreq>0

##
plt.figurefigsize = (8, 4)
plt.plot(fftFreq[i], 10*np.log10(signalPSD[i]));
#plt.xlim(0, 100);
plt.xlabel('Frequency [Hz]');
plt.ylabel('PSD [dB]')

间距正好等于 xInterp[1]-xInterp[0] 。

原文由 user3123955 发布，翻译遵循 CC BY-SA 4.0 许可协议

Answer 1

因此，我在 IPython 笔记本中运行您的代码的功能等效形式：

 %matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
import scipy.fftpack

# Number of samplepoints
N = 600
# sample spacing
T = 1.0 / 800.0
x = np.linspace(0.0, N*T, N)
y = np.sin(50.0 * 2.0*np.pi*x) + 0.5*np.sin(80.0 * 2.0*np.pi*x)
yf = scipy.fftpack.fft(y)
xf = np.linspace(0.0, 1.0/(2.0*T), N//2)

fig, ax = plt.subplots()
ax.plot(xf, 2.0/N * np.abs(yf[:N//2]))
plt.show()

我得到了我认为非常合理的输出。

自从我在工程学校考虑信号处理以来，这比我愿意承认的要长，但 50 和 80 的峰值正是我所期望的。那么问题是什么？

响应发布的原始数据和评论

这里的问题是您没有定期数据。您应该始终检查您输入 任何 算法的数据，以确保它是适当的。

 import pandas
import matplotlib.pyplot as plt
#import seaborn
%matplotlib inline

# the OP's data
x = pandas.read_csv('http://pastebin.com/raw.php?i=ksM4FvZS', skiprows=2, header=None).values
y = pandas.read_csv('http://pastebin.com/raw.php?i=0WhjjMkb', skiprows=2, header=None).values
fig, ax = plt.subplots()
ax.plot(x, y)

原文由 Paul H 发布，翻译遵循 CC BY-SA 4.0 许可协议

Answer 2

fft 的重要之处在于它只能应用于时间戳统一的数据（ 即 时间上的统一采样，如上所示）。

在非均匀采样的情况下，请使用函数来拟合数据。有几个教程和功能可供选择：

https://github.com/tiagopereira/python_tips/wiki/Scipy%3A-curve-fitting http://docs.scipy.org/doc/numpy/reference/generated/numpy.polyfit.html

如果拟合不是一种选择，您可以直接使用某种形式的插值将数据插值到均匀采样：

https://docs.scipy.org/doc/scipy-0.14.0/reference/tutorial/interpolate.html

当你有统一的样本时，你只需要担心样本的时间增量（ t[1] - t[0] ）。在这种情况下，您可以直接使用 fft 函数

Y    = numpy.fft.fft(y)
freq = numpy.fft.fftfreq(len(y), t[1] - t[0])

pylab.figure()
pylab.plot( freq, numpy.abs(Y) )
pylab.figure()
pylab.plot(freq, numpy.angle(Y) )
pylab.show()

这应该可以解决您的问题。

原文由 ssm 发布，翻译遵循 CC BY-SA 3.0 许可协议