基于圖拉普拉斯嵌入的合成孔徑雷達(dá)時變窄帶干擾抑制算法

傅東寧, 廖桂生, 黃 巖, 張邦杰, 王 幸

(1. 西安電子科技大學(xué)雷達(dá)信號處理國家重點(diǎn)實(shí)驗(yàn)室, 陜西 西安 710071;2. 東南大學(xué)毫米波國家重點(diǎn)實(shí)驗(yàn)室, 江蘇 南京 210096)

0 引 言

合成孔徑雷達(dá)(synthetic aperture radar,SAR)以其高分辨成像和全天候全天時工作的能力獲得了廣泛關(guān)注[1-2]。由于頻譜沖突問題,寬帶SAR系統(tǒng)經(jīng)常受到電磁干擾影響[3],常見的電磁干擾如通信系統(tǒng)和調(diào)頻(frequency modulation,FM)系統(tǒng)等窄帶干擾(narrow-band interference,NBI)信號[4-5],其帶寬通常比SAR信號帶寬小得多(文獻(xiàn)[6]將帶寬比率定為1%)。實(shí)際上,SAR系統(tǒng)采用二維匹配濾波處理產(chǎn)生距離和方位上的高分辨率圖像,本身具有潛在的抗干擾能力。但是,NBI影響過大,降低了系統(tǒng)的輸出信干噪比(signal to interference and noise ratio, SINR),導(dǎo)致SAR圖像質(zhì)量嚴(yán)重下降?？偟膩碚f,當(dāng)前抗NBI方法主要分為3類:參數(shù)方法、非參數(shù)方法以及半?yún)?shù)方法[7]。

非參數(shù)方法主要在其他表征域中將NBI和真實(shí)信號分離,并設(shè)計(jì)濾波器抑制干擾。究其本質(zhì),非參數(shù)方法主要利用了NBI功率遠(yuǎn)大于真實(shí)信號這一特點(diǎn)。其優(yōu)勢在于處理效率高且工程上容易實(shí)現(xiàn),能夠避免干擾模型不準(zhǔn)確所導(dǎo)致的性能損失,典型算法如陷波濾波器[8]、最小均方(least mean square, LMS)濾波器[9]和特征子空間投影(eigenva-lue subspace projection, ESP)方法[10]等。這些方法在實(shí)際系統(tǒng)中具有良好的性能與計(jì)算效率,但它們僅利用了干擾信號和有用回波之間的功率差,而沒有考慮干擾信號本身的固有特性,因此沒有實(shí)現(xiàn)性能的最大化。

參數(shù)化方法不同于非參數(shù)化方法,其基于正弦模型表征干擾,并估計(jì)模型中的頻率和復(fù)振幅。為了準(zhǔn)確估計(jì)正弦模型參數(shù),諸如迭代自適應(yīng)方法(iterative adaptive approach,IAA)[11]和RELAX(relaxation)[12-13]等多種算法被用于逐個脈沖提取NBI。但提取過程中,有用回波信號會受到破壞,且必須在多脈沖基礎(chǔ)上進(jìn)行多次迭代提取參數(shù),計(jì)算代價巨大。

近年來稀疏處理[13]和低秩恢復(fù)技術(shù)[14]在信號處理領(lǐng)域得到廣泛應(yīng)用?；谙∈枵齽t化和低秩恢復(fù)算法可以在抑制NBI的同時保護(hù)有用回波。由于這些方法的性能通過超參數(shù)優(yōu)化問題得到進(jìn)一步優(yōu)化,因此將其歸為半?yún)?shù)方法?；谡倚盘柲Ｐ?NBI在整個合成孔徑時間內(nèi)具有低秩特性,但在單次采樣快拍內(nèi)是稀疏的。魯棒主成分分析(robust principal component analysis,RPCA)[7,15]方法、無參分解方法[16]、矩陣分解(matrix factorization, MF)[17]法和復(fù)加權(quán)張量分解(complex reweighted tensor factorization, CRTF)[18]等算法都假定NBI在合成孔徑時間內(nèi)具有穩(wěn)定的頻譜特性,都屬于線性低秩恢復(fù)方法。然而實(shí)際中,對于脈間波形和頻率變化的NBI,其頻率峰值會在較小的頻帶范圍內(nèi)局部改變。

本文分析了在時變NBI情況下的信號模型和其固有的非線性低秩表達(dá)式,采用圖論[19]方法,將接收信號中位于不同方位的采樣視作頂點(diǎn),并構(gòu)造頂點(diǎn)之間的邊以挖掘其潛在的局部聯(lián)系,再利用這些聯(lián)系對局部時變干擾進(jìn)行聚類。最后,本文提出一種基于圖論的算法,即用圖拉普拉斯矩陣表示整個圖形,以提取NBI聚類,并用稀疏正則化來保護(hù)有用回波信號。本文稱之為圖拉普拉斯嵌入(graph Laplacian embedding,GLE)算法,同時為整個優(yōu)化問題的每一步嚴(yán)格推導(dǎo)出了閉式解,且具有較快的收斂速度。通過對實(shí)測的SAR和NBI數(shù)據(jù)進(jìn)行處理,證明了所提算法的有效性。

1 信號模型

NBI通常具有較高的功率,且占用相對固定的頻段,即使在很長的距離外,也會與SAR系統(tǒng)使用的頻段產(chǎn)生沖突,降低其成像質(zhì)量。SAR信號模型一般將接收到的信號堆疊到二維作用域中,例如距離-方位,或者快時間域和慢時間域[15],即

y(t,ts)=x(t,ts)+r(t,ts)+n(t,ts)

(1)

式中:t代表快時間;ts代表慢時間;x(t,ts)為真實(shí)回波信號;r(t,ts)為NBI;n(t,ts)為噪聲。通常NBI的功率很高,而參數(shù)方法將NBI建模為多個復(fù)雜正弦波的總和,即

(2)

將公式(11)中所有時域采樣整合,接收信號矩陣Y∈CM×N可以表示為

Y=X+R+N

(3)

式中:X,R,N分別為式(1)中對應(yīng)的真實(shí)回波數(shù)據(jù)矩陣、NBI矩陣和噪聲矩陣。接下來將詳細(xì)分析時變NBI的特性。

多脈沖NBI信號頻譜如圖1所示,其在距離頻域和方位時域中產(chǎn)生多個“條紋”。但實(shí)際上,對于具有脈沖間頻率變化特性的NBI信號,其實(shí)時頻率會圍繞某一中心頻率振蕩。圖2(a)為一個實(shí)測NBI信號的頻譜。圖2(b)為從圖2(a)所示的NBI中提取的頻率峰值變化情況?？梢钥闯?該NBI頻率峰值在不同脈沖的頻帶內(nèi),會在局部范圍內(nèi)呈現(xiàn)跳躍現(xiàn)象。

(4)

式中:‖·‖*表示核范數(shù);‖·‖1表示矩陣的1范數(shù);λ和δ是超參數(shù)。實(shí)際上,PCA將線性流形上的數(shù)據(jù)進(jìn)行了嵌入。從另一個角度來看,圖形結(jié)構(gòu)具有一個基本假設(shè),即局部時變干擾位于平滑的非線性低維流形上或其附近。因此,可以將圖形嵌入技術(shù)用于NBI聚類,并且通過解決以下優(yōu)化問題來獲得嵌入的干擾分量:

(5)

(6)

式中:I是單位矩陣;D是度矩陣,定義為

(7)

式中:diag代表對角化。

然后,同時利用低秩主投影和非線性GLE,將式(4)和式(5)所代表的問題組合為基于圖論的NBI抑制模型:

(8)

式中:γ是超參數(shù)。

通過GLE約束,可以同時保留干擾矩陣R的線性與非線性特征,也就是干擾的主要成分以及不同脈沖的干擾間局部變化。由此可將局部時變NBI聚類,并更好地從接收信號中提取出來。接下來,本文將集中于如何解決上述優(yōu)化問題。

2 用于NBI抑制的GLE算法

交替方向乘子法(alternating direction method of multipliers,ADMM)非常適合快速求解大規(guī)模優(yōu)化問題[7,21]。本節(jié)中,將使用ADMM算法來求解式(8)所代表的優(yōu)化問題。

首先,引入輔助變量Z來重新表述問題:

(9)

將式(9)中的增廣拉格朗日函數(shù)寫為

(10)

式中:Y1和Y2是拉格朗日乘子;μ1和μ2是懲罰系數(shù)。然后在保持其他變量不變的基礎(chǔ)上,交替更新其中一個變量。在第k+1次迭代中,干擾矩陣R的子優(yōu)化問題可以表達(dá)為

(11)

則R可以更新為

(12)

式中:UΣcVH表示矩陣C的SVD分解,矩陣C為

(13)

并且,奇異值門限可以表示為

(14)

接下來,真實(shí)回波矩陣X的子優(yōu)化問題可以表達(dá)為

(15)

式(15)所代表的問題可以通過軟閾值方法[16]求解:

(16)

其中,軟閾值可以表示為

(17)

輔助變量Z可以由以下子優(yōu)化問題表示:

(18)

對式(18)求導(dǎo)可得

(19)

最后,更新拉格朗日乘子和懲罰參數(shù),并通過迭代求解,直到整個算法的剩余誤差達(dá)到設(shè)定的門限值(如10-2)或達(dá)到最大循環(huán)次數(shù),算法工作流程結(jié)束。

以上即為所提方法的主要步驟,在進(jìn)行干擾抑制前首先需要判斷干擾存在與否,或者是否需要進(jìn)行干擾抑制(針對干擾強(qiáng)度較弱的情況)。該過程可以采用文獻(xiàn)[2]或文獻(xiàn)[22]中方法來自適應(yīng)設(shè)置門限進(jìn)行判斷,如文獻(xiàn)[2]中自適應(yīng)門限定義為

(20)

式中:μ0和σ0分別表示無干擾的脈沖回波均值和標(biāo)準(zhǔn)差;erf-1(·)表示誤差函數(shù)逆變換;α表示虛警概率。

由此,對回波信號做短時傅里葉變換,求回波信號的負(fù)熵,具體步驟可參考文獻(xiàn)[2]。當(dāng)負(fù)熵超過門限值,則需要抑制干擾。此時,可采用本文所提方法進(jìn)行干擾抑制,然后再進(jìn)行SAR成像;若未超過門限,則可以判定無需采取干擾抑制而直接進(jìn)行SAR成像。整個算法的主要流程如圖3所示。

3 實(shí)驗(yàn)結(jié)果

本節(jié)將使用實(shí)測SAR以及NBI數(shù)據(jù)分析所提GLE算法的干擾抑制性能。系統(tǒng)參數(shù)如表1所示。實(shí)驗(yàn)中的SAR圖像數(shù)據(jù)和NBI數(shù)據(jù)是分別收集的,這樣可以調(diào)整兩個數(shù)據(jù)間的比例(功率),來模擬不同輸入SINR的情況。為了更好地定量評價所提方法的有效性,采用輸出均方根誤差(root mean square error, RMSE)來量化干擾抑制算法的性能,這里RMSE定義為

表1 SAR參數(shù)

(21)

實(shí)驗(yàn)中使用距離多普勒(range-Doppler,RD)算法生成SAR圖像,根據(jù)第1節(jié)所述,將相鄰頂點(diǎn)數(shù)設(shè)為5來構(gòu)造接收信號的圖形結(jié)構(gòu),并設(shè)定輸入SINR=-20 dB,則無NBI影響的原始圖像和受NBI影響的圖像如圖4所示。為更清晰展現(xiàn)NBI影響效果,將圖4中矩形框區(qū)域放大,如圖5所示。可見整個場景被高功率NBI影響,完全無法識別地面場景的有用信息。

對受干擾數(shù)據(jù)分別采用潔化(CLEAN)算法[5]、RPCA算法、特征子空間投影(eigenvalue subspace projection, ESP)算法和所提出的GLE算法進(jìn)行干擾抑制處理，結(jié)果如圖6所示。

圖6中矩形框區(qū)域放大圖如圖7所示。其中，CLEAN算法作為參數(shù)方法的代表，ESP算法作為非參數(shù)算法的代表，而RPCA算法和本文所提算法則屬于半?yún)?shù)算法。

由圖6和圖7可以看出,本文所提的GLE算法與其他方法相比具有更好的干擾抑制性能,即圖像中殘留的干擾分量更少。而ESP算法和CLEAN算法作為非參數(shù)方法和參數(shù)方法,并沒有考慮保護(hù)有用回波,即無法在抑制干擾的同時避免對有用回波信號產(chǎn)生影響。而RPCA算法沒有對NBI的時變特性建模,與本文提出的GLE算法相比,圖像中殘留更多的噪聲。表2列出了所有方法的RMSE,這也證明了本文所提算法的出色性能。

表2 不同抗干擾算法的RMSE結(jié)果

實(shí)際中,干擾功率是變化的,因此更關(guān)心不同SINR下干擾抑制性能的穩(wěn)健性。不同SINR條件下各算法的RMSE在表3中列出。調(diào)整輸入SINR分別為-15 dB和-10 dB,并利用圖5中同樣的方法進(jìn)行了實(shí)驗(yàn),結(jié)果如圖8和圖9所示。

表3 不同輸入SINR下各抗干擾算法的RMSE結(jié)果

由表3所示不同SINR下的干擾抑制算法的RMSE結(jié)果以及圖6、圖8和圖9可以看出,所提GLE算法在不同輸入SINR條件下均能保持穩(wěn)健的抗干擾性能。原因是所提方法充分考慮了NBI的時變特性在脈沖間的聯(lián)系,并繼承了半?yún)?shù)方法的優(yōu)勢,保護(hù)了有用回波信號,因此得到了最優(yōu)的干擾抑制性能。

4 結(jié) 論

本文提出一種用于SAR系統(tǒng)抑制NBI的GLE算法,該算法基于NBI的局部時變特性,利用圖論構(gòu)造了不同脈沖之間的拉普拉斯嵌入聯(lián)系,時變干擾被聚類于非線性低維流形中,減輕其對接收信號的影響。同時，使用實(shí)際的SAR數(shù)據(jù)以及實(shí)測NBI數(shù)據(jù)，通過分析所提GLE算法的NBI抑制性能和處理結(jié)果，證明了所提算法的有效性。