陳筱月, 馬曉川,*, 李 璇
(1. 中國(guó)科學(xué)院聲學(xué)研究所水下航行器信息技術(shù)重點(diǎn)實(shí)驗(yàn)室, 北京 100190;2. 中國(guó)科學(xué)院大學(xué)電子電氣與通信工程學(xué)院, 北京 100049)
來波方位(direction of arrival, DOA)估計(jì)作為陣列信號(hào)處理的研究熱點(diǎn)之一,在目標(biāo)跟蹤,聲源定位以及無線通信等多個(gè)領(lǐng)域當(dāng)中具有廣泛的應(yīng)用[1-4]。常規(guī)波束形成(conventional beamforming, CBF)[5]是最為經(jīng)典的DOA估計(jì)算法,其通過陣列加權(quán)得到指向性波束,調(diào)整加權(quán)向量使波束指向不同方位從而得到方位譜,方位譜峰值對(duì)應(yīng)的方位即真實(shí)來波方位。基于無失真約束與波束輸出噪聲最小準(zhǔn)則,Capon[6]提出了最小方差無失真方法(minimum variance distortionless response, MVDR),利用接收信號(hào)協(xié)方差矩陣自適應(yīng)的調(diào)節(jié)波束加權(quán)向量,提高了估計(jì)結(jié)果的準(zhǔn)確性。子空間類算法[7-8]利用噪聲子空間與信號(hào)子空間正交的特性,實(shí)現(xiàn)了超分辨。近年來,一些基于稀疏優(yōu)化、解卷積以及神經(jīng)網(wǎng)絡(luò)的DOA算法被提出。稀疏優(yōu)化類算法構(gòu)造過完備基矩陣將DOA估計(jì)視作優(yōu)化問題,并通過對(duì)解向量進(jìn)行稀疏先驗(yàn)約束從而獲得稀疏解[9-11];解卷積后處理算法基于逆濾波的原理,消除點(diǎn)擴(kuò)散函數(shù)造成的方位譜擴(kuò)展,從而得到理想方位譜[12-13];神經(jīng)網(wǎng)絡(luò)類算法將各個(gè)角度分為不同的類別,然后利用深度神經(jīng)網(wǎng)絡(luò)等方法進(jìn)行分類,從而獲得估計(jì)結(jié)果[14-15]。
上述大多方法獲得高精度DOA估計(jì)結(jié)果的前提是進(jìn)行密集的方位搜索。然而,密集的方位搜索會(huì)導(dǎo)致計(jì)算量的升高,特別是對(duì)于二維DOA來說[16],計(jì)算量會(huì)按平方急劇膨脹。而對(duì)于實(shí)際的聲納雷達(dá)系統(tǒng),由于實(shí)時(shí)性要求及硬件系統(tǒng)的限制,難以通過密集的方位掃描來獲得高精度的DOA結(jié)果。因此,實(shí)際工程當(dāng)中一般會(huì)采用幅度插值[17]的方法得到來波方位,該方法基于CBF,但與傳統(tǒng)密集譜峰搜索不同的是,該方法僅在少數(shù)方位上計(jì)算波束輸出結(jié)果,然后利用不同的函數(shù)對(duì)輸出結(jié)果進(jìn)行方位譜擬合,插值得到來波方位。其中,最常用的擬合函數(shù)有線性函數(shù)、二次函數(shù)、高斯函數(shù)[18]。研究者針對(duì)幅度比較測(cè)向法的誤差成因以及估計(jì)精度的提高做了大量相關(guān)工作。顧敏劍[19]使用高斯函數(shù)對(duì)波束圖進(jìn)行近似,理論分析了系統(tǒng)噪聲,通道幅度特性不一致以及波束軸角指向?qū)烙?jì)結(jié)果的影響。在高斯噪聲下,文獻(xiàn)[20]利用泰勒級(jí)數(shù)展開得到了振幅比較單脈沖雷達(dá)的DOA估計(jì)均方誤差的解析表達(dá)式,避免了使用蒙特卡羅計(jì)算均方誤差帶來的高計(jì)算量。文獻(xiàn)[21]中自適應(yīng)地使用相鄰或交替天線進(jìn)行DOA測(cè)量,有效的降低了由于噪聲與來波方位偏離波束指向所造成的誤差。文獻(xiàn)[22]中利用神經(jīng)網(wǎng)絡(luò)模型對(duì)一階泰勒展開引起的非線性誤差進(jìn)行補(bǔ)償,得到最終的DOA結(jié)果,提高了估計(jì)精度。文獻(xiàn)[17]中考慮波束偏離0°時(shí),左右波束的不對(duì)稱,利用波束圖-3 dB點(diǎn)作為求解條件,通過非對(duì)稱函數(shù)進(jìn)行方位譜擬合,得到了更好的估計(jì)結(jié)果。文獻(xiàn)[23]中,提出了一種基于雷達(dá)模式高斯斜率修正的優(yōu)化算法,取得了比傳統(tǒng)算法更好的精度。而文獻(xiàn)[24]直接在真實(shí)方向圖基礎(chǔ)上進(jìn)行公式推導(dǎo),進(jìn)而對(duì)目標(biāo)角度進(jìn)行求解。上述工作大多均在一維進(jìn)行,針對(duì)二維情況,文獻(xiàn)[25]通過蜂窩網(wǎng)結(jié)構(gòu)分布的二維多波束,使用二元二次方程進(jìn)行曲面擬合,得到了較準(zhǔn)確的測(cè)向結(jié)果。
本文首先通過對(duì)二維平面陣方位譜的推導(dǎo),將二維幅度插值DOA估計(jì)解耦為兩個(gè)相互獨(dú)立的一維幅度插值DOA估計(jì),從而避免了文獻(xiàn)[25]中使用多元高次方程進(jìn)行曲面擬合。對(duì)于解耦后的一維幅度插值DOA估計(jì),考慮波束指向偏離中心時(shí)左右波束不對(duì)稱,若使用對(duì)稱函數(shù)進(jìn)行插值會(huì)帶來誤差,所以使用了不同的非對(duì)稱函數(shù)進(jìn)行幅度插值,并且在角度求解的過程中,區(qū)別于文獻(xiàn)[17]中使用-3 dB波束寬度作為求解條件,本文提出以左右波束斜率比作為求解條件,提高了算法的性能表現(xiàn)。針對(duì)二維幅度插值中出現(xiàn)的角度失配問題,提出了以迭代插值的方式對(duì)估計(jì)結(jié)果進(jìn)行修正,進(jìn)一步地提高了角度估計(jì)精度。最后,通過仿真實(shí)驗(yàn)分析了不同插值方法的誤差曲線以及角度失配,信噪比對(duì)插值精度的影響,并通過對(duì)湖試數(shù)據(jù)的處理對(duì)算法進(jìn)行了驗(yàn)證,本文算法能夠應(yīng)用于對(duì)實(shí)時(shí)性有要求的系統(tǒng)中進(jìn)行高精度單目標(biāo)方位估計(jì)。
本文所使用的三維坐標(biāo)系統(tǒng)[26]如圖1所示,O為三維坐標(biāo)系的原點(diǎn),在此坐標(biāo)系統(tǒng)中定義空間角為Ω=(θ,φ),其中θ為Ω與平面yOz的夾角,范圍為[-90°,90°],而φ為Ω與平面xOz的夾角,范圍為[-90°,90°]。在上述定義下,θ與φ實(shí)際上無法同時(shí)取到各自的邊界值,兩者之間存在著如下的限制關(guān)系:

圖1 三維坐標(biāo)系示意圖Fig.1 Diagram of the three-dimensional coordinate system
|θ|+|φ|≤90°
(1)
即超出式限制的空間角Ω是不存在的。
二維平面陣位于xOy平面上,陣元個(gè)數(shù)為M,陣元坐標(biāo)為pm=[pmx,pmy,pmz]T,則陣元坐標(biāo)矩陣為P=[p1,p2,…,pM],在下文中,為了仿真的計(jì)算方便,使用的為四元矩形陣,且陣元間距d等于半波長(zhǎng),但是實(shí)際上本文方法對(duì)于任意二維平面陣均適用。現(xiàn)一窄帶平面波s(t)由遠(yuǎn)場(chǎng)入射至二維平面陣上,其頻率為f,波長(zhǎng)為λ,聲速為c,入射角度為Ω0=(θ0,φ0),其入射方向的單位方向矢量v0為
(2)
則陣列接收數(shù)據(jù)x(t)可以寫為
x(t)=a(θ0,φ0)s(t)+n(t)
(3)
式中:n(t)為加性高斯白噪聲;a(θ0,φ0)為(θ0,φ0)方向上的導(dǎo)向矢量,其表達(dá)式為

(4)
式(3)為連續(xù)形式,實(shí)際中會(huì)經(jīng)過采樣得到其離散形式如下:
x(n)=a(θ0,φ0)s(n)+n(n),n=1,2,…,N
(5)
式中:N為快拍數(shù);x(n)∈CM×N;a(θ0,φ0)∈CM×1;s(n)∈C1×N;n(n)∈CM×N。
二維平面陣通過對(duì)采樣得到的陣元接收信號(hào)進(jìn)行加權(quán)求和,即可實(shí)現(xiàn)波束形成,對(duì)期望方向的信號(hào)進(jìn)行輸出。常規(guī)波束形成加權(quán)向量表達(dá)式如下:
(6)
波束輸出信號(hào)為
y(θ,φ,n)=ω(θ,φ)Hx(n)
(7)
進(jìn)而可以計(jì)算出其功率為
(8)

為了應(yīng)對(duì)上述情況,實(shí)際系統(tǒng)中多采用幅度插值的方法進(jìn)行角度估計(jì)。其原理如圖2所示,圖中的曲面是按照0.1°的間隔對(duì)整個(gè)空間進(jìn)行掃描得到的方位幅度譜,可近似認(rèn)為是一個(gè)連續(xù)的曲面(這里需要注意的是,盡管根據(jù)式(1),圖中某些角度并不存在,但是為了之后插值與理解的方便,仍然計(jì)算了這些角度的波束輸出值),紅色星號(hào)所對(duì)應(yīng)的位置(20°,10°)為通過波束密集掃描估計(jì)的來波方位。幅度插值法僅在少數(shù)的空間角度計(jì)算波束輸出,如圖2中的品紅色圓點(diǎn)對(duì)應(yīng)的角度,該操作相當(dāng)于對(duì)真實(shí)的方位幅度譜進(jìn)行空間采樣。這里的空間采樣間隔則與陣列在水平豎直方向上的波束寬度有關(guān)。對(duì)于矩形陣來講,其在水平豎直兩個(gè)方向上的分辨率相同,所以在這兩個(gè)方向上的采樣間隔也相同,如圖3(a)所示。但對(duì)于任意陣,如2×4的五元L形陣,如圖4(a)所示,其在水平豎直兩個(gè)方向上的分辨率不同,因此在這兩個(gè)方向上的采樣間隔也不一樣,在分辨率更高水平的方向上則需要更密的采樣間隔,以免波束主瓣落在兩個(gè)采樣點(diǎn)之間,一般要求采樣間隔小于主瓣寬度的一半。之后利用曲面函數(shù)對(duì)最大采樣點(diǎn)及其周圍多個(gè)采樣點(diǎn)進(jìn)行曲面擬合,插值得到的曲面最大值所對(duì)應(yīng)的角度則為估計(jì)的來波角度。該方法避免了對(duì)二維空間進(jìn)行密集掃描,有效的降低了計(jì)算量。

圖2 二維幅度插值原理示意圖Fig.2 Two-dimensional amplitude interpolation principle diagram

圖3 二維方位幅度譜及其一維切片(矩形陣)Fig.3 Two-dimensional azimuth amplitude spectrum and its one-dimensional section (rectangular array)

圖4 二維方位幅度譜及其一維切片(2×4 L形陣)Fig.4 Two-dimensional azimuth amplitude wpectrum and its one-dimensional section (2×4 L-shaped array)
然而,進(jìn)行二維曲面擬合會(huì)面臨復(fù)雜的多元高次方程組[25],求解較為困難,且面對(duì)任意陣型時(shí),并不一定能夠找到合適的擬合函數(shù),從而會(huì)導(dǎo)致估計(jì)結(jié)果誤差增加。因此,本文提出將二維幅度插值解耦為兩個(gè)不同方向(θ方向,φ方向)的一維幅度插值,下面對(duì)該方式的可行性進(jìn)行驗(yàn)證。對(duì)于式(8),在忽略噪聲項(xiàng)后得到的幅度表達(dá)式如下:
(9)
進(jìn)一步將式(6)代入式(9)可以得到
(10)

然而,θ0與φ0是所需要估計(jì)的量,實(shí)際上其具體值事先未知,僅知道在方位幅度譜上的等間隔采樣值(即圖3(a)與圖4(a)中的黑色下三角點(diǎn)與品紅色圓點(diǎn)處),而其中黑色下三角點(diǎn)的交叉處則是這些采樣點(diǎn)中的最大幅度所對(duì)應(yīng)的角度,記作(θm,φm)。因此,將二維曲面插值分解為(θ,φm)與(θm,φ)兩個(gè)切面(黑色下三角點(diǎn)指示)上的一維曲線插值。圖3(b)、圖3(c)與圖4(b)、圖4(c)分別給出了四元矩形陣與2×4的五元L形陣的方位幅度譜及其在(θ0,φ0)的兩個(gè)方向上(紅色上三角)與(θm,φm)的兩個(gè)方向上(黑色下三角)的一維切片。可以看出,對(duì)于矩形陣而言,其方位幅度譜在(θ0,φ0)的兩個(gè)方向上與(θm,φm)的兩個(gè)方向上的一維切片雖然存在幅度上的差異,但是最大值均所對(duì)應(yīng)的角度相同。并且對(duì)于矩形陣來說,在任意角度的兩個(gè)方向上的一維切片均滿足此條性質(zhì)。然而,對(duì)于其他任意二維平面陣,該條性質(zhì)并不一定滿足,如圖4所示,L形陣方位幅度譜在(θ0,φ0)的兩個(gè)方向上與(θm,φm)的兩個(gè)方向上的一維切片極值點(diǎn)并不相同,若利用(θm,φm)的兩個(gè)方向上采樣點(diǎn)進(jìn)行插值,則會(huì)產(chǎn)生誤差,在本文中這種誤差被稱作角度失配誤差。為了使算法能適用于任意二維平面陣,后續(xù)章節(jié)提出了使用多次迭代插值的方法來降低該項(xiàng)誤差。


圖5 φ=φm的一維切面Fig.5 One-dimensional section of φ=φm
3.1.1 對(duì)稱線性插值
(11)

(12)

(13)
3.1.2 對(duì)稱二次插值
二次插值則是用二次函數(shù)對(duì)方位幅度譜進(jìn)行擬合,二次函數(shù)的表達(dá)式如下:
(14)

(15)
3.1.3 對(duì)稱高斯插值
高斯插值則是用高斯函數(shù)對(duì)方位幅度譜進(jìn)行擬合,其表達(dá)式如下:
(16)
(17)
上述方法均認(rèn)為對(duì)真實(shí)二維方位幅度譜進(jìn)行切片得到的一維方位幅度譜是關(guān)于其極大值左右對(duì)稱的。然而,事實(shí)上,這一點(diǎn)僅當(dāng)一維方位幅度譜的極大值位于0°的時(shí)候滿足,當(dāng)其極大值偏離0°的時(shí)候,左右兩側(cè)將不再對(duì)稱。如圖6所示,為真實(shí)來波位于不同方位時(shí)的一維方位幅度譜,可以看出,隨著來波角度偏離0°越嚴(yán)重,方位幅度譜的左右不對(duì)稱性就越來越嚴(yán)重。因此,如果用對(duì)稱函數(shù)進(jìn)行插值擬合,會(huì)給估計(jì)結(jié)果帶來較大的誤差,為了解決這個(gè)問題,下面提出使用非對(duì)稱函數(shù)進(jìn)行插值擬合,以提高估計(jì)精度。

圖6 不同來波方位下的方位幅度譜Fig.6 Azimuth amplitude spectrum in different wave directions
3.2.1 非對(duì)稱線性插值
非對(duì)稱線性插值認(rèn)為一維方位幅度譜在其極大值左右兩側(cè)的斜率不同,因此線性擬合函數(shù)如下:
(18)
將α、β、γ與相應(yīng)的波束幅度響應(yīng)Uα、Uβ、Uγ(Uα≥Uγ)代入式(18)可以得到方程組:
(19)
(20)

(21)

(22)

(23)
3.2.2 非對(duì)稱二次插值
非對(duì)稱二次插值使用兩個(gè)不同的二次函數(shù)對(duì)一維方位幅度譜極值點(diǎn)左右兩側(cè)進(jìn)行擬合,兩個(gè)二次函數(shù)僅在二次項(xiàng)系數(shù)上存在差別,其插值函數(shù)如下:
(24)
將α、β、γ與相應(yīng)的波束幅度響應(yīng)Uα、Uβ、Uγ(Uα≥Uγ)代入式(24)可以得到方程組:
(25)

(26)
A1=ηq(β)-1-M1ηq(β)+M1
A2=1-ηq(β)+M2ηq(β)-M2
B1=2β-2ηq(β)α+2M1ηq(β)α-2M1γ
B2=2ηq(β)α-2γ-2M2ηq(β)β+2M2γ
C1=ηq(β)α2-β2-M1ηq(β)α2+M1γ2
C2=γ2-ηq(β)α2+M2ηq(β)β2-M2γ2

(27)

(28)
3.2.3 非對(duì)稱高斯插值
非對(duì)稱高斯插值使用兩個(gè)不同的高斯函數(shù)對(duì)一維方位幅度譜極值點(diǎn)左右兩側(cè)進(jìn)行擬合,使用的非對(duì)稱高斯插值函數(shù)如下:
(29)
將α、β、γ與相應(yīng)的波束幅度響應(yīng)Uα、Uβ、Uγ(Uα≥Uγ)代入式(29)可以得到方程組:
(30)

(31)
M3=[ln(Uα/Uβ)]/[ln(Uα/Uγ)]
M4=ln(Uγ/Uα)/ln(Uγ/Uβ)
A3=-ηg(β)+1+M3ηg(β)-M3
A4=-1+ηg(β)-M4ηg(β)+M4
B3=-2β+2ηg(β)α-2M3ηg(β)α+2M3γ
B4=-2ηg(β)α+2γ+2M4ηg(β)β-2M4γ
C3=-ηg(β)α2+β2+M3ηg(β)α2-M3γ2
C4=-γ2+ηg(β)α2-M4ηg(β)β2+M4γ2

(32)

(33)

(34)
下面將對(duì)二維幅度插值方法的估計(jì)誤差進(jìn)行討論,并分析不同因素對(duì)算法估計(jì)結(jié)果精度的影響。
仿真使用的陣列為四元矩形陣并按照半波長(zhǎng)布陣,來波信噪比為20 dB。在與z軸夾角小于50°的區(qū)域當(dāng)中,按照10°的間隔進(jìn)行幅度采樣,使用上面提出的6種插值方法分別進(jìn)行角度估計(jì),并通過100次蒙特卡羅實(shí)驗(yàn)得到該區(qū)域內(nèi)的二維誤差曲面。圖7給出了非對(duì)稱高斯插值的估計(jì)誤差曲面,為了方便對(duì)不同方法的估計(jì)精度進(jìn)行觀察與對(duì)比,之后不依次展示各個(gè)方法的估計(jì)誤差曲面,而是展示各方法二維誤差曲面在φ=0°處的切面,如圖8所示。

圖7 非對(duì)稱高斯插值二維誤差曲面Fig.7 Two-dimensional error surface of asymmetric Gaussian interpolation

圖8 不同插值方法方位估計(jì)誤差曲線(φ=0°切面)Fig.8 Azimuth estimation error curves of different interpolation methods (φ=0° section)
可以看出,隨著入射角度逐漸偏離0°,估計(jì)誤差呈現(xiàn)振蕩上升。對(duì)稱插值方法中,由于線性函數(shù)與方位譜的相似程度較差,其估計(jì)結(jié)果誤差最大,而二次函數(shù),高斯函數(shù)與方位譜具有較好的相似度,其擬合結(jié)果更加準(zhǔn)確。而非對(duì)稱插值相對(duì)對(duì)稱插值,線性函數(shù)的整體估計(jì)誤差并沒有明顯的改變,只是其誤差曲線的極值點(diǎn)出現(xiàn)了移動(dòng),二次函數(shù)與高斯函數(shù)的誤差出現(xiàn)了較為明顯的下降,得到了更加精確的估計(jì)結(jié)果。進(jìn)一步可以發(fā)現(xiàn),非對(duì)稱二次插值與非對(duì)稱高斯插值的誤差曲線的極小值與幅度采樣角度一致,而誤差曲線的極大值出現(xiàn)在兩個(gè)相鄰采樣角度中間,即當(dāng)入射信號(hào)方位在幅度采樣角度附近的時(shí)候可以得到更小的估計(jì)誤差。

由圖8(b)可以看出,非對(duì)稱二次插值與非對(duì)稱高斯插值的誤差曲線極大值點(diǎn)出現(xiàn)在兩個(gè)相鄰的幅度采樣角度中點(diǎn)處,極小值點(diǎn)則剛好出現(xiàn)在幅度采樣角度上。這是因?yàn)楫?dāng)真實(shí)角度剛好位于兩個(gè)相鄰幅度采樣點(diǎn)中點(diǎn)的時(shí)候,真實(shí)角度離兩側(cè)的幅度采樣點(diǎn)距離最遠(yuǎn),所產(chǎn)生的角度失配最嚴(yán)重,誤差較高;而入射角度剛好位于幅度采樣點(diǎn)處時(shí),幾乎沒有角度失配產(chǎn)生,所以誤差很小。于是本文提出多次迭代插值,在每次迭代插值中,以上一次的插值結(jié)果為基礎(chǔ)建立新的采樣網(wǎng)格,這樣,真實(shí)來波方位距離新的采樣網(wǎng)格點(diǎn)更近,因此產(chǎn)生的失配誤差就越小,下面分別對(duì)矩形陣與任意二維陣的迭代插值方式進(jìn)行介紹。


經(jīng)過多次迭代插值,即可以一步步地降低角度失配誤差,使估計(jì)結(jié)果逼近真實(shí)值。下面通過仿真對(duì)提出的方法進(jìn)行驗(yàn)證,仿真條件與第4.1節(jié)所使用的仿真條件相同,圖9給出了非對(duì)稱插值方法經(jīng)過一次迭代插值后得到的誤差曲線,*代表迭代插值結(jié)果。

圖9 單次非對(duì)稱插值誤差估計(jì)曲線與一次迭代非對(duì)稱插值誤差估計(jì)曲線對(duì)比Fig.9 Comparison of the error estimation curve of single asymmetric interpolation with that of one-iteration asymmetric interpolation
由仿真可以看出,通過以第一次插值結(jié)果為基準(zhǔn),重新進(jìn)行幅度采樣并插值得到的估計(jì)誤差有了明顯的降低。實(shí)際上,對(duì)于矩形陣而言,只需要進(jìn)行一次迭代插值就能夠得到理想的估計(jì)結(jié)果,而對(duì)于任意二維陣,其角度失配更加嚴(yán)重,需要經(jīng)過多次迭代插值才能夠得到理想的估計(jì)結(jié)果。

從圖10可以發(fā)現(xiàn),常規(guī)波束形成的方位譜僅在高信噪比的情況下與波束圖有較高的一致性,而隨著信噪比的逐漸降低,兩者之間的一致性也出現(xiàn)了明顯的下降。在第4.2節(jié)進(jìn)行非對(duì)稱幅度插值的過程當(dāng)中,為了對(duì)方程組進(jìn)行求解,文獻(xiàn)[17]引入了波束圖的-3 dB點(diǎn)作為求解條件,本文則是定義了波束圖的左右波束斜率比作為求解條件。而幅度插值實(shí)際上是針對(duì)方位譜進(jìn)行插值,在低信噪比下,由于真實(shí)方位譜會(huì)相對(duì)波束圖出現(xiàn)展寬,此時(shí)使用基于波束圖的額外求解條件進(jìn)行求解會(huì)導(dǎo)致最終的插值結(jié)果出現(xiàn)誤差。盡管如此,由于無法事先知道方位譜,所以只能使用基于波束圖的求解條件。這里想要說明的是,本文所使用的求解條件波束圖斜率比對(duì)信噪比的變化并不敏感,其相對(duì)于-3 dB波束寬度更加穩(wěn)定。

圖10 不同信噪比下二維方位譜切片(φ=0°)Fig.10 Section of two-dimensional azimuth amplitude spectrum at different signal to noise ratios (φ=0°)
圖11給出了不同來波方位下,方位譜的-3 dB寬度隨信噪比的變化,圖12給出了不同來波方位下,方位譜的斜率比隨信混比的變化曲線。

圖11 不同DOA下方位譜的-3 dB寬度隨信噪比的變化曲線Fig.11 Variation curve of -3 dB width of azimuth spectrum with signal to noise ratio under different wave directions

圖12 方位幅度譜的斜率比隨信噪比的變化曲線Fig.12 Curve of slope ratio of the azimuth spectrum varies with the signal to noise ratio
從圖11可以明顯看出,方位譜的-3 dB寬度隨著信噪比的降低逐漸變大,甚至當(dāng)信噪比低到某個(gè)程度之后,方位譜的-3 dB寬度已經(jīng)不存在了,并且對(duì)于小孔徑的陣列來說,當(dāng)來波方位偏離0°之后,其方位譜的一側(cè)也不存在-3 dB寬度。因此,使用波束圖的-3 dB點(diǎn)作為求解條件進(jìn)行非對(duì)稱插值,其估計(jì)結(jié)果受信噪比的影響較為嚴(yán)重,雖然文獻(xiàn)[17]中通過粗估信噪比,并根據(jù)信噪比對(duì)插值過程中左右波束寬度進(jìn)行修正,但這樣仍然會(huì)存在較大的誤差。而本文定義的斜率比隨信噪比的降低變化較小,相對(duì)穩(wěn)定,這代表雖然低信噪比下方位譜相對(duì)于波束圖有所展寬,但是在這個(gè)過程中斜率比幾乎不改變,因此使用波束圖的斜率比進(jìn)行非對(duì)稱插值,能夠在低信噪比下獲得較高的插值精度。
接下來以10°作為采樣間隔,在不同的信噪比下計(jì)算了不同插值方法在與圖1中z軸之間空間夾角小于50°范圍內(nèi)的平均估計(jì)誤差,經(jīng)過100次蒙特卡羅計(jì)算得到的誤差曲線如圖13所示。


圖13 不同插值方法估計(jì)誤差隨信噪比的變化曲線Fig.13 Estimation error curves of different interpolation methods varies with the signal to noise ratio
在圖13(a)中,線性插值由于其與真實(shí)方位譜擬合較差,所以具有較高的估計(jì)誤差;而二次函數(shù)插值與高斯插值方法的性能相近,在各個(gè)信噪比下均低于線性插值的誤差,且在考慮真實(shí)方位譜的非對(duì)稱性后其誤差有明顯的下降。同時(shí),對(duì)比圖13(b)可以發(fā)現(xiàn),在經(jīng)過一次迭代插值之后,對(duì)稱線性插值誤差存在略微下降,而非對(duì)稱線性插值誤差下降更加明顯;對(duì)稱高斯插值與對(duì)稱二次插值誤差反而有所升高,這是因?yàn)閳D8(a)中這兩種方式的誤差曲線的極小值點(diǎn)并不是在采樣點(diǎn)附近,而是在兩個(gè)采樣點(diǎn)之間;而非對(duì)稱二次插值與非對(duì)稱高斯插值在低信噪比下誤差相對(duì)單次插值改變不明顯,但是在較高信噪比下(10 dB以上),出現(xiàn)了下降,取得了上述所有方法中最優(yōu)秀的估計(jì)效果。
為了進(jìn)一步對(duì)本文提出的非對(duì)稱幅度插值方位估計(jì)算法的有效性與估計(jì)精度進(jìn)行驗(yàn)證,下面將對(duì)實(shí)際的湖試采集數(shù)據(jù)進(jìn)行處理。數(shù)據(jù)于杭州千島湖采集,千島湖是一個(gè)人造湖,水底存在一定起伏,水面較為平靜。實(shí)驗(yàn)區(qū)域湖深約80 m,無水底山體干擾,實(shí)驗(yàn)時(shí)接收陣列與目標(biāo)均布放在水深20 m附近。其中,接收陣列為均勻半波長(zhǎng)布陣四元矩形陣,陣列朝向斜向下,由于一些其他原因,四元矩形陣的右下角陣元沒有采集到有效數(shù)據(jù),因此實(shí)際上使用的陣列為三元的L形陣列。目標(biāo)為一人工聲源,發(fā)射信號(hào)為10 kHz的連續(xù)波(continuous wave, CW)信號(hào)。
對(duì)接收信號(hào)按照0.01°的間隔在整個(gè)空間進(jìn)行密集波束形成掃描,得到二維方位譜如圖14所示,可以發(fā)現(xiàn)峰值出現(xiàn)在 (-17.85°,43.51°)的位置,該處即為密集掃描估計(jì)得到的來波方位。而在二維方位譜的左下角也存在峰值,需要說明的是,根據(jù)圖2坐標(biāo)系的定義,該峰值所對(duì)應(yīng)角度|θ|+|φ|>90°,實(shí)際上這個(gè)空間角并不存在,因此將其排除。以二維方位譜最大值所對(duì)應(yīng)的空間角作為基準(zhǔn),采樣間隔為10°,在表1中給出了不同的幅度插值方法的估計(jì)結(jié)果進(jìn)行對(duì)比。表1中,密集波束掃描的掃描間隔為0.01°,*代表進(jìn)行了4次迭代插值。

表1 不同插值方法的估計(jì)結(jié)果Table 1 Estimated results of different interpolation methods (°)

圖14 密集掃描方位幅度譜Fig.14 Azimuth amplitude spectrum with intensive scanning
從表1的插值結(jié)果可以看出,在只進(jìn)行單次插值的結(jié)果中,非對(duì)稱插值的估計(jì)結(jié)果雖然略微優(yōu)于對(duì)稱插值,但是整體的估計(jì)結(jié)果與密集波束掃描得到的結(jié)果相差較大。這是因?yàn)椴杉瘮?shù)據(jù)使用的陣列為三元L形陣列而非矩形陣,而對(duì)于任意陣來說,在幅度采樣最大值處將其二維幅度插值分解為兩個(gè)方向上獨(dú)立的一維幅度插值會(huì)因?yàn)榻嵌仁洚a(chǎn)生誤差,這一點(diǎn)已經(jīng)在第3節(jié)中進(jìn)行了詳細(xì)論述。為了降低由角度失配造成的誤差,使用了在第4.2節(jié)中提出的多次插值的方法,經(jīng)過4次迭代插值后,可以發(fā)現(xiàn)各種插值方法的估計(jì)精度有了較為明顯的提升,其中非對(duì)稱二次插值與非對(duì)稱高斯插值的估計(jì)精度提升最為明顯,其估計(jì)結(jié)果已經(jīng)相當(dāng)?shù)慕咏芗ㄊ鴴呙璧慕Y(jié)果。而對(duì)稱插值法在經(jīng)過多次插值后,其估計(jì)結(jié)果精度依然較差,事實(shí)上繼續(xù)增加迭代插值次數(shù),對(duì)稱插值算法的估計(jì)精度也不會(huì)有較好的提升,這是由于這類插值方法沒有充分的考慮方位譜的左右不對(duì)稱性所導(dǎo)致。另外需要說明的是,雖然多次插值會(huì)給計(jì)算量帶來一定的提升,但是每次迭代插值只需要多計(jì)算5次波束輸出,而這所需的計(jì)算量是相當(dāng)小的,遠(yuǎn)小于以0.01°間隔進(jìn)行空間掃描的計(jì)算量。
本文提出了一種任意陣列的二維非對(duì)稱幅度插值方位估計(jì)算法。算法將二維幅度插值解耦為兩個(gè)獨(dú)立的一維幅度插值,然后利用左右波束斜率比作為求解條件,通過非對(duì)稱函數(shù)對(duì)少量離散波束輸出進(jìn)行幅度插值得到來波角度。進(jìn)一步的,針對(duì)插值過程中的角度失配誤差提出迭代插值的方法進(jìn)行修正。通過仿真,本文提出的算法相對(duì)其他插值方法估計(jì)誤差更小,湖試數(shù)據(jù)的處理結(jié)果說明,本文方法能夠以較小的計(jì)算量取得與密集掃描常規(guī)波束形成類似的估計(jì)精度。本文算法具有如下優(yōu)點(diǎn):一是采用插值的方法進(jìn)行方位估計(jì),計(jì)算量很小;二是對(duì)陣型沒有特殊要求,且能夠在陣列孔徑較小的情況下取得精確的估計(jì)結(jié)果。因此,本文算法可以應(yīng)用于一些對(duì)實(shí)時(shí)性有要求的小型航行器上,實(shí)現(xiàn)對(duì)單目標(biāo)的高精度實(shí)時(shí)方位估計(jì)。但由于本文算法是基于常規(guī)波束形成進(jìn)行插值,目前只局限于對(duì)單目標(biāo)的估計(jì),如何將算法擴(kuò)展至多目標(biāo),是下一步的研究目標(biāo)。