基于稀疏恢復的自適應角度多普勒補償方法

郭藝奪，宮 健

(空軍工程大學防空反導學院，陜西西安 710051)

0 引言

空時自適應處理(Space Time Adaptive Processing, STAP)[1-2]方法能夠抑制地面雜波，檢測慢速目標，其雜波抑制性能取決于由訓練樣本單元估計的雜波協方差的準確度。在實際中，待測單元的雜波協方差是由與待測單元獨立同分布(Independent Identically Distributed, IID)的訓練樣本單元估計得到的，但在機載非正側陣雷達系統中，雜波多普勒頻率和空間頻率具有嚴重的距離相關性，尤其是近程條件下，距離相關性更加顯著，使得訓練樣本單元與待測單元波分布不滿足IID條件[3-4]，進而無法以訓練樣本單元來準確估計待測單元的雜波協方差矩陣，使得STAP方法的雜波抑制性能下降，無法有效檢測慢速目標。

補償距離相關性的方法有很多，包括多普勒彎曲(Doppler Warping, DW)[5]、角度多普勒補償(Angle Doppler Compensation, ADC)[6]、空時內插補償(Space Time Interpolating Technique, STINT)[7]和基于配準補償(Registration Based Compensation, RBC)[8]等方法。自適應角度多普勒補償(Adaptive Angle Doppler Compensation, A2DC)方法[9]及其改進方法[10]是一種補償非正側陣雷達雜波距離相關性的常用方法。該方法利用雷達雜波回波數據通過子孔徑平滑估計雜波主特征向量的最小方差無失真響應(Minimum Variance Distortionless Response, MVDR)譜，基于此估計訓練樣本單元和待測單元的雜波譜中心位置，然后將訓練樣本單元的雜波譜中心分別沿角度和多普勒方向平移，使得平移后訓練樣本單元的MVDR譜中心與待檢測單元的MVDR譜中心重合。由于該方法利用雷達雜波回波數據自適應計算各訓練樣本單元的轉換矩陣，因此可有效降低了載機平臺提供的飛行配置參數誤差對補償性能造成的影響。

然而，A2DC方法估計的雜波分布特性是基于子孔徑平滑獲得的，孔徑損失使得空時譜分辨率有一定程度的損失，系統自由度(Degree of Freedom, DOF)下降，影響后續補償過程性能，導致變換后的各訓練樣本單元平穩性下降，雜波抑制性能受限。此外，A2DC方法在計算各距離單元雜波的主特征向量的MVDR譜中心位置時，通常先要利用特征分解求得各距離單元雜波協方差矩陣的主特征向量，而典型的矩陣特征分解方法如特征值分解和奇異值分解算法的運算量都非常大，這使得A2DC方法實現起來較為復雜，限制了其應用。文獻[10]利用投影逼近子空間追蹤(Projection Approximation Subspace Tracking，PAST)方法解決A2DC方法運算量大的問題，但PAST方法并不是全局收斂的，可能陷入局部極小點而無法得到全局最優解。

針對上述問題，本文提出了基于稀疏恢復(Sparse Recovery, SR)[11]的A2DC方法——SR-A2DC方法。該方法首先對不同距離單元的回波數據進行稀疏恢復，得到雜波空時分布譜，并計算雜波協方差矩陣，已避免常規A2DC方法孔徑損失帶來的影響；然后利用正交近似投影子空間追蹤(Orthogonal PAST, OPAST)[12]方法估計不同距離單元的雜波主特征向量，降低運算量；最后估計主特征向量的空間頻率和多普勒頻率，計算空時轉換矩陣，將訓練樣本單元回波數據與待測單元的回波數據進行配準。理論分析和仿真實驗表明，該方法可以實現雜波距離相關性的自適應補償，提高主瓣雜波的補償效果，減少運算量，在存在誤差時也能獲得較好的處理性能。

1 雜波模型及A2DC方法

1.1 雜波模型

機載雷達第l個距離單元的雜波數據由該距離單元上的多個離散雜波塊疊加而成：

(1)

式中，P為離散雜波塊個數，σi為第i個雜波塊對應的復幅度，wt和ws分別為雜波的歸一化多普勒頻率和空間頻率，Si為對應的空時二維導向矢量：

Si(wt,ws)=St,i(wt,i)?Ss,i(ws,i)

(2)

式中，St,i和Ss,i分別為時域導向矢量和空域導向矢量：

(3)

式中，wt,i=2πfd,i/fprf，ws,i=2πdfs,i/λ，N和K分別為天線陣元數和相干脈沖數，fd,i和fs,i分別為雜波對應的多普勒頻率和空間頻率，d和λ分別為陣元間距和波長。

非正側陣雷達雜波的多普勒頻率和空間頻率與俯仰角和方位角的耦合關系為

(4)

fs,i=cosφicosθi

(5)

式中，θi和φi分別為雜波對應的方位角和俯仰角，v0為載機速度，θp為天線陣面與載機速度之間的夾角。

1.2 A2DC方法

A2DC方法首先以空域和時域的孔徑損失為代價估計各距離單元的雜波協方差矩陣；接著對估計得到的雜波協方差矩陣進行特征分解，得到主特征向量，然后計算各個距離單元主特征向量的MVDR譜中心位置；最后以待測單元為基準，通過角度和多普勒平移，使得各訓練樣本單元主特征向量的MVDR譜中心與待測單元主特征向量的MVDR譜中心重合。其具體實現步驟如下。

第一步，利用子孔徑平滑法獲得該距離單元的雜波協方差矩陣。

子孔徑平滑法首先將第l個距離單元的雷達回波數據Xl∈CNK×1變換成矩陣Xl∈CN×K：

(6)

設空域子孔徑和時域子孔徑分別為P和Q，對Xl進行空時平滑可以得到D=(N-P+1)·(K-Q+1)個子矩陣Xu,v∈CP×Q為

(7)

式中：u=1,2,…,N+P-1；v=1,2,…,K-Q+1。

由子孔徑平滑得到的雜波回波數據子矩陣，可估計第l個距離單元回波的協方差矩陣Rsl∈CPQ×PQ為

(8)

其中，Vec(Xu,v)表示將矩陣Xu,v的第2列置于第1列的下面，第3列置于第2列的下面，以此類推，使矩陣Xu,v變為一列向量。D的值應不小于雜波協方差矩陣Rl的列數，即D≥PQ，以確保Rl為滿秩矩陣。

第二步，對雜波協方差矩陣進行特征值分解，得到相應的主特征向量。

在得到各距離單元的雜波協方差矩陣估計后，A2DC方法利用特征值分解求得雜波協方差矩陣Rl的最大特征值λl,max及其對應的最大特征向量hl,max。

第三步，根據主特征向量的相位信息估計得到主特征向量的MVDR譜中心的多普勒頻率和空間頻率。

主特征向量的MVDR譜中心的位置可以通過二維峰值搜索計算得到，也可以直接通過主特征向量的相位信息估計得到。后者需要估計雜波協方差矩陣Rl的主特征向量的空間相位斜率κl,s和時間相位斜率κl,t，利用它們與空間頻率和多普勒頻率的線性關系計算主特征向量的MVDR譜中心的多普勒頻率fl,d和空間頻率fl,s。

(9)

第四步，計算訓練單元與待測單元主特征向量MVDR譜中心的多普勒頻率和空間頻率之間的差值為

(10)

接著，根據訓練樣本單元和待測單元的MVDR譜中心的多普勒頻率和空間頻率，計算各個訓練樣本單元的空時轉換矩陣TA2DC,l，對訓練樣本單元接收回波數據進行補償，其中，第l個訓練單元的轉換矩陣為

TA2DC,l=Ttl?Tsl

(11)

式中，

(12)

(13)

經過A2DC補償后的第l個距離單元的雷達回波數據為

(14)

2 SR-A2DC方法

由以上分析可知，A2DC方法可以自適應實現雜波距離相關性補償，但該方法需要對雜波協方差矩陣Rl進行特征分解，使得計算空時轉換矩陣TA2DC,l的運算量巨大，且由于孔徑損失，Rl的估計準確度不高，以上兩點使得A2DC方法的應用受限。

2.1 基于稀疏的雜波協方差矩陣求解

由式(1)可知，機載雷達回波數據是由不同空間頻率和多普勒頻率的回波數據疊加而成的，將空間頻率和多普勒頻率分別遍歷并離散為Ns=ρsN，Nd=ρdN個分辨單元，則第l個距離單元的雜波回波數據可以表示為

(15)

式中，ρs和ρd分別表示空間頻率和多普勒頻率的離散化程度，在高分辨情況下一般遠大于1；Sq(wt,q,ws,q)為第q個空時導向矢量，γq為其對應的復幅度，wt,q和ws,q為對應的多普勒頻率和空間頻率；αl=[αl,1,αl,2,…,αl,NsNd]代表雷達回波數據在空間頻率-多普勒頻率域上的幅度分布，即雜波空時譜；Ψ為超完備基矩陣：

(16)

估計雜波空時譜等同于在方程(15)中已知Xk和Ψ而求解αl。由于Ψ的列數NsNd遠大于行數NK，因此方程(15)屬于欠定方程，存在多個可能解。根據稀疏恢復理論，當雜波空時譜αl具有稀疏性時，欠定方程可以較高的概率求解。

針對上述問題的求解方法主要有3種，即凸優化方法[13]、FOCUSS算法[14]和MP算法[15]，其中，FOCUSS算法不利用可能存在誤差的先驗知識，采用后驗知識迭代加權，使解的能量逐步集中，可以避免誤差帶來的不利影響，且FOCUSS算法在初值設置合理時通常可獲取較好的譜估計性能，同時運算量相比范數最小化方法大大降低。因此綜合考慮算法性能和運算量，本文選擇FOCUSS算法進行計算。

實際中，僅利用待測單元恢復得到的雜波空時譜與真實雜波空時譜相比存在一定差異，一方面，如果待測單元中存在目標，則后續的處理會將目標同時抑制；另一方面，稀疏恢復估計的雜波空時譜存在過于稀疏的問題。為了更準確地表示雜波空時譜分布，一般采用對臨近訓練樣本單元稀疏恢復進行平均的方法。需要特別強調的是，訓練樣本單元必須選擇與待測單元相鄰，以保證雜波近似獨立同分布。

對式(15)進行求解，得到雜波的空時譜分布αl，則可獲得對應的雜波協方差矩陣估計：

Sq(wt,q,ws,q)H

(17)

式中，αq為雜波空時譜估計在第q個位置上的雜波復幅度；Sq(wt,q,ws,q)為第q個空時導向矢量。

2.2 引入OPAST算法估計雜波主特征向量

A2DC方法的第二步是對雜波協方差矩陣進行分解得到主特征向量，但是傳統的特征分解方法如特征值分解和奇異值分解等，雖然性能優越，但運算量大，特別是對于大線陣而言，不利于實時實現。快速子空間跟蹤方法是一種新的特征值分解方法，其中PAST和PASTd算法收斂速度快，誤差小，但該方法不是全局收斂的，而正交PAST(OPAST)算法相比PAST算法，在運算量增加不多的情況下，可以保證全局收斂。本文利用OPAST算法對雜波協方差矩陣進行特征值分解，獲得主特征向量，以降低A2DC方法的運算量。

OPAST算法將特征子空間的確定轉化為一種無約束的最優化問題進行求解，定義的無約束代價函數為

(18)

式中，Rl為第l個距離單元的雜波協方差矩陣估計，矩陣變元Wl∈NK×r(r≤NK)，Tr(·)表示矩陣的跡。

在利用遞推最小二乘算法求解代價函數的全局最小值時，為保證全局收斂和矩陣變元Wl的規范正交性，在每次迭代中對矩陣變元Wl進行正交化，即引入正交化公式：

(19)

當代價函數J(Wl)取得最小值時，Wl收斂于雜波子空間的一組正交基上，特別地，當r=1，即Wl∈NK×1時，最小化J(Wl)得到Wl，即為Rl的歸一化最主要特征向量。

與常規A2DC方法相比，利用OPAST算法求解雜波協方差矩陣的主特征向量，運算量從O(NK3)減少為4NKS+O(S2)，其中，S為雜波協方差矩陣Rl的秩。

3 仿真實驗

選取斜側陣(θp=60°)機載相控陣雷達進行仿真實驗，其中實驗條件為：發射接收陣元數均為12；相干脈沖間隔內的脈沖數為12；陣元間隔0.16 m，載機高度8 km，載機速度120 m/s，波長0.32 m，脈沖重復頻率3 000 Hz；雜噪比60 dB，雷達最大作用距離800 km；雷達主波束方位角為90°；待測單元145，假設對應的距離為25 km，且不存在距離模糊。

采樣協方差求逆方法(LSMI)的訓練樣本數為2NK；稀疏恢復中，訓練樣本數為4；空間頻率和多普勒頻率離散化程度均為8，即ρs=ρd=8；子空間平滑的空域子孔徑和時域子孔徑為P=Q=6。

實驗1：雜波空時譜估計

本實驗對比真實雜波功率譜、未經補償的LSMI法的雜波功率譜、常規A2DC方法子孔徑平滑法的雜波空時譜和本文SR-A2DC方法估計的雜波空時譜，如圖1～圖4所示。

圖1 真實雜波功率譜

圖2 未經補償的LSMI法估計的雜波功率譜

圖3 A2DC方法子孔徑雜波空時譜

圖4 SR-A2DC方法全孔徑雜波空時譜

由圖1和圖2可以看出，在非正側陣的情況下，由于雜波具有距離相關性，使用LSMI方法直接利用訓練樣本單元估計待測單元的雜波協方差矩陣，不同距離單元的雜波疊加在一塊，會導致雜波譜嚴重展寬。由圖3和圖4可以看出，常規A2DC方法得到的雜波空時譜，由于子孔徑平滑造成的孔徑損失，存在能量擴散和部分雜波分布難以估計的問題；SR-A2DC方法利用FOCUSS方法估計雜波空時譜，不存在孔徑損失，因此準確度和分辨率更高。

實驗2：雜波譜中心多普勒頻率和空間頻率

本實驗對比常規A2DC方法特征值分解(Eigen Decomposition, ED)和本文SR-A2DC方法OPAST算法估計的雜波譜中心多普勒頻率和空間頻率與真實值之間的差異，如圖5、圖6所示。

圖5 歸一化多普勒頻率對比

圖6 歸一化空間頻率對比

從圖5和圖6可以看出，相比常規A2DC方法特征值分解估計的雜波譜中心多普勒頻率和空間頻率，本文SR-A2DC方法OPAST算法準確性更高。這是因為，OPAST算法能夠以較高的收斂速度保證較小的誤差，而基于稀疏恢復的雜波空時譜估計不存在孔徑損失，兩者結合性能優于常規A2DC方法。

實驗3：補償效果對比

本實驗對比常規A2DC方法和SR-A2DC方法對雜波距離相關性進行補償的補償效果，并以改善因子為基準衡量不同方法的雜波抑制性能，如圖7～圖9所示。其中，改善因子的定義為輸出信雜噪比(Signal to Clutter-plus-Noise Ratio, SCNR)與輸入SCNR的比值。

由圖7和圖8可知，A2DC方法和SR-A2DC方法均能對非正側陣雜波距離相關性進行補償，經過補償后，主瓣雜波基本配準，但由于A2DC方法是基于子孔徑平滑，存在孔徑損失，補償性能在一定程度上受限，而SR-A2DC方法是基于全孔徑的，不存在孔徑損失，可獲得更高的分辨率，設計的空時轉換矩陣更加有效，因此補償性能更優。

圖7 A2DC方法雜波距離相關性補償效果

圖8 SR-A2DC方法雜波距離相關性補償效果

由圖9可知，經過本文SR-A2DC方法進行雜波距離相關性補償后，利用LSMI方法估計雜波協方差矩陣，計算空時濾波器權值，可以在主雜波處形成深凹口，且相比A2DC方法更窄，對慢速目標的檢測能力可提高5～10 dB。而利用A2DC方法設計的空時濾波器，則存在著凹口偏移、展寬等問題。

圖9 雜波抑制性能

4 結束語

本文詳細分析了常規A2DC方法的原理和步驟，針對該方法存在孔徑損失和運算量大的問題，提出了基于稀疏恢復的SR-A2DC方法，并將OPAST算法引入到雜波協方差矩陣特征分解中。理論分析和仿真實驗表明：1)所提的SR-A2DC方法可以避免空時損失帶來的影響，減少運算量；2)SR-A2DC方法與常規A2DC方法，雜波距離相關性補償性能更優，更適于工程實現；3)由于本文使用稀疏恢復估計雜波空時譜，而稀疏恢復求解所構造的超完備矩陣存在著基失配(off-grid)等問題，因此，如何減少或者消除off-grid的影響，自適應設計超完備基矩陣，獲得更加準確的雜波空時譜分布，是下一步研究工作的重點。