SAR影像中疊掩與陰影區域的識別
——以湖北巴東為例

張同同，楊紅磊,2，李東明，李永杰，劉俊男

(1.中國地質大學(北京)，北京 100083； 2. 資源環境與災害監測山西省重點實驗室，山西 晉中 030600; 3. 水利部水工金屬結構質量檢驗測試中心，河南 鄭州 450044)

合成孔徑雷達干涉測量(interferometric synthetic aperture radar，InSAR)是一種新型主動式地表變形監測技術，廣泛應用于地質環境災害監測方面[1-3]，其通過微波成像的方式獲取影像，與可見光成像方式相比，具有全天候、無接觸、大面積、高空間分辨率、高精度等優勢[4-5]。

滑坡是坡體在自然因素或人為活動影響下無法維持自身的力學極限平衡狀態而發生的一種地質現象[6],滑坡一旦發生會造成巨大的人員傷亡和經濟損失，因此滑坡監測在國家經濟發展中十分重要。由于地形和雷達衛星運行軌跡的限制，在地形起伏較大的區域會使雷達波束從斜坡頂部返回的時間比從斜坡底部返回的時間短，導致在雷達影像中斜坡的頂部和底部顛倒顯示，即為疊掩現象[7-8]。而由于山體信號遮擋，導致雷達影像上陰影區域內沒有數據，即為陰影現象。這兩種現象都嚴重破壞了干涉相位的連續性，導致InSAR技術在相位處理工作中無法準確濾波和解纏。

近年來許多專家學者對疊掩和陰影區域的識別與提取做了相應的研究。文獻[9]利用SAR圖像中的頻譜搬移特性識別疊掩和陰影區域；文獻[10]針對單幅SAR影像中在瞬時照射到目標區域時會受到最小SAR天線面積的限制，使用多基線和多波段的方法解決疊掩問題；文獻[11]利用極大似然估計方法，融合相位圖與DEM提取疊掩和陰影區域；文獻[12]采用基于干涉信號相關矩陣的特征值提取疊掩和陰影區域；文獻[13]分析疊掩區域的干涉相位特性和幅度特性，提出了一種結合SAR幅度圖像和干涉相位對疊掩區域進行識別的算法；文獻[14]對機載SAR數據在相位域對陰影區域進行干涉處理，利用高質量DEM數據識別陰影區域；文獻[15]提出一種結合SAR圖像干涉幅度和相關系數來識別疊掩和陰影區域的算法；文獻[16]基于Tsai方法，運用連續局部閾值分割的方法對陰影區域進行檢測；文獻[17]提出一種基于局部頻率估計，幅度閾值分割以及結合形態學識別SAR影像中的疊掩和陰影區域的方法。本文利用組合坐標系統提取雷達監測過程中的幾何模型，并結合形態學方法實現疊掩和陰影區域的識別，取得了較好的試驗效果。

1 原 理

1.1 疊掩和陰影區域的檢測

根據雷達采集數據時的幾何模型獲得該時刻的瞬時局部入射角，試驗通過瞬時局部入射角判斷該區域是否為疊掩和陰影區域，如圖1所示。

局部入射角θ的具體表達式為

θ=β-α

(1)

式中，β為雷達波側視角；α為目標點的坡度角。當θ<0°時該目標區域產生頂底位移，為疊掩區域；當0°<θ<90°時該目標區域產生透視收縮；當θ>90°時該目標區域產生信號缺失，為陰影區域。

根據數學模型分析可得，在雷達影像成像時，根據信號傳輸的距離判定目標在SAR影像上的位置，進而判別是否為疊掩和陰影區域。疊掩和陰影區域分為主動區域和被動區域，如圖2所示，在檢測過程中將兩者區分檢測，首先對主動疊掩區域進行檢測，在獲得第一個主動疊掩點后，向近距端方向檢測被動疊掩區域，直到確定無被動疊掩區域。在檢測出一個完整的疊掩區域之后向遠距端方向尋找與疊掩區域初始位置距離一樣的位置，該位置與主動疊掩區域之間同樣為被動疊掩區域，同理可得，在檢測出主動陰影區域之后向該區域遠距端方向檢測被動陰影區域，直到無被動陰影區域出現。

1.2 坐標系統的組合

常用的坐標系為大地坐標系、高斯平面坐標系、以及WGS-84坐標系，在這幾種坐標系的基礎上，本文根據實際需要將高斯平面坐標系與大地高程坐標系構成組合坐標系統，獲取成像時的幾何模型。

試驗數據中衛星軌道坐標系為WGS-84坐標系，衛星軌道的參數文件中包含11組位置參數(X(t)，Y(t)，Z(t)，t)。試驗利用這些離散參數運用三項式擬合的方法，擬合出雷達衛星的運行軌跡，其擬合公式為

(2)

式中，a0、a1、a2、a3、b0、b1、b2、b3、c0、c1、c2、c3為待定參數。

衛星軌道坐標系為WGS-84坐標系，試驗需要統一坐標系，因此將衛星軌跡坐標(xyz)轉換到大地經緯度(BLH)坐標系下，其具體表達式為

(3)

式中，N為卯酉圈半徑，即

(4)

式中，a=6 378 137 m；e2=0.006 694 379 990 13。

試驗中大地坐標系為過渡坐標系，主要是將已知目標點的(B,L)轉換到高斯平面直角坐標(x,y)下，從而解算出雷達衛星與照射目標的幾何關系。

1.3 數學形態學

圖像腐蝕和膨脹算法是形態學圖像處理的基礎，眾多形態學算法也以這兩種操作為基礎[18]。

(1) 膨脹運算

x⊕B={x|X∩Bx≠φ}={x|Bx↑X}

(5)

式(5)表示數據集X用結構元素B進行膨脹時，其結果為集合x，是BX與X的交集不為空的數據集，或者說x是B集中(用↑號表示)而形成的數據集。

(2) 腐蝕運算

xΘB={x|Bx?X}

(6)

式(6)表示數據集X用結構元素B進行腐蝕時，其結果為集合x。

根據形態學開運算能夠去除孤立的小點和毛刺以及濾波靈活的特性，試驗選擇開運算去除圖像中微小錯誤雜點，即先進行腐蝕運算再進行膨脹運算。

2 數據處理與試驗分析

2.1 研究區概況

試驗區域為湖北省巴東地區，地形西高東低、南北起伏，且地質條件復雜。在新構造運動及庫區水位變化等外力作用下，使得該地區滑坡災害頻發。本文使用的試驗數據是ALOS-2衛星采集的單視復數SAR圖像和ALOS-1衛星生成的12.5 m分辨率的DEM數據。巴東地區高程介于600 m與1200 m之間，地形縱橫交錯，地貌環境復雜，使得SAR影像中容易產生疊掩和陰影現象，如圖3所示。

2.2 方法流程

本文試驗利用瞬時局部入射角的數學模型反演出SAR影像成像時的幾何模型，然后進行疊掩與陰影區域的提取，并結合形態學上的開運算去除誤差點，其具體流程如圖4所示。

2.3 試驗分析

根據強度圖像特性，在疊掩區域中強度值較大，圖像較亮；在陰影區域中強度值較小，圖像較暗。本文運用強度影像這一特性作試驗結果的對比性分析，其原始圖像如圖5所示。

由圖5可以看出，圖中疊掩和陰影區域有明顯亮度差異，因此試驗將所得結果疊加到強度圖上作為對比性試驗分析，其試驗結果如圖6所示。

根據圖5與圖6對比可知，本文試驗檢測的結果和強度圖像大致相同，但有部分微小地區檢測結果和強度圖像不一致。根據實際地形和高程數據分析主要有以下兩種原因：①有部分微小地區地形起伏較大，與周邊整體起伏變化有較大差異，而在強度圖像上信號值被被動地區掩蓋，不能被分辨出來；②試驗數據使用的是ALOS-1衛星獲取的DEM數據，時間距今相差較大，部分地區會有植被覆蓋或地形改變等因素發生。因此針對此類問題，試驗運用形態學上的開運算方法消除誤差點。檢測結果如圖7所示。

對比圖6與圖7可得，經過開運算處理的結果明顯去除錯誤的檢測區域，同時將檢測結果更加細化，完善了檢測結果的整體性，得到了良好的試驗效果。

從檢測結果可以看出，疊掩相較于陰影區域更容易產生。對比試驗幾何模型可知，試驗使用的數據日期為2017年10月29日，中心雷達波側視角為32.411 7°。通過粗略計算可得在12.5 m分辨率的DEM數據中所發生陰影現象必須保證相近兩像元之間的高差在20 m以上，而在相同情況下發生疊掩現象只需高差在8 m以上，因此產生兩種現象的高差相差巨大。依據本文試驗，在今后的雷達影像數據采集工作中根據研究區域地形地貌的不同設計相應的雷達入射角，可以提高數據采集的質量。

3 結 論

本文利用雷達衛星成像時的幾何模型與圖像處理中形態學方法綜合識別疊掩和陰影區域，并以巴東地區為例，驗證了本文試驗方法在疊掩和陰影識別領域中擁有很好的檢測能力。同時隨著DEM精度的更新與提高，會更加完善本文的試驗結果。