基于曲波變換的全極化SAR數據與TM影像融合方法

肖世忱，廖靜娟

(1.中國科學院遙感與數字地球研究所數字地球重點實驗室，北京 100094；2.中國科學院大學，北京 100049)

1 引 言

隨著傳感器技術、遙感技術等的迅猛發展，多平臺、多傳感器、多時相、多光譜和多分辨率的遙感影像構成了同一地區的多源遙感信息。在眾多的信息源中，SAR數據對目標的幾何特性非常敏感，而多光譜影像對目標的物理和化學屬性則更為敏感，主要反映不同地物的輪廓與光譜信息[1]，這些不同的信息類型涉及到完全不同的目標質量并且極大地不相關[2]，因此SAR數據與多光譜影像的融合將能充分利用它們的互補性，彌補各自的不足，這對于更好地實現對被觀測目標的理解是非常重要的。

SAR數據與多光譜影像進行融合的方法主要有IHS變換和小波變換等。其中基于IHS變換的經典融合方案是在IHS空間內，使用經過拉伸處理的單極化SAR數據替換多光譜影像的I分量的方式來產生融合影像[3]，該方式雖容易扭曲原始光譜特性，但可分離多光譜影像中的光譜信息和強度信息。小波變換因其具有良好的多尺度與時頻局部特性在影像融合領域也得到廣泛應用。近年來Candes等提出曲波變換的概念[4]，該變換不僅繼承了小波變換的優良特性，還具有多方向和各向異性的特征，在給定相同的重構精度下可以近乎最優地表示影像中具有方向性的直線及曲線等奇異性特征，對函數的稀疏表達能力較小波變換也更強。基于Curvelet變換的經典融合方案是在粗尺度層取系數平均值，細尺度層取絕對值較大的系數。許星等采用這種方案進行SAR數據與TM影像的融合時，獲得了相比于小波變換更好的融合效果[5]。

目前，研究者在進行SAR數據與多光譜影像的融合時，多采用單極化數據，主要因為SAR數據記錄的后向散射信息對地面目標的類型、方向、同質性和空間關系更為敏感，相比于光學遙感影像呈現出有望能被更好使用的結構特征[6]。相比于單極化SAR數據，多極化數據在具備目標豐富結構信息的同時，還記錄了目標在不同極化方式下的信息，尤其是全極化數據，以散射矩陣的方式記錄了地面目標的信息。隨著Radarsat-2等衛星的成功發射，全極化SAR數據逐漸豐富，其參與融合可使融合結果中包含全極化數據中蘊含的豐富地物信息。Xu等嘗試在RGB空間采用相應通道以一定的權重相加來實現SAR全極化數據與多光譜影像的融合[7]，但由于RGB顏色空間是根據紅、綠和藍3原色的混合比例來定義不同的顏色，而人眼只能通過感知顏色的色調、亮度和飽和度來識別和區分物體[8]，因此選擇在IHS空間進行融合將更加優越。

本文在IHS空間提出了一種基于第二代曲波變換融合全極化SAR數據與TM影像的方法。針對兩種數據源的影像特征，設計了一種充分保留原影像重要特征的低頻子帶系數融合規則。實驗結果表明，該方法所得融合影像較基于IHS變換的經典融合方案和基于曲波變換的經典融合方案所得影像在主觀視覺效果和客觀定量指標上均有所改善。

2 數據源與數據處理

2.1 數據源

SAR數據源是Radarsat-2衛星于2012年6月中旬在黑龍江省獲得的全極化數據。產品類型是單視復數據，空間分辨率為8m。TM影像是由Landsat-5衛星于2008年6月30日獲取的，產品等級為L4，包含7個波段，分辨率為30m。影像內大部分區域地勢平坦，包括居民區、耕地、河流和道路等典型地類。

2.2 數據處理

Radarsat-2全極化數據與TM影像融合前的數據處理主要涉及SAR數據的相干斑抑制、影像配準以及直方圖匹配等環節。由于TM影像的L4級產品是經過輻射校正、幾何校正和幾何精校正等處理后的高精度產品，因此數據處理環節中涉及到信息統一的操作將以TM影像為基準。

SAR數據的相干斑抑制步驟是先對SAR數據進行方位向：距離向等于2∶1的多視處理來解決方位向和距離向采樣大小差異較大的問題，再采用窗口大小為3×3的Gamma Map自適應濾波方法分別對4種極化數據進行降斑處理，處理后地物內部平滑均勻，地物邊緣銳化清晰。配準處理是以TM影像為基準來校正SAR數據，控制點的選取位于道路和田梗等的交叉點和拐點處，糾正方法選用的是一階多項式。配準后X、Y方向的配準誤差均控制在1個像元內，總體RMS為0.60個像元。由于從不同傳感器采集到的影像會因平臺角度、影像獲取時間以及衛星波長擴展間存在差異等造成相同的場景具有不同照明條件[9]，所以融合前還需將SAR數據各通道同TM影像相應波段進行直方圖匹配來調整SAR數據的亮度。

3 融合方法

曲波(Curvelet)變換理論可細分為兩代，其中第二代Curvelet變換[10]是一種更易于數字實現、更便于理解的快速變換，與第一代Curvelet變換[4]相比，在構造上已經完全不同，其在實現過程中沒有引入Ridgelet變換[11]，而是在頻域中直接給出了Curvelet基的具體表示形式，可以說是一種真正意義上的Curvelet變換。第二代Curvelet變換算法包括兩種具體實現：纏繞算法和非等間隔快速傅里葉算法[12]。相比于后者，纏繞算法的等間隔特性使其更容易實現且忠實于連續變換，再加上運算量和復雜度都相對較小[13]，因此本文將基于纏繞算法對數據進行Curvelet變換。SAR數據與TM影像融合方案如圖1所示，I分量Curvelet變換流程如圖2所示。

3.1 低頻子帶系數融合規則

影像經Curvelet變換后所得低頻子帶系數代表的是原影像的概貌信息，包含了影像的主要能量。由于SAR與TM影像的成像機理不同，同一目標在兩種影像中的灰度分布特性往往存在較大差異，甚至極性完全相反。因此，對低頻子帶系數沿用傳統加權平均的融合規則通常不能獲得滿意的結果。為此本文提出了一種基于SAR數據與TM影像特征的低頻子帶系數融合規則。

Radarsat-2數據的分辨率高于TM影像，在刻畫目標邊緣等細節信息方面的能力要優于TM影像，因此首先對SAR數據的低頻子帶系數求取梯度并作歸一化處理，然后適當選取梯度閾值TG，當SAR低頻子帶系數的梯度值≥TG時，融合系數F(x，y)取SAR與TM低頻子帶系數的線性加權。由于SAR數據中亮度較大的目標往往也是需要融合的目標，而TM影像在描述影像紋理特征方面要優于SAR數據，因此當SAR低頻子帶系數的梯度

影像f(x，y)在像素點(x，y)處的梯度可定義為一個向量：

(1)

梯度的模為各分量的平方和再求平方根，即為：

(2)

對于數字圖像，連續導數形式可以用求差來近似表示，則：

|G(x,y)|?|f(x,y)-f(x+1,y)|+|f(x,y)-f(x,y+1)|

(3)

低頻子帶系數的融合規則為：

(4)

3.2 高頻子帶系數融合規則

經Curvelet變換后所得高頻子帶系數描述的是影像的細節和邊緣等特征，而這些特征在曲波變換域中則表現為具有絕對值較大的高頻子帶系數。因此，高頻子帶系數采用取絕對值大的系數進行融合。

圖1 SAR與TM影像的融合方案

圖2 基于Curvelet變換的影像I分量融合流程圖

4 融合結果與效果的分析評價

融合結果的總體效果如圖3所示。由于試驗區較大，總體效果無法詳細描述細部信息，文中又選取道路和耕地結果放大圖進行具體說明，如圖4、圖5所示。每幅圖中又包括6景子圖像，其中(a)為原SAR影像，(b)為原TM影像，(c)為HV極化數據基于曲波變換經典融合方案的結果，(d)為全極化SAR數據基于IHS變換經典融合方案的結果，(e)為全極化SAR數據基于曲波變換經典融合方案的結果，(f)為文中方法的融合結果。

圖3 總體效果圖

圖4 道路效果圖

圖5 耕地效果圖

從目視判讀的效果上看，相比于影像(c)，影像(e)和(f)不僅影像內的信息更加豐富，而且影像效果更加平滑，這說明全極化SAR數據產生的I分量相比于SAR單極化通道不僅能有效綜合各極化通道的信息，還能進一步抑制單極化通道中的相干斑噪聲，從而產生更好的融合效果。對比影像(d)、(e)和(f)可以發現經典IHS變換所得結果(d)較原TM影像扭曲嚴重，光譜失真較大，而曲波變換的引入則能有效克服傳統IHS變換方式的不足，實現對原TM影像光譜信息更好的保留。對比道路、居民區等地物的光譜信息可以發現影像(f)在保持原TM影像的光譜信息以及影像清晰度方面要優于影像(e)，這說明文中方法的融合規則較經典曲波變換融合規則更加有效。

SAR數據與TM影像融合效果的定量化評價應綜合考慮對光譜信息的保持和空間細節信息的增強等兩方面。其中，反映光譜信息的指標包括扭曲程度和相關系數等，而反映空間細節信息的指標則可以用信息熵和清晰度等進行衡量。

扭曲程度直接反映了融合影像相比于原影像的光譜失真程度。扭曲程度值越小，融合影像的光譜失真程度就越少。其定義為：

(5)

式中DNF(i，j)與DN(i，j)分別為融合后和原始影像上(i，j)處的灰度值。

相關系數揭示了原始影像和融合影像之間在小尺寸結構上的相似性，反映了原始影像和融合影像的相關程度。相關系數越大，融合后的影像與原始影像的相關性越好。其定義為：

(6)

式中M與N分別代表影像的行與列，D(X)、D(Y)分別為影像X和Y的方差。

信息熵是衡量影像信息豐富程度的指標，是將影像內頻率的變化作為信息。融合影像的信息熵越大，說明融合影像的信息量增加的越多。其定義為：

(7)

式中pi為影像X中像素灰度值為i的概率。

清晰度又稱為平均梯度，是描述影像清晰程度的指標，反映了影像中微小細節反差與紋理變化特征，清晰度越大，影像越清晰，邊緣信息保留越多。

(8)

原SAR數據、原TM影像以及各方法融合結果相應指標的統計結果如表1所示。從信息熵指標上看，融合結果的信息量較原SAR數據和TM影像都有增加，這與目視判讀效果相一致，而方法間無明顯差別。從扭曲程度和相關系數兩項指標上看，融合結果較原TM影像均產生一定程度的扭曲和失相關。但在程度上文中方法所產生的扭曲和失相關均低于另外3種方法，這說明文中方法在對TM影像光譜信息的保持方面性能更優。從清晰度指標上看，融合結果的清晰度均優于原始數據。在基于全極化數據產生的融合結果中，文中方法較另外兩種融合方法所得融合結果的清晰度更高，這說明文中方法能更有效地在融合結果中注入全極化SAR數據的細節信息。所不足的是基于全極化數據的融合結果在影像清晰度方面要低于單極化SAR數據的融合結果，這主要因為全極化數據I分量的產生實質上是HH、HV和VV極化數據間的線性加權，而這種方式在一定程度上會降低原影像的對比度，從而導致融合結果影像清晰度的部分損失。從各方法執行時所耗費的時間看，IHS變換經典融合方案耗費的時間最少，原因在于沒有對數據進行曲波變換；文中方法所耗費的時間要略多于另外兩種基于曲波變換的融合方法，原因在于文中方法需要進一步計算低頻子代系數的梯度等信息。

表1 影像融合結果的定量評價

5 結束語

本文在IHS空間提出了一種基于二代Curvelet變換融合全極化SAR數據與TM影像的算法。針對全極化SAR數據與TM影像的特征，對Curvelet變換后的低頻子帶系數依據SAR數據的梯度信息與亮度信息確定其融合權值，對高頻子帶系數采用取絕對值大的系數進行融合。實驗結果表明，該方法既可以很好地保留TM影像中的光譜信息，又能有效地注入全極化SAR數據中的細節信息。在目視效果與定量化指標上均優于基于經典IHS變換和經典Curvelet變換的融合算法。但方法在綜合全極化SAR數據中的信息時存在損失影像清晰度的不足。另外，文中僅將全極化SAR數據當作彩色圖像使用，而全極化數據是可通過極化變換等方式進一步挖掘信息的，如何提取出全極化SAR數據中最有益于參與融合的I分量也是后續需要進一步研究的問題。

參考文獻：

[1] 李暉暉，郭雷，劉航.基于互補信息特征的SAR與可見光圖像融合研究[J].計算機科學，2006，33(04)：221-224.

[2] HELLWICH O，GUNZL M.Landuse classification by fusion of optical and multitemporal SAR imagery[C].Proceedings of the International Geoscience and Remote Sensing Symposium，Honolulu，Hawaii，July，2000：1921-1923.

[3] 黃登山，楊敏華，姚學恒，等.基于Contourlet與IHS變換的SAR與多光譜影像融合方法[J].宇航學報，2011，32(1)：187-192.

[4] CANDES E J，DONOHO D L.Curvelets：A surprisingly effective nonadaptive representation for objects with edges[C].Proceedings of International Conference on Curves and Surfaces，4th，Saint-Malo，France，July，1999：105-120.

[5] 許星，李映，孫瑾秋，等.基于Curvelet變換的SAR與TM圖像融合研究[J].西北工業大學學報，2008，26(03)：395-398.

[6] ZENG Y，ZHANG J X，VAN GENDEREN J L，et al.Image fusion for land cover change detection[J].International Journal of Image and Data Fusion，2010，1(2)：193-215.

[7] XU M S，XIA Z S，ZHANG F L，et al.Multi-temporal polarimetric SAR and optical data fusion for land cover mapping in southwest China[C].International Conference on Multimedia Technology，2010：1-4.

[8] 張微，丁獻文，陳漢林.基于改進的IHS方法的SAR圖像增強與融合[J].海洋學研究，2007，25(01)：73-79.

[9] GUNGOR O，SHAN J.An optimal fusion approach for optical and SAR images[C].Proceedings of the ISPRS Commission VII Mid-term Symposium “Remote Sensing：From Pixels to Processes”，Enschede，The Netherlands，May，2006：111-116.

[10] CANDES E J，DONOHO D L.New tight frames of curvelets and optimal representations of objects with piecewise C2 singularities[J].Communications on Pure and Applied Mathematics，2003，57(2)：219-266.

[11] CANDES E J.Ridgelets：Theory and applications[D].Stanford University，1998.

[12] CANDES E J，DEMANET D L，YING L.Fast discrete curvelet transforms[J].Multiscale Model & Simulation，2006，5(3)：861-899.

[13] SVEINSSON J R，BENEDIKTSSON J A.Combined wavelet and curvelet denoising of SAR images using TV segmentation[C].Proceedings of the International Geoscience and Remote Sensing Symposium，Barcelona，Spain，2007：503-506.