2016-05-14 33 views
3

我仍然试图在Python中使用FFT对此data进行频率分析。 采样率是每分钟1个数据点。使用Python的FFT - 意外的低频

我的代码是:

from scipy.fftpack import fft 
df3 = pd.read_csv('Pressure - Dates by Minute.csv', sep=",", skiprows=0) 
df3['Pressure FFT'] = df3['ATMOSPHERIC PRESSURE (hPa) mean'] - df3['ATMOSPHERIC PRESSURE (hPa) mean'].mean() 
Pressure = df3['Pressure FFT'] 
Fs = 1/60 
Ts = 1.0/Fs 
n = len(Pressure) 
k = np.arange(n) 
T = n/Fs 
t = np.arange(0,1,1/n) # time vector 
frq = k/T # two sides frequency range 
frq = frq[range(int(n/2))] # one side frequency range 

Y = np.fft.fft(Pressure)/n # fft computing and normalization 
Y = Y[range(int(n/2))] 

fig, ax = plt.subplots(2, 1) 
ax[0].plot(t,Pressure) 
ax[0].set_xlabel('Time') 
ax[0].set_ylabel('Amplitude') 
ax[1].plot(frq,abs(Y),'r') # plotting the spectrum 
ax[1].set_xlabel('Freq (Hz)') 
ax[1].set_ylabel('|Y(freq)|') 

但结果得出:

enter image description here

所以我的问题是:

1)为什么没有频率呢?数据显然是周期性的。

2)为什么频谱如此之低? (0 - 0.009)

3)也许我应该尝试不同的过滤技术?

任何见解?

谢谢!

+1

由FFT函数返回的数组中的第一项具有DC分量,即原始数组中的所有值的总和。这通常比周期性分量大几个数量级。在绘图之前尝试绘制'Y [1:]'或绘制'Y [0] = 0',您应该看到您的频率出现。 – Jaime

+0

我试图在绘制之前先做Y [0] = 0,但仍然没有频率。也许这是正常化?因为周期性行为是一天两次。 – ValientProcess

+1

通常,您要使用对数刻度绘制Y轴。您还忘了在FFT之前应用合适的窗口功能。 –

回答

2

1)为什么没有频率呢?数据显然是周期性的。

那么,有频率的内容,它只是不完全可见,因为它的结构。试着改变了标绘的频谱线,从ax[1].plot(frq,abs(Y),'r')ax[1].semilogy(frq,abs(Y),'r')

这将导致到:

semilog Spectrum

如果我们现在已经应用了简单的改造是提升低值和高限制值。欲了解更多信息,请参阅​​。当然,删除DC(就像你在代码的第3行那样)也有帮助。

这似乎仍然有点模糊,它是,但如果我们放大到频谱的下部,我们可以看到这一点:

semilog Transform

其中显示了一个峰值大约2.3E-05 Hz,相当于大约12小时。

2)为什么频谱如此之低? (0-0。009)

因为您每60秒采样一次,所以您的采样频率为(大约)0.016 Hz。您的频谱包含DC(0Hz)和0.0083Hz之间的所有内容。欲了解更多信息,请参阅this link

3)也许我应该尝试不同的过滤技术?

如果你不能解决谐波但你看起来并不需要它,你可以尝试开窗。

希望这会有所帮助。

1

这些频率看起来如此之低的部分原因是因为您的幅度图中的时间轴会被奇怪地缩放。如果您每60秒确实有一个样本,那么x轴的范围应介于0到1690260秒之间(即〜20天!)。

enter image description here

通过眼睛,你似乎有每50000秒(〜2每天)约一小峰,这将对应于大约2x10⁻⁵Hz的频率。因此,您的周期图对我而言看起来相当合理,因为x轴的尺寸有多大。