如何利用matlab處理音頻信號

1、MatIab處理音頻信號一、問題的提出：數(shù)字語音是信號的一種，我們處理數(shù)字語音信號，也就是對一 種信號的處理，那信號是什么呢？ 信號是傳遞信息的函數(shù)。一、問題的提出： 數(shù)字語音是信號的一種，我們處理數(shù)字語音信號，也就是對一種信號的處理，那信 號是什么呢？ 信號是傳遞信息的函數(shù)。離散時間信號 26mdash;%26mdas序列26mdash;%26mdas可以用圖形來表示。按信號特點的不同，信號可表示成一個或幾個獨立變量的函數(shù)。例如，圖像信號就 是空間位置（二元變量）的亮度函數(shù)。一維變量可以是時間，也可以是其他參量， 習(xí)慣上將其看成時間。信號有以下幾種：（1）連續(xù)時間信號：在連續(xù)時間范圍內(nèi)定義的

2、信號，但信號的幅值可以是連續(xù)數(shù) 值，也可以是離散數(shù)值。當(dāng)幅值為連續(xù)這一特點情況下又常稱為模擬信號。實際上 連續(xù)時間信號與模擬信號常常通用，用以說明同一信號。（2）離時間信號：時間為離散變量的信號,即獨立變量時間被量化了。而幅度仍是連續(xù)變化的（3）數(shù)字信號：時間離散而幅度量化的信號。語音信號是基于時間軸上的一維數(shù)字信號，在這里主要是對語音信號進(jìn)行頻域上的分析。在信號分析中，頻域往往包含了更多的信息。對于頻域來說，大概有 8種波 形可以讓我們分析：矩形方波，鋸齒波，梯形波，臨界阻尼指數(shù)脈沖波形，三角 波，余旋波，余旋平方波，高斯波。對于各種波形，我們都可以用一種方法來分 析，就是傅立葉變換：將時域

3、的波形轉(zhuǎn)化到頻域來分析。于是，本課題就從頻域的角度對信號進(jìn)行分析，并通過分析頻譜來設(shè)計出合適的濾 波器。當(dāng)然，這些過程的實現(xiàn)都是在 MATLAB軟件上進(jìn)行的，MATLAB軟件在數(shù) 字信號處理上發(fā)揮了相當(dāng)大的優(yōu)勢。設(shè)計方案:利用MATLAB中的WaVreaC命令來讀入（采集）語音信號，將它賦值給某一向 量。再將該向量看作一個普通的信號，對其進(jìn)行FFT變換實現(xiàn)頻譜分析，再依據(jù)實際情況對它進(jìn)行濾波。對于波形圖與頻譜圖（包括濾波前后的對比圖）都可以用 MATLAB畫出。我們還可以通過SOUnd命令來對語音信號進(jìn)行回放，以便在聽覺 上來感受聲音的變化。選擇設(shè)計此方案，是對數(shù)字信號處理的一次實踐。在數(shù)字信

4、號處理的課程學(xué)習(xí)過程 中，我們過多的是理論學(xué)習(xí)，幾乎沒有進(jìn)行實踐方面的運(yùn)用。這個課題正好是對數(shù) 字語音處理的一次有利實踐，而且語音處理也可以說是信號處理在實際應(yīng)用中很大 眾化的一方面。這個方案用到的軟件也是在數(shù)字信號處理中非常通用的一個軟件26mdash;%26mdash;MATLAB軟件。所以這個課題的設(shè)計過程也是一次數(shù)字信 號處理在MATLAB中應(yīng)用的學(xué)習(xí)過程。課題用到了較多的 MATLAB語句，而由于 課題研究范圍所限，真正與數(shù)字信號有關(guān)的命令函數(shù)卻并不多。主體部分:（一）、語音的錄入與打開:y,fs,bits=wavread（'Blip',N1 N2）用于讀取語音，采樣

5、值放在向量 y中，fs表示采樣 頻率（Hz），bits表示采樣位數(shù)。N1 N2表示讀取從N1點到N2點的值（若只有一 個N的點則表示讀取前N點的采樣值）。SOund（x,fs,bits）;用于對聲音的回放。向量y則就代表了一個信號（也即一個復(fù)雜 的%26ldqu0函數(shù)表達(dá)式%26rdquo）也就是說可以像處理一個信號表達(dá)式一樣處理 這個聲音信號。FFT的MATLAB實現(xiàn)在MATLAB的信號處理工具箱中函數(shù)FFT和IFFT用于快速傅立葉變換和逆變 換。下面介紹這些函數(shù)。函數(shù)FFT用于序列快速傅立葉變換。函數(shù)的一種調(diào)用格式為y=fft(x)其中，X是序列，y是序列的FFT, X可以為一向量或矩陣，

6、若X為一向量，y是X 的FFTO且和X相同長度。若X為一矩陣，則y是對矩陣的每一列向量進(jìn)行FFTO如果X長度是2的幕次方，函數(shù)fft執(zhí)行高速基2FFT算法；否則fft執(zhí)行一種混 合基的離散傅立葉變換算法，計算速度較慢。函數(shù)FFT的另一種調(diào)用格式為y=fft(x,N)式中，X, y意義同前，N為正整數(shù)。函數(shù)執(zhí)行N點的FFT。若X為向量且長度小于N ,貝U函數(shù)將X補(bǔ)零至長度N。若 向量X的長度大于N ,則函數(shù)截短X使之長度為N。若X為矩陣，按相同方法對X 進(jìn)行處理。經(jīng)函數(shù)fft求得的序列y 般是復(fù)序列，通常要求其幅值和相位。MATLAB提供求復(fù)數(shù)的幅值和相位函數(shù)：abs ange這些函數(shù)一般和FF

7、T同時使用。函數(shù)abs(x用于計算復(fù)向量X的幅值，函數(shù)angle(x用于計算復(fù)向量的相角，介于 和之間，以弧度表示。函數(shù)UnWraP(P用于展開弧度相位角P ,當(dāng)相位角絕對變化超過時，函數(shù)把它擴(kuò)展 至。用MATLAB工具箱函數(shù)fft進(jìn)行頻譜分析時需注意：(1) 函數(shù)fft返回值y的數(shù)據(jù)結(jié)構(gòu)對稱性若已知序列 x=4,3,2,6,7,8,9,0求 X(k)=DFTx(n)。利用函數(shù)fft計算，用MATLAB編程如下：N=8;n=0:N-1;Xn=4 3 2 6 7 8 9 0'XK=fft(x n)結(jié)果為：XK =39.0000-10.7782 + 6.2929i0 - 5.0000i4.

8、7782 - 7.7071i5.00004.7782 + 7.7071i0 + 5.0000i-10.7782 - 6.2929i由程序運(yùn)行所得結(jié)果可見，X(k)和x(n)的維數(shù)相同，共有8個元素。X(k)的第一行 元素對應(yīng)頻率值為0,第五行元素對應(yīng)頻率值為 NyquiSt頻率，即標(biāo)準(zhǔn)頻率為1.因 此第一行至第五行對應(yīng)的標(biāo)準(zhǔn)頻率為 01。而第五行至第八行對應(yīng)的是負(fù)頻率， 其X(k)值是以NyqUiSt頻率為軸對稱。(注：通常表示為NyqUiSt頻率外擴(kuò)展，標(biāo)以正值。)一般而言，對于N點的x(n)序列的FFT是N點的復(fù)數(shù)序列，其點n=N2+1對應(yīng) NyqUiSt頻率，作頻譜分析時僅取序列 X(k

9、)的前一半，即前N/2點即可。X(k)的后 一半序列和前一半序列時對稱的。(2) 頻率計算若N點序列x(n)(n=0,1,N-1)是在采樣頻率 下獲得的。它的FFT也是N點序列， 即X(k)(k=0,1,2,N-1),則第k點所對應(yīng)實際頻率值為f=k*f /N.(3) 作FFT分析時，幅值大小與FFT選擇點數(shù)有關(guān)，但不影響分析結(jié)果。2、設(shè)計內(nèi)容：(1) 下面的一段程序是語音信號在MATLAB中的最簡單表現(xiàn)，它實現(xiàn)了語音的讀 入打開，以及繪出了語音信號的波形頻譜圖。x,fs,bits=wavread('di ng.wav',1024 5120);SOun d(x,fs,bits)

10、;X=fft(x,4096);magX=abs(X);an gX=a ngle(X);SUbPlOt(221);PlOt(x);title原始信號波形');SUbPIot(222);PIot(X); title(原始信號頻譜');SUbPlOt(223);plot(magX);title原始信號幅值');SUbPIOt(224);plot(angX);title原始信號相位');程序運(yùn)行可以聽到聲音，得到的圖形為：(2) 定點分析：已知一個語音信號，數(shù)據(jù)采樣頻率為100Hz,試分別繪制N = 128 點DFT的幅頻圖和N = 1024點DFT幅頻圖。編程如下：X

11、=WaVreadcdi ng.wav');SOun d(x);fs=100;N=128;y=fft(x,N);magy=abs(y);f=(0:Ie ngth(y)-1)'*fsle ngth(y);SUbPlOt(221);plot(f,magy);xlabel(頻率(Hz) ')ylabel幅值'); title('N=128(a)');gridSUbPlOt(222);plot(f(1:N/2),magy(1:N/2); xlabel(頻率(Hz) ');ylabel幅值'); title('N=128(b)'

12、;);gridfs=100;N=1024;y=fft(x,N);magy=abs(y);f=(0:le ngth(y)-1)'*fsle ngth(y);SUbPlOt(223);plot(f,magy);XIabe1(頻率(Hz) ')ylabel幅值'); title('N=1024(c)');gridSubPlOt(224);plot(f(1:N/2),magy(1:N/2); xlabel(頻率(Hz) ')ylabel幅值'); title('N=1024(d)');grid運(yùn)行結(jié)果如圖：上圖(a) (b)為N=

13、128點幅頻譜圖，(c)、(d)為N=1024點幅頻譜圖。由于采樣頻率f =100HZ ,故 NyqUiSt頻率為 50Hz。(a)、(C)是 0100Hz頻譜圖，(b)、(d)是 050Hz 頻譜圖。由(a)或(C)可見，整個頻譜圖是以NyqUiSt頻率為軸對稱的。因此利用fft對 信號作頻譜分析，只要考察0NyqUiSt頻率(采樣頻率一半)范圍的幅頻特性。比較(a和(C)或(b)和(d)可見，幅值大小與fft選用點數(shù)N有關(guān)，但只要點數(shù)N足夠 不影響研究結(jié)果。從上圖幅頻譜可見，信號中包括15Hz和40Hz的正弦分量。(3)若信號長度T=25.6s,即抽樣后x(n)點數(shù)為TTs=256 ,所得

14、頻率分辨率為HZ ,以此觀察數(shù)據(jù)長度N的變化對DTFT分辨率的影響：編程如下：x,fs,bits=wavread('di ng.wav');N=256;f=0:fs/N:fs/2-1/N;X=fft(x);X=abs(X);subplot(211)XIabel('Hz'),ylabel('H(ejw)') %數(shù)據(jù)長度N擴(kuò)大4倍后觀察信號頻譜N=N*4;f=0:fs/N:fs/2-1/N;X=fft(x);X=abs(X);subplot(212)plot(f(45*4:4*60),X(4*45:4*60);gridXIabeI('Hz&#

15、39;),ylabel('H(ejw)')結(jié)果如圖：(三)、濾波器設(shè)計：1、相關(guān)原理：設(shè)計數(shù)字濾波器的任務(wù)就是尋求一個因果穩(wěn)定的線性時不變系統(tǒng)，并使系統(tǒng)函數(shù) H(Z)具有指定的頻率特性。數(shù)字濾波器從實現(xiàn)的網(wǎng)絡(luò)結(jié)構(gòu)或者從單位沖激響應(yīng)分類，可以分成無限長單位沖激 響應(yīng)(IIR)數(shù)字濾波器和有限長單位沖激響應(yīng)(FIR)數(shù)字濾波器。數(shù)字濾波器頻率響應(yīng)的三個參數(shù)：(1) 幅度平方響應(yīng):(2) 相位響應(yīng)其中，相位響應(yīng)(3) 群時延響應(yīng)IIR數(shù)字濾波器：IIR數(shù)字濾波器的系統(tǒng)函數(shù)為的有理分?jǐn)?shù)，即IIR數(shù)字濾波器的逼近問題就是求解濾波器的系數(shù) 和，使得在規(guī)定的物理意義上 逼近所要求的特性的問題

16、。如果是在 S平面上逼近，就得到模擬濾波器，如果是在 Z平面上逼近，則得到數(shù)字濾波器。FIR數(shù)字濾波器：設(shè)FIR的單位脈沖響應(yīng)h(n)為實數(shù)，長度為N ,則其Z變換和頻率響應(yīng)分別為按頻域采樣定理FIR數(shù)字濾波器的傳輸函數(shù) H(Z)和單位脈沖響應(yīng)h(n)可由它的N 個頻域采樣值H(k)唯一確定。MATLAB中提供了幾個函數(shù)，分別用于實現(xiàn)IIR濾波器和FIR濾波器。(1) 卷積函數(shù)COnV卷積函數(shù)COnV的調(diào)用格式為c=conv(a,b)該格式可以計算兩向量a和b的卷積，可以直接用于對有限長信號采用FIR濾波器的濾波。(2) 函數(shù) filter函數(shù)filter的調(diào)用格式為y=filter(b,a,

17、x)該格式采用數(shù)字濾波器對數(shù)據(jù)進(jìn)行濾波，既可以用于IIR濾波器，也可以用于FIR濾波器。其中向量b和a分別表示系統(tǒng)函數(shù)的分子、分母多項式的系數(shù)，若a=1此時表示FIR濾波器，否則就是IIR濾波器。該函數(shù)是利用給出的向量 b和 a,對X中的數(shù)據(jù)進(jìn)行濾波，結(jié)果放入向量y。(3)函數(shù) fftfilt函數(shù)fftfilt的調(diào)用格式為y=fftfilt(b,x)該格式是利用基于FFT的重疊相加法對數(shù)據(jù)進(jìn)行濾波，這種頻域濾波技術(shù)只對 FIR濾波器有效。該函數(shù)是通過向量b描述的濾波器對X數(shù)據(jù)進(jìn)行濾波。關(guān)于用butter函數(shù)求系統(tǒng)函數(shù)分子與分母系數(shù)的幾種形式。b,a=butter(N,wc,'high&

18、#39;)設(shè)計N階高通濾波器，WC為它的3dB邊緣頻率，以 為單 位，故。b,a=butter(N,wc):當(dāng)WC為具有兩個元素的矢量 wc=w1,w2時，它設(shè)計2N階帶通濾 波器，3dB通帶為，W的單位為。b,a=butter(N,wc,'stop')若 wc=w1,w2,則它設(shè)計 2N 階帶阻濾波器，3dB 通帶為，W的單位為。如果在這個函數(shù)輸入變元的最后，加一個變元26ldquo;s%26rdquq；表示設(shè)計的是模擬濾波器。這里不作討論。為了設(shè)計任意的選項巴特沃斯濾波器，必須知道階數(shù)N和3dB邊緣頻率矢量WCO這可以直接利用信號處理工具箱中的buttord函數(shù)來計算。如果

19、已知濾波器指標(biāo)，和，則調(diào)用格式為N,wc=buttord(wp,ws,Rp,As)對于不同類型的濾波器，參數(shù) WP和WS有一些限制：對于低通濾波器，wp%26lt;WS對于高通濾波器，wp%26gt;WS對于帶通濾波器，WP和WS分別為具 有兩個元素的矢量，wp=wp1,wp2和ws=ws1,ws2,并且 ws1%26lt;wp1%26lt;wp2%26lt;ws2 對于帶阻濾波器wp1%26lt;ws1%26lt;ws2%26lt;wp22、設(shè)計內(nèi)容:(1)濾波器示例：在這里為了說明如何用 MATLAB來實現(xiàn)濾波，特舉出一個簡單的函數(shù)信號濾波實 例(對信號x(n)=sin( n4)+5cos

20、( n2)進(jìn)行濾波，信號長度為500點)，從中了解濾 波的實現(xiàn)過程。程序如下：Wn=0.2*pi;N=5;b,a=butter(N,Wnpi);n=0:499;X=Si n(pi* n4)+5*cos(pi* n/2);X=fft(x,4096);SubPlOt(221);PlOt(x);title濾波前信號的波形');SUbPlot(222);PlOt(X);title(濾波前信號的頻譜');y=filter(b,a,x);Y=fft(y,4096);SUbPlOt(223);plot(y);title濾波后信號的波形');SUbPlOt(224);PlOt(Y);t

21、itle濾波后信號的頻譜');在這里，是采用了 butter命令，設(shè)計出一個巴特沃斯低通濾波器，從頻譜圖中可以 很明顯的看出來。下面，也就是本課題的主要內(nèi)容，也都是運(yùn)用到了butter函數(shù),以便容易的得到系統(tǒng)函數(shù)的分子與分母系數(shù)，最終以此來實現(xiàn)信號的濾波。(2)N階高通濾波器的設(shè)計(在這里，以5階為例，其中WC為其3dB邊緣頻率，以為單位)，程序設(shè)計如下：X=WaVreadcdi ng.wav');SOun d(x);N=5;wc=0.3;b,a=butter(N,wc,'high');X=fft(x);SubPlOt(321);PlOt(x);title濾波前

22、信號的波形');SUbPlot(322);PlOt(X);title(濾波前信號的頻譜');y=filter(b,a,x);Y=fft(y);SUbPlOt(323);plot(y);title('IIR濾波后信號的波形');SUbPlOt(324);PlOt(Y);title('IIR濾波后信號的頻譜');z=fftfilt(b,x);Z=fft(z);SUbPIOt(325);PIOt(Z);title('FIR濾波后信號的波形');得到結(jié)果如圖:(3) N階低通濾波器的設(shè)計(在這里，同樣以 5階為例，其中WC為其3dB邊緣

23、頻率，以為單位)，程序設(shè)計如下：X=WaVreadcdi ng.wav');SOun d(x);N=5;wc=0.3;b,a=butter(N,wc);X=fft(x);SUbPlOt(321);PlOt(x);title濾波前信號的波形');SUbPlot(322);PlOt(X);title(濾波前信號的頻譜');y=filter(b,a,x);Y=fft(y);SUbPlOt(323);plot(y);title('IIR濾波后信號的波形');SUbPlOt(324);PlOt(Y);title('IIR濾波后信號的頻譜');z=f

24、ftfilt(b,x);SUbPlOt(326);PlOt(Z);title('FIR濾波后信號的頻譜');得到結(jié)果如圖：(4) 2N階帶通濾波器的設(shè)計(在這里，以10階為例，其中WC為其3dB邊緣頻率，以 為單位，wc=w1,w2,w1 WC w2 ,程序設(shè)計如下：X=WaVread('di ng.wav');SOun d(x);N=5;wc=0.3,0.6;b,a=butter(N,wc);X=fft(x);SUbPlOt(321);PlOt(x);title濾波前信號的波形');SUbPIot(322);PIot(X);title(濾波前信號的頻譜

25、');y=filter(b,a,x);Y=fft(y);SUbPIOt(323);plot(y);title('IIR濾波后信號的波形');SUbPlOt(324);PlOt(Y);title('IIR濾波后信號的頻譜');z=fftfilt(b,x);Z=fft(z);SUbPlOt(325);PlOt(Z);title('FIR濾波后信號的波形');SUbPlOt(326);PlOt(Z);title('FIR濾波后信號的頻譜');得到結(jié)果如圖：(5) 2N階帶阻濾波器的設(shè)計(在這里，以10階為例，其中WC為其3dB邊

26、緣頻率，以 為單位，wc=w1,w2,w1 WC w2 ,程序設(shè)計如下：X=WaVread('di ng.wav');SOUn d(x);N=5;wc=0.2,0.7;b,a=butter(N,wc,'stop');X=fft(x);SUbPlOt(321);PlOt(x);title濾波前信號的波形');SUbPIot(322);PIot(X);title(濾波前信號的頻譜');y=filter(b,a,x);Y=fft(y);SUbPIOt(323);plot(y);title('IIR濾波后信號的波形');SUbPlOt(3

27、24);PlOt(Y);title('IIR濾波后信號的頻譜');z=fftfilt(b,x);Z=fft(z);SUbPIOt(325);PIOt(Z);title('FIR濾波后信號的波形');SUbPlOt(326);PlOt(Z);title('FIR濾波后信號的頻譜');得到結(jié)果如圖：(6) 小結(jié)：以上幾種濾波，我們都可以從信號濾波前后的波形圖以及頻譜圖上看 出變化。當(dāng)然，也可以用sou nd(函數(shù)來播放濾波后的語音，從聽覺上直接感受語 音信號的變化，但由于人耳聽力的限制，有些情況下我們是很難聽出異同的。同樣，通過函數(shù)的調(diào)用，也可以將信

28、號的頻譜進(jìn)行%26ldquo分離觀察26rdquo;,如顯出信號的幅值或相位。下面，通過改變系統(tǒng)函數(shù)的分子與分母系數(shù)比，來觀察 信號濾波前后的幅值與相位。并且使結(jié)果更加明顯，使人耳得以很容易的辨聽。X=WaVread('di ng.wav');soun d(x);b=100;a=5;y=filter(b,a,x);X=fft(x,4096);SubPlOt(221);PlOt(x);title濾波前信號的波形');SUbPlOt(222);plot(abs(X);title濾波前信號的幅值');Y=fft(y,4096);SubPlOt(223);plot(y)

29、;title濾波后信號的波形');SubPlOt(224);plot(abs(Y);title濾波后信號的幅值');結(jié)果如圖:%26gt%26gt; SOun d(y);可以聽到聲音明顯變得高亢了。從上面的波形與幅值(即幅頻)圖，也可看出，濾 波后的幅值變成了濾波前的20倍。%26gt;%26gt; figure,SUbPlOt(211);plot(angle(X);title濾波前信號相位');SUbPlOt(212);plot(angle(Y);title濾波后信號相位');得圖：可以看到相位譜沒什么變化。(四)、界面設(shè)計：直接用M文件編寫GUl程序很繁瑣，

30、而使用GUIDE設(shè)計工具可以大大提高工作 效率。GUIDE相當(dāng)于一個控制面板，從中可以調(diào)用各種設(shè)計工具以輔助完成界面 設(shè)計任務(wù)，例如控件的創(chuàng)建和布局、控件屬性的編輯和菜單設(shè)計等。使用GUIDE設(shè)計GUI程序的一般步驟如下：1. 將所需控件從控件面板拖拽到 GUIDE的設(shè)計區(qū)域；2. 利用工具條中的工具(或相應(yīng)的菜單和現(xiàn)場菜單)，快速完成界面布局；3. 設(shè)置控件的屬性。尤其是tag屬性，它是控件在程序內(nèi)部的唯一標(biāo)識；4. 如果需要，打開菜單編輯器為界面添加菜單或現(xiàn)場菜單；5. 保存設(shè)計。GUIDE默認(rèn)把GUI程序保存為兩個同名文件：一個是fig文件，用來保存窗體布局和所有控件的界面信息；一個是m文件，該文件的初始內(nèi)容是 GUIDE自動產(chǎn)生的程序框架，其中包括了各個控件回調(diào)函數(shù)的定義。該M文件與一般的M文件沒有本質(zhì)區(qū)別，但是鑒于它的特殊性，MATALAB把這類文件統(tǒng)稱為GU