陳 誠, 宋曉驥, 何志華, 劉 濤, 曹來保, 粟 毅

(國防科技大學(xué)電子科學(xué)學(xué)院, 湖南 長沙 410073)

0 引 言

探地雷達(ground penetrating radar,GPR)利用電磁波對介質(zhì)的穿透性,可以實現(xiàn)對地下目標的成像探測,如地雷、管線等[1-3]。相比基于聲學(xué)、光學(xué)、化學(xué)和核物理等的技術(shù),基于 GPR 的探測技術(shù)得益于電磁波對介質(zhì)的穿透能力、非接觸的測量方式和有效獲取目標散射特性的能力,能夠更好地滿足地下探測特別是地雷探測對準確性和安全性的要求,因而受到了研究者的密切關(guān)注[4-5]。

在GPR探測時,接收到的雷達回波中通常不僅包含目標響應(yīng),還包含了地表反射回波和天線間直接耦合等。尤其是地表反射回波,其往往遠強于目標信號,因此在雷達圖像中目標響應(yīng)往往被掩蓋于強雜波之中,導(dǎo)致目標檢測十分困難[6]。雜波抑制方法研究是GPR探測領(lǐng)域中的一個重要方向,傳統(tǒng)雜波抑制方法如平均值對消(mean subtraction, MS)法[7]、時間門技術(shù)[8]、奇異值分解(singular value decomposition, SVD)法[9-11]等,在地面平整均勻且參數(shù)已知的情況下雜波抑制效果較好,但當實際地表粗糙起伏時,往往由于環(huán)境復(fù)雜,難以準確估計地表回波時延或是雜波子空間,導(dǎo)致抑制效果變差甚至失效。形態(tài)學(xué)成分分析(morphological component analysis,MCA)[12-14]通過結(jié)合圖像的稀疏表示理論,利用信號組成成分的形態(tài)差異將圖像分解為目標和雜波,其性能依賴于稀疏表示字典的構(gòu)建。

低秩稀疏分解(low rank and sparse decomposition,LRSD)是指將探測數(shù)據(jù)分解為一個低秩矩陣與一個稀疏矩陣的疊加,在模式識別、圖像處理、雷達探測領(lǐng)域均有應(yīng)用,主要用于提取低秩背景或稀疏目標[15-17]。近年來,穩(wěn)健主成分分析(robust principle component analysis,RPCA)將LRSD理論應(yīng)用于雷達穿透成像領(lǐng)域以實現(xiàn)雜波抑制[16-17],并展示了其應(yīng)用潛力。由于GPR在探測時,目標響應(yīng)往往被地表強反射淹沒,因此雜波以地表背景為主,且由于目標(主要包括地雷等)數(shù)量有限,通常呈稀疏分布,因此在GPR領(lǐng)域,目標圖像可以被視為稀疏成分,而雜波可以表征為低秩矩陣[17]。由于矩陣的秩是非凸非連續(xù)的,直接計算其解,將原始圖像分解為目標和雜波信號是非確定性多項式(non-deterministic polynomial,NP)困難問題。RPCA通過求解核范數(shù)最小化,將其松弛為凸優(yōu)化問題,但在核范數(shù)最小化時引入的SVD迭代導(dǎo)致計算量增大[11],降低了求解效率。此外,將核范數(shù)作為正則化項時,大奇異值在優(yōu)化過程中起決定性作用,容易造成矩陣信息損失[18-19]。文獻[20]提出加權(quán)核范數(shù)最小化(weighted nuclear norm minimization,WNNM)的思想,即用小權(quán)值約束大奇異值,用大權(quán)值約束小奇異值,有利于更好地恢復(fù)出低秩矩陣,但該方法仍存在計算效率低的問題[21]。文獻[22]提出隨機SVD算法,將高維矩陣投影到低維后再進行分解運算,能夠在提高計算效率的同時保持精度。

針對上述問題,本文提出一種基于低秩稀疏分解的雜波抑制方法。該方法將加權(quán)核范數(shù)正則化項引入RPCA方法,通過控制核范數(shù)的權(quán)重,有利于精確估計出低秩矩陣,并利用隨機SVD代替核范數(shù)計算時的奇異值分解以減少計算量,結(jié)合交替方向乘子(alternating direction method of multipliers,ADMM)法,建立帶有增廣拉格朗日項的優(yōu)化方程求解目標函數(shù),以進一步提高求解效率[23-26]。本文實驗應(yīng)用仿真與實測數(shù)據(jù),通過與多種方法進行對比,定量和定性地分析了本文方法的性能優(yōu)勢,表明所提方法能夠在GPR應(yīng)用中有效實現(xiàn)雜波抑制。

1 LRSD及RPCA理論

根據(jù)低秩稀疏分解理論,對于GPR系統(tǒng),可以將獲取到的m×n維(假設(shè)m≥n)雷達數(shù)據(jù)矩陣D進行分解[27],即

D=L+S+N

(1)

式中:L表示含少數(shù)非零奇異值的低秩矩陣;S為含少數(shù)非零項的稀疏矩陣;N表示隨機噪聲,其Frobenius范數(shù)滿足

(2)

RPCA方法將上述分解轉(zhuǎn)換為如下易解的凸優(yōu)化問題[15]:

(3)

2 基于LRSD的雜波抑制方法

鑒于核范數(shù)在低秩稀疏分解中的局限性,將加權(quán)核范數(shù)引入到RPCA方法中,通過控制核范數(shù)的權(quán)重來約束各奇異值,提高求解的精確性。矩陣L的加權(quán)核范數(shù)定義為

(4)

式中:w=[w1,w2,…,wn]T;σi(L)為矩陣L的第i個奇異值;wi表示該奇異值對應(yīng)的權(quán)重。

因此,式(3)中的問題可轉(zhuǎn)換為

(5)

在本算法中,定義加權(quán)核范數(shù)的權(quán)重值為

wi=μ/[σi(D)+ε]

(6)

式中,μ為正則化參數(shù);ε為大于0的極小值?？梢钥闯?當ε→0時,有如下關(guān)系:

表明定義的加權(quán)核范數(shù)能有效逼近數(shù)據(jù)矩陣的秩。

為了提高算法效率,在進行矩陣分解時,引入基于隨機SVD算法[28]過程的隨機SVD算法,即對m×n維雷達數(shù)據(jù)矩陣D,生成n×l維高斯隨機矩陣X(l遠小于n)。根據(jù)圖1所示流程生成樣本矩陣并進行正交三角分解,通過正交變換Q構(gòu)建遠小于原數(shù)據(jù)矩陣D的l×n維矩陣E。由于正交變換的穩(wěn)定性,在變換過程中能夠保持原矩陣的各類性質(zhì),因此對E進行奇異值分解能夠在保持運算精度的條件下提高運算效率。

圖1 隨機SVD流程圖Fig.1 Flowchart of randomized SVD

根據(jù)ADMM算法思想求解式(5),引入等式約束,得到增廣拉格朗日函數(shù)為

(7)

式中:U表示拉格朗日乘子矩陣;(·)H表示對矩陣進行共軛轉(zhuǎn)置;β為懲罰系數(shù)。通過引入拉格朗日乘子項為優(yōu)化函數(shù)增加等式約束,能夠?qū)⒛繕撕瘮?shù)的解約束于可行域內(nèi),提升算法的求解效率。

在求解式(7)中的3個變量時,可通過分別固定其他兩個變量、更新其中一個變量的方式,具體迭代過程如下:

(8)

(9)

Uk+1=Uk+β(D-Lk+1-Sk+1)

(10)

式中:k表示雷達數(shù)據(jù)第k次的迭代值。

綜上,可提出本文基于LRSD的雜波抑制方法流程,如算法1所示。

算法 1 基于LRSD的雜波抑制方法輸入 原始雷達矩陣D,正則化系數(shù)λ,β,最大迭代次數(shù)K,精度閾值η。輸出 低秩矩陣L和稀疏矩陣S步驟 1 初始化L0=S0=U0=0。步驟 2 更新低秩雜波矩陣Lk+1=argminLLw,*+(β/2)D-L-Sk+(Uk/β)2F計算核范數(shù)最小化時,通過隨機SVD算法對待分解矩陣進行降維處理。步驟 3 更新稀疏目標矩陣Sk+1=argminSλS1+(β/2)D-Lk+1-S+(Uk/β)2F步驟 4 更新拉格朗日乘子矩陣Uk+1=Uk+β(D-Lk+1-Sk+1)步驟 5 若不滿足停止準則,繼續(xù)步驟2～步驟4,若滿足停止準則,跳出循環(huán),得到輸出低秩雜波矩陣L和稀疏目標圖像矩陣S。

3 實驗結(jié)果和分析

為了驗證所提低秩稀疏分解方法的有效性,本節(jié)利用仿真與實測數(shù)據(jù)來驗證其雜波抑制性能,并通過與MS、SVD、MCA和RPCA這4種方法處理后,將成像結(jié)果進行對比,證明算法在成像質(zhì)量和速度上的優(yōu)勢。

3.1 仿真實驗設(shè)計

仿真通過GPRMax[29]電磁仿真軟件獲取場景數(shù)據(jù),實驗場景為表面略有起伏(不超過10 mm)的土壤層,相對介電常數(shù)為5,土壤內(nèi)部埋藏有兩個相對介電常數(shù)為2.9的塑料圓柱,用于模擬埋藏在土壤中的地雷,如圖2所示。左側(cè)地雷直徑為60 mm,高為40 mm,距離地面90 mm,右側(cè)地雷直徑為60 mm,高為40 mm,距離地面100 mm,兩地雷間隔為25 mm。天線掃描面距離介質(zhì)表面150 mm,掃描步進為10 mm,發(fā)射信號為雷克子波,頻率為2.2 GHz。

圖2 仿真實驗場景示意圖Fig.2 Schematic diagram of simulation test

3.2 實測實驗設(shè)計

實驗1通過如圖3所示的雷達探測系統(tǒng)進行,實驗?zāi)繕藶槿舾?2式反步兵地雷,地雷直徑為79 mm,高度為35 mm,等間隔埋藏于沙地中,埋藏深度約為100 mm,每兩個地雷之間間隔為50 mm。天線距離地面約200 mm,掃描頻率范圍為0.75～4.75 GHz,空間采樣間隔為20 mm。

圖3 實驗1場景圖Fig.3 Schematic diagram of Experiment 1

實驗2借鑒了喬治亞理工學(xué)院出版的雷達數(shù)據(jù)[30-31],實驗頻率范圍為60 MHz～8.06 GHz,頻段內(nèi)共計401個頻點,頻點在頻率范圍內(nèi)均勻分布,天線相位中心距離地面278 mm,空間采樣間隔為20 mm。實驗場景及目標布設(shè)情況如圖4所示,這里選取x=45 cm處的數(shù)據(jù)(即圖4中紅線標注位置)進行數(shù)據(jù)處理及對比實驗。對上述數(shù)據(jù)進行成像,在此基礎(chǔ)上分別利用MS、SVD、MCA、RPCA與本文所提方法進行處理,處理結(jié)果如圖5～圖7所示。

圖4 實驗2場景圖Fig.4 Schematic diagram of Experiment 2

圖5 仿真實驗結(jié)果圖Fig.5 Simulation experiment results

圖6 實驗1結(jié)果圖Fig.6 Experiment 1 imaging results

圖7 實驗2結(jié)果圖Fig.7 Experiment 2 imaging results

3.3 實驗結(jié)果與分析

圖5所示為對仿真數(shù)據(jù)使用各算法處理前后的成像結(jié)果。其中,圖5(a)為直接成像結(jié)果,可以看出,土壤回波對目標回波干擾嚴重。圖5(c)和圖5(d)表明,由于土壤表層的起伏,平均值對消法與SVD方法在抑制雜波時幾近失效。從圖5(e)可以看出,MCA方法盡管能夠有效增強圖像中的目標,但其構(gòu)建字典的特性使得目標幾何結(jié)構(gòu)發(fā)生變化。圖5(b)和圖5(f)中目標圖像明顯,表明低秩稀疏分解方法能夠在保留目標完整性的前提下有效抑制成像結(jié)果中的雜波,且圖5(b)中殘留的雜波明顯弱于圖5(f),表明本文方法抑制雜波的效果優(yōu)于RPCA方法。

圖6為對實驗1數(shù)據(jù)處理前后的結(jié)果。其中,圖6(a)為處理前的成像結(jié)果,從中可以看到明顯的強地表回波,埋藏目標無法顯現(xiàn);從圖6(b)～圖6(f)可以看出,幾種算法均能對成像結(jié)果中的雜波進行抑制。但從雜波抑制效果來看,本文方法效果最好,處理結(jié)果中殘留的雜波相對于目標已經(jīng)十分微弱,而平均值對消法、SVD方法與RPCA方法處理結(jié)果中仍存在著明顯雜波;MCA方法處理后視覺效果較好,但在處理結(jié)果中引入了大量規(guī)則雜波。圖7為實驗2結(jié)果圖,從中可以看出原始圖像分別經(jīng)過平均值對消法、MCA方法、SVD方法及兩種低秩稀疏分解方法處理后,目標信息逐漸增強。對比圖7與圖4中的目標分布,RPCA方法處理結(jié)果中存在虛警目標,而本文方法能更好地反映出目標的真實分布情況。圖8為對3組數(shù)據(jù)使用本文方法分別計算出的低秩雜波圖像。

圖8 本文方法處理各實驗數(shù)據(jù)后的低秩雜波圖像Fig.8 Low-rank clutter images with experimental data processed by the proposed method

為了量化評估各方法的細節(jié)表現(xiàn),表1記錄了各方法處理前后結(jié)果圖的目標信雜比(signal-to-clutter ratio, SCR),即圖像目標區(qū)域平均功率與雜波區(qū)域平均功率的比值,其定義為

表1 不同方法處理后成像結(jié)果SCR參數(shù)Table 1 SCR of imaging results processed by different methods dB

(11)

式中:Rt和Rc分別表示目標和雜波區(qū)域;Nt和Nc分別表示目標和雜波區(qū)域內(nèi)的像素點個數(shù);I(p)和I(q)分別表示目標和雜波區(qū)域各點的像素值大小。

從表1中的SCR評價指標來看,采用低秩稀疏分解方法處理后,圖像的SCR明顯高于其他方法,且本文方法SCR也始終高于RPCA方法,證實了所提方法抑制雜波的有效性。

為了進一步評估本文方法的運算效率,表2列出了各方法的處理時間及計算復(fù)雜度(設(shè)圖像維度為m×n,迭代次數(shù)為r。由于MCA方法參數(shù)初始化過程復(fù)雜,未評估其計算復(fù)雜度)。由表2可以看出,對于各實驗中的單幅圖像,所提方法處理時間均不足1 s,且運算效率與RPCA方法相比提升了5倍以上。綜上分析,本文方法能夠兼顧精度與效率,在提升雜波抑制效果的同時,提高了運算效率。

表2 運算時間和計算復(fù)雜度對比Table 2 Comparison of running time and computation complexity

4 結(jié) 論

本文提出了一種改進的低秩稀疏分解雜波抑制算法。該方法在RPCA方法的基礎(chǔ)上,引入加權(quán)核范數(shù)正則化項與隨機SVD計算方法,并結(jié)合ADMM方法求解目標函數(shù),實現(xiàn)低秩稀疏矩陣的準確、高效分解。多組數(shù)據(jù)的實驗結(jié)果驗證了本文方法抑制雜波的有效性與優(yōu)越性,且證實本文算法運算效率優(yōu)于RPCA方法。但當GPR探測僅關(guān)注特定目標時,該方法無法排除場景中由其他雜波物體造成的虛警,因此結(jié)合目標散射機理研究目標與雜波物體的特性差異將是后續(xù)研究工作的重點。