劉智學(xué)，薛東劍*

（1.成都理工大學(xué) 地球科學(xué)學(xué)院，四川 成都 610059）

我國西南地區(qū)受印度板塊向亞歐板塊俯沖影響，是地震活動(dòng)發(fā)生極頻繁且強(qiáng)烈的區(qū)域，地震帶來的次生地質(zhì)災(zāi)害不可計(jì)數(shù)。在所有地質(zhì)災(zāi)害中，滑坡是最常見的類型，占地質(zhì)災(zāi)害總量的70%以上[1]。由2008年汶川地震（里氏震級(jí)8.0 Ms、矩震級(jí)8.3 Mw）引發(fā)的大光包滑坡是世界上規(guī)模最大的滑坡點(diǎn)之一，其內(nèi)潛藏的大量安全隱患足以對(duì)當(dāng)?shù)卦斐蓴?shù)十年的影響。因此，為避免震后次生災(zāi)害威脅，需要對(duì)高風(fēng)險(xiǎn)區(qū)域進(jìn)行長期監(jiān)測。

早期的滑坡形變探測方法包括大地測量法、GPS法、自動(dòng)伸縮計(jì)法和分布式光纖法等[2]；但這類以點(diǎn)探面的監(jiān)測方法很難全面反映滑坡區(qū)域的形變特征，易導(dǎo)致分析結(jié)果與實(shí)際情況間存在較大誤差。合成孔徑雷達(dá)（SAR）是一種有源成像系統(tǒng)，基于其發(fā)展的DInSAR地表形變監(jiān)測技術(shù)可大面積連續(xù)監(jiān)測地表 cm級(jí)甚至mm級(jí)的形變[3]，在地震形變場提取、火山活動(dòng)以及礦區(qū)沉降監(jiān)測等方面卓有成效，如侯麗燕[4]等基于多次DInSAR處理研究了地震序列發(fā)震斷層特征，加深了對(duì)周邊活動(dòng)斷裂性質(zhì)的認(rèn)識(shí)；Wempen J M[5]利用DInSAR技術(shù)定量分析了礦區(qū)開采初期和隨著采礦方向推進(jìn)所造成的地面沉降，為沉降動(dòng)力的深入研究提供了理論依據(jù)。然而，傳統(tǒng)的DInSAR技術(shù)更注重短期形變監(jiān)測，在解決失相干現(xiàn)象時(shí)略有不足。因此，長時(shí)間序列差分干涉測量方法應(yīng)運(yùn)而生，主要包括適用于城市、鄉(xiāng)鎮(zhèn)等相干目標(biāo)密集區(qū)域的永久散射體干涉測量（PS-InSAR）和同時(shí)適用于城市區(qū)域和自然表面的小基線干涉測量（SBAS-InSAR）方法，后者采用了更多的主影像來構(gòu)造干涉像對(duì)，有效避免了長時(shí)間基線帶來的失相干，在地表沉降與滑坡形變監(jiān)測中應(yīng)用廣泛。程海琴[6]利用時(shí)序InSAR技術(shù)分析了龍門山地區(qū)震后的滑坡災(zāi)害點(diǎn)時(shí)空分布特征，定量揭示了滑坡物質(zhì)遷移特征和震后滑坡體的地表形變特征；陸超然[7]等分別利用差分干涉和SBAS方法進(jìn)行了大范圍地災(zāi)隱患點(diǎn)的探測，并總結(jié)了二者的優(yōu)劣；莫玉娟[8]等基于SBAS技術(shù)提取了小金縣的地表沉降 趨勢。

本文分別利用DInSAR和SBAS-InSAR技術(shù)對(duì)包含研究區(qū)域的22景哨兵一號(hào)（Sentinel-1）C波段數(shù)據(jù)進(jìn)行處理，提取分析了該滑坡點(diǎn)的視線向累積形變量、時(shí)序形變數(shù)據(jù)和年平均形變速率，并對(duì)比了兩種數(shù)據(jù)提取方法的形變結(jié)果，為大光包滑坡穩(wěn)定性探測工作提供了數(shù)據(jù)參考。

1 研究區(qū)概況

大光包滑坡（31e38′N、104e07′E）位于四川省綿陽市安州區(qū)，四川盆地西北與山地區(qū)域交界處，西北方向有大面積山林覆蓋，地勢落差大，平均海拔相對(duì)較高；地質(zhì)結(jié)構(gòu)復(fù)雜多樣且降水量較大，地貌構(gòu)造包括平壩、丘陵以及低中型山脈3種類型。該區(qū)域地處龍門山斷裂帶，地質(zhì)災(zāi)害較多，尤以地震為最，中小型地震頻繁發(fā)生（圖1）。在本文選取的雷達(dá)影像數(shù)據(jù)的成像時(shí)間內(nèi)，四川地區(qū)共發(fā)生36次里氏4.0級(jí)以上的淺源地震，其中包括11次不低于里氏5.0級(jí)的中強(qiáng)地震。滑坡分別距離2019年12月9日安州區(qū)發(fā)生的里氏4.6級(jí)地震和2020年2月3日成都市青白江區(qū)發(fā)生的里氏5.1級(jí)地震只有18 km和85 km。

圖1 研究區(qū)位置示意圖 （審圖號(hào)：GS(2019)3333）

區(qū)域內(nèi)亞熱帶濕潤季風(fēng)氣候顯著，全年降雨量足，冬春與夏秋季變換時(shí)，旱、雨季交替出現(xiàn)，中小型河流交織密集，水資源充沛。為分析研究區(qū)旱、雨季降雨量變化與地表形變速率變化之間的線性關(guān)系，本文制作了滑坡周邊100 km范圍內(nèi)的月均降雨量數(shù)據(jù)柱狀圖，如圖2所示。

圖2 研究區(qū)降雨量數(shù)據(jù)

2 數(shù)據(jù)來源與研究方法

2.1 數(shù)據(jù)來源

Sentinel-1衛(wèi)星所提供的C頻段雷達(dá)遙感影像能提取大范圍、高精度的地面運(yùn)動(dòng)數(shù)據(jù)，是歐洲委員會(huì)和歐空局針對(duì)歐洲航天哥白尼計(jì)劃的第一顆對(duì)地觀測衛(wèi)星[9]。該衛(wèi)星為雙星座衛(wèi)星系統(tǒng)，成像重訪周期最短可達(dá)6 d，且所有數(shù)據(jù)均免費(fèi)提供給公眾。本文采用Sentinel-1A衛(wèi)星干涉寬幅模式數(shù)據(jù)進(jìn)行研究。

2.2 DInSAR技術(shù)

DInSAR技術(shù)基于重復(fù)獲取的同一區(qū)域后向散射信號(hào)，通過共軛相乘獲取干涉相位[10]，獲取相干性良好的雷達(dá)復(fù)數(shù)影像差分干涉兩次成像之間的相位差信息，再去除干涉相位中包含的地形起伏相位信息，減少大氣效應(yīng)、噪聲等的影響，得到最終形變結(jié)果。根據(jù)DInSAR處理時(shí)對(duì)不同時(shí)相SAR影像的數(shù)目要求，常用的 DInSAR 技術(shù)可分為兩軌法、三軌法和四軌法[11]3種，本文采用兩軌法進(jìn)行處理。該方法將研究區(qū)域前后不同時(shí)間形成的兩景雷達(dá)影像作為干涉像對(duì)來進(jìn)行差分干涉處理。由于其獲得的干涉相位中既包含形變相位，又包含地形相位信息（忽略大氣擾動(dòng)的影響），因此需以實(shí)驗(yàn)數(shù)據(jù)成像區(qū)域的先驗(yàn)DEM數(shù)據(jù)為外部參考，對(duì)相同空間幾何條件下的地形相位進(jìn)行模擬，并將其從合成干涉數(shù)據(jù)中去除，這也意味著DEM數(shù)據(jù)的準(zhǔn)確性將直接影響DInSAR的結(jié)果。其數(shù)學(xué)表達(dá) 式為：

式中，ΔΦD為差分相位，是模擬得到的地形相位。

2.3 SBAS-InSAR技術(shù)

Berardino P等為提高差分干涉測量的精度，利用SBAS-InSAR技術(shù)來解決失相干等問題。該方法基于不同時(shí)間段的單視數(shù)據(jù)提取相干目標(biāo)[12]，從而獲取大范圍時(shí)序形變規(guī)律。其優(yōu)勢在于，對(duì)數(shù)據(jù)量要求不高，通過鏈接多個(gè)長時(shí)間序列影像，構(gòu)成小基線集合，進(jìn)而增加數(shù)據(jù)利用率，最終獲得mm級(jí)形變速率結(jié)果。SBAS-InSAR的數(shù)據(jù)處理流程如圖3所示。

圖3 SBAS-InSAR的數(shù)據(jù)處理流程

在SBAS-InSAR測量方法中，將研究區(qū)某一時(shí)間段內(nèi)的N景數(shù)據(jù)按預(yù)設(shè)條件進(jìn)行配對(duì)組合，生成最多M組像對(duì)，即

為獲取合適的干涉像對(duì)進(jìn)行形變結(jié)果提取，需設(shè)置一個(gè)合適的時(shí)空基線閾值來對(duì)數(shù)據(jù)進(jìn)行選擇。由于Sentinel傳感器平臺(tái)重訪周期較短，為提高形變探測精度、避免完全失相關(guān)，本文將最大臨界基線閾值設(shè)置為45%。像對(duì)組合經(jīng)過干涉工作流處理后，將生成解纏后的干涉相位圖，則第i組像對(duì)生成的干涉圖中任一像素的干涉相位為：

式中，x、r分別為像素點(diǎn)的方位向和距離向坐標(biāo)；A、B分別為第i組像對(duì)的主輔影像。

為更準(zhǔn)確地獲取真實(shí)相位值，解算結(jié)果中(-π,π)的相位值，解纏后相位表達(dá)式為：

式中，tn、tn-1分別為該組像對(duì)兩幅影像的成像時(shí)間；vn為對(duì)應(yīng)的平均相位速率。

加入控制點(diǎn)文件對(duì)解纏后的結(jié)果進(jìn)行軌道精煉，用以校正干涉相位以及解纏后所得的相位信息，如恒定相位的去除。本文輸入雷達(dá)坐標(biāo)系下的控制點(diǎn)進(jìn)行優(yōu)化處理，將上述結(jié)果描述為矩陣形式，即

得到一個(gè)M×N矩陣，為避免矩陣秩虧，采用奇異值分解[13]最小二乘處理后，再對(duì)各成像間隔內(nèi)速度積分，即可得到研究區(qū)緩慢形變的最終形變時(shí)間序列。SBAS時(shí)序累積形變?nèi)鐖D4所示。

圖4 SBAS時(shí)序累積形變

3 研究結(jié)果與分析

為了方便對(duì)DInSAR形變數(shù)據(jù)進(jìn)行綜合分析，將每組像對(duì)的累積形變量制作成圖表（圖5），從而直觀反映形變結(jié)果。將每組影像像對(duì)間的累積形變值與地震、降水量數(shù)據(jù)進(jìn)行對(duì)比，可以看出，它們之間存在一定的相關(guān)關(guān)系，如滑坡區(qū)在11月—次年3月步入旱季，降雨量較少，暫不考慮雨水引發(fā)的滑坡體滑動(dòng)，主要查看地震數(shù)據(jù)，這期間四川省內(nèi)共發(fā)生淺源可感地震11次，其中2020年2月的5.1級(jí)地震距離滑坡點(diǎn)僅不到30 km，較之以往冬季的形變量更大，但也沒有產(chǎn)生明顯滑動(dòng)，說明中小型淺源地震的發(fā)生對(duì)滑坡體影響較小，推測地震發(fā)生頻次更高、且距離滑坡體更近時(shí)，滑坡形變數(shù)據(jù)才會(huì)產(chǎn)生更明顯的變化；滑坡點(diǎn)周邊在2018年3月—10月未發(fā)生5級(jí)以上的中型地震，可將地震對(duì)滑坡滑動(dòng)的影響忽略不計(jì)，主要考慮降水量變化，隨著雨季降水量的明顯增加，滑坡也開始表現(xiàn)出更大的形變數(shù)據(jù)波動(dòng)，受滑坡體附近水流沖刷作用影響，滑坡部分地區(qū)有一定程度的 抬升。

圖5 DInSAR形變結(jié)果的累積形變區(qū)間

SBAS-InSAR測量得到的研究區(qū)近兩年累積形變量（圖6）表明，滑坡點(diǎn)及其周邊位置雖一直處于動(dòng)態(tài)變化中，但其在固定時(shí)間段內(nèi)的總體抬升與沉降量相互平衡，因此其一維視線方向上的形變值變化并不明顯，重大地質(zhì)災(zāi)害發(fā)生的可能性較小。

圖6 SBAS累積形變量

為了更深入地了解滑坡點(diǎn)附近的地表形變，本文選取滑坡上的A、B、C、D、E五個(gè)位置，將其平均沉降速率繪制成曲線圖進(jìn)行分析，如圖7所示，可以看出，D、E點(diǎn)有少量抬升，B、C點(diǎn)有沉降現(xiàn)象， A點(diǎn)附近則近于平穩(wěn)。根據(jù)滑坡點(diǎn)周圍地質(zhì)條件推測，坡體高處的風(fēng)化剝蝕堆積物導(dǎo)致了D、E點(diǎn)的抬升，而B、C點(diǎn)受河流沖刷影響，緩慢下降。

圖7 部分點(diǎn)位平均形變速率/(mm/a)

4 結(jié) 語

針對(duì)大光包滑坡實(shí)驗(yàn)區(qū)地形地貌復(fù)雜，實(shí)地勘測地理位置信息難度較大、危險(xiǎn)性較高的特點(diǎn)，本文通過非接觸傳感器獲取的遙感影像來研究地震和降水量對(duì)滑坡體產(chǎn)生的影響。利用Sentinel-1A衛(wèi)星影像獲取實(shí)驗(yàn)區(qū)相位數(shù)據(jù)，下載并處理了2018年3月— 2020年2月的22景干涉寬幅模式的單視復(fù)數(shù)影像；利用DInSAR和SBAS-InSAR技術(shù)計(jì)算并提取了滑坡點(diǎn)在該時(shí)期內(nèi)的階段形變量、時(shí)序累積形變量和年平均形變速率；結(jié)合旱、雨季降雨量變化以及周邊地區(qū)淺源地震發(fā)生頻次對(duì)滑坡形變速率變化規(guī)律進(jìn)行了總結(jié)。研究結(jié)果表明，滑坡體在實(shí)驗(yàn)數(shù)據(jù)成像時(shí)間段內(nèi)的形變量不大，且波動(dòng)較??；滑坡區(qū)域的最大抬升約為10 cm，累積沉降量最大不超過6 cm，年平均形變速率趨于穩(wěn)定，最高不超過100 mm/a，說明大光包滑坡已處于較穩(wěn)定的狀態(tài)；在春、冬交替之際，隨著降雨量的減少，河流沖刷作用減弱，形變趨勢也稍有減緩；周邊地區(qū)中小型淺源地震的發(fā)生沒有對(duì)滑坡點(diǎn)造成顯著危害。