生物聽覺耦合機制的矢量組合陣方位估計方法

方爾正，桂晨陽，孫純，陳峰，王歡

(1.哈爾濱工程大學 水聲技術重點實驗室，黑龍江 哈爾濱 150001；2.海洋信息獲取與安全工業和信息化部重點實驗室(哈爾濱工程大學)，黑龍江 哈爾濱 150001；3.哈爾濱工程大學 水聲工程學院，黑龍江 哈爾濱 150001)

水下目標方位估計是水聲陣列信號處理的重要研究方向之一，通常在空間分布若干傳感器采集聲場中的信號信息，然后利用方位估計技術提取信號的來波方向(direction of arrival, DOA)。常用的目標測向技術有基于常規波束形成(conventional beamforming, CBF)和最小方差無畸變響應(minimum variance distortionless response, MVDR)波束形成的DOA估計算法，以及子空間類DOA估計算法[1]。

隨著海洋國防科技的快速發展，基于水下潛器的自主作戰平臺已經成為各國海軍的重要發展領域之一。為了提高水下無人機的隱蔽性和靈活性，其工作平臺尺度越來越小，所搭載的聲吶陣列孔徑也隨之變小，基于小孔徑陣列的目標探測技術已然成為未來海洋軍事發展的重要內容之一。水聲工程實踐場景中，常會遇到陣列空間有限制的情況，為了提高小孔徑陣列的DOA估計性能，通常采取的方法是虛擬陣元技術和合成孔徑技術。前者通過不同陣元數據擬合，虛擬出新陣元擴展陣列，但是虛擬陣元信號與真實陣元信號相關性較高，無法帶來較高信噪比收益，因此適用于高信噪比情況[2]。后者是在平臺緩慢地勻速運動時，虛擬不同位置陣元，提高目標探測性能，主要針對勻速運動平臺[3-4]。倘若在小平臺直接搭載遠小于信號波長的陣列，陣元間相位差較小，導致CBF和MVDR方法的波束寬度過寬，陣增益較小，無法有效分辨水下多目標。

奧米亞寄生蠅的獨特聽覺機制為解決這一問題提供了新方向。奧米亞寄生蠅的雙耳距離約為1.2 mm，而蟋蟀叫聲的波長約為7 cm遠遠大于雙耳距離，但奧米亞寄生蠅卻能依靠聽覺定位寄主蟋蟀的位置進行寄生。Miles等[5]對奧米亞寄生蠅的聽覺系統展開了研究，揭示了奧米亞寄生蠅可以精確定位蟋蟀位置的原因，建立了一個機械耦合模型模擬寄生蠅的雙耳聽覺系統結構。基于此機械耦合模型，開展了仿奧米亞寄生蠅聽覺系統的聲源定位、超短傳感器陣列制作等研究內容[6-9]。然而，該耦合機制在放大信號的同時也將噪聲同等放大，導致在低信噪比條件下該方法的DOA估計算法估計性能欠佳[6]。

矢量傳感器可同步共點的拾取聲場中一點的聲壓和2個或3個相互正交的振速分量。相比傳統的聲壓傳感器，矢量傳感器可獲得更豐富的信息。文獻[10-13]的研究表明，對于有限尺度的聲源，在其遠程聲場，矢量傳感器接收的聲壓與振速信號是相干的，而對于各向同性噪聲場互不相關，從理論上揭示了聲強方法具有抑制噪聲能力。

本文受聲強和奧米亞寄生蠅聽覺系統啟發，提出了仿生物聽覺系統機制的小孔徑矢量組合陣列的DOA估計方法。本文采用一只矢量傳感器若干聲壓傳感器設計了一小孔徑矢量組合陣列，設計了一種聲強濾波器，結合奧米亞寄生蠅聽覺結構，將其拓展至多元小孔徑陣，推導出了耦合矩陣的具體表達式，最后利用CBF和MVDR方法獲得了空間譜圖，并進行仿真驗證所提方法的有效性。

1 組合陣接收信號模型

本文研究的是小孔徑矢量組合陣的方位估計問題。如圖1所示，組合陣由2U只標量傳感器和1只矢量傳感器組成，陣元間距為d，矢量傳感器位于坐標原點為“0”號陣元，原點兩側各布放U個聲壓傳感器，分別用“-U,…，-1”和1,…,U進行標號。假設K個不相關單頻信號s1(t),…,sK(t)分別從θ1,…,θK入射至矢量線陣，以坐標原點為參考點，則在t時刻，M元組合陣的接收數據可表示為(M=2U+1)：

圖1 矢量組合陣模型

(1)

式中：s(t)=[s1(t),s2(t),…,sK(t)]T表示K×1維信號向量；p(t)表示M×1維組合陣聲壓通道接收的數據；vx，0(t)和vy，0(t)表示1×1維矢量傳感器的2個振速通道接收數據；np(t)表示M×1維零均值的高斯白噪聲，各傳感器的噪聲互不相關，且與信號也不相關；nx，0(t)和ny，0(t)表示矢量傳感器2個振速通道接收的噪聲；A(θ)=[a-U(θ),…,a-1(θ),a0(θ),a1(θ),…,aU(θ)]T表示M×K維矢量陣聲壓分量的導向矢量矩陣。則標號為m的傳感器對應的導向矢量可表示為(m=-U,-U+1,…,U)：

am(θ)=[am(θ1)…am(θk)…am(θK)]

(2)

式中：am(θk)=ej2πmdcos θk/λk；λk=c/fk表示第k個單頻信號的波長;c表示水下聲速;fk表示第k個信號的頻率。

式中:ux和uy分別表示x、y振速分量的方向矢量，可表示為：

(3)

2 矢量聲強濾波器的設計

本節基于矢量水聽器設計了矢量聲強濾波器，然后對組合陣的聲壓通道數據進行濾波。

首先，將矢量陣接收的p0(t)、vx，0(t)和vy，0(t)做時域傅里葉變換，有：

(4)

式中：P0(f)、Vx，0(f)和Vy，0(f)分別表示p0(t)、vx，0(t)和vy，0(t)的傅里葉變換；S(f)表示傅里葉變換后的信號矢量；Np(f)、Nx(f)和Nx(f)表示噪聲的傅里葉變換。

(5)

根據文獻[12]，可得到廣義互相關濾波器，：

(6)

根據式(6)，可得應用于聲壓陣的濾波器：

Hz(f)=Hs(f)IM

(7)

式中IM表示M×M維單位矩陣。

利用上式求得的濾波器對組合陣聲壓通道進行濾波，可得濾波后的聲壓陣列的頻域信號：

(8)

3 仿生物聽覺耦合算法

3.1 寄生蠅雙耳耦合機制

如圖2所示，x1(t)和x2(t)表示蟋蟀聲波達到耳朵兩側鼓膜的聲波信號；y1(t)和y2(t)是奧米亞寄生蠅聽覺耦合響應后的信號，在機械系統中模擬的是鼓膜的位移；奧米亞寄生蠅的鼓室建模為2個剛性橋；m0為質量塊模擬雙耳兩側運動部分；奧米亞寄生蠅鼓膜和周圍結構的動力學特性分別用2組彈性和阻尼系數(α1,β1)、(α2,β2)模擬；鼓膜膜間橋用彈性和阻尼系數(α3,β3)的耦合連接2個剛性橋。當兩側鼓膜的聲波信號使鼓膜振動，鼓膜將發生位移，該機械耦合系統受到了彈性力、阻尼力以及慣性力的作用，其運動方程可表示為：

圖2 奧米亞寄生蠅雙耳機械耦合系統[8]

(9)

對式(9)進行傅里葉變換，可得到：

(10)

其中：

引入變量：

(11)

將式(11)代入式(10)，可得：

(12)

其中：

通常情況下寄生蠅雙耳為對稱，即α1=α2=α0，β1=β2=β0，因此D1(f)=D2(f)、B4(f)=B2(f)、B3(f)=-B1(f)、以及G(f)=2B1(f)B2(f)。以雙耳中心為參考點，設x1(t)=ejπfτ,x2(t)=e-jπfτ(τ為兩耳接收信號的時延差)，則雙耳輸出響應可表示為：

(13)

式中：τ=dccosθs/cc;dc=1.2 mm為雙耳距離;cc≈344 m/s為空氣聲速;θc表示聲波入設角度。

圖3、4給出了耦合后(信號Y1(f)、Y2(f))和沒有耦合(即雙耳接收信號X1(f)、X2(f))的幅度以及相位隨頻率的變化關系，其中入射角度為45°。從圖3、4可以看出，經過寄生蠅雙耳耦合放大系統后的2個信號幅度和相位的差別明顯增大，說明了該耦合放大系統增大了到達兩耳信號的幅度和相位差(等效于擴大雙耳間距)。

圖3 寄生蠅雙耳耦合和非耦合歸一化幅度隨頻率的變化

圖4 寄生蠅雙耳耦合和沒有耦合的相位隨頻率的變換

3.2 多元線陣耦合放大矩陣

本文根據上述原理，將寄生蠅雙耳耦合系統推廣至多元線陣。通過寄生蠅耦合方法對廣義矢量濾波后的信號進行耦合處理，增強2個相鄰陣元的。

根據式(13)不難發現，二元傳感器陣列的耦合放大矩陣的具體表達形式為：

Tc(f)=Hc(f)Λc

(14)

式中Hc(f)、Λc為2×2維矩陣，分別為：

其中：Γ(f)=α0+2α3+j2πf(β0+2β3)-m0(2πf)2，Υ(f)=α0+j2πfβ0-m0(2πf)2。

3.3 線陣耦合DOA算法

式(14)給出了相鄰2個陣元的耦合放大矩陣的具體表達式，下面針對M元線陣，推導出相應的耦合放大矩陣。具體做法為：以相鄰的2個陣元為一組將M元線陣分成N=M-1組，根據式(14)可以得到N組耦合放大矩陣，假設M元傳感器的幅度、相位響應均一致，則N組耦合放大矩陣也是相同的，因此，將這些耦合放大矩陣合并，可以得到M元線陣的耦合放大矩陣TM為：

TM(f)=H(f)Λ

(15)

式中H(f)和Λ分別為2N×2N和2N×N矩陣，具體為：

將式(15)左乘式(8)，便可得到經過耦合放大的陣列接收信號：

(16)

(17)

根據式(17)，可以得到耦合放大的陣列信號協方差矩陣：

(18)

(19)

式中θl表示空間掃描角度。由于奧米亞寄生蠅的聽覺耦合系統不僅增加了相位差，同時也增加了幅度差，為了降低幅度影響，對導向矢量進行去模處理為：

(20)

(21)

4 小孔經空間譜算法仿真

4.1 單目標空間譜對比分析

圖5給出了不同入射角度單目標的CBF和MVDR空間譜估計結果。信號頻率為1 000 Hz。信噪比為-10 dB。從圖5可以看出，傳統CBF和MVDR方法的旁瓣級都大于-3 dB，主瓣寬度很大。相比之下，基于本文所提方法的CBF和MVDR方法的旁瓣級很低，而且主瓣寬度較窄。產生這種現象的原因主要是本文設計了矢量聲強濾波器抑制了聲壓陣列的噪聲，然后又引入了生物聽覺系統的耦合放大機制放大了相鄰陣元的相位差。

圖5 不同入射角度的單聲源空間譜結果

4.2 多目標空間譜對比分析

圖6(a)給出了2個不相關等幅信號的CBF和MVDR空間譜估計結果。信號頻率分別為1 000、980 Hz，入射角度分別為70°和110°，2個信號的信噪比均為15 dB。從圖5(a)可以看出，傳統CBF方法無法分辨2個不相關信號的來波方向，而本文所提CBF方法有2個譜峰，可分辨2個不相關信號，但方位估計誤差較大。盡管傳統MVDR估計的空間譜有2個譜峰，但譜峰高度太低，且方位估計誤差較大，相比之下本文MVDR方法估計空間譜有2個超過25 dB的譜峰，可有效分辨2個不相關信號，而且方位估計精度較高。

圖6(b)給出了3個不相關目標的CBF和MVDR空間譜估計結果。信號頻率分別為1 000、980和1 080 Hz，對應的入射角度分別為55°、80°和110°，信噪比分別為15、15和9 dB。從圖可以看出傳統CBF、MVDR以及本文CBF方法都無法分辨3個不相關信號，相比之下本文MVDR方法有3個譜峰，可以分辨3個不相關信號，而且方位估計精度較高。

圖6 多目標的空間譜結果

本節仿真結果表明本文所提方法相比傳統方法更適合應用于水下小平臺的多目標測向。

4.3 方位估計性能分析

本節通過200次蒙特卡羅仿真實驗分析傳統CBF、MVDR方法與所提CBF、MVDR方法的目標測向的均方根誤差和目標分辨能力。

4.3.1 測向精度

本節采用的均方根誤差定義為：

(22)

圖7給出了入射角度為90°和50°的單目標方位估計的RMSE隨信噪比的變化情況。從圖7可以看出，隨著信噪比的增加，CBF和MVDR方法的均方根誤差逐漸減小，而且傳統CBF、MVDR與本文所提CBF、MVDR方法在相同信噪比條件下方位估計精度幾乎相等。對比圖7(a)和圖7(b)可以發現，當入射角度偏離陣列法線方向時，CBF和MVDR的測向結果變差，但對本文MVDR和傳統CBF和MVDR的影響較小。盡管仿生學的耦合機制在放大信號的同時也放大了噪聲，本文在耦合放大技術前加入的矢量聲強濾波可以減弱耦合放大技術對噪聲的放大作用，因此，所提方法在低信噪比條件下仍可保持較好的估計精度。

圖7 不同入射角度單聲源的測向RMSE隨信噪比的變化情況

圖8給出了CBF和MVDR方法對2個等幅目標測向RMSE隨信噪比的變化情況。2個不相關信號頻率分別為1 000和980 Hz，入射角度分別為70°和110°。從圖8可以看出傳統CBF方法的均方根誤差一直為20°，說明傳統CBF方法始終只有一個譜峰；隨著信噪比的增加，傳統MVDR和本文所提CBF和MVDR方法的均方根誤差逐漸減小；而且本文CBF相比傳統CBF具有較高的DOA估計精度，類似地，本文MVDR相比傳統的MVDR具有較高的DOA估計精度；此外由于MVDR分辨力較高，因而本文MVDR方法較本文CBF方法具有更高的方位估計精度。

圖8 入射角度為70°和110°的雙聲源測向誤差隨信噪比的變化情況

本節仿真結果說明在單目標情況下本文所提CBF、MVDR方法與傳統CBF和MVDR方法的測向精度近似相同，但在多目標情況下，本文所提CBF和MVDR的測向精度要優于傳統CBF和MVDR方法。

4.3.2 目標分辨能力

圖9給出了雙目標情況下CBF和MVDR方法的目標分辨概率隨2個信號入射角度間隔的變化情況。2個不相關信號頻率分別為1 000 Hz和980 Hz，入射角度分別為60°和60°+Δθ，Δθ從0°以3°為間隔變化至60°，信噪比為15 dB。從圖9可以看出，傳統CBF方法無法分辨2個不相關信號；本文的CBF方法在角度間隔不小于42°時便可分辨2個信號；傳統MVDR在角度間隔不小于39°時才可完全分辨2個信號；相比之下本文MVDR方法在角度間隔達到15°時便可100%分辨2個目標。

圖9 雙目標分辨概率隨角度間隔的變化

本節仿真結果說明本文所提CBF、MVDR方法較傳統方法具有更高的目標分辨能力，而且本文所提MVDR相比其他方法的目標分辨能力更強。

5 結論

1)本文通過利用設計的矢量聲強濾波器抑制了各向同性噪聲，提高了組合陣聲壓通道接收信號的信噪比，使所提CBF和MVDR方法在較低信噪比條件下仍可保持較好的DOA估計精度。

2)本文所提CBF和MVDR方法有效的結合了矢量聲強濾波器的抑制噪聲能力和寄生蠅聽覺系統的耦合放大機理，與傳統的CBF和MVDR方法相比，具有更低的旁瓣，更高的測向精度和目標分辨能力。

本文僅對小平臺中的單頻信號的DOA估計問題進行了研究，而水下信號還包括寬帶信號，又因為耦合參數選取與頻率息息相關，因此小孔徑陣的寬帶信號DOA估計問題是下一步要研究的內容。