郭浩泉， 朱代柱,2

(1.上海船舶電子設備研究所， 上海 201108； 2.水聲對抗技術重點實驗室， 上海 201108)

高方位分辨力的波束形成在目標檢測、參數(shù)估計與識別中具有重要的意義，在多目標和強弱目標同時存在的情況下，有效抑制多目標相干干擾是后續(xù)信號處理的基礎。導向最小方差(steered minimum variance, STMV)波束形成算法以其高方位分辨力、無旁瓣以及良好穩(wěn)健性等優(yōu)越性能而受到廣泛關注，其中最主要的運算量取決于協(xié)方差矩陣的求逆過程，對陣元數(shù)為M的基陣進行經典STMV分析，每個導向方位的每個頻率點都需要進行O(M3)次乘法運算[1-2]，這對大規(guī)?；嚨膶拵盘柼幚硎请y以承受的，必須尋求更低復雜度的簡化算法。

為了降低STMV算法中矩陣求逆過程的復雜度，許多學者做了大量研究，基本可以歸納為降維處理、遞歸迭代、矩陣分解和矩陣分塊4大類方向：1)在降維處理方面，文獻[3-7]主要利用盲分離的自適應波束形成、對導向矢量進行處理以降低干擾的自適應波束形成、結合穩(wěn)健Capon波束形成器和Krylov子空間、基于陣形在線估計聚焦和Krylov子空間等降維技術。2)在遞歸迭代方面，文獻[8-9]主要利用基于統(tǒng)計優(yōu)化準則的迭代MVDR算法、基于矩陣求逆引理的遞歸MVDR算法進行迭代運算減少計算量，文獻[10]提出了基于一階遞歸濾波器和矩陣求逆引理的迭代STMV算法，將計算量降低至經典STMV算法的2/M。3)在矩陣分解方面，文獻[11] 主要利用基于QR分解的自適應波束形成算法將求解權向量的問題轉化為求解三角線性方程的問題，避免了自相關矩陣的估計和求逆運算；文獻[12-14]提出了子陣STMV波束形成算法，大幅度降低了運算量。4)在矩陣分塊方面，文獻[15-17]主要利用協(xié)方差矩陣為正定Hermite陣的對稱性并結合Strassen矩陣求逆的高效性降低協(xié)方差矩陣求逆的復雜度。

這些研究將全陣平均分為多個子陣的矩陣分解方法計算量最小，但由于削弱了子陣之間的信號相關性，在本質上是一種近似算法，所以與全陣STMV算法相比性能略有損失；并且，如何分配子陣數(shù)目和子陣陣元數(shù)成為一個繞不開的問題，缺乏尋求最優(yōu)解的理論依據。其他降維處理、遞歸迭代和矩陣分塊研究方面，計算量最小的是文獻[10]在保持STMV性能的同時將計算量降低至經典STMV算法的2/M。

本文針對經典STMV的大運算量問題，提出了一種低復雜度的BISTMV算法，通過對矩陣求逆公式的層層分解與分塊迭代，可以大幅度降低協(xié)方差矩陣求逆的復雜度，更加有利于工程應用。

1 導向最小方差波束形成

1.1 STMV的基本原理

基于導向協(xié)方差矩陣的導向最小方差波束形成算法的基本原理是通過波束旋轉，使指定方向的響應為1，而陣列輸出功率最小。

假設各陣元時域輸出為x(t)=[x1(t)…xM(t)]T，相應的頻域輸出為X(fk)，各陣元的權向量為W=[w1w2…wM]T，則導向矩陣為：

(1)

式中：M為陣元數(shù)；θ為方位角；d為陣元間距；c為聲速。

陣列的輸出為：

(2)

式中：h為頻率點上限；l為頻率點下限。

陣列輸出的平均功率為：

P(θ)=E[y(t,θ)y(t,θ)H]

(3)

STMV波束形成器的最佳權向量w是使指定方向的響應為1，而陣列輸出功率最小，故：

(4)

1.2 迭代STMV算法原理

STMV波束形成時，對N次快拍的頻域輸出Xn(fk)進行求和，得到平均頻域輸出為：

(5)

則導向協(xié)方差矩陣Rstmv(θ)可表示為：

(6)

(7)

(8)

式中：α為一階遞歸濾波器的系數(shù)；β=1-α。

由矩陣求逆引理：

(9)

(10)

(11)

1.3 分塊迭代STMV算法原理

1.3.1 分塊矩陣求逆算法

關于2階分塊矩陣的求逆公式，文獻[18]已詳細闡述。

(12)

(13)

把式(13)命名為“組合求逆公式”。

證明考慮到：

則

(14)

式中E-1可由式(12)求得。

(15)

1.3.2 算法原理

假設U被等分成K段，記為Ui，i=1,2，…,K，每段的長度為k=M/K。則導向協(xié)方差矩陣Rstmv(θ)被等分成K×K塊，每塊矩陣的階數(shù)為k?？紤]到Rstmv(θ)為埃爾米特矩陣，則式(7)可表示為：

(16)

(17)

2) 構造矩陣

(18)

(19)

3) 應用分塊矩陣求逆公式遞推求解E-1；

2 計算量分析

2.1 經典STMV算法計算量

假設陣元數(shù)為M，信號帶寬為B，利用N次快拍的頻域輸出估計互譜密度矩陣，注意到乘法運算的計算時間為加法運算的數(shù)十倍(約40倍)，因此本文僅考慮乘法運算的計算量。其計算量主要集中在2部分：

1)估計互譜密度矩陣時，需在每個頻點處需進行NM2次乘法運算，總共需進行NM2B次乘法運算；

2) 估計空間譜估計時，需在每個波束進行如下計算：

①估計導向協(xié)方差矩陣，在每個頻點處需進行2M2次乘法運算，總共需進行2M2B次乘法運算；

②對導向協(xié)方差矩陣求逆需進行M3次乘法運算，然后利用式(4)進行空間譜估計，需進行M2+M次乘法運算。

綜上所述，假設對L個波束進行空間譜估計，所需的計算量為：

NM2B+[M3(2B+1)+M2+M]L

(20)

2.2 迭代STMV算法計算量

假設m階矩陣求逆的計算量為C(m)，則根據式(12)、(14)、(15)，2×2、3×3、4×4分塊矩陣求逆的計算量依次為：

x2=2C(m)+4m3

(21)

x3=x2+C(m)+11m3

(22)

x4=x3+C(m)+21m3

(23)

綜上所述，考慮n×n分塊矩陣求逆的計算量的通項公式。n×n分塊矩陣求逆的計算量為：

(24)

結合式(21)～(23)，可得通項公式為：

3(n-1))m3

(25)

假設陣元數(shù)為M，信號帶寬為B。其計算量主要集中在2部分：

綜上所述，假設對L個波束進行空間譜估計，所需的計算量為：

[M2(4B+1)+M(B+1)]L

(26)

2.3 分塊迭代STMV算法計算量

假設陣元數(shù)為M，信號帶寬為B，U被等分成K段，記為Ui，i=1,2，…,K。每段的長度為k=M/K。其計算量主要集中在5部分：

4) 利用分塊矩陣求逆公式求解E-1，需進行(K-1)C(k)+((K-2)3+4(K-2)2+3(K-2))k3/2次乘法運算；

綜上所述，假設對L個波束進行空間譜估計，所需的計算量為：

(M2B+(3k2+k)B+(K(K+1)/2-1)k2B+

(K-1)C(k))L+((K-2)3+4(K-2)2+

3(K-2))k3L/2+((K-1)2k3+(K-1)k3+

M2+M)L

(27)

式中：C(k)為k階矩陣求逆的計算量，取C(k)=k3。同時考慮到k=M/K，式(27)可化簡為：

(28)

(29)

假設信號帶寬B=5 kHz，在0°～180°方位內每隔2°形成一個預成波束，即對L=91個波束進行空間譜估計，快拍數(shù)N=4，K=16，STMV算法與分塊迭代STMV算法計算量的比值隨陣元數(shù)M的變化如圖 1 所示。

圖1 2種算法計算量比值

文獻[10]提出的迭代算法的計算量大致是STMV算法計算量的2/M。為便于分析，同時考慮到L?N、B?M，式(20)和式(29)分別可近似化簡為：

(30)

當K足夠大時，本文提出算法的計算量大致是原STMV算法計算量的1/4M，約為文獻[10]迭代STMV算法計算量(2/M)的1/8。當M?5時，計算量顯著降低。

3 仿真分析

仿真條件如下：設48元均勻線陣接收，陣元間距d為信號上限頻率對應的半波長，每快拍處理時間為3 s，采樣頻率為fs，處理帶寬為0.05～0.25fs，聲速為1 500 m/s。兩目標輻射噪聲信號及背景噪聲均為相互獨立的高斯白噪聲，目標1的方位為30°保持不變，陣元輸入端信噪比為6 dB；目標2的初始方位為50°，陣元輸入端信噪比為0 dB，在0.1fs處疊加一根線譜，該線譜相對于周圍連續(xù)譜的信噪比為15 dB，對總能量的貢獻可忽略不計，目標2每快拍的目標方位增加0.4°。取α=0.5，由STMV算法和分塊迭代BISTMV算法(K=6)估計的130 s空間譜歷程圖看不出明顯差異，截取不同時刻的空間譜如圖2所示，圖3為圖2(b)的局部放大圖。

圖2 不同時刻2種算法的空間譜

圖3 線譜目標局部放大圖

4 海試數(shù)據分析

為驗證算法的穩(wěn)健性，對海上試驗數(shù)據進行處理，試驗態(tài)勢如圖4所示。以拖船的航向為0°且直線航行，其后拖著一條均勻線列陣。目標1距離拖船較近，其輻射噪聲較強，目標4的初始方位為75°且朝著拖船方向航行。試驗時，航運繁忙且目標3、4相對于目標1、2距離拖船較遠，其輻射噪聲較弱。

圖4 試驗態(tài)勢

對接收信號進行預成波束形成，處理帶寬為0.025～0.1fs，由STMV算法和BISTMV算法(K=6)估計的空間譜歷程圖也看不出明顯差異，截取第88 s和141 s時刻的空間譜，如圖5所示，當存在強目標時(目標1、3)，仍可實現(xiàn)對弱目標(目標2、4)的有效檢測。

圖5 88 s和141 s 2種算法的空間譜

由圖5(a)和5(b)可以看出，BISTMV算法在同時存在多個強弱目標的情況下，對每個強弱目標的檢測效果都與STMV算法比較貼合，對其中的線譜目標還有1～4 dB(平均2.2 dB)的檢測增益,如圖6(a)和6(b)所示。究其原因，是由于BISTMV算法在式(8)對協(xié)方差矩陣進行求逆的過程中引入了一階遞歸濾波器。一階遞歸濾波器將前后快拍之間的信息進行了關聯(lián)，由于目標輻射噪聲信號是平穩(wěn)相關的，而背景噪聲是弱相關的，因此，這種關聯(lián)在大多數(shù)情況下對目標檢測是有利的：當處理寬帶信號時幾乎沒有差異，但當處理窄帶信號時效果會有改善，尤其是當目標輻射噪聲信號中包含有窄帶線譜信息時(如圖5(a)和5(b)中75°左右的目標4)；只有當目標方位快速變化時，這種關聯(lián)才會略微展寬該目標所在波束的寬度(如圖5(a)中127°左右的目標1)，造成輕微負面影響。但總體而言，BISTMV算法在繼承STMV算法的優(yōu)良特性的同時大幅減小了計算量，對窄帶線譜信息豐富的目標檢測還可以獲得一定程度的增益提升。

圖6 線譜目標局部放大圖

5 結論

1)本文提出了基于一階遞歸濾波器、矩陣求逆引理和矩陣分塊求逆的BISTMV算法，進行了計算量分析，大致是經典STMV算法計算量的1/4M，是迭代STMV算法計算量的1/8，顯著減小了計算量，提高了其工程可實現(xiàn)性。

2)仿真和海試數(shù)據處理結果表明，該算法保持了STMV算法高方位分辨力和對微弱目標檢測的優(yōu)勢。

3)該算法對窄帶線譜信息豐富的目標檢測可以獲得平均2 dB左右的增益提升，對工程應用，特別是線譜目標的檢測具有重要的意義。