IIR滤波器的设计方法

实验六IIR滤波器的设计方法■实验学时：2学时■实验类型：设计、验证一、实验目的次数字滤波器的常用指标和设计过程的理解；响应不变法和双线性变换法设计IIR滤波器的原理，设计方法，步骤;现象，比较脉冲响应不变法和双线性变换法设计IIR滤波器的特点熟悉乂21匕山信号处理工具箱中的常用函数二、实验原理与方法1次数字滤波器设计步骤：(一)先设计模拟低通原型滤波器模拟滤波器设计的步骤♦确定滤波器的技术指标：通带截止频率、通带衰减阻带截止频率、阻带衰减♦将技术指标设计转化为低通原型滤波器幅度平方函数的参数-ButterwortH低通滤波器(MQ。-CHebysHev低通滤波器=Q「)♦构造归一化低通原型滤波器的系统函数H°(p)♦反归一，求得H0(s)=珥①心皿,簿为归一化参考频率(二)八尸数字化为DF1．脉冲响应不变法用数字滤波器的单位脉冲响应序列h(n)模仿模拟滤波器的冲激响应h(t)，让ah(n)正好等于h(t)的采样值，即h(n)=h(nT)，其中T为采样间隔，如果以H(s)及H(Z)分别表示ha(t)的拉式变换及h(n)的Z变换，则TOC\o"1-5"\h\z1f 2人、=—ZH(s+j——m)\o"CurrentDocument"z=eSTT a Tm=-82．双线性变换法S平面与z平面之间满足以下映射关系：s=2.匕£1 (***)T1+z-is平面的虚轴单值地映射于z平面的单位圆上，s平面的左半平面完全映射到z平面的单位圆内。双线性变换不存在混叠问题。双线性变换是一种非线性变换(源(万)=方)，这种非线性引起的幅频特性畸变可通过预畸而得到校正。双线性变换法特别适合用于设计常见的选频性滤波器，下面以双线性变换法设计低通数字滤波器为例，将设计步骤归纳如下：(1)根据需要确定数字滤波器的性能指标:通带截止频率了,、阻带截止频率fp>通带最大衰减ap、阻带内的最小衰减as，采样频率(；(2)确定相应的数字角频率，3=2裾T；3=2亏T；(3)计算经过预畸的相应模拟低通原型的截止频率，23 23。=—tg(k)，。=—tg(f)；pT2sT2(4)根据。p，ap和。s，as计算模拟低通原型滤波器的传递函数Ha(s)；(5)用上述(***)表示的双线性变换公式代入H(s)，求出H(z)；分析滤波器特性，检查其指标是否满足要求。a3、IIR低通、高通、带通数字滤波器设计采用双线性原型变换公式：变换类型变换关系式备注低通L3妙产2)田=2^fT高通s=w仁丁烟目”带通z2-2zcos4+1___cos4-cosmA Z-1 ”k1sin(由]十田(2):可山口= ，也/£in叼+sin田口为带通的上下边带临界频率相关mtlab函数介绍1、Matlab信号处理工具箱中提供了设计巴特沃思模拟滤波器的函数，格式如下:用于设计巴特沃思模拟滤波器H(s)设计的函数模拟滤波器传输函数表示为：H(s)=B(s)-=俨一1++久i"bMa A(S)asn-1++as+a1 …N-1 N其中M是向量B的长度，N是向量A的长度，bk和ak分别是B和A的第k个元素。(1)[N,W]=buttord(W,W,R,R,'s')用于根据设计指标计算巴特沃思:ALPF幅度平方函数的参数N和3dB截止频率Wc(即本书中的符号QJ。其中，Wp和Ws分别是滤波器的通带截止频率Qp和阻止截止频率Qs，单位为;ad/s；Rp和Rs分别是通带最大衰减系数ap和阻带最小衰减系数

，。(2)[B,A]=butter(N,W,'ftype','s')用于根计算巴特沃思ALPF系统函数Hq(s)中分子和分母多项式系数向量B和A,a其中N和7分别是滤波器的阶和3dB截止频率Q，返回向量B和A，其中ak和b分别是上面的Ha(s)表示式中的分母和分子系数。'ftype缺省时，设计低通滤波器;ftype=high时，a设计高通滤波器；ftype=stop时，设计带阻滤波器，此时/为两个元素的向量[W1,W]，且ftype缺省时，设计带通滤波器，带通的频率区间为WKQ〈W。5缺省时，设计数字滤波器。 cc例如：设计一个满足下列指标要求的巴特沃思模拟低通滤波器。指标=4rad/ss%滤波器参数%计算滤波器的阶和3dB截止频率%计算滤波器系统函数分子、分母多项式系数%在0~6rad/s频率范围内取1001个频率点%计算频率向量W上的滤波器频率响应=4rad/ss%滤波器参数%计算滤波器的阶和3dB截止频率%计算滤波器系统函数分子、分母多项式系数%在0~6rad/s频率范围内取1001个频率点%计算频率向量W上的滤波器频率响应%纵坐标的单位取dB%绘制幅频响应曲线Matlab程序:p sW=2;W=4;R=1;R=30；[N,W]=buttord(W,W,R,R/s')；[B,A]=butter(N,W:s')；W=0：0.001：6；[H，W]=freqs(B,A,W);H=20*log10(abs(H));plot(W,H),gridonxlabel('频率/(rad/s)';ylabel('幅度/dB')运行结果如下：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即：HaHa(s)=s6+8.6916s5+37.7720s4+104.0667s3+191.1447s2+222.5793s+129.5917幅频响应曲线(略)。2、Matlab中设计数字滤波器的函数都是采用双线性变换法，将模拟滤波器转换为数字滤波器。这些函数及其凋用格式如下(以巴特沃思数字滤波器为例)一 b+bz-1++bz-1+bz-n数字滤波器的系统函数为：H(z)=T—2z Nz一如二a+az-1++az-1+az-N其中，b和ak分别是向量B和A的第3元素。:N N+1(1)[N,W]=buttord(W,W,R,R)该格式用于‘计算巴特沃思氤字滤波器的阶N和3dB截止频率的归一化值W(关于C

—八7。Jl/心IZ^VMHIQ«^LriH ]]、/八IIcI-L_uih ]]、/八I的归一化值（关于冗归一化），要求0<W<1和0<W<1,其中1表示数字频率兀（对ps应模拟频率fs/2,fs为采样频率）。Rp和R分别是通带最大衰减和阻带最小衰减，单位dB。当诂<Wp时，设计高通滤波器；当Wp和W是二元向量时，设计带通（Ws<W）或带阻（Ws2W）滤波器，这时返回参数也是二元向量。（2）[B,A]=butter（N,W,'ftype'）；该格式用于计算巴特沃点模拟滤波器系统函数中分子和分母多项式系数向量B和A。调用参数N和WC分别是巴特沃思数字滤波器的阶和3dB截止频率归一化值（关于兀归一化）。当ftype缺省时，设计低通滤波器；当ftype=high时，设计高通滤波器；ftype=stop时，设计带阻滤波器，此时%为二元向量[W1,W]，W1和W分别是带阻滤波器的通带3dB下，上截止频率的归一化值（关于c兀归一化）；也为向量[W1,W]，且ftype缺省时，设计带通滤波器，带通的频率区间为W1®〈W:注意设计出的带通和带阻数字滤波器是2N阶的，这是因为带通滤波器可表示为一个N阶低通滤波器与一个N阶高通滤波器的级联。（3）设计实例例如：用双线性变换法设计巴特沃思数字低通滤波器,指标要求通带截止频率,=0.2冗rad，通带最大衰减a广1dB；阻带截止频率W=0.35冗rad，阻带最小衰减0P=10dB解设计步骤加下： s（1）给出数字滤波器的指标。W=0.2;W=0.35;R=10;（2）计算巴特反思数字滤波器的阶N和3dB截止频率。[N,W]=buttord（W,W,R,R）；（3）用双线性变换法设计巴特沃思数字低通滤波器。[B,A]=butter（N,W）；运行结果如下： CBZ=0.0335 0.1006 0.1006 0.0335AZ=1.000 -1.4245 0.8827 -1.1900所以,设计出的数字滤波器系统函数为0.0335+0.1006z-1+0.1006z-2+0.0335z-3程序如下：W=0.2;W=0.35;R=10;R=10;[N,W]=buttord程序如下：W=0.2;W=0.35;R=10;R=10;[N,W]=buttord(W,W,R,R);[BZ,AZ]=butter(N,WC);%省去绘图部分程序%滤波器参数%计算数字滤波器的阶和3dB截止频率%计算数字滤波器系统函数分子和分母多项式的系数向量B和A三、实验内容1、设采样周期T=250〃s(采样频率fs=4kHz),用脉冲响应不变法和双线性变换法设计一个三阶巴特沃兹滤波器，其3dB边界频率为fc=1kHz。参考程序如下：%脉冲响应不变法设计H⑵[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⑵[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”按降幂的次序排列，且其传递函数来自下式：_二⑸_E工店品+占〔2注品一工+…+五山+L

'-A(s)~sn+a(2)sn-1+■■■+a(n+1)(2)在MATLAB中，可分别用impinvar和历1册2「实现从模拟滤波器到数字滤波器的脉冲响应不变映射和双线性映射，调用格式为：[b,a]=impinvar(c,d,fs)；脉冲响应不变法将H(s)fH⑵[b,a]=bilinear(c,d,£$)；双线性变换法将H(s)fH⑵其中，c、d分别为模拟滤波器H(s)的分子和分母多项式系数向量；fs为采样频率。缺省值fs=1Hz；b、a分别为数字滤波器H⑵分子和分母多项式系数向量。(3)求解数字滤波器的频率响应[H,w]=freqz(B,A,N)B和A分别为离散系统的系统函数H⑵分子、分母多项式的系数向量，N为正整数，返回量H则包含了离散系统频率响应在0〜pi范围内N个频率等分点的值，向量w则包含范围内N个频率等分点。调用中若N默认，默认值为512。(4)程序说明程序中第一个butter的边界频率2nX1000rad/s,为脉冲响应不变法原型低通滤波器的边界频率；第二个butter的边界频率2/T=2/0.00025,为双线性变换法原型低通滤波器的边界频率，此值是通过预畸变得到的。f=1000Hz, Q=2对, ①=QT,经过预畸变得到的模拟角频率：Q=3tg[巴]=3tg(0.25兀)=-cT(2)T T边界频率；第二个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('')丫匕6。1('幅度分8')3、设计一巴特沃兹带通滤波器，其3dB边界频率分别为f2=110kHz和f1=90kHz,在阻带f3=80kHz和f4=120kHz处的最小衰减大于10dB,采样频率fs=400kHz。参考程序：w1=2*400*tan(2*pi