赞
踩
- import numpy as np#导入一个数据处理模块
-
- import matplotlib.pyplot as plt#导入一个绘图模块
-
- # 依据快速傅里叶算法得到信号的频域
- def test_fft():
- sampling_rate = 8192 # 采样率
- fft_size = 8192 # FFT取样长度
- t = np.arange(0, 8.192, 1.0 / sampling_rate)
- #np.arange(起点,终点,间隔)产生8.192s长的取样时间
- x=0.6*np.sin(2*np.pi*500*t)+0.6*np.sin(2*np.pi*50*t)
- # 两个正弦波叠加,500HZ和50HZ
- # N点FFT进行精确频谱分析的要求是N个取样点包含整数个取样对象的波形。
- # 因此N点FFT能够完美计算频谱对取样对象的要求是n*Fs/N(n*采样频率/FFT长度),
- # 因此对8KHZ和512点而言,完美采样对象的周期最小要求是8000/512=15.625HZ,
- # 所以156.25的n为10,234.375的n为15。
-
- xs = x[:fft_size]# 从波形数据中取样fft_size个点进行运算
-
- xf = np.fft.rfft(xs) / fft_size # 返回fft_size/2+1 个频率
- #利用np.fft.rfft()进行FFT计算,rfft()是为了更方便对实数信号进行变换,
- # 由公式可知 / fft_size为了正确显示波形能量
- # rfft函数的返回值是N/2+1个复数,分别表示从0(Hz)到sampling_rate/2(Hz)的分。
-
- # 于是可以通过下面的np.linspace计算出返回值中每个下标对应的真正的频率:
- freqs = np.linspace(0, sampling_rate*10, fft_size/2+1 ) # 表示频率
- #freqs = np.linspace(0, sampling_rate/2 , fft_size/2 + 1) # 表示频率
-
- xfp = 20 * np.log10(np.clip(np.abs(xf), 1e-20, 1e100))
- #xfp = np.abs(xf) * 2 # 代表信号的幅值,即振幅
- # 最后我们计算每个频率分量的幅值,并通过 20*np.log10()将其转换为以db单位的值。
- # 为了防止0幅值的成分造成log10无法计算,我们调用np.clip对xf的幅值进行上下限处理
-
- plt.figure(figsize=(8, 4))
- plt.subplot(211)
- plt.plot(t[:fft_size], xs)
- plt.xlabel(u"时间(秒)", fontproperties='FangSong')
- plt.title(u"500Hz和50Hz的波形和频谱", fontproperties='FangSong')
-
- plt.subplot(212)
- plt.plot(freqs, xfp)
- plt.xlabel(u"频率(Hz)", fontproperties='FangSong')
- #字体FangSong
- plt.ylabel(u'幅值', fontproperties='FangSong')
- plt.subplots_adjust(hspace=0.4)
- '''subplots_adjust(left=None, bottom=None, right=None, top=None, wspace=None, hspace=None)
- 有六个可选参数来控制子图布局。值均为0~1之间。其中left、bottom、right、top围成的区域就是子图的区域。
- wspace、hspace分别表示子图之间左右、上下的间距。实际的默认值由matplotlibrc文件控制的。
- '''
- plt.show()
-
-
-
- test_fft()
Copyright © 2003-2013 www.wpsshop.cn 版权所有,并保留所有权利。