張 錚，徐雅南，孔 強

(杭州應用聲學研究所，浙江杭州 310023)

0 引言

陣列信號處理作為現(xiàn)代信號處理的重要組成部分，被廣泛應用于聲吶、雷達中。水平線列陣作為一種基本陣型，是聲吶系統(tǒng)中的一個重要組成部分。在傳統(tǒng)的聲吶系統(tǒng)中，通過常規(guī)波束形成(Conventional BeamForming,CBF)處理來獲得陣增益，提高目標探測能力。但是，CBF主瓣波束寬度寬、旁瓣級高，降低了其空間分辨能力，在多目標時容易產(chǎn)生混疊，對有強干擾的弱目標檢測能力較弱。為了提高空間分辨能力，一批高分辨波束形成算法被提出，比如最小方差無畸變響應法(Minimum Variance Distortionless Response,MVDR)[1]、最大熵法[2]、主模式抑制法(Dominant Mode Rejection,DMR)[3]、多信號分類法(Multiple Signal Classification,MUSIC)[4]、基于旋轉(zhuǎn)不變技術的信號參數(shù)估計法(Estimating Signal Parameter Via Rotational Invariance Techniques,ESPRIT)[5]等。但是，這些方法在提高了分辨力的同時，也有其使用限制。如對陣型精度要求高、需要較多數(shù)據(jù)構建協(xié)方差矩陣、計算量大、穩(wěn)定性差[6]等。

為了獲得兼顧高分辨能力和穩(wěn)定性的波束形成算法，人們開始將其他學科的成果應用到陣列信號處理中來。反卷積算法起初多應用于圖像處理領域[7-8]，用于提高圖像分辨率，Yang將反卷積算法與常規(guī)波束形成結合，提出了一種穩(wěn)健的高分辨波束形成算法，并通過水平線列陣[9]和圓陣[10]對反卷積波束形成的穩(wěn)定性、分辨力和處理增益進行了仿真研究和實驗驗證，結果表明，該方法具有和CBF一樣的穩(wěn)定性，與MVDR等高分辨波束形成算法分辨能力相近，并且提供了一定的處理增益。

盡管文獻[9]的反卷積波束形成算法具有穩(wěn)健、分辨力高和提高處理增益的優(yōu)點，但在工程應用中還會面臨一些實際問題。以水平線列陣聲吶為例，通常工作在一個較寬的頻帶內(nèi)，而反卷積算法中用到的陣列波束圖也可稱作陣列的點擴散函數(shù)(Point Spread Function,PSF)與頻率有關，較寬的頻帶即意味著需要計算更多頻點的PSF，增加了計算量；文獻[9]中使用基于Richardson-Lucy(R-L)迭代算法的反卷積波束形成算法，因為信號截斷，在波束圖的兩端會有邊界模糊問題。反卷積波束形成可以看成是在CBF的基礎上進行反卷積處理的算法，因此其計算量高于CBF。如何降低該方法的計算量，將計算量增加部分控制在一個可接受的范圍內(nèi)，是工程實踐中需要解決的實際問題。本文對以上工程應用中的問題進行分析、討論和優(yōu)化，提出了一種利于工程應用的快速反卷積波束形成方法。

1 利于工程應用的快速反卷積波束形成方法

1.1 均勻線陣常規(guī)波束形成

以一個具有N個陣元的水平線列陣為例，陣元間距為d，聲源在遠場，聲信號以平面波形式到達各陣元，信號的入射方向為θ，如圖1所示。入射信號以R=AiPi+Ns表示，其中i為目標計數(shù)，Ai為第i個目標信號強度，為第i個目標在各陣元的相移向量，其中p定義如下pin=exp[- jk(n-1)dcosθi] ，n=1，2，...，N，符號 T 表示轉(zhuǎn)置，j表示虛數(shù)，k為波數(shù)，Ns為噪聲。

圖1 均勻線陣布陣示意圖Fig.1 Schematic diagram Uniform line array

根據(jù)常規(guī)波束形成理論，定義導向向量S=[s1s2...sN]T，其中，有K個互不相關的隨機信號，則波束形成結果為

式中：H為共軛轉(zhuǎn)置(具體推導參見文獻[9])。由式(1)可以看出，B關于cosθ具有移不變性，即B(cosθ-cosθi)=B(cosθ|cosθi)，|B(cosθ|cosθi)|即為均勻線列陣在θi角度的指向性函數(shù)[11]?？蓪⑹?1)改寫成積分形式：

1.2 反卷積波束形成

由式(1)、(3)可以看出，Bp是一個只與陣元個數(shù)、陣元間距和信號波長有關的函數(shù)(即均勻線列陣的PSF函數(shù))，即在陣型確定的情況下，是可以提前確定的[12]；S是包含了目標方位信息的δ函數(shù)的和，這方位信息即所要求解的目標信息[13]；P作為CBF輸出的功率譜，可以由陣接收到的數(shù)據(jù)獲得。因此，由Bp、S和P三者之間的卷積關系可以推出，對P求反卷積即可獲得目標的方位信息[13]。

目前，有多種算法可以進行反卷積，比如CLEAN[14]，反卷積波束形成法(Deconvolution Approach for the Mapping of Acoustic Sources Method,DAMAS)[15]，非負最小二乘法(Non-Negative Least squares,NNLS)[16]， Richardson-Lucy(R-L)迭 代 算法[17-18]以及相應的衍生算法[19]。由文獻[16]可知，對于均勻線列陣這種PSF具有移不變特性的陣列，RL法相比于DAMAS、NNLS等算法具有更低的計算量，更易于工程應用。

R-L算法是一個基于貝葉斯準則的最優(yōu)化迭代計算過程，將CBF功率譜結果作為初始值，每次迭代得到的功率譜估計值和CBF功率譜結果差值越小，則認為迭代得到的目標函數(shù)越接近真實的聲源分布，取得極小值后即可停止迭代。也可直接指定迭代次數(shù)，在完成指定次數(shù)迭代后停止迭代[20]。具體的迭代公式為

1.3 寬帶隨機信號處理

在實際的工程應用中，待分析的信號多為寬帶隨機信號，而由式(4)可知，反卷積算法中用到的PSF是一個與信號頻率相關的參數(shù)，即：

分析式(7)可知，PSF在一個較窄的頻帶內(nèi)區(qū)別不大，而且頻率越高，具有近似PSF的頻帶越寬。如圖2所示。將寬帶信號劃分成若干個不等寬的較窄頻帶，每個頻帶內(nèi)取一個頻點的PSF做近似反卷積計算，這樣可以極大地降低PSF和反卷積的計算量。

圖2 各頻率的點擴散函數(shù)Fig.2 PSF at each frequency

1.4 邊界模糊問題的解決

人們在用R-L算法對二維圖像進行反卷積處理時發(fā)現(xiàn)[9]，因為圖像數(shù)據(jù)截止問題，恢復圖像邊緣部分會出現(xiàn)模糊。同理，對一維的CBF功率譜進行R-L迭代，也會在功率譜邊界位置出現(xiàn)模糊現(xiàn)象，對于本文中的均勻線列陣來說，出現(xiàn)在cosθ=±1附近。文獻[9]利用波束指向性圖的周期性，給出了一種將功率譜積分范圍擴展的方法，即將式(3)改寫為

這樣就將邊界問題移動到cosθ=±1.5附近，在1＜|cosθ|≤1.5范圍內(nèi)會得到一個復數(shù)角度，僅參與計算，在具有實際物理意義的|cosθ|≤1范圍內(nèi)不受邊界模糊影響。但是，這種邊界擴展的方式在解決了R-L算法邊界模糊問題的同時，增加了波束形成50%的計算量，尤其當波束個數(shù)較多時，這種計算量的增加極大地限制了算法在工程上的應用。

考慮到|cosθ|＞1范圍內(nèi)的波束形成結果，獲得的θ值沒有實際物理意義，在此花費較多時間計算有些得不償失，為此，本文提出一種近似的優(yōu)化方法，在不增加波束形成計算量的前提下，解決反卷積R-L算法帶來的邊界模糊問題。該方法具體步驟如下：

(1)按式(3)進行M個波束的常規(guī)波束形成計算，獲得CBF的功率譜PM；

(2)以PM的左右端點為起點向兩邊各擴展M/4個波束，擴展波束以PM左右端點值為基準取值，獲得擴展后的功率譜PkM，k=1.5；

(3)按-1.5≤cosθ≤1.5計算kM個波束的PSF函數(shù)；

(4)由PkM和PSF函數(shù)做R-L反卷積計算，取中間M個波束作為反卷積波束形成的輸出結果。

2 仿真分析

通過計算仿真分析，首先驗證了反卷積波束形成高分辨力和提高處理增益的能力。以64元的均勻線列陣為例，按信號頻率為6 kHz進行半波長布陣，6 kHz窄帶信號在不同信噪比和角度的仿真結果如圖3所示。

圖3中橫軸為以cosθ值表示的角度，縱軸為以dB為單位表示的歸一化后波束能量值。圖3(a)為-90°方向的單目標，信噪比(Signal to Noise Ratio,SNR)為50 dB，圖3(c)為圖3(a)中虛線框局部放大，可以清楚地看出，在這種“理想”條件下，反卷積波束形成3 dB主瓣寬度約是CBF主瓣寬度的1/5。圖3(b)為31°和33°方位的等強度雙目標，SNR為-10 dB，3(d)為其虛線框局部放大，可以看出在這種仿真條件下，CBF已無目標分辨能力，而反卷積波束形成依然可以清楚分辨雙目標。對比圖中噪聲部分，反卷積波束形成噪聲能量比CBF低了約10 dB。

圖3 窄帶信號反卷積結果Fig.3 Deconvolution results of narrowband signal

以上仿真的單頻窄帶信號在實際工程應用中比較少見，對于常見的寬帶信號，可以看作是多個窄帶信號的疊加。以64元、6 kHz半波長布陣的均勻直線陣為例，工作帶寬B=6 kHz，由圖2可知，各頻點PSF不同，即按傳統(tǒng)的反卷積波束形成算法，需要對帶寬內(nèi)的所有頻點都做反卷積運算，這無疑增加了很多計算量。為了減少不同頻點的反卷積運算，利于工程應用，本文提出將寬帶拆分的方法，將整個工作帶寬B拆分成(0～200 Hz)、(200～500 Hz)、 (500～1 000 Hz)、 (1000～3 000 Hz)和(3 000～6 000 Hz)共5個相對較窄的頻帶，并取每個頻帶中心頻點的PSF做反卷積波束形成。因各頻帶仿真結果具有相似性，本文僅給出1 000～3 000 Hz頻帶隨機信號的仿真結果。

圖4中橫軸為以cosθ值表示的角度，縱軸為以dB為單位表示的歸一化后波束能量值。圖4(a)中信號SNR為0 dB，圖4(b)中信號SNR為-20 dB。對比兩圖可以看出，因為利用了不同頻點的PSF信息，全頻點反卷積在性能上優(yōu)于窄帶中心頻點反卷積，但在信號SNR較低時，兩種方法檢測能力和分辨力近似，都優(yōu)于CBF。對弱目標檢測而言，窄帶中心頻點反卷積具有計算量更小的優(yōu)勢。

圖4 全頻點和窄帶中心頻點反卷積結果Fig.4 Deconvolution results at total frequency points and at the central frequency point of each narrow band

對于R-L迭代算法的邊界模糊問題，以64元6 kHz半波長布陣的均勻直線陣進行仿真，SNR為0 dB，目標來波方向分別為90°和10°，仿真結果如圖5(a)、5(b)所示。

圖5 中橫軸為以cosθ表示的角度，縱軸為歸一化后波束能量值。圖5(c)為圖5(a)中虛線框局部放大。圖5(d)為圖5(d)中虛線框局部放大，由圖5可以看出，未作周期擴展的反卷積結果在cosθ=±1附近會出現(xiàn)波形畸變，這會給檢測帶來虛假目標；當目標出現(xiàn)在邊界位置時，會造成目標方位估計誤差。與之相對的，進行了周期擴展的反卷積相比于CBF提高了檢測能力和分辨力，并且結果準確。

圖5 R-L迭代算法的邊界模糊問題驗證Fig.5 Validation of boundary blur problem of R-L iterative algorithm

由以上仿真可知，利用R-L迭代算法的反卷積波束形成需要進行周期擴展。但是，對波束進行擴展會極大地增加CBF的計算量，為此本文提出一種優(yōu)化方法，并進行如下仿真：64元6 kHz半波長布陣的均勻直線陣，(1 000～3 000 Hz)寬帶隨機信號，對CBF的能量譜做邊界擴展和反卷積，與對波束進行擴展的反卷積結果進行對比，結果如圖6所示。

圖6(a)中目標在60°，SNR為0 dB可以看出在較高信噪比條件下，兩種反卷積方法檢測能力一致，都優(yōu)于CBF；改變信號SNR為-20 dB，仿真結果如圖6(b)所示，對比圖6(a)、6(b)可以看出，在不同的SNR條件下，當目標不在邊界位置，兩種反卷積方法對目標檢測能力和分辨力一致，在邊界處的噪聲有一點差別，但對目標檢測來說，基本沒有影響；圖6(c)為目標在10°，SNR為0 dB的仿真結果，圖中虛線框局部放大如圖6(e)所示。因為目標在邊界位置，對CBF能量譜進行邊界擴展時會將目標信號的能量擴展，在做反卷積時，這部分擴展的信號能量會泄露到噪聲中，增加了噪聲的譜級。對波束擴展，相當于增加了虛擬的角度范圍，對cosθ=±1附近的目標不會有影響；改變圖6(c)中的信號SNR為-20 dB，仿真結果如圖6(d)所示，圖中虛線框局部放大如圖6(f)所示，與圖6(c)、6(e)對比可知，在較低的SNR條件下，目標分辨能力有所降低。由圖6(a)～6(e)可知，對波束擴展的反卷積算法不受目標位置影響，能一直保持較好的檢測能力和分辨力。對CBF能量譜擴展的反卷積算法的性能會受目標位置影響，當目標不在邊界位置時，其具有與波束形成時周期擴展的反卷積算法一樣的檢測能力和分辨力；當目標在邊界位置時，算法性能會有一定下降，但依然優(yōu)于CBF。另外，考慮實際應用情況，均勻直線陣在cosθ=±1位置附近為檢測“盲區(qū)”，一般不作為基陣的主要觀測方向，以降低“盲區(qū)”的部分檢測性能換取計算量的大幅下降，CBF能量譜邊界擴展的反卷積算法在工程應用上依然有其優(yōu)勢。

圖6 波束擴展反卷積和譜擴展反卷積的處理結果對比Fig.6 Comparison of processing results between beam expansion deconvolution and spectrum expansion deconvolution

綜合以上兩種計算速度的提升方法，對一寬帶隨機信號，做129個波束全頻點的CBF，周期擴展為193個波束，每12個快拍進行一次能量積分，共1 200個快拍數(shù)據(jù)，反卷積迭代10次，對本文所提的快速方法與傳統(tǒng)反卷積波束形成算法進行了計算速度的對比，記錄仿真總時間(PSF計算采用預存方式，不計入總時間)，結果如表1所示(因仿真平臺性能差異，此計算時間僅作對比參考)。

由表1中的第2行可以看出，在單頻點情況下，兩種方法反卷積計算時間一樣，計算速度差距主要體現(xiàn)在波束個數(shù)上，因傳統(tǒng)方法要在波束域做1.5倍擴展，整體計算時間也是本文方法的約1.5倍。傳統(tǒng)方法反卷積迭代所耗時間與分析頻點個數(shù)成正比，本文方法反卷積迭代所耗時間與拆分窄帶個數(shù)成正比，隨著分析帶寬的增加，頻點個數(shù)會遠大于拆分窄帶個數(shù)，傳統(tǒng)方法所耗時間的增加也會遠大于本文方法所增加的時間。表1所示仿真結果與理論分析一致，即相比于傳統(tǒng)方法，本文方法整體計算速度在單頻點時提升了約50%。這也是本文方法速度提升的下限，隨著分析頻點的增多，本文方法的速度提升也隨之增大。

表1 傳統(tǒng)的和本文提出的反卷積波束形成方法的計算速度對比Table.1 Speed comparison between conventional and the proposed deconvolution beamforming methods

3 試驗驗證

選取一段海上試驗數(shù)據(jù)對本方法進行了驗證，試驗時兩艘船在一均勻直線陣正橫方向交叉航行，按本文方法將分析帶寬劃分成多個較窄頻帶，分別進行CBF和反卷積處理，選取低頻和高頻兩個頻段進行對比分析，具體結果如圖7所示。

圖7中橫軸為方位，以cosθ值顯示；縱軸為實驗的時間批次，1為最早時刻。對比圖7(a)、7(b)，可以直觀地看出，本文快速反卷積算法處理后的主瓣波束寬度要明顯小于CBF，圖7(a)中兩目標交叉無法分辨時間約380批次，而圖7(b)中僅約240批次。從圖7(b)中可以看出，在cosθ(0.6～0.75)范圍內(nèi)有一較弱目標，而在圖7(a)的相同位置，此目標并不明顯，大部分“淹沒”在噪聲中，由此可以看出，本文中的快速反卷積算法對弱目標的檢測能力有了一定提升。圖7(c)、7(d)中的高頻處理結果，相比于低頻，高頻時主瓣波束寬度都有明顯降低，在圖7(c)中兩目標交叉無法分辨時間約190批次，在圖7(d)中約130批次。圖7(c)中兩個目標左側都存在較高的旁瓣，與之相比，圖7(d)中的旁瓣被削弱。

圖7 兩種方法對實際海試數(shù)據(jù)的處理結果Fig.7 Processing results of real sea trial data by the two methods

4 結論

相比于CBF方法，反卷積波束形成方法具有更高的分辨力、更高的處理增益和更好的旁瓣抑制能力，但作為CBF的一種后置處理方法，計算量有所增加。對于寬帶隨機信號的處理，因為利用了不同頻率的PSF，全頻點反卷積在高SNR條件下性能優(yōu)于本文提出的窄帶中心頻點反卷積方法；但在低SNR條件下，兩者檢測性能相近，而中心頻點反卷積方法因其計算量遠小于全頻點反卷積方法，在工程應用上更具有優(yōu)勢。對于R-L算法的邊界模糊問題，本文提出了一種在CBF能量譜上進行邊界擴展的方法，相比于在波束角度上進行擴展，既解決了R-L算法的邊界模糊問題，又極大地降低了波束形成的計算量，更利于工程實現(xiàn)。綜合兩種速度提升方法后，本文方法相比于傳統(tǒng)反卷積波束形成速度至少提高了50%。

海試數(shù)據(jù)的處理結果也驗證了本方法的有效性。結果顯示，本方法具有更高的分辨力，提高了對弱目標的檢測能力，對旁瓣級的抑制在5 dB以上。