实验三参考IIR数字滤波器的设计

实验三IIR数字滤波器的设计一、实验目的(1)掌握双线性变换法及脉冲相应不变法设计IIR数字滤波器的具体设计方法;(2)熟悉用双线性变换法及脉冲响应不变法设计低通、高通和带通IIR数字滤波器的计算机编程。二、实验原理在MATLAB中，可以用下列函数辅助设计IIR数字滤波器：1)利用buttord和cheb1ord可以确定低通原型巴特沃斯和切比雪夫滤波器的阶数和截止频率；2)[num,den]=butter（N,Wn）(巴特沃斯)和[num,den]=cheby1（N,Wn）,[num,den]=cheby2（N,Wn）(切比雪夫1型和2型)可以进行滤波器的设计；3）lp2hp，lp2bp，lp2bs可以完成低通滤波器到高通、带通、带阻滤波器的转换；4)使用bilinear可以对模拟滤波器进行双线性变换，求得数字滤波器的传输函数系数；5)利用impinvar可以完成脉冲响应不变法的模拟滤波器到数字滤波器的转换。三、要求（1）在MATLAB中，熟悉函数butter、cheby1、cheby2的使用，其中：[num,den]=butter（N,Wn）巴特沃斯滤波器设计；[num,den]=cheby1（N,Wn）切比雪夫1型滤波器设计；[num,den]=cheby2（N,Wn）切比雪夫2型滤波器设计。（2）阅读附录中的实例，学习在MATLAB中进行数字滤波器的设计；（3）给出IIR数字滤波器参数和滤波器的冲激响应，绘出它们的幅度和相位频响曲线，讨论它们各自的实现形式和特点。四、实验内容利用MATLAB编程，用脉冲响应不变法和双线性变换法设计一个数字带通滤波器，指标要求如下：通带边缘频率：45.01P，65.02P，通带峰值起伏：][1dBp。阻带边缘频率：3.01S，8.02S，最小阻带衰减：][40dBS。附录：例1设采样周期T=250μs（采样频率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,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')程序中第一个butter的边界频率2π×1000，为脉冲响应不变法原型低通滤波器的边界频率；第二个butter的边界频率2/T=2/0.00025，为双线性变换法原型低通滤波器的边界频率.图1给出了这两种设计方法所得到的频响，虚线为脉冲响应不变法的结果；实线为双线性变换法的结果。脉冲响应不变法由于混叠效应，使得过渡带和阻带的衰减特性变差，并且不存在传输零点。同时，也看到双线性变换法，在z=-1即Ω=π或f=2000Hz处有一个三阶传输零点，这个三阶零点正是模拟滤波器在ω=∞处的三阶传输零点通过映射形成的。例2设计一数字高通滤波器，它的通带为400～500Hz，通带内容许有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,'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=120kHz处的最小衰减大于１０dB，采样频率fs=400kHz。w1=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([w1w2],[0wr],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')例4一数字滤波器采样频率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')参考源程序：closeall,clearall,clc;mp1=2*tan(0.2*pi/2);mp2=2*tan(0.3*pi/2);ms1=2*tan(0.1*pi/2);ms2=2*tan(0.4*pi/2);%预畸变r=mp2-mp1;%参考频率rp1=mp1/r;rp2=mp2/r;rs1=ms1/r;rs2=ms2/r;%归一化模拟带通滤波器的边界频率r0=sqrt(rp1*rp2);mpL1=(rp1^2-r0^2)/rp1;%归一化模拟低通原型滤波器mpL2=(rp2^2-r0^2)/rp2;msL1=(rs1^2-r0^2)/rs1;msL2=(rs2^2-r0^2)/rs2;mp=mpL2;ms=min(abs(msL1),abs(msL2));ap=1;as=25;%转换成模拟低通滤波器后的参数[N,mc]=buttord(mp,ms,ap,as,'s')%确定模拟低通原型的阶数和截止频率[zl,pl,kl]=buttap(N);[num,den]=zp2tf(zl,pl,kl);zplane(zl,pl);%模拟低通滤波器的零极点图figure;wo=sqrt(mp1*mp2);%模拟带通的中心频率[bt,at]=lp2bp(num,den,wo,mp2-mp1+0.1);%将模拟低通滤波器转换为模拟带通滤波器[zb,pb,kb]=tf2zp(bt,at);zplane(zb,pb);%模拟带通的零极点图figure;fs=1;%双线性变换法求得数字带通[zd,pd,kd]=bilinear(zb,pb,kb,fs);[b,a]=zp2tf(zd,pd,kd);zplane(zd,pd);%数字带通的零极点图figure;freqz(b,a);