Python scipy利用快速傅里叶变换实现滤波
作者:微小冷
这篇文章主要为大家详细介绍了Python scipy如何利用快速傅里叶变换实现滤波,文中的示例代码讲解详细,感兴趣的小伙伴可以跟随小编一起学习一下
fft模块简介
scipy官网宣称,fftpack模块将不再更新,或许不久之后将被废弃,也就是说fft将是唯一的傅里叶变换模块。
Fourier变换极其逆变换在数学上的定义如下
下表整理出一部分与Fourier变换相关的函数,其中FFT为快速Fourier变换(Fast Fourier Transform);DCT为离散余弦变换(Discrete Cosine Transform);DST为离散正弦变换(discrete sine transform),另外,函数的前缀和后缀有如下含义
- i表示逆变换;
- 2, n分别表示2维和n维
正变换 | 逆变换 | |
---|---|---|
通用 | fft, fft2, fftn | ifft, ifft2, ifftn |
实数域 | rfft, rfft2, rfftn | irfft, irfft2, irfftn |
厄米对称 | hfft, hfft2, hfftn | ihfft, ihfft2, ihfftn |
离散余弦变换 | dct, dctn | idct, idctn |
离散正弦变换 | dst, dstn | idst, idstn |
汉克尔变换 | fht | ifht |
移动零频 | fftshift | ifftshift |
DFT采样频率 | fftfreq | ifftfreq |
fft函数示例
在数值计算中,一切输入输出均为离散值,所以实际上用到的是离散Fourier变换,即DFT,其功能是将离散的采样变换为一个离散的频谱分布。
下面将手动创建一组采样点,并添加一点噪声,然后通过FFT获取其频域信息。
import numpy as np from scipy import fft PI = np.pi*2 fs = 60 #采样频率 T = 100 #采样周期数 N = fs*T #采样点 t = np.linspace(0, T, N) noise = 2 * np.random.randn(*t.shape) s = 2 * np.sin(PI * t) + 3 * np.sin(22 * PI * t) + noise F = fft.fft(s) f = fft.fftfreq(N, 1.0/fs)
其中,t
为时间序列,s
为模拟的采样点,F
是Fourier变换之后的结果。但由于fft
默认是在复数域上的,故而可以查看其实部、虚部、模和辐角的值。
下面对采样点以及Fourier变换的结果进行绘制
import matplotlib.pyplot as plt fig = plt.figure() ax = fig.add_subplot(2,2,1) ax.plot(t, s) ax.set_title("t vs s") f_abs = np.abs(F) ax = fig.add_subplot(2,2,2) ax.plot(f, f_abs) ax.set_title("fs vs |F|") xlims = [[0,2], [21,23]] for i, xlim in enumerate(xlims): ax = fig.add_subplot(2,2,3+i) ax.stem(f, f_abs) ax.set_title("fs vs |F|") ax.set_xlim(xlim) plt.show()
结果为
即f=1和f=22处被筛选了出来。
滤波
有了这个,就可以在频域上对数据进行滤波,其思路是,对f_abs
中的值进行阈值分割,例如,只筛选出低频部分,然后看一下滤波效果
fig = plt.figure(1) f_filt = F * (np.abs(f) < 2) s_filt = fft.ifft(f_filt) ax = fig.add_subplot() ax.plot(t, s, lw=0.2) ax.plot(t, s_filt.real, lw=2) ax.set_title("threshold=2") ax.set_xlim([0,10]) plt.show()
效果如下
到此这篇关于Python scipy利用快速傅里叶变换实现滤波的文章就介绍到这了,更多相关Python scipy实现滤波内容请搜索脚本之家以前的文章或继续浏览下面的相关文章希望大家以后多多支持脚本之家!