IIR滤波器的设计方法

1实验六IIR滤波器的设计方法实验学时：2学时实验类型：设计、验证一、实验目的IR数字滤波器的常用指标和设计过程的理解；响应不变法和双线性变换法设计IIR滤波器的原理，设计方法，步骤；现象，比较脉冲响应不变法和双线性变换法设计IIR滤波器的特点熟悉Matlab信号处理工具箱中的常用函数二、实验原理与方法IIR数字滤波器设计步骤：（一）先设计模拟低通原型滤波器（二）AF数字化为DF1．脉冲响应不变法用数字滤波器的单位脉冲响应序列)(nh模仿模拟滤波器的冲激响应)(tha,让)(nh正好等于)(tha的采样值，即)()(nThnha，其中T为采样间隔，如果以)(sHa及)(zH分别表示)(tha的拉式变换及)(nh的Z变换，则)2(1)(mTjsHTzHmaezsT2．双线性变换法S平面与z平面之间满足以下映射关系：11112zzTs（***）s平面的虚轴单值地映射于z平面的单位圆上，s平面的左半平面完全映射到z平面的单位圆内。双线性变换不存在混叠问题。2双线性变换是一种非线性变换，这种非线性引起的幅频特性畸变可通过预畸而得到校正。双线性变换法特别适合用于设计常见的选频性滤波器，下面以双线性变换法设计低通数字滤波器为例，将设计步骤归纳如下：（1）根据需要确定数字滤波器的性能指标：通带截止频率pf、阻带截止频率pf、通带最大衰减αp、阻带内的最小衰减αs，采样频率sf；(2)确定相应的数字角频率，Tfpp2；Tfss2；(3)计算经过预畸的相应模拟低通原型的截止频率，)2(2pptgT，)2(2sstgT；(4)根据Ωp，αp和Ωs，αs计算模拟低通原型滤波器的传递函数)(sHa；(5)用上述（***）表示的双线性变换公式代入)(sHa，求出)(zH；分析滤波器特性，检查其指标是否满足要求。3、IIR低通、高通、带通数字滤波器设计采用双线性原型变换公式：变换类型变换关系式备注低通高通带通为带通的上下边带临界频率相关matlab函数介绍1、Matlab信号处理工具箱中提供了设计巴特沃思模拟滤波器的函数，格式如下：用于设计巴特沃思模拟滤波器Ha(s)设计的函数模拟滤波器传输函数表示为：111111()()()MMMaNNNbsbsbBsHsAsasasa其中M是向量B的长度，N是向量A的长度，kkba和分别是B和A的第k个元素。（1）[,](,,,,CPSPSNWbuttordWWRRs用于根据设计指标计算巴特沃思ALPF幅度平方函数的参数N和3dB截止频率Wc(即本书中的符号c)。其中，Wp和Ws分别是滤波器的通带截止频率p和阻止截止频率s，单位为rad/s；Rp和Rs分别是通带最大衰减系数p和阻带最小衰减系数3s，单位为dB。（2）[,](,,,)CBAbutterNWftypes用于根计算巴特沃思ALPF系统函数Ha(s)中分子和分母多项式系数向量B和A，其中N和CW分别是滤波器的阶和3dB截止频率c，返回向量B和A，其中ka和kb分别是上面的()aHs表示式中的分母和分子系数。ftype缺省时，设计低通滤波器；ftype=high时，设计高通滤波器；ftype=stop时，设计带阻滤波器，此时CW为两个元素的向量[1,]ccuWW，且ftype缺省时，设计带通滤波器，带通的频率区间为1ccuWW。S缺省时，设计数字滤波器。例如：设计一个满足下列指标要求的巴特沃思模拟低通滤波器。指标sraddBsraddsspp/4,30,/2,1BMatlab程序：2;4;1;30[,](,,,,[,](,,)00.0016[](pspsCPSPSCWWRRNWbuttordWWRRsBAbutterNWsWHWfreqs；%滤波器参数；%计算滤波器的阶和3dB截止频率；%计算滤波器系统函数分子、分母多项式系数：：；%在0~6rad/s频率范围内取1001个频率点，,,);20*log10(());(,),(//);(/)BAWWHabsHBplotWHgridxlabelradsylabeldB%计算频率向量上的滤波器频率响应%纵坐标的单位取don%绘制幅频响应曲线频率（幅度运行结果如下：N=6Wc=2.2496B=000000129.5917A=1.00008.6916;37.7720;104.0667;191.1447;222.5973;129.5917即：65432129.5917()8.691637.7720104.0667191.1447222.5793129.5917aHsssssss幅频响应曲线（略）。2、Matlab中设计数字滤波器的函数都是采用双线性变换法，将模拟滤波器转换为数字滤波器。这些函数及其凋用格式如下（以巴特沃思数字滤波器为例）数字滤波器的系统函数为：1112111121()NNNNNNbbzbzbzHzaazazaz其中，kkba和分别是向量B和A的第k个元素。（1）[,](,,,)CPSPSNWbuttordWWRR该格式用于计算巴特沃思数字滤波器的阶N和3dB截止频率的归一化值CW（关于4归一化）。调用参数Wp和Ws分别是数字滤波器的通带截止频率和阻带截止频率的归一化值（关于归一化），要求01pW和01sW，其中1表示数字频率（对应模拟频率/2sf，sf为采样频率）。pR和sR分别是通带最大衰减和阻带最小衰减，单位dB。当spWW时，设计高通滤波器；当pW和sW是二元向量时，设计带通（spWW）或带阻（spWW）滤波器，这时返回参数CW也是二元向量。（2）[,](,,)CBAbutterNWftype；该格式用于计算巴特沃思模拟滤波器系统函数中分子和分母多项式系数向量B和A。调用参数N和WC分别是巴特沃思数字滤波器的阶和3dB截止频率归一化值（关于归一化）。当ftype缺省时，设计低通滤波器；当ftype=high时，设计高通滤波器；ftype=stop时，设计带阻滤波器，此时CW为二元向量[1,]ccuWW，1cW和cuW分别是带阻滤波器的通带3dB下，上截止频率的归一化值（关于归一化）；CW为向量[1,]ccuWW，且ftype缺省时，设计带通滤波器，带通的频率区间为1ccuWW。注意设计出的带通和带阻数字滤波器是2N阶的，这是因为带通滤波器可表示为一个N阶低通滤波器与一个N阶高通滤波器的级联。（3）设计实例例如：用双线性变换法设计巴特沃思数字低通滤波器，指标要求通带截止频率0.2pwrad，通带最大衰减1pdB；阻带截止频率0.35swrad，阻带最小衰减10sdB解设计步骤如下：（1）给出数字滤波器的指标。0.2;0.35;10;PSSWWR（2）计算巴特沃思数字滤波器的阶N和3dB截止频率。[,](,,,)CPSPSNWbuttordWWRR；（3）用双线性变换法设计巴特沃思数字低通滤波器。[,](,)ZZCBAbutterNW；运行结果如下：ZB=0.03350.10060.10060.0335ZA=1.000-1.42450.8827-1.1900所以，设计出的数字滤波器系统函数为1231230.03350.10060.10060.0335()1.0001.42450.88270.1900zzzHzzzz程序如下：0.2;0.35;10;10;[,](,,,);[,](,);PSPSCPSPSZZCWWRRNWbuttordWWRRBAbutterNW%滤波器参数%计算数字滤波器的阶和3dB截止频率%计算数字滤波器系统函数分子和分母多项式的系数向量B和A%省去绘图部分程序5三、实验内容1、设采样周期T=250μs（采样频率fs=4kHz），用脉冲响应不变法和双线性变换法设计一个三阶巴特沃兹滤波器，其3dB边界频率为fc=1kHz。参考程序如下：%脉冲响应不变法设计H(z)[B,A]=butter(3,2*pi*1000,'s');%N=3,截止频率wc=2*pi*1000[num1,den1]=impinvar(B,A,4000);[h1,w]=freqz(num1,den1);%双线性变换法设计H(z)[B,A]=butter(3,2/0.00025,'s');%N=3,预畸变后的截止频率wc=2/T=2/0.00025[num2,den2]=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')程序说明;(1)[B,A]=butter(n,Wn,'s')，该函数设计一个n阶的低通模拟的Butterworth滤波器，其截止角频率为Wnrad/s。它转换得到滤波器的系数A、B是n+1行的向量，其中“s”按降幂的次序排列，且其传递函数来自下式：(2)在MATLAB中，可分别用impinvar和bilinear实现从模拟滤波器到数字滤波器的脉冲响应不变映射和双线性映射，调用格式为：[b,a]=impinvar(c,d,fs)；脉冲响应不变法将H(s)→H(z)[b,a]=bilinear(c,d,fs)；双线性变换法将H(s)→H(z)其中，c、d分别为模拟滤波器H(s)的分子和分母多项式系数向量；fs为采样频率。缺省值fs=1Hz；b、a分别为数字滤波器H(z)分子和分母多项式系数向量。(3)求解数字滤波器的频率响应[H,w]=freqz(B,A,N)B和A分别为离散系统的系统函数H(z)分子、分母多项式的系数向量，N为正整数，返回量H则包含了离散系统频率响应在0~pi范围内N个频率等分点的值，向量w则包含范围内N个频率等分点。调用中若N默认，默认值为512。（4）程序说明程序中第一个butter的边界频率2π×1000rad/s，为脉冲响应不变法原型低通滤波器的边界频率；第二个butter的边界频率2/T=2/0.00025，为双线性变换法原型低通滤波器的边界频率，此值是通过预畸变得到的。T2)25.tg(0T22tgT2,,2,1000cccccccTfHzf角频率：经过预畸变得到的模拟6程序中第一个butter的边界频率2π×1000，为脉冲响应不变法原型低通滤波器的边界频率；第二个butter的边界频率2/T=2/0.00025，为双线性变换法原型低通滤波器的边界频率。2、设计一数字切比学雪夫高通滤波器，它的通带截止频率为400Hz，通带内容许有0.5dB的波动(即通带最大衰减为0.5dB)，阻带截止频率为317Hz，阻带最小衰减为19dB，采样频率为1000Hz。参考程序如下：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,'high','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]);grid;xlabel('')ylabel('幅度/dB')3、设计一巴特沃兹带通滤波器，其３dB边界频率分别为f2=110kHz和f1=90kHz，在阻带f3=80kHz和f4=120kHz处的最小衰减大于１０dB，采样频率fs=400kHz。参考程序：w1=2*400*tan(2*pi*90/(2*400));%频率预畸变w2=2*400*tan(