王 進(jìn)，牛樹華

(海裝裝備審價(jià)中心,北京 100071)

近年來,無人艇、UUV、魚雷等平臺(tái)技術(shù)發(fā)展迅猛。對(duì)于這些尺寸受限的載體,聲吶系統(tǒng)布放的空間資源相對(duì)緊張。矢量水聽器可以同步獲取聲場(chǎng)的聲壓和振速數(shù)據(jù),因此,不用成陣僅依靠單個(gè)矢量水聽器即可實(shí)現(xiàn)對(duì)目標(biāo)的方位估計(jì)[1-3]。所以,與常規(guī)標(biāo)量探測(cè)設(shè)備相比,矢量水聽器更適合應(yīng)用于此類平臺(tái)。

方位估計(jì)依賴于信號(hào)頻率這一先驗(yàn)條件,但實(shí)際應(yīng)用中信號(hào)頻率通常未知,因此，針對(duì)方位頻率聯(lián)合估計(jì)問題開展研究具有現(xiàn)實(shí)意義。利用矢量水聽器獲得的聲壓振速信息,有利于實(shí)現(xiàn)對(duì)信號(hào)方位頻率的聯(lián)合估計(jì)。文獻(xiàn)[4-7]基于單矢量水聽器分別利用ESPIRT算法、波達(dá)方向矩陣法和平行因子模型對(duì)多個(gè)信號(hào)的方位角、俯仰角和頻率進(jìn)行了聯(lián)合估計(jì)。但由于算法模型的原因,這三種方法估計(jì)的目標(biāo)數(shù)目有限。文獻(xiàn)[8-9]分別使用狀態(tài)空間模型方法和稀疏分解理論獲得了單矢量水聽器接收信號(hào)的方位角和頻率參數(shù)信息。

在文獻(xiàn)[9]的研究基礎(chǔ)上,本文提出了基于稀疏分解理論的單矢量水聽器多維參數(shù)(即頻率、方位角和俯仰角)聯(lián)合估計(jì)方法。首先根據(jù)聲源信號(hào)的特征,建立單矢量水聽器的時(shí)空接收數(shù)據(jù)模型。然后構(gòu)造相同形式的過完備原子庫,利用稀疏分解理論進(jìn)行求解,由此可以得到聲源信號(hào)頻率、方位角和俯仰角的精確估計(jì)。

1 單矢量水聽器時(shí)空接收數(shù)據(jù)模型

假設(shè)遠(yuǎn)場(chǎng)有K個(gè)窄帶聲源信號(hào)入射到空間分布的某個(gè)三維矢量水聽器上,則該單矢量水聽器t時(shí)刻的接收信號(hào)可以表示為

Y(t)=[p(t),vx(t),vy(t),vz(t)]T=

U·X(t)+N(t)

(1)

式中,p(t)表示矢量水聽器的聲壓通道數(shù)據(jù),vx(t)、vy(t)和vz(t)分別為矢量水聽器的三個(gè)振速通道數(shù)據(jù);U=[u1,…,uk,…,uK]為4×K維的方向矩陣,其中uk=[1,sinφkcosθk,sinφksinθk,cosφk]T為三維矢量水聽器的方向向量,φk和θk分別表示第k個(gè)聲源的俯仰角和方位角;X(t)為聲源信號(hào)信息,N(t)為零均值的觀測(cè)加性高斯白噪聲。

對(duì)接收信號(hào)進(jìn)行采樣,采樣點(diǎn)為T個(gè),可以得到:

Y(n)=[p(n),vx(n),vy(n),vz(n)]T=

U·X(n)+N(n)

(2)

式中,n=1,2,…,T-1。根據(jù)時(shí)空變換技術(shù)(TST),將矢量水聽器各通道的時(shí)域接收數(shù)據(jù)按每一個(gè)采樣點(diǎn)進(jìn)行排列寫為式(3)。參照式(3)將式(2)重新整理,如式(4)所示。式(4)中,I是元素全為1的K維列向量,W為觀測(cè)噪聲。為了方便,假設(shè)聲源信號(hào)均為復(fù)正弦信號(hào),如式(5)所示。

Y=[p(0),vx(0),vy(0),vz(0),…,p(T-1),vx(T-1),vy(T-1),vz(T-1)]T

(3)

(4)

xk(t)=μkej(2πfkt+σk)+nk(t)=μkejσk·ej2πfkt+nk(t)

(5)

式(5)中,μk、fk和σk為第k個(gè)聲源信號(hào)的幅度、頻率和相位,k∈[1,K];nk(t)為高斯白噪聲。令E=[μ1ejσ1,μ2ejσ2,…,μKejσK]T,則式(4)可以寫為式(6)所示。

(6)

令Fk=[1,ej2πfk,…,ej2πfk(T-1)]T,因此,式(6)可以寫為如下形式:

Y=[F1?u1,F2?u2,…,FK?uK]E+W

(7)

式中,?為Kronecker積的運(yùn)算符號(hào);令A(yù)=[F1?u1,F2?u2,…,FK?uK],A包含了聲源信號(hào)的時(shí)域和空域信息,稱為時(shí)空陣列流型。式(7)為單矢量水聽器的時(shí)空接收數(shù)據(jù)模型。

2 多參數(shù)聯(lián)合估計(jì)方法

首先確定參數(shù)的搜索空間,對(duì)搜索空間根據(jù)精度需求選擇合適的間隔進(jìn)行網(wǎng)格劃分。假設(shè)頻率的搜索點(diǎn)數(shù)為Nf個(gè),方位角的搜索點(diǎn)數(shù)為Nθ個(gè),俯仰角的搜索點(diǎn)數(shù)為Nφ個(gè),則參數(shù)搜索空間維度為Ns=Nf×Nθ×Nφ。

假設(shè)每個(gè)網(wǎng)格點(diǎn)上都存在聲源,則式(7)可以改寫

成如下形式:

Y=AsZ+W

(8)

其中,

As=[a(f1,θ1,φ1),…,a(fNs,θNs,φNs)]
Z=[z1,…,zNs]T

(9)

式(8)中,a(fi,θi,φi)=Fi?ui,i∈[1,Ns];As為4T×Ns維度的矩陣,當(dāng)4T?Ns時(shí),As被稱為過完備原子庫。

但實(shí)際上,只有聲源真實(shí)存在的位置對(duì)應(yīng)的網(wǎng)格點(diǎn)上才有信號(hào),其余網(wǎng)格點(diǎn)上的信號(hào)強(qiáng)度為0。所以,Z是個(gè)稀疏度為K的稀疏向量,其非零元素的索引包含聲源的位置信息。式(8)將單矢量水聽器多維參數(shù)估計(jì)問題轉(zhuǎn)化成欠定線性方程組的求解問題。通常,欠定線性方程組有無窮個(gè)解,因此對(duì)該方程組的求解需要引入先驗(yàn)性稀疏約束,從而構(gòu)建稀疏分解模型,如下所示:

(10)

式中,‖·‖0為l0范數(shù),ε為稀疏分解誤差容忍項(xiàng),與數(shù)據(jù)所含噪聲有關(guān)。

綜上所述,本文提出的單矢量水聽器多維參數(shù)聯(lián)合估計(jì)方法具體流程歸納如下:

步驟1:將單矢量水聽器各通道接收的時(shí)域數(shù)據(jù)選取合適的采樣點(diǎn)數(shù)參照式(3)進(jìn)行整理,得到待處理的觀測(cè)數(shù)據(jù)Y;

步驟2:確定參數(shù)的搜索空間,根據(jù)估計(jì)精度需求選取合適的間隔對(duì)搜索網(wǎng)格進(jìn)行劃分,按照a(fi,θi,φi)=Fi?ui的形式構(gòu)造過完備原子庫As;

步驟3:將觀測(cè)數(shù)據(jù)Y和過完備原子庫As代入式(10),利用正交匹配追蹤(OMP)算法進(jìn)行求解,從而得到稀疏向量Z,將Z中非零元素的索引值進(jìn)行換算即可獲得對(duì)參數(shù)的估計(jì)。

3 數(shù)值仿真分析

3.1 方法性能分析

仿真條件1:假設(shè)有2個(gè)窄帶目標(biāo)源信號(hào)入射,幅值分別為2和1,頻率分別為32 Hz和80 Hz,方位角分別為83°和227°,俯仰角分別為75°和141°。采樣頻率為200 Hz。頻率搜索范圍為1 Hz～100 Hz,搜索間隔為1 Hz;方位角搜索范圍為1°～360°,俯仰角搜索范圍為1°～180°,搜索間隔均為1°。

基于仿真條件1,使用本文方法進(jìn)行處理,處理結(jié)果如圖1所示。

圖1 單矢量水聽器頻率-方位角-俯仰角聯(lián)合估計(jì)結(jié)果(2個(gè)目標(biāo))

由圖1a)可知,當(dāng)信噪比(SNR)為20 dB、快拍數(shù)(T)為500時(shí),本文方法可以實(shí)現(xiàn)對(duì)目標(biāo)頻率、方位角和俯仰角參數(shù)的準(zhǔn)確估計(jì)(本文以單矢量水聽器各通道的一個(gè)采樣點(diǎn)的數(shù)據(jù)作為一個(gè)快拍)。同時(shí),從圖中可以得到兩個(gè)目標(biāo)信號(hào)幅值的估計(jì)值分別約等于2和1,與仿真理論值相符。當(dāng)信噪比下降為5 dB且快拍數(shù)仍為500,仿真結(jié)果如圖1b)所示,此時(shí)方位角和俯仰角的估計(jì)已出現(xiàn)偏差,但頻率估計(jì)值仍是準(zhǔn)確的。保持信噪比為5 dB,增加快拍數(shù)至3 000,仿真結(jié)果如圖1c)所示,此時(shí)方位角和俯仰角的估計(jì)值與仿真理論值相符?？梢钥闯?信噪比和快拍數(shù)對(duì)本文方法的估計(jì)性能有重要的影響。

為方便比較,在上文仿真條件1中增加一個(gè)新目標(biāo):幅值為0.8,頻率為65 Hz,方位角為110°,俯仰角為30°。圖2展示了3個(gè)目標(biāo)時(shí)不同的信噪比和快拍數(shù)條件下的仿真結(jié)果。對(duì)比圖2b)和圖1b)可知,在同等信噪比和快拍數(shù)的條件下,目標(biāo)數(shù)增加會(huì)導(dǎo)致方位角和俯仰角估計(jì)結(jié)果出現(xiàn)誤差,頻率估計(jì)值依舊準(zhǔn)確。保持信噪比不變,增加快拍數(shù)至5 000,處理結(jié)果見圖2c)。在該條件下,三個(gè)目標(biāo)源的頻率、方位角和俯仰角參數(shù)估計(jì)都精確。保持信噪比不變,減少快拍數(shù)至10,處理結(jié)果見圖2a)。在該情形下,除了方位角和俯仰角之外,頻率的估計(jì)值也出現(xiàn)了誤差。將圖2的三個(gè)子圖進(jìn)行比較,可以發(fā)現(xiàn)頻率估計(jì)受快拍數(shù)影響較大。

圖2 單矢量水聽器頻率-方位角-俯仰角聯(lián)合估計(jì)結(jié)果(3個(gè)目標(biāo))

3.2 與傳統(tǒng)方法對(duì)比

仿真條件2:假設(shè)有2個(gè)目標(biāo)源信號(hào)入射,幅值分別為2和1,頻率分別為65 Hz和80 Hz,方位角分別為33°和107°,俯仰角分別為75°和141°。采樣頻率為200 Hz。為減少計(jì)算量,采用變步長(zhǎng)的搜索網(wǎng)格。頻率搜索范圍在64 Hz～66 Hz和79 Hz～81 Hz內(nèi)時(shí)搜索間隔為0.1 Hz,其余搜索間隔為1 Hz;方位角搜索范圍在32°～34°、106°～108°,俯仰角搜索范圍在74°～76°、140°～142°時(shí),搜索間隔為0.1°,其余搜索間隔為1°。蒙特卡洛試驗(yàn)次數(shù)為500。

基于仿真條件2分別使用本文方法和文獻(xiàn)[4]中的ESPIRT方法進(jìn)行處理。使用均方根誤差(RMSE)作為衡量標(biāo)準(zhǔn)比較兩種方法的估計(jì)性能:

(11)

兩種方法處理結(jié)果如圖3和圖4所示。圖3中展示了固定快拍數(shù)為100的條件下,兩種方法在不同信噪比時(shí)參數(shù)估計(jì)的RMSE曲線。圖4中展示了固定信噪比為10 dB的條件下,兩種方法在不同快拍數(shù)時(shí)參數(shù)估計(jì)的RMSE曲線。

圖3 不同信噪比下參數(shù)估計(jì)的均方根誤差

圖4 不同快拍數(shù)下參數(shù)估計(jì)的均方根誤差

圖3、4中使用本文方法估計(jì)頻率獲得的RMSE曲線分別在15 dB和100個(gè)快拍之后截?cái)?原因是在該條件下RMSE的估計(jì)值為0,而使用對(duì)數(shù)坐標(biāo)系繪制RMSE曲線無法處理0值,因此圖中呈現(xiàn)截?cái)酄顟B(tài)。由圖3、4可知,相比于傳統(tǒng)的ESPIRT方法,本文方法對(duì)方位角、俯仰角和頻率的估計(jì)精度更高。

4 結(jié)束語

針對(duì)單矢量水聽器多維參數(shù)估計(jì)問題,本文提出一種基于稀疏分解理論的單矢量水聽器多維參數(shù)聯(lián)合估計(jì)方法。通過數(shù)值仿真驗(yàn)證了該方法的有效性,同時(shí)得出了以下結(jié)論:1)目標(biāo)數(shù)、信噪比和快拍數(shù)對(duì)本文方法的估計(jì)性能有重要的影響。目標(biāo)數(shù)越大,維持較高的估計(jì)精度需要更大的信噪比或者更多的快拍數(shù)；2)相同條件下,本文方法相較于經(jīng)典的ESPIRT方法估計(jì)精度更高。