基于短時傅里葉變換的語音信號時頻分析

基于STFT的語音信號時頻分析摘要：視頻分析是近年來信號處理的新熱點。本文首先介紹了語音信號STFT的相關知識，隨后利用MATLAB將采集到的語音信號進行處理，并進行了信號時域和頻域的相關分析。關鍵詞：語音信號STFT時頻分析語音信號的短時傅里葉變換傅里葉變換是一種信號的整體變換，要么完全在時域，要么完全在頻域進行分析處理，無法給出信號的頻譜如何隨時間變化的規(guī)律。而有些信號，例如語音信號，它具有很強的時變性，在一段時間內呈現(xiàn)出周期性信號的特點，而在另一段時間內呈現(xiàn)出隨機信號的特點，或者呈現(xiàn)出兩個混合的特性。對于頻譜隨時間變化確實定性信號以及非平穩(wěn)隨機信號，利用傅里葉變換分析方法有很大的局限性，或者說是不適宜的。傅里葉變換無法針對性的分析相應時間區(qū)域內信號的頻率特征?？梢杂靡粋€窗函數與時間信號相乘積，當該窗函數的時寬足夠窄，使取出的信號可以被看成是平穩(wěn)信號時，就可以對乘積信號進行傅里葉變換，從而反映該時寬中的信號頻譜變化規(guī)律。早在1946年，Gabor就提出了短時傅立葉變換〔ShortTimeFourierTransform，STFT〕的概念，用以測量聲音信號的頻率定位[64]。 給定一信號，其STFT定義為 式中 及 ，并且窗函數應取對稱函數。STFT的含義可解釋如下：在時域用窗函數去截〔注：將，的時間變量換成〕，對截下來的局部信號作傅立葉變換，即得在時刻得該段信號得傅立葉變換。不斷地移動，也即不斷地移動窗函數的中心位置，即可得到不同時刻的傅立葉變換。這些傅立葉變換的集合，即是，如圖1所示。顯然，是變量的二維函數。 設語音波形分幀處理后得到的第n幀語音信號為Xn(m)，那么Xn(m)滿足下式：其中，n＝0，1T，2T，…，并且N為幀長，T為幀移長度。某一幀的短時傅里葉變換的定義如下：式中w(n-m)是窗函數。不同的窗函數，可得到不同的傅里葉變換的結果?？梢钥闯龆虝r傅里葉變換有兩個變量，即離散時間n及連續(xù)頻率w。假設令，那么可得到離散的短時傅里葉變換如下：它實際上就是頻率抽樣。將上述某一幀語音信號的傅里葉變換寫為可以看出時變傅里葉變換是時間標號n的函數，當n變化時，窗函數w(n-m)沿著x(m)滑動。圖1窗函數w(n-m)沿著x(m)滑動可以得出結論：短時傅里葉變換實際就是窗選語音信號的標準傅里葉變換。這里，窗w(n-m)是一個“滑動的〞窗口，它隨n的變化而沿著序列X(n)滑動。由于窗口是有限長度的，滿足絕對可和條件，所以這個變換是存在的。當然窗口函數不同，博里葉變換的結果也將不同。對于w(n-m)窗來說，它除了具有選出x(m)序列中被分析局部作用外，其形狀對時變傅里葉變換的特性也具有重要作用，從標準傅里葉變換可以方便的解釋這種作用。如果被看成是w(n-m)x(m)序列的標準傅里葉變換，同時假設x(m)及w(m)的標準傅里葉變換存在，即：當n固定時，序列w(n-m)的傅里葉變換為根據卷積定理，有：因為上式右邊兩個卷積項均為關于角頻率w的以2π為周期的連續(xù)函數，所以也可將其寫成以下的卷積積分形式：假設x(m)的DTFT是，且的DTFT是，那么是和的周期卷積。根據信號的時寬帶寬的積為一常數這一根本性質，可知主瓣寬度與窗口寬度成反比，N越大，的主瓣越窄。為了使忠實再現(xiàn)的特性，相對于來說必須是—個沖激函數。所以為了使，需；但是N值太大時，信號的分幀又失去了意義。尤其是N大于語音的音素長度時，已不能反映該語音音素的頻譜了。因此，應折衷選擇窗的寬度N。另外，窗的形狀也對短時博氏頻譜有影響，如矩形窗，雖然頻率分辨率很高(即主辯狹窄鋒利)，但由于第一旁瓣的衰減很小，有較大的上下沖，采用矩形窗時求得的與的偏差較大，這就是Gibbs效應，所以不適合用于頻譜成分很寬的語音分析中。而漢明窗在頻率范圍中的分辨率較高，而且旁瓣的衰減大，具有頻譜泄漏少的優(yōu)點。所以在求短時頻譜時一船采用具有較小上下沖的漢明窗。 短時傅里葉變換有以下性質：(1)時移性設那么(2)頻移性設，那么短時傅里葉變換具有頻移不變性，不具有時移不變性.但是，它在某一調制范圍內即相差一相位因子的范圍內保持時移不變性.(3)假設綜合窗函數滿足“完全重構條件〞：那么可由逆變換完全重構，即語音信號的時域分析語音信號的時域分析就是分析和提取語音信號的時域參數。進行語音分析時，最先接觸到并且也是最直觀的是它的時域波形。語音信號本身就是時域信號，因而時域分析是最早使用，也是應用最廣泛的一種分析方法，這種方法直接利用語音信號的時域波形。時域分析通常用于最根本的參數分析及應用，如語音的分割、預處理、大分類等。這種分析方法的特點是：=1\*GB3①表示語音信號比擬直觀、物理意義明確。=2\*GB3②實現(xiàn)起來比擬簡單、運算量少。=3\*GB3③可以得到語音的一些重要的參數。④只使用示波器等通用設備，使用較為簡單等。MATLAB數據采集箱中提供的進行語音信號分析的函數命令如下：wavread:wavread用于讀取擴展名為“.wav〞的聲音文件。其調用形式為:y=wavread(‘filename’)。其作用是讀取wave文件，將讀取的采樣數據送到y(tǒng)中。sound:音頻信號是以向量的形式表示聲音采樣的。sound函數用于將向量轉換為聲音，其調用形式為:sound(y,fs)，作用是向揚聲器送出向量y中的音頻信號(采樣頻率為fs)。通過wavread和plot()函數即可顯示語音信號的時域波形。如圖2所示。圖2語音信號的時域波形貫穿于語音分析全過程的是“短時分析技術〞。因為，語音信號從整體來看其特性及表征其本質特征的的參數均是隨時間而變化的，所以它是一個非平穩(wěn)過程，不能用處理平穩(wěn)信號的數字信號處理技術對其進行分析處理。但是，由于不同的語音是由人的口腔肌肉運動產生的，相對于語音頻率來說是非常緩慢的，所以從另一方面來看，雖然語音信號具有時變特性，但是在一個短時間范圍內〔一般可認為在10~30ms的短時間內〕，其特性根本保持不變即相對穩(wěn)定，因而可以將其看作是一個準穩(wěn)態(tài)過程，即語音信號具有短時平穩(wěn)性。所以任何語音信號的分析和處理必須建立在“短時〞的根底上，即進行“短時分析〞，將語音信號分為一段一段來分析其特征參數，其中每一段稱為一“幀〞，幀長一般取為10~30ms。這樣，對于整體的語音信號來講，分析出的是由每一幀特征參數組成的特征參數時間序列。一般每秒的幀數約為33~100幀，視實際情況而定。分幀雖然可以采用連續(xù)分段的方法，但一般要采用交疊分段的方法，這是為了使幀與幀之間平滑過渡，保持其連續(xù)性。前一幀和后一幀的交疊局部稱為幀移。幀移與幀長的比值一般取為0~1/2。分幀是用可移動的有限長度窗口進行加權的方法來實現(xiàn)的，這就是用一定的窗函數w(n)來乘s(n)，從而形成加窗語音信號Sw(n)=s(n)w(n)。在語音信號數字處理中常用的窗函數是矩形窗和漢明窗等，它們的表達式如下〔其中N為幀長〕：矩形窗：漢明窗：這兩種窗函數都有低通特性，通過分析這兩種窗的頻率響應幅度特性可以發(fā)現(xiàn)：矩形窗的主瓣寬度小，具有較高的頻率分辨率，旁瓣峰值大，會導致泄漏現(xiàn)象；漢明窗的主瓣寬8*pi/N，旁瓣峰值低，可以有效的克服泄漏現(xiàn)象，具有更平滑的低通特性。因此在語音頻譜分析時常使用漢明窗，在計算短時能量和平均幅度時通常用矩形窗。語音信號的時域參數有短時能量、短時過零率、短時自相關函數和短時平均幅度差函數等，這是語音信號的一組最根本的短時參數，在各種語音信號數字處理技術中都要應用。短時平均過零率分析是時域分析中最為簡單的方法，即單位時間內時域波形穿過坐標軸的次數，用來區(qū)別濁音還是清音。語音中濁音過零率較低。清音的過零率較高。但由于這種上下只是相對而言，難以細化，因此用過零率來描述語言信號頻率就比擬粗略。短時平均能量分析是另一種常用的時域分析方法。它也可以區(qū)分濁音和清音，清音的能量較小，濁音的能量較大。在高信噪比條件下可以區(qū)分有無語音。短時自相關分析也是視頻分析中比擬重要的分析方法，常用語對基音周期的估計。在決定基音周期時，利用了短時自相關函數基音周期各整數倍的點上具有較高的峰起值。只要找到第一大峰值點的位置并計算它與原點的間隔，便能估計出基音周期。本作業(yè)只對短時能量進行簡單分析。由于語音信號的能量隨時間變化，清音和濁音之間的能量差異相當顯著。因此對語音的短時能量進行分析，可以描述語音的這種特征變化情況。定義n時刻某語音信號的短時平均能量En為：其中N為窗長,可見短時能量為一幀語音樣值點的加權平方和。特殊地，當窗函數為矩形窗時，可簡化為：本次語音信號在矩形窗長N=400時的短時能量如以下圖3所示：圖3N=400時的短時能量短時平均能量的主要用途如下：①可以作為區(qū)分清音和濁音的特征參數。實驗結果說明濁音的能量明顯高于清音。通過設置一個能量門限值，可以大致判定濁音變?yōu)榍逡艋蛘咔逡糇優(yōu)闈嵋舻臅r刻，同時可以大致劃分濁音區(qū)間和清音區(qū)間。②在信噪比擬高的情況下，短時能量還可以作為區(qū)分有聲和無聲的依據。③可以作為輔助的特征參數用于語音識別中。語音信號的頻域分析語音信號的頻域分析就是分析語音信號的頻域持征。從廣義上講，語音信號的頻域分析包括語音信號的頻譜、功率譜、倒頻譜、頻譜包絡分析等，而常用的頻域分析方法有帶通濾波器組法、傅里葉變換法、線件預測法等幾種。本文介紹的是語音信號的傅里葉分析法。因為語音波是一個非平穩(wěn)過程，因此適用于周期、瞬變或平穩(wěn)隨機信號的標準傅里葉變換不能用來直接表示語音信號，而應該用短時傅里葉變換對語音信號的頻譜進行分析，相應的頻譜稱為“短時譜〞。在MATLAB的信號處理工具箱中函數FFT和IFFT用于快速傅立葉變換和逆變換。函數FFT用于序列快速傅立葉變換，其調用格式為y=fft(x)，其中，x是序列，y是序列的FFT，x可以為一向量或矩陣，假設x為一向量，y是x的FFT且和x相同長度；假設x為一矩陣，那么y是對矩陣的每一列向量進行FFT。函數FFT的另一種調用格式為y=fft(x,N)，式中，x，y意義同前，N為正整數。函數執(zhí)行N點的FFT，假設x為向量且長度小于N，那么函數將x補零至長度N；假設向量x的長度大于N，那么函數截短x使之長度為N；假設x為矩陣，按相同方法對x進行處理。利用上述函數即可畫出語音信號的頻譜圖，如以下圖4所示。圖4語音信號的頻譜圖從圖中可以確定語音參數，例如共振峰頻率及基頻。條紋的起點相當于聲門脈沖的起點，條紋之間的距離表示基音周期。條紋越密表示基音頻率越高，從圖4中可以看出，當發(fā)濁音時，條紋相對密集些，其基音頻率可以到達最大，說明濁音比清音的基音頻率高。共振蜂表示基音脈沖的某些頻率成分被加強，這在頻譜圖上呈現(xiàn)為條紋區(qū)更寬更黑。參考文獻[1]韓紀慶，張磊，鄭鐵然.語音信號處理[M].清華大學出版社，2004[2]易克初，田斌，付強.語音信號處理[M].國防工業(yè)出版社，2000[3]周輝，董正宏.數字信號處理根底及MATLAB實現(xiàn)[M].北京希望電子出版社，2006[4]鄒理和.語音信號處理[M].北京：國防工業(yè)出版社，1985[5]叢玉良,王宏志.數字信號處理原理及其MATLAB實現(xiàn)[M].北京：電子工業(yè)出版社，2005[6]何強，何英.MATLAB擴展編程[M].清華大學出版社，2002[7]王世一.數字信號處理[M].北京：北京理工大學出版社，2005附錄程序closeall; clearall;[x1]=wavread('xinhao.wav');%讀取語音信號的數據，賦給變量x1figure(1);plot(x1);%做原始語音信號的時域圖形title('語音信號時域波形');xlabel('采樣點');ylabel('幅值');n=length(x1);y1=fft(x1,n);%對信號做FFT變換figure(2);plot(abs(y1));%原始語音信號的頻域圖形axis([0150000250])title('語音信號頻譜圖');xlabel('Hz');ylabel('幅值');用MATLAB計算語音信號的短時平均能量x=wavread('xinhao.wav');s=fra(400,400,x);%對輸入信號進行分幀s2=s.^2;%一幀內各樣點的能量enery=sum(s2,2);%求一幀能量plot(enery);%畫出N=50時的語音能量圖xlabel('幀數');%橫坐標ylabel('短時能量E');%縱坐標legend('N=400');%曲線標識axis([0,190,5,15])其中fra()為分幀函數，其MATLAB程序如下：%fra.mfunctionf=fra