李家強(qiáng),胡張燕,姚昌華,郭桂祥,陳金立

(南京信息工程大學(xué) 電子與信息工程學(xué)院,南京 210044)

0 引 言

合成孔徑雷達(dá)(Synthetic Aperture Radar,SAR)是最早提出并投入實用的成像雷達(dá),廣泛用于各類重點場景監(jiān)視[1-2]。由雷達(dá)分辨率理論和奈奎斯特采樣定理可知,合成孔徑雷達(dá)系統(tǒng)性能的提高通常伴隨著采樣數(shù)據(jù)量的增加,這給系統(tǒng)設(shè)計和實現(xiàn)帶來了困難。

近年來,隨著壓縮感知[3]和相關(guān)理論的發(fā)展,稀疏信號處理在SAR成像中廣泛應(yīng)用,這表明當(dāng)回波信號具有稀疏性或可壓縮性時,能夠以遠(yuǎn)低于奈奎斯特的采樣頻率,用較少的觀測數(shù)據(jù)高概率地恢復(fù)出原信號。傳統(tǒng)的壓縮感知重構(gòu)算法如貝葉斯壓縮感知(Bayesian Compressed Sensing,BCS)[4]算法,存在重構(gòu)精度不高、抗噪性能差等問題?，F(xiàn)有的壓縮感知SAR成像模型中,L1正則化是最常見的重構(gòu)模型[5]。該模型能夠有效抑制場景稀疏時的旁瓣和噪聲,但L1范數(shù)判罰函數(shù)是凸的,解析解是有偏估計[6],會導(dǎo)致重建結(jié)果的偏差效應(yīng)并最終影響重建精度。針對出現(xiàn)的問題,許多研究集中在非凸正則化和加權(quán)L1正則化算法上。如極小極大凹罰(Minimax Concave Penalty,MCP)函數(shù)[7]、Lq范數(shù)判罰(0

針對上述問題,本文提出一種基于自適應(yīng)加權(quán)的非凸正則化和全變分的稀疏SAR成像算法,將非凸判罰函數(shù)與TV-norm罰函數(shù)線性結(jié)合,在成像模型中形成一個復(fù)合正則化器,采用交替方向乘子法求解SAR成像模型,并使用方位-距離解耦算子替換測量矩陣及其厄米特轉(zhuǎn)置以減少存儲空間且加快信號處理速度。實驗結(jié)果表明,本文算法具有很好的重構(gòu)性能與抗噪能力。

1 SAR回波信號二維稀疏表示模型

圖1為SAR成像幾何模型。假設(shè)雷達(dá)平臺沿y軸以恒定的速度v前行,若雷達(dá)發(fā)射線性調(diào)頻(Linearly Frequency Modulated,LFM)信號,則理想點目標(biāo)P(x,y)的基帶回波為

圖1 合成孔徑雷達(dá)成像幾何模型

(1)

式中:D={(x,y)|xmin≤x≤xmax,ymin≤y≤ymax};tr和ta分別為距離向快時間和方位向慢時間;r為點目標(biāo)P(x,y)到雷達(dá)航線的垂直距離;g(x,y)為目標(biāo)后向散射系數(shù)矩陣;aa(·)為方位向窗函數(shù);ar(tr)=rect(tr/Tp)為發(fā)射信號包絡(luò),Tp為脈沖持續(xù)時間;c為光速;γ為LFM信號線性調(diào)頻率;λ為波長;R(ta;r)為ta時刻雷達(dá)與目標(biāo)間距離,

(2)

忽略距離徙動后理想點目標(biāo)的回波經(jīng)采樣后為

S=ΦrGΦa。

(3)

式中:S∈M×N為均勻采樣回波信號,M為距離向采樣點數(shù),N為方位向采樣點數(shù);G∈M×N為觀測場景目標(biāo)點的后向散射矩陣;矩陣Φr∈M×M為距離向卷積矩陣,其列是距離向參考函數(shù)的共軛與快時間域中采樣數(shù)的卷積,其中距離向參考函數(shù)為矩陣Φa∈N×N為方位向卷積矩陣,其行是方位向參考函數(shù)的共軛與慢時間域中采樣數(shù)的卷積,其中方位向參考函數(shù)為為多普勒調(diào)頻率。

對信號在距離向與方位向均進(jìn)行降采樣,如圖2所示,其中灰色部分代表采樣點信號。此時均勻采樣回波信號S∈M×N變?yōu)殡S機(jī)采樣的回波信號S2CS∈P×Q;P(P≤M)為距離向隨機(jī)采樣點數(shù),Q(Q≤N)方位向隨機(jī)采樣點數(shù)。距離向卷積矩陣Φr變?yōu)棣礡CS∈P×M,方位向卷積矩陣Φa變?yōu)棣礎(chǔ)CS∈N×Q,二維隨機(jī)降采樣合成孔徑雷達(dá)回波信號的欠定方程為

圖2 SAR回波數(shù)據(jù)隨機(jī)采樣

S2CS=ΦRCSGΦA(chǔ)CS。

(4)

考慮到現(xiàn)實場景中加性噪聲影響,式(4)可表示為

S2ZCS=ΦRCSGΦA(chǔ)CS+NZ。

(5)

式中:S2CS∈P×Q和S2ZCS∈P×Q為隨機(jī)采樣回波信號;ΦRCS∈P×M為距離維觀測矩陣;ΦA(chǔ)CS∈N×Q為方位維觀測矩陣;G∈M×N為待重構(gòu)目標(biāo)圖像;NZ∈P×Q為復(fù)高斯白噪聲信號。

式(5)表示的信號模型是二維卷積模型,可用兩個連續(xù)一維方程表示。首先,令X∈M×Q作為方位維壓縮數(shù)據(jù),表達(dá)式為

X=GΦA(chǔ)CS,

(6)

則式(5)可寫為

S2ZCS=ΦRCSX+NZ。

(7)

式(6)表示方位向壓縮感知模型,式(7)表示距離向壓縮感知模型。

求解式(6)和式(7)的過程如下:首先,求解欠定方程式(7)得到距離向壓縮后的信號;其次,將接收到的回波S2ZCS分解為列向量形式,其表達(dá)式為

S2ZCS=[s1,s2,…,sQ]。

(8)

同時,X分解為列向量形式:

X=[x1,x2,…,xQ]。

(9)

則式(7)表示為

sq=ΦRCSxq+nz,q=1,2,…,Q。

(10)

式中:sq∈P×1;xq∈M×1;nz∈P×1;Q為方位向接收到的脈沖數(shù)。

稀疏SAR成像過程如式(10)可視為一個線性逆問題,求解線性逆問題的常用算法是最小正則化線性最小二乘代價函數(shù),其表達(dá)式為

(11)

式中:p(·)為懲罰項或正則項;λ為正則化參數(shù)。為提高xq的稀疏性,應(yīng)選L0范數(shù)作為xq的判罰函數(shù)。但L0范數(shù)屬于非確定性多項式時間復(fù)雜度問題,需組合搜索解決,維度增加時難以實現(xiàn),需要修改判罰函數(shù)將其轉(zhuǎn)化為等價問題求解。對于稀疏信號,當(dāng)測量矩陣ΦRCS滿足有限等距性質(zhì)條件時,L0范數(shù)優(yōu)化問題等價于L1范數(shù)優(yōu)化問題,即將非凸問題轉(zhuǎn)換為凸問題。基于L1正則化的SAR成像模型為

(12)

(13)

2 基于TV和非凸正則化的合成孔徑雷達(dá)成像

2.1 基于MC函數(shù)的稀疏SAR成像

基于L1正則化的SAR成像模型會影響重建精度,而基于非凸正則化的SAR成像模型因其非凸函數(shù)的幾何性質(zhì)更接近L0范數(shù)可以在提高稀疏性的同時減少偏差。基于非凸正則化的SAR成像模型為

(14)

式中:pNC(·)表示非凸判罰。

求解非凸優(yōu)化問題關(guān)鍵在于計算近端算子問題,本文采用迭代收縮閾值(Iterative Shrinkage Threshoding,IST)[12]算法來求解近端算子問題。

(15)

(16)

λpNC(xq)。

(17)

(18)

(19)

(20)

本文采用方位-距離解耦算子代替測量矩陣及其厄米特轉(zhuǎn)置加快處理速度。首先假設(shè)G(·)是原始數(shù)據(jù)生成算子,I(·)是成像算子,方位-距離解耦方案[14〗被描述為

(21)

2.2 基于MC函數(shù)和TV的稀疏SAR成像

非凸正則化旨在增強(qiáng)SAR圖像中基于點的特征。為提高目標(biāo)重建精度,在稀疏正則化基礎(chǔ)上引入TV-norm。對于二維場景,圖像向量xq的離散各項同性TV-norm定義如下:

(22)

(23)

式中:

Du|A|i,j=|A[i+1,j]|-|A[i,j]|;

(24)

Dv|A|i,j=|A[i,j+1]|-|A[i,j]|。

(25)

將MC判罰函數(shù)和TV-norm判罰函數(shù)線性組合成一個復(fù)合正則化器,則非凸正則化的SAR成像模型式(14)表示為

(26)

在變量分裂方法的指導(dǎo)下,將無約束問題式(26)轉(zhuǎn)化為以下等價的約束優(yōu)化問題:

(27)

式中:z1和z2是兩個輔助變量。根據(jù)拉格朗日乘數(shù)法原理,上述優(yōu)化問題可由以下公式解決:

(28)

式中:l1和l2為拉格朗日乘數(shù)。相應(yīng)最優(yōu)解表達(dá)式為

(29)

為最小化L(xq,z1,z2,l1,l2),對變量xq,z1,z2采用交替最小化方法,同時隨著迭代逐漸增加拉格朗日乘數(shù)l1和l2的值。迭代過程如下:

(30)

(31)

(32)

式(30)中的最小化,由于目標(biāo)函數(shù)是二次的,產(chǎn)生了一個有解的線性方程組:

(33)

式中:I是單位矩陣。式(31)的最優(yōu)解為

(34)

式中:λMC=λ1/(2l1)。

采用Chambolle算法[15]對TV-norm正則化對偶問題求解從而得到式(32)的解。首先對式進(jìn)行參數(shù)變換,變換后為式(35),等價于標(biāo)準(zhǔn)TV-norm正則化問題。

(35)

式中:λTV=λ2/(2l2)。上述優(yōu)化問題解決方案為

(36)

式中:p=(p1,p2)∈Y是z2∈X的對偶變量,X表示為歐幾里得空間。因本文圖像大小為M×N,則X表示為RM×N,Y=X×X。對p可用梯度下降法[16]進(jìn)行求解,過程如下:

(37)

(38)

同時,將方向-距離解耦因子也應(yīng)用到MC判罰與全變分的稀疏SAR成像模型中。

2.3 基于自適應(yīng)加權(quán)MC函數(shù)和TV的稀疏SAR成像

本文將非凸判罰與全變分的正則化模型同加權(quán)正則化模型相結(jié)合,提出一種自適應(yīng)加權(quán)的極小極大凹罰函數(shù)和全變分的稀疏SAR成像模型,將非凸正則化的SAR成像模型表示為

(39)

(40)

式中:φ(xq(i,j))是按式(41)定義的標(biāo)量加權(quán)MC函數(shù),

(41)

標(biāo)量MC函數(shù)φ(xq)被表示為模函數(shù)與標(biāo)量Huber函數(shù)之間的差值,即φ(xq)=|xq|-ρ(xq),其中ρ(xq)是標(biāo)量Huber函數(shù),表達(dá)式為

(42)

經(jīng)過推導(dǎo)與整理,最終標(biāo)量MC函數(shù)φ(xq)為

(43)

式中:β≥0;λ>0;γ>1。根據(jù)上述推導(dǎo),式(39)的優(yōu)化問題表示為

(44)

基于自適應(yīng)加權(quán)MC函數(shù)的表達(dá)式如式(44)所示,再引入TV-norm,過程如2.2節(jié)所示,區(qū)別是輔助變量z1的表達(dá)式不同。

(45)

本文采用交替最小化方法解決式(44)優(yōu)化問題,首先對輔助變量β按式(46)進(jìn)行更新:

(46)

進(jìn)一步,

(47)

對xq進(jìn)行如下更新:

(48)

則式(48)所示優(yōu)化問題可以采用快速軟閾值迭代算法(Fast Iterative Shrinkage Thresholding Algorithm,FISTA)[17]來求解。之后,對Λ進(jìn)行更新:

(49)

綜上,基于自適應(yīng)加權(quán)的極小極大凹罰函數(shù)和全變分的算法偽代碼如下:

1 當(dāng)res>ε且t

5Λ(t+1)=β(t+1)

11t=t+1

12 結(jié)束循環(huán)

3 仿真與實測數(shù)據(jù)結(jié)果與分析

為綜合比較分析本文方法的性能及其工程實用性,將L1正則化算法、BCS算法、MC算法、MC-TV算法與本文算法(AWMC-TV)作為對比,同時采用相對均方根誤差(Relative Root Mean Square Error,RRMSE)[7]、目標(biāo)雜波比(Target Clutter Ratio,TCR)[18]和圖像熵(Image Entropy,IE)[18]來定量評價成像質(zhì)量。

3.1 目標(biāo)二維成像仿真實驗分析

為充分驗證算法性能,對23個點目標(biāo)仿真產(chǎn)生的數(shù)據(jù)進(jìn)行二維成像仿真實驗,仿真參數(shù)設(shè)置如表1所示。目標(biāo)散射點模型如圖3所示。將圖4中RD算法的成像圖作為真實值,對本文成像算法性能進(jìn)行評價。為便于分析,將距離向采樣率與方位向采樣率同時設(shè)為原采樣率的1/2。仿真結(jié)果如圖5所示,可見除基于L1算法成像結(jié)果中有少量假目標(biāo)點外,其余算法均可清晰實現(xiàn)目標(biāo)點SAR成像。

表1 雷達(dá)仿真參數(shù)

圖3 目標(biāo)散射點模型

圖4 RD算法SAR圖像

圖5 各種算法二維成像結(jié)果

為進(jìn)一步比較上述算法的成像質(zhì)量,利用目標(biāo)雜波比和圖像熵對成像質(zhì)量進(jìn)行定量評價。表2給出了在采樣率為原采樣率1/2時幾種算法的TCR和IE值,顯然基于AWMC-TV算法的圖像TCR值大于其他算法,IE值小于其他算法,相比于其他算法有更好的聚焦性能和重建精度。

表2 目標(biāo)的成像結(jié)果評價參數(shù)

采用相對均方根誤差作為評估標(biāo)準(zhǔn),反應(yīng)重建圖像中每個散射點的信號強(qiáng)度與真實值間誤差,該值越小說明重建效果越好。表3顯示了本文算法相比其他算法RRMSE較小,重構(gòu)效果更好。

表3 幾種算法RRMSE

分別對各算法進(jìn)行計算復(fù)雜度對比,如表4所示。表4中,K1和K2分別表示各算法的循環(huán)迭代次數(shù),Ti和Ji(i=1,2,3)分別表示各算法的外部和內(nèi)部迭代次數(shù)。

表4 各類算法的計算復(fù)雜度比較

綜上,在相同的實驗條件下,相比其他幾種算法,本文算法雖然復(fù)雜度有所增加,但目標(biāo)的重建精度顯著提高,表明該方法在SAR成像領(lǐng)域具有一定的潛力和優(yōu)勢。

3.2 實測數(shù)據(jù)處理

為充分驗證本文算法在實際應(yīng)用中的有效性,采用加拿大溫哥華的RADARSAT-1精細(xì)模式2(加拿大航天局版權(quán)所有)的實測數(shù)據(jù)進(jìn)行實驗,參數(shù)設(shè)置如表5所示。圖6(a)為基于RD算法的英國灣附近溫哥華地區(qū)的實際場景,選取紅色矩形區(qū)域內(nèi)的6個目標(biāo)點(6艘貨船)進(jìn)行實驗。這6個目標(biāo)點相比海平面是強(qiáng)散射點,相比英吉利海灣是稀疏的。圖6(b)～(f)為幾種算法在距離采樣率為原采樣率的30%、方位采樣率為原采樣率20%的成像結(jié)果。顯然,在低采樣率情況下,AWMC-TV算法依然能夠準(zhǔn)確重建目標(biāo)。

表5 溫哥華場景RADARSAT-1參數(shù)

圖6 各算法成像結(jié)果

為比較幾種算法在有噪聲的真實情況下的性能,對圖6(a)紅框中綠色標(biāo)記的點目標(biāo)沿距離向進(jìn)行切片,計算該點目標(biāo)幅值在所在距離單元內(nèi)的比值。表6給出了該目標(biāo)點振幅在各算法下占該距離單元內(nèi)的百分比,可見AWMC-TV算法下目標(biāo)點幅值在該距離向上占比最大,表明AWMC-TV算法能夠消除圖像中的模糊和拖尾,在重建性能和抗噪性能上具有優(yōu)越性。

表6 目標(biāo)點幅度在該距離單元內(nèi)占比

4 結(jié)束語

針對L1正則化存在的偏差效應(yīng)及傳統(tǒng)壓縮感知重構(gòu)算法重構(gòu)精度低的問題,本文將隸屬于非凸函數(shù)簇中的極小極大凹罰函數(shù)應(yīng)用到SAR成像中,同時將非凸正則化與全變分判罰線性組合作為成像模型的復(fù)合正則化器,提出一種基于自適應(yīng)加權(quán)的極小極大凹罰函數(shù)和全變分的稀疏SAR成像算法,并采用FISTA算法與交替最小化算法求解該模型。仿真與實測數(shù)據(jù)成像結(jié)果驗證了該算法的有效性,在抗噪性能及重構(gòu)精度等方面比其他算法更具優(yōu)勢。

由于本文中壓縮感知模型只考慮到噪聲的影響,因此后續(xù)還需研究雜波對算法的影響。此外,能否通過多種變換找到復(fù)雜場景的稀疏表示是壓縮感知理論進(jìn)一步在SAR成像領(lǐng)域應(yīng)用的關(guān)鍵,這也是下一步研究的內(nèi)容。