張宇軒 金禹希 陳世進 吳永清 郝程鵬

（1.中國科學院聲學研究所，北京 100190；2.中國科學院大學，北京 100049）

1 引言

當目標具有一定運動速度時，目標信號與雜波信號在角度-多普勒平面上是可分離的，而空時自適應處理（Space-Time Adaptive Processing，STAP）技術正是利用這種可分離性對雜波加以抑制，從而提升機載雷達的運動目標檢測性能［1-2］。在實際應用中，STAP 濾波器的設計需要通過訓練樣本來估計雜波協(xié)方差矩陣（Clutter Covariance Matrix，CCM）。傳統(tǒng)STAP 方法通過極大似然準則來估計CCM，然而根據(jù)Reed-Mallett-Brennan 準則，最大似然估計方法需要二倍于系統(tǒng)自由度的獨立同分布訓練數(shù)據(jù)，才能使輸出信雜噪比損失小于3 dB。由于在非均勻環(huán)境中樣本量過少，STAP 雜波抑制能力將顯著下降［3-4］。為了緩解非均勻雜波環(huán)境中訓練樣本數(shù)量不足時的STAP 性能下降問題，研究者提出了降維STAP 和知識輔助STAP（Knowledge-Aided STAP，KA-STAP）等方法來減少對訓練數(shù)據(jù)量的需求［5-8］。降維STAP 方法理論上可以將數(shù)據(jù)需求量從系統(tǒng)自由度的兩倍降至CCM 秩的兩倍，此類方法本質(zhì)上是DOF 與數(shù)據(jù)需求量的一種權衡，在減少數(shù)據(jù)需求量的同時也會損失系統(tǒng)自由度，造成最小可檢測速度增加、系統(tǒng)誤差容忍能力降低、雜波抑制性能下降等問題［6-7］。KA-STAP 技術旨在利用各種不同形式的環(huán)境知識提高CCM 的估計精度，但一般情況下環(huán)境知識獲取難度很高，大幅增加了系統(tǒng)設計成本，并且當系統(tǒng)獲取的環(huán)境知識存在誤差時會進一步導致STAP性能下降［8］。

近期國內(nèi)外學者提出了稀疏恢復STAP 方法［9-13］（Sparse Recovery STAP，SR-STAP），此類方法通過將角度-多普勒平面網(wǎng)格化構造過完備空時導向向量字典，然后利用該字典對雜波功率進行稀疏恢復從而獲得CCM，從而在訓練樣本數(shù)量受限情況下實現(xiàn)了比傳統(tǒng)STAP 算法更好的雜波抑制性能。最常見的SR-STAP算法主要為正交匹配追蹤STAP算法［10］（Orthogonal Matching Pursuit STAP，OMP-STAP）和欠定系統(tǒng) 聚焦解STAP 算法［11］（FOcal Underdeter?mined System Solution STAP，F(xiàn)OCUSS-STAP）。這類算法簡單、高效，其缺點是性能十分依賴于超參數(shù)的設置，理論上這些超參數(shù)可根據(jù)真實環(huán)境的信噪比進行調(diào)整，然而在實際應用中真實信噪比的獲取是難以實現(xiàn)的。為拓展算法的實用性，文獻［12］提出一種迭代自適應STAP算法（Iterative Adaptive Ap?proach STAP，IAA-STAP），這種算法是一種非參數(shù)化稀疏恢復算法，對稀疏度不敏感，可以避免OMPSTAP 和FOCUSS-STAP 算法中超參數(shù)難以設置的問題。上述SR-STAP 算法均只從待檢測單元（Cell Under Test，CUT）中恢復雜波能量，是一種單觀測樣本（Single Measurement Vector，SMV）方法，因為單個樣本包含的信息有限，SMV 方法在信噪比較低或者目標移速較慢時，雜波抑制效果嚴重下降。為了充分應用多個樣本所包含的信息，文獻［13］提出了一種多快拍自適應迭代STAP 算法（Mul?tiple snapshots Iterative Adaptive Approach STAP，MIAA-STAP），將SR-STAP 擴展至多觀測樣本（Multiple Measurement Vector，MMV）的應用背景中。MIAA-STAP 先從多個訓練樣本中恢復雜波譜，再根據(jù)雜波譜恢復結(jié)果設計相應的STAP 濾波器，獲得了比IAA-STAP 更好的雜波抑制性能。

另一方面，當陣列為均勻線陣時，空域協(xié)方差矩陣可具有斜對稱特性，文獻［14］和［15］均利用了該先驗知識來提高訓練數(shù)據(jù)數(shù)量不足時陣列協(xié)方差矩陣的估計精度。受到這些工作的啟發(fā)，本文注意到當機載雷達的接收陣列為等間隔均勻線陣，并且系統(tǒng)在一個相干處理間隔（Coherent Processing Interval，CPI）中脈沖重復頻率恒定時，CCM 同樣具有斜對稱特性。因此可以考慮將該特性推廣到STAP 中，并與SR-STAP 方法相結(jié)合，以聯(lián)合利用CCM 的斜對稱特性和雜波譜稀疏性兩項先驗知識，進一步提升CCM 估計精度和雜波抑制性能。該算法首先使用斜對稱變換矩陣對CUT 中的數(shù)據(jù)和訓練樣本進行預處理，從而完成了對等效訓練樣本的擴展；隨后，將擴展后的訓練樣本與協(xié)方差稀疏迭代算法［16］（Sparse Iterative Covariance-based Estima?tion，SPICE）相結(jié)合，實現(xiàn)雜波譜的稀疏恢復和CCM 的重構。仿真結(jié)果表明，提出的算法能夠在訓練樣本量不足時表現(xiàn)出穩(wěn)健的雜波譜恢復和CCM重構能力，相應的STAP 濾波器能夠有效地抑制雜波。

2 雜波譜稀疏恢復問題的建模

本節(jié)首先對機載雷達的雜波信號模型進行介紹，隨后根據(jù)該模型構建空時導向向量字典，將雷達地雜波信號和CCM進行稀疏表示。

2.1 機載雷達地雜波信號模型

考慮如圖1 所示的工作環(huán)境，機載雷達飛行速度為v，運動方向與x軸的正向夾角為θv，平臺高度為H，θ、φ和?分別表示雜波塊與平臺的方位角、俯仰角和空間角。

假設在一個相干處理間隔（Coherent Processing Interval，CPI）內(nèi)發(fā)射機發(fā)射M個脈沖，脈沖重復頻率為fr；接收陣列為等間隔均勻線陣，其中陣元個數(shù)為N，陣元間距為d=λ/2，λ表示工作波長，則雷達接收的一個距離環(huán)內(nèi)雜波信號xc∈CMN×1可以表示為［17-18］

式中，[?]T表示轉(zhuǎn)置運算。

2.2 雜波譜的稀疏恢復模型

SR-STAP 所需要的超完備字典由多個具有不同多普勒頻率和空域頻率的空時導向向量所構成，若將整個角度-多普勒平面離散為Ns×Nd個網(wǎng)格點，則所有網(wǎng)格點對應的空時導向向量就可以組成一個字典，其中Nd=ρdN和Ns=ρsM分別表示多普勒軸和空間軸上的網(wǎng)格點數(shù)量，參數(shù)ρs、ρd分別為網(wǎng)格的空域、多普勒域分辨率。利用超完備字典，可將雜波信號近似表示為［19］

SR-STAP 方法關鍵在于估計雜波空時功率譜p，空時功率譜的求解方法可以分為SMV 方法和MMV方法，其中SMV方法僅從CUT中的數(shù)據(jù)恢復γ并計算出空時功率譜p，即

式中[?]*表示共軛運算，⊙為Hadamard 積，和分別表示向量的l0范數(shù)和l2范數(shù)，參數(shù)?為噪聲的允許誤差水平。因為單個樣本包含的信息有限，這種方法在信噪比較低或目標運動速度較慢時，雜波譜恢復效果和雜波抑制效果下降。相對的，MMV 方法旨在從多個訓練樣本中恢復空時功率譜p并估計CCM。具體來說，假設有K個獨立同分布的訓練樣本y1，y2，…，yK，其表達式分別為

其中nk為第k個距離單元的噪聲向量，滿足與目標距離單元的噪聲向量獨立同分布條件。此時，MMV方法的實現(xiàn)可描述如下

3 聯(lián)合利用雜波雙先驗知識的穩(wěn)健STAP算法

由于MMV 類SR-STAP 算法的雜波抑制效果和樣本數(shù)量呈正相關，當訓練樣本不足時，CCM 的最大似然估計估計精度下降，從而導致式（10）中雜波譜的恢復性能下降。針對這一問題，本節(jié)提出一種聯(lián)合利用CCM 斜對稱特性和雜波譜稀疏性的STAP 算法以減少訓練樣本需求量，進一步提升雜波譜恢復性能。算法通過引入斜對稱變換矩陣對訓練樣本進行預處理以提升CCM 的最大似然估計準確度，隨后將處理后的樣本與SPICE算法相結(jié)合，對CCM進行重構。

3.1 構造基于CCM斜對稱特性的稀疏恢復模型

首先對斜對稱特性進行介紹，對于任意一個N維向量a=，若該向量滿足［20］

則稱該向量具有斜對稱特性。容易驗證當時域?qū)蛳蛄亢涂沼驅(qū)蛳蛄肯辔恢行娜∮谥虚g位置時，兩者均滿足斜對稱特性。

當一個N維矩陣A∈CN×N滿足［20］

則稱該矩陣具有斜對稱特性，其中Aij表示矩陣A第i+1 行、第j+1 列的元素。更一般的，具有斜對稱特性的向量和矩陣滿足a=JNa*，A=JN A*JN，其中，JN為N維變換矩陣，表達式為：

根據(jù)斜對稱特性的定義和性質(zhì)，可以得到如下3條推論：

推論1：單位矩陣具有斜對稱特性。

推論2：若N維向量c=αa，其中α為一復數(shù)域上的數(shù)且α滿足斜對稱特性，則N維矩陣C=ccH也具有斜對稱特性。

推論3：當N維向量a和b滿足斜對稱特性時，c=a?b也具有斜對稱特性。

推論1、2 可通過將單位陣、c=αa和C=ccH代入式（13）進行矩陣運算得到，推論3 可見附錄1。由于均勻線陣的空域?qū)蛳蛄亢屠走_脈沖間隔相等時的時域?qū)蛳蛄烤鶟M足式（11），根據(jù)上述三條推論可以將斜對稱特性推廣至空時導向向量和CCM中。

通過引入MN×MN維斜對稱變換矩陣T，可以將斜對稱特性應用到CCM 的最大似然估計中。矩陣T的表達式如下所示

變換矩陣T滿足以下性質(zhì)［21］：

性質(zhì)1：對于任意斜對稱向量h∈CMN×1，Th是實向量；

性質(zhì)2：對于任意斜對稱矩陣H∈CMN×MN，THTH是實矩陣。

利用上述性質(zhì)，可以對CUT 中的數(shù)據(jù)和訓練樣本進行以下預處理：

利用是實矩陣的特性，可以證明其最大似然估計為［22］

其中，Re[·]表示取實部運算。與不使用斜對稱特性相比，式（18）的估計精度等效于將原訓練樣本數(shù)量擴充為原來的兩倍［22］。此外，變換后的CCM具有如下稀疏表示

3.2 算法求解和STAP濾波器設計

將式（18）與式（19）代入式（10）中，可以獲得如下優(yōu)化問題

該優(yōu)化問題的求解面臨兩方面的困難：一方面是參數(shù)?受環(huán)境影響較大，很難給出一個合適的值進行優(yōu)化；在另一方面，上式為一個NP 難問題，往往難以求解。因此本文引入SPICE 算法［16］，該算法是一種無參數(shù)稀疏恢復算法，同時具有較高的運算效率。SPICE 算法通過以下優(yōu)化問題求解空時功率譜p

文獻［23］證明上式對協(xié)方差中參數(shù)的估計結(jié)果可以漸進達到克拉美勞下界。為了描述方便，將變換后的CCM重新表示為

又因為最大似然估計為一致估計，則有：

根據(jù)上式結(jié)果，可以將（23）轉(zhuǎn)化為以下漸進等價的優(yōu)化問題

該問題具有半正定規(guī)劃形式，為凸優(yōu)化問題。另一方面，注意到約束項是累加求和的形式，這種形式正是加權l(xiāng)1范數(shù)約束，因此事實上該優(yōu)化問題為一個稀疏恢復問題。為了求解該優(yōu)化問題，引入矩陣，將優(yōu)化問題（27）等價變化為以下形式

注意該優(yōu)化問題中的待優(yōu)化變量包括bn以及矩陣G?？梢允褂醚h(huán)優(yōu)化算法［24］分別對兩者進行優(yōu)化：先固定矩陣G，對bn進行優(yōu)化；再固定bn，對G進行優(yōu)化，二者交替循環(huán)迭代。假設在第i次迭代中獲得了矩陣G(i)，需要求解。不妨記Q(i)=，則有：

且不等式的等號成立的條件為

不等式（30）中，左邊為原優(yōu)化問題的代價函數(shù)，右邊是與bn無關的常數(shù)，因此該不等式取等號時即為優(yōu)化問題（28）的解。根據(jù)等號成立的條件可知原優(yōu)化問題中關于bn的優(yōu)化結(jié)果為

將式（32）代入原優(yōu)化問題（28）中的第一個約束項，即可獲得

至此為止獲得了矩陣G被固定時，（28）對bn進行優(yōu)化的結(jié)果。接下來考慮固定bn，求解對G的優(yōu)化。假設在第i+1 次迭代中獲得了bn(i+1)，此 時（28）中第一個約束項不再有意義，該優(yōu)化問題可改寫為以下形式

結(jié)合式（32）、（33）和（35），可以獲得求解優(yōu)化問題（28）的循環(huán)優(yōu)化算法。

算法的初始值可設置如下

在獲得變換后的CCM 后，需根據(jù)線性約束最小方差方法準則設計STAP 濾波器，由于CUT 中的數(shù)據(jù)受到了變換矩陣T的預處理，為與之匹配，目標空時導向向量也要進行相應處理，最后獲得的STAP濾波器設計準則為

其中=Tvt，vt表示目標的空時導向向量，計算得到的最優(yōu)權向量為

綜上，提出的穩(wěn)健STAP算法實現(xiàn)流程如表1所示。

表1 算法步驟Tab.1 Procedure of the proposed algorithm

4 仿真實驗和計算復雜度分析

通過實驗評估所提出的STAP 算法，檢驗所提算法的雜波譜恢復性能和雜波抑制性能，并與IAASTAP、MIAA-STAP算法進行比較。為了分析引入CCM斜對稱特性對算法的影響，也將提出的算法與未應用斜對稱變換矩陣T預處理時的實驗結(jié)果進行比較，后者稱為SPICE-STAP。最后，對算法的復雜度進行分析。仿真設置的機載脈沖多普勒雷達系統(tǒng)參數(shù)如表2 所示。

表2 系統(tǒng)仿真參數(shù)Tab.2 System simulation parameters

從距離20 km 處選擇雜波樣本，雜波環(huán)寬度為100 m，設置每個雜波環(huán)劃分為360 個雜波塊并選擇網(wǎng)格分辨率參數(shù)ρd=ρs=10，收斂閾值設為κ=10-3。仿真結(jié)果通過100 次蒙特卡洛實驗的均值獲得。

4.1 雜波譜稀疏恢復性能

第一個實驗比較本文算法與上述各算法恢復出的雜波Capon譜，雜波Capon譜的表達式為

圖2（a）是仿真原始雜波譜，其為當R已知時計算出的雜波Capon 譜，是對仿真數(shù)據(jù)的超分辨結(jié)果。圖2（b）至（e）分別展示了IAA-STAP、MIAASTAP、SPICE-STAP 和本文算法估計得到的雜波譜，其中除了IAA-STAP 算法為SMV 類法，訓練樣本被設置為1 外，其他算法訓練樣本量均設置為3。在各距離單元統(tǒng)計特性相同的前提下，由于這四種算法都是建立在恢復雜波譜的基礎上估計CCM，而圖2（a）中假設CCM 已知，因此仿真原始雜波譜可以視為這些算法的性能上限。根據(jù)雜波譜恢復結(jié)果可以看出，由于單個樣本包含的信息有限，IAA-STAP 算法恢復出的雜波譜在整個平面上有較為嚴重的擴散。SPICE 算法相較于MIAA算法能夠有效抑制雜波譜的展寬，恢復出的雜波譜更加集中在雜波脊線上。本文算法在SPICESTAP 的基礎上利用了CCM 斜對稱結(jié)構，將等效訓練樣本的數(shù)量擴展成原來的兩倍，因此形成的空時雜波譜更加清晰，更接近圖2（a）中的仿真原始雜波譜。

4.2 雜波抑制性能分析

第二個實驗采用改善因子（Improvement Factor，IF）來評估各算法的雜波抑制能力，其定義式為［25］

圖3 與圖4 分別展示了理想情況與存在雜波內(nèi)部起伏（Intrinsic Clutter Motion，ICM）時的實驗結(jié)果，其中雜波信號表達式為［26］

式中α=[α0，α1，…，αM-1]∈CM×1為雜波信號的時域復起伏，其時間自相關函數(shù)為

式中σv為速度標準偏移量，設置為0.3。

圖3 和圖4 分別給出理想條件下和存在ICM 時的IF性能變化曲線，其中圖3（a）和圖4（a）展示了各算法在樣本量為3 且目標歸一化空間頻率為0 時IF隨歸一化多普勒頻率的變化曲線，圖3（b）和圖4（b）展示了各算法在目標歸一化多普勒頻率為0 時IF隨歸一化空間頻率的變化曲線，其樣本量與圖3（a）、4（a）相同。由圖3（a）和圖4（a）可以觀察到，IAA-STAP算法抑制效果最差，SPICE-STAP算法和MIAA-STAP 算法在歸一化多普勒頻率處于-1 至0.4 部分時雜波抑制效果近似，但是SPICE-STAP 算法的凹口展寬程度相較于后者更小，這是因為其恢復出的雜波譜更加集中在真實位置處。其中本文算法的IF 曲線凹口最窄，在各多普勒頻率上均優(yōu)于另外3 個算法。圖3、4（b）中各算法的IF 性能曲線趨勢與（a）中類似，在理想條件下和存在ICM 時本文算法的凹口展寬程度均為最小，在各空間頻率上都表現(xiàn)出更好的雜波抑制性能。圖3（c）與圖4（c）為目標歸一化多普勒頻率為0.5，歸一化空間頻率為0 時IF 隨樣本量的變化曲線。其中，當樣本量為1 時，MIAA-STAP 算法退化為IAA-STAP 算法。從圖3（c）與圖4（c）可以看出，樣本量增多時，SPICESTAP 算法的雜波抑制性能提升更大。由于利用了CCM 的斜對稱結(jié)構，本文算法的IF 性能與SPICESTAP取兩倍樣本量時的IF性能近似。

4.3 各距離單元輸出功率分析

本小節(jié)使用各算法計算出的空時濾波器對不同距離單元回波信號進行濾波，其中假設目標位于第200 個距離單元，且歸一化多普勒頻率為0.4，歸一化空間頻率為0。圖5（a）和（b）分別給出在理想條件下和存在ICM 時第180 個距離單元至第220 個距離單元的輸出功率曲線。由實驗結(jié)果可以看出IAA-STAP、MIAA-STAP、SPICE-STAP 和本文算法均能探測出目標，其中本文算法雜波抑制性能最優(yōu)，在不含目標的距離單元中輸出功率相比SPICESTAP 算法減少3 dB 左右，而SPICE-STAP 性能略微優(yōu)于MIAA-STAP，當存在ICM 時，這種優(yōu)勢進一步降低。IAA 算法雖然也能檢測出目標，由于其僅僅利用了一個距離單元中的數(shù)據(jù)，信息量有限，因此其雜波抑制性能最差。

4.4 算法復雜度分析

為比較各算法的計算復雜度，以算法所需復數(shù)乘法的次數(shù)作為計算量的衡量標準。算法復雜度的分析結(jié)果由表3給出，其中KIAA、KMIAA、KSPICE和KPerSPICE分別表示IAA-STAP、MIAA-STAP、SPICE-STAP 和本文方法所需要的迭代次數(shù)。因為NdNs遠大于MN，故每次迭代中IAA-STAP算法和MIAA-STAP算法的計算復雜度都可以視為O(2NdNs(MN)2)，而SPICESTAP 和本文算法的計算復雜度可以近似為O(4NdNs(MN)2)，因此在每次迭代中計算量約為前者的兩倍。為了更加客觀地評估各算法的計算量，需要分析其收斂速度，實驗結(jié)果如圖6 所示。可以看到，本文算法的收斂速度稍慢于MIAA-STAP，在設置收斂閾值κ=10-3的條件下，MIAA-STAP 算法和本文算法的平均迭代次數(shù)分別為17次和21次，算法所需收斂次數(shù)較為相近。綜上可知，本文提出的算法能夠提高精度，但是缺點是算法的計算復雜度更高，約為MIAA-STAP的兩倍。

表3 計算復雜度比較Tab.3 Comparison of computational complexity

5 結(jié)論

為進一步提升稀疏STAP算法的協(xié)方差估計精度和雜波抑制性能，本文提出一種聯(lián)合利用CCM斜對稱特性和雜波譜稀疏性的穩(wěn)健STAP算法。該算法利用CCM的斜對稱特性將等效訓練樣本量擴展至原來的兩倍，然后結(jié)合該特性和SPICE算法實現(xiàn)雜波譜的稀疏恢復并設計相應STAP濾波器。算法不需要人為設置超參數(shù)，實際應用中易于操作，仿真實驗表明，本文算法能夠準確地恢復出雜波譜，在理想條件和存在ICM條件下均保持較好的雜波抑制效果。

附錄1：推論3

附錄2：優(yōu)化問題（34）求解過程

下面將證明(34)的解為(35)。首先，將(35)代入(34)中的約束項，易于驗證等式成立。另一方面，將(35)代入(34)的待優(yōu)化函數(shù)中，有

如果可以證明對于任意矩陣G，以下不等式恒成立(X≥Y表示X-Y為半正定矩陣)，則可以證明(34)的解為(35)。

上式可以被等價寫為

其中，0表示相應維度的全零元素矩陣，且

綜上，證明了對于任意矩陣G，(52)恒成立，因此證明了固定時，優(yōu)化問題(34)關于G的解為(35)。