淺海環(huán)境中被動定位方法研究及軟件實現(xiàn)

劉亞琴，劉 巖，祝 獻

（杭州應用聲學研究所，浙江 杭州 310023）

0 引言

目標被動定位是水聲學研究的一個重要研究領域，其在海防力量以及海洋強國建設中起著重要作用?，F(xiàn)有的定位方法很多，如三元陣被動測距、匹配場定位、匹配模定位、基于波導不變量和基于warping 變換以及?；盘柼幚硭暥ㄎ患夹g等，在利用水聲信號進行定位方面已經(jīng)取得顯著成果。然而，涉及具有彈性海底環(huán)境時，常規(guī)定位方法應用較少，如利用直達波與地聲界面波的傳播時延差對水中目標進行距離估計[1]，利用地震表面波垂直振動與水平振動之間的關系對目標進行方位角估計[2]，利用Scholte 波頻散特性對距離進行估計[3]，利用簡正波階數(shù)和地聲參數(shù)聯(lián)合起來作為反演參數(shù)進行反演，確定地聲參數(shù)以及聲源的距離等[4]。文獻[3]利用垂直波阻抗作為測量向量，基于粒子濾波算法，對具有彈性海底環(huán)境中的目標聲源進行了定位研究。在此基礎上，本文將該方法應用于液態(tài)Pekeris 波導，形成一套完整的、適合Pekeris 波導以及具有彈性海底的Pekeris 波導環(huán)境的定位軟件。該定位軟件將目標定位過程所涉及的水文環(huán)境參數(shù)、海底參數(shù)以及聲源信息、粒子濾波參數(shù)進行分區(qū)設置，操作簡單。

1 垂直波阻抗

對于水平分層介質波導，水中聲壓、垂直振速及垂直波阻抗可表示為[4-5]：

式中：zs表示聲源深度，z表示接收深度，r表示接收距離，S(ω)表示聲源譜，ω表示角頻率，ρ1表示水的密度，φn和ξn分別表示第n階本征函數(shù)和本征值，為第一類第0 階漢克爾函數(shù)，表示本征函數(shù)相對于深度z 的導數(shù)。

由于垂直波阻抗抵消了聲源譜信息，因此，當使用垂直波阻抗進行聲源定位，不需要考慮聲源頻譜信息。這對于被動定位來說是極其有利的。

2 粒子濾波

粒子濾波是一種基于蒙特卡羅方法（Monte Carlo Method）和遞推貝葉斯估計的統(tǒng)計濾波方法[6-9]。其核心思想是利用一些離散隨機采樣點（這些采樣點稱為粒子）來近似隨機變量的概率密度函數(shù)，以樣本均值代替貝葉斯估計中的積分運算，從而獲得狀態(tài)的估計值。粒子濾波算法在非線性、非高斯系統(tǒng)表現(xiàn)優(yōu)異，在處理海洋環(huán)境這一復雜系統(tǒng)時能夠提高估計性能。為了利用粒子濾波進行聲源定位，需要構建狀態(tài)方程和測量方程。其中，狀態(tài)方程描述狀態(tài)隨時間的變化過程，測量方程描述測量值與狀態(tài)向量之間的關系。狀態(tài)方程和測量方程分別為：

式中：狀態(tài)向量x=[r,zs]T，r和zs分別表示接收器與聲源的距離以及聲源深度；測量向量y=[Zz1,…,Zzm,…,ZzM]T，Zzm為角頻率ωm時垂直波阻抗；wk和vk分別為狀態(tài)噪聲向量和測量噪聲向量，對任意i,k滿足E{wkwiT}=Qkδki，E{vkviT}=Rkδki，E{vkviT}=0，這里δ(·)是Dirac delta 函數(shù)，Qk和Rk分別為k時刻狀態(tài)噪聲和測量噪聲的方差矩陣；d(·)為描述測量值與狀態(tài)向量之間的函數(shù)，d(xk)=[d1(x),…,dm(x),…,dM(x)]T，對于液態(tài)Pekiers 波導，如圖1（a）所示，=β1nm/(Hβ1nm-sinβ1nmHcosβ1nmHb2tanβ1nmHsin2β1nmH)，這里，ξnm為角頻率ωm時的本征值，H為海深，為水中垂直波數(shù)，k1nm=ωm/c1為水中波數(shù)，c1為水中聲速。

對于具有彈性海底Pekeris 波導，如圖1（b）所示，，

圖1 環(huán)境波導示意圖

這里，b=ρ1/ρ2為水中密度與海底密度之比，k2m=ωm/cp、χm=ωm/cs分別為半無限彈性海底縱波、橫波波數(shù)，cp、cs分別為彈性海底縱波、橫波聲速，分別為彈性海底縱波、橫波垂直波數(shù)，。

給出狀態(tài)方程和測量方程后，定位問題即可用數(shù)學模型描述為：已知測量值垂直波阻抗和狀態(tài)方程、測量方程，找到最優(yōu)估計參量值(r,zs)即可實現(xiàn)對聲源的定位。

假設狀態(tài)向量的初始概率密度函數(shù)Pd(x0)已知，它可以被假設為均勻分布或高斯分布等形式，或用數(shù)值技術估算。記Yk=[y1,y2,…,yk]為1 到k時刻測量向量組成的集合，Xk=[x1,x2,…,xk]。從貝葉斯理論角度出發(fā)，定位問題關鍵在于求解狀態(tài)向量x的后驗概率密度Pd(xk|Yk)，然后依據(jù)最小均方誤差或后驗概率最大準則估計目標參量最優(yōu)值。

基于蒙特卡羅積分和序貫重要性采樣理論，后驗概率可從概率密度分布q(Xk|Yk)采樣，得到：

3 定位軟件設計實現(xiàn)

定位軟件利用接收點處聲壓、垂直振速時域波形，獲取接收點處垂直波阻抗，將接收點處垂直波阻抗作為測量向量，基于粒子濾波算法實現(xiàn)定位。定位軟件實現(xiàn)流程如圖2 所示。

圖2 定位軟件實現(xiàn)流程圖

其中，接收點處聲壓時域波形、垂直振速時域波形可表示為：

則，接收點處垂直波阻抗可由聲壓、垂直振速時域信號獲得，如下：

將接收點處垂直波阻抗作為測量值，基于粒子濾波算法，即可完成定位估計。

定位軟件界面主要包括參數(shù)設置區(qū)和結果顯示區(qū)。參數(shù)設置區(qū)主要涉及水文環(huán)境及工作參數(shù)配置區(qū)、海底參數(shù)設置、聲源信號設置以及PF（Particle Filtering）參數(shù)設置。參數(shù)設置完成后，可通過“時域波形預報”“粒子濾波定位”按鈕分別完成聲壓、振速時域波形預報計算以及定位算法的運行，預報及定位結果在顯示區(qū)展示。初始界面如圖3 所示。

圖3 定位軟件初始界面

其中，水文環(huán)境及工作參數(shù)配置區(qū)主要設置海深、聲速剖面（Sound Speed Profile，SSP）、聲源深度、接收深度以及接收距離信息。聲速剖面內置良好水文、中等水文以及惡劣水文三種水文條件；海底參數(shù)設置主要包括縱波聲速、橫波聲速、縱波衰減、橫波衰減以及密度；聲源信號設置主要包括信號中心頻率、采樣頻率、信號幅度以及脈沖周期，聲源信號形式默認為Hanning 加權的四周期正弦波；PF 參數(shù)設置主要包括狀態(tài)向量初始值、粒子數(shù)目以及迭代次數(shù)。

所有參數(shù)設置完成后，點擊“時域波形預報”“粒子濾波定位”按鈕，完成時域波形預報與定位結果估計。圖4、圖5 分別為Pekeris 波導、具有彈性海底的Pekeris 波導環(huán)境下的計算結果。從圖4、圖5可以看出，兩種環(huán)境下，距離估計曲線、聲源深度估計曲線均收斂于設置值附近，表明定位軟件能夠有效實現(xiàn)目標的被動定位。

圖4 液態(tài)Pekeris 波導環(huán)境中定位結果

圖5 具有彈性海底Pekeris 波導環(huán)境中定位結果

4 結語

定位軟件基于粒子濾波定位算法，以頻率域上垂直波阻抗作為測量向量，構建粒子濾波測量方程，實現(xiàn)淺海環(huán)境聲源被動定位。定位軟件預報了聲源信號傳播到接收位置處的聲壓、垂直振速時域波形，在此基礎上，獲取頻率域上垂直波阻抗，基于粒子濾波算法完成距離與聲源深度估計。在兩種不同波導環(huán)境中進行了軟件的測試，測試結果表明在兩種不同的波導環(huán)境中，定位軟件計算的距離估計曲線、聲源深度估計曲線均收斂于真值附近，表明定位軟件能有效完成距離和聲源深度的估計。