鄭恩明， 陳新華， 李嶷

(中國科學(xué)院聲學(xué)研究所, 北京 100190)

0 引言

在線列陣接收數(shù)據(jù)中，真實(shí)目標(biāo)信號往往會被強(qiáng)干擾所掩蓋，使得目標(biāo)探測與跟蹤非常困難。研究如何阻塞強(qiáng)干擾，提高對弱目標(biāo)的探測能力顯得尤為重要，尤其是當(dāng)存在寬帶干擾時(shí)[1-9]。

阻塞矩陣方法因運(yùn)算量小，常被工程應(yīng)用[7-9]，該方法最早被用于解決波束形成中協(xié)方差矩陣估計(jì)時(shí)的信號混入問題，然后被用于干擾阻塞方面研究[10-13]，被稱為干擾阻塞方法(本文在時(shí)域上對干擾進(jìn)行處理，稱其為“時(shí)域干擾阻塞方法”)。干擾阻塞方法是通過阻塞矩陣實(shí)現(xiàn)對預(yù)知方位的干擾抑制，會對線列陣接收數(shù)據(jù)造成一定自由度損失，改變了線列陣接收信號原始形態(tài)，在部分區(qū)間降低了對目標(biāo)的探測性能[14]。為了進(jìn)一步闡述干擾阻塞方法對目標(biāo)信號探測性能的影響，許稼等[15]和高陽等[16-17]利用輸出信干噪比增益、波束圖對該方法抗干擾性能進(jìn)行衡量，分析了干擾阻塞方法存在的問題及適用性，但并未給出改善方法。為了改善干擾阻塞方法對干擾鄰近弱目標(biāo)的探測性能，葛士斌等[18-19]從協(xié)方差矩陣上進(jìn)行分析論證，通過加大干擾阻塞所用陣元間距，減小了干擾附近的衰減區(qū)間，提高了對干擾鄰近弱目標(biāo)的探測性能。但是該方法并未解決干擾阻塞方法形成的“寬凹”或“凸”字型波束對目標(biāo)探測結(jié)果的影響，且通過加大干擾阻塞所用陣元間距也會造成類似柵瓣現(xiàn)象的存在，影響對弱真目標(biāo)的探測。

針對時(shí)域干擾阻塞方法導(dǎo)致的波束畸變和探測盲區(qū)問題，本文依據(jù)干擾阻塞方法形成的“寬凹”或“凸”字型波束與干擾形式無關(guān)，只與處理數(shù)據(jù)頻帶、干擾角度有關(guān)(在線列陣確知情況下)[20]，提出一種基于引導(dǎo)信號修正的時(shí)域干擾阻塞方法。該方法采用由引導(dǎo)信號形成的衰減曲線對時(shí)域干擾阻塞方法進(jìn)行非線性化處理，降低了時(shí)域干擾阻塞方法對線列陣波束形成輸出信號能量造成的影響，在時(shí)域上對干擾阻塞方法輸出波束實(shí)現(xiàn)了校正，降低了“寬凹”或“凸”字型波束對目標(biāo)探測結(jié)果的影響。通過數(shù)值仿真和海試數(shù)據(jù)處理進(jìn)行進(jìn)一步驗(yàn)證，本文方法通過引導(dǎo)信號有效地減小了時(shí)域干擾阻塞方法對波束形成輸出信號能量造成的衰減區(qū)間，對原衰減區(qū)間內(nèi)的弱目標(biāo)實(shí)現(xiàn)了有效探測。

1 時(shí)域干擾阻塞方法

以一個(gè)干擾和一個(gè)目標(biāo)情況為例，對其進(jìn)行分析論述，其他情況可由該情況進(jìn)行遍歷。在t時(shí)刻，相對等間距線列陣(接收陣)首端方向，有弱目標(biāo)s(t)從方位角θ1入射，干擾角i(t)從方位角θ0入射，則線列陣第k個(gè)陣元拾取的數(shù)據(jù)xk(t)可表示為

(1)

式中：1≤k≤K為線列陣陣元數(shù)；d為陣元間距；c為聲速；nk(t)為第k個(gè)陣元拾取的加性高斯白噪聲。

第k+1個(gè)陣元在t時(shí)刻拾取的數(shù)據(jù)xk+1(t)為

(2)

將第k+1個(gè)陣元拾取數(shù)據(jù)相對第k個(gè)陣元拾取數(shù)據(jù)延遲時(shí)間τk+1,k與第k個(gè)陣元拾取數(shù)據(jù)相減可得新數(shù)據(jù)yk(t)為

(3)

當(dāng)τk+1,k=dcosθ0/c時(shí)，由(3)式得到的新數(shù)據(jù)可表示為

(4)

將yk(t)按頻域形式可表示為

(5)

式中：ω=2πf為數(shù)據(jù)角頻率，f為頻率；Nk(ω)為第k個(gè)陣元拾取的背景噪聲功率譜。

由(5)式可知，yk(t)中的目標(biāo)信號相對于真實(shí)目標(biāo)信號s(t)也發(fā)生了失真。下面將進(jìn)一步討論s(t)中目標(biāo)信號的變化情況。

2 目標(biāo)信號能量變化分析

2.1 理論分析

由理論分析可知，干擾阻塞后，目標(biāo)信號相對于真實(shí)目標(biāo)信號s(ω)發(fā)生了(1-ejωd(cos θ1-cos θ0)/c)的變化。為了更直觀地分析時(shí)域干擾阻塞方法對波束形成輸出造成的影響，接下來從線列陣波束形成輸出信號和噪聲變化上進(jìn)行分析論述，以便得到改善時(shí)域干擾阻塞方法。

令線列陣各陣元拾取干擾、信號和背景噪聲，背景噪聲和背景噪聲之間相關(guān)性為0，干擾阻塞前，波束形成在掃描角度為θ，θ∈[0°,180°]上輸出信號、噪聲能量分別為

(6)

式中：τk=(k-1)dcosθ/c為波束形成第k個(gè)陣元相對參考陣元(本文選第1個(gè)陣元為參考陣元)進(jìn)行的延時(shí)，T為一次處理數(shù)據(jù)采樣長度。

為了更好說明，將(6)式按頻域形式表示為

(7)

式中：ωl為數(shù)據(jù)處理所用濾波器下限；ωh為數(shù)據(jù)處理所用濾波器上限；(·)*為共軛函數(shù)；M為快速傅里葉變換長度；同一背景下，可令N(ω)=Nk(ω)為各陣元拾取的背景噪聲功率譜。

同理，干擾阻塞后，波束形成在掃描角度θ上輸出信號、噪聲能量分別為

(8)

同樣，將(8)式按頻域形式表示為

(9)

對比(7)式和(9)式可知，干擾阻塞前后，波束形成輸出信號能量變化與γ2(ω)有關(guān)，噪聲能量變化與β2(ω)有關(guān)。

當(dāng)掃描角度θ=θ1時(shí)，(7)式和(9)式可進(jìn)一步表示為

(10)

(11)

對比(10)式和(11)式可知，干擾阻塞前后，波束形成在目標(biāo)方向θ1上輸出信號、噪聲能量變化可表示為

(12)

由(12)式所示函數(shù)可知，時(shí)域干擾阻塞方法對波束形成輸出信號能量、噪聲能量造成的變化主要與γ2(ω)有關(guān)，其變化趨勢一致。對于確知線列陣，γ2(ω)主要與數(shù)據(jù)處理頻帶[ωl,ωh]、目標(biāo)方位角θ1和干擾方位角θ0有關(guān)。

2.2 數(shù)值仿真分析

為了進(jìn)一步驗(yàn)證時(shí)域干擾阻塞方法對波束形成輸出信號能量、噪聲能量造成的變化主要與γ2(ω)有關(guān)，且其變化趨勢一致，進(jìn)行如下仿真分析。仿真分析中陣元數(shù)K=64，陣元間距d=4 m，聲速c=1 500 m/s，目標(biāo)信號為寬帶線性調(diào)頻信號，背景噪聲為高斯白噪聲。接下來分析處理數(shù)據(jù)頻帶[ωl,ωh]、目標(biāo)方位角θ1和干擾方位角θ0對波束形成輸出信號能量、噪聲能量造成的影響及其變化趨勢。

1)數(shù)據(jù)處理頻帶為[100 Hz,200 Hz]。不同干擾方位角θ0情況下，波束形成輸出信號能量、噪聲能量在不同方位角θ1上的仿真結(jié)果如圖1所示。

圖1 輸出信號和噪聲量衰減區(qū)間(干擾方位變化)Fig.1 Attenuation intervals of output signal and noise energy (interference bearing change)

2)干擾方位角θ0=30°. 不同處理頻帶情況下，波束形成輸出信號能量、噪聲能量在不同方位角θ1上的仿真結(jié)果如圖2所示。

圖2 輸出信號和噪聲能量衰減區(qū)間(處理頻帶變化)Fig.2 Attenuation intervals of output signal and noise energy (processing frequency band change)

由圖1和圖2可知，時(shí)域干擾阻塞方法會對波束形成輸出信號能量和噪聲能量造成不同程度的衰減，輸出信號能量和噪聲能量在一部分方位區(qū)間被增強(qiáng)，而在另一部分方位區(qū)間被降低。由于輸出信號能量的較大差異，將影響衰減區(qū)間對弱目標(biāo)的探測，降低了弱目標(biāo)探測性能。

3 基于引導(dǎo)信號修正的時(shí)域干擾阻塞方法

3.1 理論分析

由于實(shí)際數(shù)據(jù)處理中目標(biāo)信號形式和方位角θ1未知，無法通過模擬處理頻帶目標(biāo)信號作為引導(dǎo)信號修正干擾阻塞對波束形成輸出信號能量的影響。但是由第2節(jié)理論分析和數(shù)值仿真結(jié)果可知，干擾阻塞前后線列陣波束形成輸出信號能量、噪聲能量變化趨勢一致。對此，可通過模擬處理頻帶背景噪聲作為引導(dǎo)信號修正干擾阻塞對波束形成輸出信號能量的影響，在干擾角度處，通過干擾阻塞方法對噪聲引導(dǎo)信號進(jìn)行處理，得到干擾阻塞后輸出能量變化趨勢，然后依據(jù)輸出能量變化趨勢修正干擾阻塞對實(shí)際數(shù)據(jù)處理的影響。

為了能夠真實(shí)模擬干擾阻塞對線列陣波束形成輸出噪聲能量的影響規(guī)律，本文采用事先錄取的海洋背景噪聲作為引導(dǎo)信號(也可通過模擬高斯白噪聲作為引導(dǎo)信號)。按照第2節(jié)所示干擾阻塞方法，可得到阻塞θ0方位數(shù)據(jù)后的新數(shù)據(jù)為

(13)

干擾阻塞后，波束形成在掃描角度θ上輸出波束為

(14)

Kalman濾波是一組以遞推關(guān)系給出隨機(jī)系統(tǒng)狀態(tài)的最優(yōu)線性濾波算法，根據(jù)觀測數(shù)據(jù)自動改變系數(shù)并保證下一時(shí)刻真實(shí)值均方誤差最小[21]，本文接下來采用Kalman濾波對Py′,n(θ)進(jìn)行光滑處理。

由于本文中所需濾波數(shù)據(jù)為一維非相關(guān)數(shù)據(jù)，其狀態(tài)轉(zhuǎn)移和觀測方程模型可表示為

(15)

式中：S(θ)=Py′,n(θ)；v(θ)和w(θ)分別為狀態(tài)噪聲和觀測噪聲，其方差分別為Q和R.

Kalman濾波方程可表示為

(16)

在(16)式Kalman濾波方程中，由于干擾阻塞引導(dǎo)信號后輸出波束Py′,n(θ)存在一定波動，且波動帶寬較寬，求取觀測噪聲方差R時(shí)，所用數(shù)據(jù)為[0°,180°]內(nèi)所有數(shù)據(jù)；在求取狀態(tài)噪聲方差Q時(shí)，所用數(shù)據(jù)為該濾波點(diǎn)前后L個(gè)點(diǎn)數(shù)據(jù)(L根據(jù)實(shí)際干擾阻塞后波束“寬凹”或“凸”變化程度設(shè)置)；預(yù)測方差P(θ)初始值所用數(shù)據(jù)為干擾阻塞引導(dǎo)信號后輸出波束Py′,n(θ)前L個(gè)點(diǎn)數(shù)據(jù)；狀態(tài)初始更新值Y(θ)為干擾阻塞引導(dǎo)信號后輸出波束Py′,n(θ)第1個(gè)點(diǎn)值。

利用Y(θ)作為衰減曲線，按(17)式對(8)式輸出能量進(jìn)行非線性化處理，可使干擾阻塞前后，線列陣波束形成輸出信號能量比值近似為1，即經(jīng)非線性化處理可有效降低時(shí)域干擾阻塞方法對線列陣輸出信號能量造成的影響。

(17)

3.2 實(shí)現(xiàn)流程及運(yùn)算量分析

由第3.1節(jié)理論分析可知，基于引導(dǎo)信號修正的時(shí)域干擾阻塞方法可分為如下步驟實(shí)現(xiàn)：

步驟2根據(jù)預(yù)設(shè)干擾阻塞角度θ0，按(14)式得到引導(dǎo)信號時(shí)域干擾阻塞后波束形成在掃描角度θ∈[0°,180°]上輸出波束Py′,n(θ)；

步驟3利用第3.1節(jié)所述Kalman濾波方程和參數(shù)設(shè)置方法，對步驟2輸出波束進(jìn)行光滑處理，得到波束衰減曲線Y(θ)；

步驟4根據(jù)處理頻帶[ωl,ωh]對線列陣采集數(shù)據(jù)進(jìn)行濾波處理，得到處理頻帶內(nèi)數(shù)據(jù)xk(t)；

步驟5根據(jù)預(yù)設(shè)干擾阻塞角度θ0，按(8)式得到線列陣時(shí)域干擾阻塞后波束形成在掃描角度θ∈[0°,180°]上輸出波束Py(θ)=Py,s(θ)+Py,n(θ)；

步驟7在處理頻帶[ωl,ωh]和干擾阻塞角度θ0未發(fā)生變化情況下，更新線列陣采集數(shù)據(jù)，重復(fù)步驟4～步驟6得到最新數(shù)據(jù)經(jīng)干擾阻塞后波束形成輸出結(jié)果；在處理頻帶[ωl,ωh]發(fā)生變化情況下，更新線列陣采集數(shù)據(jù)，重復(fù)步驟1～步驟6得到最新數(shù)據(jù)經(jīng)干擾阻塞后波束形成輸出結(jié)果；在處理頻帶[ωl,ωh]未發(fā)生變化，干擾角度θ0發(fā)生變化情況下，更新線列陣采集數(shù)據(jù)，重復(fù)步驟2～步驟6得到最新數(shù)據(jù)經(jīng)干擾阻塞后波束形成輸出結(jié)果。

由以上實(shí)現(xiàn)過程可知，在處理頻帶[ωl,ωh]和干擾阻塞角度θ0未發(fā)生變化情況下，本文方法相比時(shí)域干擾阻塞方法，只需在統(tǒng)計(jì)時(shí)間內(nèi)對引導(dǎo)信號進(jìn)行一次時(shí)域?yàn)V波、波束形成和Kalman濾波就行，在統(tǒng)計(jì)時(shí)間內(nèi)增加運(yùn)算量有限，可忽略不記；在處理頻帶[ωl,ωh]發(fā)生變化情況下，本文方法相比時(shí)域干擾阻塞方法，在統(tǒng)計(jì)時(shí)間內(nèi)每次更新數(shù)據(jù)需要對引導(dǎo)信號進(jìn)行一次時(shí)域?yàn)V波、波束形成和Kalman濾波，在統(tǒng)計(jì)時(shí)間內(nèi)增加近1倍運(yùn)算量；在處理頻帶[ωl,ωh]未發(fā)生變化，干擾角度θ0發(fā)生變化情況下，本文方法相比時(shí)域干擾阻塞方法，在統(tǒng)計(jì)時(shí)間內(nèi)只需對引導(dǎo)信號進(jìn)行一次時(shí)域?yàn)V波，每次更新干擾阻塞角度時(shí)，需要對濾波后引導(dǎo)信號進(jìn)行波束形成和Kalman濾波，在統(tǒng)計(jì)時(shí)間內(nèi)增加運(yùn)算量1倍波束形成運(yùn)算量和181個(gè)數(shù)據(jù)的Kalman濾波。

3.3 數(shù)值仿真分析

為了進(jìn)一步驗(yàn)證本文方法可以有效地降低時(shí)域干擾阻塞方法對線列陣波束形成輸出信號造成的影響，進(jìn)行如下數(shù)值仿真分析。

令干擾信號頻帶、目標(biāo)信號頻帶均為[100 Hz,200 Hz]，干擾與目標(biāo)方位角分別為θ0=60°和θ1=63°，背景噪聲為帶寬[100 Hz,200 Hz]的高斯白噪聲，聲速c=1 500 m/s，陣元間距d=4 m，陣元數(shù)K=64. 目標(biāo)與干擾譜級比為-20 dB，目標(biāo)與背景噪聲譜級比為-15 dB. 圖3～圖6為不同方法在處理頻帶[100 Hz,200 Hz]內(nèi)對θ0=60°干擾阻塞前后所得結(jié)果，在時(shí)域?yàn)V波處理中濾波器階數(shù)為256，Kalman濾波方程中L=10，其他相關(guān)參數(shù)根據(jù)(16)式描述所得。

圖3 方位歷程圖(干擾阻塞前)Fig.3 Bearing-time history (before interference blocking)

圖4 方位歷程圖(時(shí)域干擾阻塞方法)Fig.4 Bearing-time history (interference blocking in time domain)

圖6 波束圖(t=1 s)Fig.6 Beam pattern (t=1 s)

由圖3～圖6可知：在阻塞干擾時(shí)，由于時(shí)域干擾阻塞方法會對線列陣波束形成輸出信號能量造成不同程度的衰減，干擾阻塞后弱目標(biāo)不能被有效探測；而本文方法采用了由引導(dǎo)信號形成的衰減曲線對時(shí)域干擾阻塞方法進(jìn)行了非線性化處理，降低了時(shí)域干擾阻塞方法對線列陣波束形成輸出信號能量造成的影響，在阻塞干擾時(shí)，原衰減區(qū)間內(nèi)的弱目標(biāo)可被有效探測。數(shù)值仿真結(jié)果驗(yàn)證了本文方法通過引導(dǎo)信號可有效地減小干擾阻塞衰減區(qū)間，實(shí)現(xiàn)波束校正，對原衰減區(qū)間內(nèi)弱目標(biāo)實(shí)現(xiàn)有效探測。

為了進(jìn)一步驗(yàn)證2種方法檢測不同位置處弱目標(biāo)情況，進(jìn)行如下數(shù)值仿真。仿真中令干擾、目標(biāo)信號頻帶均為[100 Hz,200 Hz]，干擾方位角θ0=60°，目標(biāo)方位角θ1為0°～180°，聲速c=1 500 m/s，陣元間距d=4 m，陣元數(shù)K=64. 目標(biāo)與干擾譜級比為-20 dB，目標(biāo)與背景噪聲譜級比為-15 dB，在時(shí)域?yàn)V波處理中濾波器階數(shù)為256，Kalman濾波方程中L=10，其他相關(guān)參數(shù)根據(jù)(16)式描述所得。圖7為θ1∈[1°,180°]時(shí)，2種方法進(jìn)行200次獨(dú)立統(tǒng)計(jì)所得檢測弱目標(biāo)成功概率。

圖7 檢測弱目標(biāo)成功概率Fig.7 Success probability of weak target detection

由圖7可知，在本文仿真條件下，相比時(shí)域干擾阻塞方法，本文方法未成功檢測目標(biāo)(成功概率小于50%)對應(yīng)方位區(qū)間寬度由40°縮小到5°.

3.4 海試數(shù)據(jù)處理

試驗(yàn)數(shù)據(jù)為2016年在南海進(jìn)行目標(biāo)探測試驗(yàn)所得。試驗(yàn)采用64元水平線陣接收信號，陣元間隔為4 m，水平線陣尾端方向設(shè)為180°.

數(shù)據(jù)處理實(shí)例1：處理頻帶為[40 Hz,50 Hz]，在時(shí)域?yàn)V波處理中濾波器階數(shù)為256，Kalman濾波方程中L=10，其他相關(guān)參數(shù)根據(jù)(16)式描述所得。圖8～圖10分別為3種方法所得方位歷程圖，干擾阻塞角度有2個(gè)，分別為圖8中20°～40°能量最大值對應(yīng)角度θ0和θ0=160°.

圖8 方位歷程圖(干擾阻塞前)Fig.8 Bearing-time history (before interference blocking)

圖9 方位歷程圖(時(shí)域干擾阻塞方法)Fig.9 Bearing-time history (interference blocking in time domain)

圖10 方位歷程圖(本文方法)Fig.10 Bearing-time history (the proposed method)

由圖8～圖10可知：在阻塞θ0為20°～40°和θ0=160°干擾時(shí)，位于衰減區(qū)域內(nèi)的θ1為100°～140°處目標(biāo)不能被時(shí)域干擾阻塞方法實(shí)現(xiàn)有效探測；而本文方法利用了由引導(dǎo)信號形成的衰減曲線對時(shí)域干擾阻塞方法進(jìn)行了非線性化處理，降低了時(shí)域干擾阻塞方法對線列陣波束形成輸出信號能量造成的影響，位于原衰減區(qū)域內(nèi)的θ1為100°～140°處目標(biāo)被有效探測。試驗(yàn)數(shù)據(jù)1處理結(jié)果驗(yàn)證了本文方法通過引導(dǎo)信號減小了時(shí)域干擾阻塞衰減區(qū)間，實(shí)現(xiàn)了波束校正，對原衰減區(qū)域內(nèi)的弱目標(biāo)實(shí)現(xiàn)有效探測。

數(shù)據(jù)處理實(shí)例2：處理頻帶為[80 Hz,250 Hz]，在時(shí)域?yàn)V波處理中濾波器階數(shù)為256，Kalman濾波方程中L=10，其他相關(guān)參數(shù)根據(jù)(16)式描述所得。圖11～圖13分別為3種方法所得方位歷程圖，干擾阻塞角度同樣有2個(gè)，分別為圖11中20°～40°能量最大值對應(yīng)角度θ0和θ0=160°.

圖11 方位歷程圖(干擾阻塞前)Fig.11 Bearing-time history (before interference block)

圖12 方位歷程圖(時(shí)域干擾阻塞方法)Fig.12 Bearing-time history (interference block in time domain)

圖13 方位歷程圖(本文方法)Fig.13 Bearing-time history (the proposed method)

由圖11～圖13可知，在阻塞θ0為20°～40°和θ0=160°附近干擾時(shí)，相比時(shí)域干擾阻塞方法，本文方法可對120°～140°處運(yùn)動弱目標(biāo)實(shí)現(xiàn)有效探測。試驗(yàn)數(shù)據(jù)2處理結(jié)果進(jìn)一步驗(yàn)證了本文方法通過引導(dǎo)信號減小了時(shí)域干擾阻塞衰減區(qū)間，實(shí)現(xiàn)了波束校正，對原衰減區(qū)間內(nèi)的弱目標(biāo)實(shí)現(xiàn)有效探測。

4 結(jié)論

1) 本文論述了時(shí)域干擾阻塞方法實(shí)現(xiàn)原理，分析了時(shí)域干擾阻塞前后線列陣波束形成輸出信號與噪聲能量變化函數(shù)，并通過數(shù)值仿真驗(yàn)證了波束形成輸出信號與噪聲能量變化趨勢一致性。根據(jù)理論和數(shù)值仿真分析結(jié)果，提出了基于引導(dǎo)信號修正的時(shí)域干擾阻塞方法，修正方法通過模擬噪聲能量變化趨勢，減小了時(shí)域干擾阻塞方法對線列陣波束形成輸出信號能量造成的衰減區(qū)間，實(shí)現(xiàn)了波束校正。

2) 在仿真條件下，相比時(shí)域干擾阻塞方法，本文方法未成功檢測目標(biāo)(成功概率小于50%)對應(yīng)方位區(qū)間寬度由40°縮小到5°. 數(shù)值仿真和海試數(shù)據(jù)處理結(jié)果進(jìn)一步驗(yàn)證，本文方法利用引導(dǎo)信號形成的衰減曲線對輸出波束進(jìn)行非線性化處理，有效地縮小了線列陣波束形成輸出信號能量差異，降低了原時(shí)域干擾阻塞方法對衰減方位區(qū)間內(nèi)目標(biāo)探測性能影響，提高了時(shí)域干擾阻塞方法抗干擾性能和探測性能，為進(jìn)一步提高時(shí)域干擾阻塞方法的性能提供了一種可行方案。