基于短時(shí)自相關(guān)函數(shù)的基因周期檢測

1、基音周期的檢測什么是基音周期檢測 人在發(fā)音時(shí)，根據(jù)聲帶是否震動(dòng)可以將語音信號分為清音跟濁音兩種。濁音又稱有聲語言，攜帶者語言中大部分的能量，濁音在時(shí)域上呈現(xiàn)出明顯的周期性；而清音類似于白噪聲，沒有明顯的周期性。發(fā)濁音時(shí)，氣流通過聲門使聲帶產(chǎn)生張弛震蕩式振動(dòng)，產(chǎn)生準(zhǔn)周期的激勵(lì)脈沖串。這種聲帶振動(dòng)的頻率稱為基音頻率，相應(yīng)的周期就成為基音周期。 通常，基音頻率與個(gè)人聲帶的長短、薄厚、韌性、勁度和發(fā)音習(xí)慣等有關(guān)系，在很大程度上反應(yīng)了個(gè)人的特征。此外，基音頻率還跟隨著人的性別、年齡不同而有所不同。一般來說，男性說話者的基音頻率較低，大部分在70200Hz的范圍內(nèi)，而女性說話者和小孩的基音頻率相對較高，在

2、200450Hz之間。 目前，基音檢測的算法有很多種，常用的檢測算法有自相關(guān)法、平均幅度差函數(shù)法、并行處理法、倒譜法、簡化逆濾波法等。 本次基音周期的檢測，我們選用的是短時(shí)自相關(guān)函數(shù)法，包括四個(gè)模塊。 第一個(gè)模塊為基音的端點(diǎn)檢測，主要為了區(qū)分濁音和清音。第二個(gè)模塊為基音檢測中的帶通濾波器，主要為了減少共振峰的干擾。第三個(gè)模塊為短時(shí)自相關(guān)函數(shù)法做基音檢測，主要為了計(jì)算出基音周期。第四個(gè)模塊為平滑處理，主要為了消除偏離值點(diǎn)?；糁芷诘臋z測（端點(diǎn)檢測）1.設(shè)計(jì)數(shù)據(jù)文件名稱:audioread功能：讀入以.wav為擴(kuò)展名的音頻文件。調(diào)用格式：y,fs,nbits=audioread(wavFilena

3、me)說明：輸出參數(shù)Y是數(shù)字化的音頻；fs是信號采樣頻率；nbits是每個(gè)樣點(diǎn)在編碼時(shí)的二進(jìn)制位數(shù)。2.分幀處理先定義一個(gè)幀數(shù)*幀長的數(shù)組，然后計(jì)算幀的偏移量來，按語音數(shù)據(jù)的位置把x放入y中。名稱：enframe功能：把語音信號按幀長和幀移進(jìn)行分幀調(diào)用格式：f=enframe(x,win,inc)3.端點(diǎn)檢測名稱：pitch_vad1功能：用能熵比法進(jìn)行端點(diǎn)檢測。說明：輸入?yún)?shù)y是分幀后的數(shù)組，一般是一列表示一幀數(shù)據(jù);fn是信號的總幀數(shù)；T1是一個(gè)閾值；輸出參數(shù)voiceseg是語音端點(diǎn)的信息；vosl是幾個(gè)有話段，當(dāng)該幀為有話段時(shí)SF=1；否則為0，Ef是能熵比值。pitch_vad1for

4、i=1:fnSp=abs(fft(y(:,i);%FFT取幅值(:表示所有y(:,1)表示第一列)Sp=Sp(1:wlen/2+1); %只取正頻率部分Esum(i)=sum(Sp.*Sp);%計(jì)算能量值（能量放入Esum里）prob=Sp/(sum(Sp); %計(jì)算概率H(i)=-sum(prob.*log(prob+eps);%求譜熵值（eps表示很小的數(shù)，避免為零）（譜熵放入H里）end譜熵值很類似與過零率值，在有話段內(nèi)的譜熵值小于噪聲段的譜熵值，所以如同能零比一樣，用能熵比更能突出噪聲段的差別?；糁芷诘臋z測（基音檢測中的帶通濾波器）什么是帶通濾波器 帶通濾波器（band-passfi

5、lter）是一個(gè)允許特定頻段的波通過同時(shí)屏蔽其他頻段的濾波器。帶通濾波器的作用 在使用相關(guān)法和AMDF法的基音檢測之前常用到低通濾波器和帶通濾波器，其主要目的是減少共振峰的干擾。我用了500Hz作為濾波器的上限頻率，選擇截止頻率高不利于減少噪聲和共振峰的影響。所以在本次的基音檢測中的預(yù)濾波器選擇的帶寬為60500Hz，高頻截止頻率選擇500Hz，是因?yàn)榛l區(qū)間的高端就在這個(gè)區(qū)域中，低頻截止頻率選擇60Hz是為了減少工頻和低頻噪聲的干擾。什么是共振峰 共振峰是指在聲音的頻譜中能量相對集中的一些區(qū)域，共振峰不但是音質(zhì)的決定因素，而且反映了聲道（共振腔）的物理特征。 聲音在經(jīng)過共振腔時(shí)，受到腔體的濾

6、波作用，使得頻域中不同頻率的能量重新分配，一部分因?yàn)楣舱袂坏墓舱褡饔玫玫綇?qiáng)化，另一部分則受到衰減。由于能量分布不均勻，強(qiáng)的部分猶如山峰一般，故而稱之為共振峰。為什么要減小共振峰 聲道特性與基音周期估計(jì)有關(guān)，會對基音周期檢測造成一定干擾。從語音信號中去除聲道的影響，直接取出僅與聲帶振動(dòng)有關(guān)的聲源信息并非易事。聲道共振峰可能強(qiáng)烈改變聲門波形的結(jié)構(gòu)，從而嚴(yán)重影響激勵(lì)信號的諧波結(jié)構(gòu)，給基因檢測造成困難。我選用IIR濾波器中的橢圓濾波器，因?yàn)镮IR濾波器的運(yùn)算量比FIR少，當(dāng)然IIR濾波器會帶來延遲，也就是相位的變化，但語音信號是對相位不敏感的信號；又選用橢圓濾波器，因?yàn)樗诮?jīng)典濾波器設(shè)計(jì)中相同過渡帶和

7、帶寬條件下，需要的階數(shù)比較小。橢圓濾波器（Ellipticfilter）是在通帶和阻帶等波紋的一種濾波器。橢圓濾波器在通帶和阻帶內(nèi)都具有等波紋幅頻響應(yīng)特性。由于其極點(diǎn)位置與經(jīng)典場論中的橢圓函數(shù)有關(guān)，所以由此取名為橢圓濾波器。對于給定的階數(shù)和波紋要求，橢圓濾波器能得到較其它濾波器更窄的過渡帶寬，可以獲得對理想濾波器幅頻響應(yīng)的最好逼近，是一種性價(jià)比很高的濾波器，橢圓濾波器相比其他類型的濾波器，在階數(shù)相同的條件下有著最小的通帶和阻帶波動(dòng)。濾波器的要求為采樣頻率8000Hz，通帶是60500Hz，通帶波紋為1dB，阻帶分別為30Hz和2000Hz，阻帶衰減為40dB。 濾波器的階數(shù)，就是指過濾諧波的次

8、數(shù)，一般來講，同樣的濾波器，其階數(shù)越高，濾波效果就越好，但是，階數(shù)越高，成本也就越高，因此，選擇合適的階數(shù)是非常重要的。所以本次實(shí)驗(yàn)選用橢圓6階帶通濾波器。實(shí)現(xiàn)代碼%pr8_1_1clearall;clc;closeall;fs=8000;fs2=fs/2;%采樣頻率Wp=60500/fs2;%濾波器通帶Ws=202000/fs2;%濾波器阻帶Rp=1;Rs=40;%通帶的波紋和阻帶的衰減n,Wn=ellipord(Wp,Ws,Rp,Rs);%計(jì)算濾波器的階數(shù)b,a=ellip(n,Rp,Rs,Wn);%計(jì)算濾波器的系數(shù)fprintf(b=%5.6f%5.6f%5.6f%5.6f%5.6f%5

9、.6f%5.6fn,b)fprintf(a=%5.6f%5.6f%5.6f%5.6f%5.6f%5.6f%5.6fn,a)代碼理解 Ellipord函數(shù)的功能是求濾波器的最小階數(shù)，其調(diào)用格式為n,Wn=ellipord(Wp,Ws,Rp,Rs,s) 其中各參量分別為：n-橢圓濾波器最小階數(shù)；Wp-橢圓濾波器通帶；Ws-橢圓濾波器阻帶；Rp-通帶波紋（dB）；Rs-阻帶衰減(dB)； Ellip函數(shù)的功能是用來設(shè)計(jì)橢圓濾波器，其調(diào)用格式： b,a=ellip(n,Rp,Rs,Wp) 由程序計(jì)算出濾波器系數(shù)為： b=0.012280-0.0395080.0421770.000000-0.04217

10、70.039508 a=1.000000-5.52714612.854342-16.11030711.479789-4.4101790.713507基音周期的檢測（短時(shí)自相關(guān)函數(shù)法一 算法原理 設(shè)語音信號的時(shí)間序列是x(n)，它的自相關(guān)函數(shù)定義為： R(k)=x(n)x(n+k)，K是時(shí)間的延遲量。 短時(shí)自相關(guān)函數(shù)有以下重要性質(zhì) 如果信號x(n)具有周期性，周期是P，那么它的自相關(guān)函數(shù)R(k) 也具有周期性，而且周期與信號x(n)的周期性相同。 當(dāng)k=0時(shí)，短時(shí)自相關(guān)函數(shù)具有最大值，即在延遲量為0，P，2P，時(shí)，周期信號的自相關(guān)函數(shù)也達(dá)到最大值。 短時(shí)自相關(guān)函數(shù)是偶函數(shù)，即R(k) =R(-k

11、)。短時(shí)自相關(guān)函數(shù)基因檢測的主要原理大都是利用它的這些性質(zhì)，通過比較原始信號和它的延遲后的信號之間的類似性質(zhì)來確定基音周期。如果延遲量等于基音周期，那么兩個(gè)信號具有最大類似性，或是直接找出短時(shí)自相關(guān)函數(shù)的兩個(gè)最大值間的距離，即為基音周期的初估值。在用短時(shí)自相關(guān)函數(shù)檢測基音時(shí)，常用歸一化的自相關(guān)函數(shù)，表達(dá)式為r(k)=R(k)/R(0)以上性質(zhì)中已指出，k=0，R(0)為最大值。所以r(k)的模值永遠(yuǎn)小于或等于1。已知基音頻率范圍為60Hz500Hz之間，采樣率為fs時(shí)，則基音周期（樣點(diǎn)值）范圍在fs/500fs/60之間，然后在這個(gè)范圍之間尋找歸一化相關(guān)函數(shù)的最大值，對應(yīng)的延遲量就是基音周期。

12、二代碼編寫用短時(shí)自相關(guān)函數(shù)法得基音周期的程序：functionperiod=ACF_corr(y,fn,vseg,vsl,lmax,lmin)pn=size(y,2);ifpn=fn,y=y;end%把y轉(zhuǎn)換為每列數(shù)據(jù)表示一幀語音信號wlen=size(y,1);%取得幀長period=zeros(1,fn);%初始化fori=1:vsl%只對有話段數(shù)據(jù)處理ixb=vseg(i).begin;ixe=vseg(i).end;ixd=ixe-ixb+1;%求取一段有話段的幀數(shù)fork=1:ixd%對該段有話段數(shù)據(jù)處理u=y(:,k+ixb-1);%取來一幀數(shù)據(jù)ru=xcorr(u,coeff);

13、%計(jì)算歸一化自相關(guān)函數(shù)ru=ru(wlen:end);%取延遲量為正值的部分tmax,tloc=max(ru(lmin:lmax);%在PminPmax范圍內(nèi)尋找最大值period(k+ixb-1)=lmin+tloc-1;%給出對應(yīng)最大值的延遲量endend三實(shí)驗(yàn)結(jié)果四結(jié)果分析在端點(diǎn)檢測和濾波后，再經(jīng)過短時(shí)自相關(guān)函數(shù)法處理得到上面的實(shí)驗(yàn)結(jié)果圖，可看出基音周期有很多的野點(diǎn)，需要進(jìn)行平滑處理，也可以看出基音周期大約為40個(gè)樣點(diǎn)數(shù)，即基音周期為40/fs,采樣率為8000Hz，得周期為0.005s,基音頻率為200Hz.基音周期的檢測（平滑處理） 線性濾波算法 讓基音軌跡序列通過一個(gè)“低通”FIR

14、濾波器。由于基音軌跡中的野點(diǎn)都是一些突發(fā)跳動(dòng)點(diǎn)，它們具有很高的“頻率分量”，采用“低通”濾波即可以將其去除。 線性平滑是用滑動(dòng)窗口進(jìn)行線性濾波處理，即 式中，w(m)，m=-L,-L+1,0,1,2,L為（2L+1）點(diǎn)平滑窗，滿足平滑處理中值濾波算法：一種非線性濾波在被平滑點(diǎn)的左右各取L個(gè)樣點(diǎn)。連同被平滑點(diǎn)共同構(gòu)成一組信號采樣值（共（2L+1）個(gè)樣值），然后將這（2L+1）個(gè)樣值按大小次序排成一對，取此隊(duì)列中間者作為平滑器的輸出。L值一般取為1或2，即中值平滑的“窗口”一般套住3或5個(gè)樣值。中值平滑的優(yōu)點(diǎn)是既可以有效地去除少量的野點(diǎn)，又不會破壞基音周期軌跡中兩個(gè)平滑段之間的階躍性變化。matl

15、ab函數(shù)y=medfilt1(x,n)x為輸入序列；k為窗長，即套住的樣點(diǎn)數(shù)，一般取3或5。y是中值濾波后的輸出序列。設(shè)有一個(gè)一維序列f1,f2,fn,取窗口長度（點(diǎn)數(shù)）為m（m為奇數(shù)）,對其進(jìn)行中值濾波，就是從輸入序列中相繼抽出m個(gè)數(shù)fi-v,fi-1,fi,fi+1,fi+v(其中fi為窗口中心值,v=(m-1)/2),再將這m個(gè)點(diǎn)按其數(shù)值大小順序排序，取其序號的中心點(diǎn)的那個(gè)數(shù)作為濾波輸出。數(shù)學(xué)公式表示為：Yi=Medfi-v,fi-1,fi,fi+1,fi+viNv=(m-1)/2；Yi稱為序列fi-v,fi-1,fi,fi+1,fi+v的中值舉例來說，輸入：Y1-10:1,2,3,4,

16、5,6,7,8,9,10.取區(qū)間2k=4,所以k=2;執(zhí)行中值濾波K=中值濾波(Y)由x-k+1=1，所以當(dāng)k=2時(shí)，x=2濾波時(shí)：K1=Y1K2=(Y1、Y2、Y3、Y4)的中間值，即為2或3 T0=pitfilterm1(period,voiceseg,vosl); functiony=pitfilterm1(x,vseg,vsl) y=zeros(size(x);%初始化 fori=1:vsl%有段數(shù)據(jù) ixb=vseg(i).begin;%該段的開始位置 ixe=vseg(i).end;%該段的結(jié)束位置 u0=x(ixb:ixe);%取來一段數(shù)據(jù) y0=medfilt1(u0,5);%5點(diǎn)的中值濾波 v0=linsmoothm(y0,5);%線性平滑 y(ixb:ixe)=v0;%賦值給y end程序functiony=linsmoothm(x,n)ifnargin2n=3;endwin=han