河海大学-通信工程-DSP实验4

1、实验四无限长单位脉冲响应滤波器设计班级：姓名：学号：一、实验目的1、掌握双线形变换法及脉冲响应不变法设计IIR 数字滤波器的具体设计方法及其原理，熟悉用双线形变换法及脉冲响应不变法设计低通，高通和带通IIR数字滤波器的计算机编程。2、观察双线形变换法及脉冲响应不变法设计的滤波器的频域特性，了解双线形变换法及脉冲响应不变法的特点。3、熟悉巴特沃思滤波器，切比雪夫滤波器和椭圆滤波器的频率特性。；二、实验原理1、将模拟滤波器的参数指标通过频率变换转换为原型模拟低通滤波器的参数指标，设计满足指标要求的原型模拟低通滤波器，通过频率变换将原型模拟低通滤波器的系统函数变换为其他类型的模拟滤波器2 、脉冲响应

2、不变法是使数字滤波器的单位脉冲响应序列 h（ n）逼近模拟滤波器的冲激响应 hs（ t ），让 h（n）正好等于 hs（t ）的采样值。3 、双线性变换法使得 s平面与 z平面是一一映射关系，消除了多值变换性，克服了脉冲响应不变法产生的频响的混叠失真。4 、巴特沃思滤波器通带内有最大平坦的幅度特性且随着频率的升高而单调的下降和切比雪夫滤波器逼近误差峰值在一个规定的频段上位最小。三、实验内容1 、fc=0.3kHz, =0.8dB,fr=0.2kHz,At=20dB,T=1ms; 设计一切比雪夫高通滤波器， 观察其通带损耗和阻带衰减是否满足要求。2 、 fc=0.2kHz, =1dB,fr=0.

3、3kHz,At=25dB,T=1ms; 分别用脉冲响应不变法及双线形变换法设计一巴特沃思数字低通滤波器， 观察所设计数字滤波器的幅频特性曲线， 记录带宽和衰减量，检查是否满足要求。比较这两种方法的优缺点。3 、利用双线性变换法分别设计满足下列指标的巴特沃思滤波器，切比雪夫滤波器和椭圆滤波器，并作图验证设计结果：fc=1.2kHz, <=0.5dB,fr=2kHz,At>=40dB,fs=8kHz。4 、分别用脉冲响应不变法及双线形变换法设计一巴特沃思数字带通滤波器，已知fs=3kHz ， 其 等 效 的 模 拟 滤 波 器 指 标 为 <3kHz, 2kHz<f<

4、=3kHz, At>=5dB, f>=6kHz, At>=20dB,f<=1.5kHz 。5 、利用双线性变换法设计满足下列指标的切比雪夫型数字带阻滤波器，并作图验证设计结果：当 1kHz<=f<=2kHz 时， At>=18dB；当 f<=500Hz 以及 f>=3kHz 时， <=3dB ；采样频率 fs=10kHz 。四、实验结果与分析1、解：实验程序：wc=2*1000*tan(2*pi*300/(2*1000);wr=2*1000*tan(2*pi*200/(2*1000);N,wn=cheb1ord(wc,wr,0.8,2

5、0,'s');B,A=cheby1(N,0.8,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 on;xlabel('频率 /Hz');ylabel('幅度 /dB');实验图形：100-10-20B-30d/度幅 -40 -50 -60 -70 -80050100150200250300350400450500频 率 /Hz理

6、论分析：f=200Hz 时阻带衰减大于30dB，通过修改axis(0,fs/2,-80,10) 为 axis(200,fs/2,-1,1)，发现通带波动 rs 满足 <0.8。10.80.60.4B0.2d/0度幅-0.2-0.4-0.6-0.8-1250300350400450500200频 率 /Hzbz =0.0262-0.10470.1570-0.10470.0262az =1.00001.52891.65370.94520.2796系统函数为：0.0262 - 0.1047 z 10.1570 z 20.1047 z 30.0262 z 4H ( z)1.5289 z 11.6

7、537 z 20.9452 z 30.2796z 412、解：实验程序：fs=1000;fc=200;fr=300;rp=1;rs=25;wp=2*pi*fc;ws=2*pi*fr;N, wn = buttord(wp, ws, rp, rs, 's');b1 a1=butter(N,wn,'s');bz1,az1=impinvar(b1,a1,fs);h1,w=freqz(bz1,az1);wp=2*fs*tan(2*pi*fc/fs/2);ws=2*fs*tan(2*pi*fr/fs/2);N, wn = buttord(wp, ws, rp, rs, &#

8、39;s');b2 a2=butter(N,wn,'s');bz2,az2=bilinear(b2,a2,fs);h2,w=freqz(bz2,az2);f=w/(2*pi)*fs;figure; plot(f,abs(h1),'-.r',f,abs(h2),'-b');grid; xlabel('频率 /Hz'); ylabel('幅度 ');legend('脉冲响应不变法','双线性变换法');title('巴特沃思低通滤波器幅频特性');实验图形：巴特沃

9、思低通滤波器幅频特性1.4脉冲响应不变法双线性变换法1.210.8度幅0.60.40.20050100150200250300350400450500频 率 /Hz理论分析：bz1 =0.00000.00020.01530.09950.14440.06110.00750.00020.00000az1 =1.0000 -1.9199 2.5324 -2.2053 1.3868 -0.6309 0.2045 -0.0450 0.0060 -0.0004脉冲响应不变法的系统函数为:H imp( z)0.0002 z 10.0153 z 20.0995 z 30.1444 z 40.0611z 50.

10、0075 z 60.0002 z 71 2.5324 z 22.2053 z 31.3868 z 40.6309 z 50.2045 z 60.0450 z 70.0060 z 8 0.0004 z 91 1.9199 zbz2 =0.01790.10720.26810.35750.26810.10720.0179az2 =1.0000-0.60190.9130-0.29890.1501-0.02080.0025双线性变换法的系统函数为:Hbil (z)0.01790.1072z 10.2681z 20.3575z 30.2681z 40.1072z 50.0179z 61 0.6019z 1

11、 0.9130z 20.2989z 30.1501z 40.0208z 50.0025z 6脉冲响应不变法的N=9，双线性变换法的N=6，由图知它们都满足要求，但脉冲响应的衰减较快，双线性变换的过渡带窄一些，且阶数比脉冲小，容易实现。3、解：实验程序：fs=8000;fc=1200;fr=2000;rp=0.5;rs=40;wp=2*fs*tan(2*pi*fc/fs/2);ws=2*fs*tan(2*pi*fr/fs/2);N, wn = buttord(wp, ws, rp, rs, 's');b1a1=butter(N,wn,'s');bz1,az1=bi

12、linear(b1,a1,fs);h1,w=freqz(bz1,az1);H1=20*log10(abs(h1);f=w/(2*pi)*fs;figure; plot(f,H1);axis(0,fs/2,-100,10);grid; xlabel('频率 /Hz'); ylabel('幅度 ');title('巴特沃思低通滤波器对数幅度谱');wc=2*fs*tan(2*pi*fc/(2*fs);wt=2*fs*tan(2*pi*fr/(2*fs);N,wn=cheb1ord(wc,wt,rp,rs,'s');b2,a2=cheb

13、y1(N,rp,wn,'low','s');bz2,az2=bilinear(b2,a2,fs);h2,w=freqz(bz2,az2);H2=20*log10(abs(h2);f=w*fs/(2*pi);figure;plot(f,H2);axis(0,fs/2,-100,10);grid;xlabel('频率 /Hz');ylabel('幅度 /dB');title('切比雪夫低通滤波器对数幅度谱');wp=2*fs*tan(2*pi*fc/fs/2);N,wp=ellipord(wp,ws,rp,rs,

14、9;s');b3,a3=ellip(N,rp,rs,wp,'low','s');bz3,az3=bilinear(b3,a3,fs);h3,w=freqz(bz3,az3);H3=20*log10(abs(h3);f=w/(2*pi)*fs;figure;plot(f,H3);axis(0,fs/2,-100,10);grid; xlabel('频率 /Hz'); ylabel('幅度 /dB');title('椭圆型数字低通滤波器对数幅度谱');实验图形（ 1）：巴特沃思低通滤波器对数幅度谱100-10-

15、20-30-40度幅 -50 -60-70 -80 -90-10005001000150020002500300035004000频 率 /Hz理论分析：bz1 =0.00040.00320.01290.0302 0.04530.04530.03020.01290.00320.0004az1 =1.0000-2.79964.4582-4.54123.2404-1.63300.5780 -0.13700.0197-0.0013系统函数为 :H10.00040.0032z 10.0129z 20.0302z 30.0453z 40.0453z 50.0302z 60.0129z 70.0032z

16、80.0004z 9(z)2.7996z 14.4582z 24.5412z 33.2404z 41.6330z 50.5680z 60.1370z 70.0197z 80.0013z 91N=9，为九阶巴特沃思低通滤波器，从图中可以看出通带波动和阻带衰减都满足设计要求。实验图形（ 2）：切比雪夫低通滤波器对数幅度谱100-10-20-30Bd-40/度幅 -50 -60 -70 -80-90-10005001000150020002500300035004000频 率 /Hz理论分析：bz2=0.00260.01320.02640.02640.01320.0026az2= 1.0000-2.

17、97754.2932-3.51241.6145-0.3334系统函数为：0. 0026 0. 0132z 10. 0264z 20. 0264z 30. 0132z 40. 0026z 5H2( z )14. 2932z 23. 5124z 31. 6154z 40. 3334z 51 2. 9775zN=5，为五阶切比雪夫低通滤波器，从图中可以看出通带波动和阻带衰减都满足设计要求。实验图形（ 3）：椭圆型数字低通滤波器对数幅度谱100-10-20-30Bd-40/度幅 -50 -60 -70 -80 -90-10005001000150020002500300035004000频 率 /Hz

18、理论分析：bz =0.03890.03630.06650.03630.0389az =1.0000-2.14442.3658-1.32500.3332系统函数为：0. 0389 - 0. 036310. 066520. 0363z30. 0389z4H3( z )zz12. 3658z 2- 1.3250 z 30. 3332z 41 - 2. 1444zN=4，为四阶椭圆型数字低通滤波器，从图中可以看出通带波动和阻带衰减都满足设计要求。4、解：实验程序：fs=30000;fc=2000,3000;fr=1500,6000;rp=3;rs=20;wp=2*pi*fc;ws=2*pi*fr;N,

19、wn=buttord(wp, ws, rp, rs, 's');b1,a1=butter(N,wn,'s');bz1,az1=impinvar(b1,a1,fs);h1,w=freqz(bz1,az1);wp=2*fs*tan(2*pi*fc/fs/2);ws=2*fs*tan(2*pi*fr/fs/2);N,wn=buttord(wp, ws, rp, rs, 's');b2,a2=butter(N,wn,'s');bz2,az2=bilinear(b2,a2,fs);h2,w=freqz(bz2,az2);f=w/(2*pi)

20、*fs;figure;plot(f,abs(h1),'-.r',f,abs(h2),'-b');grid;xlabel(' 频率 /Hz');ylabel(' 幅度 ');legend(' 脉冲响应不变法',' 双线性变换法 ');实验图形：10.9脉冲响应不变法双线性变换法0.80.70.6度 0.5幅0.40.30.20.10500010000150000频 率 /Hz理论分析：bz1 =-0.00000.0057-0.01220.00250.0089-0.00490az1 =1.0000-4.

21、805610.2376-12.26258.7012-3.47190.6145脉冲响应不变法的系统函数为:0.0057 z 10 .0122 z 20 .0 025 z30 .0089 z40 .0049 z5H ( z) 1110 .2376 z 212 .2625 z 38 . 7012 z 43 .4719 z 50. 6145 z 61 4. 8056 zbz2 =0.0014-0.0000-0.0042-0.00000.0042-0.0000-0.0014az2 =1.0000-4.807110.2473-12.28388.7245-3.48490.6176双线性变换法的系统函数为:H ( z) 20 .00140 .0042 z 20. 0042 z 40 .0014 z 61234561 4 .8071 z10 .2473 z12 .2838 z8.7245 z3 .4849 z0 .6176 z5、解：实验程序：f3=1000;f4=2000;ar=18;tao=3;f1=500;f2=3000;fs=10000;fp1=pi*f1/fs;fp2=pi*f2/fs;fc1=pi*f3/fs;fc2=pi*f4/fs;w1=2*fs*tan(fp1) 2*fs*tan(