集成InSAR與高分影像的地質災害隱患早期識別與應用研究

張文龍，張童康

(1.自然資源陜西衛星技術應用中心，陜西 西安 710100； 2.中國煤炭地質總局航測遙感局，陜西 西安 710100)

我國山地丘陵區約占國土面積的65%，地質條件復雜，是世界上地質災害最嚴重、受威脅人口最多的國家之一[1-2]。近年來，我國地質災害處于多發態勢，主要的地質災害有崩塌[3]、滑坡[4]、泥石流和地面塌陷[5]。地質災害多發地區主要有陜北黃土高原地區[6]、山西采煤沉陷區[7]、三峽庫區、滇西北和藏東南等地區。地質災害防治形勢嚴峻，給人民生命安全和經濟帶來嚴重威脅[8]。目前，利用InSAR技術對地質災害早期識別的研究取得了一些成功范例。張路等[9]采用相干散射體時序InSAR方法，識別出了17處持續變形中的不穩定坡體；戴可人等[10]利用InSAR成功對潛在的滑坡進行了早期識別；Grabriel等[11]利用D-InSAR技術探測厘米級的地表形變，測量了美國加利福尼亞州東南部的英佩瑞爾河谷灌溉區的地表形變；黃潔慧等[12]采用二軌D-InSAR技術對米林滑坡進行形變監測，得到米林滑坡災前形變趨勢；張曉倫等[13]將D-InSAR技術和光學遙感技術結合對東川城區周邊滑坡災害進行了有效識別。傳統的D-InSAR技術能夠有效探測地表微小變形，但是容易受空間和時間失相干的影響，大多數干涉圖中都會存在大區域信號不相關或沒有變形信息的現象，而且在最終的變形相位中往往包含明顯的大氣延遲相位，嚴重影響地表變形測量的精度和效率[14]。為了剔除這些誤差，提高監測精度，國外學者Ferreti等[15-16]提出了PS-InSAR和SBAS-InSAR技術為代表的時間序列InSAR技術。

代志宏等[14]采用PS-InSAR技術對南寧地表沉降進行了監測與分析，該研究可為地質災害預測和防治工作提供新思路。牛全福等[17]將PS-InSAR技術成功應用于蘭州市地面形變的監測。馮文凱等[18]利用SBAS-InSAR技術對金沙江流域沃達村滑坡進行地表形變監測，為類似老滑坡監測預警提供了新的思路與借鑒。Liu Peng等[19]利用SBAS-InSAR技術處理了三峽庫區2003—2010年間EnvisatSAR圖像，識別出巴東縣2個山體滑坡，認為季節滑坡運動與水位變化之間存在明顯的相關性。張詩茄等[20]嘗試采用小基線方法對一定區域范圍內的歷史滑坡進行監測分析，為區域監測提供了有益借鑒。可見，SBAS-InSAR方法不僅可以在區域尺度上獲取長時序的緩慢地表形變信息，而且在潛在滑坡識別、已知滑坡監測等方面具有較大的應用潛力[21]，但是其在地質災害隱患大范圍早期識別中應用較少。本文采用SBAS-InSAR技術結合高分影像對大范圍地質災害易發區進行早期識別，通過SBAS-InSAR技術提取形變區，然后利用高分影像進行遙感解譯，剔除不可靠的形變區，有效提高了地質災害早期隱患識別的精度。

1 SBAS-InSAR技術

SBAS-InSAR技術主要去除或減弱了時空失相干等誤差的影響，假設在時間段內，獲取S景SAR影像，選取主影像并與從影像進行配準。然后進行差分干涉處理，一共得到N幅差分干涉圖，且滿足：

(1)

對于tA和(tA

(2)

假定不同干涉圖間的形變速率為υk，k+1，則tA～tB間的累積形變可表示為：

(3)

對N幅干涉條紋圖進行三維時空相位解纏，即可求出不同SAR獲取時間的形變速率。

時間序列InSAR技術流程[23]如圖1所示。

圖1 SBAS-InSAR數據處理流程Fig.1 Data processing flow of SBAS-InSAR

2 研究區概況與實驗數據

2.1 研究區概況

子長市位于黃土高原中部，總面積為2 405 km2[24]，具體位置如圖2所示。子長市屬典型的黃土高原丘陵溝壑區，基巖產狀近于水平，無明顯褶皺和斷裂，而節理裂隙構造較為發育。境內峁梁起伏，溝壑縱橫，地面支離破碎，崩塌、滑坡等地質災害頻繁發生[25]。

圖2 研究區位置Fig.2 Location of study area

2.2 實驗數據

研究采用Sentinel-1數據作為InSAR監測主要數據源，研究區范圍為Path/frame：11/116，共獲取39景Sentinel-1數據，部分數據參數如下。

數據獲取日期共有10日，激化方式VV，入射角39.02°，分辨率20 m，拍攝模式為升軌。

在形變監測中，為了去除干涉相位中的地形相位，在處理過程中還需要高精度的DEM數據作為輔助數據進行配準、去除地形相位等操作。本次選用最新的AW3D DEM作為輔助數據。

InSAR數據處理中，在影像數據的圖像配準、干涉圖生成、相位解纏、形變圖生成等步驟中均涉及到軌道信息。此次工作利用ENVI中的SARscape拓展模塊自動讀取精密軌道數據，并添加到干涉處理過程中。

此次遙感監測工作主要采用國產高分系列(GF-1、GF-2)衛星數據，研究區涉及GF-2數據和少量GF1數據共19景，覆蓋解譯區全區。

3 結果分析

3.1 時空基線分析

將軌道Path/Frame為11/116的Sentinel-1數據源進行格式轉化，采用時序處理方法中的小基線集技術進行時間序列InSAR技術處理。時間基線閾值設置為90 d，空間基線閾值設置為5%，構建干涉像對組合，用于后續的差分干涉處理。

為提高解纏的精確性，選取可用于3D解纏的連接圖。選取39景Sentinel-1數據共計134個干涉像對組合，干涉圖組合圖如圖3—圖5所示。

圖3 時空基線連接Fig.3 Space-time baseline connection

圖4 時間和位置的連接Fig.4 Connection of time and position

圖5 3D解纏連接Fig.5 3D untangled connection

3.2 差分干涉圖

對差分干涉圖進行空間濾波和相位解纏，獲取差分干涉圖(圖6)。研究中，為避免基線誤差的影響，采用GCP點結合地面高程值和對應的解纏干涉相位進行基線優化，優化后再次進行差分干涉、濾波和相位解纏，獲取解纏后的差分干涉相位。本次選取其中一景延安市子長市瓦窯堡鎮(Path/frame：11/116)差分干涉圖和相干系數圖，如圖6所示(雷達坐標系下)。

圖6 SAR干涉對查分結果Fig.6 SAR interference matching results

從差分干涉圖上能夠直接發現由形變引起的干涉條紋，見圖6中的虛線框位置，同時在相干系數圖(圖6 (c))可以看出這些形變區域明顯。為進一步分析潛在的形變區域，將所有建立連接的影像對進行差分干涉處理，選取干涉質量較優的差分干涉像對，高相干點識別方法選擇時序相干均值大于0.35相干目標點進行SBAS-InSAR反演，對去除差分干涉相位進行時間序列分析，對殘差相位采用空間域濾波方法進行大氣延遲相位提取，獲取大氣延遲相位后，再次進行差分干涉得到差分干涉圖。此時仍然存在恒定相位和越變相位，進行軌道精煉和重去平，消除上述相位引起的誤差，保證干涉圖質量，提高SBAS-InSAR反演結果的精度。對比時序SAR像對獲取的相干系數和解纏結果，找到多數干涉圖中相干性高的區域并在此選擇GCP點，進行重去平處理。選擇線性模型，通過Goldstein濾波消除大氣效應引起的相位誤差，得到最終的時序形變量，并將得到的平均形變速率圖及時序形變量進行地理編碼，轉換到地理坐標系下。

3.3 SBAS-InSAR監測結果分析

基于SARscape軟件，進行了基線估算、配準，Sentinel-1數據配準精度均優于0.001個像元，滿足干涉的要求；通過差分干涉、濾波、解纏等，利用39景Sentinel-1數據獲取了2018年1月—2020年2月陜北子長地區的地表形變信息(沿雷達視線向，即LOS向)。監測結果包括平均形變速率、形變區分布、累計形變量等信息(圖7、圖8、圖9)。利用視坡夾角將沿雷達視線方向的形變速率轉化成沿斜坡坡度方向的形變速率，以表達沿斜坡面的地表形變信息，形變速率圖和累計形變量圖中紅色代表遠離衛星視線方向，即代表沉降；藍色代表靠近衛星視線方向，即代表抬升。

圖7 地表形變速率Fig.7 Surface deformation rate

圖8 地表累計形變Fig.8 Cumulative deformation diagramof surface deformation

圖9 形變區分布Fig.9 Distribution of deformation areas

由形變速率圖和累計形變圖可知，本次監測子長地表最大形變速率為-76 mm/a，累計最大形變量為-152 mm。形變區域主要分布在子長市瓦窯堡鎮和余家坪鎮，其中子長市南部的形變區域相較于北部較多，主要集中于南部的瓦窯堡鎮、余家坪鎮和楊家園子鎮；北部的澗峪岔鎮和南溝岔鎮的形變區域分布較少。

子長市地表形變主要為煤礦開采導致地表塌陷變形，其變形區形態在InSAR結果上較為規則，形態多呈長條狀或類似矩形，與煤礦開采區形態相關。

3.4 典型點分析

子長市新磘上滑坡(ZC1),累積形變如圖10所示。

圖10 滑坡點ZC1的累積形變(沿LOS向形變)Fig.10 Cumulative deformation of ZC1 at landslidepoint (deformation along LOS direction)

從圖10中可以看出，滑坡點ZC1有明顯的形變信號，結合光學影像和地形分析，此處滑坡在雷達視線的背坡向，其從高向低滑動，在InSAR監測結果中表現出負的形變(遠離衛星)。滑坡區域中某特征點的形變曲線顯示(圖11)，此處滑坡在2018年1月—2019年11月間不同時刻有明顯的形變趨勢。

圖11 ZC1滑坡區中某特征點的時序形變曲線Fig.1 Time-series deformation curve of a characteristic point in ZC1 landslide area

4 結論

基于SBAS-InSAR技術開展陜北子長地區地表形變InSAR監測，獲取了2018年1月—2020年2月子長地區的地表沉降信息(沿雷達視線向，即LOS向)。子長地表最大形變速率為-76 mm/a，累計最大形變量為-152 mm。形變區域主要分布在子長市瓦窯堡鎮和余家坪鎮。

通過基于InSAR地表形變監測、綜合遙感技術監測、野外調查等工作，排除偽隱患點，形成了對區域地質災害發育及嚴重程度的正確認識。初步總結了一套InSAR形變分類方法，由于地表形變的形成因素較多以及InSAR數據獲取、處理過程中的誤差影響，InSAR監測獲取的形變到地質災害隱患確定過程尤為重要，對形變點(區)進行分類，有益于地質災害隱患點的獲取，為隨后數據的處理優化也提供了借鑒。初步總結了一套適用于子長地區地質災害InSAR監測技術體系，本次基于InSAR技術開展區域性地質災害監測預警，綜合應用高分辨率多光譜光學影像遙感、合成孔徑雷達干涉測量、三維傾斜攝影測量等對地觀測技術，獲取區域性地質災害的分布特征、變形狀況，結合實地調查核實工作，圈定并監測重點地表形變區，確定地質災害隱患點，編制地質災害隱患分布圖。