有限脈沖響應數(shù)字濾波器設計實驗報告

1、成績：數(shù)字信號處理作業(yè)與上機實驗（第二章）班 級：學 號：姓 名：任課老師：完成時間：信息與通信工程學院20142015學年第1學期第7章 有限脈沖響應數(shù)字濾波器設計1、教材 p238:19 .設信號x(t) = s(t) + v(t),其中v(t)是干擾，s與v(t)的頻譜不混 疊，其幅度譜如題19圖所示。要求設計數(shù)字濾波器，將干擾濾除， 指標是允許|s|在0Wfwi5 kHz頻率范圍中幅度失真為± 2%( 5 1 = 0.02) ;f > 20 kHz,衰減大于40 dB ( 8 2=0.01);希望分別設計 性價比最高的FIR和IIR兩種濾波器進行濾除干擾。 請選擇合適

2、的濾波器類型和設計方法進行設計，最后比較兩種濾波器的幅頻特性、相頻4,叔郴數(shù)。_ y/kHz題19圖(1)matlab 代碼：%基于雙線性變換法直接設計IIR數(shù)字濾波器Fs=80000;fp=15000;fs=20000;rs=40;wp=2*pi*fp/Fs;ws=2*pi*fs/Fs;Rp=-20*log10(1-0.02);As=40;N1,wp1=ellipord(wp/pi,ws/pi,Rp,As);B,A=ellip(N1,Rp,As,wp1);Hk,wk1=freqz(B,A,1000);mag=abs(Hk);pah=angle(Hk);%窗函數(shù)法設計FIR數(shù)字濾波器Bt=ws

3、-wp;alph=0.5842*(rs-21)A0.4+0.07886*(rs-21);N=ceil(rs-8)/2.285/Bt);wc=(wp+ws)/2/pi;hn=fir1(N,wc,kaiser(N+1,alph);M=1024;Hk=fft(hn,M);k=0:M/2-1;wk=(2*pi/M)*k;%畫出各種比較結果圖 figure(2);plot(wk/pi,20*log10(abs(Hk(k+1), hold on,'linewidth',2.5)；plot(wk1/pi,20*log10(mag), hold off'linewidth',2

4、);legend( 'FIR 濾波器，IIR axis(0,1,-80,5);xlabel( title('損耗函數(shù)');figure(3) plot(wk/pi,angle(Hk(k+1)/pi, hold on濾波器');'w/pi');ylabel( '幅度/dB','linewidth',2.5););圖1損耗函數(shù)比較圖圖2相頻特性比較圖plot(wk1/pi,pah/pi, hold off'linewidth' ,2);legend( 'FIR 濾波器，IIR 濾波器);xlab

5、el( 'w/pi' );ylabel( '相位八pi');title('相頻特性曲線);1、2(2)兩種數(shù)字濾波器的損耗函數(shù)和相頻特性的比較分別如圖所示:損耗函數(shù)0-10-20FIR濾波器一 IIR濾波器4610.80.60.4O 1,20 2 4 6s-1 o-0-0-0-0/.位相FIR濾波器IIR濾波器(3)IIR數(shù)字濾波器階數(shù)：N=5FIR數(shù)字濾波器階數(shù)：N=36(4)運行結果分析：由圖2及階數(shù)可見，IIR階數(shù)低得多，但相 位特性存在非線性失真，F(xiàn)IR具有線性相位特性。20 .調用MATLAB工具箱函數(shù)萬門設計線性相位低通FIR濾波器， 要求希

6、望逼近的理想低通濾波器通帶截止頻率3 c=兀/4 rad,濾波器長度N=21。分別選用矩形窗、Hanning窗、Hamming窗和Blackman窗進行設計，繪制用每種窗函數(shù)設計的單位脈沖響應h(n)及其損耗函數(shù)曲線，并進行比較,觀察各種窗函數(shù)的設計性臺匕 目匕。(1)matlab 代碼:wc=pi/4;N=21;hn_boxcar=fir1(N-1,wc/pi,boxcar(N);hn_hanning=fir1(N-1,wc/pi,hanning(N);hn_hamming=fir1(N-1,wc/pi,hamming(N);hn_blackman=fir1(N-1,wc/pi,blackm

7、an(N);n=0:N-1; plot(n,hn_boxcar);hold on plot(n,hn_hanning, plot(n,hn_hamming, plot(n,hn_blackman, hold offxlabel( 'n' );ylabel(':','linewidth''+' , 'linewidth' 'o');'h(n)');,2);,2);legend('矩形窗，漢寧窗，哈明窗，布萊克曼窗);title( '單位沖激響應');M=1024

8、;Hk=fft(hn_boxcar,M);k=0:M/2-1;wk=(2*pi/M)*k;'linewidth',2);figure();plot(wk/pi,20*log10(abs(Hk(k+1), Hk=fft(hn_hanning,M)hold onplot(wk/pi,20*log10(abs(Hk(k+1), Hk=fft(hn_hamming,M)plot(wk/pi,20*log10(abs(Hk(k+1), Hk=fft(hn_blackman,M)plot(wk/pi,20*log10(abs(Hk(k+1), hold off,'linewidth

9、''o'););,3)；圖4損耗函數(shù)比較圖legend('矩形窗，漢寧窗，哈明窗，布萊克曼窗)；axis(0,1,-80,5);xlabel('w/pi' );ylabel( '幅度/dB');title('損耗函數(shù)');(2)四種窗函數(shù)設計的單位脈沖響應的比較如圖3所示:單位沖激響應0.3矩形窗-0.05-0.1101214161820圖3 單位脈沖響應比較圖四種窗函數(shù)設計的損耗函數(shù)的比較如圖4所示:損耗函數(shù)0-10-20-30-40-50-60-70-8000.10.20.30.40.50.60.70.80.9

10、1w/ 7：.運行結果分析：由圖4可見，當濾波器長度N不變時，矩形窗設計的濾波器的過渡帶最窄，阻帶最小衰減最??；布萊克曼窗設計的濾波器的過渡帶最寬，同時阻帶最小衰減最大。21 .將要求改成設計線性相位高通 FIR濾波器，重作題20。(1)matlab 代碼:wc=pi/4;N=21;':','linewidth''+' , 'linewidth' 'o');,2);,2);hn_boxcar=fir1(N-1,wc/pi, hn_hanning=fir1(N-1,wc/pi, hn_hamming=fir1(N-1

11、,wc/pi, hn_blackman=fir1(N-1,wc/pi, n=0:N-1; plot(n,hn_boxcar); hold on plot(n,hn_hanning, plot(n,hn_hamming, plot(n,hn_blackman, hold offhigh' ,boxcar(N);'high' ,hanning(N);'high' ,hamming(N);'high' ,blackman(N)xlabel( 'n' );ylabel( 'h(n)');legend（'矩形

12、窗，漢寧窗，哈明窗，布萊克曼窗）;title（ '單位沖激響應'）;M=1024;Hk=fft（hn_boxcar,M）;k=0:M/2-1;wk=(2*pi/M)*k;'linewidth' ,2);':','linewidth' ,3);'o');figure();plot(wk/pi,20*log10(abs(Hk(k+1), Hk=fft(hn_hanning,M)hold onplot(wk/pi,20*log10(abs(Hk(k+1), Hk=fft(hn_hamming,M)plot(wk/pi,2

13、0*log10(abs(Hk(k+1), Hk=fft(hn_blackman,M)plot(wk/pi,20*log10(abs(Hk(k+1), hold off legend（'矩形窗，漢寧窗，哈明窗，布萊克曼窗）；axis（0,1,-80,5）;xlabel（'w八pi'）;ylabel（ '幅度/dB'）;title（'損耗函數(shù)'）;(2)四種窗函數(shù)設計的單位脈沖響應的比較如圖5所示:單位沖激響應0.80.60.40.20-0.2*矩形窗漢寧窗*哈明窗16182002468101214n圖5 單位脈沖響應比較圖四種窗函數(shù)設計的損

14、耗函數(shù)的比較如圖 6所示:圖6損耗函數(shù)比較圖布萊克曼窗(5)運行結果分析：由圖6可見，當濾波器長度N不變時，矩形窗設計的濾波器的過渡帶最窄，阻帶最小衰減最??；布萊克曼窗 設計的濾波器的過渡帶最寬，同時阻帶最小衰減最大。25 .調用MATLAB工具箱函數(shù)fi門設計線性相位高通FIR濾波器。要求通帶截止頻率為0.6兀rad,阻帶截止頻率為0.45兀，通帶最大衰減為0.2 dB,阻帶最小衰減為45 dBo顯示所設計的單位脈沖響應h(n)的數(shù)據(jù)，弁畫出損耗函數(shù)曲線(1)matlab 代碼：wp=0.6*pi;ws=0.45*pi;Bt=wp-ws;N0=ceil(6.6*pi/Bt);N=N0+mod

15、(N0+1,2);wc=(wp+ws)/2/pi;hn=fir1(N-1,wc, 'high' ,hamming(N);M=1024;Hk=fft(hn,M);n=0:N-1;stem(n,hn);xlabel( 'n' );ylabel( 'h(n)');title( '單位沖激響應');k=0:M/2-1;wk=(2*pi/M)*k;figure(2);plot(wk/pi,20*log10(abs(Hk(k+1););ylabel( '幅度/dB');損耗函數(shù)如圖 7、8所示:axis(0,1,-80,5);

16、xlabel('w/pi'title( '損耗函數(shù)');grid on(2)高通FIR濾波器的單位脈沖響應、0.5單位沖激響應0.40.30.20.1I I-0.1-0.2-0.3-0.4圖8損耗函數(shù)圖7單位脈沖響應26 .調用MATLAB工具箱函數(shù)firl設計線性相位帶通FIR濾波器。要求通帶截止頻率為 0.55兀rad和0.7兀rad,阻帶截止頻率為 0.45兀rad和0.8兀rad,通帶最大衰減為 0.15 dB,阻帶最小衰 減為40 dBo顯示所設計的單位脈沖響應 h(n)的數(shù)據(jù)，弁畫出損 耗函數(shù)曲線。(1)matlab 代碼:wp1=0.55*pi;w

17、p2=0.7*pi;ws1=0.45*pi;ws2=0.8*pi;Bt=wp2-wp1;N=ceil(6.2*pi/Bt);wc=(wp1+ws1)/2/pi,(ws2+wp2)/2/pi;hn=fir1(N-1,wc,hanning(N);M=1024;Hk=fft(hn,M);n=0:N-1;stem(n,hn);xlabel( 'n' );ylabel( 'h(n)');title( '單位沖激響應');k=0:M/2-1;wk=(2*pi/M)*k;figure(2);plot(wk/pi,20*log10(abs(Hk(k+1););y

18、label( '幅度/dB');axis(0,1,-80,5);xlabel('w/pi'title( '損耗函數(shù)');grid on(2)帶通FIR濾波器的單位脈沖響應、損耗函數(shù)如圖9、10所示:1圖10損耗函數(shù)圖9單位脈沖響應2、 某信號 x(t)為：x(t) =0.5cos(2n f1t) +0.7cos(2兀 f2t +0.1兀)+ 0.4cos(2n f3t), 其中fi =100Hz, f2 =130Hz, fi =600Hz.設計最低階FIR數(shù)字濾波器， 按下圖所示對x(t)進行數(shù)字濾波處理，實現(xiàn)：1)將f3頻率分量以高于50dB的

19、衰減抑制，同時以低于 2dB的衰減通過i和f2頻率分量；一、基于窗函數(shù)法設計FIR數(shù)字濾波器：(1)matlab 代碼：Fs=3800;fp=130;fs=600;rs=50;wp=2*pi*fp/Fs;ws=2*pi*fs/Fs;Bt=ws-wp;alph=0.5842*(rs-21)A0.4+0.07886*(rs-21);N=ceil(rs-8)/2.285/Bt);wc=(wp+ws)/2/pi;hn=fir1(N,wc,kaiser(N+1,alph);M=1024;Hk=fft(hn,M);k=0:M/2-1;wk=(2*pi/M)*k;figure(2);plot(wk/pi,2

20、0*log10(abs(Hk(k+1);axis(0,1,-80,5);xlabel('w八pi');ylabel('幅度 /dB');title( '損耗函數(shù)');grid onfigure(3)plot(wk/pi,angle(Hk(k+1)/pi);gridonxlabel('w八pi');ylabel('相位八pi'); title('相頻特性曲線);(2)數(shù)字濾波器的損耗函數(shù)和相頻特性分別如圖11、12所示：損耗函數(shù)相頻特性曲線圖11損耗函數(shù)曲線圖12相頻特性曲線二、按直接型網(wǎng)絡結構編程編寫濾波程

21、序(1)matlab 代碼：N=500;n=0:N-1;f=2800;T=1/f;t=n*T;x=0.5*cos(2*pi*100*t)+0.7*cos(2*pi*130*t+0.1*pi)+0.4*cos(2*pi*600*t)；m1=0;m2=0;m3=0;m4=0;m5=0;m6=0;m7=0;m8=0;m9=0;m10=0;m11=0;m12=0;m1 3=0;m14=0;m15=0;m16=0;m17=0;m18=0;m19=0;m20=0;m21=0;m22=0;m23=0;m24=0;for m=1:length(x) y(m)=0.0012*x(m)+m1*0.0011-m2*

22、0.0014-m3*0.0072-m4*0.0147-m5*0.019 3-.m6*0.0145+m7*0.0055+m8*0.0423+m9*0.0910+m10*0.1410+m11*0.1786+m12*0.1926+.m13*0.1786+m14*0.1410+m15*0.0910+m16*0.0423+m17*0.0055-m18*0.0145 m19*0.0193-m20*0.0147-m21*0.0072-m22*0.0014+m23*0.0011+m24*0.0012;m24=m23;m23=m22;m22=m21;m21=m20;m20=m19;m19=m18;m18=m17

23、;m17=m16;m 16=m15;m15=m14;m14=m13;m13=m12;m12=m11;m11=m10;m10=m9;m9=m8;m8=m7;m7=m6;m6=m5;m5=m4;m4=m3;m3=m2;m2=m1;m1=x(m);endplot(n,x);title( '信號x(n)');ylabel( '幅值');xlabel( 'n');S=fft(x,N);fs=n/(N*T);figure(2)plot(fs,abs(S);axis(0,1500,0,180);title('原信號幅度頻譜(采樣點數(shù)為500)'

24、);xlabel('頻率/Hz' );ylabel( '幅值');figure(3)plot(n,y);title('信號 y(n)');ylabel('幅值');xlabel( 'n');S=fft(y,N);fs=n/(N*T);figure(4)plot(fs,abs(S);axis(0,1500,0,160);title( '幅度頻譜');xlabel( '頻率/Hz' );ylabel( '幅值);(2)原信號及其幅度頻譜分別如圖13、14所示:信號x(n)21.5

25、-2，050100150200250300350400450500n圖13信號x(n)波形圖14幅度頻譜(3)濾波后信號y(n)及其幅度頻譜分別如圖15、16所示:圖15信號y(n)波形圖16幅度頻譜2)將fi和f2頻率分量以高于50dB的衰減抑制，同時以低于2dB的衰減通過f3頻率分量；一、基于頻率采樣法設計 FIR數(shù)字濾波器：(1)matlab 代碼:T=0.48;Fs=3800;fp=600;fs=100;wp=2*pi*fp/Fs;ws=2*pi*fs/Fs; datB=wp-ws;wc=wp; m=1;N=ceil(m+1)*2*pi/datB+1); N=N+mod(N+1,2);

26、Np=fix(wc/(2*pi/N);Ns=N-2*Np-1;Ak=zeros(1,Np+1),ones(1,Ns),zeros(1,Np); Ak(Np+2)=T;Ak(N-Np尸T; thetak=-pi*(N-1)*(0:N-1)/N;Hk=Ak.*exp(1j*thetak); hn=real(ifft(Hk); M=1024; Hk=fft(hn,M); k=0:M/2-1; wk=(2*pi/M)*k; figure(2);plot(wk/pi,20*log10(abs(Hk(k+1);axis(0,1,-80,5);xlabel('w/pi'title( '

27、;損耗函數(shù)');grid on figure(3) plot(wk/pi,angle(Hk(k+1)/pi);grid xlabel('w八pi');ylabel( '相位八pi'title('相頻特性曲線););ylabel( '幅度/dB');on)；(2)數(shù)字濾波器的損耗函數(shù)和相頻特性分別如圖17、18所示:圖17損耗函數(shù)曲線圖18相頻特性曲線二、按直接型網(wǎng)絡結構編程編寫濾波程序：(1)matlab 代碼：N=500;n=0:N-1;f=2800;T=1/f;t=n*T;x=0.5*cos(2*pi*100*t)+0.7*c

28、os(2*pi*130*t+0.1*pi)+0.4*cos(2*pi*600*t)；m1=0;m2=0;m3=0;m4=0;m5=0;m6=0;m7=0;m8=0;m9=0;m10=0;m11=0;m12=0;m13=0;m14=0;m15=0;m16=0;for m=1:length(x)y(m)=-0.0009*x(m)-m1*0.0169-m2*0.0128+m3*0.0282+m4*0.0627+m5*0.0198- .m6*0.1198-m7*0.2827+m8*0.6447-m9*0.2827-m10*0.1198+m11*0.0198+m12*0.0627+.m13*0.0282

29、-m14*0.0128-m15*0.0169-m16*0.0009;m16=m15;m15=m14;m14=m13;m13=m12;m12=m11;m11=m10;m10=m9;m9=m8;m8=m 7;m7=m6;m6=m5;m5=m4;m4=m3;m3=m2;m2=m1;m1=x(m);endplot(n,x);title( '信號x(n)' );ylabel( '幅值');xlabel( 'n');S=fft(x,N);fs=n/(N*T);figure(2)plot(fs,abs(S);axis(0,1500,0,180);title(&

30、#39;原信號幅度頻譜(采樣點數(shù)為500)');xlabel('頻率/Hz' );ylabel( '幅值');figure(3)plot(n,y);title( '信號 y(n)');ylabel('幅值');xlabel( 'n');S=fft(y,N);fs=n/(N*T);figure(4)plot(fs,abs(S);axis(0,1500,0,160);title( '幅度頻譜)；xlabel('頻率/Hz' );ylabel( '幅值');(4)原信號及其

31、幅度頻譜分別如圖19、20所示：信號x(n)原信號幅度頻譜(采樣點數(shù)為500)2180 值幅001500圖19信號x(n)波形圖20幅度頻譜(5)濾波后信號y(n)及其幅度頻譜分別如圖21、22所示:幅度頻譜160信號y(n)0.8-0.61111,111050100150200250300350400450500n值幅5001000頻率/Hz1500圖21信號y(n)波形圖22幅度頻譜要求：按數(shù)字濾波器直接型結構圖編寫濾波程序，求得 y(n)； 1)中的FIR濾波器采用窗函數(shù)法設計；2)中的FIR濾波器采用頻率采樣法設計。畫出所設計的濾波器頻率特性圖、信號時域圖；給出濾波器設計的 MATLA

32、B代碼與濾波器實現(xiàn)的代碼；選擇合適的信號采樣周期To3)與第6章作業(yè)2的IIR濾波方法進行比較研究。一、低通濾波器部分：(1)matlab 代碼：Fs=3800;fp=130;fs=600;rs=50;wp=2*pi*fp/Fs;ws=2*pi*fs/Fs;Bt=ws-wp;alph=0.5842*(rs-21)A0.4+0.07886*(rs-21);N=ceil(rs-8)/2.285/Bt);wc=(wp+ws)/2/pi;hn=fir1(N,wc,kaiser(N+1,alph);M=1024;Hk=fft(hn,M);k=0:M/2-1;wk=(2*pi/M)*k;wp2=2*fp/

33、Fs;ws2=2*fs/Fs;Rp=2;As=50;N1,wp1=ellipord(wp2,ws2,Rp,As);B,A=ellip(N1,Rp,As,wp1)Hk1,wk1=freqz(B,A);mag=abs(Hk1);pah=angle(Hk1);on':','linewidth',3);gridonplot(wk1/pi,20*log10(mag);grid hold onplot(wk/pi,20*log10(abs(Hk(k+1), onhold offxlabel('w八pi');ylabel(title('損耗函數(shù)曲線&#

34、39;);legend( 'IIR' , 'FIR');figure(2)plot(wk1/pi,pah/pi);gridhold onplot(wk/pi,angle(Hk(k+1)/pi, hold offxlabel('w八pi');ylabel(title('相頻特性曲線)； legend( 'IIR' , 'FIR');幅度/dB')on':','linewidth' ,3);grid相位八pi');(2)兩種濾波器的損耗函數(shù)、相頻特性的比較圖見圖2

35、3、 24:損耗函數(shù)曲線加. * itI*iIIR11. VV #i.*t T* 事1* a r .K4, * .sS . : a 'tI>14-iS” ： , n-f* « H"*I* 4ft,: ii.i胃 ll.I a t-| -ft *«> 1 1 .» W以L*i1事; n =; * * «；»i-*A 1* *1 P HM n k 1 _ , k 11. * « I* i 事. . * *1 I il* r h； .* *W 4 i- r * *B .；» ,* K:JV 1I； 41-一同,! «1 p： 1 J_L二 u1 , Iw*m v1 ;* «相頻特性曲線0.80.60.40.20-0.2-0.4-0.6-0.8/二 位相00.10.20.30.40.50.60.70.80.91w/ ：-1圖24相頻特性比較圖圖23損耗函數(shù)比較圖IIR濾波器的階數(shù)：N1=3FIR濾波器的階數(shù)：N=17(1)matlab 代碼：T=0.48;Fs=3800;fp=600;fs=100;wp=2*pi*fp/Fs;ws=2*pi*fs/Fs;datB=wp-ws;wc=wp;m=1;N=ceil(m+1)