集成InSAR與高分影像的地質(zhì)災(zāi)害隱患早期識別與應(yīng)用研究

張文龍，張童康

(1.自然資源陜西衛(wèi)星技術(shù)應(yīng)用中心，陜西 西安 710100； 2.中國煤炭地質(zhì)總局航測遙感局，陜西 西安 710100)

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

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

1 SBAS-InSAR技術(shù)

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

(1)

對于tA和(tA

(2)

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

(3)

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

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

圖1 SBAS-InSAR數(shù)據(jù)處理流程Fig.1 Data processing flow of SBAS-InSAR

2 研究區(qū)概況與實驗數(shù)據(jù)

2.1 研究區(qū)概況

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

圖2 研究區(qū)位置Fig.2 Location of study area

2.2 實驗數(shù)據(jù)

研究采用Sentinel-1數(shù)據(jù)作為InSAR監(jiān)測主要數(shù)據(jù)源，研究區(qū)范圍為Path/frame：11/116，共獲取39景Sentinel-1數(shù)據(jù)，部分?jǐn)?shù)據(jù)參數(shù)如下。

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

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

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

此次遙感監(jiān)測工作主要采用國產(chǎn)高分系列(GF-1、GF-2)衛(wèi)星數(shù)據(jù)，研究區(qū)涉及GF-2數(shù)據(jù)和少量GF1數(shù)據(jù)共19景，覆蓋解譯區(qū)全區(qū)。

3 結(jié)果分析

3.1 時空基線分析

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

為提高解纏的精確性，選取可用于3D解纏的連接圖。選取39景Sentinel-1數(shù)據(jù)共計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 差分干涉圖

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

圖6 SAR干涉對查分結(jié)果Fig.6 SAR interference matching results

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

3.3 SBAS-InSAR監(jiān)測結(jié)果分析

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

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

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

圖9 形變區(qū)分布Fig.9 Distribution of deformation areas

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

子長市地表形變主要為煤礦開采導(dǎo)致地表塌陷變形，其變形區(qū)形態(tài)在InSAR結(jié)果上較為規(guī)則，形態(tài)多呈長條狀或類似矩形，與煤礦開采區(qū)形態(tài)相關(guān)。

3.4 典型點分析

子長市新磘上滑坡(ZC1),累積形變?nèi)鐖D10所示。

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

從圖10中可以看出，滑坡點ZC1有明顯的形變信號，結(jié)合光學(xué)影像和地形分析，此處滑坡在雷達(dá)視線的背坡向，其從高向低滑動，在InSAR監(jiān)測結(jié)果中表現(xiàn)出負(fù)的形變(遠(yuǎn)離衛(wèi)星)?；聟^(qū)域中某特征點的形變曲線顯示(圖11)，此處滑坡在2018年1月—2019年11月間不同時刻有明顯的形變趨勢。

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

4 結(jié)論

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

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