- # -*- coding: utf-8 -*-
- import numpy as np
- import matplotlib.pyplot as plt
- from scipy import signal
- #iirfilter函数参数说明,详见[1]
- #第一个参数是滤波器的阶数
- #第二个参数是滤波器的通带频率
- #rp是通带纹波最大波动1dB
- #rs是阻带最小衰减
- #btype是滤波器的滤波方式,除了band以外,还有lowpass、highpass等
- #analog是选择该滤波器为模拟滤波器(True)还是数字滤波器(False)
- #ftype是滤波器的种类,除了ellip,还有bessel、cheby2等
- #fs为数字信号的采样率
- #output是滤波器函数的返回值类型,这也是官方推荐的返回值类型
- sos = signal.iirfilter(10,[20,100], rp = 1, rs=60, btype='bandpass',
- analog=False, ftype='ellip',fs=1000,output='sos')
- #通过iirfilter得到了滤波器,我们接着使用sosfreqz获得该滤波器的频率响应
上述滤波器的幅频特性如下图所示
得到了这个椭圆滤波器,我们接下来看看,这个滤波器的过滤效果,看下面代码
- # -*- coding: utf-8 -*-
- import matplotlib.pyplot as plt
- from scipy import signal
- import numpy as np
- fs1 = 1000 #采样率1000Hz
- t = np.arange(fs1)/1000
- y1 = np.sin(2*np.pi*50*t)
- y2 = np.sin(2*np.pi*110*t)
- s = y1 + y2
- sos = signal.iirfilter(10,[20,100], rp = 1, rs=60, btype='bandpass',
- analog=False, ftype='ellip',fs=1000,output='sos')
- #使用滤波器sos过滤信号s,详见[3]
- y = signal.sosfilt(sos, s)
- plt.subplot(3, 1, 1)
- plt.plot(t,s)
- plt.ylabel('50Hz+110Hz')
- plt.subplot(3, 1, 2)
- plt.plot(t,y1)
- plt.ylabel('50Hz')
- plt.subplot(3, 1, 3)
- plt.plot(t,y)
- plt.ylabel('filtered waveform')
- plt.show()
过滤结果对比如下如所示,可以看到,在第三个图中很好的过滤掉了混合信号s中110Hz的信号
最后,我们来对比一下,采用上述生成的滤波器,仅仅更换不同的滤波器种类,对比一下几种不同滤波器的幅频特性。这里,我直接给出对比图。
可以看出,不论从通带、阻带、过渡带来看,椭圆滤波器的滤波器性能是最优的,其次是切比雪夫II型滤波器,然后是巴特沃斯滤波器,最后才是贝塞尔滤波器。当然,要注意的是,虽然椭圆滤波器的滤波性能优良,但是就一个滤波器而言,可能还需考虑计算成本,以及相频特性,结合三者才能合理的选择一个合适的滤波器。尽管如此,但对于研究算法而言,先在计算机上仿真,用最佳性能的滤波器处理信号以让我们看到或者接近算法的边界依然是很有参考价值的。
对于更多更多的滤波器设计,请参考scipy的signal模块。
[1]https://docs.scipy.org/doc/scipy/reference/reference/generated/scipy.signal.iirfilter.html#scipy.signal.iirfilter
[2]https://docs.scipy.org/doc/scipy/reference/reference/generated/scipy.signal.sosfreqz.html#scipy.signal.sosfreqz
[3]https://docs.scipy.org/doc/scipy/reference/reference/generated/scipy.signal.sosfilt.html#scipy.signal.sosfilt