數(shù)字信號處理期末綜合實驗報告.doc

數(shù)字信號處理綜合實驗報告實驗題目：基于Matlab的語音信號去噪及仿真專業(yè)名稱：學 號：姓 名：日 期：報告內(nèi)容：一、 實驗原理1、 去噪的原理1.1 采樣定理在進行模擬/數(shù)字信號的轉(zhuǎn)換過程中，當采樣頻率fs.max大于信號中，最高頻率fmax的2倍時，即：fs.max=2fmax,則采樣之后的數(shù)字信號完整地保留了原始信號中的信息，一般實際應用中保證采樣頻率為信號最高頻率的510倍；采樣定理又稱奈奎斯特定理。 1924年奈奎斯特(Nyquist)就推導出在理想低通信道的最高大碼元傳輸速率的公式: 理想低通信道的最高大碼元傳輸速率=2W*log2 N (其中W是理想低通信道的帶寬,N是電平強度)為什么把采樣頻率設(shè)為8kHz?在數(shù)字通信中，根據(jù)采樣定理, 最小采樣頻率為語音信號最高頻率的2倍頻帶為F的連續(xù)信號 f(t)可用一系列離散的采樣值f(t1),f(t1t)，f(t12t)，.來表示,只要這些采樣點的時間間隔t1/2F，便可根據(jù)各采樣值完全恢復原來的信號f(t)。 這是時域采樣定理的一種表述方式。 時域采樣定理的另一種表述方式是：當時間信號函數(shù)f(t)的最高頻率分量為fM時,f(t)的值可由一系列采樣間隔小于或等于1/2fM的采樣值來確定,即采樣點的重復頻率f2fM。圖為模擬信號和采樣樣本的示意圖。 時域采樣定理是采樣誤差理論、隨機變量采樣理論和多變量采樣理論的基礎(chǔ)。對于時間上受限制的連續(xù)信號f(t)（即當tT時,f(t)=0,這里T=T2-T1是信號的持續(xù)時間），若其頻譜為F（）,則可在頻域上用一系列離散的采樣值 （1-1）采樣值來表示,只要這些采樣點的頻率間隔 （1-2）。1.2 采樣頻率采樣頻率，也稱為采樣速度或者采樣率，定義了每秒從連續(xù)信號中提取并組成離散信號的采樣個數(shù)，它用赫茲（Hz）來表示。采樣頻率的倒數(shù)是采樣周期或者叫作采樣時間，它是采樣之間的時間間隔。通俗的講采樣頻率是指計算機每秒鐘采集多少個聲音樣本，是描述聲音文件的音質(zhì)、音調(diào)，衡量聲卡、聲音文件的質(zhì)量標準。采樣頻率只能用于周期性采樣的采樣器，對于非周期性采樣的采樣器沒有規(guī)則限制。 采樣頻率的常用的表示符號是 fs。 通俗的講采樣頻率是指計算機每秒鐘采集多少個聲音樣本，是描述聲音文件的音質(zhì)、音調(diào)，衡量聲卡、聲音文件的質(zhì)量標準。采樣頻率越高，即采樣的間隔時間越短，則在單位時間內(nèi)計算機得到的聲音樣本數(shù)據(jù)就越多，對聲音波形的表示也越精確。采樣頻率與聲音頻率之間有一定的關(guān)系，根據(jù)采樣定理，只有采樣頻率高于聲音信號最高頻率的兩倍時，才能把數(shù)字信號表示的聲音還原成為原來的聲音。這就是說采樣頻率是衡量聲卡采集、記錄和還原聲音文件的質(zhì)量標準。采樣位數(shù)和采樣率對于音頻接口來說是最為重要的兩個指標，也是選擇音頻接口的兩個重要標準。無論采樣頻率如何，理論上來說采樣的位數(shù)決定了音頻數(shù)據(jù)最大的力度范圍。每增加一個采樣位數(shù)相當于力度范圍增加了6dB。采樣位數(shù)越多則捕捉到的信號越精確。對于采樣率來說你可以想象它類似于一個照相機，44.1kHz意味著音頻流進入計算機時計算機每秒會對其拍照達441000次。顯然采樣率越高，計算機攝取的圖片越多，對于原始音頻的還原也越加精確2、 數(shù)字濾波器設(shè)計的基本原理1.1 FIR數(shù)字濾波器的設(shè)計FIR：有限脈沖響應濾波器。有限說明其脈沖響應是有限的。與IIR相比，它具有線性相位、容易設(shè)計的優(yōu)點。這也就說明，IIR濾波器具有相位不線性，不容易設(shè)計的缺點。而另一方面，IIR卻擁有FIR所不具有的缺點，那就是設(shè)計同樣參數(shù)的濾波器，F(xiàn)IR比IIR需要更多的參數(shù)。這也就說明，要增加DSP的計算量。DSP需要更多的計算時間，對DSP的實時性有影響。FIR濾波器的設(shè)計比較簡單，就是要設(shè)計一個數(shù)字濾波器去逼近一個理想的低通濾波器。通常這個理想的低通濾波器在頻域上是一個矩形窗。根據(jù)傅里葉變換我們可以知道，此函數(shù)在時域上是一個采樣函數(shù)。通常此函數(shù)的表達式為：sa（n）sin（n）/n（1-3）但是這個采樣序列是無限的，計算機是無法對它進行計算的。故我們需要對此采樣函數(shù)進行截斷處理。也就是加一個窗函數(shù)。就是傳說中的加窗。也就是把這個時域采樣序列去乘一個窗函數(shù)，就把這個無限的時域采樣序列截成了有限個序列值。但是加窗后對此采樣序列的頻域也產(chǎn)生了影響：此時的頻域便不在是一個理想的矩形窗，而是成了一個有過渡帶，阻帶有波動的低通濾波器。通常根據(jù)所加的窗函數(shù)的不同，對采樣信號加窗后，在頻域所得的低通濾波器的阻帶衰減也不同。通常我們就是根據(jù)此阻帶衰減去選擇一個合適的窗函數(shù)。如矩形窗、漢寧窗、漢明窗、BLACKMAN窗、凱撒窗等。1.2 IIR數(shù)字濾波器的設(shè)計對于數(shù)字高通、帶通濾波器的設(shè)計，通用方法為雙線性變換法?？梢越柚谀M濾波器的頻率轉(zhuǎn)換設(shè)計一個所需類型的過渡模擬濾波器，再經(jīng)過雙線性變換將其轉(zhuǎn)換策劃那個所需的數(shù)字濾波器。具體設(shè)計步驟如下：（1）確定所需類型數(shù)字濾波器的技術(shù)指標。（2）將所需類型數(shù)字濾波器的邊界頻率轉(zhuǎn)換成相應的模擬濾波器的邊界頻率，轉(zhuǎn)換公式為=2/T tan(0.5) （1-4）(3)將相應類型的模擬濾波器技術(shù)指標轉(zhuǎn)換成模擬低通濾波器技術(shù)指標。（4）設(shè)計模擬低通濾波器。（5）通過頻率變換將模擬低通轉(zhuǎn)換成相應類型的過渡模擬濾波器。（6）采用雙線性變換法將相應類型的過渡模擬濾波器轉(zhuǎn)換成所需類型的數(shù)字濾波器。我們知道，脈沖響應不變法的主要缺點是會產(chǎn)生頻譜混疊現(xiàn)象，使數(shù)字濾波器的頻響偏離模擬濾波器的頻響特性。為了克服之一缺點，可以采用雙線性變換法。二、 實驗步驟 其大概流程框圖可如下表示：（圖2-1）語音信號采集效果顯示、對比語音信號錄入語音信號變換信號加噪語音信號濾波 圖2-11. 濾波器的設(shè)計。 利用窗函數(shù)法設(shè)計FIR濾波器的步驟。如下：（1）根據(jù)對阻帶衰減及過渡帶的指標要求，選擇串窗數(shù)類型（矩形窗、三角窗、漢寧窗、哈明窗、凱塞窗等），并估計窗口長度N。先按照阻帶衰減選擇窗函數(shù)類型。原則是在保證阻帶衰減滿足要求的情況下，盡量選擇主瓣的窗函數(shù)。（2）構(gòu)造希望逼近的頻率響應函數(shù)。（3）計算h(n).。（4）加窗得到設(shè)計結(jié)果。接下來，我們根據(jù)語音信號的特點給出有關(guān)濾波器的技術(shù)指標：低通濾波器的性能指標：fp=1000Hz，fc=1200Hz，As=50db ,Ap=1dB高通濾波器的性能指標：fp=3500Hz，fc=4000Hz，As=50dB，Ap=1dB在Matlab中，可以利用函數(shù)fir1設(shè)計FIR濾波器，利用Matlab中的函數(shù)freqz畫出各步步器的頻率響應。MATLAB信號處理工具箱函數(shù)cheblap,cheblord和cheeby1是切比雪夫I型濾波器設(shè)計函數(shù)。我們用到的是cheeby1函數(shù)，其調(diào)用格式如下：B,A=cheby1(N,Rp,wpo,ftypr)B,A=cheby1(N,Rp,wpo,ftypr,s)利用模擬濾波器設(shè)計IIR數(shù)字低通濾波器的步驟：如下：（1）確定數(shù)字低通濾波器的技術(shù)指標：通帶邊界頻率、通帶最大衰減，阻帶截止頻率、阻帶最小衰減。（2）將數(shù)字低通濾波器的技術(shù)指標轉(zhuǎn)換成相應的模擬低通濾波器的技術(shù)指標。（3）按照模擬低通濾波器的技術(shù)指標設(shè)計及過渡模擬低通濾波器。（4）用雙線性變換法，模擬濾波器系統(tǒng)函數(shù)轉(zhuǎn)換成數(shù)字低通濾波器系統(tǒng)函數(shù)。MATLAB信號處理工具箱函數(shù)cheblap,cheblord和cheeby1是切比雪夫I型濾波器設(shè)計函數(shù)。我們用到的是cheeby1函數(shù)，其調(diào)用格式如下：B,A=cheby1(N,Rp,wpo,ftypr)B,A=cheby1(N,Rp,wpo,ftypr,s)函數(shù)butter,cheby1和ellip設(shè)計IIR濾波器時都是默認的雙線性變換法，所以在設(shè)計濾波器時只需要代入相應的實現(xiàn)函數(shù)即可。2. 語音信號的錄制。單擊自己的電腦開始程序，選擇所有程序，接著選擇附件，再選擇錄音。自己錄入語音信號，然后保存在MATLAB文件夾里面，命名為“chenghaijie1.wav”。3. 在MATLAB平臺上讀入語音信號、繪制頻譜圖并回放原始語音信號。利用MATLAB中的wavread命令來讀入（采集）語音信號，將它賦值給某一向量。y,fs,bits=wavread( N1 N2);用于讀取語音，采樣值放在向量y中，fs表示采樣頻率(Hz)，bits表示采樣位數(shù)。N1 N2表示讀取從N1點到N2點的值（若只有一個N的點則表示讀取前N點的采樣值。 4. 利用matlab函數(shù)randn編程加入一段隨機噪音信號，再利用設(shè)計的FIR和IIR濾波器去噪，分別繪制去噪后的頻譜圖、回放語音信號與原始信號的頻譜圖、原始語音信號比較，并且比較兩種濾波器的優(yōu)缺點和得出語音信號的頻段。三、 實驗結(jié)果下面我們將給出設(shè)計FIR數(shù)字濾波器的主要程序和圖像。 FIR低通濾波器程序見附錄1； FIR低通濾波器圖像：（圖22）圖22 FIR低通濾波器 FIR高通濾波程序見附錄2； FIR高通濾波圖像：（圖23）圖23 FIR高通濾波器下面我們將給出IIR數(shù)字濾波器的主要程序。IIR低通濾波器程序見附錄3；IIR低通濾波器圖像：（圖24）圖24 IIR低通濾波器IIR濾波器高通程序見附錄4；IIR濾波器高通圖像：（圖25）圖25 IIR高通濾波器原始語音信號的回放及頻譜分析程序：見附錄5；原始語音信號及其頻譜分析圖像、運行程序可以聽到原始信號的回放：（圖26） 圖26 原始語音信號及其頻譜分析圖像利用randn函數(shù)把一段隨機噪音信號加入原始語音信號的信號處理過程程序：見附錄6；加噪后語音信號的時域波形、頻譜圖：（圖27） （圖27）加噪后語音信號的時域波形、頻譜圖利用FIR低通濾波器法去噪程序：見附錄7；FIR濾波器去噪前與去噪后語音信號的時域波形、頻譜圖：（圖28） （圖28）FIR濾波器去噪前與去噪后語音信號的時域波形、頻譜圖利用IIR低通濾波器法去噪程序：見附錄8；IIR濾波器去噪前與去噪后語音信號的時域波形、頻譜圖：（圖29） （圖29）IIR濾波器去噪前與去噪后語音信號的時域波形、頻譜圖四、 實驗分析1. 由（圖28）FIR濾波器去噪前與去噪后語音信號的時域波形、頻譜圖和（圖29）IIR濾波器去噪前與去噪后語音信號的時域波形、頻譜圖可看出原始語音信號和加噪語音信號時域波形和頻譜圖的區(qū)別。加噪后的語音信號的時域波形比原始語音信號要模糊得多，頻譜圖則是在頻率5000Hz以后出現(xiàn)了明顯的變化。再通過濾波前的信號波形和頻譜圖的對比，可以明顯看出濾波后的波形開始變得清晰了，有點接近原始信號的波形圖了。濾波后信號的頻譜圖也在5000Hz以后開始逐漸接近原始語音信號的頻譜圖。再從對語音信號的回放，人耳可以明顯辨別出加噪后的語音信號比較渾濁，還有很明顯嘎吱嘎吱的雜音在里面。濾波后，語音信號較加噪后的信號有了明顯的改善，基本可以聽清楚了，而且雜音也沒有那么強烈，但是聲音依然沒有原始語音信號那么清晰脆耳。2. 結(jié)合去噪后的頻譜圖對比可知FIR濾波器和IIR濾波器的優(yōu)缺點IIR數(shù)字濾波器采用遞歸型結(jié)構(gòu)，即結(jié)構(gòu)上帶有反饋環(huán)路。IIR濾波器運算結(jié)構(gòu)通常由延時、乘以系數(shù)和相加等基本運算組成，可以組合成直接型、正準型、級聯(lián)型、并聯(lián)型四種結(jié)構(gòu)形式，都具有反饋回路。由于運算中的舍入處理，使誤差不斷累積，有時會產(chǎn)生微弱的寄生振蕩。 （1）IIR數(shù)字濾波器的相位特性不好控制，對相位要求較高時，需加相位校準網(wǎng)絡(luò)。FIR濾波器則要求較低。（2）IIR濾波器運算誤差大，有可能出現(xiàn)極限環(huán)振蕩，F(xiàn)IR相比之下運算誤差較小，不會出現(xiàn)極限環(huán)振蕩。（3）IIR幅頻特性精度很高，不是線性相位的，可以應用于對相位信息不敏感的音頻信號上； （4）與FIR濾波器的設(shè)計不同，IIR濾波器設(shè)計時的階數(shù)不是由設(shè)計者指定，而是根據(jù)設(shè)計者輸入的各個濾波器參數(shù)（截止頻率、通帶濾紋、阻帶衰減等），由軟件設(shè)計出滿足這些參數(shù)的最低濾波器階數(shù)。在MATLAB下設(shè)計不同類型IIR濾波器均有與之對應的函數(shù)用于階數(shù)的選擇。 （5）IIR單位響應為無限脈沖序列FIR單位響應為有限的 （6）FIR幅頻特性精度較之于iir低，但是線性相位，就是不同頻率分量的信號經(jīng)過FIR濾波器后他們的時間差不變。這是很好的性質(zhì)。 （7）IIR濾波器有噪聲反饋，而且噪聲較大，F(xiàn)IR濾波器噪聲較小。FIR幅頻特性精度較之于IIR低，但是線性相位，就是不同頻率分量的信號經(jīng)過FIR濾波器后他們的時間差不變，這是很好的性質(zhì)。3. 通過對語音信號做頻譜分析和通過比較加噪前后，語音的頻譜和語音回放，能明顯的感覺到加入噪聲后回放的聲音與原始的語音信號有很大的不同,前者隨較尖銳的干擾嘯叫聲。從含噪語音信號的頻譜圖中可以看出含噪聲的語音信號頻譜,在整個頻域范圍內(nèi)分是布均勻。其實，這正是干擾所造成的。通過濾波前后的對比，低通濾波后效果最好，高通濾波后的效果最差。由此可見，語音信號主要分布在低頻段，而噪聲主要分布在高頻段。五、 實驗總結(jié)基于Matlab的語音信號去噪及仿真課題的特色在于它將語音信號看作一個向量，于是就把語音數(shù)字化了。那么，就可以完全利用數(shù)字信號處理的知識來解決語音及加噪處理問題。我們可以像給一般信號做頻譜分析一樣，來對語音信號做頻譜分析，也可以較容易的用數(shù)字濾波器來對語音進行濾波處理。通過比較加噪前后，語音的頻譜和語音回放，能明顯的感覺到加入噪聲后回放的聲音與原始的語音信號有很大的不同,前者隨較尖銳的干擾嘯叫聲。從含噪語音信號的頻譜圖中可以看出含噪聲的語音信號頻譜,在整個頻域范圍內(nèi)分是布均勻。其實，這正是干擾所造成的。通過兩種濾波器濾波前后的對比，可知FIR濾波器和IIR濾波器的優(yōu)缺點，低通濾波后效果最好，高通濾波后的效果最差。由此可見，語音信號主要分布在低頻段，而噪聲主要分布在高頻段。六、 參考文獻1 Sanjit K.Mitra，數(shù)字信號處理實驗指導書，電子工業(yè)出版社，2005 年 1月 2 胡航，語音信號處理，哈爾濱工業(yè)大學出版社，2000 年 5 月 3姚天任.數(shù)字語音處理M.武漢：華中科技大學出版社，2005 年 3 月七、 附錄附錄 1Ft=8000;%設(shè)置抽樣頻率為8000HzFp=1000;%設(shè)置通帶邊界為1000HzFs=1200;%設(shè)置阻帶邊界為1200Hzwp=2*Fp/Ft;%計算歸一化通帶邊界角頻率ws=2*Fs/Ft;%計算歸一化阻帶邊界角頻率rp=1;%設(shè)置通帶紋波為1dBrs=50;%設(shè)置阻帶衰減為50dBp=1-10.(-rp/20);%通帶阻帶波紋s=10.(-rs/20);fpts=wp ws;mag=1 0;dev=p s;n21,wn21,beta,ftype=kaiserord(fpts,mag,dev);%kaiserord求階數(shù)截止頻率b21=fir1(n21,wn21,Kaiser(n21+1,beta);%由fir1設(shè)計濾波器h,w=freqz(b21,1);%得到頻率響應plot(w/pi,abs(h);%繪制FIR低通濾波器title(FIR低通濾波器);%標題附錄2Ft=8001;%設(shè)置抽樣頻率為8001HzFp=4000;%設(shè)置通帶邊界為4000HzFs=3500;%設(shè)置阻帶邊界為3500Hzwp=2*Fp/Ft;%計算歸一化通帶邊界角頻率ws=2*Fs/Ft;%計算歸一化阻帶邊界角頻率rp=1;%設(shè)置通帶紋波為1dBrs=50;%設(shè)置阻帶衰減為50dBp=1-10.(-rp/20);%通帶阻帶波紋s=10.(-rs/20);fpts=ws wp;mag=0 1;dev=p s;n23,wn23,beta,ftype=kaiserord(fpts,mag,dev);%kaiserord求階數(shù)截止頻率b23=fir1(n23,wn23,high,Kaiser(n23+1,beta);%由fir1設(shè)計濾波器h,w=freqz(b23,1);%得到頻率響應plot(w*12000*0.5/pi,abs(h);%繪制FIR高通濾波器title(FIR高通濾波器);%標題axis(3000 6000 0 1.2);%確定橫縱坐標的范圍附錄3Ft=8000;%設(shè)置抽樣頻率為8000HzFp=1000;%設(shè)置通帶邊界為1000HzFs=1200;%設(shè)置阻帶邊界為1200Hzwp=2*pi*Fp/Ft;%計算歸一化通帶邊界角頻率ws=2*pi*Fs/Ft;%計算歸一化阻帶邊界角頻率fp=2*Ft*tan(wp/2);%計算通帶邊界頻率fs=2*Fs*tan(wp/2);%計算阻帶邊界頻率n11,wn11=buttord(wp,ws,1,50,s);%求低通濾波器的階數(shù)和截止頻率b11,a11=butter(n11,wn11,s);%利用butter求S域的頻率響應的參數(shù)num11,den11=bilinear(b11,a11,0.5);%雙線性變換實現(xiàn)S域到Z域的變換h,w=freqz(num11,den11);%根據(jù)參數(shù)求出頻率響應plot(w*8000*0.5/pi,abs(h);%繪制IIR低通濾波器title(IIR低通濾波器);%標題legend(用butter設(shè)計);%添加圖例的標注grid;%使繪制的圖有網(wǎng)格附錄4Ft=8001;%設(shè)置抽樣頻率為8001HzFp=4000;%設(shè)置通帶邊界為4000HzFs=3500;%設(shè)置阻帶邊界為3500Hzwp1=tan(pi*Fp/Ft);%通帶到低通濾波器參數(shù)轉(zhuǎn)換ws1=tan(pi*Fs/Ft);%阻帶到低通濾波器參數(shù)轉(zhuǎn)換wp=1;%計算歸一化通帶邊界角頻率ws=wp1*wp/ws1;%計算歸一化阻帶邊界角頻率n13,wn13=cheb1ord(wp,ws,1,50,s);%求模擬的低通濾波器階數(shù)和截止頻b13,a13=cheby1(n13,1,wn13,s);%利用cheby1求S域的頻率響應的參數(shù)num,den=lp2hp(b13,a13,wn13); %將S域低通參數(shù)轉(zhuǎn)為高通的num13,den13=bilinear(num,den,0.5); %利用雙線性變換實現(xiàn)S域到Z域轉(zhuǎn)h,w=freqz(num13,den13);%得到頻率響應plot(w*21000*0.5/pi,abs(h);%繪制IIR高通濾波器title(IIR高通濾波器);%標題legend(用cheby1設(shè)計);%添加圖例的標注附錄5x,fs,bits=wavread(chenghaijie1.wav);%利用MATLAB中的wavread命令來讀入（采集）語音信號，將它賦值給某一向量sound(x,fs,bits);%播放語音信號X=fft(x,100000);%對信號做50000點FFT變換magX=abs(X);%求幅度angX=angle(X);%求相位plot(x);title(原始信號波形);%繪制原始信號的波形plot(X); title(原始語音信號采樣后的頻譜圖)%繪制原始語音信號采樣后的頻譜圖附錄6y,fs,bits=wavread(chenghaijie1.wav);%利用MATLAB中的wavread命令來讀入（采集）語音信號，將它賦值給某一向量sound(y,fs)%播放語音信號n=length(y)%計算向量y的長度y_p=fft(y,n);%對向量y進行fft變換f=fs*(0:n/2-1)/n;%對應點的頻率figure(1)%建立第一幅圖subplot(2,1,1);%分配繪制的位置為第一塊位置plot(y);%繪制原始語音信號采樣后時域波形title(原始語音信號采樣后時域波形);%標題xlabel(時間軸)%標注X軸ylabel(幅值 A)%標注Y軸subplot(2,1,2);%分配繪制的位置為第二塊位置plot(f,abs(y_p(1:n/2);%繪制原始語音信號采樣后的頻譜圖title(原始語音信號采樣后的頻譜圖);xlabel(頻率Hz);%標注X軸ylabel(頻率幅值);%標注Y軸L=length(y)%計算向量y的長度noise=0.1*randn(L,2);%設(shè)置加入的隨機噪聲函數(shù)y_z=y+noise;%將原始信號與隨機噪聲函數(shù)相疊加sound(y_z,fs)%試聽加入隨機噪聲函數(shù)的信號n=length(y);%計算向量y的長度y_zp=fft(y_z,L);%對加入正弦噪聲向量y進行fft變換f=fs*(0:L/2-1)/L;%計算對應點的頻率figure(2)%建立第二幅圖subplot(2,1,1);%分配繪制的位置為第一塊位置plot(y_z);%繪制加噪語音信號時域波形title(加噪語音信號時域波形);%標題xlabel(時間軸)%標注X軸ylabel(幅值 A)%標注Y軸subplot(2,1,2);%分配繪制的位置為第二塊位置plot(f,abs(y_zp(1:L/2);%繪制加噪語音信號頻譜圖title(加噪語音信號頻譜圖);%標題xlabel(頻率Hz);%標注X軸ylabel(頻率幅值);%標注Y軸附錄7y,fs,bits=wavread(chenghaijie1.wav);%利用MATLAB中的wavread命令來讀入（采集）語音信號，將它賦值給某一向量sound(y,fs)%播放語音信號n=length(y)%計算向量y的長度y_p=fft(y,n);%對向量y進行fft變換f=fs*(0:n/2-1)/n;%對應點的頻率figure(1)%建立第一幅圖subplot(2,1,1);%分配繪制的位置為第一塊位置plot(y);%繪制原始語音信號采樣后時域波形title(原始語音信號采樣后時域波形);%標題xlabel(時間軸)%標注X軸ylabel(幅值 A)%標注Y軸subplot(2,1,2);%分配繪制的位置為第二塊位置plot(f,abs(y_p(1:n/2);%繪制原始語音信號采樣后的頻譜圖title(原始語音信號采樣后的頻譜圖);xlabel(頻率Hz);%標注X軸ylabel(頻率幅值);%標注Y軸L=length(y)%計算向量y的長度noise=0.1*randn(L,2);%設(shè)置加入的隨機噪聲函數(shù)y_z=y+noise;%將原始信號與隨機噪聲函數(shù)相疊加sound(y_z,fs)%試聽加入隨機噪聲函數(shù)的信號n=length(y);%計算向量y的長度y_zp=fft(y_z,L);%對加入正弦噪聲向量y進行fft變換f=fs*(0:L/2-1)/L;%計算對應點的頻率figure(2)%建立第二幅圖subplot(2,1,1);%分配繪制的位置為第一塊位置plot(y_z);%繪制加噪語音信號時域波形title(加噪語音信號時域波形);%標題xlabel(時間軸)%標注X軸ylabel(幅值 A)%標注Y軸subplot(2,1,2);%分配繪制的位置為第二塊位置plot(f,abs(y_zp(1:L/2);%繪制加噪語音信號頻譜圖title(加噪語音信號頻譜圖);%標題xlabel(頻率Hz);%標注X軸ylabel(頻率幅值);%標注Y軸Ft=5000;%設(shè)置抽樣頻率為5000HzFp=1000;%設(shè)置通帶邊界為1000HzFs=1200;%設(shè)置阻帶邊界為1200Hzwp=2*Fp/Ft;%計算歸一化通帶邊界角頻率ws=2*Fs/Ft;%計算歸一化阻帶邊界角頻率rp=1;%設(shè)置通帶紋波為1dBrs=50;%設(shè)置阻帶衰減為50dBp=1-10.(-rp/20);%通帶阻帶波紋s=10.(-rs/20);fpts=wp ws;mag=1 0;dev=p s;n21,wn21,beta,ftype=kaiserord(fpts,mag,dev);%kaiserord求階數(shù)截止頻率b21=fir1(n21,wn21,Kaiser(n21+1,beta);%由fir1設(shè)計濾波器h,w=freqz(b21,1);%得到頻率響應plot(w/pi,abs(h);%繪制FIR低通濾波器title(FIR低通濾波器);%標題x=fftfilt(b21,y_z);%重疊相加法實現(xiàn)線性卷積的計算X=fft(x,n);%計算x的fftfigure(3);%建立第三幅圖subplot(2,2,1);plot(f,abs(y_zp(1:n/2);%繪制濾波前信號的頻譜title(濾波前信號的頻譜);%標題subplot(2,2,2);plot(f,abs(X(1:n/2);%繪制濾波后信號的頻譜title(濾波后信號的頻譜);%標題subplot(2,2,3);plot