信号与信号处理实验八IIR数字滤波器的设计

实验八IIR数字滤波器的设计按频率选择特性不同，滤波器可以分为低通、高通、带通、带阻等类型。|)(|jH|)(|jH|)(|jH|)(|jH|)(|jeH000002|)(|jeH02|)(|jeH02|)(|jeH02低通高通带通带阻由于抽样过程造成频率特性延拓，数字滤波器所有的频率特性位于数字角频率轴上的0－2π区间。其低频频率成分位于π的偶数倍附近，而高频成分位于π的奇数倍附近。理想滤波器对应的时域响应为非因果的，因而实际的滤波器频率特性是用一个具有因果冲激响应的系统函数对理想滤波器频率特性的逼近。|)(|jeH通带阻带过渡带1112cts11|()|1jHe2|)(|jeHtc||s过渡带c||通带c通带截止频率1通带容限||ts阻带ts阻带截止频率2阻带容限低通滤波器幅频特性IIR滤波器设计过程按设计任务，确定滤波器性能要求，制定技术指标。将数字滤波器指标转换成模拟滤波器指标。设计一个满足设计指标的模拟滤波器系统函数H(s)逼近此性能指标。将设计的模拟滤波器系统函数H(s)转换成数字滤波器系统函数H(z)。在MATLAB中，经典设计IIR数字滤波器的步骤：（1）根据给定的性能指标和方法，首先对设计性能指标中的频率指标进行转换，转换后的频率指标作为模拟滤波器原型设计指标；（2）估计模拟滤波器最小阶数和边界频率，可利用MATLAB工具函数buttord，cheb1ord,cheb2ord等。格式：[N，Wc]=buttord(wp,ws,Rp,Rs,’s’)wp,ws,Wc均以弧度/秒为单位。（3）设计模拟低通滤波器原型，利用函数buttap,cheb1ap等。格式：[z,p,k]=buttap(N)功能：返回N阶归一化原型Butterworth模拟滤波器的零极点增益模型。利用zp2tf函数可求出滤波器的传递函数模型。[B，A]=zp2tf[z,p,k](4)由模拟低通原型经频率变换得到模拟滤波器（低通、高通，带通、带阻），函数为lp2lp,lp2hp,lp2bp,lp2bs.格式：[Bt，At]=lp2lp(B,A,Wc)(5)将模拟滤波器离散化得到IIR数字滤波器，函数为bilinear,impinvar.格式：[num2,den2]=bilinear(Bt,At,fs)例：试设计一个巴特沃斯数字带通滤波器，要求性能指标为：通带内的衰减≥-3dB,0.45π≤w≤0.55π,阻带内的衰减≤-15dB,0≤w≤0.3π,0.7π≤w≤π,求此滤波器的系统函数和频率特性。wp1=0.45*pi;wp2=0.55*pi;Rp=3;T=1ws1=0.3*pi;ws2=0.7*pi;Rs=15;%数字滤波器指标Wp1=2/T*tan(wp1/2);Wp2=2/T*tan(wp2/2);Ws1=2/T*tan(ws1/2);Ws2=2/T*tan(ws2/2);%模拟带通滤波器指标Wpp1=(Wp1*Wp1-Wp1*Wp2)/Wp1/(Wp2-Wp1)Wpp2=(Wp2*Wp2-Wp1*Wp2)/Wp2/(Wp2-Wp1)Wp=max(abs(Wpp1),abs(Wpp2))Wss1=(Ws1*Ws1-Wp1*Wp2)/Ws1/(Wp2-Wp1)Wss2=(Ws2*Ws2-Wp1*Wp2)/Ws2/(Wp2-Wp1)Ws=min(abs(Wss1),abs(Wss2))%模拟低通滤波器的截止频率[N,Wc]=buttord(Wp,Ws,Rp,Rs,'s');%确定模拟滤波器的最小阶数和截止频率[z,p,k]=buttap(N);%模拟原型低通滤波器零极点增益型[ba,aa]=zp2tf(z,p,k)Wo=sqrt(Wp1*Wp2);B=Wp2-Wp1;[b,a]=lp2bp(ba,aa,Wo,B)%模拟原型低通转化为带通滤波器[bz,az]=bilinear(b,a,1/T)%采用双线性变换转化为数字滤波器freqz(bz,az,512)00.10.20.30.40.50.60.70.80.91-200-1000100200NormalizedFrequency(rad/sample)Phase(degrees)00.10.20.30.40.50.60.70.80.91-150-100-500NormalizedFrequency(rad/sample)Magnitude(dB)直接设计数字滤波器的MATLAB函数：[N,wn]=buttord(wp,ws,Rp,Rs)%数字频率采用标准化频率，取值范围为0~1之间，标准化频率1对应的数字频率为π，对应的模拟频率为采样频率的一半。设计带通滤波器时，wp=[wp1,wp2];ws=[ws1,ws2][b,a]=butter(N,wn)[b,a]=butter(N,wn,’ftype’)%N为滤波器的阶数，wn为滤波器的截止频率（0~1），“ftype”为滤波器的类型：‘high’为高通，‘stop’为带阻，截止频率为wn=[w1,w2];带通滤波时，wn=[w1,w2];缺省时为低通和带通滤波器1．wp=0.2*pi;%digitalPassbandfreqinradws=0.3*pi;%digitalStopbandfreqinradFs=1000;wp1=2*Fs*tan(wp/2);ws1=2*Fs*tan(ws/2);Rp=1;Rs=15;[N,Wn]=buttord(wp1,ws1,Rp,Rs,'s');%返回模拟滤波器的最小阶数和截止频率[Z,P,K]=buttap(N);%模拟低通滤波器原型[Bap,Aap]=zp2tf(Z,P,K);[b,a]=lp2lp(Bap,Aap,Wn);[bz,az]=bilinear(b,a,Fs);[H,W]=freqz(bz,az);plot(W,abs(H));gridxlabel('频率/弧度')ylabel('频率响应幅度')[H,W]=freqz(bz,az,50)figurestem(W,abs(H))第一题直接设计数字滤波器程序•wp=0.2;%归一化通带截止频率•ws=0.3;%归一化阻带截止频率•Rp=1;Rs=15;•[N,Wn]=buttord(wp,ws,Rp,Rs);%返回模拟滤波器的最小阶数和截止频率•[bz,az]=butter(N,Wn);•[H,W]=freqz(bz,az);•plot(W,abs(H))•grid•xlabel('频率/弧度')•ylabel('频率响应幅度')•[H,W]=freqz(bz,az,50)•figure•stem(W,abs(H))例：设计一个butterworth高通数字滤波器，满足：通带边界频率为400Hz，阻带边界频率为200Hz；通带纹波小于3dB；阻带衰减大于15dB；采样频率为1000Hz。MATLAB程序如下：fs=1000;wp=400*2/fs;ws=200*2/fs;rp=3;rs=15;Nn=128;[N,wn]=buttord(wp,ws,rp,rs);[ba,az]=butter(N,wn,'high')freqz(b,a,Nn,fs)050100150200250300350400450500050100150200Frequency(Hz)Phase(degrees)050100150200250300350400450500-100-80-60-40-200Frequency(Hz)Magnitude(dB)•fs=1000;•wp=[100250]*2/fs;•ws=[50300]*2/fs;•rp=3;•rs=30;•Nn=128;•[N,wn]=cheb2ord(wp,ws,rp,rs)•[b,a]=cheby2(N,rp,wn,'bandpass')•freqz(b,a,Nn,fs)•运行结果4、带通滤波器程序