范文濤 唐 波 李朋偉 王英志

（中國(guó)人民解放軍92578部隊(duì)，北京 100161）

1 引言

水下被動(dòng)聲檢測(cè)核心問(wèn)題是如何從復(fù)雜時(shí)變、充滿干擾的環(huán)境中檢測(cè)出感興趣的弱目標(biāo)信號(hào)［1］。常用處理手段是通過(guò)波束形成將無(wú)源傳感器陣列接收到的被動(dòng)聲信號(hào)進(jìn)行時(shí)間域和空間域處理，從而挖掘出目標(biāo)的寬帶和窄帶信號(hào)特征，實(shí)現(xiàn)后續(xù)目標(biāo)跟蹤和識(shí)別。類似于雷達(dá)和語(yǔ)音信號(hào)處理領(lǐng)域［2］，應(yīng)用于水下被動(dòng)聲檢測(cè)的波束形成方法可分為兩類：常規(guī)和自適應(yīng)。其中最常用的檢測(cè)算法是常規(guī)波束形成（Conventional Beamforming，CBF）［3］，該方法可以適配各種水下感知系統(tǒng)和海洋環(huán)境，因此也非常受到水下信號(hào)處理工程師的青睞。但是，隨著安靜型潛艇的出現(xiàn)，所輻射的噪聲級(jí)越來(lái)越低，導(dǎo)致傳統(tǒng)CBF很難檢測(cè)到目標(biāo)。一些研究者考慮使用最優(yōu)陣處理來(lái)解決這個(gè)問(wèn)題，盡量壓低波束圖的旁瓣，使得系統(tǒng)獲得更大的陣增益，但是這種方法存在計(jì)算量大的缺點(diǎn)。其中一個(gè)著名的檢測(cè)算法叫最小方差無(wú)失真響應(yīng)（Minimum Variance Distortionless Response，MVDR）波束形成器［4］，這種方法主要是基于信號(hào)和噪聲協(xié)方差矩陣來(lái)求解最優(yōu)加權(quán)函數(shù)，使得波束形成器能夠自適應(yīng)的將波束圖凹槽對(duì)準(zhǔn)干擾方向，將主瓣對(duì)準(zhǔn)目標(biāo)。MVDR 波束形成器比起CBF 方法大大增強(qiáng)了在多目標(biāo)分辨率和干擾抑制方面的性能，但是MVDR 波束形成器在實(shí)際系統(tǒng)中應(yīng)用時(shí)不夠穩(wěn)健，尤其是快拍數(shù)較少時(shí)，需要對(duì)接收信號(hào)協(xié)方差矩陣進(jìn)行對(duì)角加載，不同的加載量會(huì)造成一定的檢測(cè)增益損失，到目前為止，還未出現(xiàn)一種既能保住MVDR 檢測(cè)性能又能增強(qiáng)其穩(wěn)健度的對(duì)角加載方法。

最近，一些研究者嘗試將CBF 或者M(jìn)VDR 波束形成器與其他的非線性處理方法融合用以提高傳統(tǒng)方法在水下被動(dòng)聲檢測(cè)中的檢測(cè)性能。盲源分離（Blind Source Separation，BSS）［5］就是一種在語(yǔ)音信號(hào)處理領(lǐng)域中較著名的處理方法，該方法主要是通過(guò)從接收的陣列混合信號(hào)矩陣中提取出感興趣的信號(hào)。小波變換（wavelet transform，WT）［5］和經(jīng)驗(yàn)?zāi)B(tài)分解（Empirical Mode Decomposition，EMD）方法［6］也被考慮用以改善聽(tīng)音波束、LOFAR 譜［2］、DEMON 譜［2］的信噪比。還有二階錐處理［4］和遺傳算法［3］也被用來(lái)幫助波束形成器找到適合自己的最優(yōu)權(quán)向量。但是遺憾的是，所有這些非線性處理技術(shù)都存在兩個(gè)問(wèn)題，一是需要對(duì)環(huán)境和基陣進(jìn)行假設(shè)，二是工程化實(shí)現(xiàn)起來(lái)非常復(fù)雜。

在反潛戰(zhàn)領(lǐng)域，大孔徑線列陣通常被用來(lái)對(duì)安靜型潛艇目標(biāo)搜索和跟蹤的［2］，例如，舷側(cè)陣、拖曳陣、海底固定陣［4］、高級(jí)分布式陣［2］等等。但是這些陣都面臨一個(gè)相同的問(wèn)題，就是當(dāng)目標(biāo)運(yùn)動(dòng)到陣列端首附近時(shí)，波束形成器的檢測(cè)能力就極大地下降。一般稱基陣端首方向?yàn)槊^(qū)（Blind Zone，BZ）［4］。但是，這個(gè)問(wèn)題又不能回避，因?yàn)锽Z區(qū)域太重要了，從戰(zhàn)術(shù)上講，魚雷就喜歡從這個(gè)方向?qū)ε灤l(fā)起攻擊。一些學(xué)者提出用逆波束形成來(lái)抑制盲區(qū)的本艦自噪聲強(qiáng)干擾［4］，但是該方法是利用語(yǔ)音中譜減的思路來(lái)將所有波束中的干擾剪掉。要知道譜減的過(guò)程中，信號(hào)也會(huì)被減掉，這樣會(huì)引起檢測(cè)背景的提升。

本文提出一種新的波束輸出譜加權(quán)方法，步驟為第一步，從波束輸出中估計(jì)背景噪聲級(jí)；第二步，根據(jù)背景噪聲級(jí)從波束譜中提取感興趣的頻譜；第三步，基于ECKART 最優(yōu)檢測(cè)濾波器估計(jì)波束最優(yōu)權(quán)值；最后，權(quán)值被加到波束形成器輸出波束上完成弱目標(biāo)檢測(cè)。通過(guò)使用實(shí)際拖曳陣和舷側(cè)陣的數(shù)據(jù)進(jìn)行驗(yàn)證，從實(shí)驗(yàn)結(jié)果來(lái)看，所提方法在淺海多干擾環(huán)境下，對(duì)弱目標(biāo)能夠獲得4 dB 左右的檢測(cè)增益。

2 寬帶模型

考慮M個(gè)水聽(tīng)器組成的均勻線列陣，假設(shè)遠(yuǎn)場(chǎng)有D個(gè)寬帶聲源（D＜M）以平面波的形式到達(dá)基陣，波達(dá)方向?yàn)棣?[θ1，...，θd，...，θD]。那么第m個(gè)水聽(tīng)器所接收到的時(shí)域信號(hào)為

其中sd(t)表示第d個(gè)聲源信號(hào)，ξn(θd)表示第d個(gè)聲源相對(duì)于第n個(gè)水聽(tīng)器的傳播延遲，vn(t)表示疊加的高斯白噪聲。由于被動(dòng)聲納基陣接收的信號(hào)為寬頻帶信號(hào)，而大部分方位估計(jì)算法都是基于窄帶假設(shè)，所以處理時(shí)一般通過(guò)短時(shí)傅里葉變換（Shorttime Fourier Transform，STFT），將寬帶信號(hào)分段，并劃分為互不重疊的多個(gè)窄帶，然后再在每個(gè)子帶應(yīng)用相應(yīng)的信號(hào)處理算法［7］。那么，基陣輸出的M×1維快拍向量數(shù)據(jù)的頻域響應(yīng)可表示為

其中fj表示頻點(diǎn)（frequency bin），頻點(diǎn)數(shù)就是STFT對(duì)陣元域數(shù)據(jù)分段后，每個(gè)段做FFT 的點(diǎn)數(shù)，頻點(diǎn)數(shù)就是子帶數(shù)。s(fj)=[s1(fj)，...，sD(fj)]T為D×1維的源信號(hào)向量（(·)T表示轉(zhuǎn)置）。A(Θ，fj)=[a(θ1，fj)，...，a(θD，fj)]是M×D維的陣列流形矩陣。導(dǎo)引矢量a(θi，fj)［8］是：

3 所提出的方法

圖1 是本文所提方法的框架結(jié)構(gòu)。在圖中，陣列接收的時(shí)域信號(hào)被STFT 變換成頻域多子帶信號(hào)，以保證窄帶信號(hào)處理算法能夠工作。公式（2）中的X(fj，m，k)表示第m時(shí)間段下第j個(gè)頻點(diǎn)下的第k根譜線：

第二步是自適應(yīng)波束形成，本文中使用MVDR算法［4］。MVDR 方法的權(quán)值通過(guò)對(duì)如下公式的最小化來(lái)實(shí)現(xiàn)：

其中R(ω)是協(xié)方差矩陣［7］。在實(shí)際應(yīng)用中，一般對(duì)協(xié)方差矩陣采取時(shí)間平滑處理措施，則協(xié)方差矩陣為［7］：

公式（5）是一個(gè)趨向于單位1的約束優(yōu)化過(guò)程；

公式（5）和（6）可以通過(guò)拉格朗日乘子法來(lái)求解，所得到的結(jié)果為［7］：

頻域MVDR波束輸出為：

因此，頻率-方位矩陣（Frequency Azimuth，F(xiàn)RAZ）可表示為：

ZFRAZ(fj)是J×B的矩陣，J是代表頻點(diǎn)數(shù)和B是波束數(shù)。

利用分裂雙通窗（Two-Pass Split Window，TPSW）［8］算法來(lái)估計(jì)每個(gè)波束的噪聲級(jí)（Noise Level，NL）并且以此提取出背景噪聲譜。TPSW 是一種常用的估計(jì)噪聲譜的處理技術(shù)，計(jì)算量少，實(shí)現(xiàn)簡(jiǎn)單，其窗函數(shù)的公式為［8］：

TPSW 將沿著信號(hào)滑動(dòng)來(lái)獲得噪聲平均值。TPSW最后輸出背景噪聲譜矩陣NFRAZ(fj)：

在寬帶檢測(cè)理論中，ECKART 濾波器［8，10］被認(rèn)為是一種能夠最大化輸出信噪比的最優(yōu)檢測(cè)器。但是，ECKART 需要對(duì)信號(hào)和噪聲具有一定的先驗(yàn)信息，實(shí)際上信號(hào)和噪聲譜是未知的，需要估計(jì)信號(hào)和噪聲功率譜。實(shí)際被動(dòng)陣列信號(hào)的FRAZ 矩陣中，信號(hào)譜分量往往表現(xiàn)出一定的空間相關(guān)性，而噪聲譜分量往往空間不相關(guān)，本文使用矩形窗在FRAZ 矩陣中沿著每三個(gè)波束頻譜分量中表現(xiàn)為“峰值”的值：

公式（13）尋找到的“峰值”zpeak(fj，θb)被認(rèn)為候選信 號(hào)譜分量。反之，被定義為“谷值”：

由此，本文設(shè)定規(guī)則：當(dāng)公式（13）估計(jì)得到的“峰值”大于對(duì)應(yīng)的公式（12）背景噪聲譜分量值，ECKART 濾波器的信號(hào)分量即為該“峰值”，噪聲分量為常數(shù)1；當(dāng)“峰值”小于對(duì)應(yīng)的背景噪聲譜分量值時(shí)，ECKART 濾波器的信號(hào)分量即為該“峰值”，噪聲分量為背景噪聲譜分量的平方；當(dāng)公式（14）估計(jì)得到的“谷值”大于對(duì)應(yīng)的公式（12）背景噪聲譜分量值，即認(rèn)為該“谷值”為噪聲或者干擾，對(duì)該值疊加背景噪聲譜值后賦為負(fù)值，作為抑制；當(dāng)“谷值”小于對(duì)應(yīng)的噪聲譜分量值時(shí)，即直接對(duì)該值賦負(fù)值。因此，波束譜特征加權(quán)值為：

因此，波束形成器的輸出波束可以被加權(quán)為［11］：

最后，寬帶方位歷程譜能通過(guò)非相干累加Tnew(fj)的方式得到。所提方法處理框架如圖1所示。

4 數(shù)值試驗(yàn)

仿真采用陣列信號(hào)模擬器產(chǎn)生數(shù)據(jù)，場(chǎng)景設(shè)置為：拖曳線列陣為41元水聽(tīng)器均勻線列陣，其陣元間距滿足半波長(zhǎng)的假設(shè)；拖船輻射噪聲源級(jí)設(shè)為SL拖船=120dB，拖速v拖速=14節(jié)，拖纜長(zhǎng)度R拖纜=400米，陣深為30米；目標(biāo)A輻射噪聲源級(jí)設(shè)為SL目標(biāo)A=120dB，距離拖曳陣R目標(biāo)A=5.4海里，航速v目標(biāo)A=18節(jié)；目標(biāo)B 輻射噪聲源級(jí)設(shè)為SL目標(biāo)B=105 dB，距離拖曳陣R目標(biāo)B=2.7 海里，航速v目標(biāo)B=8節(jié)；目標(biāo)C輻射噪聲源級(jí)設(shè)為SL目標(biāo)C=115 dB，距離拖曳陣R目標(biāo)C=3.7 海里，航速v目標(biāo)C=14節(jié)；目標(biāo)A、B、C初始狀態(tài)分別位拖曳陣舷角是95°、115°和135°，向拖曳陣端射方向運(yùn)動(dòng)。其中陣列信號(hào)數(shù)據(jù)采樣率fs=48 kHz，F(xiàn)FT 長(zhǎng)度N=1024，頻率分辨率Δf=fs/N=46.875 Hz，所用的數(shù)據(jù)段總點(diǎn)數(shù)為48000 采樣點(diǎn)，將該數(shù)據(jù)段進(jìn)行分段，每段長(zhǎng)度為1024，相鄰兩段數(shù)據(jù)重疊率為50%，處理頻段為100～600 Hz，相應(yīng)的頻點(diǎn)數(shù)為93。MVDR 方法和本文所提方法均采用“白噪聲增益”對(duì)角加載，加載量取0.092。

圖2 是MVDR 與本文所提方法的寬帶警戒方位-時(shí)間-幅度瀑布顯示圖。對(duì)比（a）、（b）子圖可以看出，本文所提方法對(duì)三個(gè)目標(biāo)的檢測(cè)增益要高于MVDR 方法，尤其是過(guò)了40度接近拖曳線列陣聲納檢測(cè)盲區(qū)時(shí)，本文所提方法對(duì)目標(biāo)C依然保持接觸；本文所提方法的檢測(cè)背景起伏要大于MVDR 方法，原因在于本文所提方法在提取譜特征過(guò)程中，會(huì)存在將部分噪聲譜分量作為信號(hào)提取的問(wèn)題，導(dǎo)致檢測(cè)背景出現(xiàn)起伏。

圖3 為對(duì)弱目標(biāo)利用蒙特卡洛仿真方法，比較不同輸入信干噪比、距離和信噪比條件下，本文所提方法與MVDR 方法的陣增益性能，其中干噪比為10 dB，設(shè)定-3 dB為檢測(cè)閾值。從子圖（a）中可以看出，本文所提方法在輸入信干噪比為-34 dB 時(shí)就能檢測(cè)目標(biāo)，比MVDR 方法小了9 dB；子圖（b）可以看出，本文所提方法的檢測(cè)距離比MVDR 方法要遠(yuǎn)5 公里左右；子圖（c）可以看出，MVDR 方法陣增益信噪比門限比本文所提方法要低3 dB，意味著本文所提方法的抗噪性能不如MVDR方法。

可以得到以下結(jié)論：

（1）從方位歷程圖的檢測(cè)結(jié)果來(lái)看，本文所提方法增大了背景的起伏，但是對(duì)于弱目標(biāo)尤其是盲區(qū)附近目標(biāo)的檢測(cè)性能要好于MVDR方法。

（2）表1 為總結(jié)圖3 得出的有關(guān)兩種方法的弱目標(biāo)最小可檢測(cè)性能，包括最小輸入信干噪比、最遠(yuǎn)檢測(cè)距離和最小輸入信噪比三個(gè)度量指標(biāo)。從表中可以看出，本文所提方法比MVDR 方法的最小輸入信干噪比低了9 dB；本文所提方法比MVDR 方法的最遠(yuǎn)檢測(cè)距離多了4 公里；本文所提方法比MVDR方法的信噪比門限高了3 dB。

表1 弱目標(biāo)最小可檢測(cè)性能表Tab.1 Weak target minimum detectable performance table

5 實(shí)際海試結(jié)果

本文利用舷側(cè)陣數(shù)據(jù)驗(yàn)證算法性能。水聽(tīng)器陣為32 陣元等間隔直線陣，間距1 米。試驗(yàn)中目標(biāo)以6 節(jié)航速定深航行，有合作目標(biāo)。根據(jù)多次測(cè)量中比較典型的不同陣元數(shù)和不同頻率上的噪聲能量分布，設(shè)定檢測(cè)處理頻段為310～630 Hz，采樣頻率為6000 Hz。處理數(shù)據(jù)時(shí)間長(zhǎng)度為3000 s，該段時(shí)間內(nèi)共有6 批目標(biāo)。將該數(shù)據(jù)段進(jìn)行分段，每段長(zhǎng)度為2048，相鄰兩段數(shù)據(jù)重疊率為75%。MVDR 方法和本文所提方法均采用“白噪聲增益”對(duì)角加載，加載量取0.092，積分時(shí)間為2 s。本次試驗(yàn)中沒(méi)有合作目標(biāo)，通過(guò)目力觀測(cè)和導(dǎo)航雷達(dá)記錄，觀察到了6批次水面漁船和商船目標(biāo)，其中起始于170度舷角的是大型商船，其余為漁船。

圖4 子圖（a）是MVDR 方法的寬帶警戒方位-時(shí)間-幅度瀑布顯示二維、三維對(duì)比圖。從子圖（a）檢測(cè)結(jié)果來(lái)看，目標(biāo)A、目標(biāo)B、目標(biāo)C、目標(biāo)D、目標(biāo)E和目標(biāo)F度附近有存在目標(biāo)，其中目標(biāo)C 能量較強(qiáng)，目標(biāo)A、目標(biāo)B、目標(biāo)D 和目標(biāo)E 的前1500 s 歷程非常弱，目標(biāo)F 從舷側(cè)陣端首往另一側(cè)端首（0 度）運(yùn)動(dòng)，靠近端首盲區(qū)的歷程部分MVDR 檢測(cè)效果不佳，出現(xiàn)了方位模糊現(xiàn)象。

圖4 子圖（b）為本文所提方法得到的三維和二維方位歷程顯示，可以看出本文所提方法比MVDR方法對(duì)六批目標(biāo)的歷程跟蹤較為清晰，其中對(duì)于起始于目標(biāo)A和目標(biāo)B的兩個(gè)漁船目標(biāo)檢測(cè)增益明顯提升，全程歷程顯示清晰；對(duì)于起始于目標(biāo)D 和目標(biāo)E 的兩個(gè)漁船目標(biāo)，也較為明顯的將其前1500 s的歷程呈現(xiàn)出來(lái)；對(duì)于170 度目標(biāo)靠近兩端端首部分的歷程顯示最清晰完整。本文所提方法在100度至120度方位扇區(qū)內(nèi)的檢測(cè)背景起伏要高于MVDR方法，原因在于受到正橫目標(biāo)C 附近強(qiáng)干擾能量泄露的影響，本文所提方法在提取信號(hào)譜特征分量過(guò)程中，將泄露至該方位扇區(qū)內(nèi)的干擾噪聲分量作為信號(hào)分量提取，導(dǎo)致最后在該方位扇區(qū)內(nèi)的背景起伏大于MVDR方法。

表2為舷側(cè)陣數(shù)據(jù)陣增益的統(tǒng)計(jì)結(jié)果。由于該實(shí)錄數(shù)據(jù)無(wú)法預(yù)知準(zhǔn)確的目標(biāo)態(tài)勢(shì)、信干噪比、干噪比和信信比的信息，因此本文考慮選取整個(gè)數(shù)據(jù)分析歷程中的四秒（10 s、500 s、1000 s 和2000 s）時(shí)間，針對(duì)目標(biāo)A 附近的弱目標(biāo)方位譜峰值與平均背景的對(duì)數(shù)比值作為檢測(cè)陣增益，經(jīng)過(guò)比較本文所提方法比MVDR方法平均陣增益提高了2.785 dB。

表2 舷側(cè)陣數(shù)據(jù)陣增益統(tǒng)計(jì)表Tab.2 Flank array data array gain statistics table

由于本文所提方法，為了更加直觀地說(shuō)明本文所提方法的復(fù)雜度，將本文所提方法與MVDR 方法處理同一批試驗(yàn)數(shù)據(jù)，統(tǒng)計(jì)觀察時(shí)間2 s內(nèi)方位譜輸出所需的計(jì)算時(shí)間。計(jì)算機(jī)硬件采用Intel i7-8500CPU，主頻3.1 Hz，內(nèi)存8.0 GB。計(jì)算軟件采用MATLAB2019a版本。如表3所示，MVDR，本文所提方法所需計(jì)算時(shí)間分別為0.83 s和0.97 s。對(duì)比可以看出本文所提方法基于MVDR方法波束輸出譜增加了特征提取和ECKART 濾波兩個(gè)步驟，相比MVDR方法慢了0.14 s，計(jì)算復(fù)雜度沒(méi)有明顯提升。

表3 方法性能比較表Tab.3 Method performance comparison table

6 結(jié)論

本文提出一種波束譜特征加權(quán)的水下弱目標(biāo)檢測(cè)方法，通過(guò)對(duì)波束輸出信號(hào)進(jìn)行背景級(jí)估計(jì)，從而利用ECKART最優(yōu)濾波器從波束輸出譜中提取感興趣的目標(biāo)譜信號(hào)，實(shí)現(xiàn)對(duì)弱目標(biāo)的檢測(cè)。計(jì)算機(jī)仿真和實(shí)際海試數(shù)據(jù)證明該方法對(duì)于水下被動(dòng)弱目標(biāo)檢測(cè)具有一定的效果，同時(shí)所提方法基于MVDR 方法波束輸出譜進(jìn)行特征處理，計(jì)算復(fù)雜度未明顯提升，利于在實(shí)際水下傳感器處理系統(tǒng)上工程實(shí)現(xiàn)。