語音信號的數(shù)字濾波處理(五)

1、 語音信號的數(shù)字濾波處理目 錄1 緒論1 1.1 課程設(shè)計背景11.2 課程設(shè)計目的12 課程設(shè)計預(yù)習(xí)2 2.1 卷積運算的演示22.2 採樣定理的演示43 程式設(shè)計的基本原理6 3.1 IIR濾波器的設(shè)計原理73.2 FIR濾波器的設(shè)計原理74 設(shè)計程式的調(diào)試和運行結(jié)果84.1 濾波器的設(shè)計9 4.1.1設(shè)計方案10 4.1.2設(shè)計步驟11 4.2 IIR濾波器的設(shè)計12 4.2.1巴特沃斯低通程式的設(shè)計13 4.2.2巴特沃斯高通程式的設(shè)計14 4.2.3巴特沃斯帶通程式的設(shè)計15 4.3 FIR濾波器的設(shè)計16 4.3.1 Kaiser窗低通程式的設(shè)計17 4.3.2 Kaiser窗高通

2、程式的設(shè)計19 4.3.3 Kaiser窗帶通程式的設(shè)計20 4.3 回放語音信號225 心得與體會23參考文獻(xiàn)24附錄251 緒論1.1課程設(shè)計背景數(shù)字濾波器是數(shù)字信號處理的基礎(chǔ)，用來對信號進(jìn)行過濾、檢測與參數(shù)估計等處理，在通信、圖像、語音、雷達(dá)等許多領(lǐng)域都有著十分廣泛的應(yīng)用。尤其在圖像處理、數(shù)據(jù)壓縮等方面取得了令人矚目的進(jìn)展和成就。由此可見，數(shù)字濾波器的設(shè)計尤為重要。目前數(shù)字濾波器的設(shè)計有許多現(xiàn)成的高級語言設(shè)計程式，但他們都存在設(shè)計效率較低，不具有可視圖形，不便於修改參數(shù)等缺點，而Matlab為數(shù)字濾波的研究和應(yīng)用提供了一個直觀、高效、便捷的利器。他以矩陣運算為基礎(chǔ)，把計算、可視化、程式設(shè)

3、計融合到了一個互動式的工作環(huán)境中。尤其是Matlab工具箱使各個領(lǐng)域的研究人員可以直觀方便地進(jìn)行科學(xué)研究與工程應(yīng)用。其中的信號處理工具箱、圖像處理工具箱、小波工具箱等更是為數(shù)字濾波研究的蓬勃發(fā)展提供了可能。本次課程設(shè)計主要介紹基於Matlab環(huán)境的IIR數(shù)字濾波器和FIR數(shù)字濾波器的設(shè)計，然後將設(shè)計結(jié)果應(yīng)用於語音信號的處理，針對處理結(jié)果進(jìn)一步分析了濾波器性能的。1.2課程設(shè)計目的 (1)掌握數(shù)字信號處理的基本概念，基本理論和基本方法。(2)熟悉離散信號和系統(tǒng)的時域特性。(3)掌握序列快速傅裏葉變換方法。(4)學(xué)會MATLAB的使用，掌握MATLAB的程式設(shè)計方法。(5)掌握利用MATLAB對語

4、音信號進(jìn)行頻譜分析。(6)掌握濾波器的網(wǎng)路結(jié)構(gòu)。(7)掌握MATLAB設(shè)計IIR、FIR數(shù)字濾波器的方法和對信號進(jìn)行2 課程設(shè)計預(yù)習(xí)2.1卷積運算的演示1. 線性卷積序列x1(n)=2 0 1 2 5 7 0 5 0 2 1 3，序列x2(n)= 2 0 1 1 9 9。動態(tài)演示兩個序列進(jìn)行線性卷積x1(n)x2(n)的翻轉(zhuǎn)、移位、乘積、求和的過程。其中翻轉(zhuǎn)採用fliplr1，程式見附錄A。圖2.1 線性卷積圖示圖2.2 動態(tài)演示N=2和N=72.1.2 迴圈卷積序列x1(n)= 2 0 1 2 5 7 0 5 0 2 1 3，序列x2(n)= 2 0 1 1 9 9 ，N=12。動態(tài)演示兩個

5、序列進(jìn)行迴圈卷積x1(n)x2 (n)的翻轉(zhuǎn)、移位、乘積、求和的過程,程式見附錄B。圖2.3 迴圈卷積圖示當(dāng)N<12時，迴圈卷積發(fā)生混疊；當(dāng)N20時，迴圈卷積等於線性卷積。2.2 採樣定理的演示信號進(jìn)行採樣，A為幅度因數(shù)，a為衰減因數(shù)，為模擬角頻率，調(diào)試結(jié)果如Error! Reference source not found.所示，程式見附錄C圖2.4 採樣定理圖示3 課程設(shè)計的基本原理3.1 IIR濾波器的設(shè)計原理IIR數(shù)字濾波器的設(shè)計一般是利用目前已經(jīng)很成熟的模擬濾波器的設(shè)計方法來進(jìn)行設(shè)計，通常採用模擬濾波器原型有巴特沃斯濾波器、切比雪夫濾波器、橢圓濾波器等。IIR數(shù)字濾波

6、器單位取樣回應(yīng)是時寬無限的h(n)，即n0n。數(shù)字濾波器按照實現(xiàn)的方法和結(jié)構(gòu)形式分為遞歸型或非遞歸型兩類。遞歸型數(shù)字濾波器的當(dāng)前輸出y(n)是輸入x(n)的當(dāng)前值和以前各輸入值x(n)，x(n1)，.，及以前各輸出值y(n)，y(n1)，.的函數(shù)。一個N階遞歸型數(shù)字濾波器(IIR濾波器)的差分方程如公式2-1所示： 2-1其中，式中的係數(shù)至少有一項不為零。說明必須將延時的輸出序列進(jìn)行回饋。遞歸系統(tǒng)的傳統(tǒng)函數(shù)定義如公式2-2所示： 2-2遞歸系統(tǒng)的傳遞函數(shù)H(z)在Z平面上不僅有零點，而且有極點。非遞歸型數(shù)字濾波器當(dāng)前的輸出值y(n)僅為當(dāng)前及以前的輸入序列的函數(shù)，而與以前的各個輸出值無關(guān)，因此

7、從結(jié)構(gòu)上看非遞歸系統(tǒng)沒有回饋環(huán)路。3.2 FIR濾波器的設(shè)計原理FIR數(shù)字濾波器的單位取樣回應(yīng)是時寬有限的h(n)，即n1nn2 。一個N階的非遞歸型數(shù)字濾波器(FIR濾波器)的差分方程如公式2-3所示： 2-3差分方程式中的係數(shù)等於單位取樣回應(yīng)的序列值h(n)，其系統(tǒng)函數(shù)H(z)可以表示為公式2-4所示的形式： 2-4H(z)是的多項式，因此它的極點只能在Z平面的原點上。FIR濾波器通常採用窗函數(shù)方法來設(shè)計。窗設(shè)計的基本思想是，首先選擇一個適當(dāng)?shù)睦硐脒x頻濾波器（它總是具有一個非因果，無限持續(xù)時間脈衝回應(yīng)），然後街區(qū)（加窗）它的脈衝回應(yīng)得到線性相位和因果FIR濾波器。我們用表示理想的選頻濾波器

8、，它在通帶上具有單位增益和線性相位，在阻帶上具有零回應(yīng)。一個帶寬wc<pi的低通濾波器由下式給定：為了從hd(n)得到一個FIR濾波器，必須同時在兩邊截取hd（n）。而要得到一個因果的線性相位濾波器，它的h(n)長度為N，必須滿足公式2-5： 2-5h(n)可以看做是hd(n)與窗函數(shù)w（n）的乘積，如公式2-6所示： 2-6這種操作叫做加窗，其中根據(jù)w（n）的不同定義，可以得到不同的窗結(jié)構(gòu)。在頻域中，因果FIR濾波器回應(yīng)由和窗回應(yīng)的週期卷積得到，如公式2-7所示: 2-7常用的窗函數(shù)有矩形窗、漢寧（HANNING）窗、海明（矩形）窗、布萊克曼（BLACKMAN）窗、凱澤（KAISER）

9、窗等。四 設(shè)計方案及步驟4.1 程式設(shè)計4.1.1 設(shè)計方案1.要求利用windows下的答錄機(jī)（開始程式附件娛樂答錄機(jī)，檔屬性立即轉(zhuǎn)換8000Hz，8位，單聲道）錄製一段自己的話音，或者採用Windows自帶的聲音檔（默認(rèn)為22050Hz），時間控制在幾秒左右。然後在MATLAB軟體平臺下，利用函數(shù)wavread對語音信號進(jìn)行採樣，記住採樣頻率和採樣點數(shù)。通過wavread函數(shù)的使用，要求理解採樣頻率、採樣位數(shù)等概念。wavread函數(shù)調(diào)用格式：y=wavread(file)，讀取file所規(guī)定的wav檔，返回採樣值放在向量y中。y,fs,nbits=wavread(file)，採樣值放在向

10、量y中，fs表示採樣頻率（Hz），nbits表示採樣位數(shù)。y=wavread(file,N)，讀取前N點的採樣值放在向量y中。y=wavread(file,N1,N2)，讀取從N1點到N2點的採樣值放在向量y中。2語音信號的頻譜分析要求首先畫出語音信號的時域波形；然後對語音信號進(jìn)行頻譜分析，在MATLAB中，可以利用函數(shù)fft對信號進(jìn)行快速付立葉變換，得到信號的頻譜特性；從而加深對頻譜特性的理解。3設(shè)計數(shù)字濾波器和畫出頻率回應(yīng)根據(jù)語音信號的特點給出有關(guān)濾波器的性能指標(biāo)，例如：1）低通濾波器性能指標(biāo)，fp=1000Hz，fc=1200 Hz， As=100dB，Ap=1dB；2）高通濾波器性能指

11、標(biāo)，fc=2800 Hz，fp=3000 Hz，As=100dB，Ap=1dB；3）帶通濾波器性能指標(biāo)，fp1=1200 Hz，fp2=3000 Hz，fc1=1000 Hz，fc2=3200 Hz，As=100dB，Ap=1dB。4.回放語音信號在Matlab中,函數(shù)sound可以對聲音進(jìn)行回放。其調(diào)用格式：sound(x,fs,bits);可以感覺濾波前後的聲音有變化。5.方案設(shè)計流程圖如圖3.5.1所示。（圖為巴特沃思濾波器，Kaiser窗濾波器與巴特沃思類似）圖4.1 方案設(shè)計流程圖4.1.2設(shè)計步驟（1）語音信號的採集Matlab程式見附錄D:運行結(jié)果如圖U/V|X(eiw)|圖4.

12、2 語音信號（2）對語音信號加噪U/VMatlab程式見附錄E：運行結(jié)果如圖|X(eiw)|U/V|X(eiw)|圖4.3 對語音信號加噪4.2 用巴特沃斯設(shè)計IIR濾波器4.2.1巴特沃斯低通濾波器Matlab程式見附錄F巴特沃斯低通濾波器程式所示，運行程式波形圖分析如下：U/V|X(eiw)|圖4.44.2.2 巴特沃斯帶通濾波器Matlab程式見附錄G切比雪夫帶通濾波器程式所示，運行程式波形圖分析如下：U/V|X(eiw)|U/V|X(eiw)|U/V|X(eiw)|圖4.54.2.3 巴特沃斯高通濾波器Matlab程式見附錄H比雪夫高通濾波器程式所示，運行程式波形圖分析如下：|X(ei

13、w)|U/V|X(eiw)|U/V|X(eiw)|U/V圖4.64.3 Kaiser窗設(shè)計FIR濾波器4.3.1 Kaiser高通濾波器調(diào)試程式見附錄IKaiser窗高通濾波器程式所示，運行程式波形圖分析如下：U/V|X(eiw)|U/V|X(eiw)|U/V|X(eiw)|X(eiw)|U/V圖4.74.3.2 Kaiser低通濾波器U/V調(diào)試程式見附錄JKaiser窗低通濾波器程式所示，運行程式波形圖分析如下：|X(eiw)|U/V|X(eiw)|U/V|X(eiw)|U/V|X(eiw)|圖4.84.3.3 Kaiser帶通濾波器調(diào)試程式見附錄KKaiser窗帶通濾波器程式所示，運行程式

14、波形圖分析如下：U/V|X(eiw)|U/V|X(eiw)|U/V|X(eiw)|U/V|X(eiw)|圖4.95 心得與體會在本次課程設(shè)計的過程中，主要遇到的問題是原始音樂信號的選取。有些音樂信號在MATLAB中運行時出錯，開始以為自己的操作問題，在跟同學(xué)商討多次更換音樂信號後得到瞭解決。其次就是在濾波器參數(shù)的選取上，開始時不清楚怎麼選取，又在課本中找了好久才找到。學(xué)習(xí)的過程是相互討論共同進(jìn)步的，多多討論課題中遇到的問題，可以鞏固我們的知識掌握能力，增加熟練運用度。從本次課程設(shè)計的中心來看，課題是希望將數(shù)字信號處理技術(shù)應(yīng)用於某一實際領(lǐng)域，這裏就是指對音樂的處理。作為存儲於電腦中的音樂信號，其

15、本身就是離散化了的向量，我們只需將這些離散的量提取出來，就可以對其進(jìn)行處理了。在這裏，用到了處理數(shù)字信號的強有力工具M(jìn)ATLAB，通過MATLAB裏幾個命令函數(shù)的調(diào)用，很輕易的在實際化音樂與數(shù)字信號的理論之間搭了一座橋。課題的特色在於它將音樂看作了一個向量，於是音樂數(shù)位化了，則可以完全利用數(shù)字信號處理的知識來解決。我們可以像給一般信號做頻譜分析一樣，來給音樂信號做頻譜分析，也可以較容易的用數(shù)字濾波器來對音樂進(jìn)行濾波處理。改變參數(shù)，理論結(jié)合實際，分析各參數(shù)對圖形的影響，從而加深對各個參數(shù)的理解。在完成這次課程設(shè)計過程中學(xué)到了許多東西，進(jìn)一步理解了濾波器設(shè)計方法和各參數(shù)意義，通過分析信號時域和頻域

16、的關(guān)係等，加深了對濾波性能的理解，而且學(xué)會了使用Matlab一些基本函數(shù)，增加了進(jìn)一步學(xué)習(xí)Matlab軟體的興趣。同時，通過本次課程設(shè)計，鍛煉了我的動手能力，和提高了我分析問題，解決問題的能力。課程設(shè)計是我們專業(yè)課程知識綜合應(yīng)用的實踐訓(xùn)練，著是我們邁向社會，從事職業(yè)工作前一個必不少的過程”千里之行始於足下”，通過這次課程設(shè)計，我深深體會到這句千古名言的真正含義我今天認(rèn)真的進(jìn)行課程設(shè)計，學(xué)會腳踏實地邁開這一步，就是為明天能穩(wěn)健地在社會大潮中奔跑打下堅實的基礎(chǔ)在這次的課程設(shè)計中，犯了些錯誤，發(fā)現(xiàn)自己還不夠認(rèn)真不夠細(xì)心以前的知識學(xué)習(xí)的不是很牢固，因為數(shù)字信號處理學(xué)的不是很好，很多相關(guān)知識的運用不是很

17、熟悉，還有MATLAB 的很多相關(guān)指令的運用不熟練，導(dǎo)致課程設(shè)計做起來有點吃力，但是只要努力去認(rèn)真做就會得到你想要的，通過這次DSP課程設(shè)計，不但鞏固了以前學(xué)的知識，還學(xué)習(xí)了新內(nèi)容，收穫頗豐。原來濾去語音中的雜訊我們運用自己學(xué)到的知識能做到，並更加鞏固了MATLAB軟體的使用，在從這次的課程設(shè)計中體會到，我覺得什麼事都要我們真真正正用心去做了，才會使自己更加成長，沒有學(xué)習(xí)就不可能有實踐的機(jī)會，沒有自己的實踐就不會有所突破，希望這次的經(jīng)歷能讓我們在以後的學(xué)習(xí)生活中不斷成長，在以後的學(xué)習(xí)中更能與實際相結(jié)合，為自己今後在社會上能立足打好堅實的基礎(chǔ)。想起大二的實習(xí)和這次的課程設(shè)計，至今我仍感慨頗多，從

18、理論到實踐，在這段日子裏，可以說得是苦多於甜，但是可以學(xué)到很多很多的東西，同時不僅可以鞏固了以前所學(xué)過的知識，而且學(xué)到了很多在書本上所沒有學(xué)到過的知識。通過這次課程設(shè)計使我懂得了理論與實際相結(jié)合是很重要的，只有理論知識是遠(yuǎn)遠(yuǎn)不夠的，只有把所學(xué)的理論知識與實踐相結(jié)合起來，從理論中得出結(jié)論，才能真正為社會服務(wù)，從而提高自己的實際動手能力和獨立思考的能力。在設(shè)計的過程中遇到問題，可以說得是困難重重，但可喜的是最終都得到瞭解決。最後在此感謝我們的袁志輝老師，老師嚴(yán)謹(jǐn)細(xì)緻的作風(fēng)一直是我學(xué)習(xí)中的榜樣。這次課程設(shè)計的每個實驗細(xì)節(jié)和每個數(shù)據(jù)，都離不開老師您的細(xì)心指導(dǎo)。而您開朗的個性和寬容的態(tài)度，幫助我能夠很順

19、利的完成了這次課程設(shè)計。同時感謝對我?guī)椭^的同學(xué)們，謝謝你們對我的幫助和支持，讓我感受到同學(xué)的友誼參考文獻(xiàn)1 丁玉美.數(shù)字信號處理M.西安：西安電子科技大學(xué)出版社，2001.2 王創(chuàng)新，文弅.數(shù)字信號處理試驗指導(dǎo)書M.長沙：長沙理工大學(xué)出版 社,2010.3 陳懷琛.數(shù)字信號處理及其MATLAB實現(xiàn)M.北京:電子工業(yè)出版社,2005.4 陳懷琛.MATLAB及在電子資訊課程中的應(yīng)用M，北京:電子工業(yè)出版 社,2005.5 A.V.奧本海姆,R.W.謝弗.數(shù)字信號處理M.北京：科學(xué)出版社,2010.6 胡廣書.數(shù)字信號處理理論、演算法與實現(xiàn)M.北京：電子工業(yè)出版 社,2010.附錄附錄A 線性卷

20、積程式%線性卷積xn=2,0,1,2,5,7,0,5,0,2,1,3;nx=0:11;hn=2,0,1,1,9,9;nh=0:5;nys=nh(1)+nx(1);nyf=nh(end)+nx(end);ny=nys:nyf;yn=conv(hn,xn);subplot(221)stem(ny,yn);xlabel('n');ylabel('幅度')xlabel('ny');ylabel('yn'); %卷積xm=xn;subplot(222)stem(nh,hn);xlabel('n');ylabel('幅

21、度')xlabel('nh');ylabel('hn');subplot(223)hm=fliplr(hn);mh=-fliplr(nh);stem(mh,hm); xlabel('n');ylabel('幅度')xlabel('mh');ylabel('h(-m)'); %翻轉(zhuǎn)for j=0:11 xmm=xm,zeros(1,6); hmm=zeros(1,j) hm zeros(1,12-j); d=sum(xmm.*hmm)'); subplot(224) stem(j,d)

22、;xlabel('n');ylabel('幅度') hold onend %乘積與求和for i=0:17 figure()stem(xn,'r',':');hold on;stem(mh+i,hn,'g')end%移位for j=0:17 xmm=xm,zeros(1,6); hmm=zeros(1,j) hm zeros(1,12-j); d=sum(xmm.*hmm)'); subplot(224) stem(j,d);xlabel('n');ylabel('幅度')

23、hold onend 附錄B迴圈卷積程式%迴圈卷積與其動態(tài)演示n=-8:15;x=2,0,1,1,9,9;x2=zeros(1,8) x zeros(1,10);subplot(321);stem(n,x2);title('初始序列 x2')xlabel('n');ylabel('幅度')axis(-8,20,0,10)xx=x,zeros(1,2);x21 =xx, xx, xx;subplot(324); stem(n,x21); title('週期延伸')xlabel('n');ylabel('幅度&

24、#39;)axis(-8,20,0,10)x22 = x21(3:24), xx(1:2);subplot(3,2,2); stem(n,x22); title('迴圈移位')xlabel('n');ylabel('幅度')axis(-8,20,0,10)N=12;x1=2,0,0,8,5,7,0,5,0,2,1,5;x11=fliplr(x1);for i=0:11; x12=x11(12-i:12),x11(1:11-i); x23=x,zeros(1,6); y=sum(x12.*x23); subplot(323) stem(i,y);t

25、itle('N=12的迴圈卷積') xlabel('n');ylabel('幅度') hold onend subplot(325)f=0:16;yn=conv(x1,x);stem(f,yn);title('線性卷積')xlabel('n');ylabel('幅度')NN=17;xx1=x1,zeros(1,5);xx11=fliplr(xx1);xx23=x,zeros(1,11);for j=0:16; xx12=xx11(17-j:17),xx11(1:16-j); yy=sum(xx12.

26、*xx23); subplot(326) stem(j,yy);title('NN=17的迴圈卷積') xlabel('n');ylabel('幅度') hold onend 附錄C 採樣定理程式n=0:50; A=10;a=10*sqrt(2)*pi;w0=10*sqrt(2)*pi;t0=0.001;t1=0.020;t2=0.005;x=A*exp(-a*n*t0).*sin(w0*n*t0); y1=A*exp(-a*n*t1).*sin(w0*n*t1); y2=A*exp(-a*n*t2).*sin(w0*n*t2); subplot

27、(321);stem(n,x,'r'),grid on title('離散時間信號')xlabel('n');ylabel('幅度')subplot(322);plot(n,x,'g'),grid on title('連續(xù)時間信號')xlabel('t');ylabel('幅度')subplot(323);stem(n,y1),grid ontitle('50Hz連續(xù)時間信號')xlabel('n');ylabel('幅度

28、9;)subplot(324);stem(n,y2,'m'),grid ontitle('200Hz理想採樣信號序列');xlabel('n');ylabel('幅度')subplot(325);stem(n,x,'k'),grid ontitle('1000Hz連續(xù)時間信號')xlabel('n');ylabel('幅度') t=0:0.001:1; xx=A*exp(-a*n*t0).*sin(w0*n*t0); X=fft(xx);subplot(326)plo

29、t(abs(X);xlabel('f/Hz');ylabel('|X(ejw)|') grid on；附錄D 語音信號的採集x,FS,bits=wavread('C:WINDOWSMediatada.wav');x=x(:,1);figure(1);subplot(2,1,1);plot(x);sound(x,FS,bits); title('語音信號時域波形圖')y=fft(x,3260);f=(FS/1630)*1:1630; subplot(2,1,2);plot(f(1:1630),abs(y(1:1630);title(

30、'語音信號頻譜圖');附錄E 對語音信號加噪t=0:length(x)-1;zs=0.05*cos(2*pi*10000*t/22050);zs0=0.05*cos(2*pi*10000*t/22050000);figure(2);subplot(2,1,1)plot(zs0)title('雜訊信號波形');zs1=fft(zs,1200);sound(zs,FS,bits); subplot(2,1,2)plot(f(1:600),abs(zs1(1:600);title('雜訊信號頻譜');x1=x+zs'sound(x1,FS,bi

31、ts); y1=fft(x1,1200);figure(3);subplot(2,1,1);plot(x1);title('加入雜訊後的信號波形');subplot(2,1,2);plot(f(1:600),abs(y1(1:600);title('加入雜訊後的信號頻譜');附錄F 巴特沃斯低通濾波器程式設(shè)計Fs=22050;x,FS,bits=wavread('C:WINDOWSMediatada.wav');x=x(:,1); figure(1); subplot(2,1,1); plot(x); sound(x,FS,bits); %回放語

32、音 title('語音信號時域波形圖')y=fft(x,3260); f=(FS/1630)*1:1630; subplot(2,1,2);plot(f(1:1630),abs(y(1:1630);title('語音信號頻譜圖'); %產(chǎn)生雜訊信號並加到語音信號 t=0:length(x)-1; zs=0.002*cos(2*pi*10000*t/22050); zs0=0.002*cos(2*pi*10000*t/22050000); figure(2); subplot(2,1,1); plot(zs0) title('雜訊信號波形');zs

33、1=fft(zs,1200); sound(zs,FS,bits); %回放噪音 subplot(2,1,2);plot(f(1:600),abs(zs1(1:600);title('雜訊信號頻譜'); x1=x+zs' sound(x1,FS,bits); %回放加入雜訊後的語音 y1=fft(x1,1200); figure(3); subplot(2,1,1);plot(x1); title('加入雜訊後的信號波形'); subplot(2,1,2); plot(f(1:600),abs(y1(1:600); title('加入雜訊後的信號

34、頻譜'); %低通濾波 fp=3000;fs=3500;Fs=22050; rp=1;rs=10; wp=2*pi*fp/Fs; ws=2*pi*fs/Fs; Fs1=1; wap=2*tan(wp/2); was=2*tan(ws/2); N,wc=buttord(wap,was,rp,rs,'s'); B,A=butter(N,wc,'s'); Bz,Az=bilinear(B,A,Fs1); figure(4); h,w=freqz(Bz,Az,512,Fs1*22050); plot(w,abs(h); title('巴特沃斯低通濾波器&

35、#39;);xlabel('頻率HZ');ylabel('耗損dB'); grid on; yd=filter(Bz,Az,x1); figure(5); subplot(2,1,1);plot(yd); title('濾波後信號波形'); ydd=fft(yd,1200); subplot(2,1,2);plot(f(1:600),abs(ydd(1:600); title('濾波後信號頻譜'); sound(yd,FS,bits);附錄G 巴特沃斯帶通濾波器程式設(shè)計Fs=22050; x,FS,bits=wavread(

36、9;C:WINDOWSMediatada.wav'); x=x(:,1); figure(1); subplot(2,1,1); plot(x); sound(x,FS,bits); %回放語音 title('語音信號時域波形圖'); y=fft(x,3260); f=(FS/1630)*1:1630; subplot(2,1,2); plot(f(1:1630),abs(y(1:1630); title('語音信號頻譜圖');%產(chǎn)生雜訊信號並加到語音信號 t=0:length(x)-1; zs0=0.002*cos(2*pi*100*t/22050);

37、 figure(2); subplot(2,1,1);plot(zs0) title('雜訊信號波形'); zs1=fft(zs0,1200); sound(zs,FS,bits); %回放噪音 subplot(2,1,2);plot(f(1:600),abs(zs1(1:600); title('雜訊信號頻譜'); x1=x+zs0' sound(x1,FS,bits); %回放加入雜訊後的語音 y1=fft(x1,1200); figure(3); subplot(2,1,1);plot(x1); title('加入雜訊後的信號波形'

38、); subplot(2,1,2); plot(f(1:600),abs(y1(1:600); title('加入雜訊後的信號頻譜'); %帶通濾波 fp=600,6000;fs=400,7000;Fs=22050; rp=1;rs=10; wp=2*pi*fp/Fs; ws=2*pi*fs/Fs; T=1;Fs1=1; wap=2*tan(wp/2); was=2*tan(ws/2); N,wc=buttord(wap,was,rp,rs,'s'); B,A=butter(N,wc,'s'); Bz,Az=bilinear(B,A,Fs1);

39、figure(4); h,w=freqz(Bz,Az,512,Fs1*22050); plot(w,abs(h); title('巴特沃斯帶通濾波器'); xlabel('頻率HZ');ylabel('耗損dB'); grid on; yd=filter(Bz,Az,x1); figure(5); subplot(2,1,1);plot(yd); ydd=fft(yd,1200); subplot(2,1,2);plot(f(1:600),abs(ydd(1:600);sound(yd,FS,bits);附錄H 巴特沃斯高通濾波器程式設(shè)計Fs=2

40、2050; x,FS,bits=wavread('C:WINDOWSMediatada.wav'); x=x(:,1); figure(1); subplot(2,1,1); plot(x);sound(x,FS,bits); %回放語音 title('語音信號時域波形圖') y=fft(x,3260); f=(FS/1630)*1:1630; subplot(2,1,2); plot(f(1:1630),abs(y(1:1630); title('語音信號頻譜圖'); %產(chǎn)生雜訊信號並加到語音信號 t=0:length(x)-1; zs0=0.

41、002*cos(2*pi*100*t/22050); figure(2); subplot(2,1,1); plot(zs0) title('雜訊信號波形'); zs1=fft(zs0,1200); sound(zs,FS,bits); %回放噪音 subplot(2,1,2); plot(f(1:600),abs(zs1(1:600); title('雜訊信號頻譜'); x1=x+zs0' sound(x1,FS,bits); %回放加入雜訊後的語音 y1=fft(x1,1200); figure(3); subplot(2,1,1);plot(x1)

42、; title('加入雜訊後的信號波形'); subplot(2,1,2); plot(f(1:600),abs(y1(1:600); title('加入雜訊後的信號頻譜'); %高通濾波 fp=600;fs=400;Fs=22050; rp=1;rs=10; wp=2*pi*fp/Fs; ws=2*pi*fs/Fs; T=1;Fs1=1; wap=2*tan(wp/2); was=2*tan(ws/2); N,wc=buttord(wap,was,rp,rs,'s'); B,A=butter(N,wc,'high','s

43、'); Bz,Az=bilinear(B,A,Fs1); figure(4); h,w=freqz(Bz,Az,512,Fs1*22050); plot(w,abs(h); title('巴特沃斯高通濾波器'); xlabel('頻率HZ');ylabel('耗損dB'); grid on; yd=filter(Bz,Az,x1); figure(5); subplot(2,1,1);plot(yd); title('濾波後信號波形'); ydd=fft(yd,1200); subplot(2,1,2);plot(f(1:

44、600),abs(ydd(1:600); title('濾波後信號頻譜'); sound(yd,FS,bits);附錄I Kaiser窗高通濾波器程式設(shè)計x,FS,bits=wavread('C:WINDOWSMediatada.wav');x=x(:,1);figure(1);subplot(2,1,1);plot(x);sound(x,FS,bits)title('語音信號時域波形圖')y=fft(x,3260);f=(FS/1630)*1:1630; subplot(2,1,2);plot(f(1:1630),abs(y(1:1630);t

45、itle('語音信號頻譜圖'); t=0:length(x)-1;zs0=0.005*cos(2*pi*100*t/22050);figure(2);subplot(2,1,1)plot(zs0)title('雜訊信號波形');zs1=fft(zs0,1200); subplot(2,1,2)plot(f(1:600),abs(zs1(1:600);title('雜訊信號頻譜');x1=x+zs0'sound(x1,FS,bits)y1=fft(x1,1200);figure(3);subplot(2,1,1);plot(x1);titl

46、e('加入雜訊後的信號波形');subplot(2,1,2);plot(f(1:600),abs(y1(1:600);title('加入雜訊後的信號頻譜');fp=600,fc=400; wp=2*pi*fp/FS;ws=2*pi*fc/FS;Bt=wp-ws; N0=ceil(11*pi/Bt); N=N0+mod(N0+1,2);wc=(wp+ws)/2/pi; hn=fir1(N-1,wc,'high',kaiser(N); X=conv(hn,x); sound(X,FS,bits);X1=fft(X,1200);figure(5);su

47、bplot(211);plot(X);title('濾波後的信號波形');subplot(212);plot(f(1:600),abs(X1(1:600); title('濾波後的信號頻譜')%顯示幅頻曲線fp=600;fs=400;Fs=22050; rp=1;rs=10; wp=2*pi*fp/Fs; ws=2*pi*fs/Fs; T=1;Fs1=1; wap=2*tan(wp/2); was=2*tan(ws/2); N,wc=buttord(wap,was,rp,rs,'s'); B,A=butter(N,wc,'high'

48、;,'s'); Bz,Az=bilinear(B,A,Fs1); figure(4); h,w=freqz(Bz,Az,512,Fs1*22050); plot(w,abs(h); title('巴特沃斯高通濾波器'); xlabel('頻率HZ');ylabel('耗損dB'); grid on;附錄J Kaiser窗低通濾波器程式設(shè)計x,FS,bits=wavread('C:WINDOWSMediatada.wav');x=x(:,1);figure(1);subplot(2,1,1);plot(x);soun

49、d(x,FS,bits)title('語音信號時域波形圖')y=fft(x,3260);f=(FS/1630)*1:1630; subplot(2,1,2);plot(f(1:1630),abs(y(1:1630);title('語音信號頻譜圖'); t=0:length(x)-1;zs=0.005*cos(2*pi*10000*t/22050);zs0=0.005*cos(2*pi*10000*t/22050000);figure(2);subplot(2,1,1)plot(zs0)title('雜訊信號波形');zs1=fft(zs,1200

50、);subplot(2,1,2)plot(f(1:600),abs(zs1(1:600);title('雜訊信號頻譜');x1=x+zs'sound(x1,FS,bits)y1=fft(x1,1200);figure(3);subplot(2,1,1);plot(x1);title('加入雜訊後的信號波形');subplot(2,1,2);plot(f(1:600),abs(y1(1:600);title('加入雜訊後的信號頻譜'); fp=3000;fs=3500;Fs=22050; rp=1;rs=10; wp=2*pi*fp/Fs; ws=2*pi*fs/Fs; Fs1=1; wap=2*tan(wp/2); was=2*tan(ws/2); N,wc=buttord(wap,was,rp,rs,'s'); B,A=butter(N,wc,'s'); Bz,Az=bilinear(B,A,Fs1); figure(4); h,w=freqz(Bz,Az,512,Fs1*22050); plot(w,abs(h); title('巴特沃斯低通濾波器');xlabel('頻率HZ');ylabel(