参考自:https://blog.csdn.net/qq_27825451/article/details/88553441
核心是 numpy.fft.fft 以及 numpy.fft.ifft
我对原文的代码略作了一点修改,如下
import numpy as np
import matplotlib.pyplot as plt
#采样点选择1400个,因为设置的信号频率分量最高为600赫兹,根据采样定理知采样频率要大于信号频率2倍,所以这里设置采样频率为1400赫兹(即一秒内有1400个采样点,一样意思的)
x=np.linspace(0,1,1400)
#设置需要采样的信号,频率分量有200,400和600
y=7*np.sin(2*np.pi*200*x) + 5*np.sin(2*np.pi*400*x)+3*np.sin(2*np.pi*600*x)
fft_y=np.fft.fft(y) #快速傅里叶变换
N=1400
x = np.arange(N) # 频率个数
half_x = x[range(int(N/2))] #取一半区间
abs_y=np.abs(fft_y) # 取复数的绝对值,即复数的模(双边频谱)
angle_y=np.angle(fft_y) #取复数的角度
normalization_y=abs_y/N #归一化处理(双边频谱)
normalization_half_y = normalization_y[range(int(N/2))] #由于对称性,只取一半区间(单边频谱)
plt.subplot(121)
plt.plot(x,y)
plt.title('original')
plt.subplot(122)
plt.plot(half_x,normalization_half_y,'blue')
plt.title('one side, normalized',fontsize=9,color='blue')
plt.show()
运行代码会看到原信号,和频谱图。
- numpy.fft.fft对信号进行操作,会返还一个复数数组,每个复数的振幅是频率,辐角是相位。
- 代码中只取频谱的一半,这个我没理解,原文也没交待,应该是快速傅里叶变换算法里来的。
- 经过实验发现,x 的区间长度范围只能是 1 ,否则会得到错误结果。