如何利用matlab處理音頻信號

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

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

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

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

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

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

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

8、7782 - 7.7071i5.00004.7782 + 7.7071i0 + 5.0000i-10.7782 - 6.2929i由程序運行所得結果可見，X(k)和x(n)的維數相同，共有8個元素。X(k)的第一行 元素對應頻率值為0,第五行元素對應頻率值為 NyquiSt頻率，即標準頻率為1.因 此第一行至第五行對應的標準頻率為 01。而第五行至第八行對應的是負頻率， 其X(k)值是以NyqUiSt頻率為軸對稱。(注：通常表示為NyqUiSt頻率外擴展，標以正值。)一般而言，對于N點的x(n)序列的FFT是N點的復數序列，其點n=N2+1對應 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點所對應實際頻率值為f=k*f /N.(3) 作FFT分析時，幅值大小與FFT選擇點數有關，但不影響分析結果。2、設計內容：(1) 下面的一段程序是語音信號在MATLAB中的最簡單表現，它實現了語音的讀 入打開，以及繪出了語音信號的波形頻譜圖。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原始信號相位');程序運行可以聽到聲音，得到的圖形為：(2) 定點分析：已知一個語音信號，數據采樣頻率為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運行結果如圖：上圖(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選用點數N有關，但只要點數N足夠 不影響研究結果。從上圖幅頻譜可見，信號中包括15Hz和40Hz的正弦分量。(3)若信號長度T=25.6s,即抽樣后x(n)點數為TTs=256 ,所得

14、頻率分辨率為HZ ,以此觀察數據長度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)') %數據長度N擴大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)')結果如圖：(三)、濾波器設計：1、相關原理：設計數字濾波器的任務就是尋求一個因果穩(wěn)定的線性時不變系統(tǒng)，并使系統(tǒng)函數 H(Z)具有指定的頻率特性。數字濾波器從實現的網絡結構或者從單位沖激響應分類，可以分成無限長單位沖激 響應(IIR)數字濾波器和有限長單位沖激響應(FIR)數字濾波器。數字濾波器頻率響應的三個參數：(1) 幅度平方響應:(2) 相位響應其中，相位響應(3) 群時延響應IIR數字濾波器：IIR數字濾波器的系統(tǒng)函數為的有理分數，即IIR數字濾波器的逼近問題就是求解濾波器的系數 和，使得在規(guī)定的物理意義上 逼近所要求的特性的問題

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

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

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

19、已知濾波器指標，和，則調用格式為N,wc=buttord(wp,ws,Rp,As)對于不同類型的濾波器，參數 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、設計內容:(1)濾波器示例：在這里為了說明如何用 MATLAB來實現濾波，特舉出一個簡單的函數信號濾波實 例(對信號x(n)=sin( n4)+5cos

20、( n2)進行濾波，信號長度為500點)，從中了解濾 波的實現過程。程序如下：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命令，設計出一個巴特沃斯低通濾波器，從頻譜圖中可以 很明顯的看出來。下面，也就是本課題的主要內容，也都是運用到了butter函數,以便容易的得到系統(tǒng)函數的分子與分母系數，最終以此來實現信號的濾波。(2)N階高通濾波器的設計(在這里，以5階為例，其中WC為其3dB邊緣頻率，以為單位)，程序設計如下：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濾波后信號的波形');得到結果如圖:(3) N階低通濾波器的設計(在這里，同樣以 5階為例，其中WC為其3dB邊緣

23、頻率，以為單位)，程序設計如下：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濾波后信號的頻譜');得到結果如圖：(4) 2N階帶通濾波器的設計(在這里，以10階為例，其中WC為其3dB邊緣頻率，以 為單位，wc=w1,w2,w1 WC w2 ,程序設計如下：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濾波后信號的頻譜');得到結果如圖：(5) 2N階帶阻濾波器的設計(在這里，以10階為例，其中WC為其3dB邊

26、緣頻率，以 為單位，wc=w1,w2,w1 WC w2 ,程序設計如下：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濾波后信號的頻譜');得到結果如圖：(6) 小結：以上幾種濾波，我們都可以從信號濾波前后的波形圖以及頻譜圖上看 出變化。當然，也可以用sou nd(函數來播放濾波后的語音，從聽覺上直接感受語 音信號的變化，但由于人耳聽力的限制，有些情況下我們是很難聽出異同的。同樣，通過函數的調用，也可以將信

28、號的頻譜進行%26ldquo分離觀察26rdquo;,如顯出信號的幅值或相位。下面，通過改變系統(tǒng)函數的分子與分母系數比，來觀察 信號濾波前后的幅值與相位。并且使結果更加明顯，使人耳得以很容易的辨聽。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濾波后信號的幅值');結果如圖:%26gt%26gt; SOun d(y);可以聽到聲音明顯變得高亢了。從上面的波形與幅值(即幅頻)圖，也可看出，濾 波后的幅值變成了濾波前的20倍。%26gt;%26gt; figure,SUbPlOt(211);plot(angle(X);title濾波前信號相位');SUbPlOt(212);plot(angle(Y);title濾波后信號相位');得圖：可以看到相位譜沒什么變化。(四)、界面設計：直接用M文件編寫GUl程序很繁瑣，

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