基于matlab的布萊克曼窗函數法設計的低通濾波器
第一章設計要求 1.1基本要求 (1)理解FIR數字低通濾波器的作用及應用領域,掌握布菜克曼窗函數法設計 FIR數字低通濾波器的原理及其在 matlab仿真中的實現方法 (2)掌握 matlab的編程方法; (3)通過腳本編程或 SIMULINK實現FIR數字LPF;使用布萊克曼窗函數;M-1, n=[0:1:M-1;Wc=0.2米pi; (4)完成課程設計報告 1.2提高要求 (1)實現Wc和M均可交的布萊克曼窗 (2)用設計岀的濾波器對聲音信號加噪聲后進行濾波,對濾波前后信號進行對 比,并在時域和頻域進行分別進行分析。 第二章設計原理 2.1FIR濾波器 FIR( Finite Impulse response)濾波器:有限長單位沖激響應濾波器,又稱 為非遞歸型濾波器,是數字信號處理系統中最基本的元件,它可以在保證幅頻特 性的冋時具有嚴格的線性相頻特性,同時其單位抽樣響應是有限長的,因而濾波 器是穩定的系統。因此,FIR濾波器在通信、圖像處理、模式識別等領域都有著 廣泛的應用。 有限長單位沖激響應(FIR)濾波器有以下特點 (1)、系統的單位沖激響應h(n)在有限個n值處不為零: (2)、系統函數H(z)在z>0處收斂,極點全部在z=0處(因果系統) (3)、結構上主要是非遞歸結構,沒有輸出到輸入的反饋,但有些結構中(例 如頻率抽樣結構)也包含有反饋的遞歸部分。 FIR低通濾波器的設計方法主要分為兩類:第一類是基于逼近理想濾波器特 性的方法,包括窗函數法、頻夲采樣法和等波紋最佳逼近法:第二類是最優設計 法。實驗采用布萊克曼窗函數法設計。 2.2窗函數設計法 窗函數設計法的基木思路是用FIR低通濾波器逼近希望的濾波特性。設希 望逼近的濾波器的頻率響應函數為Hd(e"),其單位脈沖響應用hd(n)表示。其 中 Hd(e/w)=y hd(n e/ wn hd()=1∫ Hd(elwewndw WC hd(n)是無限時寬的,是非因果序列。窗函數設計法就是截取hdan)為有限 長的一段因果序列,并用合適的窗函數進行加權做為FIR低通濾波器的單位脈 沖響應h(n)。 窗口法的設計步驟如下: (1)、通過傅里葉變換后的理想濾波器的單位脈沖響應hd(n)。 (2)、根據指標選擇窗口截止頻率Wc和階數M。 (3)、給定理想頻響由Hd(e")和hd(n),加窗得h(n)=w(n)hd(n) 窗函數又包含有多種,包括矩形窗,三角窗,漢明窗,布萊克曼窗等,實驗 采用布萊克曼窗 2.3布萊克曼窗函數 布萊克曼窗的時域形式可表示為: wvn()=(042-05c0s2m+0.08c0sm)R2n(n) 它的頻域特性為: Wa(ey)=042(y)-02WR(e0x-)+Wk{(e)+ 0.04[WR(ejw-w-1)+WR(e(wt N-1) (2-4) 其幅度函數為 Wn(o)=042Wn()-0.25(0-2x)+W/(+2)1+ 00(-)+Wm(+) (2-5) 增加一個二次諧波余弦分量,可進一步降低旁瓣,但主瓣寬度進一步增加,為 12π∧N。加N可減少過渡帶。布萊克曼窗函數的最大旁瓣之比主鐮值低57db, 但是主瓣寬度是矩形窗函數的主瓣寬度的三倍。布萊克曼窗主瓣寬,旁瓣小,頻 率識別精度最低,但幅值識別精度最高 第三章設計步驟及結果分析 實驗分基礎和提髙兩部分,基礎部分實現窗函數階數與截止頻率均固定的布 萊克曼窗,提高部分實現窗函數階數與截止頻率均可變的布萊克曼窗并測試通過 不同頻率雜音信號,實現濾除高頻分量,并且顯示其時域及頻域圖 3.1基礎部分 3.1.1設計流程圖 輸入M,Wc 產生單位康 向 生布菜克曼 加窗了 頻率響應 繪圈 結束 圖3.1基礎部分流程圖 由流程圖,基礎部分通過輸入固定M,Wc,調用已有函數進行布菜克曼濾波 器設計,然后繪制波形圖進行觀察。 3.1.2程序設計 (1) i deal_1p函數 該函數用于產生理想低通脈沖響應 function hd=ideal lp(wc, N) alpha(N-1)/2 n=[0:1:N-1] m=n-alpha hd-sin(wc*m). /(pi*m) en (2)窗函數階數與截止頻率均固定的布萊克曼窗程序 clc: clear all 4 M-11; n=[0:l:M-1]; Wc=0.2*p1; hd= ideal_lp(We,M);%產生單位脈沖響應 Wbla=( blackman(M)’;%產生布萊克曼窗 h-hd.*Wbla;%加窗 [H,W]= frege(h,I);%頻率響應 If 1figure: set(Hf 1, NumberTitle', off, Name, M-16, Wc-0 2*pi 1i3k 克曼窗’);%設置標題 subplot(2,2,1);stem(n,hd); title('理想脈沖響應’);axis(0,M1,0.1,0 3]);xlabel(M); ylabel(hd (M)); grid subplot(2,2,2);stem(n,wbla); title('布萊克曼窗’);axis([0,M-1,0,1,1]); xlabel('w): ylabel(wBl(M)); grid subplot(2,2,3);stem(m,h); Litle(實際脈沖響應’);axis([0,M-1,-0.1,0.3]); xlabel(M); ylabel(h(M))); grid subplot(2, 2,1); plot(w/pi, 20*log10(abs(h))); axis([o, 1,-200, 10]); title (幅度響應/dB'); xlabe1(o/π’); y labe(20lg|Ig(o)');grid 3.1.3實驗結果及其分析 M=11c=02pi菜克窗 File Edit View Insert Tools Desktop Window Help 固出冷氣國函鳥日圖日回 NaN+0000000coo000000i Nan NeNi 理想脈沖響應 布菜克曼窗 Nan+ NaN Nan+ NaNi NaN十 NaNI Nan+ NaNI m0.5 Nan Nani Nan+ NaNi NaN+ NaN↓ NaNI NaN十 NeNi 實際尿沖響應 幅度響應d NaN十 Neni 33 Nan+ NaNi Nan+ NaN Nan+ NaNi Nan NeNi Nan+ NaN Nan+ NaN NaN+ NaN WTT NaN十 NENI 圖3.2M=11時的窗函數設計 由圖,布萊克曼窗已經顯示出來,但由于階數太小,幅度響應并不能顯示(見 圖32左),因為經過 fretz函數后,幅度是非數(見圖32右)。經測試,當M>=16 時,幅度響應才能夠正確顯示。可以看出,效果并不是很好(見圖3.3)。 M-16W-02*pi萊克曼窗 File Edit View Insert Tools Desktop Window Hell 理想脈 布菜克曼窗 C.3 .2 .1 0.5 5 15 15 實際脈沖響應 幅度響應dB .1 15 10 15 圖3.3M=16時的窗函數設計 3.2提高部分 3.2.1設計流程圖 本課程設計是對錄制的語音信號進行加噪處理并分析加噪前后語音信號的 時域圖與頻域圖,再用布萊克曼窗設計一個FIR濾波器,而濾波器的設計必須 符合其相應的指標,否則將不能濾掉加入的噪聲。最后將濾波前后的波形圖進行 比較看是否相同。下面是整個設計的流程圖如圖3.4所示: 開始 原始信號時域圖←讀取語音信號}原始信號頻譜 加噪信號時城圖←指音信號加噪}“加信號譜 否 是否輸 入MWc 是 默以。個用讀器 匚濾波信號時域圖←「濾波器濾除噪聲濾波信號煩譜 繪圖 結束 圖3.4語音信號的整個流程圖 由圖知,系統首先通過手機錄制聲音信號,然后用 matlab軟件讀取,再對 語音信號加噪聲,經過設計好的濾波器進行濾波,最后比鉸通過濾波器前后的聲 音信號的時域及頻域的區別。 由于整個系統過于復雜,因此分為以下幾個子系統分別進行實現 3.2.2讀取語音信號 打開 MATLAB軟件平臺,將用手機錄制好的語音信號(為了區別低頻高頻, 這里通過講話聲和口哨聲區分)利用函數 audiohead對語音信號進行采樣,記住 采樣頻率的和釆樣點數。再凋用函數 sound函數聽錄制的語音。采樣完后再語音 信號中加入一個單頻噪聲,單頻的噪聲的頻率可以自己設置。按照加入噪聲后的 采樣頻率再次調用 sound函數,這時可以明顯的聽見播放的語音信號中有尖銳的 噪聲。 fs=44100 x1- audiohead(C: USers、Al1 bymyself\ Desktop\Fir.aac;%輸入參數為文件的 全路徑和文件名,輸出參數是每個樣木的值 t-0:1/44100:(1 ength(x1)-1)/44100;%計算時間范圍,樣本數除以采樣頻率 yl=fft(x1,441000); f=fs米(0:220499)/441000 % sound(x1,44100);%原始聲音信號 Au=0.01; d=LAU*COsS(2*p*50000kt)”;%單頻噪聲,頻率參數可改 x2=x1+d; % sound(x2,44100);%明顯聽出有尖銳的單頻嘯叫聲 y2=fft(x2,111000 3.2.3低通濾波器設計 (1) black man函數設計 該函數用于產生布萊克曼窗,需要注意的是這個函數與自帶函數b1 ackman產生 列向量不同,它產生的是行向量 function[x]= black man (N) n-0:N-1; x=(0.42-0.5米COS((2來pi*n)/(N-1)+0.08*COs(4**n)/(N-1)) end 7 (2) blackman2 kebian數設計 該函數為頻率和階數均可變的布萊克曼窗函數此函數設計了缺省參數,當不 輸入參數,或只輸入一個參數時,函數會設置默認參數,并且把Ⅵ和Wc顯示在 左上角,方便觀察 function [ blackman2kebian(, Wc) %設置缺省參數,當輸入一個參數,默認Wc=0.2*pi,當輸入零個參數,默認M-=11 Wc=0. 2*pi if nargin--1 Wc=0. 2*pi; end if nargin==0 M=ll; Wc=0. 2*pi; end %設置標題 a b=num2str(m) c=strcat(a, b) d=,wC b=num2str(Wc) d=strcat(a, b) e=strcat(c, d) e-strcat( blackman', e) e= sTreat(e,”)’) Hf 1=figure; set(Hf 1, NumberTitle', off, Name', e n-[0:1:M-1; hd= ideal_lp(We,M;%理想低通濾波器脈沖響應 wbla=( black man(M));%調用自編 black man窗 hn=hd.水wbla; H,W]= fretz(hn,1);%頻率參數 繪圖 subplot(2,2,1);stem(n,hd); title('理想脈沖響應’);axis(0,M1,0.1,0. 3); xlabel(M); ylabel(hd(M))); grid subplot(2,2,2);stem(m,wbla); lille(布萊克曼窗’);axis([0,M-1,0,1.1]) xlabel(M); ylabel(’wB1(M’);grid subplot(2,2,3);stem(m,hn); title('實際脈沖響應');axis(0,M-1,-0.1,0