基于高分辨率SAR的三峽庫(kù)區(qū)地表形變分析

婁連惠 劉強(qiáng) 譚玉敏

摘要：長(zhǎng)江三峽地區(qū)位于川鄂交接的山地峽谷地區(qū)，山多坡陡，一旦遇上暴雨天氣或地震災(zāi)害極易發(fā)生滑坡、泥石流或巖崩等災(zāi)害。為了更好地掌握三峽庫(kù)區(qū)地面沉降變形情況，以庫(kù)區(qū)內(nèi)地質(zhì)災(zāi)害風(fēng)險(xiǎn)較高的巴東縣城區(qū)為例，選取該區(qū)域20景2016～2017年間的PALSAR2數(shù)據(jù)，基于時(shí)間序列InSAR方法進(jìn)行了地表形變分析;并利用高切坡位移監(jiān)測(cè)數(shù)據(jù)進(jìn)行了驗(yàn)證，以更好地了解該區(qū)域地形穩(wěn)定情況。結(jié)果表明：InSAR監(jiān)測(cè)結(jié)果與該區(qū)域同時(shí)期高切坡監(jiān)測(cè)實(shí)際位移變化情況一致;巴東城區(qū)在此期間整體形變速率較小，絕大部分區(qū)域年均形變速率在20 mm/a以下;部分變形較大區(qū)域與已查明的滑坡區(qū)域和土地利用類(lèi)型變化區(qū)域一致。

關(guān) 鍵 詞：

地表變形; 形變監(jiān)測(cè); InSAR; 巴東城區(qū); 三峽庫(kù)區(qū)

中圖法分類(lèi)號(hào)： TP79

文獻(xiàn)標(biāo)志碼： A

DOI：10.16232/j.cnki.1001-4179.2021.05.018

1 研究背景

長(zhǎng)江三峽地區(qū)位于川鄂交接的山地峽谷地區(qū)，山多坡陡，一旦遇上暴雨天氣或地震災(zāi)害極易發(fā)生滑坡、泥石流或巖崩等災(zāi)害。監(jiān)測(cè)分析三峽庫(kù)區(qū)沿岸地表形變以有效識(shí)別地質(zhì)災(zāi)害隱患對(duì)三峽水庫(kù)的安全運(yùn)行和區(qū)域地質(zhì)災(zāi)害的監(jiān)測(cè)預(yù)防非常重要[1]。巴東屬于三峽庫(kù)區(qū)建成后的遷移建設(shè)城區(qū)，是長(zhǎng)江三峽沿岸地質(zhì)條件最為復(fù)雜的地區(qū)之一，區(qū)內(nèi)地質(zhì)構(gòu)造復(fù)雜，是地質(zhì)災(zāi)害高風(fēng)險(xiǎn)區(qū)域[2-3]，歷史上曾發(fā)生黃土坡滑坡等大型滑坡，城區(qū)位置受到地質(zhì)災(zāi)害的影響曾發(fā)生過(guò)三次選址和兩次搬遷。

相對(duì)于水準(zhǔn)測(cè)量、GPS等常規(guī)地面沉降測(cè)量方法，InSAR方法具有全天候、大范圍同步的優(yōu)點(diǎn)，其在大區(qū)域高精度地表形變和滑坡監(jiān)測(cè)應(yīng)用中的優(yōu)勢(shì)已經(jīng)被廣泛認(rèn)可，已成為國(guó)際上地質(zhì)災(zāi)害、地質(zhì)環(huán)境變化監(jiān)測(cè)與防治工作的重要技術(shù)手段[4-10]。在三峽庫(kù)區(qū)滑坡監(jiān)測(cè)方面已有多項(xiàng)關(guān)于InSAR應(yīng)用的研究[11-14]，以及時(shí)間序列InSAR在巴東滑坡監(jiān)測(cè)中的應(yīng)用研究[15-16]，尚未有研究從InSAR用于庫(kù)岸區(qū)域整體穩(wěn)定性的角度來(lái)對(duì)周期性蓄退水狀態(tài)下的庫(kù)區(qū)地質(zhì)災(zāi)害風(fēng)險(xiǎn)進(jìn)行研究。

本文選用高分辨率的PALSAR2數(shù)據(jù)，基于時(shí)間序列InSAR方法進(jìn)行地表形變分析，并利用高切坡位移監(jiān)測(cè)數(shù)據(jù)進(jìn)行驗(yàn)證，以更好地了解該區(qū)域地形穩(wěn)定情況。

2 數(shù)據(jù)來(lái)源及處理

2.1 研究區(qū)數(shù)據(jù)來(lái)源

本文用到的數(shù)據(jù)包括高分辨率SAR數(shù)據(jù)PALSAR2和高切坡監(jiān)測(cè)點(diǎn)位移數(shù)據(jù)。ALOS-PALSAR數(shù)據(jù)是繼JERS-1和ADEOS之后的新一代L波段的SAR數(shù)據(jù)，PALSAR數(shù)據(jù)分為PALSAR1和PALSAR2兩種，PALSAR1數(shù)據(jù)只到2011年，PALSAR2數(shù)據(jù)從2014年開(kāi)始，數(shù)據(jù)的時(shí)間分辨率為14 d，數(shù)據(jù)分辨率遠(yuǎn)高于PALSAR1，各種極化數(shù)據(jù)均有，但數(shù)據(jù)價(jià)格較高，可用數(shù)量受限。本文研究中選用了20景巴東地區(qū)2016年8月到2017年10月間的PALSAR2數(shù)據(jù)，數(shù)據(jù)分辨率為3 m。

2.2 PALSAR數(shù)據(jù)處理

SBAS作為較為典型的時(shí)間序列InSAR形變監(jiān)測(cè)分析技術(shù)已經(jīng)得到廣泛應(yīng)用[13，17]。本文中SAR數(shù)據(jù)的處理主要利用Sarscape5.3軟件工具，采用基于SBAS的時(shí)間序列InSAR數(shù)據(jù)處理方法，步驟如下。

（1） 生成連接圖。根據(jù)時(shí)間和空間基線(xiàn)情況進(jìn)行干涉像對(duì)配對(duì)，將會(huì)對(duì)配對(duì)的像對(duì)進(jìn)行干涉處理。

（2） 干涉工作流處理。對(duì)步驟（1）中生成的所有干涉像對(duì)進(jìn)行干涉處理，包括生成相干性系數(shù)圖、去平、濾波和相位解纏。

（3） 軌道精煉和重去平?？梢怨浪愫腿コ龤堄嘞辔灰约敖饫p后依然存在的相位坡道。

（4） 形變速率反演。該步驟是SBAS技術(shù)的核心步驟之一，主要是估算形變速率和殘余相位并進(jìn)行大氣濾波，通過(guò)時(shí)間序列形變情況估算和去除大氣相位的影響，進(jìn)而得到時(shí)間序列形變圖。

（5） 地理編碼。對(duì)所獲得的雷達(dá)坐標(biāo)系下的速率圖進(jìn)行地理編碼。

首先是連接圖生成，本文設(shè)置了2016年12月19日的1景PALSAR數(shù)據(jù)作為超級(jí)主影像，時(shí)間基線(xiàn)閾值設(shè)置為120 d，空間基線(xiàn)閾值設(shè)置為垂直基線(xiàn)的1.2%，共組成了74個(gè)干涉對(duì)，如圖1所示。

連接圖生成后進(jìn)行干涉工作流處理，PALSAR數(shù)據(jù)雖然分辨率較高，但是考慮到巴東地區(qū)的復(fù)雜地表情況，本文選取了Goldstein濾波方法，濾波效果較為明顯。對(duì)于軌道精煉和重去平步驟，PALSAR數(shù)據(jù)分辨率高，較長(zhǎng)的波長(zhǎng)也反映出了更高的相干性，故在選取控制點(diǎn)時(shí)選擇了地表比較穩(wěn)定且相干性高的一些點(diǎn)。重去平之后的相位圖如圖2所示（以2016年8月15日與2016年10月24日為例）。

相位解纏是SBAS技術(shù)中最重要的步驟之一，解纏結(jié)果直接決定了形變結(jié)果的精度。目前存在多種解纏算法，最小費(fèi)用流法通過(guò)先求解纏繞相位與解纏相位的差的最小值，進(jìn)而將求解轉(zhuǎn)換為最小費(fèi)用流的問(wèn)題，該算法在計(jì)算過(guò)程中會(huì)根據(jù)流的大小和方向?qū)Ω鱾€(gè)相位矩陣進(jìn)行積分，得到高質(zhì)量的解纏結(jié)果，再向低質(zhì)量區(qū)域進(jìn)行積分，得到所有的結(jié)果，這樣能夠很有效地避免相位誤差從低相干性區(qū)域傳到高相干性區(qū)域。巴東地區(qū)地表情況復(fù)雜，城區(qū)以外范圍植被茂密，存在大范圍的低相干性區(qū)域，采用最小費(fèi)用流算法進(jìn)行解纏可以得到高精確度的城區(qū)形變結(jié)果。

地表形變反演分為2次進(jìn)行。第一次反演用時(shí)較長(zhǎng)，主要進(jìn)行地表形變速率的估算和殘余地形相位的估算，同時(shí)會(huì)對(duì)干涉圖進(jìn)行二次解纏并優(yōu)化。第二次反演主要是根據(jù)所有數(shù)據(jù)進(jìn)行大氣的濾波，估算和去除大氣引起的相位變化，進(jìn)一步明確區(qū)域?qū)嶋H地形變相位信息，反演得到的地形變速率如圖3所示。

從圖3可以看出：巴東地區(qū)在2016～2017年間整體地形變速率不高，絕大部分區(qū)域年均地表形變速率在20 mm/a以下，居民區(qū)穩(wěn)定性較好。對(duì)比同時(shí)期光學(xué)遙感影像可以知道：在2016年8月到2017年6月的10個(gè)月間，巴東地區(qū)大部分地區(qū)土地利用類(lèi)型沒(méi)有發(fā)生變化，整體較為穩(wěn)定;江兩岸有一塊緩沖帶的土地類(lèi)型由水體變?yōu)榫用竦睾吐阃?，地表形變較為明顯;北岸非居民地地區(qū)，植被覆蓋有所增加，居民區(qū)裸土的面積有所增加，地表形變速率較高;南岸主要變化發(fā)生在中間的道路兩側(cè)，有兩塊相對(duì)較大的新增裸土地塊;遙感影像與圖3中信息較為吻合。

3 形變精度驗(yàn)證與分析

對(duì)于InSAR應(yīng)用中的精度驗(yàn)證問(wèn)題，有很多學(xué)者做了研究與對(duì)比[18-19]，通常采用將InSAR結(jié)果與GPS或者精密水準(zhǔn)測(cè)量數(shù)據(jù)進(jìn)行對(duì)比的方法?？赡苡绊懢缺容^的因素主要有以下幾點(diǎn)：

（1） InSAR形變結(jié)果為面狀數(shù)據(jù)集，在相干性較差的區(qū)域會(huì)考慮利用周?chē)喔尚院玫膮^(qū)域進(jìn)行輔助解纏，因此在這些區(qū)域InSAR結(jié)果會(huì)存在一些誤差，而實(shí)測(cè)數(shù)據(jù)為點(diǎn)狀數(shù)據(jù)，測(cè)量過(guò)程很少考慮地形以及地表覆蓋的情況，因此把實(shí)測(cè)數(shù)據(jù)和InSAR監(jiān)測(cè)數(shù)據(jù)統(tǒng)一在高相干性區(qū)域是得到精確驗(yàn)證結(jié)果的前提。

（2） InSAR結(jié)果和實(shí)測(cè)數(shù)據(jù)往往處于不同坐標(biāo)系下，坐標(biāo)系的轉(zhuǎn)換會(huì)造成位置的偏差。針對(duì)這些位置偏差，在整體形變連續(xù)且穩(wěn)定的區(qū)域可以采用臨近點(diǎn)搜索的方式進(jìn)行驗(yàn)證[19]。

（3） InSAR反演的地表形變往往為沿衛(wèi)星視線(xiàn)向（LOS）方向的位移，由于不同衛(wèi)星的視線(xiàn)方向?qū)τ诓煌挠^(guān)測(cè)點(diǎn)都不盡相同，想要得到真實(shí)的地表形變位移需要結(jié)合升降軌數(shù)據(jù)進(jìn)行解算，但是由于數(shù)據(jù)具體情況在有些位置難以實(shí)現(xiàn)，可將實(shí)測(cè)的三維位移數(shù)據(jù)投影到衛(wèi)星視線(xiàn)向，與衛(wèi)星視線(xiàn)向的位移進(jìn)行比較，但是該方法對(duì)于實(shí)測(cè)數(shù)據(jù)的三維精度要求較高，在推算過(guò)程涉及到復(fù)雜的投影計(jì)算。

在投影方面，本文給出了詳細(xì)的投影公式和原理。對(duì)于本文實(shí)驗(yàn)區(qū)，通過(guò)實(shí)測(cè)數(shù)據(jù)得出地表的整體位移為垂直方向，水平向的位移較小，故采用了將InSAR觀(guān)測(cè)值往垂直方向進(jìn)行投影計(jì)算再與實(shí)測(cè)數(shù)據(jù)進(jìn)行比對(duì)驗(yàn)證的方法，驗(yàn)證流程如圖4所示。首先對(duì)SBAS結(jié)果與高切坡的實(shí)際測(cè)量值進(jìn)行處理，即通過(guò)坐標(biāo)轉(zhuǎn)換將高切坡的實(shí)測(cè)值轉(zhuǎn)換到WGS84坐標(biāo)系下，與SBAS結(jié)果進(jìn)行統(tǒng)一。本文采用了63個(gè)高切坡監(jiān)測(cè)點(diǎn)的同期位移監(jiān)測(cè)數(shù)據(jù)，點(diǎn)位分布如圖5所示（地圖來(lái)自谷歌地球）。從圖5中可以看出，高切坡監(jiān)測(cè)點(diǎn)主要集中在城市區(qū)域，一般城市區(qū)域人工地物較多，SAR圖像可表現(xiàn)出較好的相干性。

由圖6可以看出：SBAS監(jiān)測(cè)結(jié)果和高切坡實(shí)測(cè)結(jié)果在大部分點(diǎn)位都較為吻合，但是也存在一些誤差較大的點(diǎn)，比如序號(hào)為33，35，41，42，45，48的誤差均超過(guò)了20 mm，最大誤差達(dá)到85 mm，導(dǎo)致所有高切坡監(jiān)測(cè)點(diǎn)位移變化與對(duì)應(yīng)SBAS提取值的差值中誤差為±15.50 mm，大誤差存在的原因包括：因InSAR監(jiān)測(cè)范圍為整個(gè)面，在一些植被覆蓋茂密的區(qū)域往往會(huì)造成數(shù)據(jù)相干性較低，過(guò)低的相干性導(dǎo)致出現(xiàn)一些誤差較大點(diǎn);其次，對(duì)于一些形變量較大的區(qū)域，存在較大的形變梯度，在應(yīng)用SBAS過(guò)程中對(duì)于大的形變梯度同樣可能會(huì)出現(xiàn)失相干現(xiàn)象[20]。

圖6中誤差較大的點(diǎn)主要位于巴東城區(qū)以外，多為植被覆蓋茂密或者地形坡度較大的區(qū)域，這些區(qū)域居住人口往往較少，對(duì)災(zāi)害預(yù)防和工程建設(shè)的需求較小。對(duì)于復(fù)雜地區(qū)的地表形變監(jiān)測(cè)也要保證在一些失相干區(qū)域產(chǎn)生的誤差不會(huì)傳遞到相干性較好的區(qū)域內(nèi)。故本文接下來(lái)選取了城市區(qū)域等多為人工地物附近的高切坡實(shí)測(cè)數(shù)據(jù)進(jìn)一步進(jìn)行了精度評(píng)定，排除6個(gè)差值異常點(diǎn)后的中誤差為±3.72mm，結(jié)果如圖7所示。

由城區(qū)的高相干性區(qū)域的擬合結(jié)果可以看出，高切坡的實(shí)測(cè)值和SBAS監(jiān)測(cè)結(jié)果存在著高度的一致性。本文進(jìn)一步針對(duì)城區(qū)的高相干性區(qū)域進(jìn)行精度評(píng)定，線(xiàn)性回歸圖如圖8所示。

由擬合結(jié)果可以看出，在城市人工地物較多的的區(qū)域，設(shè)置臨近點(diǎn)搜索閾值為9 m情況下，高切坡實(shí)測(cè)值和SBAS檢測(cè)結(jié)果表現(xiàn)出了高度一致性。

4 區(qū)域地表形變與江邊距離關(guān)聯(lián)性分析

本文中還將PALSAR數(shù)據(jù)得到的形變分析結(jié)果與與江邊距離分區(qū)進(jìn)行了關(guān)聯(lián)性分析，分析過(guò)程如下：

（1） 通過(guò)設(shè)置高相干性閾值提取巴東城區(qū)高精度形變監(jiān)測(cè)點(diǎn)的時(shí)間序列信息;

（2） 將研究區(qū)域按距離岸邊遠(yuǎn)近情況，每隔185 m設(shè)置緩沖區(qū);

（3） 根據(jù)緩沖區(qū)分割步驟（1）中獲取的高精度點(diǎn)位，進(jìn)行統(tǒng)計(jì)并分析。

本文以0.90的高相干性系數(shù)閾值對(duì)獲取的形變圖進(jìn)行點(diǎn)群提取，點(diǎn)群分割示意圖如圖9所示。

按照距離江邊遠(yuǎn)近順序?qū)Ⅻc(diǎn)群分割成了點(diǎn)群1～5。為了進(jìn)一步統(tǒng)計(jì)不同緩沖區(qū)內(nèi)的整體形變趨勢(shì)，同樣將緩沖區(qū)1～5內(nèi)的點(diǎn)分別進(jìn)行統(tǒng)計(jì)，其中每一個(gè)點(diǎn)都包含了2016年8月15日到2017年8月15日之間的時(shí)間序列形變信息，時(shí)間序列形變情況的時(shí)間分辨率由相同時(shí)間內(nèi)所輸入的PALSAR的數(shù)據(jù)量決定。然后，通過(guò)求得每一個(gè)時(shí)期衛(wèi)星視線(xiàn)向的形變量的平均值反映該緩沖內(nèi)的整體形變趨勢(shì)，形變序列折線(xiàn)圖如圖10所示。

從圖10可以看出：各個(gè)緩沖區(qū)的形變量均在-10～10 mm之間，形變較小，形變的整體趨勢(shì)比較明顯，隨著點(diǎn)群1～4與江邊的距離變遠(yuǎn)，緩沖區(qū)1～4內(nèi)的形變趨勢(shì)呈現(xiàn)出逐步減小的態(tài)勢(shì)，并且緩沖區(qū)4和緩沖區(qū)5內(nèi)的形變趨勢(shì)逐步失去規(guī)律性。

5 結(jié) 論

研究結(jié)果表明：巴東城區(qū)在2016～2017期間的整體形變速率較小，已查明局部滑坡區(qū)域和地表土地利用類(lèi)型發(fā)生變化區(qū)域有明顯地表形變。PALSAR數(shù)據(jù)在對(duì)樹(shù)木等地表覆蓋物的穿透性上表現(xiàn)出了良好的性能，整體表現(xiàn)出良好的相干性，在地形變監(jiān)測(cè)中可以達(dá)到很高的精度，但價(jià)格較高，為了控制成本可能會(huì)導(dǎo)致在長(zhǎng)時(shí)間序列的分析過(guò)程中時(shí)間分辨率較低。另外，本文所用PALSAR數(shù)據(jù)均為升軌數(shù)據(jù)，無(wú)法結(jié)合同期降軌數(shù)據(jù)進(jìn)行分析，因此文中所得到的形變速率以及形變量方向均為衛(wèi)星視線(xiàn)向。

參考文獻(xiàn)：

[1] 黃波林，殷躍平，李濱，等.三峽工程庫(kù)區(qū)巖溶岸坡巖體劣化及其災(zāi)變效應(yīng)[J].水文地質(zhì)工程地質(zhì)，2020，47（4）：51-61.

[2] 李浩，涂鵬飛，范國(guó)榮，等.城鎮(zhèn)地質(zhì)災(zāi)害危險(xiǎn)性評(píng)估的研究：以巴東沿渡河鎮(zhèn)為例[J].中國(guó)水運(yùn)（下半月），2017，17（10）：219-221，225.

[3] 陳彪，張錦高，陳彪，等.三峽庫(kù)區(qū)巴東新城區(qū)滑坡地質(zhì)災(zāi)害危險(xiǎn)性分區(qū)評(píng)價(jià)[J].科技管理研究，2009，29（1）：256-259.

[4] 謝謨文，王增福，胡嫚，等.高山峽谷區(qū)D-InSAR滑坡監(jiān)測(cè)數(shù)據(jù)特征分析[J].測(cè)繪通報(bào)，2012（4）：18-21.

[5] 葛大慶.區(qū)域性地面沉降InSAR監(jiān)測(cè)關(guān)鍵技術(shù)研究[D].北京：中國(guó)地質(zhì)大學(xué)，2013.

[6] 陳德良，陸燕燕，賈東振.基于StaMPS-InSAR的常州地區(qū)地表形變監(jiān)測(cè)研究[J].人民長(zhǎng)江，2018，49（12）：59-65.

[7] 劉朋俊，張璐，陳元申，等.南水北調(diào)中線(xiàn)濕陷性黃土區(qū)InSAR時(shí)序分析[J].人民長(zhǎng)江，2020，51（6）：125-128.

[8] JIE D，LIAO M，QIANG X，et al.Detection and displacement characterization of landslides using multi-temporal satellite SAR interferometry：a case study of Danba County in the Dadu River Basin[J].Engineering Geology，2018（240）：95-109.

[9] 白林，江利明，汪漢勝.利用高分辨率TerraSAR-X數(shù)據(jù)監(jiān)測(cè)武漢地區(qū)2013～2015年地面沉降[J].大地測(cè)量與地球動(dòng)力學(xué)，2019，39（8）：832-836.

[10] 王國(guó)強(qiáng)，陳仁琛.基于InSAR技術(shù)的流域庫(kù)岸地質(zhì)災(zāi)害監(jiān)測(cè)研究[J].中國(guó)農(nóng)村水利水電，2020（12）：206-210.

[11] 范景輝，邱闊天，夏耶，等.三峽庫(kù)區(qū)范家坪滑坡地表形變InSAR監(jiān)測(cè)與綜合分析[J].地質(zhì)通報(bào)，2017，36（9）：1665-1673.

[12] SHI X G，ZHANG L，LIAO M S，ET AL.Deformation monitoring of slow-moving landslide with L- and C-band SAR interferometry[J].Remote Sensing Letters，2014，5（11）：951-960.

[13] 鄭岳澤.InSAR技術(shù)在三峽庫(kù)區(qū)滑坡監(jiān)測(cè)中的應(yīng)用研究[D].北京：中國(guó)地質(zhì)大學(xué)（北京），2019.

[14] 徐帥，王尚曉，牛瑞卿.基于InSAR技術(shù)的三峽庫(kù)區(qū)巫山—奉節(jié)段潛在滑坡識(shí)別[J].安全與環(huán)境工程，2020，27（1）：32-38.

[15] LIU P，LI Z H，HOEY T，et al.Using advanced InSAR time series techniques to monitor landslide movements in Badong of the Three Gorges region，China[J].International Journal of Applied Earth Observation and Geoinformation，2013（21）：253-264.

[16] 楊長(zhǎng)江，易祎，趙蓉.基于Sentinel-1數(shù)據(jù)的時(shí)序InSAR技術(shù)在滑坡監(jiān)測(cè)方面的應(yīng)用：以巴東地區(qū)為例[J].科技創(chuàng)新與生產(chǎn)力，2017（4）：55-58，63.

[17] 黃潔慧，謝謨文，王立偉.基于SBAS-InSAR技術(shù)的白格滑坡形變監(jiān)測(cè)研究[J].人民長(zhǎng)江，2019，50（12）：101-105.

[18] BERARDINO P，COSTANTINI M，F(xiàn)RANCESCHTTI G，et al.Use of differential SAR interferometry in monitoring and modelling large slope instability at Maratea （Basilicata，Italy）[J].Engineering Geology，2003，68（1-2）：31-51.

[19] 陳強(qiáng)，劉國(guó)祥，丁曉利，等.永久散射體雷達(dá)差分干涉應(yīng)用于區(qū)域地表沉降探測(cè)[J].地球物理學(xué)報(bào)，2007，7（3）：87-93.

[20] 田馨，廖明生.InSAR技術(shù)在監(jiān)測(cè)形變中的干涉條件分析[J].地球物理學(xué)報(bào)，2013，56（3）：812-823.

（編輯：劉 媛）