一種基于聯(lián)合稀疏恢復(fù)的空時自適應(yīng)處理方法

段克清 王澤濤 謝文沖 高 飛 王永良

①(空軍預(yù)警學(xué)院 武漢 430019)

②(國防科學(xué)技術(shù)大學(xué)電子科學(xué)與工程學(xué)院 長沙 410073)

1 引言

機載相控陣?yán)走_體制下，空時自適應(yīng)處理(Space-Time Adaptive Processing, STAP)可實現(xiàn)對強雜波的有效抑制以及慢速運動目標(biāo)的檢測[1,2]。該技術(shù)需利用相鄰多個距離門數(shù)據(jù)(訓(xùn)練樣本)估計待檢測距離門雜波協(xié)方差矩陣(Clutter Covariance Matrix, CCM)，進而計算空時自適應(yīng)濾波器權(quán)系數(shù)。根據(jù)RMB (Reed-Mallett-Brennan)準(zhǔn)則，為使得由CCM估計誤差所帶來的雜波抑制性能損失小于3 dB，所需獨立同分布訓(xùn)練樣本數(shù)至少為2倍系統(tǒng)自由度[3]。在實際應(yīng)用中，由于機載雷達往往面臨復(fù)雜多變的下視場景，常常難以獲取足夠的均勻訓(xùn)練樣本，因此小樣本、高性能STAP算法是目前該領(lǐng)域中的熱點研究方向。

近年來，壓縮感知或稀疏恢復(fù)技術(shù)被廣泛應(yīng)用于信號處理領(lǐng)域，并初步應(yīng)用于機載雷達STAP技術(shù)中[4－10]。由于天線的方向性，雜波在空時 2維譜上的分布是稀疏的，因此可利用少量甚至單個訓(xùn)練樣本基于超完備字典來獲取空時譜上強雜波點的位置，即稀疏系數(shù)中的非零元素位置。超完備字典可利用超分辨空時2維導(dǎo)向矢量構(gòu)造，求解稀疏系數(shù)的問題實則為稀疏約束下欠定方程求解問題。針對該問題的求解目前主要有3種途徑，即凸優(yōu)化方法[11－13]、貪婪算法[14－16]和 FOCUSS 算法[17]。其中，貪婪算法運算量較小，算法性能比較穩(wěn)定，應(yīng)用較為廣泛，因此本文重點討論貪婪算法在STAP技術(shù)中的應(yīng)用。目前，文獻[6-8]主要探討了利用單樣本恢復(fù)或多個樣本分別恢復(fù)后再聯(lián)合處理的稀疏恢復(fù)STAP算法，其本質(zhì)上仍為單觀測矢量情況(Single Measurement Vector, SMV)；文獻[9,10]針對稀疏恢復(fù)STAP問題進行了深入且系統(tǒng)的研究，一方面研究了基于稀疏濾波器的STAP技術(shù)，另一方面針對空時譜稀疏問題重點研究了非參數(shù)化迭代稀疏恢復(fù)算法，旨在增強算法對正則化參數(shù)設(shè)置的穩(wěn)健性和運算量的降低。然而，上述研究均未深入涉及多觀測矢量(Multiple Measurement Vector, MMV)聯(lián)合稀疏問題。

實際上，近年來在腦成像等領(lǐng)域的技術(shù)需求，激發(fā)了人們基于MMV進行聯(lián)合稀疏恢復(fù)的興趣[18]，也提出了一些較為有效的算法[18－23]。其中同時正交匹配追蹤算法(Simutaneous Orthogonal Matching Pursuit, SOMP)[18－21]是貪婪算法由SMV到MMV的直接拓展，具有較低的運算量，同時在無噪聲背景下具有理想的信號恢復(fù)性能，是目前最為高效和典型的聯(lián)合稀疏恢復(fù)算法。然而，在實際應(yīng)用背景下，機載雷達回波信號中除目標(biāo)和強雜波外不可避免地存在各種噪聲，因此，如果簡單將SOMP算法引入到STAP技術(shù)中，受噪聲影響，雜波抑制性能的穩(wěn)定性將得不到保證。

針對上述問題，本文提出了一種基于雜波子空間的貪婪算法進行聯(lián)合稀疏恢復(fù)，并應(yīng)用于STAP進行雜波抑制處理。首先，該方法充分利用了多個訓(xùn)練樣本中的雜波信息，因此較現(xiàn)有稀疏恢復(fù)STAP方法具備更強的信號恢復(fù)能力；其次，不同于 SOMP算法中利用剩余雜波矩陣來尋找最大相關(guān)原子的策略，該方法采用了剩余雜波子空間投影來尋找最大相關(guān)原子，并利用2l范數(shù)作為準(zhǔn)則來衡量相關(guān)性的大小，從而有效回避了雜波矩陣中的噪聲分量對稀疏恢復(fù)算法的不利影響，同時可高效地尋找到過完備字典中相關(guān)性最強的原子。本文最后利用仿真實驗驗證了所提STAP算法的有效性，并給出了相應(yīng)結(jié)論。

2 機載雷達信號模型及空時譜估計

不考慮距離模糊情況下，機載雷達在某個距離門上的回波信號主要由該距離環(huán)上多個離散雜波塊的回波疊加而成：

其中p為距離環(huán)上雜波塊個數(shù)，γp為第p個雜波塊的散射系數(shù)，s為目標(biāo)信號，n為噪聲，?p為第p個雜波塊對應(yīng)空時2維導(dǎo)向矢量：

其中, St,p和 Ss,p分別為第p個雜波塊對應(yīng)的時間導(dǎo)向列矢量和空間導(dǎo)向列矢量，且有

其中，θs,p為第p個雜波塊對應(yīng)的空間錐角，fd,p=2V cos( θs,p)/λ 為第p個雜波塊對應(yīng)多普勒頻率，V為載機速度，λ為波長，N為天線陣元數(shù)，d為陣元間距，M為相干脈沖數(shù)，fr為脈沖重復(fù)頻率。

由式(1)可知，機載雷達回波信號實際上由不同來向不同多普勒頻率的回波信號疊加而成，如果將所有空間頻率和多普勒頻率均遍歷并分別離散化為Ns=ρsN和Nd=ρdM個頻率點，且只考慮訓(xùn)練樣本中回波信號，則式(1)中信號還可表示為：

其中，ρs和ρd分別表示離散化程度，在高分辨情況下一般均遠大于1, αq為第q個空時網(wǎng)格單元幅度，?q為第q個空時網(wǎng)格單元對應(yīng)的空時2維導(dǎo)向矢量，。由于字典ψ的行數(shù)NM遠小于列數(shù) NsNd，因此式(5)為病態(tài)方程或欠定方程。

即從式(6)中求解稀疏系數(shù)矩陣A中的非零元素位置。需要注意的是，本文研究背景為正側(cè)視機載雷達，因此各距離門中雜波分布特性一致，也就是說各距離門中雜波的稀疏結(jié)構(gòu)可認(rèn)為是相同的。式(6)最終需恢復(fù)得到稀疏系數(shù)矩陣A中非零行的坐標(biāo)，即支撐集。

一旦求得A，則可計算出角-多普勒網(wǎng)格上雜波所在位置網(wǎng)格點的幅度，得到雜波空時2維譜的高分辨估計，進而可構(gòu)造空時自適應(yīng)權(quán)值進行雜波抑制處理。

3 基于子空間正交匹配追蹤的STAP方法(Subspace-based OMP, SbOMP)

3.1 算法原理

無噪聲情況下，貪婪算法可有效地尋找到支撐集。然而，當(dāng)觀測矢量中含有噪聲分量時，現(xiàn)有貪婪算法的稀疏恢復(fù)性能會出現(xiàn)較大波動。針對該問題，本文提出了一種基于子空間的貪婪算法進行聯(lián)合稀疏恢復(fù)，并利用其估計空時2維譜構(gòu)造自適應(yīng)權(quán)值進行雜波抑制處理。該方法對剩余雜波矩陣進行SVD或QR分解，得到剩余雜波子空間，再構(gòu)造投影矩陣基于l2范數(shù)匹配尋找字典中最大相關(guān)原子，然后將該原子對應(yīng)觀測矩陣(訓(xùn)練樣本)中雜波分量剔除，得到新的剩余雜波矩陣，再返回第1步，不斷迭代，直到滿足終止條件。在通過聯(lián)合稀疏恢復(fù)得到空時譜中強雜波分量位置后，基于最小二乘可計算出空時譜中雜波分量幅度，進而估計得到雜波協(xié)方差矩陣，最終形成空時2維自適應(yīng)權(quán)值進行雜波抑制處理。該方法具體表述如下：

輸入：

?NM× L 維訓(xùn)練樣本矩陣X。

?稀疏度K。

輸出：

?含有K個坐標(biāo)的集合ΛK。

?NM× NM維投影矩陣 PK。

?NM× L維近似矩陣 ZK。

?NM× L維剩余雜波矩陣 YK。

?NM× NM維雜波協(xié)方差矩陣R。

? NM× 1維空時導(dǎo)向矢量W。

程序：

第 1步：初始化剩余雜波矩陣Y0=X，坐標(biāo)集Λ0=?和迭代計數(shù)k=1。

第2步：計算第k次迭代后剩余雜波矩陣對應(yīng)的雜波子空間Sk?1以及該子空間對應(yīng)的投影算子Pk?1。

第3步：通過求解以下最優(yōu)化問題求解第k次迭代時的坐標(biāo)kλ

第4步：將第k次迭代時的坐標(biāo)λk并入坐標(biāo)集Λ

第5步：根據(jù)當(dāng)前坐標(biāo)集Λk對應(yīng)字典中的原子計算投影算子 Pk。

第6步：獲取新的雜波分量 Zk并根據(jù)該分量求得剩余雜波分量 Yk

第7步：令 k = k+ 1并返回第2步直到k = K。

第8步：計算雜波協(xié)方差矩陣R和空時自適應(yīng)權(quán)值W

第2步中可通過SVD分解或QR分解算法[24]來計算雜波子空間；第3步中，? Λk表示坐標(biāo)集?中除Λk內(nèi)所含坐標(biāo)以外的所有坐標(biāo)的集合；第7步中，迭代次數(shù)K與雜波秩有關(guān)，根據(jù)經(jīng)驗一般選取K為1.5～2倍雜波秩r；第8步中，ψΛ為由ψ中相對于Λ的原子組成的矩陣，AΛ為相對于Λ的稀疏系數(shù)矩陣，σ2為噪聲方差，目標(biāo)空時導(dǎo)向矢量ψ0可通過式(2)求取。

3.2 方法分析

已有貪婪算法均利用剩余雜波矩陣與各原子進行相關(guān)，依據(jù)相關(guān)性大小來篩選出與雜波分量最大相關(guān)的原子。無噪聲情況下，該類方法可高效篩選出r個與雜波相關(guān)性最強的原子。然而，存在噪聲時受噪聲分量影響，該類方法每次迭代并非恰好篩選出與雜波相關(guān)的原子，如果迭代次數(shù)依舊為r個，則篩選出與雜波相關(guān)的原子往往是欠缺的，因此一般迭代次數(shù)選取比r稍大來確保所有與雜波相關(guān)的原子均被搜索到。同時，與噪聲分量相關(guān)的部分原子不可避免地也被搜索到，表現(xiàn)為雜波譜上出現(xiàn)部分偽峰，導(dǎo)致雜波抑制性能下降。

本文方法利用剩余雜波子空間構(gòu)造投影算子來篩選與雜波分量相關(guān)的原子，保證了每次與原子進行相關(guān)運算的分量為純雜波分量，因此不會受到噪聲等其他因素的干擾。采用該算法每次篩選到的原子必然分布在理想雜波脊上，并隨著迭代次數(shù)的增加連續(xù)性增強。如果迭代次數(shù)太少，則體現(xiàn)在雜波譜上較為離散；而如果迭代次數(shù)太多，一方面運算量較大；另一方面會導(dǎo)致自動搜索到雜波脊附近相關(guān)性相對較強的原子；因此，根據(jù)經(jīng)驗值，一般選取迭代次數(shù)為1.5～2倍雜波秩r。

4 仿真實驗

選取正側(cè)視均勻線陣機載相控陣?yán)走_進行仿真實驗，其中載機高度為 9000 m，載機速度為 150 m/s，波長為0.3 m，重頻為2000 Hz，主波束俯仰角為3o，天線陣元個數(shù)N=8，相干脈沖數(shù)M=8，空間頻率和多普勒頻率超分辨因子ρs和ρd均設(shè)為4，待檢測距離門為300。共選取4種STAP方法進行對比，包括傳統(tǒng)的對角加載協(xié)方差求逆法[1,2]，文獻[6]中所提稀疏恢復(fù)方法，文獻[20]中所提聯(lián)合稀疏恢復(fù)算法和本文算法。方便起見，上述4種STAP算法在以下仿真實驗中分別簡稱為 LSMI, SR,SOMP和 SbOMP。需要注意的是，在實驗 1～3中雜噪比(Clutter-to-Noise Ratio, CNR)均設(shè)定為35 dB，實驗1～2中訓(xùn)練樣本數(shù)均選取為12(依照RMB準(zhǔn)則，樣本數(shù)應(yīng)不小于30，因此該情況為樣本嚴(yán)重不足情況)。

4.1 實驗1：雜波空時譜估計對比

本部分主要對比采用不同算法得到的空時譜分布情況，其中LSMI法采用傳統(tǒng)的MVDR譜。由圖1可以看出，采用稀疏恢復(fù)算法估計得到雜波譜的分辨率均顯著高于傳統(tǒng) MVDR法。其中，由圖1(b)可以看出，采用 SR法得到的譜分布離散化較為嚴(yán)重，且雜波脊附近區(qū)域存在偽峰；由圖1(c)可以看出，采用SOMP法在雜波脊線上可恢復(fù)得到約r個強雜波點，但在其他區(qū)域也存在大量離散偽峰，這主要是由于噪聲的影響所導(dǎo)致；由圖1(d)可以看出，采用所提SbOMP法可沿雜波脊精確恢復(fù)出連續(xù)的強雜波點，且在其他區(qū)域無偽峰出現(xiàn)。

圖1 采用不同算法估計空時譜分布情況Fig.1 Space-time spectra corresponding to different estimation algorithm

4.2 實驗2：訓(xùn)練樣本不足情況下雜波抑制性能對比

本節(jié)研究在訓(xùn)練樣本不足情況下上述4種算法的雜波抑制性能對比。本文采用改善因子(Improvement Factor, IF)作為基準(zhǔn)來衡量雜波抑制性能。其中，改善因子定義為輸出信雜噪比(Signal-to-Clutter-plus-Noise Ratio, SCNR)與輸入SCNR的比值，具體表述為：

由圖2可以得出以下3點結(jié)論：一是在訓(xùn)練樣本不足的情況下，采用稀疏恢復(fù)算法的性能要顯著優(yōu)于傳統(tǒng)STAP方法；二是采用聯(lián)合稀疏恢復(fù)算法的主雜波區(qū)的抑制性能要顯著優(yōu)于已有稀疏恢復(fù)算法，這主要是已有稀疏恢復(fù)算法沒有充分利用多個訓(xùn)練樣本中所蘊含的雜波信息；三是噪聲環(huán)境下本文所提SbOMB算法的雜波抑制性能要優(yōu)于經(jīng)典聯(lián)合稀疏恢復(fù)算法，這驗證了所提算法在噪聲環(huán)境下具有更好的穩(wěn)健性。

圖2 雜波抑制性能對比Fig.2 Clutter suppression performance comparison

4.3 實驗3：不同輸入CNR情況下雜波抑制性能對比

本節(jié)主要研究在不同噪聲背景下，所提基于子空間的聯(lián)合稀疏恢復(fù)算法與經(jīng)典聯(lián)合稀疏恢復(fù)算法的雜波抑制性能對比。其中，Known Cov.算法為已知雜波協(xié)方差矩陣下協(xié)方差矩陣求逆方法，相應(yīng)STAP處理器為最優(yōu)處理器，提供算法性能上限，以方便對比。本實驗中，采用各多普勒通道對應(yīng)改善因子的平均值，即平均改善因子(Average Improvement Factor, AIF)作為標(biāo)準(zhǔn)來衡量雜波抑制性能。共進行50次蒙特卡諾仿真。由圖3可以看出，在CNR低于35 dB的情況下，所提SbOMP算法的雜波抑制性能比SOMP算法可改善約3～10 dB；而在噪聲相對較弱(雜噪比較大)情況下，相較于SOMP算法，SbOMP算法的雜波抑制性能也有約1 dB的改善。

5 結(jié)論

本文將經(jīng)典聯(lián)合稀疏恢復(fù)算法應(yīng)用于機載雷達STAP處理，并針對其所存在的問題提出了一種基于子空間的聯(lián)合稀疏恢復(fù)STAP方法。研究表明，相比于已有稀疏恢復(fù)算法，采用聯(lián)合稀疏恢復(fù)STAP算法后，對主瓣雜波的抑制性能有了顯著提升；同時，所提基于子空間聯(lián)合稀疏恢復(fù) STAP算法解決了經(jīng)典聯(lián)合稀疏恢復(fù)算法在噪聲背景下性能下降的問題。實際上，本文算法由于涉及子空間處理，因此運算量較經(jīng)典聯(lián)合稀疏恢復(fù)算法有所提升；此外，在極少訓(xùn)練樣本情況下，本文算法的旁瓣抑制性能并不理想，下一步重點針對這些問題開展進一步研究。

圖3 不同輸入CNR情況下平均IF對比Fig.3 Average IF comparison of different input CNRs

[1]Klemm R.Principles of Space-Time Adaptive Processing[M].London: Institute of Electical Engineering, 2006.

[2]王永良, 彭應(yīng)寧.空時自適應(yīng)信號處理[M].北京: 清華大學(xué)出版社, 2000.Wang Y L and Peng Y N.Space-Time Adaptive Processing[M].Beijing: Tsinghua University Press, 2000.

[3]Reed I S and Mallett J D.Rapid convergence rate in adaptive arrays[J].IEEE Transactions on Aerospace and Electronics Systems, 1974, 10(6): 853-863.

[4]Selesnick I W, Pillai S U, Li K Y, et al..Angle-Doppler processing using sparse regularization[C]. IEEE International Conference on in Acoustics, Speech and Signal Processing, Dallas, TX, USA, 2010: 2750-2753.

[5]Li J, Zhu X, Stoica P, et al..Iterative space-time adaptive processing[C].IEEE Digital Signal Processing Workshop,Marco Island, FL, 2009: 440-445.

[6]Sun K, Zhang H, Li G, et al..A novel STAP algorithm using sparse recovery technique[C].IEEE International Conference on Geoscience & Remote Sensing Symposium, Cape Town,2009: 336-339.

[7]Sun K, Meng H, Wang Y, et al..Direct data domain STAP using sparse representation of clutter spectrum[J].Signal Processing, 2011, 91(9): 2222-2236.

[8]Sun K, Meng H, Lapierre F D, et al..Registration-based compensation using sparse representation in conformal-array STAP[J].Signal Processing, 2011, 91(10): 2268-2276.

[9]Yang Z C, de Lamare R C, and Li X.l1-Regularized STAP algorithms with a generalized sidelobe canceler architecture for airborne radar[J].IE EE Transactions on Signal Processing, 2012, 60(2): 674-686.

[10]Yang Z C, de Lamare R C, and Li X.Sparsity-aware Space-Time Adaptive Processing algorithms with l1-norm regularization for airborne radar[J].IET Signal Processing,2012, 6(5): 413-423.

[11]Chen S, Donoho D, and Saunders M.Atomic decomposition by basis pursuit[J].SIAM Review, 2001, 43(1): 129-159.

[12]Tropp J A.Just relax: convex programming methods for identifying sparse signals in noise[J].IEEE Transactions on Information Theory, 2006, 53(3): 1030-1051.

[13]Wright S J, Nowak R D, and Figueiredo M A T.Sparse reconstruction by separable approximation[J].IEEE Transactions on Signal Processing, 2009, 57(7): 2479-2493.

[14]Tropp J A.Greed is good: algorithmic results for sparse approximation[J].IEE E Transactions on Information Theory, 2004, 50(10): 2231-2242.

[15]Tropp J A and Gilbert A C.Signal recovery from random measurements via orthogonal matching pursuit[J].IEEE Transactions on Information Theory, 2007, 53(12): 4655-4666.

[16]Needell D and Vershynin R.Signal recovery from incomplete and inaccurate measurements via regularized orthogonal matching pursuit[J].IEEE Journal of Selected Topics on Signal Processing, 2009, 4(2): 310-316.

[17]Gorodnitsky I F and Rao B D.Sparse signal reconstruction from limited data using FOCUSS: a re-weighted minimum norm algorithm[J].IEEE Transactions on Signal Processing,1997, 45(3): 600-616.

[18]Cotter S, Rao B, Engan K, et al..Sparse solutions to linear inverse problems with multiple measurement vectors[J].IEEE Transactions on Signal Processing, 2005, 53(7):2477-2488.

[19]Tropp J, Gilbert A, and Strauss M.Algorithms for simultaneous sparse approximation.Part I: Greedy pursuit[J].Signal Processing, 2006, 86(3): 572-588.

[20]Chen J and Huo X.Theoretical results on sparse representations of multiple-measurement vectors[J].IEEE Transactions on Signal Processing, 2006, 54(12): 4634-4643.

[21]Gribonval R, Rauhut R, Schnass K, et al..Atomsof all channels, unite! average case analysis of multi-channel sparse recovery using greedy algorithms[J].Journal of Fourier Analysis Applications, 2008, 14(5/6): 655-687.

[22]Tropp J A.Algorithms for simultaneous sparse approximation.Part II: Convex relaxation[J].Signal Processing, 2006, 86(3): 589-602.

[23]Wipf D P and Rao B D.An empirical Bayesian strategy for solving the simultaneous sparse approximation problem[J].IEEE Transactions on Signal Processing, 2007, 55(7):3704-3716.

[24]Goldstein J S and Reed I S.Theory of partially adaptive radar[J].IEEE Transactions on Aerospace and Electronic Systems, 1997, 33(4): 1309-1325.