数字信号处理实验四

西南大学工程技术学院学生实验报告姓名：钟昌余学号：222013322270082班级：2专业：自动化实验日期：2015.11.01实验学时：3同组人：无00.10.20.30.40.50.60.70.80.91-35-30-25-20-15-10-5000.10.20.30.40.50.60.70.80.91-120-100-80-60-40-20020实验四IIR和FIR数字滤波器设计一、实验目的1、掌握IIR数字滤波器的设计原理、设计方法和设计步骤；2、掌握IIR数字滤波器的设计原理、设计方法和设计步骤；3、学会根据给定的滤波器指标进行滤波器设计；4、掌握数字巴特沃斯滤波器、数字切比雪夫滤波器的设计原理和步骤；5、掌握窗函数法、频域抽样法和最优化方法设计FIR数字滤波器。二、实验内容1、%巴特沃斯+冲激响应变化+低通（IIR）DFfs=10000;wp=1000*2*pi;ws=3000*2*pi;rp=1;rs=20;[N,Wn]=buttord(wp,ws,rp,rs,'s');[b,a]=butter(N,Wn,'s');[bz,az]=impinvar(b,a,fs);[h,w]=freqz(bz,az);plot(w/pi,20*log10(abs(h)));grid2、%巴特沃斯+冲激响应变化+带通（IIR）DFfs=2000;rp=2;rs=40;wp1=2*pi*400;wp2=2*pi*600;ws1=2*pi*200;ws2=2*pi*800;wp=[wp1,wp2];ws=[ws1,ws2];[N,Wn]=buttord(wp,ws,rp,rs,'s');[b,a]=butter(N,Wn,'s');[bz,az]=impinvar(b,a,fs);[h,w]=freqz(bz,az);plot(w/pi,20*log10(abs(h)));grid3、%巴特沃斯+双线性映射+带阻（IIR）DFfs=2000;rp=2;rs=40;wp1=2*pi*200;wp2=2*pi*800;ws1=2*pi*400;ws2=2*pi*600;西南大学工程技术学院学生实验报告姓名：钟昌余学号：222013322270082班级：2专业：自动化实验日期：2015.11.01实验学时：3同组人：无00.10.20.30.40.50.60.70.80.91-350-300-250-200-150-100-50000.10.20.30.40.50.60.70.80.91-350-300-250-200-150-100-50000.10.20.30.40.50.60.70.80.91-300-250-200-150-100-500Wp1=2*fs*tan(wp1/2/fs);Wp2=2*fs*tan(wp2/2/fs);Ws1=2*fs*tan(ws1/2/fs);Ws2=2*fs*tan(ws2/2/fs);wp=[Wp1,Wp2];ws=[Ws1,Ws2];[N,Wn]=buttord(wp,ws,rp,rs,'s');[b,a]=butter(N,Wn,'stop','s');[bz,az]=bilinear(b,a,fs);[h,w]=freqz(bz,az);plot(w/pi,20*log10(abs(h)));grid4、%切比雪夫+双线性++带阻（IIR）DFfs=2000;rp=2;rs=40;wp1=2*pi*200;wp2=2*pi*800;ws1=2*pi*400;ws2=2*pi*600;Wp1=2*fs*tan(wp1/2/fs);Wp2=2*fs*tan(wp2/2/fs);Ws1=2*fs*tan(ws1/2/fs);Ws2=2*fs*tan(ws2/2/fs);wp=[Wp1,Wp2];ws=[Ws1,Ws2];[N,Wn]=cheb1ord(wp,ws,rp,rs,'s');[b,a]=cheby1(N,rp,Wn,'stop','s');[bz,az]=bilinear(b,a,fs);[h,w]=freqz(bz,az);plot(w/pi,20*log10(abs(h)));grid4、%直接+带阻（IIR）DFfs=2000;rp=2;rs=40;ws1=2*400/fs;ws2=2*600/fs;wp1=2*200/fs;wp2=2*800/fs;wp=[wp1,wp2];ws=[ws1,ws2];[N,Wn]=buttord(wp,ws,rp,rs);[b,a]=butter(N,Wn,'stop');[h,w]=freqz(b,a);plot(w/pi,20*log10(abs(h)));grid西南大学工程技术学院学生实验报告姓名：钟昌余学号：222013322270082班级：2专业：自动化实验日期：2015.11.01实验学时：3同组人：无0500100015002000250030003500-300-250-200-150-100-50000.10.20.30.40.50.60.70.80.91-1500-1000-5000NormalizedFrequency(rad/sample)Phase(degrees)00.10.20.30.40.50.60.70.80.91-100-50050NormalizedFrequency(rad/sample)Magnitude(dB)00.10.20.30.40.50.60.70.80.91-3000-2000-100001000NormalizedFrequency(rad/sample)Phase(degrees)00.10.20.30.40.50.60.70.80.91-100-50050NormalizedFrequency(rad/sample)Magnitude(dB)5、%直接+带阻（IIR）模拟滤波器fs=2000;rp=2;rs=40;ws1=2*400*pi;ws2=2*600*pi;wp1=2*200*pi;wp2=2*800*pi;wp=[wp1,wp2];ws=[ws1,ws2];[N,Wn]=buttord(wp,ws,rp,rs,'s');[b,a]=butter(N,Wn,'stop','s');[h,w]=freqs(b,a);plot(w/pi,20*log10(abs(h)));grid6、窗函数设计法+低通+FIR+DF已知：wp=0.2*pi，ws=0.4*pi，最小衰减为as=50dB解：wc=0.3*pi，hd=sin(wc*(n-a))/pi./(n-a);a=(N-1)/2;Δw=6.6*pi/N,得N=33（先由衰减确定海明窗，再确定点数）clearall;closeall;clcN=35;%33当不满足性能要求是增加点数wc=0.3*pi;a=(N-1)/2;n=0:N-1;hd=sin(wc*(n-a))/pi./(n-a);hd(a+1)=wc/pi;wn=hamming(N);hn=hd.*wn';freqz(hn,1)总结：根据最小衰减As（dB）确定窗函数类型；根据过渡带宽度确定窗函数的点数。wc=1/2（wp+ws）As=-20log[δ/(1+δ)]。高通、带通、带阻都可由低通滤波器相加减得到。7、窗函数设计法+高通+FIR+DF（利用全通减去低通得到）N=35;wc=0.3*pi;a=(N-1)/2;n=0:N-1;hd=sin(wc*(n-a))/pi./(n-a);hd(a+1)=wc/pi;wn=hamming(N);hn=hd.*wn';西南大学工程技术学院学生实验报告姓名：钟昌余学号：222013322270082班级：2专业：自动化实验日期：2015.11.01实验学时：3同组人：无00.10.20.30.40.50.60.70.80.91-3000-2000-10000NormalizedFrequency(rad/sample)Phase(degrees)00.10.20.30.40.50.60.70.80.91-100-50050NormalizedFrequency(rad/sample)Magnitude(dB)00.10.20.30.40.50.60.70.80.91-1500-1000-5000NormalizedFrequency(rad/sample)Phase(degrees)00.10.20.30.40.50.60.70.80.91-150-100-50050NormalizedFrequency(rad/sample)Magnitude(dB)wc=0.6*pi;%wc=pi（高通）a=(N-1)/2;n=0:N-1;hd=sin(wc*(n-a))/pi./(n-a);hd(a+1)=wc/pi;wn=hamming(N);hn1=hd.*wn';freqz(hn1-hn,1)8、窗函数设计法+带阻+FIR+DFN=35;wc=0.3*pi;a=(N-1)/2;n=0:N-1;hd=sin(wc*(n-a))/pi./(n-a);hd(a+1)=wc/pi;wn=hamming(N);hn1=hd.*wn';%低通1wc=0.6*pi;hd=sin(wc*(n-a))/pi./(n-a);hd(a+1)=wc/pi;wn=hamming(N);hn1=hd.*wn';%低通2wc=pi;hd=sin(wc*(n-a))/pi./(n-a);hd(a+1)=wc/pi;wn=hamming(N);hn2=hd.*wn';freqz(hn2-hn1+hn,1)%全通-低通2+低通1=带阻9、频率抽样法设计法+低通+FIR+DF已知：Wp=0.2*pi,Ws=0.4*pi;δ1=1dB,δ2=50dB求解：由δ2=50dB知过渡带点数m=1，Δw=（m+1）*2*pi/N,得N=20，截止频率点：N1=floor(wc*N/（2*pi）);N=31;wc=0.2*pi;k=0:N-1;theta=-pi*k*(1-1/N);N1=floor(wc*N/2/pi);hk=[ones(1,N1+1),zeros(1,N-2*N1-1),ones(1,N1)];hk(N1+2)=0.588;hk(N1+3)=0.106;hk(N-N1)=0.588;hk(N-N1-1)=0.106;西南大学工程技术学院学生实验报告姓名：钟昌余学号：222013322270082班级：2专业：自动化实验日期：2015.11.01实验学时：3同组人：无00.10.20.30.40.50.60.70.80.91-3000-2000-100001000NormalizedFrequency(rad/sample)Phase(degrees)00.10.20.30.40.50.60.70.80.91-100-500NormalizedFrequency(rad/sample)Magnitude(dB)00.10.20.30.40.50.60.70.80.91-80-70-60-50-40-30-20-10010Hk=hk.*exp(1j*theta);hn=real(ifft(Hk));freqz(hn,1)10、频率抽样法设计法+高通+FIR+DFN=45;wc=0.44*pi;k=0:N-1;d1=0.334;theta=-pi*k*(1-1/N);N1=floor(wc*N/2/pi);hk=[zeros(1,N1-1),d1,ones(1,N-2*N1+1),d1,zeros(1,N1-2)];Hk=hk.*exp(1j*theta);hn=real(ifft(Hk));freqz(hn,1);axis([0,1,-100,10])11、最优化设计法+带阻+FIR+DFclc;clearall;rp=3;rs=40;fs=2000;f=[200,300,500,600];a=[1,0,1