通信工程數(shù)字信號(hào)實(shí)驗(yàn)5-有限長(zhǎng)單位脈沖響應(yīng)濾波器設(shè)計(jì)

實(shí)驗(yàn)五有限長(zhǎng)單位脈沖響應(yīng)濾波器設(shè)計(jì)實(shí)驗(yàn)?zāi)康恼莆沼么昂瘮?shù)法、頻率采樣法及優(yōu)化設(shè)計(jì)法設(shè)計(jì)FIR濾波器的原理及方法，熟悉相應(yīng)的MATLAB編程。熟悉線性相位FIR濾波器的幅頻特性和相頻特性。了解各種不同窗函數(shù)對(duì)濾波器性能的影響。二、實(shí)驗(yàn)原理window=ones(1,N):產(chǎn)生N點(diǎn)矩形窗，行向量。window=hann(N):產(chǎn)生N點(diǎn)漢寧窗，列向量。window=hanning(N):產(chǎn)生N點(diǎn)非零漢寧窗，列向量。等價(jià)于去除hann(N+2)的第一個(gè)零元素和最后一個(gè)零元素，得到的N點(diǎn)非零窗函數(shù)。window=hamming(N):產(chǎn)生N點(diǎn)海明窗，列向量。window=blackman(N):產(chǎn)生N點(diǎn)布萊克曼窗，列向量。window=kaiser(N,beta):產(chǎn)生參數(shù)為beta的N點(diǎn)凱塞窗，列向量。[M,Wd,beta,ftype]=kaiserord(f,a,dev,fs):凱塞窗參數(shù)估計(jì)。f為一組邊界頻率，最高頻率為fs/2。a為f中各個(gè)頻帶的幅度值，通帶取1，阻帶取0。如果f中有2個(gè)元素，則形成3個(gè)頻帶，其中第1個(gè)和第3個(gè)是通帶或阻帶，第2個(gè)是過(guò)渡帶，a中也有2個(gè)元素，指明第1個(gè)和第3個(gè)頻帶是通帶還是阻帶；如果f中有4個(gè)元素，則形成5個(gè)頻帶，其中1，3和5是通帶或阻帶，2和4是過(guò)渡帶，a中有3個(gè)元素，指明1，3和5是通帶還是阻帶。dev的維數(shù)與a相同，指明每個(gè)頻帶上的波動(dòng)值。fs為采樣頻率。M為FIR濾波器的階數(shù)，M=N-1。Wd為歸一化邊界頻率，等于數(shù)字邊界角頻率除以π，或者邊界頻率除以fs/2。beta就是凱塞窗的參數(shù)β。ftype為濾波器的類型。fs=8000;f=[10001500];a=[10];dev=[0.050.01];[M,Wd,beta,ftype]=kaiserord(f,a,dev,fs);h=fir1(M,Wd,ftype,kaiser(M+1,beta));freqz(h);fs=8000;f=[1000130022102410];a=[010];dev=[0.010.050.01];[M,Wd,beta,ftype]=kaiserord(f,a,dev,fs);hh=fir1(M,Wd,ftype,kaiser(M+1,beta));[H,f]=freqz(hh,1,1024,fs);figure;plot(f,abs(H));gridon;b=fir1(M,Wd,'ftype',window):用窗函數(shù)法求FIR濾波器的系數(shù)b(單位脈沖響應(yīng))。M為濾波器的階數(shù)，M=N-1。Wd為一組歸一化邊界頻率，通帶和阻帶間隔分布，無(wú)過(guò)渡帶；只有一個(gè)元素，表示低通或高通濾波器；有兩個(gè)元素表示帶通和帶阻濾波器；有三個(gè)及以上元素，表示多帶濾波器。'ftype'表示濾波器類型，'high'表示高通濾波器，'stop'表示帶阻濾波器，'DC-0'表示多帶濾波器的第一個(gè)頻帶為阻帶，'DC-1'表示多帶濾波器的第一個(gè)頻帶為通帶。window為窗口類型，缺省為海明窗。b=fir2(M,f,m,window):用頻率采樣法求FIR濾波器的系數(shù)b。M為濾波器的階數(shù)，M=N-1。f為一組歸一化頻率，第一個(gè)元素必須為0，最后一個(gè)元素必須為1(對(duì)應(yīng)奈奎斯特頻率，即采樣頻率的一半)，中間的元素按升序排列。m的維數(shù)與f相同，指明f中每個(gè)頻率上的理想幅度。window為窗口類型，缺省為海明窗。Fir2可以實(shí)現(xiàn)任意幅度特性的濾波器。三、實(shí)驗(yàn)內(nèi)容1、用窗函數(shù)法設(shè)計(jì)一個(gè)線性相位FIR高通濾波器，通帶邊界頻率為0.6π，阻帶邊界頻率為0.4π，要求阻帶衰減不小于40dB，窗函數(shù)從矩形窗、漢寧窗、海明窗和布萊克曼窗中選取，且要求濾波器的階數(shù)最小。給出該濾波器的單位脈沖響應(yīng)h(n)的解析式，并用MATLAB繪出時(shí)域波形和幅頻特性。2、用窗函數(shù)法設(shè)計(jì)一個(gè)線性相位FIR帶通濾波器，通帶邊界頻率為4kHz和6kHz，阻帶邊界頻率為2kHz和8kHz，采樣頻率為20kHz，要求阻帶衰減不小于50dB，窗函數(shù)從矩形窗、漢寧窗、海明窗和布萊克曼窗中選取，且要求濾波器的階數(shù)最小。給出該濾波器的單位脈沖響應(yīng)h(n)的解析式，并用MATLAB繪出時(shí)域波形和幅頻特性。3、用窗函數(shù)法設(shè)計(jì)一線性相位FIR帶阻濾波器，通帶邊界頻率為600Hz和1400Hz，阻帶邊界頻率為800Hz和1200Hz，采樣頻率為4000Hz，要求阻帶衰減不小于50dB，窗函數(shù)從矩形窗、漢寧窗、海明窗和布萊克曼窗中選取，且要求濾波器的階數(shù)最小。給出該濾波器的單位脈沖響應(yīng)h(n)的解析式，并用MATLAB繪出時(shí)域波形和幅頻特性。4、用凱塞窗設(shè)計(jì)一個(gè)多帶線性相位濾波器，幅頻特性如下圖所示，N=40，β分別取4、6、10，比較不同β值時(shí)的幅頻特性和相頻特性。5、用頻率采樣法設(shè)計(jì)第(4)題的多帶濾波器，比較過(guò)渡帶不設(shè)采樣點(diǎn)和設(shè)置一個(gè)采樣點(diǎn)的幅頻特性，過(guò)渡點(diǎn)的值設(shè)為0.5，窗函數(shù)選矩形窗。四、實(shí)驗(yàn)結(jié)果與分析1、用窗函數(shù)法設(shè)計(jì)一個(gè)線性相位FIR高通濾波器，通帶邊界頻率為0.6π，阻帶邊界頻率為0.4π，要求阻帶衰減不小于40dB，窗函數(shù)從矩形窗、漢寧窗、海明窗和布萊克曼窗中選取，且要求濾波器的階數(shù)最小。給出該濾波器的單位脈沖響應(yīng)h(n)的解析式，并用MATLAB繪出時(shí)域波形和幅頻特性。解：由題可得:求理想低通濾波器的邊界頻率。理想低通濾波器的邊界頻率是實(shí)際濾波器幅度函數(shù)正負(fù)肩峰頻率的中心，但這里不知道正負(fù)肩峰的頻率，因此用通帶和阻帶邊界頻率的中心近似:選擇窗函數(shù)和窗口長(zhǎng)度。阻帶衰減不小于40dB，因此選擇海明窗(其阻帶最小衰減為-44dB)；濾波器的過(guò)渡帶寬為0.6π-0.4π=0.2π，因此窗口長(zhǎng)度為:線性相位延遲常數(shù)為:根據(jù)理想邊界頻率和線性相位延遲常數(shù)，求理想單位脈沖響應(yīng):窗函數(shù)與理想單位脈沖響應(yīng)相乘，得到線性相位FIR低通濾波器的單位脈沖響應(yīng):matlab程序見(jiàn)附錄程序一運(yùn)行結(jié)果如下：分析：由圖知阻帶衰減最小值大于40，滿足要求。2、用窗函數(shù)法設(shè)計(jì)一個(gè)線性相位FIR帶通濾波器，通帶邊界頻率為4kHz和6kHz，阻帶邊界頻率為2kHz和8kHz，采樣頻率為20kHz，要求阻帶衰減不小于50dB，窗函數(shù)從矩形窗、漢寧窗、海明窗和布萊克曼窗中選取，且要求濾波器的階數(shù)最小。給出該濾波器的單位脈沖響應(yīng)h(n)的解析式，并用MATLAB繪出時(shí)域波形和幅頻特性。解：求數(shù)字邊界頻率:求理想濾波器的邊界頻率:求理想單位脈沖響應(yīng):選擇窗函數(shù)。阻帶最小衰減為-50dB，因此選擇漢明窗(其阻帶最小衰減為-53dB)；濾波器的過(guò)渡帶寬為0.4π-0.2π=0.8π-0.6π=0.2π，因此求FIR濾波器的單位脈沖響應(yīng)h(n):matlab程序見(jiàn)附錄程序二運(yùn)行結(jié)果如下：分析：由圖可知，阻帶衰減不小于50dB滿足要求。3、用窗函數(shù)法設(shè)計(jì)一線性相位FIR帶阻濾波器，通帶邊界頻率為600Hz和1400Hz，阻帶邊界頻率為800Hz和1200Hz，采樣頻率為4000Hz，要求阻帶衰減不小于50dB，窗函數(shù)從矩形窗、漢寧窗、海明窗和布萊克曼窗中選取，且要求濾波器的階數(shù)最小。給出該濾波器的單位脈沖響應(yīng)h(n)的解析式，并用MATLAB繪出時(shí)域波形和幅頻特性。解：求數(shù)字邊界頻率:求理想濾波器的邊界頻率:求理想單位脈沖響應(yīng):選擇窗函數(shù)。阻帶最小衰減為-50dB，因此選擇漢明窗(其阻帶最小衰減為-53dB)；濾波器的過(guò)渡帶寬為0.4π-0.3π=0.7π-0.6π=0.1π，因此求FIR濾波器的單位脈沖響應(yīng)h(n):matlab程序見(jiàn)附錄程序三運(yùn)行結(jié)果如下：分析：由圖可知，阻帶衰減不小于50dB，滿足要求。4、用凱塞窗設(shè)計(jì)一個(gè)多帶線性相位濾波器，幅頻特性如下圖所示，N=40，β分別取4、6、10，比較不同β值時(shí)的幅頻特性和相頻特性。解：matlab程序見(jiàn)附錄程序四運(yùn)行結(jié)果如下：分析：β=4時(shí)幅頻圖比較接近理想濾波器，相頻特性不太明顯，隨著β的增大，幅頻圖過(guò)渡帶開(kāi)始增加，相頻特性越來(lái)越明顯。這可能是因?yàn)樵龃螃聦?dǎo)致頻譜的主瓣寬度增加造成的。5、用頻率采樣法設(shè)計(jì)第(4)題的多帶濾波器，比較過(guò)渡帶不設(shè)采樣點(diǎn)和設(shè)置一個(gè)采樣點(diǎn)的幅頻特性，過(guò)渡點(diǎn)的值設(shè)為0.5，窗函數(shù)選矩形窗。解：由于(2π/40)×4=0.2π，(2π/32)×8=0.4π，(2π/40)×12=0.6π，(2π/40)×16=0.8π，又因?yàn)镹為偶數(shù)，屬于第二類線性相位濾波器，幅度函數(shù)在0~2π上以π為中心呈奇數(shù)對(duì)成分布。因此頻率采樣值的幅度為:頻響采樣的相位為:matlab程序見(jiàn)附錄程序五運(yùn)行結(jié)果如下：分析：在理想頻率響應(yīng)的間斷點(diǎn)的邊緣加上一個(gè)過(guò)渡采樣點(diǎn)，阻帶衰減變大，但過(guò)渡帶也就越寬。五、實(shí)驗(yàn)小結(jié)通過(guò)此次實(shí)驗(yàn)，掌握了用窗函數(shù)法、頻率采樣法及優(yōu)化設(shè)計(jì)法設(shè)計(jì)FIR濾波器的原理及方法，熟悉了相應(yīng)的MATLAB編程。熟悉了線性相位FIR濾波器的幅頻特性和相頻特性。了解了各種不同窗函數(shù)對(duì)濾波器性能的影響。在實(shí)驗(yàn)中，明白了matlab在數(shù)字信號(hào)處理中的重要性，學(xué)好matlab對(duì)以后工作學(xué)習(xí)有極大的幫助。六、附錄程序一clear;N=31;n=0:N-1;hd=(sin(pi*(n-15))-sin(0.5*pi*(n-15)))./(pi*(n-15));hd(16)=0.5;win=hanning(N);h=win'.*hd;figure;stem(n,h);xlabel('n');ylabel('h(n)');grid;title('FIR高通濾波單位脈沖響應(yīng)h(n)');[H,w]=freqz(h,1);H=20*log10(abs(H));figure;plot(w/pi,H);axis([01-10010]);xlabel('\omega/\pi');ylabel('幅度/dB');grid;title('FIR高通濾波器，hanning窗，N=31');程序二clear;N=33;n=0:N-1;hd=(sin(0.7*pi*(n-16))-sin(0.3*pi*(n-16)))./(pi*(n-16));hd(17)=0.4;win=hamming(N);h=win'.*hd;figure;stem(n,h);xlabel('n');ylabel('h(n)');grid;title('FIR帶通濾波單位脈沖響應(yīng)h(n)');[H,w]=freqz(h,1);H=20*log10(abs(H));figure;plot(w/pi,H);axis([01-10010]);xlabel('\omega/\pi');ylabel('幅度/dB');grid;title('FIR帶通濾波器，hamming窗，N=33');程序三clear;N=66;n=0:N-1;hd=(sin(pi*(n-32.5))+sin(0.35*pi*(n-32.5))-sin(0.65*pi*(n-32.5)))./(pi*(n-32.5));win=hamming(N);h=win'.*hd;figure;stem(n,h);xlabel('n');ylabel('h(n)');grid;title('FIR帶通濾波單位脈沖響應(yīng)h(n)');[H,w]=freqz(h,1);H=20*log10(abs(H));figure;plot(w/pi,H);axis([01-10010]);xlabel('\omega/\pi');ylabel('幅度/dB');grid;title('FIR帶阻濾波器，hamming窗，N=66');程序四clear;%β=4Wd=[0.20.40.60.8];M=39;beta=4;hh=fir1(M,Wd,'DC-0',kaiser(M+1,beta));[H,w]=freqz(hh,1);figure;subplot(2,1,1);plot(w/pi,abs(H));xlabel('\omega/\pi');ylabel('幅度');grid;title('幅頻特性,beta=4');subplot(2,1,2);plot(w/pi,angle(H));xlabel('\omega/\pi');ylabel('弧度');grid;title('相頻特性');%β=6;Wd=[0.20.40.60.8];M=39;beta=6;hh=fir1(M,Wd,'DC-0',kaiser(M+1,beta));[H,w]=freqz(hh,1);figure;subplot(2,1,1);plot(w/pi,abs(H));xlabel('\omega/\pi');ylabel('幅度');grid;title('幅頻特性,beta=6');subplot(2,1,2);plot(w/pi,angle(H));xlabel('\omega/\pi');ylabel('弧度');grid;title('相頻特性');%β=10;Wd=[0.20.40.60.8];M=39;beta=10;hh=fir1(M,Wd,'DC-0',kaiser(M+1,beta));[H,w]=freqz(hh,1);figure;subplot(2,1,1);plot(w/pi,abs(H));xlabel('\omega/\pi');ylabel('幅度');grid;title('幅頻特性,beta=10');subplot(2,1,2);plot(w/pi,angle(H));xlabel('\omega/\pi');ylabel('弧度');grid;title('相頻特性');程序五c