1、吹管音樂濾波去噪吹管音樂濾波去噪使用使用 TUKEYWIN 設計的設計的 FIR 濾波器濾波器 學生姓名：學生姓名： 指導老師：指導老師： 摘 要 本課程設計主要內容是設計利用窗口設計法選擇 TUKEYWIN 窗設計一個 FIR濾波器，對一段含噪吹管音樂信號進行濾波去噪處理并根據濾波前后的波形和頻譜分析濾波性能。本課程設計仿真平臺為 MATLAB7.0，開發(fā)工具是 M 語言編程。首先下載一段吹管音樂信號，并人為加入一單頻噪聲，然后對信號進行頻譜分析以確定所加噪聲頻率，并設計濾波器進行濾波去噪處理，最后比較濾波前后的波形和頻譜并進行分析。由分析結果可知，濾波器后的語音信號與原始信號基本一致，即設

2、計的 FIR 濾波器能夠去除信號中所加單頻噪聲，達到了設計目的。關鍵詞 濾波去噪；FIR 濾波器；TUKEYWIN 窗；MATLAB1 引言引言本課程設計主要解決在含噪情況下對吹管音樂信號的濾波去噪處理，處理時采用的是利用窗口設計法選擇 Tukeywin 窗設計的 FIR 濾波器。通過課程設計了解 FIR 濾波器設計的原理和步驟，掌握用 Matlab 語言設計濾波器的方法及編程。通過觀察音樂信號濾波前后的時域波形的比較，加深對濾波器作用的理解。1.1 課程設計目的課程設計目的熟悉 Matlab 語言環(huán)境，掌握 Matlab 語言的編程規(guī)則，利用 TUKEYWIN 窗函數設計法來設計符合要求的

3、FIR 濾波器來實現音樂信號的濾波去噪。并繪制濾波前后的時域波形和頻譜圖。根據圖形分析判斷濾波器設計的正確性。通過本次課程設計熟悉利用 TUKEYWIN 窗函數法設計 FIR 濾波器的過程。增強自己獨立解決問題的能力，提高自己的動手能力。加深對理論知識聯系實際問題的理解。為以后的工作奠定堅實的基礎。1.2 課程設計要求課程設計要求（1）濾波器指標必須符合工程實際。（2）設計完后應檢查其頻率響應曲線是否滿足指標。（3）處理結果和分析結論應該一致，而且應符合理論。（4）獨立完成課程設計并按要求編寫課程設計報告書。MATLAB 軟件包括五大通用功能：數值計算功能（Nemeric）；符號運算功能（Sy

4、mbolic）；數據可視化功能（Graphic）；數據圖形文字統一處理功能（Notebook）和建模仿真可視化功能（Simulink）。其中，符號運算功能的實現是通過請求 MAPLE 內核計算并將結果返回到 MATLAB 命令窗口。該軟件有三大特點：一是功能強大；二是界面友善、語言自然；三是開放性強。目前，Mathworks 公司已推出 30 多個應用工具箱。MATLAB 在線性代數、矩陣分析、數值及優(yōu)化、數理統計和隨機信號分析、電路與系統、系統動力學、信號和圖像處理、控制理論分析和系統設計、過程控制、建模和仿真、通信系統、以及財政金融等眾多領域的理論研究和工程設計中得到了廣泛應用。MATLA

5、B 在信號與系統中的應用主要包括符號運算和數值計算仿真分析。由于信號與系統課程的許多內容都是基于公式演算，而 MATLAB 借助符號數學工具箱提供的符號運算功能能基本滿足信號與系統課程的需求。例如，解微分方程、傅里葉正反變換、拉普拉斯正反變換、z 正反變換等。MATLAB 在信號與系統中的另一主要應用是數值計算與仿真分析，主要包括函數波形繪制、函數運算、沖激響應與階躍響應仿真分析、信號的時域分析、信號的頻譜分析、系統的 S 域分析、零極點圖繪制等內容。數值計算仿真分析可以幫助學生更深入理解信號與系統的理論知識，并為將來使用MATLAB 進行信號處理領域的各種分析和實際應用打下基礎1。2 設計原

6、理設計原理窗函數設計法的基本思想是用 FIRDF 逼近希望的濾波特性。設希望逼近的濾波器的頻率響應函數為，其單位脈沖響應用表示。為了設計簡單方便，通常()jwdHe( )dh n選擇為具有片段常數特性的理想濾波器。因此是無限長非因果序列，不()jwdHe( )dh n能直接作為 FIRDF 的單位脈沖響應。窗函數設計法就是截取為有限長的一段因果( )dh n序列，并用合適的窗函數進行加權作為 FIRDF 的單位脈沖響應2。( )dh n2.1 FIR 濾波器濾波器FIR 濾波器是有限長單位脈沖響應數字濾波器的簡稱，FIR 濾波器的最大優(yōu)點是可以實現線性相位濾波。FIR 濾波器滿足了在數字通信和

7、圖像傳輸與處理等應用場合對線性相位的要求，FIR 濾波器是全零點濾波器，硬件和軟件實現結構簡單，不用考慮穩(wěn)定性問題。當幅頻特性指標相同（不考慮相位特性）時，FIR 濾波器的階數比 IIRDF的階數高得多。所以在要求線性相位濾波的應用場合，一般都用 FIR 濾波器。FIR 數字濾波器設計的基本步驟如下：（1）確定指標：在設計一個濾波器之前，必須首先根據工程實際的需要確定濾波器的技術指標。在很多實際應用中，數字濾波器常常被用來實現選頻操作。因此，指標的形式一般在頻域中給出幅度和相位響應。幅度指標主要以兩種方式給出。第一種是絕對指標。它提供對幅度響應函數的要求，一般應用于 FIR 濾波器的設計。第二

8、種指標是相對指標。它以分貝值的形式給出要求。 （2）逼近：確定了技術指標后，就可以建立一個目標的數字濾波器模型。通常用理想的數字濾波器模型。之后，利用數字濾波器的設計方法，設計出一個實際濾波器模型來逼近給定的目標。 （3）性能分析和計算機仿真：上兩步的結果是得到以差分或沖激響應描述的濾波器。根據這個描述就可以分析其頻率特性和相位特性，以驗證設計結果是否滿足指標要求；或者利用計算機仿真實現設計的濾波器，再分析濾波結果來判斷。FIR 濾波器的設計問題實際上是確定能滿足所要求的轉移序列或脈沖響應的常數的問題，設計方法主要有窗函數、頻率響應法和等波紋最佳逼近法等。窗口設計法是一種通過截斷和計權的方法使

9、無限長非因果序列成為有限長脈沖響應序列的設計方法。通常在設計濾波器之前，應該先根據具體的工程應用確定濾波器的技術指標。在大多數實際應用中，數字濾波器常常被用來實現選頻操作，所以指標的形式一般為在頻域中以分貝值給出的相對幅度響應和相位響應。窗口設計法步驟如下：（1）根據過渡帶寬及阻帶衰減要求，選擇窗函數的類型并估計窗口長度 N。窗函數的類型可根據最小阻帶衰減 AS 獨立選擇。（2）根據待求濾波器的理想頻率響應求出理想單位脈沖響應 hd(n)。（3）由性能指標確定窗函數 W(n)和長度 N。（4）求得實際濾波器的單位脈沖響應 h(n)， h(n)即為所設計 FIR 濾波器系數向量b(n)。( )(

10、 ).( )dh nhn W n (2.1)常見的窗函數性能表如下圖2.1所示：圖2.1 常見窗函數性能表名稱濾波器過渡帶寬最小阻帶衰減名稱濾波器過渡帶寬最小阻帶衰減矩形21dBPARZENWIN/M56dB巴特利特25dBFLATTOPWIN108dB漢寧44dBGAUSSWIN/M60dB漢明51dBBARTHANNWIN/M40dB布萊克曼11/M74dBBLACKMANHARRIS109dBBOHMANWIN/MBCHEBWIN113dBNUTTALLWIN108dBTUKEYWIN22dB2.3 Tukeywin 窗窗w =tukeywin(M)在數組 w 中產生 M 點的 Tuke

11、ywin 窗函數。由上面查表可知：Tukeywin/M，最小阻帶衰減 22db。利用可繪制出 Tukeywin 窗函數波形如下：wvtoolN=128;t=tukeywin(N);wvtool(t);圖 2-2 Tukeywin 窗函數波形3 設計步驟設計步驟3.1 設計流程圖設計流程圖根據設計的要求，首先從網上下載一段音樂信號，對音樂信號加入噪聲干擾，再利用 Tukeywin 窗設計合理的 FIR 濾波器。最后用濾波器對干擾后的音樂信號進行濾波去噪。具體設計流程圖如下圖 3.1 所示：開始下載一段吹管音樂信號，命名為“” ，繪制出其時域波形和頻譜圖。 NO YES 圖 3-1 音樂信號濾波去

12、噪流程圖 錄制語音信號錄制語音信號在網上下載一段音樂，再利用視頻轉換器將其轉換成.格式，且為單聲道，3mpwav再將此.格式音樂控制在 3 秒內，以減少設計中的誤差。然后在 Matlab 軟件平臺下，wav利用函數 wavread 對語音信號進行采樣，記住采樣頻率和采樣點數，采集完成后在信號中加入一個單頻噪聲，繪制原音樂信號和加噪后的音樂信號的時域和頻域的波形圖。具體調用如下：x,fs,bits=wavread(g:數字信號meilideshenhua.wav); % 輸入參數為文件的全路徑和文件名，輸出的第一個參數是每個樣本的值，fs 是生成該波形文件時的采樣率，bits 是波形文件每樣本的

13、編碼位數。sound(x,fs,bits); % 按指定的采樣率和每樣本編碼位數回放加入噪聲干擾，繪制干擾后的時域波形和頻譜圖。利用 Tukeywin 窗設計合理的 FIR 濾波器。將干擾后的音樂信號通過自己設計的濾波器，進行濾波去噪。比較濾波前后的時域波形和頻譜圖，并回放音樂信號，驗證是否達到去噪效果。結束是否達到去噪效果？N=length(x); % 計算信號 x 的長度fn=2100; % 單頻噪聲頻率，此參數可改t=0:1/fs:(N-1)/fs; % 計算時間范圍，樣本數除以采樣頻率x=x(:,1); y=x+0.1*sin(fn*2*pi*t); sound(y,fs,bits);

14、 % 應該可以明顯聽出有尖銳的單頻嘯叫聲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; % 計算頻譜頻率范圍得到的波形圖如圖 3.2 所示：圖 3-2 原始音樂信號與加噪后的音樂信號時域、頻域波形圖3.3 濾波器設計濾波器設計在此用窗函數法設計一個帶阻的濾波器。FIR此處給出各濾波器的性能指標： ，As20dB，Ap1 dB，其截止頻率在單頻噪聲干擾1122200020502150

15、2200psspfHzfHzfHzfHz附近。再進行加窗設計，利用 Matlab 中的函數 freqz 畫出各濾波器的頻率響應。具體調用如下：fpd=2000;fsd=2050;fsu=2150;fpu=2200;Rp=1;As=20; % 帶阻濾波器設計指標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 為單位的數字頻率wsd=fsd/fs*2*pi

16、;wsu=fsu/fs*2*pi;M=ceil(2.4*pi/dw)+1; % 計算 TUKEYWIN 窗設計該濾波器時需要的階數n=0:M-1; % 定義時間范圍w_tuk=tukeywin(M); % 產生 M 階的 TUKEYWIN 窗 hd_bs=ideal_lp(wcd,M)+ideal_lp(pi,M)-ideal_lp(wcu,M); % 調用自編函數計算理想帶阻濾波器的脈沖響應h_bs=w_tuk.*hd_bs; % 用窗口法計算實際濾波器脈沖響應db,mag,pha,grd,w=freqz_m(h_bs,1); % 調用自編函數計算濾波器的頻率特性得到的波形圖如圖 3.3 所

17、示：圖 3-3 濾波器波形圖注：濾波器的設計指標，圖中可以明顯地看到利用窗函數所20sdBTUKEYWIN設計的濾波器基本上可以達到濾波效果。 信號濾波處理信號濾波處理用自己設計的各濾波器分別對采集的信號進行濾波，在 Matlab 中，FIR 濾波器利用函數 fftfilt 對信號進行濾波，IIR 濾波器利用函數 filter 對信號進行濾波。在同一張大圖里分別繪制原始信號 x，加噪信號 x，濾波去噪信號 y 的時域波形和頻譜，以便比較和分析。濾波程序如下：y_fil=fftfilt(h_bs,y); % 用設計好的濾波器對 y 進行濾波Y_fil=fft(y_fil);Y_fil=Y_fil

18、(1:N/2); % 計算頻譜取前一半得到的波形圖如圖 3.4 所示：圖 3.4 濾波前后音樂信號的時域波形圖和頻譜圖對比圖注：從圖中可以看出濾波后的信號的幅度有所減小，說明加噪信號通過濾波器之后不可能毫無變化，也就是說參數的設置以及窗函數的選擇都會影響濾波后信號的質量。3.5 結果分析結果分析在 Matlab 中，函數 sound 可以對經過 tukeywin 窗設計的 FIR 濾波器之后的音樂信號進行回放。其調用格式：sound (y_fil,fs,bits)；可以感覺濾波前后的聲音有變化，那尖銳的噪聲明顯減小。但是所設計的濾波器并不是理想的濾波器，并且由于 tukeywin 窗所設計的濾

19、波器本身的性能局限，所以在回放濾波后的音樂的時候還可以聽到一些噪聲。但從濾波器的波形圖中也可以明顯地看到濾波器的設計滿足所定義的指標，從另一個方面也可以斷言濾波器的設計成功。與相同信源不同濾波方法（Nuttallwin 窗）的同學的濾波器相對比，可以看出由于其衰減偏小，所以性能也較差。與不相同信源（拉弦音樂）相同濾波方法的同學的濾波器相對比，可以看出其濾波效果基本相近。4 出現的問題及解決方法出現的問題及解決方法1、 在錄音時，沒有將錄音的屬性改為 wav 的格式，當在 Matlab 軟件平臺上調用時，出現無法調用的提示，最后通過與同學討論發(fā)現，不能直接更改文件后綴，應用轉換器將其轉換為 wa

20、v 格式。2、在開始錄制音樂信號并將其導入 MATLAB 中的時候出現過錯誤，原因是我所下載的音樂信號是雙音頻信號，不符合要求，在老師的指導下我將音樂信號變成了單聲道信號，再次導入的時候 MATLAB 不再報錯。 3、在最后對加噪音樂信號進行濾波并繪圖的時也出現了錯誤，繪制出最后的濾波信號的頻譜圖時縱坐標出現了負值，但是原信號的頻譜圖中卻沒有，回頭檢查了一遍程序想到幅度值只有正值，想到是不是忘記加絕對值符號了，在對濾波的信號進行 fft 變換之后未取絕對值，改正之后再運行程序就不會再出錯。4、在濾波處理時，在加入單頻噪聲部分噪聲頻率設置為 2200，但濾波器參數設置部分最高帶阻為 2200，導

21、致濾波處理后不能濾除噪聲信號，而后來將其改為 2100，才達到了效果。5 結束語結束語本次的課程設計，我的任務是利用 Tukeywin 窗函數設計 FIR 濾波器對音樂信號濾波去噪。在本次課程設計之前，我對 Tukeywin 窗函數完全沒有了解，因此在看到這個題目時，我是一頭霧水。但是通過自己翻閱資料和詢問同學，我掌握了用 Tukeywin 窗函數設計 FIR 濾波器的方法步驟，了解了窗函數的基本設計流程。經過這兩周的課程設計經過為期兩周的課程設計，這次課程設計比較的簡單，根據老師給出的模板，再結合以前學過的知識，并 在同學和老師的幫助下我順利的完成了任務。不同于在教室里上的理論，這次的課程設

22、計需要將我們平時所學習的知識運用到實踐之中，將知識學以致用。因為是以所學理論為基礎，所以在課程設計的過程中，我又重溫濾波器、窗函數等知識，更加熟悉了 MATLAB 的操作。雖說總體上是比較的簡單，但是也會出現一些錯誤，分析并總結這些錯誤主要是語句上的一些用法出現錯誤，以及沒有及時地掌握老師所提醒的部分。但是最后還是在老師的幫助下，解決了這些問題。課程設計是我們運用所學知識，動手實踐的一個很好的機會。它既可以幫助我們加深對所學知識的理解，又能提高我們運用知識，聯系實際，動手實踐的能力。而且在設計過程中可能用到我們沒學過的知識，需要我們去查閱資料獲取相關信息，這又提高了我們查找信息和學習新知識的能

23、力。在實物的調試與檢測過程中，又會遇到許多意想不到的問題，需要我們去分析原因和解決問題。也體會到真正的去獨立地完成一件事情是很困難的，同學以及老師的幫助和提醒是必須的。通過這次課程設計，我拓寬了知識面，鍛煉了實際操作能力，綜合素質也得到了提高，進一步加深了了我們對專業(yè)的認識和激發(fā)了我們對專業(yè)的興趣。雖然課程設計結束了，但是我們的學習還沒結束，對知識的進一步學習還需要繼續(xù)，很開心成功地完成了這次設計。參考文獻參考文獻1 維納.K 英格爾，約翰.G 普羅克斯.數字信號處理（MATLAB 版）.西安：西安交通大學出版社，20072 程佩青.數字信號處理教程.北京：清華大學出版社，20023 張圣勤.

24、MATLAB7.0 實用教程.北京：機械工業(yè)出版社，20084 張志涌.M北京：北京航空航天大學出版社，2003M.北京：科學出版社，2005M.第版.西安：西安電子科技出版社，2002附錄 1：吹管音樂濾波去噪吹管音樂濾波去噪使用使用 TUKEYWIN 窗窗設計的設計的 FIR 濾波器濾波器程序一：%放出原音頻x,fs,bits=wavread(h:數字信號meilideshenhua.wav); % 輸入參數為文件的全路徑和文件名，%輸出的第一個參數是每個樣本的值，%fs 是生成該波形文件時的采樣率，%bits 是波形文件每樣本的編碼位數。sound(x,fs,bits); % 按指定的采

25、樣率和每樣本編碼位數回放程序二：%原音樂信號和加噪后的音樂信號的時域和頻域的波形圖繪制N=length(x); % 計算信號 x 的長度fn=2100; % 單頻噪聲頻率，此參數可改t=0:1/fs:(N-1)/fs; % 計算時間范圍，樣本數除以采樣頻率x=x(:,1); y=x+0.1*sin(fn*2*pi*t); sound(y,fs,bits); % 應該可以明顯聽出有尖銳的單頻嘯叫聲X=abs(fft(x);Y=abs(fft(y); % 對原始信號和加噪信號進行 fft 變換,取幅度譜X=X(1:N/2); Y=Y(1:N/2);% 截取前半部分deltaf=fs/N; % 計算

26、頻譜的譜線間隔f=0:deltaf:fs/2-deltaf; % 計算頻譜頻率范圍subplot(2,2,1);plot(t,x);xlabel(時間(t);ylabel(幅度); title(原始音樂信號);axis(0,2.5,-1.5,1.5);grid on;subplot(2,2,2);plot(f,X);xlabel(頻率(f);ylabel(幅度譜);title(原始音樂信號幅度譜);axis(0,3000,0,3000);grid on;subplot(2,2,3);plot(t,y);xlabel(時間(t);ylabel(幅度);title(加干擾后的音樂信號);axis(

27、0,2.5,-1.5,1.5);grid on;subplot(2,2,4);plot(f,Y);xlabel(頻率(f);ylabel(幅度譜);title(加干擾后的音樂信號幅度譜);axis(0,3000,0,3000);grid on;程序三：%帶阻濾波器設計fpd=2000;fsd=2050;fsu=2150;fpu=2200;Rp=1;As=20; % 帶阻濾波器設計指標fcd=(fpd+fsd)/2;fcu=(fpu+fsu)/2;df=min(fsd-fpd),(fpu-fsu); % 計算上下邊帶中心頻率，和頻率間隔wcd=fcd/fs*2*pi;wcu=fcu/fs*2*p

28、i;dw=df/fs*2*pi; % 將 Hz 為單位的模擬頻率換算為 rad 為單位的數字頻率wsd=fsd/fs*2*pi;wsu=fsu/fs*2*pi;M=ceil(2.4*pi/dw)+1; % 計算 Tukeywin 窗設計該濾波器時需要的階數n=0:M-1; % 定義時間范圍w_tuk=tukeywin(M); % 產生 M 階的 Tukeywin 窗 hd_bs=ideal_lp(wcd,M)+ideal_lp(pi,M)-ideal_lp(wcu,M); % 調用自編函數計算理想帶阻濾波器的脈沖響應h_bs=w_tuk.*hd_bs; % 用窗口法計算實際濾波器脈沖響應db,

29、mag,pha,grd,w=freqz_m(h_bs,1); % 調用自編函數計算濾波器的頻率特性subplot(2,2,1);plot(w/pi,db);title(FIR 濾波器的幅度響應圖);xlabel(w/pi);ylabel(db);line(0,0.5,-1,-1,Color,r,LineWidth,2,LineStyle,-);line(0,0.5,-20,-20,Color,r,LineWidth,2,LineStyle,-);line(wsd/pi,wsd/pi,-30,15,Color,r,LineWidth,2,LineStyle,-);line(wsu/pi,wsu/

30、pi,-30,15,Color,r,LineWidth,2,LineStyle,-);axis(0.05 0.15 -30 15);grid on;subplot(2,2,2);plot(w/pi,mag);title(FIR 濾波器的幅度響應圖);xlabel(w/pi);ylabel(幅度 mag);axis(0 0.3 0 1.5);grid on;subplot(2,2,3);plot(w/pi,pha);title(濾波器相位響應圖);xlabel(w/pi);ylabel(相位 pha);axis(0 1 -4 4);grid on;subplot(2,2,4);stem(n,h_

31、bs);title(濾波器脈沖響應圖);xlabel(n);ylabel(h(n);axis(0 1000 0 0.7);grid on;程序四：%對信號進行濾波y_fil=fftfilt(h_bs,y); % 用設計好的濾波器對 y 進行濾波Y_fil=abs(fft(y_fil);Y_fil=Y_fil(1:N/2); % 計算頻譜取前一半subplot(3,2,1);plot(t,x);title(原始音樂信號時間 x);grid on;xlabel(時間:t);ylabel(幅度);axis(0,2.5,-1.5,1.5);subplot(3,2,2);plot(f,X);title(原始音樂信號幅度譜 X);xlabel(頻率:f);ylabel(幅度);axis(0 2500 0 2500);grid on;subplot(3,2,3);plot(t,y);title(加干擾音樂信號時間 x1);grid on;xlabel(時間 t);ylabel(幅度);axis(0,2.5,-1.5,1.5);subplot(3,2,4);plot(f,Y);title(加干擾后的音樂信號幅度譜 X1);xlabel(頻率 f);ylabel(幅度);ax