李家強(qiáng),胡張燕,姚昌華,郭桂祥,陳金立
(南京信息工程大學(xué) 電子與信息工程學(xué)院,南京 210044)
合成孔徑雷達(dá)(Synthetic Aperture Radar,SAR)是最早提出并投入實(shí)用的成像雷達(dá),廣泛用于各類重點(diǎn)場(chǎng)景監(jiān)視[1-2]。由雷達(dá)分辨率理論和奈奎斯特采樣定理可知,合成孔徑雷達(dá)系統(tǒng)性能的提高通常伴隨著采樣數(shù)據(jù)量的增加,這給系統(tǒng)設(shè)計(jì)和實(shí)現(xiàn)帶來了困難。
近年來,隨著壓縮感知[3]和相關(guān)理論的發(fā)展,稀疏信號(hào)處理在SAR成像中廣泛應(yīng)用,這表明當(dāng)回波信號(hào)具有稀疏性或可壓縮性時(shí),能夠以遠(yuǎn)低于奈奎斯特的采樣頻率,用較少的觀測(cè)數(shù)據(jù)高概率地恢復(fù)出原信號(hào)。傳統(tǒng)的壓縮感知重構(gòu)算法如貝葉斯壓縮感知(Bayesian Compressed Sensing,BCS)[4]算法,存在重構(gòu)精度不高、抗噪性能差等問題?,F(xiàn)有的壓縮感知SAR成像模型中,L1正則化是最常見的重構(gòu)模型[5]。該模型能夠有效抑制場(chǎng)景稀疏時(shí)的旁瓣和噪聲,但L1范數(shù)判罰函數(shù)是凸的,解析解是有偏估計(jì)[6],會(huì)導(dǎo)致重建結(jié)果的偏差效應(yīng)并最終影響重建精度。針對(duì)出現(xiàn)的問題,許多研究集中在非凸正則化和加權(quán)L1正則化算法上。如極小極大凹罰(Minimax Concave Penalty,MCP)函數(shù)[7]、Lq范數(shù)判罰(0 針對(duì)上述問題,本文提出一種基于自適應(yīng)加權(quán)的非凸正則化和全變分的稀疏SAR成像算法,將非凸判罰函數(shù)與TV-norm罰函數(shù)線性結(jié)合,在成像模型中形成一個(gè)復(fù)合正則化器,采用交替方向乘子法求解SAR成像模型,并使用方位-距離解耦算子替換測(cè)量矩陣及其厄米特轉(zhuǎn)置以減少存儲(chǔ)空間且加快信號(hào)處理速度。實(shí)驗(yàn)結(jié)果表明,本文算法具有很好的重構(gòu)性能與抗噪能力。 圖1為SAR成像幾何模型。假設(shè)雷達(dá)平臺(tái)沿y軸以恒定的速度v前行,若雷達(dá)發(fā)射線性調(diào)頻(Linearly Frequency Modulated,LFM)信號(hào),則理想點(diǎn)目標(biāo)P(x,y)的基帶回波為 圖1 合成孔徑雷達(dá)成像幾何模型 (1) 式中:D={(x,y)|xmin≤x≤xmax,ymin≤y≤ymax};tr和ta分別為距離向快時(shí)間和方位向慢時(shí)間;r為點(diǎn)目標(biāo)P(x,y)到雷達(dá)航線的垂直距離;g(x,y)為目標(biāo)后向散射系數(shù)矩陣;aa(·)為方位向窗函數(shù);ar(tr)=rect(tr/Tp)為發(fā)射信號(hào)包絡(luò),Tp為脈沖持續(xù)時(shí)間;c為光速;γ為LFM信號(hào)線性調(diào)頻率;λ為波長;R(ta;r)為ta時(shí)刻雷達(dá)與目標(biāo)間距離, (2) 忽略距離徙動(dòng)后理想點(diǎn)目標(biāo)的回波經(jīng)采樣后為 S=ΦrGΦa。 (3) 式中:S∈M×N為均勻采樣回波信號(hào),M為距離向采樣點(diǎn)數(shù),N為方位向采樣點(diǎn)數(shù);G∈M×N為觀測(cè)場(chǎng)景目標(biāo)點(diǎn)的后向散射矩陣;矩陣Φr∈M×M為距離向卷積矩陣,其列是距離向參考函數(shù)的共軛與快時(shí)間域中采樣數(shù)的卷積,其中距離向參考函數(shù)為矩陣Φa∈N×N為方位向卷積矩陣,其行是方位向參考函數(shù)的共軛與慢時(shí)間域中采樣數(shù)的卷積,其中方位向參考函數(shù)為為多普勒調(diào)頻率。 對(duì)信號(hào)在距離向與方位向均進(jìn)行降采樣,如圖2所示,其中灰色部分代表采樣點(diǎn)信號(hào)。此時(shí)均勻采樣回波信號(hào)S∈M×N變?yōu)殡S機(jī)采樣的回波信號(hào)S2CS∈P×Q;P(P≤M)為距離向隨機(jī)采樣點(diǎn)數(shù),Q(Q≤N)方位向隨機(jī)采樣點(diǎn)數(shù)。距離向卷積矩陣Φr變?yōu)棣礡CS∈P×M,方位向卷積矩陣Φa變?yōu)棣礎(chǔ)CS∈N×Q,二維隨機(jī)降采樣合成孔徑雷達(dá)回波信號(hào)的欠定方程為 圖2 SAR回波數(shù)據(jù)隨機(jī)采樣 S2CS=ΦRCSGΦA(chǔ)CS。 (4) 考慮到現(xiàn)實(shí)場(chǎng)景中加性噪聲影響,式(4)可表示為 S2ZCS=ΦRCSGΦA(chǔ)CS+NZ。 (5) 式中:S2CS∈P×Q和S2ZCS∈P×Q為隨機(jī)采樣回波信號(hào);ΦRCS∈P×M為距離維觀測(cè)矩陣;ΦA(chǔ)CS∈N×Q為方位維觀測(cè)矩陣;G∈M×N為待重構(gòu)目標(biāo)圖像;NZ∈P×Q為復(fù)高斯白噪聲信號(hào)。 式(5)表示的信號(hào)模型是二維卷積模型,可用兩個(gè)連續(xù)一維方程表示。首先,令X∈M×Q作為方位維壓縮數(shù)據(jù),表達(dá)式為 X=GΦA(chǔ)CS, (6) 則式(5)可寫為 S2ZCS=ΦRCSX+NZ。 (7) 式(6)表示方位向壓縮感知模型,式(7)表示距離向壓縮感知模型。 求解式(6)和式(7)的過程如下:首先,求解欠定方程式(7)得到距離向壓縮后的信號(hào);其次,將接收到的回波S2ZCS分解為列向量形式,其表達(dá)式為 S2ZCS=[s1,s2,…,sQ]。 (8) 同時(shí),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)可視為一個(gè)線性逆問題,求解線性逆問題的常用算法是最小正則化線性最小二乘代價(jià)函數(shù),其表達(dá)式為 (11) 式中:p(·)為懲罰項(xiàng)或正則項(xiàng);λ為正則化參數(shù)。為提高xq的稀疏性,應(yīng)選L0范數(shù)作為xq的判罰函數(shù)。但L0范數(shù)屬于非確定性多項(xiàng)式時(shí)間復(fù)雜度問題,需組合搜索解決,維度增加時(shí)難以實(shí)現(xiàn),需要修改判罰函數(shù)將其轉(zhuǎn)化為等價(jià)問題求解。對(duì)于稀疏信號(hào),當(dāng)測(cè)量矩陣ΦRCS滿足有限等距性質(zhì)條件時(shí),L0范數(shù)優(yōu)化問題等價(jià)于L1范數(shù)優(yōu)化問題,即將非凸問題轉(zhuǎn)換為凸問題。基于L1正則化的SAR成像模型為 (12) (13) 基于L1正則化的SAR成像模型會(huì)影響重建精度,而基于非凸正則化的SAR成像模型因其非凸函數(shù)的幾何性質(zhì)更接近L0范數(shù)可以在提高稀疏性的同時(shí)減少偏差。基于非凸正則化的SAR成像模型為 (14) 式中:pNC(·)表示非凸判罰。 求解非凸優(yōu)化問題關(guān)鍵在于計(jì)算近端算子問題,本文采用迭代收縮閾值(Iterative Shrinkage Threshoding,IST)[12]算法來求解近端算子問題。 (15) (16) λpNC(xq)。 (17) (18) (19) (20) 本文采用方位-距離解耦算子代替測(cè)量矩陣及其厄米特轉(zhuǎn)置加快處理速度。首先假設(shè)G(·)是原始數(shù)據(jù)生成算子,I(·)是成像算子,方位-距離解耦方案[14〗被描述為 (21) 非凸正則化旨在增強(qiáng)SAR圖像中基于點(diǎn)的特征。為提高目標(biāo)重建精度,在稀疏正則化基礎(chǔ)上引入TV-norm。對(duì)于二維場(chǎng)景,圖像向量xq的離散各項(xiàng)同性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ù)線性組合成一個(gè)復(fù)合正則化器,則非凸正則化的SAR成像模型式(14)表示為 (26) 在變量分裂方法的指導(dǎo)下,將無約束問題式(26)轉(zhuǎn)化為以下等價(jià)的約束優(yōu)化問題: (27) 式中:z1和z2是兩個(gè)輔助變量。根據(jù)拉格朗日乘數(shù)法原理,上述優(yōu)化問題可由以下公式解決: (28) 式中:l1和l2為拉格朗日乘數(shù)。相應(yīng)最優(yōu)解表達(dá)式為 (29) 為最小化L(xq,z1,z2,l1,l2),對(duì)變量xq,z1,z2采用交替最小化方法,同時(shí)隨著迭代逐漸增加拉格朗日乘數(shù)l1和l2的值。迭代過程如下: (30) (31) (32) 式(30)中的最小化,由于目標(biāo)函數(shù)是二次的,產(chǎn)生了一個(gè)有解的線性方程組: (33) 式中:I是單位矩陣。式(31)的最優(yōu)解為 (34) 式中:λMC=λ1/(2l1)。 采用Chambolle算法[15]對(duì)TV-norm正則化對(duì)偶問題求解從而得到式(32)的解。首先對(duì)式進(jìn)行參數(shù)變換,變換后為式(35),等價(jià)于標(biāo)準(zhǔn)TV-norm正則化問題。 (35) 式中:λTV=λ2/(2l2)。上述優(yōu)化問題解決方案為 (36) 式中:p=(p1,p2)∈Y是z2∈X的對(duì)偶變量,X表示為歐幾里得空間。因本文圖像大小為M×N,則X表示為RM×N,Y=X×X。對(duì)p可用梯度下降法[16]進(jìn)行求解,過程如下: (37) (38) 同時(shí),將方向-距離解耦因子也應(yīng)用到MC判罰與全變分的稀疏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)化問題,首先對(duì)輔助變量β按式(46)進(jìn)行更新: (46) 進(jìn)一步, (47) 對(duì)xq進(jìn)行如下更新: (48) 則式(48)所示優(yōu)化問題可以采用快速軟閾值迭代算法(Fast Iterative Shrinkage Thresholding Algorithm,FISTA)[17]來求解。之后,對(duì)Λ進(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) 為綜合比較分析本文方法的性能及其工程實(shí)用性,將L1正則化算法、BCS算法、MC算法、MC-TV算法與本文算法(AWMC-TV)作為對(duì)比,同時(shí)采用相對(duì)均方根誤差(Relative Root Mean Square Error,RRMSE)[7]、目標(biāo)雜波比(Target Clutter Ratio,TCR)[18]和圖像熵(Image Entropy,IE)[18]來定量評(píng)價(jià)成像質(zhì)量。 為充分驗(yàn)證算法性能,對(duì)23個(gè)點(diǎn)目標(biāo)仿真產(chǎn)生的數(shù)據(jù)進(jìn)行二維成像仿真實(shí)驗(yàn),仿真參數(shù)設(shè)置如表1所示。目標(biāo)散射點(diǎn)模型如圖3所示。將圖4中RD算法的成像圖作為真實(shí)值,對(duì)本文成像算法性能進(jìn)行評(píng)價(jià)。為便于分析,將距離向采樣率與方位向采樣率同時(shí)設(shè)為原采樣率的1/2。仿真結(jié)果如圖5所示,可見除基于L1算法成像結(jié)果中有少量假目標(biāo)點(diǎn)外,其余算法均可清晰實(shí)現(xiàn)目標(biāo)點(diǎn)SAR成像。 表1 雷達(dá)仿真參數(shù) 圖3 目標(biāo)散射點(diǎn)模型 圖4 RD算法SAR圖像 圖5 各種算法二維成像結(jié)果 為進(jìn)一步比較上述算法的成像質(zhì)量,利用目標(biāo)雜波比和圖像熵對(duì)成像質(zhì)量進(jìn)行定量評(píng)價(jià)。表2給出了在采樣率為原采樣率1/2時(shí)幾種算法的TCR和IE值,顯然基于AWMC-TV算法的圖像TCR值大于其他算法,IE值小于其他算法,相比于其他算法有更好的聚焦性能和重建精度。 表2 目標(biāo)的成像結(jié)果評(píng)價(jià)參數(shù) 采用相對(duì)均方根誤差作為評(píng)估標(biāo)準(zhǔn),反應(yīng)重建圖像中每個(gè)散射點(diǎn)的信號(hào)強(qiáng)度與真實(shí)值間誤差,該值越小說明重建效果越好。表3顯示了本文算法相比其他算法RRMSE較小,重構(gòu)效果更好。 表3 幾種算法RRMSE 分別對(duì)各算法進(jìn)行計(jì)算復(fù)雜度對(duì)比,如表4所示。表4中,K1和K2分別表示各算法的循環(huán)迭代次數(shù),Ti和Ji(i=1,2,3)分別表示各算法的外部和內(nèi)部迭代次數(shù)。 表4 各類算法的計(jì)算復(fù)雜度比較 綜上,在相同的實(shí)驗(yàn)條件下,相比其他幾種算法,本文算法雖然復(fù)雜度有所增加,但目標(biāo)的重建精度顯著提高,表明該方法在SAR成像領(lǐng)域具有一定的潛力和優(yōu)勢(shì)。 為充分驗(yàn)證本文算法在實(shí)際應(yīng)用中的有效性,采用加拿大溫哥華的RADARSAT-1精細(xì)模式2(加拿大航天局版權(quán)所有)的實(shí)測(cè)數(shù)據(jù)進(jìn)行實(shí)驗(yàn),參數(shù)設(shè)置如表5所示。圖6(a)為基于RD算法的英國灣附近溫哥華地區(qū)的實(shí)際場(chǎng)景,選取紅色矩形區(qū)域內(nèi)的6個(gè)目標(biāo)點(diǎn)(6艘貨船)進(jìn)行實(shí)驗(yàn)。這6個(gè)目標(biāo)點(diǎn)相比海平面是強(qiáng)散射點(diǎn),相比英吉利海灣是稀疏的。圖6(b)~(f)為幾種算法在距離采樣率為原采樣率的30%、方位采樣率為原采樣率20%的成像結(jié)果。顯然,在低采樣率情況下,AWMC-TV算法依然能夠準(zhǔn)確重建目標(biāo)。 表5 溫哥華場(chǎng)景RADARSAT-1參數(shù) 圖6 各算法成像結(jié)果 為比較幾種算法在有噪聲的真實(shí)情況下的性能,對(duì)圖6(a)紅框中綠色標(biāo)記的點(diǎn)目標(biāo)沿距離向進(jìn)行切片,計(jì)算該點(diǎn)目標(biāo)幅值在所在距離單元內(nèi)的比值。表6給出了該目標(biāo)點(diǎn)振幅在各算法下占該距離單元內(nèi)的百分比,可見AWMC-TV算法下目標(biāo)點(diǎn)幅值在該距離向上占比最大,表明AWMC-TV算法能夠消除圖像中的模糊和拖尾,在重建性能和抗噪性能上具有優(yōu)越性。 表6 目標(biāo)點(diǎn)幅度在該距離單元內(nèi)占比 針對(duì)L1正則化存在的偏差效應(yīng)及傳統(tǒng)壓縮感知重構(gòu)算法重構(gòu)精度低的問題,本文將隸屬于非凸函數(shù)簇中的極小極大凹罰函數(shù)應(yīng)用到SAR成像中,同時(shí)將非凸正則化與全變分判罰線性組合作為成像模型的復(fù)合正則化器,提出一種基于自適應(yīng)加權(quán)的極小極大凹罰函數(shù)和全變分的稀疏SAR成像算法,并采用FISTA算法與交替最小化算法求解該模型。仿真與實(shí)測(cè)數(shù)據(jù)成像結(jié)果驗(yàn)證了該算法的有效性,在抗噪性能及重構(gòu)精度等方面比其他算法更具優(yōu)勢(shì)。 由于本文中壓縮感知模型只考慮到噪聲的影響,因此后續(xù)還需研究雜波對(duì)算法的影響。此外,能否通過多種變換找到復(fù)雜場(chǎng)景的稀疏表示是壓縮感知理論進(jìn)一步在SAR成像領(lǐng)域應(yīng)用的關(guān)鍵,這也是下一步研究的內(nèi)容。1 SAR回波信號(hào)二維稀疏表示模型
2 基于TV和非凸正則化的合成孔徑雷達(dá)成像
2.1 基于MC函數(shù)的稀疏SAR成像
2.2 基于MC函數(shù)和TV的稀疏SAR成像
2.3 基于自適應(yīng)加權(quán)MC函數(shù)和TV的稀疏SAR成像
3 仿真與實(shí)測(cè)數(shù)據(jù)結(jié)果與分析
3.1 目標(biāo)二維成像仿真實(shí)驗(yàn)分析
3.2 實(shí)測(cè)數(shù)據(jù)處理
4 結(jié)束語