滤波器的直接设计

4.3　IIR滤波器的直接设计

前面所给出的两种数字滤波器的设计事先都需要知道模拟滤波器的系统函数，下面来介绍一种直接进行滤波器设计的方法。常用的滤波器有Butterworth(巴特沃斯)滤波器、ChebyshevⅠ(切比雪夫Ⅰ型)、ChebyshevⅡ(切比雪夫Ⅱ型)滤波器和椭圆型滤波器。幅度归一化幅度响应的参数模型如图4-10所示(［15］中p．257)。

图4-10　幅度归一化幅度响应的参数模型图

4.3.1　巴特沃斯数字低通滤波器的设计

巴特沃斯逼近又称最平幅度逼近，一个N阶巴特沃斯低通滤波器的幅度平方函数可以定义为

其中N为整数，表示滤波器的阶次，Ω_c表示截止频率(3 dB截止频率)。

巴特沃斯低通滤波器的特点如下:

(1)当Ω =0时，，即在Ω =0处无衰减。

(2)当Ω = Ω_c时，

则其增益为3 dB，即Ω_c被称为3 dB截止频率。

(3)当Ω＜Ω_c时，通带内有最大平坦的幅度特性，随着Ω由 0变到Ω_c，单调减小，N越大，减小得越慢，也就是说通带内特性越平坦。

(4)当Ω＞Ω_c时，即在过渡带及阻带中，也随着Ω的单调增加而减少，但是因为＞1，故比通带内衰减的速度要快得多，N越大，衰减得越快。当Ω = Ω_s时，即为阻带的截止频率时，衰减为－20 log₁₀ ，它是阻带的最小衰减。

在Butterworth滤波器的设计中，对于已知通带截止频率Ω_p，阻带的截止频率Ω_s，最小通带幅度以及最大阻带波纹的情况，我们可以通过公式(4．18)对Butterworth滤波器的阶次进行估计，从而得到Butterworth的阶次N，

再利用冲激响应不变法或双线性变换法设计滤波器。

下面介绍MATLAB中的几个相关函数:

［N， Wn］= buttord(Wp，Ws，Rp，Rs， ’ s’)

说明:Wp代表通带截止频率，Ws代表阻带截止频率，Rp代表通带峰波，Rs代表最小阻带衰减。注意’ s’不能掉，’ s’代表设计的是模拟滤波器。

［b， a］= butter(N， Wn， ’ s’)

说明:N， Wn为Butterworth模拟低通滤波器的阶数和截止频率，b， a分别为此滤波器拉普拉斯变换的分子的系数向量和分母的系数向量，’ s’代表设计的是模拟滤波器。

［h， omega］= freqz(b1， a1， N)

说明:此函数返回以b1， a1为分子和分母系数向量的z变换的频率响应，h为响应，omega为对应的角频率值，N为omega的长度，默认值是N =512。

返回值N表示Butterworth滤波器的阶数和此阶数所对应的3 dB截止频率。

例4．3　用Butterworth滤波器法设计低通滤波器，给定的抽样频率为10 kHz，要求在频率小于1 kHz的通带内，幅度特性下降且小于1 dB，在频率大于1．5 kHz的阻带内，衰减大于15 dB。

可以算出

程序运行结果如图4-11所示。

图4-11　例4．3的Butterworth低通滤波器

这样b1， a1就是我们所需要的z变换的分子和分母的系数序列。

4.3.2　切比雪夫I数字低通滤波器的设计方法

一个N阶模拟低通Ⅰ型Chebyshev滤波器H_d(s)的幅度平方函数定义为

其中ε为小于1的正数，它是表示通带波纹带大小的一个参数，ε越大，波

纹也就越大，是Ω对Ω_p的归一化频率，Ω_p为通带截止频率，T_N(Ω)为N阶Chebyshev多项式，其定义为

上述多项式也可以通过下面的递归循环关系得到:

T_r(Ω)=2ΩT_r－1(Ω)－T_r－2(Ω)， r≥2，

其中T₀(Ω)=1，T₁(Ω)= Ω。

ChebyshevⅠ型多项式有如下的特性:

(1)当Ω =0时，

(2)当Ω = Ω_p时，

即所有阶次的ChebyshevⅠ型幅度函数都经过点，Ω_p为此滤波器的通带截止频率。

(3)在通带0≤Ω＜Ω_p内，即＜1时，在1～是等波纹起伏的。

(4)在通带之外，即Ω＞Ω_p时，随着Ω的增大，迅速单调趋向于0。

ChebyshevⅠ型的阶次可以由下式确定^［15］:

下面介绍MATLAB中的相关函数:

［N， Wn］= cheb1ord(Wp， Ws， Rp， Rs， ’ s’)

说明:N， Wn， Wp， Ws， Rp， Rs所代表的意思与上述情况相同，即此函数返回ChebyshevⅠ型的阶数和其截止频率值。 ’ s’代表设计的是模拟滤波器。

［b， a］= cheby1(N， Rp， Wn， ’ s’)

说明:此函数得到ChebyshevⅠ型系统函数拉普拉斯变换的分子系数向量和分母系数向量，Rp为最大通带衰减。 ’ s’代表设计的是模拟滤波器。

同样我们用上述例子:

例4．4　用ChebyshevⅠ型滤波器法设计低通滤波器，给定的抽样频率为10 kHz，要求在频率小于1 kHz的通带内，幅度特性下降且小于1 dB，在频率大于1．5 kHz的阻带内，衰减大于15 dB。

可以算出

程序运行结果如图4-12所示。

num和den便是所求的ChebyshevⅠ型滤波器Z变换的分子和分母系数。

图4-12　例4．4设计的ChebyshevⅠ型滤波器

4.3.3　切比雪夫II数字低通滤波器的设计方法

ChebyshevⅡ型模拟低通滤波器的幅度平方函数称为反Chebyshev逼近，一个N阶模拟低通ChebyshevⅡ型低通滤波器的幅度平方函数定义为

ChebyshevⅡ型函数不再像ChebyshevⅠ和Butterworth函数一样是全极点函数，它还拥有零点。

ChebyshevⅡ型函数的阶次同ChebyshevⅠ型一样可以由下式得出:

下面介绍MATLAB中的相关函数:

［N， Wn］= cheb2ord(Wp， Ws， Rp， Rs， ’ s’)

说明:N， Wn， Wp， Ws， Rp， Rs所代表的意思与上述情况相同，此函数返回ChebyshevⅡ型的阶数和其截止频率值。 ’ s’代表设计的是模拟滤波器。

［b， a］= cheby2(N， Rp， Wn， ’ s’)

说明:此函数得到ChebyshevⅡ型系统函数拉普拉斯变换的分子系数向量和分母系数向量，Rp为最大通带衰减。 ’ s’代表设计的是模拟滤波器。

同样我们用上述例子:

例4．5　用ChebyshevⅡ型滤波器法设计低通滤波器，给定的抽样频率为10 kHz，要求在频率小于1 kHz的通带内，幅度特性下降且小于1 dB，在频率大于1．5 kHz的阻带内，衰减大于15 dB。

可以算出

程序运行结果如图4-13所示。

图4-13　例4．5设计的ChebyshevⅡ型滤波器

num和den便是我们所求的ChebyshevⅠ型滤波器Z变换的分子和分母系数。

4.3.4　数字椭圆低通滤波器的设计

此外，还有一种常用的滤波器设计方法，即椭圆型滤波器设计方法。由于篇幅有限，此处不给出设计原理，直接给出其MATLAB函数:

［N， Wn］= ellipord(Wp， Ws， Rp， Rs， ’ s’);

［b， a］= ellip(N， Rp， Rs， Wn， ’ s’);

说明:函数ellipord是计算设计椭圆低通滤波器的一些参数，其中，Wp，Ws， Rp， Rs的含义如上所述，N为椭圆低通滤波器合适的阶数，而Wn返回的是此椭圆滤波器的通带截止频率。同样，函数ellip的输入参数的意义如上所述，b， a分别为所设计滤波器系统函数的分子、分母系数向量。

例4．6　设计椭圆数字滤波器。

通用MATLAB程序如下:

输入:Wp = 0．3，Ws = 0．5，Rp = 0．01，Rs =0．1，程序运行结果如图4-14所示。

b， a为其分子和分母的系数向量:

＞＞b， a

图4-14　例4．6的椭圆数字滤波器

4.3.5　数字高通、带通、带阻滤波器的设计

在进行数字高通、带通、带阻滤波器的设计时，通常有两种方法:一种是先将模拟低通滤波器变成模拟高通、模拟带通、模拟带阻滤波器，然后用冲激响应不变法或双线性变换法将其变成数字滤波器;另一种方法是先将模拟低通滤波器变成数字低通滤波器，然后再用频谱转换的方法将其转换成数字高通、数字带通、数字带阻滤波器。如图4-15所示。

1．先模拟频带变换，再数字化

(1)低通到高通的变换

模拟低通和模拟高通的变换关系为

其中s为模拟低通原型的拉普拉斯变量，Ω_p为此低通原型的通带截止频率，

图4-15　设计IIR数字滤波器的频率变换法

为模拟高通拉普拉斯变量，为对应于Ω_p的高通滤波器的截止频率，于是可以得到模拟低通到模拟高通的映射为

在模拟高通到数字高通的过程中，可以利用双线性变换，由，得

MATLAB中相关函数介绍:

［numt， dent］= lp2hp(num， den， Wo)

说明:此函数中num和den为原低通滤波器拉普拉斯变换的分子和分母的系数向量，Wo为原低通的通带截止频率，返回值numt和dent为对应高通滤波器的分子、分母的系数向量。

例4．7　设计高通滤波器，其参数如下:通带截止频率Fp =700 Hz，阻带Fs =500 Hz，通带最大衰减Rp为1 dB，阻带最小衰减为Rs =32 dB，抽样频率是2 kHz。

经过计算我们可以得到:

MATLAB程序如下:

得到的滤波器如图4-16所示。

图4-16　例4．7设计的高通滤波器

num和den就是其对应的Z变换的分子和分母系数向量，得到的值如下:

(2)低通到带通的变换

模拟低通到模拟带通的变换关系为

其中s是模拟低通拉普拉斯变量，

是模拟带通拉普拉斯变量，是模拟带通的几何中心频率。在虚轴上，令，s = iΩ，代入上式，可以得到

其中带宽。(4．26)和(4．27)两式化简可得

同时可以得到对应的阻带为

因此所求带通滤波器的系统函数为

从模拟低通到数字带通的转换过程中，仍然可以利用双线性变换法得到数字带通系统函数:

其中

注意，这里采用的是归一化的低通原型(截止频率为1)，而不需要做去归一化这一步骤。这是因为对于某一个具体的Ω_p(不一定为1)，设计滤波器的时候要用代替归一化原型中的s(即去归一化)，其中中的Ω_p就与D中的Ω_p相互抵消，因此，只需要用Ω_p=1来设计归一化原型即可。

MATLAB中相关函数介绍:

［numt， dent］= lp2bp(num， den， Wo， Bw)

说明:num和den是原模拟低通滤波器的拉普拉斯变换的分子和分母的系数向量。Wo和Bw分别是带通滤波器的中心频率和带宽，numt和dent表示转换后的带通滤波器的拉普拉斯变换的分子系数向量和分母系数向量。

例4．8　设计数字带通滤波器，其参数如下所示:

ω_p1=0．45π， ω_p2=0．65π， ω_s1=0．3π， ω_s2=0．75π，

R_p=1 dB， R_s=40 dB．

MATLAB程序如下:

得到的滤波器如图4-17所示。

图4-17　例4．8的数字带通滤波器

num和den为滤波器Z变换的分子、分母系数向量，其值如下:

(3)低通到带阻的变换

模拟低通到模拟带阻的变换关系为

或者

由

其中

可得

这样可以得到模拟带阻滤波器的系统函数:

从模拟带阻到数字带阻的变换中，用双线性变换法，同样可以得到

其中

原理同上，这里取Ω_s=1即可。

MATLAB中相关函数介绍:

［numt， dent］= lp2bs(num， den， Wo， Bw)

说明:num和den是原模拟低通滤波器拉普拉斯变换的分子、分母系数向量。Wo和Bw是带阻滤波器的中心频率和阻带宽，numt和dent表示转换后带阻滤波器的拉普拉斯变换的分子系数向量和分母系数向量。

例4．9　设计带阻滤波器，其参数如下:

ω_s1=0．45π， ω_s2=0．65π， ω_p1=0．3π， ω_p2=0．75π，

R_p=1 dB， R_s=40 dB．

MATLAB程序如下:

得到的滤波器如图4-18所示。

图4-18　例4．9设计的数字带阻滤波器

num和den为其滤波器Z变换的分子、分母系数向量，其值如下:

2．先数字化，再进行数字频带变换

在实际应用中，有时需要在不重新设计滤波器的情况下，改变滤波器的一些特性以满足新的指标。如果给定数字滤波器的低通原型系统函数为H_L(z)，则可以通过一定的变换来设计其他各种不同类型的数字滤波器系统函数 ，这种变换是将的z平面变换到平面。

(1)数字低通到数字低通

两者的变换函数为一个全通函数(详细步骤参见［1］ ，第297～300页):

其中

由此，通过(4．44)和(4．45)两式，可以由已有的数字低通滤波器H_L(z)(截止频率为θ_c)变换得到新的系统函数为H_D(z)的低通滤波器(截止频率为ω_c)，即

例4．10　考虑如下的系统函数:

在通带0～0．25π内的最大衰减为0．5 dB，而在阻带0．55π～π内的最小衰减为15 dB，要求设计一个通带截止频率为0．35π而其他性能指标不变的新的数字低通滤波器。

MATLAB程序如下:

这样，新的低通滤波器的Z变换为

(2)数字低通到数字高通、带通、带阻滤波器

数字低通变换成其他各种类型滤波器(如数字高通、数字带通、数字带阻滤波器)的变换关系式如表4-1所示。

表4-1　由截止频率为θ_c的低通滤波器H_L(z)变换成各种类型滤波器H_D(z)

例4．11　利用频谱变换法将例4．10所表示的三阶低通滤波器转换成一个数字高通滤波器。要求数字滤波器通带截止频率ω_c=0．55π，其他性能指标不变。

MATLAB程序如下:

这样，H1就是所求新的滤波器，即