基于麥克風陣列的室內三維聲源定位優化算法①

黃 靜,胡馨月

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

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

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

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

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

1 陣列信號模型的建立

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

圖1 近場聲音傳播示意圖

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

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

其中,r表示麥克風i的位置坐標;rs表示聲源位置坐標;t表示時間;c為聲速;β表示一個房間6 面墻的反射系數;R表示虛源的位置坐標;δ{·}表示虛源R到麥克風i的直達聲壓,可替代為使用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表示虛源的級數;R=((1-2q)xs+2λxLx-x,(1-2j)ys+2λyLy-y,(1-2k)zs+2λzLz-z).

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

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

本文使用一個五元全向型麥克風陣列進行仿真,以一個中型大小的會議室為場景,房間大小為9 m×8 m×3 m;5 個麥克風的位置為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),陣列結構如圖2所示.聲速c設定為340 m/s;聲源使用Noizeus 中一段3 s的純凈男聲,采樣頻率為8 kHz;聲源坐標設為(1.6,2.3,1.7),在Matlab 平臺根據式(1)模擬麥克風接收到的信號,以加性高斯噪聲作為環境噪聲.對麥克風模擬信號進行帶通濾波和分幀加窗預處理,窗函數選擇典型的Hamming 窗,幀長取256 點,幀移為50%.

圖2 麥克風陣列結構模型

2 定位算法的分析與改進

2.1 TDOA 定位算法

2.1.1 廣義互相關法時延估計

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

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

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

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

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

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

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

2.1.2 球形插值法定位估計

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

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

式(9)可簡化為:

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

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

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

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

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

2.2 SRP-PHAT 定位算法

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

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

那么,估計聲源的位置就在可控響應功率值最大的空間網格:

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

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

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

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

2.3 搜索空間收縮聚類算法

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

SRP-PHAT-SSC 算法需要對室內空間的每個集群質心進行可控功率計算,數量相比SRP-PHAT的全網格搜索雖然已經減少,但算法的實時性仍然不高,可進一步縮減搜索網格數目以減少可控響應功率的計算量,提高算法實時性.

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

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

3 仿真實驗與結果分析

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

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

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

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

圖6 不同信噪比下3 種算法的定位準確率

圖7 不同混響時間下3 種算法的定位準確率

4 結論

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