基于小波包算法的壓縮傳感SAR成像方法

時 燕 陳迪榮

(北京航空航天大學數學與系統科學學院 北京 100191)

1 引言

合成孔徑雷達(Synthetic Aperture Radar,SAR)具有全天候、全天時、高分辨率等特點，目前已經得到了廣泛的應用。合成孔徑雷達的高分辨率必然帶來高數據率和高數據量，這對數據傳送帶寬和數據存儲帶來較大壓力。因此，減少雷達數據的數據量是很有必要的。

近年來，Donoho,Emmanuel Candès等人[1－5]提出了壓縮傳感理論，可以利用低于Nyquist采樣率的數據精確重構稀疏信號。近年來已引起圖像處理、雷達應用等方面的關注及應用。Baraniuk等人[6]首次提出將壓縮傳感應用于雷達成像中，分析了其應用的可行性。將壓縮傳感應用于雷達成像中能夠有效地減少數據量，但壓縮傳感應用的前提是稀疏信號(或有稀疏表示的信號)。文獻[7]提出基于小波稀疏表示的壓縮感知 SAR成像算法(Compressive Sensing Imaging Approach Based on Wavelet Sparse Representation,WCS)，該算法利用小波的稀疏表示功能，利用壓縮傳感模型進行SAR成像。該算法在嚴重降采樣情況下對具有較多細節的 SAR場景成像效果不是很好。這很可能源于小波分解只對低頻部分分解，而不對高頻部分分解。而小波包分解是對低頻和高頻部分都進行分解。所以小波包分解可能可以更好地表達較復雜的SAR場景。

本文提出基于小波包算法的壓縮傳感 SAR成像方法，通過對 SAR圖像進行分類，對多幅同類SAR圖像進行小波包分解，獲得一個綜合小波包樹，對綜合小波包樹應用小波包算法，得到稀疏表示該類SAR場景的小波包最優基，并通過求解l1范數最小化問題，得到SAR圖像。仿真實驗結果表明，在嚴重降采樣下，本文方法仍能夠得到效果較好的SAR圖像。本文的算法雖然比WCS算法的計算復雜度大一些，但在降采樣情況下成像效果好于WCS算法。

2 傳統的SAR成像原理

SAR成像就是從回波信號中提取地表后向散射系數的2維分布。SAR回波信號可以看作是發射脈沖與地表的后向散射系數卷積的結果，成像處理可看成一個解卷積的過程。

回波信號的模型為

其中σk為散射點k的后向散射系數，N表示散射點的個數。f0為線性調頻脈沖信號的載頻，rk(η)為散射點k在脈沖發射時刻η與雷達平臺間的距離，c為光速，Tp為脈沖重復間隔，Kr為線性調頻脈沖信號的時間寬度，nx(τ)為復基帶回波信號中由弱散射中心回波信號合成的等效加性噪聲。

傳統的 SAR成像算法主要是應用匹配濾波相位相乘等來解卷積，比如距離多普勒(Range Doppler,RD)算法[8]。

3 基于小波包最優基稀疏表示的壓縮傳感SAR成像

3.1 壓縮傳感SAR成像原理

將回波信號的卷積形式寫成如下矩陣形式：

式(4)中Φ是NrNa×N的chirp矩陣(Nr表示距離向采樣數，Na表示方位向采樣數)，其顯示表示如下：

其中

SAR成像的過程可以看作一個線性系統，通過回波信號sR和 chirp矩陣Φ來恢復場景反射系數σ。

稀疏度的定義：若向量s中只有K個元素非零，其它元素均為零，則稱向量s的稀疏度為K。

稀疏性的定義：信號中大部系數能量較小，只有少量系數的能量較大。稀疏度越小，稀疏性越大。

若σ的稀疏度為K，當K足夠小時，且l0范數最小化問題式(5)具有唯一最優解，則恢復σ就是求解l0范數最小化問題。表示零范數，表示向量σ中非零元素的個數。

l0范數最小化問題需要窮舉σ中所有非零值的所有種排列可能，是一個NP難問題。

Candès等人[2]指出：如果要精確重構稀疏度為K的信號σ，測量矩陣Φ必須滿足2K階限制等距(Restricted Isometry Property,RIP)條件[4,5]。當測量矩陣滿足RIP性質時，l0范數最小化問題式(5)可以與l1范數最小化問題式(6)等價。l1范數最小化問題是一個凸優化問題。

若σ在某個正交基Ψ下的展開系數為x,x具有很好的稀疏性(即稀疏度小，只有少量元素為零)。則上述問題可轉化為優化問題式(7)。

根據壓縮傳感理論，大多數隨機矩陣以高概率滿足RIP條件，可以精確重構原稀疏信號。RIP性質并不容易驗證，通常用不相干性來間接驗證。

文獻[6]中提到雷達發射的 chirp波形經過時間平移和頻率調制而形成的字典與場景反射系數的稀疏表示基(時間、頻率、時頻基)是不相干的。文獻[9]通過分析 chirp矩陣和隨機高斯矩陣的性質，表明chirp矩陣與隨機高斯矩陣的M-RIP性質類似。文獻[10]也通過大量試驗證明，chirp矩陣可以以高概率重構K稀疏信號。所以chirp矩陣Φ可以以高概率恢復滿足稀疏條件的場景反射系數σ。

當場景反射系數σ具有較好的稀疏性時，可直接用 chirp矩陣Φ來恢復場景。當場景較為復雜，細節較多時，場景反射系數不具有滿足重構條件的稀疏性時，需要一個能夠稀疏表示該復雜場景的正交基Ψ使之滿足重構條件，應用ΦΨ來恢復場景。

3.2 小波包算法選擇稀疏表示基

正交小波變換只對信號的低頻部分做分解，而對高頻部分，即信號的細節部分不再繼續分解，所以小波變換能夠表征一大類以低頻信息為主要成分的信號，但它不能很好地分解和表示包含大量細節信息的信號，如遙感圖像。與小波變換不同的是，小波包變換可以對高頻部分提供更精細的分解，對包含大量中、高頻信息的信號能夠更好的時頻局部化分析[11]。SAR場景主要表現地物地貌，包含大量細節信息。所以小波包能夠很好地表示包含大量細節信息的SAR復雜場景。

小波包分解樹是一個二叉樹(如圖1)，其中小波包樹中的子二叉樹是一組小波包基。一個小波包樹相當于一個小波包庫。壓縮傳感需要的是稀疏的信號，需要從中選擇能夠稀疏表示的小波包基。

圖1 小波包樹(3層小波包分解)Fig.1 Wavelet package decomposition tree (3-layer decomposition)

本文對SAR場景圖像進行分類，將具有類似特征的場景圖像歸類(根據圖像內容分類，例如陸地類SAR圖像、海岸類SAR圖像等)。用小波包算法找出稀疏表示基。主要步驟如下：

第1步：選擇L張同類場景圖像作為樣本，圖像展開為1維向量，長度為N。對每個樣本進行J層小波包分解，得到L個小波包分解樹；

第2步：計算每個樣本的小波包樹的第k層，第m個節點 β=(m,k),0≤m＜2k,0≤k≤J的系數的稀疏度；

第 4步：應用小波包算法[11]對綜合小波包樹以平均稀疏度sβ作為信息代價函數值，求出最優樹。具體算法為：從綜合小波包樹的最底層開始，標記最底層的節點，將兩個孩子節點的稀疏度之和與它們的父節點的稀疏度進行比較。如果比父節點的稀疏度大，則標記父節點，去掉孩子節點的標記。反之，用兩個孩子節點的稀疏度之和代替父節點的稀疏度，依次自底向頂。最終做標記的節點即為最優樹，最優樹對應的即是小波包最優基。

本文選擇了多張128×128陸地類SAR圖像作為樣本，由于db3小波的壓縮性能較好，所以本文選擇db3小波3層小波包分解，根據上述算法，得到小波包最優樹，如圖2所示。該樹對應的該類的小波包最優基。仿真用SAR場景圖像如圖3所示。

圖2 陸地類小波包最優樹Fig.2 Best wavelet package tree for land scene

圖3 SAR場景圖像(tier)Fig.3 SAR scene (tier)

圖4是SAR場景圖像(圖3)的灰度值的直方圖，該場景的圖像的灰度值相對較為平均，稀疏度大，稀疏性小。該場景在小波基下(如圖5)和小波包最優基下(如圖6)的系數都主要集中在0附近，稀疏度小，稀疏性較強。從圖中可以看出，小波包最優基下的系數分布更“陡”一些，系數更集中，稀疏性更好。

從表1中的數據也可以看出，小波系數和小波包稀疏的分布都集中在(0,0.1]范圍內，在(0,0.1]范圍內，小波包最優基下的系數分布比小波系數多，而在(0.1,1.0]范圍內，小波包最優基下的系數分布比小波系數少。同樣體現了小波包最優基下的系數比小波基下的系數更集中，稀疏度小，稀疏性更好。

圖4 SAR場景圖像的灰度值直方圖Fig.4 Gray histogram for SAR scene

圖5 小波系數分布圖Fig.5 Distribution of the wavelet coefficients for the SAR scene

圖6 小波包最優基下的系數分布圖Fig.6 Distribution of the wavelet package best basis coefficients for the SAR scene

表1 tier場景灰度值、小波系數、小波包最優基系數的稀疏度對比(%)Tab.1 Sparsity Compared from scene gray,wavelet coefficients and wavelet package best basis coefficients for tier scene (%)

3.3 小波包最優基壓縮傳感SAR成像

小波包最優基壓縮傳感(Compressive Sensing Based on Wavelet Package Best Basis,WPCS)SAR成像的過程為：將同類型SAR圖像根據3.2節中的算法訓練得到的小波包最優基作為該類 SAR圖像的稀疏表示基Ψ，以ΦΨ為觀測矩陣，回波信號sR為觀測值，通過壓縮傳感恢復算法求得場景σ。

WPCS是以小波包最優基作為稀疏表示基Ψ,WCS是以小波基作為稀疏表示基Ψ。雖然小波包最優基的計算量要比直接小波基的計算量大，但是，小波包最優基是對一類圖像通用的，只要對每類圖像訓練一次將該類的小波包最優基存儲于計算機內，在成像時，直接使用存儲的合適的最優基即可，并不比使用小波基復雜。

4 仿真結果

本文通過對SAR圖像用chirp信號進行該場景的回波仿真，用仿真回波數據進行成像處理，比較RD算法、WCS算法[7]、WPCS算法的成像效果。表2是仿真雷達參數。本文中的降采樣是對方位向進行降采樣，距離向采樣不變。方位向采樣率決定方位向分辨率。

圖7(a)、圖7(b)、圖7(c)分別為距離多普勒算法全采樣、50%采樣、30%采樣處理回波數據的結果。圖7(d)、圖7(e)、圖7(f)為WCS全采樣、50%采樣、30%采樣處理成像結果。圖 7(g)、圖 7(h)、圖 7(i)是WPCS全采樣、50%采樣、30%采樣處理成像的結果。從圖 7中可以看出，RD算法的成像結果圖像質量依次下降，30%采樣圖像發生嚴重混疊，已模糊不清。可以很明顯地看出WCS和WPCS在3種采樣情況下均比RD算法成像的效果好。

表2 仿真雷達參數Tab.2 The parameters of radar simulation

在50%采樣情況下，RD算法成像有輕微混疊，RD算法圖像質量下降較為明顯，WCS算法和WPCS算法的圖像質量基本沒有下降，幾乎不影響視覺效果。30%采樣的RD算法成像由于混疊嚴重，已經很模糊了，而WCS和WPCS成像，仍能保持圖像的主要特征，無模糊成像。從圖像中可以看出30%采樣的成像結果中，比較圖7(f)和圖7(i)，兩幅圖中從紅色方框內的細節中可以看出，圖 7(f)中紅色方框內的細節已經模糊，看不清物體的輪廓了，圖 7(i)中紅色方框內仍能看出物體的大致輪廓。圖7(i)比圖7(f)細節保留得更完整。

表3是圖7中RD成像、WCS成像，WPCS成像 3種成像結果的相關度、信噪比(SNR)、峰值信噪比(PSNR)。

圖7 SAR圖像(tier)成像處理(陸地類)Fig.7 SAR scene (tier) imaging processing (land scene)

表3 圖7中各種成像算法的相關度、SNR和PSNRTab.3 The correlation,SNR and PSNR of various imaging algorithms in Fig.7

最小化平均平方誤差(Mean Squared Error,MSE)的定義如下：

式中：xi為原始圖像的第 i個像素；為處理后的圖像的第i個像素。

信噪比(Signal to Noise Ratio,SNR)的定義如下：

信噪比的單位是dB。

圖像的峰值信噪比(PSNR)的定義：該指標一般應用在灰度在 0～255的單色圖像的降噪評估中。

從表3中可以看出，全采樣情況下，WPCS算法的相關度(0.8686),SNR(5.7422),PSNR(22.0874)和 WCS算法的相關度(0.8676),SNR(5.7289),PSNR(22.0741)均比 RD 算法的相關度(0.8153),SNR(4.4745),PSNR(20.8197)高；WPCS算法的相關度比RD算法的高6.54%，比WCS算法高0.12%；WPCS算法的SNR,PSNR比RD算法的高28.33%,6.09%；比WCS算法高0.23%,0.06%。50%采樣情況下，WPCS算法的相關度比RD算法的高6.35%,比 WCS算法的高 4.12%；WPCS算法的 SNR,PSNR比RD算法的高22.31%,4.76%；比WCS算法高18.83%,4.11%。30%采樣情況下，WPCS算法的相關度比RD算法的高8.64%，比WCS算法高5.89%，WPCS算法的SNR,PSNR比RD算法的高36.24%,7.07%，比WCS算法的高25.23%,5.26%。

圖8是海岸類SAR圖像成像處理的仿真結果。對于海岸類圖像經過 3.2節中的算法得到的小波包最優基與陸地類小波包最優基一樣，是如圖2所示的小波包最優基樹。圖8(a)-圖8(c)成像質量依次下降，50%采樣相干斑較多，海上目標已經被混淆了；30%采樣時混疊嚴重，輪廓模糊。圖8(d)-圖8(f)和圖8(g)-圖8(i)圖像的大致輪廓都能夠很好地保持，但圖 8(f)中存在一塊較大斑點(紅色方框中)，而圖8(i)圖像的質量效果要好于圖8(f)。

在降采樣情況下 WPCS算法的成像效果好于WCS算法。

5 結論

本文通過對不同場景類型的 SAR圖像在小波基下和小波包最優基下的稀疏度進行分析比較，小波包最優基下的系數比小波基下的系數更具有稀疏性。

本文通過對SAR圖像仿真成像，分別用傳統的成像算法(距離多普勒算法，RD算法)、WCS算法、WPCS算法進行仿真成像。成像結果表明WCS算法和WPCS算法比RD算法成出的圖像質量高。而在降采樣情況下，WPCS算法比WCS算法成出的圖像效果更好一點。對于較復雜的場景，WPCS算法的優勢更明顯一些。

[1]Donoho D L.Compressed sensing[J].IEEE Transactions on Information Theory,2006,52(4): 1289-1306.

[2]Candès E,Romberg J,and Tao T.Robust uncertainty principles: exact signal reconstruction from highly incomplete frequency information[J].IEEE Transactions on Information Theory,2006,52(2): 489-509.

[3]Candès E.Compressive sampling[C].In: Proceedings of International Congress of Mathematicians,Zürich,Switzerland,European Mathematical Society Publishing House,2006:1433-1452.

[4]Candès E and Tao T.Near optimal signal recovery from random projections: universal encoding strategies?[J].IEEE Transactions on Information Theory,2006,52(12):5406-5425.

[5]Candès E and Tao T.Decoding by linear programming[J].IEEE Transactions on Information Theory,2005,51(12):4203-4215.

[6]Baraniukc Richard C and Steeghs Philippe.Compressive radar imaging[C].IEEE Radar Conference,Waltham,MA,April 17-20,2007: 128-133.

[7]王偉偉,廖桂生,吳孫勇,等.基于小波稀疏表示的壓縮感知SAR成像算法研究[J].電子與信息學報,2011,33(6):1440-1446.Wang Wei-wei,Liao Gui-sheng,Wu Sun-yong,et al..A compressive sensing imaging approach based on wavelet sparse representation[J].Journal of Electronics & Information Technology,2011,33(6): 1440-1446.

[8]Ian G Cumming,Frank H Wong,洪文,et al..合成孔徑雷達成像—算法與實現[M].北京: 電子工業出版社,2007:154-169.Ian G Cumming,Frank H Wong,Hong W,et al..Digital Processing of Synthetic Aperture Radar Data: Algorithms and Implementation[M].Beijing: Electronic Industry Press,2007: 154-169.

[9]Lorne Applebaum,Stephen D Howard,Stephen Searle,et al..Chirp sensing codes: deterministic compressed sensing measurements for fast recovery[J].Applied and Computational Harmonic Analysis,2009,26(2): 283-290.

[10]Tropp J,Wakin M,Duarte M,et al..Random filters for compressive sampling and reconstruction[C].Proc.IEEE ICASSP,2005.

[11]孫延奎.小波分析及其應用[M].北京: 機械工業出版社,2005:245-259.Sun Yan-kui.Wavelet Analysis and Applications[M].Beijing:Machinery Industry Press,2005: 245-259.