1、1数字信号处理实验指导实验四、 无限长单位脉冲响应(IIR)滤波器的设计方法（一） 实验目的加深对无限冲激响应（ IIR ）数字滤波器的常用指标和设计过程的理解。（二） 实验内容常用函数介绍：1、Matlab 信号处理工具箱中提供了设计巴特沃思模拟滤波器的函数buttord、buttap 和 butter，格式如下：（1） ,(,CPSNWbutordRs用于计算巴特沃思模拟低通滤波器的阶 N 和 3dB 截止频率 Wc (即本书中的符号 )。其中，Wp 和 Ws 分别是滤波器的通带截止频率 和阻止截止频率c p，单位为 rad/s；Rp 和 Rs 分别是通带最大衰减系数 和阻带最小衰减系数s
2、 ，单位为 dB。s（2） ,()zpGbutaN用于计算 N 阶巴特沃思归一化（ =1）模拟低通滤波器系统函数的零、极c点和增益因子，返回长度为 N 的向量 z 和 p 分别给出 N 个零点和极点，G 是滤波器增益。得到的滤波器系统函数形式如下： 12()()a NaQszsHGPp其中， 和 分别是向量 z 和 p 的第 k 个元素。如果要从零、极点得到系kz统函数的分子和分母多项式系数向量 B 和 A，可以调用结构转换函数（3） ，结构转换后系统函数的形式为,2(,)BAptf11)(MMaNNbssbHsaa其中，M 是向量 B 的长度， N 是向量 A 的长度， 分别是向量 B 和
3、A 的ka和2第 k 个元素。（3） ,(,)CBAbuterNWftypes用于计算巴特沃思模拟滤波器系统函数中分子和分母多项式系数向量 B 和 A，其中 N 和 分别是滤波器的阶和 3dB 截止频率 ，返回向量 B 和 A 中的元素C c和 分别是上面的 表示式中的分母和分子系数。ftype 缺省时，设计kab()aHs低通滤波器;ftype=high 时，设计高通滤波器；ftype=stop 时，设计带阻滤波器，此时 为向量 ，且 ftype 缺省时，设计带通滤波器，带通的频率CW1,cu区间为 。S 缺省时，设计数字滤波器。1ccu例如：设计一个满足下列指标要求的巴特沃思模拟低通滤波器
4、。指标 sradBsraddspp /4,30,/2, BMatlab 程序：;41;,(,)0.6(pspsCPSCWRNbutorWsBAeHfrqs； %滤 波 器 参 数； 计 算 滤 波 器 的 阶 和 3dB截 止 频 率； 计 算 滤 波 器 系 统 函 数 分 子 、 分 母 多 项 式 系 数： ： ； 在 06rad/s频 率 范 围 内 取 10个 频 率 点， ,;2*log1)(,)/;(/)BAWabHptidxlabersyledB 计 算 频 率 向 量 上 的 滤 波 器 频 率 响 应 纵 坐 标 的 单 位 取 n %绘 制 幅 频 响 应 曲 线频 率
5、（ 幅 度运行结果如下：N=6Wc=2.2496B= 0 0 0 0 0 0 129.5917A=1.0000 8.6916; 37.7720; 104.0667;191.1447 ;222.5973 ; 129.5917即： 65432129.57()8.9137.26.4.579312.7aHssssss幅频响应曲线如图 6.3.4 所示。2、Matlab 中设计数字滤波器的函数都是采用双线性变换法，将模拟滤波器转换为数字滤波器。这些函数及其凋用格式如下（巴特沃思数字滤波器）：（1） ,(,)CPSNWbutordR该格式用于计算巴特沃思数字滤波器的阶 N 和 3dB 截止频率的归一化值3
6、（关于 归一化） 。调用参数 Wp 和 Ws 分别是数字滤波器的通带截止频率和CW阻带截止频率的归一化值（关于 归一化） ，要求 和 ，其中01pWs1 表示数字频率 （对应模拟频率 ， 为采样频率） 。 和 分别是通带/2sfsfRs最大衰减和阻带最小衰减，单位 dB。当 时，设计高通滤波器；当 和sp p是二元向量时，设计带通（ ）或带阻（ ）滤波器，这时返回sWspWsp参数 也是二元向量。C（2） ,(,)CBAbuterNftye； 该格式用于计算巴特沃思模拟滤波器系统函数中分子和分母多项式系数向量B 和 A。调用参数 N 和 WC分别是巴特沃思数字滤波器的阶和 3dB 截止频率归一
7、化值（关于 归一化） 。当 ftype 缺省时，设计低通滤波器；当 ftype=high 时，设计高通滤波器；ftype=stop 时，设计带阻滤波器，此时 为二元向量CW， 和 分别是带阻滤波器的通带 3dB 下，上截止频率的归一化值1,cuWccu（关于 归一化） ； 为向量 ，且 ftype 缺省时，设计带通滤波器，C1,cuW带通的频率区间为 。注意设计出的带通和带阻数字滤波器是 2N 阶的，cu这是因为带通滤波器可表示为一个 N 阶低通滤波器与一个 N 阶高通滤波器的级联。由函数的返回向量 B 和 A 可写出数字滤波器的系统函数为：112() NNbzbzHzaa其中， 分别是向量
8、B 和 A 的第 k 个元素。k和例如：用双线性变换法设计巴特沃思数字低通滤波器，指标要求通带截止频率 ，通带最大衰减 ；0.2pwrad1pdB阻带截止频率 ，阻带最小衰减35s 0s解 设计步骤如下：（1）给出数字滤波器的指标。40.2;.35;10;PSSWR（2）计算巴特沃思数字滤波器的阶 N 和 3dB 截止频率。；,(,)CPSNbutordW（3）用双线性变换法设计巴特沃思数字低通滤波器。；,(,)ZCBAte运行结果如下：=0.0335 0.1006 0.1006 0.0335Z=1.000 -1.4245 0.8827 -1.1900A设计出的数字滤波器系统函数为 1230.
9、35.06.0.5()1428719zzzHz程序如下： 0.;.35;10;,(,);PSPSCZCWRNbutordWBAe %滤 波 器 参 数计 算 数 字 滤 波 器 的 阶 和 3dB截 止 频 率 计 算 数 字 滤 波 器 系 统 函 数 分 子 和 分 母 多 项 式 的 系 数 向 量 和 A%省 去 绘 图 部 分 程 序 例：利用脉冲响应不变法，把系统函数为 的模拟滤波器651)(2ssHa变换成等价的数字滤波器，采样间隔 ST1.0手算： 23651)(2sssHa2111213 605.59.80)( zzzeTzezTMatlab:B=1,1;A=1,5,6;T=
10、0.1;5Fs=1/T;Bz,Az=impinvar(B,A,Fs);%用脉冲响应不变法将模拟滤波器变换成数字滤波器运行结果：Bz=0.1000 -0.0897Az=1.000 -1.5595 0.60651、 设采样周期 T=250s（采样频率 fs =4kHz），用脉冲响应不变法和双线性变换法设计一个三阶巴特沃兹滤波器，其 3dB 边界频率为 fc =1kHz。 B,A=butter(3,2*pi*1000,s); num1,den1=impinvar(B,A,4000); h1,w=freqz(num1,den1); B,A=butter(3,2/0.00025,s); num2,den
11、2=bilinear(B,A,4000); h2,w=freqz(num2,den2); f=w/pi*2000; plot(f,abs(h1),-.,f,abs(h2),-); grid; xlabel(频率 /Hz ) ylabel(幅值 /dB) 程序中第一个 butter 的边界频率 21000，为脉冲响应不变法原型低通滤波器的边界频率；第二个 butter 的边界频率 2/T=2/0.00025，为双线性变换法原型低通滤波器的边界频率.图 3.1 给出了这两种设计方法所得到的频响，虚线为脉冲响应不变法的结果；实线为双线性变换法的结果。脉冲响应不变法由于混叠效应，使得过渡带和阻带的衰减
12、特性变差，并且不存在传输零点。同时，也看到双线性变换法，在 z=-1 即 = 或 f=2000Hz 处有一个三阶传输零点，这个三阶零点正是模拟滤波器在 =处的三阶传输零点通过映射形成的。62、 设计一数字高通滤波器，它的通带为 400500Hz，通带内容许有0.5dB 的波动，阻带内衰减在小于 317Hz 的频带内至少为 19dB，采样频率为1,000Hz。 wc=2*1000*tan(2*pi*400/(2*1000); wt=2*1000*tan(2*pi*317/(2*1000); N,wn=cheb1ord(wc,wt,0.5,19,s); B,A=cheby1(N,0.5,wn,hi
13、gh,s); num,den=bilinear(B,A,1000); h,w=freqz(num,den); f=w/pi*500; plot(f,20*log10(abs(h); axis(0,500,-80,10); 7grid; xlabel() ylabel(幅度/dB) 图 3.2 给出了 MATLAB 计算的结果，可以看到模拟滤波器在 =处的三阶零点通过高通变换后出现在 =0（z=1）处，这正是高通滤波器所希望得到的。 3、 设计一巴特沃兹带通滤波器，其dB 边界频率分别为 f2=110kHz 和f1=90kHz，在阻带 f3 = 120kHz 处的最小衰减大于dB，采样频率fs=
14、400kHz。8w1=2*400*tan(2*pi*90/(2*400); w2=2*400*tan(2*pi*110/(2*400); wr=2*400*tan(2*pi*120/(2*400); N,wn=buttord(w1 w2,0 wr,3,10,s); B,A=butter(N,wn,s); num,den=bilinear(B,A,400); h,w=freqz(num,den); f=w/pi*200; plot(f,20*log10(abs(h); axis(40,160,-30,10); grid; xlabel(频率/kHz) ylabel(幅度/dB) 图 3.3 给出
15、了 MATLAB 计算的结果，可以看出数字滤波器将无穷远点的二阶零点映射为 z=1 的二阶零点，数字带通滤波器的极点数是模拟低通滤波器的极点数的两倍。 94、 一数字滤波器采样频率 fs = 1kHz，要求滤除 100Hz 的干扰，其dB 的边界频率为 95Hz 和 105Hz，原型归一化低通滤波器为 w1=95/500; w2=105/500; B,A=butter(1,w1, w2,stop); h,w=freqz(B,A); f=w/pi*500; plot(f,20*log10(abs(h); axis(50,150,-30,10); grid; xlabel(频率/Hz) ylabel(幅度/dB)图 3.4 为 MATLAB 的计算结果 （三） 实验报告要求101、简述实验目的和实验原理，用笔算求出设计的滤波器，并用计算机计算结果验证，总结实验中的主要结论，写出收获和体会。2、编写教材 P193(1t)和 P194(11t)的 MATLAB 计算程序。
