布萊克曼(Blackman)FIR濾波器課程設計

1、成績 布萊克曼（Blackman）FIR濾波器課程設計目錄一、設計任務與要求2二、總體方案設計2三、單元電路設計與參數(shù)計算33.1 FIR濾波器33.2 窗口設計法33.3 布萊克曼窗4四、性能測試與分析54.1語音信號錄制54.2 濾波器設計74.3信號濾波處理84.4結(jié)果分析9五、結(jié)論與心得10六、參考文獻10七、附錄117.1 布萊克曼窗設計FIR濾波器的MATLAB主程序117.2 調(diào)用函數(shù)ideal_lp.m的MATLAB源代碼137.3 調(diào)用函數(shù)freqz_m.m的MATLAB源代碼13布萊克曼窗FIR濾波器的設計一、設計任務與要求1濾波器指標必須符合工程實際。； 2設計完后應檢查

2、其頻率響應曲線是否滿足指標。；3. 獨立完成課程設計并按要求編寫課程設計報告書。二、總體方案設計本課程設計是對錄制的語音信號進行加噪處理并分析加噪前后語音信號的時域圖與頻域圖，再用布萊克曼窗設計一個FIR濾波器，而濾波器的設計必須符合其相應的指標，否則將不能濾掉加入的噪聲。最后將濾波前后的波形圖進行比較看是否相同。課程設計流程圖如圖1所示。開始錄制語音信號分析語音信號的頻譜用布萊克曼窗設計FIR濾波器在語音信號中加入噪聲用濾波器對信號進行濾波比較濾波前后的波型及頻譜回放語音信號結(jié)束圖1：課程設計流程圖三、單元電路設計與參數(shù)計算3.1 FIR濾波器FIR(Finite Impulse Respo

3、nse)濾波器：有限長單位沖激響應濾波器，是數(shù)字信號處理系統(tǒng)中最基本的元件，它可以在保證任意幅頻特性的同時具有嚴格的線性相頻特性，同時其單位抽樣響應是有限長的，因而濾波器是穩(wěn)定的系統(tǒng)。FIR數(shù)字濾波器以其良好的線性特性被廣泛應用于現(xiàn)代電子通信系統(tǒng)中，是數(shù)字信號處理的重要內(nèi)容之一。在實際信號處理中，往往要求系統(tǒng)兼具實時性和靈活性，而已有的一些軟件或硬件實現(xiàn)方案(如DSP)則難以同時達到這兩方面的要求。使用具有并行處理特性的FPGA來實現(xiàn)FIR濾波器，既有很強的實時性，又兼顧了靈活性，為數(shù)字信號處理提供了一種很好的解決方案。設FIR濾波器的單位沖激響應h(n)為一個N點序列，0nN1，則濾波器的系

4、統(tǒng)函數(shù)為： (2-1)就是說，它有階極點在處，有個零點位于有限z平面的任何位置因此是永遠穩(wěn)定的。穩(wěn)定和相位特性是FIR濾波器突出的優(yōu)點。FIR濾波器有以下幾種基本結(jié)構(gòu)：橫截型（卷積型、直接型）、級聯(lián)型、頻率抽樣型、快速卷積結(jié)構(gòu)。 3.2 窗口設計法 窗函數(shù)設計法的基本思路是用FIRDF逼近希望的濾波特性。設希望逼近的濾波器的頻率響應函數(shù)為，其單位脈沖響應用表示。為了設計簡單方便，通常選擇為具有片段常數(shù)特性的理想濾波器。因此是無限長非因果序列，不能直接作為FIRDF的單位脈沖響應。窗函數(shù)設計法就是截取為有限長的一段因果序列，并用合適的窗函數(shù)進行加權(quán)做為FIRDF的單位脈沖響應。下面介紹窗函數(shù)設計

5、法的基本設計過程。窗口設計法的主要工作是計算和，但當較為復雜時，就不容易由反傅里葉變換求得。這時一般可用離散傅里葉變換代替連續(xù)傅里葉變換，求得近似值。窗口法的設計步驟如下：（1）通過傅里葉變換換的理想濾波器的單位脈沖響應。 （2）根據(jù)指標選擇窗口形狀、大小和位置。確定窗口類型的主要依據(jù)是過渡帶寬和阻帶最小衰耗的指標。 （3）給定理想頻響由和，加窗得。 （4）檢驗濾波器性能。由求是否在誤差容限之內(nèi)。窗口函數(shù)對理想特性的影響：改變了理想頻響的邊沿特性，形成過渡帶，寬為,等于的主瓣寬度；過渡帶兩旁產(chǎn)生肩峰和余振（帶內(nèi)、帶外起伏），取決的旁瓣，旁瓣多，余振多；旁瓣相對值大，肩峰強，與N無關(guān)；N增加,過

6、渡帶寬減小,肩峰值不變。因主瓣附近(2-2)其中,所以N的改變不能改變主瓣與旁瓣的比例關(guān)系，只能改變的絕對值大小和起伏的密度，當N增加時，幅值變大，頻率軸變密，而最大肩峰永遠為8.95%，這種現(xiàn)象稱為吉布斯（Gibbs）效應。3.3 布萊克曼窗布萊克曼窗的時域形式可表示為： (2-3)它的頻域特性為： (2-4)其中為矩形窗函數(shù)的幅度頻率特性。增加一個二次諧波余弦分量，可進一步降低旁瓣，但主瓣寬度進一步增加，為。加N可減少過渡帶。布萊克曼窗主瓣寬，旁瓣小，頻率識別精度最低，但幅值識別精度最高。四、性能測試與分析4.1語音信號錄制 通過手機錄制一段長度24秒的語音文件，如圖2所示。圖2：語音信號

7、錄制錄制好語音信號后，打開MATLAB軟件平臺，利用函數(shù)audioread對語音信號進行采樣，再調(diào)用函數(shù)sound此時可以聽見錄制的語音。采樣完后再語音信號中加入一個單頻噪聲，單頻的噪聲的頻率可以自己設置。按照加入噪聲后的采樣頻率調(diào)用sound函數(shù)，這時可以明顯的聽見播放的語音信號中有尖銳的單頻嘯叫聲。下面是調(diào)用該語言信號以及加入噪聲的程序：x,fs=audioread('test_yfw.mp3');% 使用audioread讀取音頻文件的采樣率fssound(x,fs); % 按指定的采樣率和每樣本編碼位數(shù)回放N=length(x); % 計算信號x的長度fn=2100;

8、% 單頻噪聲頻率，此參數(shù)可改t=0:1/fs:(N-1)/fs; % 計算時間范圍，樣本數(shù)除以采樣頻率x=x(:,1)'y=x+0.1*sin(fn*2*pi*t);sound (y,fs); %明顯聽出有尖銳的單頻嘯叫聲現(xiàn)在是對加入噪聲前后的語音信號進行頻譜分析，先對原始和加噪后的語音信號進行傅里葉變換，再計算頻譜的頻率范圍和譜線間隔。最后就可以畫出未加入噪聲和加入噪聲后的時域圖和頻域圖。所有未加和加入噪聲的時域圖和頻域圖如圖3。下面是對未加和加入噪聲的頻譜分析的程序：X=abs(fft(x);Y=abs(fft(y); % 對原始信號和加噪信號進行fft變換X=X(1:N/2);Y

9、=Y(1:N/2); % 截取前半部分deltaf=fs/N; % 計算頻譜的譜線間隔f=0:deltaf:fs/2-deltaf; % 計算頻譜頻率范圍 圖3：語音信號的時域圖和頻率圖4.2 濾波器設計本課程設計就是要設計一個濾波器慮掉加入的噪聲，使其恢復原始的語音信號。而設計濾波器的方法有很多，例如：窗函數(shù)法、頻率采樣法、脈沖響應不變法和雙線性變換法。而本課程設計采用的是窗函數(shù)法設計FIR濾波器。而FIR濾波器的設計也有很多方法。在Matlab中，可以利用矩形窗、三角窗、漢寧窗、漢明窗、布萊克曼窗、凱塞窗等設計FIR濾波器。而本次采用的是布萊克曼窗來設計濾波器。在用布萊克曼窗設計濾波器的時

10、候，首先要確定濾波器的性能指標。從六種窗函數(shù)的基本參數(shù)中我們可以得到旁瓣峰值n=-57，過度帶寬,最小阻帶衰減，這就表明在設置這些值時其參數(shù)必須不大于這些值。而其它帶阻濾波器的設計指標則要根據(jù)加入噪聲的頻率來確定。若不能按照這些來設計濾波器則不可能慮掉噪聲。當所有的指標都設置完后，可以用這些數(shù)字來計算上下邊帶的中心頻率和頻率間隔，并計算布萊克曼窗設計該濾波器所需要的階數(shù)和產(chǎn)生幾階的布萊克曼窗。當所有的準備工作完成后就可以調(diào)用自編的函數(shù)計算理想帶阻濾波器的脈沖響應和用窗函數(shù)法計算實際的濾波器的脈沖響應。最后調(diào)用freqz函數(shù)得到濾波器的頻率特性。從畫出的圖中可以清楚的看見濾波器的幅頻和相頻特性。

11、下面是用布萊克曼窗設計濾波器的整個程序：fpd=1800;fsd=2050;fsu=1950;fpu=2000;Rp=1;As=70; % 帶阻濾波器設計指標fcd=(fpd+fsd)/2;fcu=(fpu+fsu)/2;df=min(fsd-fpd),(fpu-fsu); % 計算上下邊帶中心頻率，和頻率間隔wcd=fcd/fs*2*pi;wcu=fcu/fs*2*pi;dw=df/fs*2*pi; % 將Hz為單位的模擬頻率換算為rad為單位的數(shù)字頻率wsd=fsd/fs*2*pi;wsu=fsu/fs*2*pi; M=ceil(10*pi/dw)+1; % 計算布萊克曼窗設計該濾波器時需

12、要的階數(shù)n=0:M-1; % 定義時間范圍w_black=blackman(M); % 產(chǎn)生M階的布萊克曼窗hd_bs=ideal_lp(wcd,M)+ideal_lp(pi,M)-ideal_lp(wcu,M); % 調(diào)用自編函數(shù)計算理想帶阻濾波器的脈沖響應h_bs=w_black'.*hd_bs; % 用窗口法計算實際濾波器脈沖響應db,mag,pha,grd,w=freqz_m(h_bs,1); % 調(diào)用自編函數(shù)計算濾波器的頻率特性圖4：濾波器的頻率特性4.3信號濾波處理對語音信號信號進行濾波處理主要是濾掉加入的噪聲。不同的濾波器利用不同的函數(shù)對語音信號進行濾波，F(xiàn)IR濾波器利用

13、函數(shù)fftfilt對信號進行濾波。對信號進行濾波處理要先利用函數(shù)filter對y進行濾波，然后對y進行傅里葉變換。而畫頻譜時只取前面一半。最后在同一個圖中畫出原始信號的、加入噪聲的語音信號以及濾波后語音信號的頻域圖和時域圖。這樣便于將所有的圖進行對比和分析，而且還可以直觀的觀察該課程設計是否成功。當將設計好的濾波器濾掉噪聲后我們也可以再一次調(diào)用函數(shù)sound，聽此時的聲音是否與原始語音信號基本一樣，若沒有單頻嘯叫聲則說明此次設計是成功的，否則應重新設置指標。下面是對語音信號進行濾波的程序：y_fil=filter(h_bs,1,y); %用設計好的濾波器對y進行濾波Y_fil=fft(y_fi

14、l); %對y進行傅里葉變化Y_fil=Y_fil(1:length(Y_fil)/2); % 計算頻譜取前一半圖5：濾波前后語音信號對比4.4結(jié)果分析要確定本課程設計是否成功就得看原始信號的頻域圖和時域圖與經(jīng)過濾波器后的語音信號的頻域圖和時域圖是否一樣，若一樣則表示該設計是成功的，否則是不成功的。在第一個圖中：第一幅圖和第二幅圖是原始語音信號的時域圖和頻域圖,第三幅圖和第四幅圖是加入頻率為2100的噪聲。從圖中可以看出，第一圖和第三圖相比因為加入噪聲的緣故所以第三圖y軸的幅度要比第一圖要大，但其形狀還是基本沒有改變。而第二圖與第四圖相比較在頻率f=2100時多了一個尖銳脈沖。說明原始語音信號

15、加入噪聲是成功的。在濾波器頻率特性的圖中可以看到：第一個圖是以db為單位的幅頻特性，第二圖是幅頻特性，第三個圖是濾波器的相頻特性，最后一個圖是濾波器的脈沖響應。從圖中可以清楚的了解濾波器的幅頻和相頻特性。在濾波前后信號比較的圖中我們可以得到：原始的語音信號與濾波后的信號的圖基本一樣，只是濾波后的圖在原始信號的基礎上有所延遲。所以用布萊克曼窗設計的濾波器是符合要求的，也就是說該課程設計是成功的。五、結(jié)論與心得為期兩周的數(shù)字信號處理課程設計已經(jīng)結(jié)束了，但在這次設計中我學到了許多的東西。通過這次的設計，不僅加深了我對課本基礎理論知識的理解，而且增強了我的實踐能力，同時更加認識到理論知識和實踐結(jié)合的重

16、要性。首先，更加深入理解了濾波器設計的各個關(guān)鍵環(huán)節(jié)，包括在什么情況下使用哪種方法設計FIR濾波器最好以及在選擇特定的窗函數(shù)進行濾波器的設計時我們應該怎樣確定其性能指標；其次，更加深刻的認識了語音原始信號與加噪后語音信號的波形及頻譜；再次，較大地提高了綜合運用專業(yè)基礎知識及軟件設計能力，在一定程度上對自己的動手能力有很大的幫助。雖然這次課程設計已經(jīng)完成了，但是遇到的困難也是很多的。其中最主要的問題要屬怎樣設置濾波器的指標問題，如果指標的設置有問題那么后續(xù)的工作就不可能得到原始的語音信號。在設置過程中有很多次因為設置的參數(shù)不合適而導致設計的濾波器不能慮出單頻噪聲信號。所以在設計指標問題時一定要結(jié)合

17、布萊克曼本身的特點還要考慮加入噪聲的頻率。其次就是一些函數(shù)的細節(jié)問題。雖然在這次課程設計中遇到很多的困難，但通過自己查找有關(guān)資料以及老師和同學的幫助下都一一解決了，而且在與同學交流的過程中使同學之間的感情更進一步。在此向幫助我的老師及熱心同學表示忠心的感謝！希望今后還能參加更多的課程設計，以鍛煉自己在各個方面的能力，尤其是綜合運用專業(yè)基礎知識和實踐結(jié)合的能力。設計的過程中，我通過查閱大量有關(guān)資料，與同學交流經(jīng)驗和自學，并向老師請教等方式，使我學到了不少的東西，雖然有許多的辛酸，但是看到自己課程設計完成后心中的那份激動是無法用言語來形容的。六、參考文獻1 康華光模擬電子技術(shù)基礎（第五版）M.北京

18、：高等教育出版社，20062 童詩白，華成英，模擬電子技術(shù)基礎（第三版）M.高等教育出版社，20063 何方白，張德民，陽莉，李強.數(shù)字信號處理M.高等教育出版社，20114 王靖,李永全.橢圓濾波器Matlb設計與實現(xiàn)J.現(xiàn)代電子技術(shù),2007七、附錄7.1 布萊克曼窗設計FIR濾波器的MATLAB主程序x,fs=audioread('test_yfw.mp3'); % 使用audioread讀取音頻文件的采樣率fssound(x,fs); % 按指定的采樣率和每樣本編碼位數(shù)回放N=length(x); % 計算信號x的長度fn=2100; % 單頻噪聲頻率，此參數(shù)可改t=0

19、:1/fs:(N-1)/fs; % 計算時間范圍，樣本數(shù)除以采樣頻率x=x(:,1)'y=x+0.1*sin(fn*2*pi*t);sound (y,fs); %明顯聽出有尖銳的單頻嘯叫聲X=abs(fft(x);Y=abs(fft(y); % 對原始信號和加噪信號進行fft變換X=X(1:N/2);Y=Y(1:N/2); % 截取前半部分deltaf=fs/N; % 計算頻譜的譜線間隔f=0:deltaf:fs/2-deltaf; % 計算頻譜頻率范圍figuresubplot(2,2,1);plot(t,x);xlabel('時間（單位：s）');ylabel(&#

20、39;幅度');title('原始語音信號');subplot(2,2,2);plot(f,X);xlabel('頻率（單位：Hz）');ylabel('幅度譜');title('語音信號幅度譜');subplot(2,2,3);plot(t,y);xlabel('時間（單位：s）');ylabel('幅度');title(' 加入單頻干擾后的語音信號');subplot(2,2,4);plot(f,Y);xlabel('頻率（單位：Hz）');ylabel(

21、'幅度譜');title('加入單頻干擾后的語音信號幅度譜');fpd=1800;fsd=2050;fsu=1950;fpu=2000;Rp=1;As=70; % 帶阻濾波器設計指標fcd=(fpd+fsd)/2;fcu=(fpu+fsu)/2;df=min(fsd-fpd),(fpu-fsu); % 計算上下邊帶中心頻率，和頻率間隔wcd=fcd/fs*2*pi;wcu=fcu/fs*2*pi;dw=df/fs*2*pi; % 將Hz為單位的模擬頻率換算為rad為單位的數(shù)字頻率wsd=fsd/fs*2*pi;wsu=fsu/fs*2*pi;M=ceil(10*

22、pi/dw)+1; % 計算布萊克曼窗設計該濾波器時需要的階數(shù)n=0:M-1; % 定義時間范圍w_black=blackman(M); % 產(chǎn)生M階的布萊克曼窗hd_bs=ideal_lp(wcd,M)+ideal_lp(pi,M)-ideal_lp(wcu,M); % 調(diào)用自編函數(shù)計算理想帶阻濾波器的脈沖響應h_bs=w_black'.*hd_bs; % 用窗口法計算實際濾波器脈沖響應db,mag,pha,grd,w=freqz_m(h_bs,1); % 調(diào)用自編函數(shù)計算濾波器的頻率特性figuresubplot(2,2,1);plot(w/pi,db);title('濾波

23、器的db');xlabel('w/pi');ylabel('db');subplot(2,2,2);plot(w/pi,mag);title('濾波器的幅頻特性');xlabel('w/pi');ylabel('mag');subplot(2,2,3);plot(w/pi,pha);title('濾波器的相頻特性');xlabel('w/pi');ylabel('pha');subplot(2,2,4);plot(h_bs);title('濾波器的脈沖

24、響應 ');xlabel('w/pi');ylabel('h_bs');y_fil=filter(h_bs,1,y); %用設計好的濾波器對y進行濾波Y_fil=fft(y_fil); %對y進行傅里葉變化Y_fil=Y_fil(1:length(Y_fil)/2); % 計算頻譜取前一半figuresubplot(3,2,1);plot(t,x);title('未加噪聲的時域圖');xlabel('t');ylabel('x');subplot(3,2,2);plot(f,X);title('未加

25、噪聲的頻域圖');xlabel('f');ylabel('X');subplot(3,2,3);plot(t,y);title('加噪聲后的時域圖');xlabel('t');ylabel('y');subplot(3,2,4);plot(f,Y);title('加噪聲后的頻域圖');xlabel('f');ylabel('Y');subplot(3,2,5);plot(t,y_fil);title('加噪聲后的頻域圖');xlabel(

26、9;f');ylabel('Y');title('濾去噪聲后的時域圖');xlabel('t');ylabel('y_fil');subplot(3,2,6);plot(f,Y_fil);title('濾去噪聲后的頻域圖');xlabel('f');ylabel('Y_fil');7.2 調(diào)用函數(shù)ideal_lp.m的MATLAB源代碼function hd=ideal_lp(wc,M);% Ideal LowPass filter computation% -% hd = ideal_lp(wc,M);% hd = ideal impulse response between 0 to M-1% wc = cutoff frequency in radians% M = length of the ideal filter%alpha = (M-1)/2;n = 0:1:(M-1);m = n -