王緒虎，白浩東，張群飛，田 雨

(1.青島理工大學(xué) 信息與控制工程學(xué)院，青島 266520；2.西北工業(yè)大學(xué) 航海學(xué)院，西安 710072)

波達(dá)方向(direction of arrival,DOA)估計(jì)問(wèn)題(在很多領(lǐng)域受到極大的關(guān)注，例如，雷達(dá)[1]、聲吶[2-4]，生物醫(yī)學(xué)[5]等。在過(guò)去幾十年已經(jīng)提出了許多高分辨的DOA估計(jì)算法，例如，多重信號(hào)分類算法[6](multiple signal classification,MUSIC)、求根MUSIC算法[7](Root-MUSIC)、旋轉(zhuǎn)不變子空間算法[8](estimating signal parameters via rotational invariance techniques,ESPRIT)等，但上述高分辨算法在低信噪比、少快拍數(shù)的情況下，估計(jì)性能會(huì)嚴(yán)重下降，甚至失效。

隨著壓縮感知(compress sensing,CS)理論和稀疏重構(gòu)技術(shù)的不斷發(fā)展和成熟，許多學(xué)者將其與DOA估計(jì)聯(lián)系起來(lái)，使DOA估計(jì)技術(shù)進(jìn)入了新的發(fā)展階段。相較于傳統(tǒng)的高分辨算法，基于CS和稀疏重構(gòu)技術(shù)的DOA估計(jì)方法，在低信噪比、少快拍的情況下表現(xiàn)出良好的估計(jì)性能。該種方法大致分為三類：凸優(yōu)化方法、貪婪方法[9]、稀疏貝葉斯學(xué)習(xí)方法[10-11]。例如，文獻(xiàn)[12]研究稀疏信號(hào)的表示以及DOA估計(jì)問(wèn)題，提出L1范數(shù)奇異值分解(L1norm-singular value decomposition,L1-SVD)方法，把測(cè)向問(wèn)題轉(zhuǎn)化為求解L1范數(shù)的問(wèn)題，但是其優(yōu)化問(wèn)題受到模型正則化參數(shù)的影響，影響估計(jì)的精度。王偉東等提出的基于稀疏信號(hào)功率迭代補(bǔ)償?shù)腄OA估計(jì)算法，該算法基于補(bǔ)償原理，使信號(hào)的稀疏表示近似于L0范數(shù)，把DOA估計(jì)問(wèn)題轉(zhuǎn)化為求解L0范數(shù)問(wèn)題。相比于L1-SVD算法，該算法不需要設(shè)置正則化參數(shù)，對(duì)非相干源具有較高的估計(jì)精度。然而，凸優(yōu)化方法的計(jì)算效率限制了其進(jìn)一步發(fā)展。

隨著CS理論研究的不斷深入，稀疏貝葉斯學(xué)習(xí)方法(sparse Bayesian learning,SBL)被認(rèn)為與凸優(yōu)化方法具有相同的全局收斂性，而它計(jì)算效率又遠(yuǎn)遠(yuǎn)優(yōu)于后者。Ji等將SBL引入CS領(lǐng)域中，提出貝葉斯壓縮感知算法。該算法通常需要滿足聲源的來(lái)波方向位于網(wǎng)格點(diǎn)上，才能實(shí)現(xiàn)較高的方位估計(jì)精度。當(dāng)入射聲源的波達(dá)方向偏離預(yù)先劃分好的網(wǎng)格，就導(dǎo)致方位估計(jì)精度的降低。因此，Yang等[13]提出了離網(wǎng)格的稀疏貝葉斯學(xué)習(xí)算法 (off-grid sparse Bayesian learning,OGSBL)，在離網(wǎng)格模型中，在信號(hào)真實(shí)到達(dá)角處采用一階泰勒展開式近似表示，使得估計(jì)性能進(jìn)一步提高。文獻(xiàn)[14]中，提出了求根離格稀疏貝葉斯學(xué)習(xí)算法(root-off grid sparse Bayesian learning,Root-OGSBL)，降低了OGSBL方法的計(jì)算復(fù)雜性，在網(wǎng)格間距較小的情況下，提高了計(jì)算效率的同時(shí)保證了估計(jì)的精度。文獻(xiàn)[15]中，使用平均場(chǎng)理論中的變分理論來(lái)估計(jì)參數(shù)的后驗(yàn)分布，并且提出了三層信號(hào)先驗(yàn)分布，這種分層先驗(yàn)進(jìn)一步提升了信號(hào)的稀疏性，提升了方位估計(jì)的精度。文獻(xiàn)[16]中，提出了稀疏信號(hào)新的先驗(yàn)分布框架—分層合成Lasso先驗(yàn)，相比于假設(shè)伽馬先驗(yàn)分布，具有更高的稀疏性、更低的重建誤差，提升了方位估計(jì)精度。

在實(shí)際的聲吶測(cè)向系統(tǒng)中，不可避免的存在陣列誤差。例如，水聽(tīng)器陣元間的耦合、水聽(tīng)器位置的偏差、水聽(tīng)器陣元的幅相通道不一致。導(dǎo)致現(xiàn)有常規(guī)算法的方位估計(jì)精度下降，甚至失效。在文獻(xiàn)[17-18]中，提出了陣元間存在未知互耦下波達(dá)方向估計(jì)算法，解決了陣元間存在未知互耦情況下，稀疏信號(hào)的恢復(fù)問(wèn)題。針對(duì)稀疏模型中入射信號(hào)方向與離散網(wǎng)格存在誤差，且陣元間的未知互耦沒(méi)有被同時(shí)考慮的問(wèn)題，本文提出了一種水聽(tīng)器陣元間存在不確定互耦，基于稀疏貝葉斯學(xué)習(xí)的離網(wǎng)格DOA估計(jì)方法。

1 信號(hào)模型

考慮K個(gè)不同方向的遠(yuǎn)場(chǎng)窄帶信號(hào)入射到如圖1所示的均勻線列陣，該陣列包含M個(gè)全向陣元，陣元之間存在互耦，陣元間距為d，因此，該接收模型可以表示為

圖1 M個(gè)水聽(tīng)器組成的均勻線列陣Fig.1 Uniform linear array composed of M-element hydrophones

n(t)

(1)

A=[a(θ1),a(θ2),…,a(θK)]

(2)

式中，a(θk)=[a1(θk),…,am(θk),…,aM(θk)]T為第k個(gè)入射信號(hào)的導(dǎo)向矢量，am(θk)=ej2π(m-1)d/λsin θk，λ為信號(hào)的波長(zhǎng)。

式(1)為單快拍觀測(cè)矢量模型，當(dāng)水聽(tīng)器陣元接收多快拍數(shù)據(jù)時(shí)，多快拍觀測(cè)矢量模型可以表示為

Y=CAS+N

(3)

式中，Y=[y(1),y(2),…,y(T)]∈CM×T為多快拍的接收數(shù)據(jù)矩陣，S=[s(1),s(2),…,s(T)]∈CK×T為多快拍的入射信號(hào)數(shù)據(jù)矩陣，N=[n(1),n(2),…,n(T)]∈CM×T為多快拍的陣列接收噪聲矩陣，T為快拍數(shù)。

在本文中，我們研究水聽(tīng)器陣元間存在未知耦合情況下的DOA估計(jì)問(wèn)題。由文獻(xiàn)[19]可知，水聽(tīng)器陣元間的互耦可以用矩陣C表示

(4)

該矩陣為對(duì)稱的托普利茲矩陣。根據(jù)文獻(xiàn)[20]的引理三，我們將上述模型進(jìn)行處理。首先，定義向量c=[1,c1,c2,…,cM-1]T，且滿足C=toeplitz(c)，toeplitz(c)為使用向量c生成的循環(huán)對(duì)稱矩陣 ，在t時(shí)刻的接收信號(hào)可以寫為

yt=CAst+nt=Q(st?c)+nt

(5)

則對(duì)應(yīng)水聽(tīng)器陣列的多快拍數(shù)據(jù)可表示為

Y=Q(S?c)+N

(6)

2 稀疏貝葉斯學(xué)習(xí)模型

2.1 稀疏信號(hào)模型

(7)

式(7)中的系統(tǒng)模型就可以轉(zhuǎn)化為過(guò)完備的稀疏模型

(8)

式中，矩陣X是矩陣S的稀疏擴(kuò)展

X=[x1,x2,…,xT]∈CN×T

(9)

圖2 信號(hào)X結(jié)構(gòu)稀疏化模型Fig.2 The structure of sparse matrix X

(10)

在實(shí)際的陣列測(cè)向系統(tǒng)中，接收信號(hào)并不是落在預(yù)設(shè)的網(wǎng)格點(diǎn)上，此時(shí)，采用離網(wǎng)格信號(hào)模型，使用離接收信號(hào)最近的網(wǎng)格矢量一階泰勒展開式近似表示

(11)

圖3 離格陣列空域模型Fig.3 Off-grid array spatial domain model

因此，由式(8)和(11)可以表示陣元間存在互耦下的離網(wǎng)格模型

Y≈Ω(δ)(X?c)+N

(12)

最終，水聽(tīng)器陣元間存在未知耦合的離網(wǎng)格稀疏模型建立如式(12)，接下來(lái)通過(guò)重建稀疏矩陣X來(lái)估計(jì)入射信號(hào)的波達(dá)方向，重建稀疏矩陣X中的非零行代表入射信號(hào)的實(shí)際波達(dá)方向。

2.2 模型參數(shù)分布假設(shè)

2.2.1 噪聲及其精度分布

假設(shè)接收的噪聲為獨(dú)立同分布的復(fù)對(duì)稱高斯白噪聲，其分布為

(13)

p(αn)=Γ(αn|a,b)

(14)

稀疏模型的似然函數(shù)可表示為

(15)

2.2.2 信號(hào)及其精度分布

假設(shè)稀疏信號(hào)X的多快拍樣本相互獨(dú)立，各列服從均值為0、方差為，

(16)

(17)

其中，c和d為伽馬分布的超參數(shù)。

2.2.3 互耦矢量及其精度分布

(18)

(19)

式中，e和f為伽馬分布的超參數(shù)。

2.2.4 離網(wǎng)格矢量分布

p(δ)=U(-r/2,r/2)

(20)

通過(guò)組合式(14)～(20)，可得模型的聯(lián)合概率密度

p(X,Y,δ,c,αx,αc,αn)=

p(Y|X,δ,c,αn)p(X|αx)p(c|αc)p(αn)

p(αx)p(αc)p(δ)

(21)

3 稀疏貝葉斯參數(shù)學(xué)習(xí)

利用上述的分布假設(shè)，我們提出了一種水聽(tīng)器陣元間存在互耦，基于稀疏貝葉斯學(xué)習(xí)的離網(wǎng)格波達(dá)方向估計(jì)算法。該算法中，我們使用EM算法學(xué)習(xí)更新離網(wǎng)格誤差矢量和互耦矢量，直至收斂。

首先，利用上述接收的數(shù)據(jù)Y和參數(shù)先驗(yàn)分布，可以得出信號(hào)X的后驗(yàn)分布

p(X|Y,δ,c,αx,αc,αn)=

αn)p(X|αx)

(22)

式中，p(X|αx)的概率密度函數(shù)為

(23)

X后驗(yàn)分布也為高斯分布，滿足

(24)

式中，μt為均值向量，Σx為協(xié)方差矩陣

(25)

Σx=[αnΥH(δ,c)Υ(δ,c)+diag(αx)]-1

(26)

其中，Υ(δ,c)=Ω(δ)(IN?c)，另外，定義所有快拍下均值矢量組成的矩陣

μ[μ1,…,μt,…,μT]

(27)

其次，我們對(duì)超參數(shù)進(jìn)行學(xué)習(xí)，超參數(shù)的學(xué)習(xí)過(guò)程等價(jià)于最大化所有待估計(jì)參數(shù)的后驗(yàn)概率分布

(28)

最大化待估計(jì)參數(shù)的后驗(yàn)概率分布又等價(jià)于最大化似然函數(shù)

L(δ,c,αx,αc,αn)=

ε{lnp(Y|X,δ,c,αn)p(X|αx)p(c|αc)

p(αn)p(αx)p(αc)p(δ)}

(29)

式中，ε{·}表示條件后驗(yàn)期望EX|Y,δ,c,αx,αc,αn{·}。接下來(lái)，詳細(xì)推導(dǎo)超參數(shù)迭代更新的表達(dá)式，下述推導(dǎo)中涉及到矩陣對(duì)向量變量的求導(dǎo)，參考附錄B。

忽略式(29)中與互耦矢量無(wú)關(guān)的概率項(xiàng)，可得到關(guān)于互耦矢量的似然函數(shù)

L(c)=ε{lnp(Y|X,δ,c,αn)p(c|αc)}=

(30)

因?yàn)镋M算法是不斷收斂遞增的，所以估計(jì)參數(shù)就等價(jià)于求似然函數(shù)的極值，有如下關(guān)系

(31)

因此，對(duì)式(31)求導(dǎo)，可得互耦矢量的迭代表達(dá)式(32)

(32)

令式(32)等于零，可得

(33)

式中，yt表示陣列接收數(shù)據(jù)Y的第t列，μt表示稀疏信號(hào)的均值矩陣μ的第t列。

同理，忽略與信號(hào)精度無(wú)關(guān)的概率分布，可得到如下似然函數(shù)表達(dá)式

L(αx)=ε{lnp(X|αx)p(αx)}

(34)

對(duì)似然函數(shù)求導(dǎo)，令其等于零，可得信號(hào)精度的迭代表達(dá)式

(35)

Σx(n,n)表示協(xié)方差矩陣的第n行n列的元素，μn,t表示稀疏信號(hào)的均值矩陣μ的第n行t列的元素。

同理，對(duì)于噪聲精度，我們推導(dǎo)其似然函數(shù)和迭代表達(dá)式，如下

L(αn)=ε{lnp(Y|X,δ,c,αn)p(αn)}

(36)

(37)

對(duì)于互耦精度，給出似然函數(shù)和迭代表達(dá)式

L(αc)=ε{p(c|αc)p(αc)}

(38)

(39)

對(duì)于離網(wǎng)格間距，給出似然函數(shù)和迭代表達(dá)式

L(δ)=ε{lnp(Y|X,δ,c,αn)p(δ)}

(40)

(41)

其中，P∈RN×N，v∈RN×1下式(42)(43)給出

(42)

(43)

最后，基于上述分析算法總結(jié)如下：

第一步：輸入水聽(tīng)器陣元接收數(shù)據(jù)Y；

第二步：初始化需要更新的參數(shù)X,δ,c,αx,αc,αn=mean(var(Y))，超參數(shù)a,b,c,d,e,f的設(shè)置[15]；

第三步：根據(jù)(25)(26)更新后驗(yàn)期望矢量和協(xié)方差矩陣；

第四步：根據(jù)(35)更新信號(hào)精度αx；

第五步：根據(jù)(37)更新噪聲精度αn；

第六步：根據(jù)(33)和(39)更新互耦矢量及其精度；

第七步：根據(jù)(41)更新離網(wǎng)格間隔矢量；

4 仿真分析

4.1 空間譜分析

根據(jù)文獻(xiàn)[21]，水聽(tīng)器間的耦合效應(yīng)可以表示為

(44)

其中，互耦系數(shù)中的幅度ξ服從均勻分布，即ξ～U[-0.05,0.05]，互耦系數(shù)的相位φ也服從均勻分布，即φ～U[0,2π]，參數(shù)αc為耦合系數(shù)，表征相鄰陣元之間的耦合效應(yīng)。理論上，相鄰陣元間的耦合效應(yīng)是相同的，仿真參數(shù)的設(shè)置如表1所示。

表1 仿真參數(shù)Tab.1 Simulation parameters

陣元間的耦合系數(shù)αc為-15 dB，正則化參數(shù)為2時(shí)，四種方法歸一化空間譜如圖4所示；陣元間的耦合系數(shù)αc為-15 dB，正則化參數(shù)為3時(shí)，四種方法歸一化空間譜如圖5所示；陣元間的耦合系數(shù)αc為-5 dB，正則化參數(shù)為3時(shí)，四種方法歸一化空間譜如圖6所示。

從圖4至圖6可以看出MUSIC方法的空間譜主瓣寬度較寬，旁瓣較高，當(dāng)增大耦合系數(shù)時(shí)，兩個(gè)主瓣的峰值和真實(shí)的角度值相差較大；對(duì)比圖4和圖5，可以看出當(dāng)L1-SVD方法正則化參數(shù)為2時(shí)，主瓣寬度較窄，旁瓣較低；改變正則化參數(shù)為3時(shí)，主瓣寬度變寬，且出現(xiàn)了許多旁瓣；對(duì)比圖5和圖6可以看出，當(dāng)耦合系數(shù)增大時(shí)，OGSBL方法旁瓣數(shù)量增加，且旁瓣譜級(jí)升高；對(duì)比圖4～圖6，本文提出的方法主瓣寬度和旁瓣譜級(jí)基本不變。

圖4 空間譜(正則化參數(shù)為2，αc=-15 dB)Fig.4 Spatial spectrum(Regularization parameters 2, αc=-15 dB)

圖5 空間譜(正則化參數(shù)為3，αc=-15 dB)Fig.5 Spatial spectrum(Regularization parameters 3， αc=-15 dB)

圖6 空間譜(正則化參數(shù)為3，αc=-5 dB)Fig.6 Spatial spectrum(Regularization parameters 3， αc=-5 dB)

四種方法的蒙特卡洛仿真結(jié)果如圖7所示，其中L1-SVD方法的正則化參數(shù)設(shè)置為3。圖7(a)是耦合系數(shù)為-15 dB時(shí)、圖7(b)是耦合系數(shù)為-10 dB時(shí)、圖7(c)是耦合系數(shù)為-5 dB時(shí)50次蒙特卡洛仿真的統(tǒng)計(jì)結(jié)果。從圖7可以看出，水聽(tīng)器陣元間存在耦合時(shí)，本文方法的估計(jì)性能和穩(wěn)健性是最優(yōu)的。

(a) 耦合系數(shù)αc=-15 dB

L1-SVD方法在不同正則化參數(shù)時(shí)的蒙特卡洛仿真結(jié)果如圖8所示。圖8(a)為正則化參數(shù)為2時(shí)、圖8(b)為正則化參數(shù)為3時(shí)，50次蒙特卡洛仿真統(tǒng)計(jì)結(jié)果。從圖8可以看出，改變正則化參數(shù)對(duì)L1-SVD方法的估計(jì)性能影響較大，因而該方法的實(shí)用性較差。

(a) 正則化參數(shù)為2

綜上所述，當(dāng)水聽(tīng)器陣元間存在互耦時(shí)，相比于MUSIC、L1-SVD、OGSBL方法，本文提出方法的估計(jì)性能最優(yōu)，穩(wěn)健性最高。

陣元間耦合系數(shù)為-5 dB，正則化參數(shù)為3時(shí)，四種方法信號(hào)估計(jì)結(jié)果如表2所示，表中的估計(jì)結(jié)果為50次蒙特卡洛仿真試驗(yàn)的統(tǒng)計(jì)均值。從表中可以看出，本文方法的估計(jì)結(jié)果與真實(shí)角度的偏差最小，表明本文方法的估計(jì)性能優(yōu)于MUSIC、L1-SVD、OGSBL方法。

表2 真實(shí)角度與本文方法估計(jì)DOA結(jié)果對(duì)比Tab.2 Comparison of the real angle and the DOA estimation result of the algorithm proposed in this paper

4.2 分辨力分析

不同方法的分辨能力不同，本小節(jié)從信噪比的高低、陣元間的耦合系數(shù)的大小研究本文提出方法的角度分辨能力。當(dāng)滿足下述的條件時(shí)，則認(rèn)為能夠正確分辨兩個(gè)角度[22]。

(45)

仿真試驗(yàn)條件如表1所示，設(shè)定耦合系數(shù)αc=-7 dB，圖9給出了分辨成功率隨角度間隔的變化規(guī)律，其中(a)和(b)子圖分別為信噪比為0和5 dB時(shí)的統(tǒng)計(jì)結(jié)果。從圖中可以看出，信噪比為0時(shí)，三種方法成功估計(jì)出兩個(gè)角度的間隔分別為9°(本文方法)，10°(OGSBL)，12°(MUSIC)。本文方法的角度分辨力優(yōu)于OGSBL和MUSIC方法。從圖9中可以看出，信噪比為0和5 dB時(shí)，可成功估計(jì)兩個(gè)角度的間隔為9°和7°。隨著信噪比的增加，本文方法的分辨力逐漸增強(qiáng)。

(a) SNR=0 dB 分辨成功率和角度間隔的關(guān)系

仿真試驗(yàn)條件如表1所示，設(shè)定信噪比為5 dB，圖10給出了分辨成功率隨角度間隔的變化規(guī)律，其中(a)和(b)子圖分別為耦合系數(shù)為-12 dB和-3 dB時(shí)的統(tǒng)計(jì)結(jié)果。從圖中可以看出耦合系數(shù)為-12 dB時(shí)，三種方法成功估計(jì)出兩個(gè)角度的間隔分別為6°(本文方法)，7°(OGSBL)，8°(MUSIC)。本文方法的角度分辨力優(yōu)于OGSBL和MUSIC方法。從圖10中可以看出，耦合系數(shù)為-12 dB和-3 dB時(shí)，可成功估計(jì)兩個(gè)角度的間隔為6°和11°。隨著耦合系數(shù)的增加，本文方法的分辨力逐漸降低。

(a) 耦合系數(shù)αc=-12 dB，分辨成功率和角度間隔的關(guān)系

仿真試驗(yàn)條件如表1所示，設(shè)定入射角度間隔為8°和耦合系數(shù)為-7 dB時(shí)，圖11給出了分辨成功率隨信噪比的變化規(guī)律；設(shè)定入射角度間隔為11°和信噪比為5 dB時(shí)，圖12給出了分辨成功率隨耦合系數(shù)的變化規(guī)律。

圖11 分辨成功率和信噪比的關(guān)系Fig.11 The relationship between resolution probability and SNR

從圖11可以看出隨著信噪比的增加，三種方法的分辨成功率逐漸增加。從圖中可以看出，角度間隔為8°時(shí)，本文方法和OGSBL方法可成功估計(jì)出兩個(gè)角度的信噪比分別為1 dB和4 dB，而MUSIC算法無(wú)法完全正確估計(jì)出兩個(gè)真實(shí)的入射角度。本文方法的可成功估計(jì)兩個(gè)信號(hào)的信噪比低于OGSBL方法，本文方法的分辨力優(yōu)于OGSBL方法和MUSIC方法。

從圖12可以看出隨著耦合系數(shù)的增加，三種方法的分辨成功率逐漸降低。從圖中可以看出，角度間隔為11°時(shí)，三種方法成功估計(jì)出兩個(gè)角度的耦合系數(shù)分別為-3 dB(本文方法)，-7 dB(OGSBL)，-7 dB(MUSIC)。本文方法的可成功估計(jì)兩個(gè)信號(hào)的耦合系數(shù)高于OGSBL方法和MUSIC方法，本文方法的分辨力優(yōu)于OGSBL方法和MUSIC方法。

圖12 分辨成功率和耦合系數(shù)的關(guān)系Fig.12 The relationship between resolution probability and coupling coefficient

綜上所述，本文方法的角度分辨力優(yōu)于OGSBL方法和MUSIC方法。當(dāng)耦合系數(shù)一定時(shí)，分辨成功率隨著信噪比的增加逐漸提升；當(dāng)信噪比一定時(shí)，分辨成功率隨著耦合系數(shù)的增加逐漸降低。

4.3 估計(jì)精度分析

本小節(jié)從信噪比、快拍數(shù)、網(wǎng)格間隔、陣元間的耦合系數(shù)研究不同方法的方位估計(jì)精度。把均方根誤差(Root Mean Square Error，RMSE)作為衡量方位估計(jì)精度的標(biāo)準(zhǔn)，可表示為

(46)

仿真試驗(yàn)條件如表1所示，設(shè)定信噪比為-2 dB到10 dB變化，步長(zhǎng)為2 dB，正則化參數(shù)為2，圖13給出了均方根誤差隨信噪比的變化規(guī)律。從圖中可以看出，當(dāng)信噪比為10 dB時(shí)，本文方法的估計(jì)誤差為0.4°低于0.45°(OGSBL)、0.6°(L1-SVD)、0.7°(MUSIC)。隨著信噪比的增大均方根誤差逐漸減小，且本文方法的估計(jì)精度高于OGSBL 方法、L1-SVD方法和MUSIC方法。

圖13 不同信噪比，四種方法的方位估計(jì)精度Fig.13 The DOA estimation performance of the four algorithms with different SNR

仿真試驗(yàn)條件如表1所示，設(shè)定快拍數(shù)為50到300變化，步長(zhǎng)為25，正則化參數(shù)為2，圖14給出了均方根誤差隨快拍數(shù)的變化規(guī)律。從圖中可以看出，當(dāng)快拍數(shù)為300時(shí)，本文方法的估計(jì)誤差為0.35°低于0.57°(OGSBL)、0.62°(L1-SVD)、0.9°(MUSIC)。隨著快拍數(shù)的增大均方根誤差逐漸減小，且本文方法的估計(jì)精度高于OGSBL 方法、L1-SVD方法和MUSIC方法。

圖14 不同快拍數(shù)下，四種方法的方位估計(jì)精度Fig.14 The DOA estimation performance of the four algorithms with different snapshot numbers

仿真試驗(yàn)條件如表1所示，設(shè)定網(wǎng)格間隔為{1°,1.5°,2°,2.5°,3°,4°,5°,6°}變化，圖15給均方根誤差隨網(wǎng)格間隔的變化規(guī)律。當(dāng)網(wǎng)格間隔小于2°時(shí)，本文提出的方法與OGSBL方法相比，估計(jì)精度相差不大，當(dāng)網(wǎng)格間隔增大到5°時(shí)，OGSBL估計(jì)誤差將增大，高于本文方法。隨著網(wǎng)格間隔增大，均方根誤差增大，但本文方法的估計(jì)精度要優(yōu)于OGSBL方法。

圖15 不同信噪比，網(wǎng)格間隔與方位估計(jì)精度的關(guān)系Fig.15 Comparison of the relationship between grid spacing and DOA estimation performance with different SNR

仿真試驗(yàn)條件如表1所示，設(shè)定耦合系數(shù)為-15 dB到-3 dB，步長(zhǎng)為2 dB，圖16給出了均方根誤差隨耦合系數(shù)變化的規(guī)律。當(dāng)信噪比一定時(shí)，增大耦合系數(shù)，估計(jì)誤差逐漸增大，估計(jì)精度逐漸降低；當(dāng)耦合系數(shù)為-3 dB時(shí)，增大信噪比為5 dB時(shí)，估計(jì)誤差減小3°，繼續(xù)增大信噪比為10 dB時(shí)，估計(jì)誤差減小1°。當(dāng)耦合系數(shù)一定時(shí)，隨著信噪比的增大，本文方法估計(jì)誤差逐漸減小，估計(jì)精度逐漸提升。

圖16 不同信噪比條件下，耦合系數(shù)與方位估計(jì)精度的關(guān)系Fig.16 Compare the relationship between the adjacent coupling coefficient and the DOA estimation performance with different SNR

5 結(jié) 論

為了解決水聽(tīng)器陣元間存在未知耦合，現(xiàn)有方法波達(dá)方向估計(jì)精度較低、角度分辨力差的問(wèn)題，本文提出了一種基于稀疏貝葉斯學(xué)習(xí)的DOA離網(wǎng)格估計(jì)方法。該方法通過(guò)建立新的貝葉斯學(xué)習(xí)模型，使用EM算法，推導(dǎo)離網(wǎng)格矢量和互耦矢量的分布，重新表示空間譜函數(shù)，最后通過(guò)譜峰搜索估計(jì)出波達(dá)方向。仿真試驗(yàn)結(jié)果表明：當(dāng)水聽(tīng)器陣元間存在互耦時(shí)，本文方法的方位估計(jì)性能和穩(wěn)健性優(yōu)于MUSIC、OGSBL和L1-SVD方法；本文方法的角度分辨力優(yōu)于MUSIC和OGSBL方法。鑒于上述試驗(yàn)使用的是兩個(gè)不相關(guān)信號(hào)，在后續(xù)的研究中，將會(huì)對(duì)相關(guān)信號(hào)做進(jìn)一步研究。

附錄A

引理1：矩陣C和向量c定義如正文，對(duì)于任意向量a(θk)，我們有：

Ca=Qk(a(θk))c

式中，Qk(a(θk))∈CM×M是兩個(gè)矩陣Qk1(a(θk))、Qk2(a(θk))的和：

因此，Qk(a(θk))=Qk1(a(θk))+Qk2(a(θk))，Qk1(a(θk))和Qk2(a(θk))分別為

根據(jù)上述引理1可得CAst=Q(st?c)，即

附錄B

有復(fù)向量u∈CP×1，v∈CP×1，復(fù)矩陣A∈CM×P是復(fù)向量x∈CN×1的函數(shù)，可得