時延估計(jì)的聲源定位算法及MATLAB實(shí)現(xiàn)

陳向陽，劉 松，劉 培，徐 瑩，沈 超

1．武漢工程大學(xué)計(jì)算機(jī)科學(xué)與工程學(xué)院；2．郵電與信息工程學(xué)院，湖北 武漢 430074；3．湖北省檔案局(館)科技處， 湖北 武漢 430071

0 引 言

自20世紀(jì)70年代起傳聲器陣列開始用于語音信號處理以來，到目前為止已經(jīng)設(shè)計(jì)出各種各樣的聲源定位技術(shù)[1].其中尤其是基于到達(dá)時間延遲的聲源定位技術(shù)因其良好的實(shí)時性、計(jì)算量也不大而成為目前研究的熱點(diǎn)之一.基于到達(dá)時間延遲[2]的聲源定位方法需要先估算出聲源與傳聲器陣列各陣元之間的相對時差，然后根據(jù)得到的時差推導(dǎo)計(jì)算聲源至各陣元之間的距離之差，最后再利用搜索算法或幾何算法得到聲源所處的位置.故其實(shí)現(xiàn)可分為兩步進(jìn)行，即第一步估算出聲音信號至傳聲器陣元之間的相對時間延遲，候選實(shí)現(xiàn)算法有：廣義互相關(guān)函數(shù)法[3]，最小均方(LMS)自適應(yīng)濾波算法[4]；第二步再利用所得的延遲值進(jìn)行數(shù)學(xué)建模并利用陣元空間位置信息進(jìn)行定位，具體實(shí)現(xiàn)算法有角度距離定位算法[5]、球形差值定位算法[6].通過對已知各種算法研究，結(jié)果發(fā)現(xiàn)選擇改進(jìn)的廣義互相關(guān)算法并且與球形差值定位技術(shù)結(jié)合使用能得到更好的定位效果.

1 時延估計(jì)

1.1 廣義互相關(guān)函數(shù)法

假定聲源信號表達(dá)式如式(1)、(2)

x1(n)=a1s(n-τ1)+w1(n)

(1)

x2(n)=a2s(n-τ2)+w2(n)

(2)

這里使用s(n)來表示聲源信號，利用a1、a2來表達(dá)聲波從聲源到傳聲器的傳播衰減系數(shù).而利用τ1、τ2來表達(dá)聲源傳播到兩部傳聲器時所需要的時間，而w1(n)、w2(n)則代表高斯白噪聲[7]，并且s(n)與w1(n)、w2(n)是不相關(guān)[8-9]的.

利用求兩信號之間的互相關(guān)函數(shù)來實(shí)現(xiàn)廣義互相關(guān)函數(shù).兩信號之間的相對時延就是在這種互相關(guān)函數(shù)的峰值位置處.從聲源至兩臺麥克風(fēng)的信號x1(n)和x2(n) 的互相關(guān)函數(shù)R12(τ)可表示為

R12(τ)=E(x1(n)x2(n-τ))

(3)

將式(1)和(2)代入(3)式后便得到

R12(τ)=a1a2E(s(n-τ1)s(n-τ2-τ))+a1E(s(n-τ1)w2(n-τ))+

a2E(s(n-τ2-τ)w1(n))+E(w2(n)w2(n-τ))

(4)

因?yàn)閣1(n)和w2(n)是不相關(guān)的高斯白噪聲，s(n)和w(n)也是不相關(guān)的隨機(jī)信號，故由式(4)得

R12(τ)=E(a1a2s(n-τ1)s(n-τ2-τ))=a1a2Rs(τ-(τ1-τ2))

(5)

利用相關(guān)函數(shù)的性質(zhì)就可以得到，當(dāng)τ=τ12=τ1-τ2時R12(τ)便能夠取得最大值.因此求得R12(τ)的最大值處所對應(yīng)的τ值，其實(shí)就是兩個麥克風(fēng)之間的時延τ12.

利用互功率譜和互相關(guān)函數(shù)之間的關(guān)系就可以得到式(6):

(6)

其中G12(ω)則代表麥克信號x1(n)和x2(n)之間的互功率譜.

但是由于語音與噪聲的短時間處理可能導(dǎo)致出現(xiàn)R12(τ)的峰值不夠突出，由此導(dǎo)致時延τ估計(jì)精確度的下降.為了能夠銳化得到R12(τ)的峰值，利用信號和噪聲的先驗(yàn)知識，可以在頻域內(nèi)對互功率譜G12(ω)進(jìn)行加權(quán)處理來抑制噪聲與混響部分的不利影響.然后再進(jìn)行反變換到時域中，如此就能夠得到廣義互相關(guān)函數(shù)(GCC)，即

(7)

1.2 最小均方(LMS)自適應(yīng)濾波方法

假設(shè)已知在自適應(yīng)濾波器[10]中所包含的輸入信號分別是X(n)和d(n)以及輸出信號分別為y(n)和e(n)，其中X(n)能夠利用參數(shù)可調(diào)整數(shù)字濾波器中wn(n)的加權(quán)分量，以便逐漸逼近所期望的信號d(n)，從而得到輸出響應(yīng)；并且進(jìn)一步將它與參考信號d(n)進(jìn)行對比，最終產(chǎn)生誤差信號e(n)，e(n)能夠利用某種自適應(yīng)算法對濾波器參數(shù)[11]作出改變，以便使e(n)的均方值能夠達(dá)到最小化.

e(n)=d(n-ρ)-y(n)

(8)

(9)

hm(n+1)=hm(n)+μe(n)x(n-m)

(10)

其中m=-p，-p+1，…p-1，p.

2 定位算法

2.1 角度距離定位算法

角度距離定位算法[12]最重要內(nèi)容的一部分就是：利用聲源與兩個麥克風(fēng)之間的距離差，以及兩麥克風(fēng)之間的距離估計(jì)出聲源與這兩麥克風(fēng)所在直線的夾角[13]由式(11)表示

(11)

然后，由于麥克風(fēng)擺放成一定的幾何形狀，根據(jù)幾何形狀中的幾何公式，以聲源到麥克風(fēng)之間的距離差以及麥克風(fēng)之間的距離作為輸入變量，就可以求出聲源到坐標(biāo)原點(diǎn)之間的距離.最后再根據(jù)所得出的距離與角度即可以求出聲源在坐標(biāo)軸中的坐標(biāo).模擬實(shí)驗(yàn)中至少需要有4個麥克風(fēng)，坐標(biāo)原點(diǎn)放置1個，至少有1個坐標(biāo)軸上除了原點(diǎn)的麥克風(fēng)之外還應(yīng)該放置有2個相對原點(diǎn)對稱的麥克風(fēng)，這樣做是為了能夠計(jì)算出聲源到原點(diǎn)之間的距離，而且計(jì)算量很小.除了麥克風(fēng)的擺放問題，模擬試驗(yàn)中聲源到麥克風(fēng)之間的距離差是有正負(fù)之分的，所以在仿真實(shí)驗(yàn)中一定要使之一一對應(yīng)，不然就會得到錯誤的結(jié)果.

2.2 球形插值法

球形插值法[14]中的主要計(jì)算公式如式(12):

(12)

其中Ri=|ri|，Rs=|rs|，利用ri表示mj至mi之間的矢量，利用rs表示mj至聲音源S之間的矢量，利用dij表示聲音源S至傳聲器mi、mj之間的距離差值.在該算法中使用了9個麥克風(fēng)，在坐標(biāo)原點(diǎn)處放置的麥克風(fēng)作為參照麥克風(fēng)，在前邊得出了另外8個麥克風(fēng)與參照麥克風(fēng)之間的聲源距離差之后，就會有8個參數(shù)不同的公式(12)組成的方程組.在解方程組時需要使用最小二乘法，通過它可以直接計(jì)算出聲源的三維坐標(biāo).該算法得到的精確度非常高，在1萬米的數(shù)量級上的相對誤差是百分之零點(diǎn)零幾，在10萬米時產(chǎn)生的誤差也只有1.4%左右，這個距離和誤差程度是可以接受的.此外，在仿真中針對不同的麥克風(fēng)數(shù)量進(jìn)行了大量實(shí)驗(yàn)，結(jié)果發(fā)現(xiàn)至少要使用5個麥克風(fēng)才會有相對精確的測量值.總之仿真結(jié)果說明了該方法對麥克風(fēng)的數(shù)量是具有一定限制的.

3 廣義互相關(guān)函數(shù)時延估計(jì)法改進(jìn)后的相位變換加權(quán)函數(shù)

針對前述難以選擇加權(quán)函數(shù)問題，下面對廣義互相關(guān)函數(shù)(GCC)時延估計(jì)法提出改進(jìn)方法.首先考慮在Gw1w2(ω)不為零時,應(yīng)用中就能夠估計(jì)出無音階段的噪聲互相關(guān)函數(shù)R(w1w2)(τ)，也就是R12(τ)-Rw1w2(τ)，然后再將噪聲的互相關(guān)函數(shù)從信號的互相關(guān)函數(shù)中去除掉，由此削弱因噪聲帶來的不利影響.其次針對在房間中出現(xiàn)的混響現(xiàn)象，可以將G12(ω)分解為兩個成分: 直接成分和反射成分[9].參見式(13)

G12(ω)=a1a2Gss(ω)e-iωτ12+

(13)

由于PHAT加權(quán)函數(shù)的分母|G12(ω)|不等于a1a2Gss(ω)，因此不能確保對e-jωτ12分量的白化作用.因此，考慮修正PHAT加權(quán)函數(shù)，即乘以加權(quán)因子γ，以保證|γ|G12(ω)|=a1a2Gss(ω)，由此一來便能起到較好的抗混響效應(yīng).加權(quán)因子γ實(shí)質(zhì)上就是直接分量占功率譜的比例.可以利用房間的混響時間來估計(jì)出直接分量占據(jù)整個功率譜的比例，由此便得到加權(quán)因子γ.

(14)

即在信噪比較低時，就對頻域給予一定的信噪比加權(quán).在利用兩路信號求解某個時延時，PHAT加權(quán)函數(shù)互功率譜法中需計(jì)算3次快速傅里葉變換.因此其運(yùn)算量相對還是比較低的.但具有抑制噪聲作用，因此在實(shí)際定位系統(tǒng)中，得到改進(jìn)后的PHAT加權(quán)函數(shù)法得到了很好的應(yīng)用.

由于在仿真模擬實(shí)驗(yàn)中針對聲源信號的噪音是人為加上的，所以要得出只有噪聲的信號是非常容易的，那么在得到了噪音的自相關(guān)函數(shù)[13]之后，就可以在信號的互相關(guān)函數(shù)中減掉該部分，經(jīng)過仿真之后發(fā)現(xiàn)結(jié)果的確有所改進(jìn)，但是效果還不夠明顯.考慮到信號與噪音可能并不是完全的不相關(guān)，也就是說信號和噪音的互相關(guān)性可能也會產(chǎn)生一定的影響，于是可用信號互相關(guān)函數(shù)減去一個分量，即

(15)

其中rnoise是噪音的相關(guān)函數(shù)，rxn與ryn是噪音和兩個信號的互相關(guān)函數(shù)， 是一個比例系數(shù)體現(xiàn)的是這三個分量之間的比例.經(jīng)過多次仿真試驗(yàn)發(fā)現(xiàn)當(dāng)k為0.3時消除噪音的效果最佳.

4 模擬仿真結(jié)果

改進(jìn)前廣義互相關(guān)函數(shù)算法的仿真結(jié)果如圖1所示.改進(jìn)前10db的結(jié)果如圖2所示.

圖1信噪比為20db時的互相關(guān)函數(shù)圖

Fig.1 The cross-correlation function diagram with signal-to-noise ratio=20dB

圖2 信噪比為10db時的互相關(guān)函數(shù)圖

在10db時算法經(jīng)改進(jìn)后的結(jié)果如圖3所示.由圖2與圖3的比較可以得知改進(jìn)算法后的方法對噪音的抑制效果更好.

圖3 改進(jìn)后信噪比為10db互相關(guān)函數(shù)圖

在20db時使用最小均方(LMS)自適應(yīng)算法得到的結(jié)果如圖3～4所示.從圖3～4中可以看出在20db時，時延5處出現(xiàn)明顯的峰值.上述仿真結(jié)果驗(yàn)證了該方法的正確性.

圖4 信噪比為20的權(quán)值時域圖

在10db時該方法的結(jié)果如圖5所示.由圖5可以看出，自適應(yīng)算法在信噪比為10db時已經(jīng)出現(xiàn)較多的峰值，而且在正確的時延處峰值處略有偏移.所以通過以上實(shí)驗(yàn)仿真的結(jié)果可以看出利用改進(jìn)后的廣義互相關(guān)函數(shù)算法要比自適應(yīng)算法具有更好的抑制噪音效果.

圖5 信噪比為10的權(quán)值時域圖

5 總結(jié)與展望

通過對不同聲源定位中的到達(dá)時間延遲算法以及定位算法的對比，以及在Matlab仿真軟件上的模擬，說明改進(jìn)的互相關(guān)函數(shù)法相對于自適應(yīng)算法具有更好的抗噪能力，球形插值法也具有更精確的定位能力，所以建議實(shí)際應(yīng)用中使用這兩種算法能夠?qū)崿F(xiàn)更精確的聲源定位.同時由于諸多限制，仿真模擬試驗(yàn)中還有諸多問題有待解決，如在互相關(guān)函數(shù)中減去噪音成分的比例具體應(yīng)該是多少，傅里葉變換對互相關(guān)函數(shù)的影響到應(yīng)該是什么樣的，傅里葉變化的點(diǎn)數(shù)與信號的采樣點(diǎn)數(shù)的取值是否有一定的要求，才能使得最后的峰值更加明顯等問題.

致 謝

本文得到了湖北省教育科學(xué)"十二五"規(guī)劃2013年度立項(xiàng)課題：地方性高校新興交叉學(xué)科建設(shè)的研究——以網(wǎng)絡(luò)科學(xué)為例(項(xiàng)目編號2013B060)的資助，在此表示衷心的感謝.

[1] 李承智，曲天書，吳璽宏.一種改進(jìn)的AEDA聲源定位及跟蹤算法[J].北京大學(xué)學(xué)報：自然科學(xué)版，2005,41(5)：809-814.

LI Chengzhi,Qu Tianshu,WU Xihong.A modified AEDA algorithm for sound source localization and tracking[J]. Acta Scientiarum Naturalium Universitatis Pekinensis, 2005,41(5):809-814.(in chineses)

[2] 郭蘭申. 面向反恐領(lǐng)域基于多源信息融合的機(jī)器人感知系統(tǒng)研究[D].天津：河北工業(yè)大學(xué)，2006.

GUO Lanshen. Research on robot perception system for anti-terrorism based on multi-source information fusion[D]. Tianjin：Hebei University of Technology，2006. (in chineses)

[3] 金乃高，陳喆. 基于分層采樣粒子濾波的麥克風(fēng)陣列說話人跟蹤法[J]. 電子學(xué)報，2008,36(1)： 194-198.

JIN Nai-gao，CHEN Zhe. Microphone array speaker tracking using layered sampling particle filtering[J]. ACTA Electronica Sinica,2008,36(1):194-198. (in chineses)

[4] 殷福亮，陳喆. 基于加權(quán)子空間擬合的聲源定位與跟蹤方法[J]. 電子與信息學(xué)報，2008, 30(9):2134-2137.

YIN Fu-liang, CHEN Zhe. Weighted subspace fitting sound source localization method based on particle filtering[J]. Journal of Electronics & Information Technology, 2008, 30(9):2134-2137.(in chineses)

[5] 蔡衛(wèi)平. 基于麥克風(fēng)陣列的聲源定位與跟蹤算法研究[D].南京：東南大學(xué)，2010.

CAI Weiping. Research and tracking algorithm for sound source localization based on microphone array[D].Nanjing: Southeast University,2010. (in chineses)

[6] 馮道寧，王浩. 基于麥克風(fēng)陣列和虛擬儀器的聲源定位[J]. 計(jì)算機(jī)工程與應(yīng)用，2012, 48(21) : 130-132.

FENG Dao-ning, WANG Hao. Speech source localization based on microphone arrays and labview[J]. CEA, 2012, 48(21): 130-132. (in chineses)

[7] 邵懷宗，林靜然，彭啟琮,等. 基于麥克風(fēng)陣列的聲源定位研究[J]. 云南民族大學(xué)學(xué)報，2004，13(4):256-259.

SHAO Huaizong,LIN jingrui,PENG Qicong ,et al. Study on localization for speech based on microphone array[J]. Journal of Yunnan Nationalities University,2004, 13(4):256-259. (in chineses)

[8] 陸曉燕，基于麥克風(fēng)陣列實(shí)現(xiàn)聲源定位[D]. 大連:大連理工大學(xué)，2003.

LU Xiaoyan, Implementation of acoustic source localization based on microphone array[D].Dalian: Dalian University of Technology,2003. (in chineses)

[9] 周浩祥.基于麥克風(fēng)陣列的聲源定位方法研究[D] . 大連:大連理工大學(xué)，2002.

ZHOU Haoxiang. Research on acoustic source localization methods based on microphone array[D].Dalian:Dalian University of Technology,2003. (in chineses)

[10] 嚴(yán)素清，黃冰. 傳聲器陣列的聲源定位研究[J]. 電聲技術(shù)，2004(12)：27-30.

YAN Suqing,HUANG Bing. Study on acoustic source localization based on microphone arrays[J].Loudspeaker And Microphone, 2004(12)：27-30. (in chineses)

[11] YOUN D H, AHMED N, CARTER G C. On using the LMS algorithm for time delay estimation[J]. IEEE Transactions on Acoustics, Speech and Signal Processing, 1992,30(5):798-801.

[12] JIE Huang, Kume K, Saji A. et al. Robotic spatial sound socalization and its 3D sound human Interface[C]//First International Symposium on Cyber Worlds (CW2002),Tokyo:2002:191-197.

[13] 王宏禹，邱天爽. 自適應(yīng)噪聲抵消和時間延遲估計(jì)[D]. 大連:大連理工大學(xué)出版社，1999.

WANG Hongyu,QIU Tianshuang. Adaptive noise cancellation and time delay estimation[D].Dalian: Dalian University of Technology press,1999. (in chineses)

[14] JIE Huang, OHNISHI N, SUGIE N. A biomimetic system for localization and separation of multiple sound sources instrumentation and measurement[J]. IEEE Transactions on Volume,1995,44(3):733-738.

[15] SMITH J O, ABEL J S. Closed-form least-squares source location estimation from range difference measurements[J]. IEEE Transactions on Acoustics, Speech and Signal Processing,1987,35(12):1661-1669.