黃 靜,胡馨月

(浙江理工大學(xué) 信息學(xué)院,杭州 310018)

信號源的定位是傳感器陣列信號處理領(lǐng)域重要部分,曾應(yīng)用于雷達(dá)、聲吶等大型通信設(shè)備中,到20世紀(jì)80年代,隨著語音通信發(fā)展,基于麥克風(fēng)陣列的語音信號處理也逐步受到學(xué)者們關(guān)注,在會議發(fā)言人定位、機器人聽覺系統(tǒng)、安防等場景都得到了廣泛的研究與應(yīng)用[1].目前基于小型麥克風(fēng)陣列的聲源定位大多是二維層面,要實現(xiàn)三維定位則需要增加麥克風(fēng)組成大規(guī)模麥克風(fēng)陣列采集聲源信息,增加了設(shè)備成本和計算成本.因此,如何利用小型麥克風(fēng)陣列進(jìn)行聲源的三維實時定位成了重要的研究課題.

在目前的研究中,近場聲源定位算法主要分為3 類:(1)基于高分辨譜估計的定位算法.該算法主要是針對窄帶信號,若要運用到語音信號中,需要將寬帶的語音信號先分解為窄帶信號,因轉(zhuǎn)換損失會造成偏差[2,3].(2)基于可控響應(yīng)功率的定位算法.該算法對寬帶信號和窄帶信號都適用,但需要對整個空間的網(wǎng)格點進(jìn)行搜索和功率計算,計算量大導(dǎo)致響應(yīng)速度慢,不適用于實時的聲源定位系統(tǒng)[4].(3)基于到達(dá)時間差的定位算法.該算法分為兩個部分,第1 步是時延估計,第2 步是定位估計,計算量小,實時性高,在時延估計階段有較高要求,否則會在定位估計中會放大時延估計的誤差[5].

針對SRP-PHAT 算法響應(yīng)速度慢的缺陷,很多學(xué)者提出了改進(jìn)方法.Cho 等[6]提出了針對小型麥克風(fēng)陣列的SSC 加速算法,利用相同到達(dá)時間差聚類取質(zhì)心的方法減少搜索網(wǎng)格,但僅對俯仰角和方向角進(jìn)行了計算,缺少測距,且聚類后的質(zhì)心數(shù)目仍然;Yook等[7]在SSC 優(yōu)化算法的基礎(chǔ)上,提出了針對大型麥克風(fēng)陣列的兩級搜索空間聚類(Two-Level Search Space Clustering,TL-SSC)算法,將聲源的候選位置劃分成組,找到少數(shù)可能包含最大功率的組,進(jìn)一步減少SSC算法的搜索網(wǎng)格,但對于小型麥克風(fēng)陣列的計算量減小效果不明顯;Salvati 等[8]利用搜索網(wǎng)格上TDOA 信息的非均勻積累提出了基于靈敏度的區(qū)域選擇SRP(Region selection SRP,R-SRP)算法以減少聲源定位誤差,但計算量較多,實時性差.

本文介紹了陣列信號模型和傳統(tǒng)定位算法,改進(jìn)了TDOA 定位算法以提高定位準(zhǔn)確性,改進(jìn)了SRPPHAT-SSC 定位算法以簡化聚類邏輯,同時針對小型麥克風(fēng)陣列三維定位算法實時性不足的問題,提出了TDOA 定位算法與SRP-PHAT-SSC 算法相結(jié)合的搜索空間收縮聚類算法,并在Matlab 平臺對3 種SRPPHAT 搜索方式進(jìn)行了仿真比較.算法擬針對五元小型麥克風(fēng)陣列,在SRP-PHAT-SSC 加速算法的基礎(chǔ)上進(jìn)一步減少搜索網(wǎng)格數(shù),提高SRP-PHAT 算法的實時性和三維定位精度.

1 陣列信號模型的建立

陣列模型的建立是聲學(xué)研究的基礎(chǔ).本文針對室內(nèi)的聲源定位,聲場是近場模型,聲波以球面波的形式擴(kuò)散,圖1為近場中一個聲源到兩個麥克風(fēng)的傳播示意圖[9],其中,Source 表示聲源;Mic1和Mic2表示兩個麥克風(fēng);τ為聲源信號到達(dá)兩個麥克風(fēng)的時延,單位為s;c為聲速,單位為m/s;聲源到兩個麥克風(fēng)的距離差d為τ×c.

圖1 近場聲音傳播示意圖

根據(jù)圖1,時延為τ的兩個麥克風(fēng)接收到的信號模型xi(t)(i=1,2)可表示為[10]:

其中,s(t)為聲源信號;hi為聲源到麥克風(fēng)i之間的脈沖響應(yīng)函數(shù),由式(2)表示;τi為聲源到達(dá)麥克風(fēng)i的時間,τ=τ1-τ2;ni(t)為第i個麥克風(fēng)的環(huán)境噪聲信號,s(t)、n1(t)、n2(t)之間互不相關(guān).

其中,r表示麥克風(fēng)i的位置坐標(biāo);rs表示聲源位置坐標(biāo);t表示時間;c為聲速;β表示一個房間6 面墻的反射系數(shù);R表示虛源的位置坐標(biāo);δ{·}表示虛源R到麥克風(fēng)i的直達(dá)聲壓,可替代為使用Hanning 窗的理想低通濾波器,窗長為4 ms,截止頻率為奈奎斯特頻率;P={(q,j,k) |q,j,k?{0,1}};Λ={(λx,λy,λz) |λx,λy,λz∈Z:-N≤λx,λy,λz≤N},N表示虛源的級數(shù);R=((1-2q)xs+2λxLx-x,(1-2j)ys+2λyLy-y,(1-2k)zs+2λzLz-z).

在實際使用中,一般通過設(shè)定混響時間RT60來改變聲源和麥克風(fēng)i之間的脈沖響應(yīng),根據(jù)Sabin 經(jīng)驗公式可以由混響時間求出式(2)中的墻壁反射系數(shù)β.混響時間和墻壁反射系數(shù)的關(guān)系為:

其中,V表示房間的體積,Si表示第i面墻的面積.

本文使用一個五元全向型麥克風(fēng)陣列進(jìn)行仿真,以一個中型大小的會議室為場景,房間大小為9 m×8 m×3 m;5 個麥克風(fēng)的位置為M1(2,2,1.6)、M2(2,1.9,1.5)、M3(2,2.1,1.5)、M4(1.9,2,1.5)、M5(2.1,2,1.5),陣列結(jié)構(gòu)如圖2所示.聲速c設(shè)定為340 m/s;聲源使用Noizeus 中一段3 s的純凈男聲,采樣頻率為8 kHz;聲源坐標(biāo)設(shè)為(1.6,2.3,1.7),在Matlab 平臺根據(jù)式(1)模擬麥克風(fēng)接收到的信號,以加性高斯噪聲作為環(huán)境噪聲.對麥克風(fēng)模擬信號進(jìn)行帶通濾波和分幀加窗預(yù)處理,窗函數(shù)選擇典型的Hamming 窗,幀長取256 點,幀移為50%.

圖2 麥克風(fēng)陣列結(jié)構(gòu)模型

2 定位算法的分析與改進(jìn)

2.1 TDOA 定位算法

2.1.1 廣義互相關(guān)法時延估計

TDOA 定位算法的第一步是時延估計,廣義互相關(guān)法(Generalized Cross-Correlation,GCC)是一種最常用的時延估計算法[11,12].通過最大化一對麥克風(fēng)信號在頻域上的互相關(guān)函數(shù)可估計出聲音信號到達(dá)這對麥克風(fēng)的時間差[5].

信號x1與x2之間的時延τx1x2定義為:

信號x1與x2之間的廣義互相關(guān)函數(shù)Rx1x2(τ)為:

其中,ω為頻率;j=-2πi;Gx1x2(ω)為x1與x2的互功率譜;X1(ω)與X2(ω)分別是x1與x2經(jīng)過快速傅里葉變換(FFT)后的頻域信號;(*)表示復(fù)共軛運算.

由于混響和噪聲容易對廣義互相關(guān)函數(shù)產(chǎn)生影響,為了減少混響和噪聲的影響,本文引入PHAT(PHAse Transform,相位變換)加權(quán)函數(shù)對互功率譜函數(shù)作白化處理,PHAT 加權(quán)對混響有一定的抑制效果[13,14].采用一對麥克風(fēng)進(jìn)行時延估計仿真,GCC-PHAT 函數(shù)與GCC 函數(shù)的對比如圖3所示,GCC-PHAT 函數(shù)只保留了相位信息,從而銳化互相關(guān)函數(shù)的峰值.PHAT 加權(quán)函數(shù)的公式為:

圖3 信噪比0 dB、混響時間0 s 時,GCC 與GCC-PHAT的比較

則式(5)引入PHAT 加權(quán)可寫為:

2.1.2 球形插值法定位估計

定位估計是TDOA 定位算法的第2 步,為了保證TDOA 定位算法的實時性,本文選用球形插值法(Spherical Interpolation,SI)用于定位估計.球形插值法屬于最小二乘(Least Square,LS)算法,復(fù)雜度低,基本思想是以麥克風(fēng)陣列中的一個麥克風(fēng)為參考麥克風(fēng),通過最小化聲源到參考麥克風(fēng)與其他麥克風(fēng)之間距離差的估計值與實際值的誤差平方和,求出聲源位置[15,16].

假設(shè)有N個麥克風(fēng)的坐標(biāo)(M1,M2,…,MN),若以M1為參考麥克風(fēng),并作為空間坐標(biāo)原點,此時,Mi的坐標(biāo)為qi,聲源s到Mi與M1的估計距離差為,是根據(jù)第一步的估計時延τ得到的值,Mi到坐標(biāo)原點的距離為Ri,聲源位置坐標(biāo)為qs,則Mi到原點距離的平方可表示為:

式(9)可簡化為:

則式(10)的最小二乘解可表示為:

若估計聲源的位置與實際位置中的方向角φ相差Δφ度、俯仰角θ相差Δθ度、徑向距離r相差Δr米,則記為異常值,否則記為正確值.本文使用準(zhǔn)確率作為評判準(zhǔn)則:準(zhǔn)確率=正確值的幀數(shù)和/總幀數(shù).

由于GCC-PHAT 算法的時延估計誤差會在球形插值法中放大,定位錯誤率較高,所以本文引入離群值校正,對基于TDOA 算法的定位結(jié)果進(jìn)行校正處理,將離群值校正為中位數(shù),其中,離群值為中位數(shù)絕對偏差(Median Absolute Deviation,MAD)與MAD 中位數(shù)相差超過3 倍的值.

設(shè)置方向角、俯仰角和徑向距離的容錯范圍分別是Δφ=15、Δθ=15、Δr=1.在5 元錐形麥克風(fēng)陣列的情況下,TDOA 定位算法無法準(zhǔn)確計算出俯仰角,經(jīng)過仿真實驗得知,俯仰角的準(zhǔn)確率均在20%以下.當(dāng)無混響影響時,如圖4所示,方向角和徑向距離準(zhǔn)確率隨信噪比的減少而下降,定位結(jié)果容易受到噪聲影響,做了離群值校正后,TDOA 定位算法的魯棒性增強,準(zhǔn)確率都在85%以上;當(dāng)只有混響影響時,由于PHAT 加權(quán)函數(shù)對混響的抑制作用,所以定位結(jié)果都在容錯范圍內(nèi),TDOA在方向角和徑向距離兩個方面的準(zhǔn)確率均為100%,不受混響時間的影響.

圖4 信噪比對TDOA 定位算法的影響

2.2 SRP-PHAT 定位算法

SRP-PHAT 算法屬于一步定位的算法,其PHAT濾波被用于提高SRP 算法魯棒性[17].通過對定位空間劃分格子,每一個格子進(jìn)行所有麥克風(fēng)對之間互相關(guān)函數(shù)的求和來精確功率的計算,從而確定最大功率的空間點作為估計聲源的位置[18].

定義q為空間點位置坐標(biāo),Rlm[τlm(q)]為麥克風(fēng)l和麥克風(fēng)m所接收信號的GCC-PHAT 函數(shù),則可控響應(yīng)功率P(q)為:

那么,估計聲源的位置就在可控響應(yīng)功率值最大的空間網(wǎng)格:

可控響應(yīng)功率取決于陣列中每個麥克風(fēng)對的接收信號在頻域上的相位差,而頻譜的相位差在時域上是離散的,所以如果有兩個空間點的位置很接近,則這兩個空間點到各個麥克風(fēng)對的時間差相同,那么就無法依靠式(14)求出聲源位置.

SSC 算法是一種針對基于SRP-PHAT的小型麥克風(fēng)陣列聲源定位而提出的加速算法[6].SSC 算法將具有相同TDOA的空間點集合在一起形成集群,然后求出每個集群的質(zhì)心,將每個集群的質(zhì)心作為代表坐標(biāo)存儲在一個查找表中,從而避免了上述問題的出現(xiàn),同時進(jìn)行聲源定位時大大減少了功率計算的時間,因為對于小型麥克風(fēng)陣列,最終組成的集群的數(shù)目遠(yuǎn)小于初始的空間網(wǎng)格總數(shù).

文獻(xiàn)[6]中使用的是自頂向下的聚類方法構(gòu)造查找表,本文使用邏輯與實現(xiàn)更為簡單的自底向上聚類方法,其算法流程如算法1.

算法1.自底向上聚類算法1)以球面坐標(biāo)(φ,θ,r)劃分空間網(wǎng)格,形成網(wǎng)格點坐標(biāo)集合;2) 遍歷,計算每個網(wǎng)格點中每組麥克風(fēng)對的TDOA,存儲到TDOA 查找表;3)將相同TDOA 值的網(wǎng)格點劃分成簇;4)計算每個簇的質(zhì)心,形成質(zhì)心集合{c};5)遍歷{c},計算每個集群質(zhì)心的TDOA,存儲到TDOA_SSC 查找表;6)計算每個質(zhì)心的可控響應(yīng)功率;7)功率最大的質(zhì)心位置為估計聲源位置.

2.3 搜索空間收縮聚類算法

本文改進(jìn)后的TDOA 算法提高了定位計算的魯棒性,在五元錐形麥克風(fēng)陣列且信噪比大于10 dB的情況下,可在一定容錯范圍內(nèi)得到正確的方向角和徑向距離,但無法計算得到俯仰角.

SRP-PHAT-SSC 算法需要對室內(nèi)空間的每個集群質(zhì)心進(jìn)行可控功率計算,數(shù)量相比SRP-PHAT的全網(wǎng)格搜索雖然已經(jīng)減少,但算法的實時性仍然不高,可進(jìn)一步縮減搜索網(wǎng)格數(shù)目以減少可控響應(yīng)功率的計算量,提高算法實時性.

因此,針對小型麥克風(fēng)陣列的三維定位,本文提出一種基于TDOA和SRP-PHAT-SSC的組合算法-搜索空間收縮聚類算法(Search Space Shrinking Clustering,SSSC),通過引入改進(jìn)的TDOA 算法以進(jìn)一步減少SRP-PHAT-SSC 算法中的候選聲源位置,增加算法的魯棒性并提高實時性.該算法主要分為兩個階段,粗定位階段和細(xì)定位階段,整體流程如圖5所示.首先使用自底向上聚類的SSC 算法得到整個空間的集群質(zhì)心位置集合B并保存,用于提供網(wǎng)格查詢的功能,采取5 cm的細(xì)粒度網(wǎng)格以提高SRP-PHAT-SSC算法的定位精度.首先是粗定位階段,對音頻信號幀進(jìn)行基于TDOA的定位計算,先利用式(8)求得估計時延,再利用式(12)求得估計聲源坐標(biāo),得到估計聲源位置集合Qs,對Qs進(jìn)行離群值校正以提高TDOA 定位的魯棒性,由于基于五元錐形麥克風(fēng)陣列的TDOA 定位算法無法準(zhǔn)確判斷俯仰角,所以取其方向角和徑向距離作為空間收縮的參考.TDOA 定位的容錯范圍設(shè)為Δφ=15和Δr=1,根據(jù)設(shè)置的TDOA 定位容錯范圍進(jìn)行搜索空間收縮以減少可控功率的計算量,收縮后的空間Bs為之后進(jìn)入細(xì)定位階段,根據(jù)Bs對集群質(zhì)心坐標(biāo)集合B進(jìn)行篩選,僅對收縮空間Bs內(nèi)的集群質(zhì)心使用式(13)進(jìn)行可控響應(yīng)功率計算,之后根據(jù)式(14)取功率最大的質(zhì)心位置作為聲源.

圖5 搜索空間收縮聚類算法流程

3 仿真實驗與結(jié)果分析

本文仿真實驗中,SRP-PHAT 算法的全網(wǎng)格搜索方法和搜索空間聚類方法的網(wǎng)格尺寸分別取徑向距離的精度為20 cm、10 cm和5 cm,搜索空間收縮聚類算法取5 cm的細(xì)粒度網(wǎng)格尺寸,俯仰角和方向角的精度為1°,設(shè)置Δφ=10、Δθ=10、Δr=1,若最終估計聲源位置的球坐標(biāo)在該設(shè)定范圍內(nèi),則標(biāo)記為正確定位.

表1以信噪比為30 dB,混響時間為0 s的情況,對3 種算法進(jìn)行了對比,根據(jù)表1可以看出:本文提出的SSSC 算法搜索的格子數(shù)相比SSC 算法減少了96%,從而使得計算時間比SSC 算法減少了97%,具有較好的實時性,且定位準(zhǔn)確率達(dá)95%.

表1 不同SRP-PHAT 算法的仿真結(jié)果對比

圖6和圖7是對3 種算法在球坐標(biāo)系的3 個維度上分別進(jìn)行準(zhǔn)確率判斷,其中,FGS和SSC 算法選取與搜索空間收縮聚類算法相同網(wǎng)格粒度的數(shù)據(jù)作為展示.可以看出,由于PHAT 加權(quán)函數(shù)對混響的抑制作用,3 個算法的準(zhǔn)確率在僅有混響的情況下幾乎不變.SSC 算法無法較好的估計出聲源的徑向距離,改進(jìn)后的算法通過引入TDOA 算法計算出準(zhǔn)確的徑向距離,從而彌補了SSC 算法的缺點.在低信噪比的情況下,FGS和SSC 算法的準(zhǔn)確率降低,而改進(jìn)后的算法在3 個維度上的準(zhǔn)確率都保持在80%以上,這是因為其在TDOA定位階段進(jìn)行了離群值校正,從而提高了魯棒性.

圖6 不同信噪比下3 種算法的定位準(zhǔn)確率

圖7 不同混響時間下3 種算法的定位準(zhǔn)確率

4 結(jié)論

本文分析了基于TDOA的聲源定位算法和SRPPHAT 聲源定位算法,在二者的基礎(chǔ)上提出了基于SRPPHAT的搜索空間收縮聚類優(yōu)化算法,本文選取了不同的網(wǎng)格步長、信噪比和混響時間,對全空間搜索、搜索空間聚類和搜索空間聚類算法進(jìn)行了對比分析,從實驗結(jié)果可以看出搜索空間聚類算法減少了需要進(jìn)行功率計算的網(wǎng)格數(shù),從而減少了運行時間,具有較好的實時性,并且在球坐標(biāo)的3 個維度都有著高準(zhǔn)確率,在實際應(yīng)用中有一定的價值.