侯東升,崔遜學(xué)
(陸軍炮兵防空兵學(xué)院 a.研究生管理大隊(duì); b.六系,合肥 230031)
近年來,隨著傳感器與通信領(lǐng)域的發(fā)展,基于到達(dá)時(shí)差(Time Difference of Arrival,TDOA)的聲陣列定位算法引起了廣泛關(guān)注。當(dāng)輻射源距離陣列很遠(yuǎn)時(shí),人們難以直接通過TDOA雙曲線交叉方法精確估計(jì)輻射源坐標(biāo)。但如果在不同位置分散部署多個(gè)陣列,便可以利用各陣列所接收到的信號(hào)平面波之間的時(shí)差來估計(jì)波達(dá)方向(Direction of Arrival,DOA)[1-2]。該過程計(jì)算量小、易于實(shí)現(xiàn),且只要求傳感器之間滿足嚴(yán)格時(shí)間同步,因此,便于在實(shí)際中操作與應(yīng)用[3]。
文獻(xiàn)[4]較早提出基于時(shí)差的輻射源方向問題,其設(shè)計(jì)一種線性最小二乘(Linear Least Squares,LLS)測(cè)向算法,由3個(gè)分量構(gòu)成的DOA矢量估計(jì)可計(jì)算出輻射源的方位角和俯仰角。LLS算法具有簡(jiǎn)捷、易操作的特點(diǎn),且其提供的解方案是閉式的。
基于時(shí)差的測(cè)向問題最近得到較多關(guān)注,特別是在槍聲定位等公眾安全領(lǐng)域。如果測(cè)向陣列是固定結(jié)構(gòu),則可通過聯(lián)合建立若干個(gè)TDOA測(cè)量方程并利用三角函數(shù)關(guān)系來求解三維空間中輻射源的方位角和俯仰角[5]。當(dāng)陣列結(jié)構(gòu)屬于不規(guī)則拓?fù)湫螤顣r(shí),可采用文獻(xiàn)[4]提出的LLS算法來獲得DOA估計(jì)的閉式解。但 LLS算法的不足之處在于當(dāng)TDOA測(cè)量噪聲較大時(shí),會(huì)干擾DOA的估計(jì)結(jié)果,使算法性能下降。
如果知道輻射源現(xiàn)場(chǎng)傳播信號(hào)的速度參量,例如對(duì)于聲信號(hào)而言,通過測(cè)量現(xiàn)場(chǎng)的氣象參數(shù),則可計(jì)算出實(shí)際聲速的估計(jì)值[6],再結(jié)合LLS算法的DOA估計(jì)結(jié)果,即通過采用一些優(yōu)化方法可以直接獲得輻射源方位角和俯仰角這2個(gè)變量的估計(jì)值。文獻(xiàn)[7]提出采用智能優(yōu)化方法來搜索得到這2個(gè)參變量的估計(jì)值,但需人工設(shè)置較多參數(shù),不利于工程實(shí)現(xiàn)。
有研究結(jié)果表明,由多個(gè)到達(dá)時(shí)差測(cè)量值構(gòu)成的測(cè)量方程具有非線性的特點(diǎn),通常需要將其轉(zhuǎn)換為線性方程組后再進(jìn)行求解。文獻(xiàn)[8]提出的泰勒級(jí)數(shù)展開法(以下簡(jiǎn)稱Taylor算法)對(duì)處理非線性方程組十分有效,且操作步驟較易設(shè)計(jì),文獻(xiàn)[9-10]即提供了Taylor算法的具體過程和步驟。但由文獻(xiàn)[11]可知,Taylor算法為局部收斂迭代算法,對(duì)初始值要求較高,迭代初始值必須接近于真實(shí)值,否則會(huì)發(fā)生局部最小化現(xiàn)象或者不能保證結(jié)果收斂[12]。
針對(duì)上述方法存在的問題,本文以聲音在空氣中傳播的測(cè)向問題作為應(yīng)用背景,提出一種基于萊溫伯格-馬夸特(Levenberg-Marquardt,LM)[13]的測(cè)向算法。用LLS算法求出的閉式解作為L(zhǎng)M算法的初始方位估計(jì),通過LM算法直接求得輻射源的方位,以進(jìn)行對(duì)聲源方位的高精度估計(jì)。在實(shí)驗(yàn)環(huán)節(jié)將LLS算法、Taylor算法和本文提出的LM算法進(jìn)行性能比較,以驗(yàn)證本文LM算法的定位性能。
三維空間測(cè)向問題示意圖如圖1所示。
圖1 三維空間輻射源測(cè)向示意圖
(1)
其中,3個(gè)分量kx、ky、kz分別為輻射源DOA在3個(gè)坐標(biāo)軸上的單位投影。在本文中,符號(hào)上方有“^”的值為估計(jì)值或測(cè)量值,否則為真實(shí)值。
通常假定第0個(gè)傳感器為TDOA的測(cè)量基準(zhǔn),其位于坐標(biāo)系原點(diǎn),則存在如式(2)所示的TDOA測(cè)量方程。
(2)
則式(2)變?yōu)槿缦戮仃囆问?
cτ=SK
(3)
由于真實(shí)TDOA值無法獲得,實(shí)際測(cè)量結(jié)果均存在測(cè)量誤差,用Δτ表示TDOA的測(cè)量誤差矢量,則:
(4)
傳統(tǒng)的LLS算法通過對(duì)式(3)進(jìn)行求解,得到如下的估計(jì)結(jié)果[4]:
(5)
(6)
根據(jù)DOA在三維空間的三角函數(shù)關(guān)系,方位角和俯仰角的估計(jì)值可分別表示為:
(7)
本文LM算法對(duì)初值的要求較低。將LLS算法的估計(jì)結(jié)果作為L(zhǎng)M算法的初值,完全可以滿足要求。在此基礎(chǔ)上進(jìn)行迭代運(yùn)算,可以提高方位估計(jì)的精度。
將式(2)變換為:
(8)
則目標(biāo)函數(shù)定義為:
(9)
(10)
將式(10)右邊取最小化,則:
(11)
為避免計(jì)算海森函數(shù)的二階信息,將式(11)化簡(jiǎn),結(jié)果為:
(12)
(13)
λk=λk-1/1+α
阻尼因數(shù)更新后,根據(jù)式(12)更新增長(zhǎng)矢量,然后根據(jù)式(13)更新方位估計(jì),循環(huán)運(yùn)行該過程,直至方位估計(jì)增加矢量足夠小。
總結(jié)上述過程,可將LM算法的步驟表示為:
算法1LM測(cè)向算法
輸入傳感器數(shù)目及其坐標(biāo)、TDOA測(cè)量值、迭代精度EPS
由式(8)可得測(cè)量方程為:
(14)
(15)
通過一些數(shù)學(xué)處理可以得出費(fèi)歇爾信息矩陣:
(16)
其中:
則克拉姆-拉奧下界(Cramer-Rao Lower Bound,CRLB)為:
CRLB=FIM-1(γ)=((F)TW-1(F))-1
(17)
命題LM算法估計(jì)結(jié)果的均方差為:
MSE(Δγ)=((h(γ)))-1
(18)
(19)
由式(17)可得:
(20)
cov(Δγ)=(h(γ))-1
(21)
將式(21)與克拉姆-拉奧下界進(jìn)行比對(duì):
cov(Δγ)= (h(γ))-1=
(c(F)T(c2W)-1(c(F)))-1=
4.1.1 參數(shù)設(shè)置
本次實(shí)驗(yàn)的性能評(píng)估指標(biāo)是測(cè)向誤差的累積分布函數(shù)(Cumulative Distribution Function,CDF)。CDF描述實(shí)數(shù)隨機(jī)變量的概率分布情況,其值越大,表示算法性能越好。輻射源測(cè)向誤差為輻射源的估計(jì)方向與真實(shí)方向之間的差值。模擬實(shí)驗(yàn)以聲源測(cè)向?yàn)楸尘?具體參數(shù)設(shè)置如下:
1)采用Matlab工具,每輪蒙特卡洛模擬均運(yùn)行10 000次實(shí)驗(yàn),LM算法和Taylor算法的迭代精度閾值為10-5。
2)傳感器的位置隨機(jī)部署,其x軸和y軸坐標(biāo)在15 m半徑的圓形內(nèi),z軸坐標(biāo)在 [-5 m,5 m]的范圍內(nèi)。這與通常的聲探測(cè)野外部署條件相似,且DOA估計(jì)結(jié)果位于合理范圍內(nèi)。
3)輻射源方位角與俯仰角的取值范圍均為[30°,60°]。
4)考慮傳感器數(shù)目變化,實(shí)驗(yàn)分2種情形:傳感器數(shù)目較少,在[4,8]范圍內(nèi)隨機(jī)選取;傳感器數(shù)目較多,在[8,12]范圍內(nèi)隨機(jī)選取。
5)TDOA測(cè)量誤差服從高斯分布,不考慮測(cè)量的系統(tǒng)偏差問題。標(biāo)準(zhǔn)偏差在[0,20]范圍內(nèi)隨機(jī)選取(單位為ms)。
6)由于LM算法和Taylor算法需要聲音的速度參量,模擬現(xiàn)場(chǎng)測(cè)量的氣象條件為:氣溫20 ℃,風(fēng)速10 m/s,風(fēng)向?yàn)棣?3 。野外環(huán)境下氣象參數(shù)通常會(huì)隨著時(shí)間而改變,假定氣溫變化范圍為±5°,風(fēng)速變化范圍為±6 m/s,風(fēng)向變化范圍為±π/8 。實(shí)驗(yàn)中根據(jù)文獻(xiàn)[6]方法計(jì)算現(xiàn)場(chǎng)聲速的估計(jì)值。
7)如果迭代后輸出的結(jié)果中包含如下情況:方位角超出其周期2π的整數(shù)倍或俯仰角超出其周期π的整數(shù)倍,則相應(yīng)加上或減去若干個(gè)各自的周期值,使得估計(jì)結(jié)果位于合理的角度范圍內(nèi)。否則,在統(tǒng)計(jì)測(cè)向誤差和精度時(shí)會(huì)產(chǎn)生很大的不合理偏差。
4.1.2 實(shí)驗(yàn)結(jié)果
傳感器數(shù)目在[4,8]范圍內(nèi)隨機(jī)選取時(shí)實(shí)驗(yàn)結(jié)果如圖2所示,傳感器數(shù)目在[8,12]范圍內(nèi)隨機(jī)選取時(shí)實(shí)驗(yàn)結(jié)果如圖3所示。由圖2、圖3可以看出,LM算法具有最佳的性能,Taylor算法次之,LLS算法性能表現(xiàn)最差。另外,當(dāng)傳感器數(shù)目較多時(shí),3種算法對(duì)方位角的估計(jì)性能幾乎相同,原因是傳感器的x軸和y軸坐標(biāo)均在半徑為15 m的圓內(nèi)取值,較多的傳感器數(shù)目使得各算法的方位角估計(jì)性能接近一致。
圖2 傳感器數(shù)目較少時(shí)算法DOA估計(jì)誤差的CDF值
圖3 傳感器數(shù)目較多時(shí)算法DOA估計(jì)誤差的CDF值
本次實(shí)驗(yàn)對(duì)Taylor算法和本文LM算法的不收斂率與平均迭代次數(shù)進(jìn)行統(tǒng)計(jì),以比較算法的性能。
傳感器數(shù)目較少時(shí)的實(shí)驗(yàn)結(jié)果如表1所示,傳感器數(shù)目較多時(shí)的實(shí)驗(yàn)結(jié)果如表2所示。由表1、表2可以看出,相對(duì)Taylor算法,LM算法在不收斂率及平均迭代次數(shù)上性能優(yōu)越且穩(wěn)定。另外,在實(shí)驗(yàn)中測(cè)得傳感器數(shù)目為6時(shí),Taylor算法不收斂率很高,主要原因是傳感器數(shù)目少,LLS算法的DOA估計(jì)結(jié)果與真實(shí)值并不很接近,而Taylor算法對(duì)初值有很高的要求,當(dāng)LLS算法的DOA估計(jì)結(jié)果作為Taylor算法迭代的初值時(shí),必然會(huì)導(dǎo)致Taylor算法的不收斂率增高。
表1 傳感器數(shù)目較少時(shí)算法性能比較
表2 傳感器數(shù)目較多時(shí)算法性能比較
4.3.1 參數(shù)設(shè)置
本次實(shí)驗(yàn)以測(cè)向的均方根誤差(Root Mean Square Error,RMSE)作為性能評(píng)估指標(biāo)。RMSE描述方位估計(jì)偏差[16],能夠反映測(cè)向的精度。實(shí)驗(yàn)主要參數(shù)設(shè)置如下:
1)TDOA測(cè)量誤差服從高斯分布,不考慮測(cè)量的系統(tǒng)偏差問題。標(biāo)準(zhǔn)偏差取值范圍為[10,20](單位為ms)。
2)傳感器數(shù)目固定為8個(gè),當(dāng)?shù)谝粋€(gè)傳感器坐標(biāo)平移后設(shè)置為坐標(biāo)原點(diǎn)時(shí),一組隨機(jī)生成的陣列坐標(biāo)如下:(0,0,0),(0.63,1.70,-3.45),(6.56,9.85,-6.38),(1.60,15.33,0.45),(7.66,5.28,1.55),(18.22,8.28,-5.55),(-1.96,24.40,-3.21),(14.37,11.49,1.54) (單位為m)。
3)輻射源方位角與俯仰角取值均為45°。
4.3.2 實(shí)驗(yàn)結(jié)果
各算法RMSE實(shí)驗(yàn)結(jié)果如圖4所示,其中,CRLB線為克拉姆-拉奧下界曲線。由圖4可以看出,在測(cè)量噪聲較高時(shí),LM算法的性能能夠達(dá)到測(cè)向誤差的CRLB,表明其具有較高的有效性和魯棒性。
圖4 不同算法DOA估計(jì)的RMSE值
本文研究聲源測(cè)向問題,提出一種基于萊溫伯格-馬夸特的測(cè)向算法LM。該算法能夠較好地解決TDOA測(cè)量噪聲干擾測(cè)向精度的問題,并且避免Taylor算法對(duì)初值要求較高及結(jié)果不收斂的現(xiàn)象。實(shí)驗(yàn)結(jié)果表明,相對(duì)LLS算法、Taylor算法,該算法性能優(yōu)越,能夠達(dá)到克拉姆-拉奧下界。但在某些環(huán)境下,若聲速偏差較大會(huì)嚴(yán)重干擾本文算法的測(cè)向精度,解決該問題將是今后的研究重點(diǎn)。