基于matlab的語音信號的時、頻域分析課程設計

1、<p><b>　　摘要</b></p><p>　　用MATLAB對于語音信號進行分析和處理，采集語音信號后，在MATLAB軟件平臺進行頻譜分析；并對所采集的語音信號加入干擾噪聲，對加入噪聲的信號進行頻譜分析，設計合適的濾波器濾除噪聲，恢復原信號。本次設計介紹了基于MATLAB的對語音信號的采集，處理及濾波器的設計，并使之實現(xiàn)的過程。</p><p>　　

2、關鍵詞：MATLAB;語音信號；濾波器；加噪；除噪</p><p><b>　　目錄</b></p><p><b>　　摘要1</b></p><p><b>　　1 設計原理1</b></p><p>　　1.1 理論原理1</p><p>

3、;　　1.1.1采樣頻率1</p><p>　　1.1.2采樣位數1</p><p>　　1.1.3采樣定理1</p><p>　　1.1.4時域信號的FFT分析2</p><p>　　1.1.5數字信號的濾波器原理和方法2</p><p>　　1.1.6 各種不同類型濾波器的性能比較3</p>

4、;<p>　　1.1.7離散傅立葉變換3</p><p>　　2 信號采集及讀取4</p><p>　　3 構造受干擾信號并對其進行FFT頻譜分析........................…...6</p><p>　　4 數字濾波器設計8</p><p>　　4.1 用窗函數法設計IIR帶通濾波器8<

5、/p><p>　　4.2 用窗函數法設計FIR低通濾波器10</p><p><b>　　5信號處理12</b></p><p>　　5.1 IIR帶通濾波12</p><p>　　5.2 FIR低通濾波14</p><p><b>　　6心得體會16</b></

6、p><p>　　7 參考文獻：16</p><p><b>　　1設計原理</b></p><p><b>　　1.1 理論原理</b></p><p>　　1.1.1采樣頻率</p><p>　　也稱為采樣速度或者采樣率，定義了每秒從連續(xù)信號中提取并組成離散信號的采樣個數，

7、它用赫茲（Hz）來表示。采樣頻率只能用于周期性采樣的采樣器，對于非周期性采樣的采樣器沒有規(guī)則限制。通俗的講采樣頻率是指計算機每秒鐘采集多少個聲音樣本，是描述聲音文件的音質、音調，衡量聲卡、聲音文件的質量標準。采樣頻率越高，即采樣的間隔時間越短，則在單位時間內計算機得到的聲音樣本數據就越多，對聲音波形的表示也越精確。</p><p><b>　　1.1.2采樣位數</b></p>

8、<p>　　即采樣值或取樣值，用來衡量聲音波動變化的參數，是指聲卡在采集和播放聲音文件時所使用數字聲音信號的二進制位數。</p><p><b>　　1.1.3采樣定理</b></p><p>　　在進行模擬/數字信號的轉換過程中，當采樣頻率fs.max大于信號中最高頻率fmax的2倍時(fs.max>=2fmax)，采樣之后的數字信號完整地保留了原

9、始信號中的信息，一般實際應用中保證采樣頻率為信號最高頻率的5～10倍；采樣定理又稱奈奎斯特定理。</p><p>　　1.1.4時域信號的FFT分析</p><p>　　信號的頻譜分析就是計算信號的傅立葉變換。連續(xù)信號與系統(tǒng)的傅立葉分析顯然不便于用計算機進行計算，使其應用受到限制。而FFT是一種時域和頻域均離散化的變換，適合數值運算，成為用計算機分析離散信號和系統(tǒng)的有力工具。對連續(xù)信號和系

10、統(tǒng)，可以通過時域采樣，應用DFT進行近似譜分析。</p><p>　　1.1.5數字信號的濾波器原理和方法</p><p>　　IIR數字濾波器系統(tǒng)函數：</p><p>　　H(z)=···············&

11、#183;（1.1）</p><p>　　其中H(z)稱為N階IIR數字濾波器系統(tǒng)函數。IIR濾波器設計方法有間接和直接法，間接法是借助于模擬濾波器的設計方法進行的。其步驟是：先設計過度模擬濾波器得到系統(tǒng)函數Ha(s)，然后將Ha（s）按某種方法轉換成數字濾波器的系統(tǒng)函數H(z)。</p><p>　　利用有限脈沖響應(FIR)設計濾波器。有限脈沖響應濾波器在保證幅度特性滿足技術要求的同時

12、，很容易做到有嚴格的線性相位特性。用N表示FIR濾波器單位脈沖響應h(n)的長度，其系統(tǒng)函數H（z）為</p><p>　　··················（1.2）</p><p>　　H（z）是z-1的N

13、-1次多項式，它在z平面上有N-1個零點，在原點z=0處有一個N-1重極點。因此，H(z)永遠穩(wěn)定。穩(wěn)定和線性相位是FIR濾波器最突出的優(yōu)點。</p><p>　　1.1.6 各種不同類型濾波器的性能比較</p><p>　　巴特沃斯濾波器具有單調下降的幅頻特性；切比雪夫濾波器的幅頻特性在通帶或阻帶有波紋特性，可以提高選擇性；貝塞爾濾波器通帶內有較好的線性相位特性；橢圓濾波器的選擇性相對前

14、三種是最好的，但通帶和阻帶內均呈現(xiàn)等波紋幅頻特性，線性特性的非線性稍重。</p><p>　　IIR數字濾波器最大的優(yōu)點是給定一組指標時，它的階數要比相同組的FIR濾波器低得多。IIR數字濾波器的設計方法是利用模擬濾波器成熟的理論及設計圖進行的，因而保留了一些典型模擬濾波器優(yōu)良的幅度特性。但設計中只考慮了幅度特性，沒考慮相位特性，所設計的濾波器一般是某種確定的非線性相位特性。為了得到線性相位特性，對IIR濾波器必

15、須另外增加相位校正網絡。FIR濾波器在保證幅度特性滿足技術要求的同時，很容易做到有嚴格的線性相位特性。穩(wěn)定和線性相位是FIR濾波器最突出的優(yōu)點。</p><p>　　1.1.7離散傅立葉變換</p><p><b>　　1.3</b></p><p>　　其中,N為DFT變換的區(qū)間長度。</p><p>　　2 信號采

16、集及讀取</p><p>　　用錄音工具錄一段自己的聲音，轉化為wav格式的文件，然后保存到MATLAB的文件夾下的work文件夾中。在MATLAB軟件平臺下，利用wavread函數對語音信號進行采樣，記住采樣頻率和采樣點數。</p><p>　　Wavread函數調用格式：</p><p>　　y=wavread(file),讀取file所規(guī)定的wav文件，返回采

17、樣值放在向量y中。</p><p>　　[y,fs,nbits]=wavread(‘1.wav’)；把語音信號進行載MATLAB仿真軟件平臺中，采樣值放在向量y中，fs表示采樣頻率（hz），nbits表示采樣位數。</p><p>　　sound（y,fs,nbits）；回放語音信號</p><p>　　N=length(y)；求出語音信號的長度</p>

18、<p>　　Y=fft(y,N)；傅立葉變換</p><p>　　調用函數y為被變換的時域系列向量，N是DFT變換區(qū)間長度，當N大于y的長度時，fft函數自動在y后面補零，函數返回y的N點DFT變換結果向量Y。當N小于y的長度時，fft函數計算y的前面N個元素構成的N點DFT，忽略y后面的元素。首先畫出語音信號的時域波形，然后對語音信號進行頻譜分析。在MATLAB中利用fft對新號進行快速傅立葉變換

19、，得到信號的頻譜特性。</p><p>　　subplot(2,1,1);plot(y);title(‘原始信號波形’);</p><p>　　subplot(2,1,2);plot(abs(Y));title(‘原始信號頻譜’)；</p><p><b>　　程序如下：</b></p><p>　　[y,fs,nbit

20、s]=wavread('1.wav');</p><p>　　sound(y,fs,nbits);</p><p>　　N=length(y);</p><p>　　Y=fft(y,N);</p><p>　　subplot(2,1,1);</p><p><b>　　plot(y);<

21、/b></p><p>　　title('原始信號波形');xlabel('時間');ylabel('幅度');</p><p>　　subplot(2,1,2);</p><p>　　plot(abs(Y));</p><p>　　title('原始信號頻譜');xlab

22、el('時間');ylabel('幅度');</p><p><b>　　程序結果如下圖：</b></p><p><b>　　程序如下：</b></p><p>　　[y,fs,nbits]=wavread ('1.wav');</p><p>　　

23、N = length (y) ; %求出語音信號的長度</p><p>　　noise=0.2*randn(N,1); %隨機函數產生噪聲</p><p>　　Si=y+noise; %語音信號加入噪聲 </p><p>　　sound(Si);</p>&

24、lt;p>　　subplot(2,1,1);</p><p><b>　　plot(Si);</b></p><p>　　title('加噪語音信號的時域波形');xlabel('時間');ylabel('幅度');</p><p>　　S=fft(Si);

25、 %傅里葉變換</p><p>　　subplot(2,1,2);</p><p>　　plot(abs(S)); </p><p>　　title('加噪語音信號的頻域波形'); xlabel('時間');ylabel('幅度');</p><p><b>　　程序結

26、果如下圖：</b></p><p><b>　　4 數字濾波器設計</b></p><p>　　4.1 用窗函數法設計IIR帶通濾波器</p><p><b>　　程序如下：</b></p><p>　　wp=0.2*pi; ws=0.3*pi; %性能指標</p>

27、<p>　　wdelta=ws-wp; %過渡帶寬度</p><p>　　M=ceil(3.32*pi/wdelta); %濾波器長度，朝正無窮方向舍入</p><p>　　N=2*M+1; %窗口長度</p><p>　　wc=(ws+wp)/2; %截止頻率</p>&

28、lt;p>　　win=hamming(N); %因為衰減不小于50dB，所以選擇海明窗，這里得到海明窗的時域響應</p><p>　　b=fir1(N-1,wc/pi,win);</p><p><b>　　n=0:1:N;</b></p><p>　　[hi t]=impz(b,1,n);%得到脈沖響應</p>

29、<p>　　[hf w]=freqz(b,1,512); %得到頻率響應</p><p>　　subplot(3,1,1); stem(n,hi);</p><p>　　xlabel('n'); ylabel('幅度'); title('脈沖響應');</p><p>　　subplot(3,1,2);

30、 plot(w/pi,20*log10(abs(hf)));</p><p>　　xlabel('頻率(Hz)'); ylabel('幅度(dB)');</p><p>　　title('頻率響應');</p><p>　　subplot(3,1,3); plot(w/pi,180/pi*unwrap(angle(hf

31、)));</p><p>　　xlabel('頻率(Hz)'); ylabel('相位(degrees)');</p><p>　　title('頻率響應');</p><p><b>　　程序結果如下圖：</b></p><p>　　4.2 用窗函數法設計FIR低通濾波器

32、</p><p><b>　　程序如下：</b></p><p><b>　　Ft=8000;</b></p><p><b>　　Fp=1000;</b></p><p><b>　　Fs=1200;</b></p><p>

33、　　wp=2*Fp/Ft;</p><p>　　ws=2*Fs/Ft;</p><p><b>　　rp=1;</b></p><p><b>　　rs=50;</b></p><p>　　p=1-10.^(-rp/20);</p><p>　　s=10.^(-rs/20);

34、</p><p>　　fpts=[wp,ws];</p><p>　　mag=[1,0];</p><p>　　dev=[p,s];</p><p>　　[n21,wn21,beta,ftype]=kaiserord(fpts,mag,dev);</p><p>　　b21=fir1(n21,wn21,kaiser(n

35、21+1,beta));</p><p>　　[h,w]=freqz(b21,1);</p><p>　　plot(w/pi,abs(h));</p><p>　　title('FIR低通濾波器');</p><p><b>　　grid;</b></p><p><b>

36、;　　其程序結果如下圖：</b></p><p><b>　　5信號處理</b></p><p>　　用設計好的數字濾波器對含噪聲的語音信號進行濾波，在MATLAB中FIR濾波器利用函數fftfilt對信號進行濾波，IIR濾波器利用函數filter對信號進行濾波，在一個窗口同時畫出濾波前后的波形及頻譜如下圖。</p><p>　　5

37、.1 IIR帶通濾波</p><p><b>　　程序如下：</b></p><p>　　[y,fs,nbits]=wavread('1.wav');</p><p>　　N=length(y);</p><p>　　Noise=0.2*randn(N,1);</p><p>　　

38、Si=y+Noise; %語音信號加入噪聲 </p><p>　　nn=0:1:50;</p><p>　　y=filter(b,1,Si);</p><p>　　subplot(2,1,1); stem(Si);</p><p>　　xlabel('n'); ylabel('幅度');

39、</p><p>　　title('濾波前信號頻譜'); axis([0 50 -2 2]);</p><p>　　subplot(2,1,2); stem(out);</p><p>　　xlabel('n'); ylabel('幅度');</p><p>　　title('濾波后信號

40、頻譜'); axis([0 50 -1 1]);</p><p><b>　　sound(y);</b></p><p><b>　　程序結果如下圖：</b></p><p>　　5.2 FIR低通濾波</p><p><b>　　程序如下：</b></p>

41、<p>　　[y,fs,nbits]=wavread('1.wav');</p><p>　　N = length (y) ; </p><p>　　noise=0.2*randn(N,1); </p><p>　　Si=y+noise; </p><p>　　S=

42、fft(Si);</p><p>　　z21=fftfilt(b21,Si);</p><p>　　sound(z21);</p><p>　　m21=fft(z21);</p><p>　　subplot(2,2,1);</p><p>　　plot(abs(S),'g');</p>&

43、lt;p>　　title('濾波前信號的頻譜');xlabel('時間');ylabel('幅度');</p><p><b>　　grid;</b></p><p>　　subplot(2,2,2);</p><p>　　plot(abs(m21),'r');</p&

44、gt;<p>　　title('濾波后信號頻譜');xlabel('時間');ylabel('幅度');</p><p><b>　　grid;</b></p><p>　　subplot(2,2,3);</p><p><b>　　plot(S);</b>&l

45、t;/p><p>　　title('濾波前信號波形');xlabel('時間');ylabel('幅度');</p><p><b>　　grid;</b></p><p>　　subplot(2,2,4);</p><p>　　plot(z21);</p>&l

46、t;p>　　title('濾波后信號波形');xlabel('時間');ylabel('幅度');</p><p><b>　　grid;</b></p><p><b>　　程序結果如下圖：</b></p><p><b>　　6心得體會</b>

47、</p><p>　　通過這次課程設計，我學到了很多東西。進一步加深了對MATLAB的了解。在設計的過程中我遇到了不少的問題，比如剛開始要處理的是WAV格式的語音信號，可是我發(fā)現(xiàn)不是所有的WAV文件都可以用WAVREAD函數讀取，通過軟件轉化了才可以使用。在編寫程序的時候遇到了很多錯誤，最后經過查找一些資料和自己的理解都一一糾正了。總的來說，通過這次的課程設計我對語音信號有的更多的認識，對數字信號處理的知識有了更

48、多的理解，鍛煉了我的理論與實際相結合的能力，并且熟練的應用MATLAB也加深了我對課程的理解。這次設計使我更深刻的了解MATLAB的使用方法，同時我相信，進一步加強對MATLAB的學習對我今后的學習將會起到很大的幫助。</p><p><b>　　7 參考文獻：</b></p><p>　　高西全，丁玉美編著?！稊底中盘柼幚怼肺靼?，西安電子科技大學出版社，2008&l