Python 中的快速傅立叶变换
Fast Fourier Transform in Python
我是傅里叶理论的新手,我看过很好的教程,介绍如何将 fft 应用于信号并绘制它以查看它包含的频率。不知何故,他们所有人都创建了正弦混合作为他们的数据,我无法将其适应我的实际问题。
我有 242 小时的观测值和每天的周期,这意味着我的周期是 24。所以我希望在我的 fft 图上有一个峰值在 24 左右。
我的 data.csv 的样本在这里:
https://pastebin.com/1srKFpJQ
绘制的数据:
我的代码:
data = pd.read_csv('data.csv',index_col=0)
data.index = pd.to_datetime(data.index)
data = data['max_open_files'].astype(float).values
N = data.shape[0] #number of elements
t = np.linspace(0, N * 3600, N) #converting hours to seconds
s = data
fft = np.fft.fft(s)
T = t[1] - t[0]
f = np.linspace(0, 1 / T, N)
plt.ylabel("Amplitude")
plt.xlabel("Frequency [Hz]")
plt.bar(f[:N // 2], np.abs(fft)[:N // 2] * 1 / N, width=1.5) # 1 / N is a normalization factor
plt.show()
这输出了一个非常奇怪的结果,似乎我在每个频率上都得到了相同的值。
我想问题出在 N、t 和 T 的定义上,但我在网上找不到任何帮助我清楚地理解这一点的东西。请帮助:)
编辑 1:
使用 charles answer 提供的代码,我在 0 附近有一个尖峰,这看起来很奇怪。我已经使用 rfft
和 rfftfreq
来避免频率过高。
我读到这可能是因为该系列的直流分量,所以在减去平均值后我得到:
我无法解释这一点,尖峰似乎周期性地发生,但以 Hz 为单位的值不允许我获得 24 值(总频率)。任何人都知道如何解释这个?我错过了什么?
您看到的问题是因为条形太宽,而您只看到一个条形。您必须将条形的宽度更改为 0.00001 或更小才能看到它们。
不要使用条形图,而是使用 fftfreq = np.fft.fftfreq(len(s))
制作 x 轴,然后使用绘图函数,plt.plot(fftfreq, fft)
:
import matplotlib.pyplot as plt
import pandas as pd
import numpy as np
data = pd.read_csv('data.csv',index_col=0)
data.index = pd.to_datetime(data.index)
data = data['max_open_files'].astype(float).values
N = data.shape[0] #number of elements
t = np.linspace(0, N * 3600, N) #converting hours to seconds
s = data
fft = np.fft.fft(s)
fftfreq = np.fft.fftfreq(len(s))
T = t[1] - t[0]
f = np.linspace(0, 1 / T, N)
plt.ylabel("Amplitude")
plt.xlabel("Frequency [Hz]")
plt.plot(fftfreq,fft)
plt.show()
我是傅里叶理论的新手,我看过很好的教程,介绍如何将 fft 应用于信号并绘制它以查看它包含的频率。不知何故,他们所有人都创建了正弦混合作为他们的数据,我无法将其适应我的实际问题。
我有 242 小时的观测值和每天的周期,这意味着我的周期是 24。所以我希望在我的 fft 图上有一个峰值在 24 左右。
我的 data.csv 的样本在这里: https://pastebin.com/1srKFpJQ
绘制的数据:
我的代码:
data = pd.read_csv('data.csv',index_col=0)
data.index = pd.to_datetime(data.index)
data = data['max_open_files'].astype(float).values
N = data.shape[0] #number of elements
t = np.linspace(0, N * 3600, N) #converting hours to seconds
s = data
fft = np.fft.fft(s)
T = t[1] - t[0]
f = np.linspace(0, 1 / T, N)
plt.ylabel("Amplitude")
plt.xlabel("Frequency [Hz]")
plt.bar(f[:N // 2], np.abs(fft)[:N // 2] * 1 / N, width=1.5) # 1 / N is a normalization factor
plt.show()
这输出了一个非常奇怪的结果,似乎我在每个频率上都得到了相同的值。
我想问题出在 N、t 和 T 的定义上,但我在网上找不到任何帮助我清楚地理解这一点的东西。请帮助:)
编辑 1:
使用 charles answer 提供的代码,我在 0 附近有一个尖峰,这看起来很奇怪。我已经使用 rfft
和 rfftfreq
来避免频率过高。
我读到这可能是因为该系列的直流分量,所以在减去平均值后我得到:
我无法解释这一点,尖峰似乎周期性地发生,但以 Hz 为单位的值不允许我获得 24 值(总频率)。任何人都知道如何解释这个?我错过了什么?
您看到的问题是因为条形太宽,而您只看到一个条形。您必须将条形的宽度更改为 0.00001 或更小才能看到它们。
不要使用条形图,而是使用 fftfreq = np.fft.fftfreq(len(s))
制作 x 轴,然后使用绘图函数,plt.plot(fftfreq, fft)
:
import matplotlib.pyplot as plt
import pandas as pd
import numpy as np
data = pd.read_csv('data.csv',index_col=0)
data.index = pd.to_datetime(data.index)
data = data['max_open_files'].astype(float).values
N = data.shape[0] #number of elements
t = np.linspace(0, N * 3600, N) #converting hours to seconds
s = data
fft = np.fft.fft(s)
fftfreq = np.fft.fftfreq(len(s))
T = t[1] - t[0]
f = np.linspace(0, 1 / T, N)
plt.ylabel("Amplitude")
plt.xlabel("Frequency [Hz]")
plt.plot(fftfreq,fft)
plt.show()