音樂(lè)信號(hào)濾波去噪-基于Gaussian的頻率采樣型FIR濾波器設(shè)計(jì)

1、 王杰 音樂(lè)信號(hào)濾波去噪-使用GAUSS窗設(shè)計(jì)的頻率采樣型FIR濾波器 第24頁(yè) 共23頁(yè)音樂(lè)信號(hào)濾波去噪基于Gausswin設(shè)計(jì)的FIR濾波器學(xué)生姓名：王杰 指導(dǎo)老師：高明摘 要 本次課程設(shè)計(jì)主要內(nèi)容是設(shè)計(jì)利用Gausssian窗口函數(shù)設(shè)計(jì)法設(shè)計(jì)一個(gè)FIR濾波器，將一段加入噪聲干擾的音樂(lè)信號(hào)進(jìn)行濾波去噪處理并根據(jù)濾波前后的波形和頻譜分析濾波性能。本課程設(shè)計(jì)仿真平臺(tái)為MATLAB7.0，開(kāi)發(fā)工具是M語(yǔ)言編程。首先在因特網(wǎng)上下載工具一段音樂(lè)信號(hào)，并加入預(yù)設(shè)好頻率的單頻噪聲，然后對(duì)信號(hào)進(jìn)行頻譜分析得到加噪前后的信號(hào)對(duì)比圖，然后設(shè)計(jì)濾波器進(jìn)行濾波去噪處理，最后比較濾波前后的波形和頻譜并進(jìn)行分析。由分

2、析結(jié)果得到，濾波器后的語(yǔ)音信號(hào)與原始信號(hào)基本一致，即設(shè)計(jì)的FIR濾波器能夠去除信號(hào)中所加單頻噪聲，達(dá)到了設(shè)計(jì)目的。關(guān)鍵詞 MATLAB；濾波去噪；FIR濾波器；Gauss窗1引 言 本次課程設(shè)計(jì)主要是將一段音樂(lè)信號(hào)加入噪聲，然后用某種函數(shù)法設(shè)計(jì)出的FIR濾波器對(duì)加入噪聲后的音樂(lè)信號(hào)進(jìn)行濾波去噪處理，處理時(shí)采用的是利用窗口設(shè)計(jì)法選擇GAUSSWIN 設(shè)計(jì)的FIR 濾波器， 通過(guò)課程設(shè)計(jì)了解FIR 濾波器設(shè)計(jì)的原理和步驟，掌握用Matlab語(yǔ)言設(shè)計(jì)濾波器的方法，觀察音樂(lè)信號(hào)濾波前后的時(shí)域波形的比較，加深對(duì)濾波器作用的理解。1.1 課程設(shè)計(jì)的目的數(shù)字信號(hào)處理（Digital Signal Procc

3、essing，簡(jiǎn)稱DSP）是一門涉及許多學(xué)科而廣泛應(yīng)用于許多領(lǐng)域的新型學(xué)科。20世紀(jì)60年代后，隨著計(jì)算機(jī)和信息技術(shù)的飛速發(fā)展，數(shù)字信號(hào)處理技術(shù)應(yīng)運(yùn)而生并并得到迅速的發(fā)展。在過(guò)去的二十多年時(shí)間里，數(shù)字信號(hào)處理已經(jīng)在通信領(lǐng)域得到極為廣泛的應(yīng)用。數(shù)字信號(hào)處理是利用計(jì)算機(jī)或?qū)Ｓ锰幚碓O(shè)備，以數(shù)字形式對(duì)信號(hào)進(jìn)行采集、變換、濾波、估值、增強(qiáng)、壓縮、識(shí)別等處理，以得到符合人們的需要的信號(hào)形式。在本次課程設(shè)計(jì)中，最主要的設(shè)計(jì)是設(shè)計(jì)FIR 濾波器， FIR 濾波器的設(shè)計(jì)方法主要分為兩類，第一類是基于逼近理想濾波器特性的方法，包括窗函數(shù)法、頻域采樣法和等波紋最佳逼近法；第二類是最優(yōu)設(shè)計(jì)法。本次的課程設(shè)計(jì)主要采用的

4、是第一類設(shè)計(jì)方法，是利用GAUSSWIN法設(shè)計(jì)FIR濾波器對(duì)一段音樂(lè)信號(hào)進(jìn)行濾波去噪，通過(guò)這一過(guò)程，對(duì)濾波前后波形進(jìn)行對(duì)比分析得到結(jié)論。此課程設(shè)計(jì)比較簡(jiǎn)單，主要是將書(shū)本中的知識(shí)應(yīng)用到現(xiàn)實(shí)中，并且根據(jù)自己對(duì)設(shè)計(jì)題目的理解，運(yùn)用軟件編寫(xiě)出程序?qū)崿F(xiàn)這一設(shè)計(jì)，也是我們對(duì)數(shù)字信號(hào)處理的原理進(jìn)行驗(yàn)證的一個(gè)過(guò)程。對(duì)此，也可以加深我們對(duì)所學(xué)知識(shí)的理解，培養(yǎng)我們的動(dòng)手能力。1.2 課程設(shè)計(jì)的要求（1）濾波器指標(biāo)必須符合工程實(shí)際。（2）設(shè)計(jì)完后應(yīng)檢查其頻率響應(yīng)曲線是否滿足指標(biāo)。（3）處理結(jié)果和分析結(jié)論應(yīng)該一致，而且應(yīng)符合理論。（4）獨(dú)立完成課程設(shè)計(jì)并按要求編寫(xiě)課程設(shè)計(jì)報(bào)告書(shū)。1.3 設(shè)計(jì)平臺(tái)MATLABMATLA

5、B是由美國(guó)Math Works公司20世紀(jì)80年代中期推出的數(shù)學(xué)軟件。MATLAB是“Matric Laboratory”的縮寫(xiě)，意及“矩陣實(shí)驗(yàn)室”，優(yōu)秀的數(shù)值計(jì)算能力和卓越的數(shù)據(jù)可視化能力使其很快在數(shù)學(xué)軟件中脫穎而出。Matlab已經(jīng)發(fā)展成為多學(xué)科、多種工作平臺(tái)的功能強(qiáng)大的大型軟件。在歐美的高校和研究機(jī)構(gòu)中，MATLAB是一種非常流行的計(jì)算機(jī)語(yǔ)言，許多重要的學(xué)術(shù)刊物上發(fā)表的論文均是用MATLAB來(lái)分析計(jì)算以及繪制出各種圖形。MATLAB是一完整的并可擴(kuò)展的計(jì)算機(jī)環(huán)境，是一種進(jìn)行科學(xué)和工程計(jì)算的交互式程序語(yǔ)言。它的基本數(shù)據(jù)單元是不需要指定維數(shù)的矩陣，它可直接用于表達(dá)數(shù)學(xué)的算式和技術(shù)概念，而普通

6、的高級(jí)語(yǔ)言只能對(duì)一個(gè)個(gè)具體的數(shù)據(jù)單元進(jìn)行操作。因此，解決同樣的數(shù)值計(jì)算問(wèn)題，使用MATLAB 要比使用Basic、Fortran 和C 語(yǔ)言等提高效率許多倍。許多人贊譽(yù)它為萬(wàn)能的數(shù)學(xué)“演算紙”。MATLAB 采用開(kāi)放式的環(huán)境，你可以讀到它的算法，并能改變當(dāng)前的函數(shù)或增添你自己編寫(xiě)的函數(shù)MATLAB 包含的內(nèi)容非常豐富，功能強(qiáng)大，可以概括為以下幾個(gè)方面：(1)可以在多種操作系統(tǒng)下運(yùn)行，如DOS、Windows 95/98/2000/2000/NT、CompaqAlpha、LinuxSun Solaris 等。(2)有超過(guò)500 種的數(shù)學(xué)、統(tǒng)計(jì)、科學(xué)及工程方面的函數(shù)，使用簡(jiǎn)單快捷，并且有很強(qiáng)的用戶

7、自定義函數(shù)的能力。(3)有強(qiáng)大的圖形繪制和可視化功能， 可以進(jìn)行視覺(jué)數(shù)據(jù)處理和分析，進(jìn)行圖形、圖像的顯示及編輯，能夠繪制二維、三維圖形，使用戶可以制作高質(zhì)量的圖形，從而寫(xiě)出圖文并茂的文章。(4)有從外部文件及外部硬件設(shè)備讀入數(shù)據(jù)的能力。(5)有豐富的工具箱toolbox。各個(gè)領(lǐng)域的專家學(xué)者將眾多學(xué)科領(lǐng)域中常用的算法編寫(xiě)為一個(gè)個(gè)子程序，即m 文件，這些m 文件包含在一個(gè)個(gè)工具箱中。其工具箱可以分為兩大類，即功能性工具箱和科學(xué)性工具箱。功能性工具箱主要用來(lái)擴(kuò)充MATLAB的符號(hào)計(jì)算、圖形可視化、建模仿真、文字處理等功能以及與硬件實(shí)時(shí)交互的功能。學(xué)科性工具箱是按學(xué)科領(lǐng)域來(lái)分類的，如信號(hào)處理、控制、通

8、信、神經(jīng)網(wǎng)絡(luò)圖像處理、系統(tǒng)辨識(shí)、魯棒控制、模糊邏輯、小波等工具箱。MATLAB 中的信號(hào)處理工具箱內(nèi)容豐富，使用簡(jiǎn)便。在數(shù)字信號(hào)處理中常用的算法，如FFT，卷積，相關(guān)，濾波器設(shè)計(jì)，參數(shù)模型等，幾乎都只用一條語(yǔ)句即可以調(diào)用。數(shù)字信號(hào)處理所常用的函數(shù)有波形的產(chǎn)生、濾波器的分析和設(shè)計(jì)、傅里葉變換、Z 變換等。2 設(shè)計(jì)原理2.1 FIR濾波器濾波器根據(jù)其沖激響應(yīng)函數(shù)的時(shí)域特性，可分為2 種，即無(wú)限長(zhǎng)沖激響應(yīng)（IIR）濾波器和有限長(zhǎng)沖激響應(yīng)（FIR）濾波器。FIR 和IIR 的濾波原理都是進(jìn)行卷積，就是對(duì)輸入信號(hào)進(jìn)行某種計(jì)算。FIR 用處就在于對(duì)數(shù)字信號(hào)進(jìn)行必要的處理，得到所需的輸出信號(hào)。FIR 系統(tǒng)有

9、自己突出的優(yōu)點(diǎn)：系統(tǒng)總是穩(wěn)定的；易實(shí)現(xiàn)線性相位；允許設(shè)計(jì)多通帶（或多阻帶）濾波器，后兩項(xiàng)都IIR系統(tǒng)不易實(shí)現(xiàn)的。FIR 數(shù)字濾波器的設(shè)計(jì)方法有多種，如窗函數(shù)設(shè)計(jì)法、頻率采樣法和Chebyshev逼近法等。隨著Matlab 軟件尤其是Matlab的信號(hào)處理工作箱的不斷完善，不僅數(shù)字濾波器的計(jì)算機(jī)輔助設(shè)計(jì)有了可能，而且還可以使設(shè)計(jì)達(dá)到最優(yōu)化。FIR 數(shù)字濾波器設(shè)計(jì)的基本步驟如下：（1）確定技術(shù)指標(biāo)在設(shè)計(jì)一個(gè)濾波器之前，必須首先根據(jù)工程實(shí)際的需要確定濾波器的技術(shù)指標(biāo)。在很多實(shí)際應(yīng)用中，數(shù)字濾波器常被用來(lái)實(shí)現(xiàn)選頻操作。因此，指標(biāo)的形式一般在頻域中給出幅度和相位響應(yīng)。幅度指標(biāo)主要以2種方式給出。第一種是

10、絕對(duì)指標(biāo)。他提供對(duì)幅度響應(yīng)函數(shù)的要求，一般應(yīng)用于FIR濾波器的設(shè)計(jì)。第二種指標(biāo)是相對(duì)指標(biāo)。他以分貝值的形式給出要求。本文中濾波器的設(shè)計(jì)就以線性相位FIR 濾波器的設(shè)計(jì)為例。（2）逼近確定了技術(shù)指標(biāo)后，就可以建立一個(gè)目標(biāo)的數(shù)字濾波器模型（通常采用理想的數(shù)字濾波器模型）。之后，利用數(shù)字濾波器的設(shè)計(jì)方法（窗函數(shù)法、頻率采樣法等） ，設(shè)計(jì)出一個(gè)實(shí)際濾波器模型來(lái)逼近給定的目標(biāo)。（3）性能分析和計(jì)算機(jī)仿真上兩步的結(jié)果是得到以差分或系統(tǒng)函數(shù)或沖激響應(yīng)描述的濾波器。根據(jù)這個(gè)描述就可以分析其頻率特性和相位特性，以驗(yàn)證設(shè)計(jì)結(jié)果是否滿足指標(biāo)要求；或者利用計(jì)算機(jī)仿真實(shí)現(xiàn)設(shè)計(jì)的濾波器，再分析濾波結(jié)果來(lái)判斷。2.2 窗口

11、設(shè)計(jì)法窗口法設(shè)計(jì)的基本想法是要選取某一種合適的理想頻率選擇性濾波器（這種濾波器總是有一個(gè)非因果，無(wú)限長(zhǎng)的脈沖響應(yīng)） ，然后將它的脈沖響應(yīng)截?cái)啵ɑ蚣哟埃┮缘玫揭粋€(gè)線性相位和因果的FIR 濾波器。因此，這種方法的重點(diǎn)在于選擇某種恰當(dāng)?shù)拇昂瘮?shù)和一個(gè)合適的理想濾波器?，F(xiàn)用 代表一理想頻率選擇性濾波器，它在整個(gè)通帶內(nèi)有單位幅度增益和線性相位特性，而阻帶內(nèi)有零響應(yīng)。用窗口設(shè)計(jì)法基本步驟如下1：(1) 給定要求的理想頻率響應(yīng),一般給定分段常數(shù)的理想頻率特性。(2) 由于是在時(shí)域設(shè)計(jì)故必須求出 (2-1)(3) 由于 是無(wú)限時(shí)長(zhǎng)的，故要用一個(gè)有限時(shí)長(zhǎng)的“窗函數(shù)”序列 將加以截?cái)?相乘)，窗的點(diǎn)數(shù)是N點(diǎn)。截?cái)嗪?/p>

12、的序列為h(n)： (2-2)窗的點(diǎn)數(shù)N及窗的形狀是兩個(gè)極重要的參數(shù)。(4) 求出加窗后實(shí)際的頻率響應(yīng) (2-3)(5) 檢驗(yàn)是否滿足的要求，不滿足，則需考慮改變窗形狀或改變 窗長(zhǎng)的點(diǎn)數(shù)N，重復(fù)第(3)、(4)兩步，到滿足要求為止。常見(jiàn)的窗函數(shù)性能表如下表2-1所示。表2-1常見(jiàn)的窗函數(shù)性能表名稱濾波器過(guò)渡帶寬最小阻帶衰減名稱濾波器過(guò)渡帶寬最小阻帶衰減矩形1.8/M21dBPARZENWIN6.6/M56db巴特利特6.1/M25dBFLATTOPWIN19.6/M108db漢寧6.2/M44dBGAUSSWIN5.8/M30db漢明6.6/M51dBBARTHANNWIN3.6/M21db布

13、萊克曼11/M74dBBLACKMANHARRIS16.1/M109dbBOHMANWIN10/M48dbCHEBWIN15.2/M113dbNUTTALLWIN15.4/M108dbTUKEYWIN2.4/M22db2.3 GAUSSIAN窗高斯窗是一種指數(shù)窗。主瓣較寬，故而頻率分辨力低；無(wú)負(fù)的旁瓣，第一旁瓣衰減達(dá)一55dB。常被用來(lái)截短一些非周期信號(hào)，如指數(shù)衰減信號(hào)等。對(duì)于隨時(shí)間按指數(shù)衰減的函數(shù)，可采用指數(shù)窗來(lái)提高信噪比。Gausswin的時(shí)域表達(dá)式可以表示為： (2-4)其中k與 的取值范圍必須滿足22.4濾波器結(jié)構(gòu)本次課程設(shè)計(jì)的濾波器采用的是頻率采樣型結(jié)構(gòu)。(1) 頻率采樣型結(jié)構(gòu)的導(dǎo)入

14、若FIR DF的沖激響應(yīng)為有限長(zhǎng)(N點(diǎn))序列h(n)，則有：圖2-1 關(guān)于h(n)的推導(dǎo)所以，對(duì)h(n)可以利用DFT得到H(k),再利用內(nèi)插公式： (2-5)來(lái)表示系統(tǒng)函數(shù)。(2) 頻率采樣型濾波器結(jié)構(gòu)由式2-5得到FIR濾波器的另外一種結(jié)構(gòu)：頻率采樣型結(jié)構(gòu)。它是由兩部分級(jí)聯(lián)而成。 (2-6)其中：級(jí)聯(lián)中的第一部分為梳狀濾波器 ，第二部分由N個(gè)諧振器組成的諧振貴。 (2-7)(3) 頻率采樣型結(jié)構(gòu)流圖如圖2-2所示圖2-2 頻率采樣型結(jié)構(gòu)流程圖(4) 頻率采樣型結(jié)構(gòu)特點(diǎn)它的系數(shù)H(k)直接就是在 處的頻率響應(yīng)。因此，控制濾波器的頻率響應(yīng)是很直接的。但所有的相乘系數(shù)及H(k)都是復(fù)數(shù)，應(yīng)將它們

15、先化成二階的實(shí)數(shù)，這樣乘起來(lái)比較復(fù)雜，增加了乘法次數(shù)及儲(chǔ)存量；且所有的諧振器的極點(diǎn)都是在單位圓上，由 決定?？紤]到系數(shù)量化的影響，當(dāng)系數(shù)量化時(shí)，極點(diǎn)會(huì)移動(dòng)，有些極點(diǎn)就不能被梳狀濾波器的零點(diǎn)鎖抵消。(零點(diǎn)由延時(shí)單元絕對(duì)，不受量化的影響)系統(tǒng)就會(huì)變得不穩(wěn)定1。3設(shè)計(jì)步驟3.1設(shè)計(jì)流程圖濾波器設(shè)計(jì)流程圖如圖3-1所示。圖3-1 濾波器設(shè)計(jì)流程圖3.2下載并截取音樂(lè)信號(hào)從網(wǎng)上下載一段音樂(lè)，從中截取一段格式為.wav 的語(yǔ)音信號(hào)，時(shí)間為23秒，并加入噪聲。在MALAB 平臺(tái)上，觀察原始語(yǔ)音信號(hào)與加入噪聲后的時(shí)域和頻譜圖。源程序如下所示：x,fs,bits=wavread('G:數(shù)字信號(hào)處理課程設(shè)

16、計(jì)_rEfunk.wav');%讀取音樂(lè)信號(hào)的數(shù)據(jù)sound(x,fs,bits);N=length(x); %計(jì)算信號(hào)x的長(zhǎng)度t=0:1/fs:(N-1)/fs; %計(jì)算時(shí)間范圍，樣本除以采樣頻率x=x'fn=2436; %單頻噪聲功率y=x+0.1*sin(fn*2*pi*t); %加單頻噪聲X=abs(fft(x);Y=abs(fft(y); %對(duì)原始信號(hào)和加噪信號(hào)進(jìn)行fft變換X=X(1:N/2);Y=Y(1:N/2); %截取前半部分deltaf=fs/N; %計(jì)算頻譜的譜線間隔f=0:deltaf:fs/2-deltaf; %計(jì)算頻譜的頻率范圍原始信號(hào)與含噪信號(hào)的時(shí)

17、域圖和頻域圖如圖3-2所示，由時(shí)頻域分析圖可以看出，音樂(lè)信號(hào)能量主要集中在1000Hz以內(nèi)。圖3-2 原始信號(hào)與含噪信號(hào)的時(shí)頻域分析圖3.3濾波器設(shè)計(jì)截取好原始信號(hào)，接下來(lái)的工作是設(shè)計(jì)一個(gè)Gaussian窗濾波器。利用公式1： RP=-min(db(1:1:ceil(wp1/dww)+1,ceil(wp2/dww)+1:1:end) (3-1)AS=-max(db(wsd/dww+1:wsu/dww+1) (3-2)分別得到通帶波紋及通帶最小阻帶衰減，并輸入相應(yīng)濾波器的頻率參數(shù)，接著用freqz_m求得濾波器的響應(yīng)頻率特性。源程序如下所示:fp1=1836;fs1=2406;fs2=2476;

18、fp2=3036; %帶阻濾波器設(shè)置指標(biāo)Rp=6;As=28;fcd=(fp1+fs1)/2; %計(jì)算上下邊帶中心頻率fcu=(fp2+fp2)/2;df=min(fs1-fp1),(fp2-fs2); wcd=fcd/fs*2*pi;wcu=fcu/fs*2*pi;dw=df/fs*2*pi; % 將Hz為單位的模擬頻率換算為rad為單位的數(shù)字頻率wsd=fs1/fs*2*pi;wsu=fs2/fs*2*pi;M=ceil(5.8*pi/dw)+1; %計(jì)算Gauss窗設(shè)計(jì)該濾波器時(shí)需要的階數(shù)n=0:M-1; % 定義時(shí)間范圍gausswf=gausswin(M); %利用gausswin函

19、數(shù)產(chǎn)生一個(gè)Gaussian窗函數(shù)hd_bs=ideal_lp(wcd,M)+ideal_lp(pi,M)-ideal_lp(wcu,M); %計(jì)算理想阻帶濾波器的脈沖響應(yīng)h_bs=gausswf'.*hd_bs; %窗口法計(jì)算實(shí)際濾波脈沖響應(yīng)db,mag,pha,grd,w=freqz_m(h_bs,1); %計(jì)算濾波器的頻率特性wp1=fp1/fs*2*pi;wp2=fp2/fs*2*pi;dww=2*pi/1000;RP=-min(db(1:1:ceil(wp1/dww)+1,ceil(wp2/dww)+1:1:end) %計(jì)算通帶最大衰減AS=-max(db(wsd/dww+1:

20、wsu/dww+1) %計(jì)算出阻帶最小衰減濾波器幅度響應(yīng)圖、相位響應(yīng)圖、脈沖響應(yīng)圖如圖3-3所示。由圖可知，濾波器的實(shí)際通帶最大衰減高于設(shè)定的通帶最大衰減，而實(shí)際的阻帶最小衰減低于設(shè)定的阻帶最小衰減，性能達(dá)標(biāo)。圖3-3濾波器性能分析圖3.4信號(hào)濾波處理設(shè)計(jì)好濾波器后，我們要對(duì)語(yǔ)音信號(hào)進(jìn)行濾波，對(duì)比濾波前和濾波后的音樂(lè)信號(hào)。源程序如下所示：y_fil=fftfilt(h_bs,y); %用設(shè)計(jì)好的濾波器對(duì)y進(jìn)行濾波Y_fil=abs(fft(y_fil);Y_fil=Y_fil(1:length(Y_fil)/2);% %計(jì)算頻譜取前一半sound (y_fil,fs,bits); %音樂(lè)信號(hào)回

21、放wavwrite(y_fil,fs,'d:funk_filted.wav'); %產(chǎn)生去噪后的音樂(lè)文件原始音樂(lè)信號(hào)、濾波后信號(hào)的時(shí)域、頻域分析圖，縱坐標(biāo)為對(duì)數(shù)坐標(biāo)的去噪前后的幅度頻譜圖對(duì)比圖如圖3-4所。圖3-4 音樂(lè)信號(hào)時(shí)域頻域分析圖3.5結(jié)果分析：通過(guò)觀察上圖，濾波后的語(yǔ)音信號(hào)發(fā)生了衰減，說(shuō)明濾波器起到了濾波作用，同時(shí)通過(guò)頻譜對(duì)比，可以看出濾波器濾掉了一部分頻率范圍內(nèi)的信號(hào)。分別聽(tīng)原始語(yǔ)音和濾波后的語(yǔ)音信號(hào)，發(fā)現(xiàn)濾波后的語(yǔ)音信號(hào)噪聲減小了，同時(shí)原始信號(hào)強(qiáng)度稍有減弱，基本達(dá)到了濾波的效果。3.6濾波器結(jié)構(gòu)設(shè)計(jì)本次課程設(shè)計(jì)采用的是頻率采樣型的濾波器結(jié)構(gòu)，將按照之前步驟設(shè)計(jì)出的

22、直接型濾波器通過(guò)tf2fs函數(shù)轉(zhuǎn)換成頻率采樣型1。源程序如下所示：function C,B,A=tf2fs(h)%直接型到頻率抽樣型的轉(zhuǎn)換%-%C,B,A=tf2fs(h)%C=包含各并聯(lián)部分增益的行向量%B=包含按行排列的分子系數(shù)矩陣%A=包含按行排列的分母系數(shù)矩陣%h=FIR濾波器單位沖激響應(yīng)向量%-N=length(h);H=fft(h,N);r=input('r=');magH=abs(H);phaH=angle(H)'%check even or odd Nif (N=2*floor(N/2) L=N/2-1; A1=1,-r,0;1,r,0; C1=real

23、(H(1),real(H(L+2);else L=(N-1)/2; A1=1,-r,0; C1=real(H(1);endk=1:L'%初始數(shù)組B和數(shù)組AB=zeros(L,2);A=zeros(L,3);%計(jì)算分母系數(shù)A(1:L,1)=1;A(1:L,2)=-2*r*cos(2*pi*k/N);A(1:L,3)=r.2A=A;A1;%計(jì)算分子系數(shù)B(1:L,1)=cos(phaH(2:L+1);B(1:L,2)=-cos(phaH(2:L+1)-(2*pi*k/N);%計(jì)算增益系數(shù)C=2*magH(2:L+1),C1'根據(jù)計(jì)算得出的C(各并聯(lián)部分增益的行向量),B(分子系數(shù)向

24、量),A(分母系數(shù)向量)畫(huà)出頻率采樣型濾波器結(jié)構(gòu)圖，如圖3-4所示，C,B,A參數(shù)如圖3-5所示。圖3-4 頻率采樣型濾波器結(jié)構(gòu)圖圖3-5 濾波器的C，B，A參數(shù)表4出現(xiàn)的問(wèn)題及解決方法1、 下載得到的音樂(lè)文件為.mp3文件，用錄音機(jī)軟件轉(zhuǎn)換格式時(shí)無(wú)法識(shí)別，后來(lái)用Adobe Audition 3.0軟件將下載得到的音樂(lè)文件轉(zhuǎn)換成.wav文件就能順利轉(zhuǎn)換文件格式了。2、 在畫(huà)原始信號(hào)，加噪后的信號(hào)及濾波信號(hào)的時(shí)頻域分析圖時(shí)，發(fā)現(xiàn)畫(huà)出的圖形坐標(biāo)系上下限不一致，導(dǎo)致不能產(chǎn)生直觀的對(duì)比效果。經(jīng)過(guò)參考matlab的help文檔之后得知可以采用axis函數(shù)手動(dòng)確定坐標(biāo)系上下限，直觀地畫(huà)出了對(duì)比圖。3、 在

25、生成加噪之后的音樂(lè)文件時(shí)，按照老師的提示使用wavwrite()函數(shù)生成，但是在括號(hào)里填寫(xiě)好目標(biāo)存儲(chǔ)目錄之后發(fā)現(xiàn)matlab報(bào)錯(cuò)，不能成功的生成。后來(lái)經(jīng)常參考matlab的help文件，得知wavwrite()函數(shù)需要填寫(xiě)源數(shù)據(jù)作為參數(shù)，于是將加噪之后的音頻源數(shù)據(jù)y作為參數(shù)，填寫(xiě)到wavwrite()成功生成了加了噪聲之后的文件。4、 在調(diào)用自編函數(shù)計(jì)算理想阻帶濾波器的脈沖響應(yīng)及計(jì)算濾波器的頻率特性時(shí)，發(fā)現(xiàn)ideal_lp函數(shù)以及freqz_m函數(shù)報(bào)錯(cuò)，仔細(xì)檢查才發(fā)現(xiàn)沒(méi)有將已經(jīng)編譯好的m文件放入matlab的m文件夾內(nèi)，導(dǎo)致無(wú)法調(diào)用自編函數(shù)。將編譯好的函數(shù)m文件復(fù)制進(jìn)目標(biāo)文件夾后，順利地進(jìn)行了

26、計(jì)算。5、 設(shè)計(jì)好濾波器后，對(duì)濾波器的性能進(jìn)行驗(yàn)證時(shí)，不知道該如何檢驗(yàn)。后來(lái)按照老師所給的提示，通過(guò)公式計(jì)算出濾波器的實(shí)際Rp以及As與預(yù)設(shè)的Rp和As進(jìn)行比較并且在性能曲線圖中畫(huà)出作為參考的預(yù)設(shè)的Rp,As與實(shí)際的性能曲線相比較判斷性能是否達(dá)標(biāo)，經(jīng)對(duì)比發(fā)現(xiàn)濾波器符合設(shè)計(jì)標(biāo)準(zhǔn)，順利的完成了濾波器性能的檢驗(yàn)6、 進(jìn)行的最后的濾波器結(jié)構(gòu)圖設(shè)計(jì)時(shí)，發(fā)現(xiàn)由tf2fs函數(shù)計(jì)算得出的C,B,A矩陣長(zhǎng)度過(guò)長(zhǎng)，達(dá)到了幾百個(gè)，將會(huì)使得結(jié)構(gòu)圖的繪制非常困難。經(jīng)過(guò)仔細(xì)調(diào)整通帶和阻帶截止頻率將高斯窗的階數(shù)設(shè)定在了43階，將C,B,A矩陣減少到了22行，簡(jiǎn)化了濾波器的結(jié)構(gòu)圖，使得繪制順利的完成了。5結(jié)束語(yǔ)這次課程設(shè)計(jì)，

27、是運(yùn)用上學(xué)期的所學(xué)數(shù)字信號(hào)處理理論知識(shí)的一次實(shí)踐，加深了對(duì)所學(xué)知識(shí)的理解，所以這次課程設(shè)計(jì)對(duì)我來(lái)說(shuō)是收獲頗多的。本次的課程設(shè)計(jì)，我選到的課題主要內(nèi)容是通過(guò)數(shù)字信號(hào)處理技術(shù)實(shí)際處理信號(hào)，這里就是指對(duì)音樂(lè)信號(hào)的處理。作為存儲(chǔ)于計(jì)算機(jī)中的音樂(lè)信號(hào)，其本身就是離散化了的向量，我們只需將這些離散的量提取出來(lái)，就可以對(duì)其進(jìn)行處理了。本次課程設(shè)計(jì)的主要設(shè)計(jì)平臺(tái)是MATLAB，通過(guò)MATLAB 里調(diào)用幾個(gè)命令函數(shù)，很輕易的在實(shí)際化音樂(lè)與數(shù)字信號(hào)的理論之間搭了一座橋。課題的特色在于它將音樂(lè)看作了一個(gè)向量，于是音樂(lè)數(shù)字化了，則可以完全利用數(shù)字信號(hào)處理的知識(shí)來(lái)解決。我們可以像給一般信號(hào)做頻譜分析一樣，來(lái)給音樂(lè)信號(hào)做

28、頻譜分析，也可比較容易地用數(shù)字濾波器來(lái)對(duì)音樂(lè)進(jìn)行濾波處理。改變參數(shù)，理論結(jié)合實(shí)際，分析各參數(shù)對(duì)圖形的影響，從而加深對(duì)各個(gè)參數(shù)的理解。在完成這次課程設(shè)計(jì)過(guò)程中學(xué)到了許多東西，進(jìn)一步理解了濾波器設(shè)計(jì)方法和各參數(shù)意義，通過(guò)分析信號(hào)時(shí)域和頻域的關(guān)系等，加深了對(duì)濾波性能的理解，而且學(xué)會(huì)了使用Matlab 一些基本函數(shù)，增加了進(jìn)一步學(xué)習(xí)Matlab 軟件的興趣。同時(shí)，通過(guò)本次課程設(shè)計(jì)，鍛煉了我的動(dòng)手能力，和提高了我分析問(wèn)題，解決問(wèn)題的能力。同學(xué)們?cè)谝黄鹚伎紗?wèn)題，通過(guò)親自動(dòng)手設(shè)計(jì)一個(gè)FIR 濾波器，發(fā)現(xiàn)了許多潛在的問(wèn)題，而這些問(wèn)題是在平常被我們所忽略的，甚至認(rèn)為是不成問(wèn)題的。讓我再次感受到我們應(yīng)當(dāng)把所學(xué)知識(shí)

29、和實(shí)踐相結(jié)合，才能夠提高自己的知識(shí)水平，取得更大的發(fā)展。參考文獻(xiàn)1程佩青.數(shù)字信號(hào)處理教程M.北京：清華大學(xué)出版社，2013:307,476,278-283.2yangpan011.幾種常見(jiàn)的窗函數(shù)以及其Matlab實(shí)現(xiàn)EB/OL. 3林愛(ài)英, 谷小青, 鄭寶周. 用頻率采樣法設(shè)計(jì)FIR濾波器J. 現(xiàn)代電子技術(shù), 2010, 33(17):85-87.附錄1：音樂(lè)信號(hào)濾波去噪設(shè)計(jì)源程序清單% 程序名稱： funk_filter.m% 程序功能：采用基于gausswin的窗口設(shè)計(jì)法，設(shè)計(jì)FIR 濾波器對(duì)含噪音樂(lè)信號(hào)進(jìn)行濾波去噪處理。% 程序作者：王杰% 最后修改日期： 2017-2-28x,fs

30、,bits=wavread('G:數(shù)字信號(hào)處理課程設(shè)計(jì)_rEfunk.wav'); %讀取音樂(lè)信號(hào)的數(shù)據(jù)sound(x,fs,bits);N=length(x); %計(jì)算信號(hào)x的長(zhǎng)度t=0:1/fs:(N-1)/fs; %計(jì)算時(shí)間范圍，樣本除以采樣頻率x=x'fn=2436; %單頻噪聲功率y=x+0.1*sin(fn*2*pi*t); %加單頻噪聲X=abs(fft(x);Y=abs(fft(y); %對(duì)原始信號(hào)和加噪信號(hào)進(jìn)行fft變換X=X(1:N/2);Y=Y(1:N/2); %截取前半部分deltaf=fs/N; %計(jì)算頻譜的譜線間隔f=0:deltaf:fs/

31、2-deltaf; %計(jì)算頻譜的頻率范圍figure(1);subplot(2,2,1),plot(t,x);grid on;axis(0,2.5,-2,2);title('原始音樂(lè)信號(hào)');xlabel('時(shí)間（單位：s）');ylabel('幅度');subplot(2,2,2),plot(f,X);grid on;title('音樂(lè)信號(hào)幅度譜圖');xlabel('頻率（單位：Hz）');ylabel('幅度譜');subplot(2,2,3),plot(t,y);grid on;axis(0

32、,2.5,-2,2);title('加入噪聲干擾后的音樂(lè)信號(hào)');xlabel('時(shí)間（單位：s）');ylabel('幅度');subplot(2,2,4),plot(f,Y);grid on;title('加入噪聲干擾后的音樂(lè)信號(hào)幅度譜圖');xlabel('頻率（單位：Hz）');ylabel('幅度譜');wavwrite(y,fs,'d:funk_converted.wav');fp1=1836;fs1=2406;fs2=2476;fp2=3036; %帶阻濾波器設(shè)置指標(biāo)R

33、p=6;As=28;fcd=(fp1+fs1)/2; %計(jì)算上下邊帶中心頻率fcu=(fp2+fp2)/2;df=min(fs1-fp1),(fp2-fs2); wcd=fcd/fs*2*pi;wcu=fcu/fs*2*pi;dw=df/fs*2*pi; %將Hz為單位的模擬頻率換算為rad為單位的數(shù)字頻率wsd=fs1/fs*2*pi;wsu=fs2/fs*2*pi;M=ceil(5.8*pi/dw)+1; %計(jì)算Gauss窗設(shè)計(jì)該濾波器時(shí)需要的階數(shù)n=0:M-1; % 定義時(shí)間范圍gausswf=gausswin(M); %利用gausswin函數(shù)產(chǎn)生一個(gè)Gaussian窗函數(shù)hd_bs=

34、ideal_lp(wcd,M)+ideal_lp(pi,M)-ideal_lp(wcu,M); %計(jì)算理想阻帶濾波器的脈沖響應(yīng)h_bs=gausswf'.*hd_bs; %窗口法計(jì)算實(shí)際濾波脈沖響應(yīng)db,mag,pha,grd,w=freqz_m(h_bs,1); %計(jì)算濾波器的頻率特性wp1=fp1/fs*2*pi;wp2=fp2/fs*2*pi;dww=2*pi/1000;RP=-min(db(1:1:ceil(wp1/dww)+1,ceil(wp2/dww)+1:1:end) %計(jì)算通帶最大衰減AS=-max(db(wsd/dww+1:wsu/dww+1) %計(jì)算出阻帶最小衰減f

35、igure(2);subplot(2,2,1);plot(w/pi,db);xlabel('w/pi');ylabel('db');title('FIR濾波器的幅度響應(yīng)圖');axis(0.3 0.9 -80 5);line(0,1.5,-RP,-RP,'Color','r','LineWidth',2,'LineStyle','-'); line(0,1.5,-AS,-AS,'Color','r','LineWidth'

36、;,2,'LineStyle','-'); line(wp1/pi,wp1/pi,-100,10,'Color','r','LineWidth',2,'LineStyle','-'); line(wp2/pi,wp2/pi,-100,10,'Color','r','LineWidth',2,'LineStyle','-');line(wsd/pi,wsd/pi,-100,10,'Color'

37、;,'r','LineWidth',2,'LineStyle','-'); line(wsu/pi,wsu/pi,-100,10,'Color','r','LineWidth',2,'LineStyle','-'); grid on; subplot(2,2,2);plot(w/pi,mag);xlabel('w/pi');ylabel('幅度mag');title('FIR濾波器的幅度響應(yīng)圖');axi

38、s(0 1 -0.2 1.2);grid on; subplot(2,2,3);plot(w/pi,pha);xlabel('w/pi');ylabel('相位pha');title('濾波器相位響應(yīng)圖');axis(0 1 -4 4);grid on;subplot(2,2,4);stem(n,h_bs);xlabel('n');ylabel('h(n)');title('濾波器脈沖響應(yīng)圖');grid on;y_fil=fftfilt(h_bs,y); %用設(shè)計(jì)好的濾波器對(duì)y進(jìn)行濾波 Y_fil

39、=abs(fft(y_fil);Y_fil=Y_fil(1:length(Y_fil)/2); %計(jì)算頻譜取前一半 sound (y_fil,fs,bits); %音樂(lè)信號(hào)回放wavwrite(y_fil,fs,'d:funk_filted.wav'); %產(chǎn)生去噪后的音樂(lè)文件figure(3); %作圖subplot(4,2,1);plot(t,x);grid on; title('原始音樂(lè)信號(hào)');xlabel('時(shí)間（t）');ylabel('幅度');axis(0,2.5,-2,2);subplot(4,2,2);plot(f,X);grid on; axis(0,3000,0,1500);title('原始音樂(lè)信號(hào)幅度譜');xlabel('頻率（f）');ylabel('幅度譜'); subplot(4,2,5);plot(t,y_fil);grid on; axis(0,2.5,-2,2);title('濾波后的音樂(lè)信號(hào)');xlabel('時(shí)間（t）');ylabel('幅度');subplot(4,2,6);plot