鄭萬基, 孫 倩, 劉小鴿
(1. 中南大學 地球科學與信息物理學院,湖南 長沙 410083;2. 湖南師范大學 資源與環(huán)境科學學院,湖南 長沙 410081;3. 湖南師范大學 地理空間大數(shù)據(jù)挖掘與應用湖南省重點實驗室,湖南 長沙 410081)
金沙江位于長江上游,流經(jīng)中國西藏、云南及四川等地,流域內(nèi)坡陡水急,落差高達5100 m,蘊含豐富的水能資源,是中國具有重要戰(zhàn)略地位的水電基地。然而,金沙江流域內(nèi)地形情況極為復雜,氣候多變,易發(fā)生滑坡及泥石流等自然災害。伴隨著流域內(nèi)水庫蓄水,庫區(qū)內(nèi)的地下水會發(fā)生相應的變化,使得穩(wěn)定的堆積體有再次發(fā)生滑動的可能。庫區(qū)內(nèi)發(fā)生滑坡,可能導致建筑物發(fā)生損毀、水庫造成淤積,形成涌流,給經(jīng)濟生產(chǎn)建設(shè)帶來極大的安全隱患,因此需要進行長期的地表形變監(jiān)測。目前常用的監(jiān)測技術(shù)主要包括GPS、水準等地面測量手段,然而這些傳統(tǒng)手段已經(jīng)難以滿足大范圍、高精度、高時空分辨率的監(jiān)測需求,并且需要投入大量人力物力,作業(yè)人員的人身安全難以得到保障。近些年發(fā)展起來的時序合成孔徑雷達干涉測量(InSAR)技術(shù),在地震[1]、火山[2]、滑坡[3]及城市基礎(chǔ)設(shè)施監(jiān)測[4]等方面已經(jīng)得到較為廣泛的應用,該技術(shù)彌補了傳統(tǒng)監(jiān)測手段空間分辨率低的不足,無需投入大量的人力物力,兼具較高的監(jiān)測精度,其cm到mm級的監(jiān)測精度已得到廣泛的驗證[5-6]。
目前使用較為廣泛的時序InSAR監(jiān)測方法主要是SBAS(小基線集)技術(shù)[7]及PSI(永久散射體InSAR)技術(shù)[8]。然而SBAS技術(shù)難以解決大氣延遲等與空間相關(guān)的誤差,且需要進行相位解纏,在面對一些“孤島”區(qū)域時,易產(chǎn)生解纏誤差;而PSI技術(shù)則需要監(jiān)測對象在時間序列上保持較高的相干性,否則將難以選取較多的形變監(jiān)測點,降低InSAR監(jiān)測的空間分辨率。TCP-InSAR技術(shù)是近年發(fā)展起來的一種新型時序InSAR技術(shù)[9],該技術(shù)融合PS技術(shù)與SBAS技術(shù)的特點,在無需解纏的情況下,可以有效地克服與空間相關(guān)的誤差,并確保可以選取較多的監(jiān)測點,保證InSAR監(jiān)測的空間分辨率,其監(jiān)測精度在滑坡[3]、沉降監(jiān)測[10]及火山活動監(jiān)測[11]等領(lǐng)域已得到驗證。
本文基于歐空局的新一代中高分辨率C波段SAR衛(wèi)星——哨兵一號A星(Sentinel-1A),使用TCP-InSAR技術(shù),獲取金沙江中游興培當至草可都段的滑坡結(jié)果。最后根據(jù)監(jiān)測結(jié)果,結(jié)合從NOAA獲取的氣象資料及當?shù)氐刭|(zhì)情況,分析6處滑坡的影響因素。
本文的研究區(qū)域如圖1所示,地處云南省西北部,距離麗江市以北約69 km處的金沙江中游地區(qū),該河段距離下游的梨園水電站約28 km。如圖1中(c)所示,金沙江連通區(qū)域南北,格吉河與白水河于此處匯入金沙江,右岸有草可都、興培當?shù)染用顸c,左岸則為下只恩水電站和白水河水電站,即從北緯27°27′37.2″至27°30′15.4″,東經(jīng)100°8′45.1″至100°12′22.7″的范圍。該地區(qū)處于滇西北高原地區(qū),地勢起伏較大,山川走向以北北西和北北東為主。研究區(qū)內(nèi)氣候分為干、濕兩季,濕季主要集中于5—10月,年降水量約為900 mm,在雨季集中降雨時,常發(fā)生山洪及泥石流等災害。該地區(qū)地層巖性以玄武巖及灰?guī)r為主,存在有多處活動中的大型堆積體滑坡,滑坡體上分布有少量的耕地和居民,滑坡體物質(zhì)以碎塊石土、殘坡積層、粘土及砂土為主[12]。本文通過TCP-InSAR技術(shù)對該區(qū)域的滑坡體進行監(jiān)測,而InSAR技術(shù)通常會受地表覆蓋植被的影響,影響大小與植被高度、SAR信號頻率有關(guān)[13]。研究區(qū)海拔均處于2 500 m以上,沒有高大的樹林遮擋,而由圖1(c)可知,在第一景SAR影像時間點前后,該區(qū)域植被覆蓋較少,個別有植被覆蓋的區(qū)域植被也較為稀疏,不會對C波段的哨兵一號衛(wèi)星造成失相關(guān)。此外,本文采用的TCP-InSAR技術(shù),將根據(jù)時間序列上干涉圖的相干性進行選點,若存在低相干性的植被觀測點,將不會被選中用于形變估計。
實驗采用的數(shù)據(jù)為歐洲空間局的新一代中高分辨率衛(wèi)星哨兵一號A星(Sentinel-1A)數(shù)據(jù)。該衛(wèi)星為C波段的合成孔徑雷達衛(wèi)星,采用全新的TOPSAR成像模式,實現(xiàn)在中等分辨率(距離向3.7 m,方位向14 m)下,具備12 d的重返周期及250 km×250 km的大范圍觀測能力。在A星與B星串聯(lián)模式下,更可實現(xiàn)6 d的重返觀測。此外,基于歐盟的哥白尼計劃,該衛(wèi)星目前可免費提供全球最新的實時觀測數(shù)據(jù),這為實現(xiàn)高時空分辨率的滑坡監(jiān)測提供重要的數(shù)據(jù)保障。實驗從歐洲空間局的哨兵科學數(shù)據(jù)中心處獲取17景Sentinel-1A數(shù)據(jù),時間跨度為2015-12-18—2017-02-10,覆蓋區(qū)域如圖1(b)中紅色方框所示。由于哨兵衛(wèi)星具備較好的軌道控制能力,因此此次獲取的哨兵數(shù)據(jù)其垂直基線皆控制在150 m以內(nèi),具體數(shù)據(jù)參數(shù)如表1所示。
注:圖(a)(b)底圖為當?shù)谼EM數(shù)據(jù),圖(c)底圖為google獲取的2015-11-22影像圖1 研究區(qū)位置及數(shù)據(jù)覆蓋情況
序號時間升/降軌入射角/(°)垂直基線/m時間基線/d12015-12-18升軌33.90022016-01-11升軌33.926.92432016-02-04升軌33.94.84842016-02-28升軌33.9-67.57252016-04-16升軌33.9123.712062016-05-10升軌33.9-71.114472016-06-03升軌33.9-18.216882016-07-21升軌33.950.821692016-08-14升軌33.9-44.4240102016-09-07升軌33.9-10.5264112016-10-01升軌33.913.6288122016-10-25升軌33.9-50.3312132016-11-18升軌33.9-42.3336142016-12-12升軌33.951.8260152017-01-05升軌33.98.3384162017-01-29升軌33.9-14.7408172017-02-10升軌33.916.6420
TCP-InSAR技術(shù)是近幾年發(fā)展起來的一種新興的時序InSAR技術(shù)[9],該技術(shù)是在早期PSI和SBAS等常用時序InSAR技術(shù)的基礎(chǔ)上,充分結(jié)合二者的優(yōu)勢發(fā)展出來的。在干涉圖組合上,TCP-InSAR技術(shù)使用SBAS技術(shù)中的多主影像、短時空基線的干涉數(shù)據(jù)組合方法,有效減少干涉圖中的時空失相干噪聲;而在觀測量的選取上,該技術(shù)則使用類似于PSI技術(shù)的以構(gòu)網(wǎng)弧段作為觀測量的方法,通過相鄰點互相差分的方式,有效地抑制了大氣等空間相關(guān)的噪聲。而與二者不同的是,該技術(shù)在參數(shù)估計時,除估計線性形變速率及高程殘余外,還同時估計由于基線不準確導致的軌道趨勢項[14];并且,該技術(shù)選取的相干點,并不需要在整個時間序列上保持高相干特征,該方法根據(jù)數(shù)據(jù)配準時的偏移量信息或相干性信息,選取出在整個時間序列及在部分序列上保持高相干的點,極大擴展觀測點的數(shù)量[15]。此外,TCP-InSAR技術(shù)的另一大優(yōu)勢是,該方法不需要對差分干涉圖進行解纏處理[16]?;诖蟛糠窒噜廡CP點在短時空基線的差分干涉對上不包含整周模糊度的思想,TCP-InSAR技術(shù)在進行第一次參數(shù)估計之后,將根據(jù)觀測殘差判斷出具有整周模糊度的弧段,并將其剔除。TCP-InSAR技術(shù)的技術(shù)流程如圖2所示。
圖2 TCP-InSAR數(shù)據(jù)處理流程
實驗首先使用GAMMA軟件將所有的SAR數(shù)據(jù)配準于同一空間坐標系下,為了避免TOPSAR模式的Sentinel-1A數(shù)據(jù)在配準時在不同的burst上產(chǎn)生的相位跳變現(xiàn)象,實驗過程中利用USGS提供的30 m分辨率的SRTM數(shù)據(jù)用于輔助SAR數(shù)據(jù)配準。在完成配準后,將試驗區(qū)域單獨裁出,進行差分干涉處理,以減少數(shù)據(jù)的運算量,在處理過程中使用5×1的多視視數(shù)(最終影像分辨率約為14 m×14 m)。差分過程使用USGS提供的30 m分辨率的SRTM數(shù)據(jù)作為外部DEM去除地形相位。為了避免產(chǎn)生孤立的小基線集,降低解算結(jié)果的穩(wěn)健性,同時還需保持較高的相干性以抑制時空失相干噪聲,以80 m作為空間基線閾值;此外,對于滑坡而言,當滑坡的累積形變超出InSAR的監(jiān)測梯度后,會造成失相關(guān),影響InSAR監(jiān)測結(jié)果精度,因此以90 d作為時間基線的閾值,由此最終選出如圖3所示的24對干涉基線構(gòu)成短基線集。此外,在滑坡地區(qū),短時間內(nèi)滑坡可能會因為崩塌而對InSAR信號造成嚴重的失相關(guān),影響TCP-InSAR的監(jiān)測質(zhì)量,為此,以相干性0.4作為閾值選取出滑坡上未崩塌的形變點,最終在該地區(qū)沿金沙江兩岸共選出33 987個TCP觀測點構(gòu)成觀測網(wǎng)絡。由于該試驗區(qū)域選取的TCP點位在金沙江兩岸,中間有金沙江隔斷,為了使結(jié)果具有較高的穩(wěn)健性,將TCP點分成東岸TCP點及西岸TCP點單獨進行處理。在西岸,通過第一次參數(shù)估計,獲得每個弧段的線性形變速率及高程殘差,基于第一次參數(shù)估計的結(jié)果,以1.5作為TCP-InSAR模糊度探測的閾值[15],剔除具有模糊度的觀測弧段,同時將軌道趨勢誤差從每個觀測弧段中去除,并進行第二次的參數(shù)再估計,最終得到線性形變速率及高程殘余,將線性形變速率加到每個SLC的殘余相位中,得到時間序列上的形變結(jié)果。求得每個弧段的線性形變速率及時序形變序列后,選取如圖1(c)所示金沙江西岸的下之恩水電站作為穩(wěn)定的解算起始點,進行網(wǎng)絡平差,求得每個TCP點的年均形變速率及時序形變結(jié)果。在金沙江東岸,以同樣的參數(shù)進行處理,即得到東岸TCP點的線性形變速率、高程殘差及形變時間序列,在東岸選取的解算起始點為遠離形變區(qū)、地形較高且地勢平坦的高相干點,如圖1(c)金沙江東岸所示。
圖3 TCP-InSAR構(gòu)建的小基線集
圖4 研究區(qū)雷達LOS向年均形變速率場,底圖為google earth 2015-11-22影像
如圖4所示,通過TCP-InSAR技術(shù)最終獲取金沙江中游興培當至草可都段的雷達視線向(LOS向)的年均地表形變信息。LOS向形變?yōu)榈乇碚鎸嵭巫冊诶走_實現(xiàn)方向上的投影,靠近衛(wèi)星為正值,遠離衛(wèi)星為負值。由于目前多數(shù)星載SAR衛(wèi)星均以極地軌道的方式運行,飛行方向近似南北向,因此LOS向形變中以東西向和垂直向投影結(jié)果的貢獻最多。從圖4中可以發(fā)現(xiàn),該地區(qū)有多處滑坡區(qū)域,其中對金沙江會產(chǎn)生直接影響的滑坡則主要有4處,如圖4紅色虛線區(qū)域所示,形變方向主要以由兩岸向金沙江滑動為主,按照形變量的大小,分別為位于格吉河口與白水河口之間的滑坡A、草可都滑坡B、白水河口右岸滑坡C、興培當滑坡D。由于本次實驗采用的是升軌SAR數(shù)據(jù),因此,在金沙江西側(cè)形變?yōu)樨撝?,在金沙江東側(cè)形變?yōu)檎?。位于金沙江東岸的滑坡體為興培當滑坡B和草可都滑坡D,草可都滑坡B的形變量約為25~50 mm/a,由于形變較大的緣故,在滑坡的中后部出現(xiàn)形變超出監(jiān)測梯度的情況,因此,在草可都滑坡的中后部無法選出相干性較高的點;興培當滑坡D的形變量約為15~27 mm/a,二者均為厚層大型堆積體滑坡。在金沙江西岸的主要滑坡體則為白水河口滑坡C,及位于格吉河口與白水河口之間的滑坡A。其中,滑坡A的形變量級為-40~-65 mm/a,而滑坡C為厚層中型堆積體滑坡,形變量級在-35~-47 mm/a之間。
根據(jù)Brabb等人的實驗結(jié)果,地質(zhì)、土壤及坡度是影響滑坡穩(wěn)定性的主要原因[17]。草可都滑坡B的前緣坡度為50°~65°,中后部約為10°~30°,而興培當滑坡D的前緣坡度為45°~75°,中后部為5°~25°,因此滑坡坡度較大的草可都滑坡的形變量比興培當滑坡的形變量要大。此外,根據(jù)現(xiàn)場考察資料[12],在草可都滑坡體的兩側(cè)有發(fā)育小沖溝,前緣有泉水點出露,從而導致該地區(qū)土壤含水量增加,孔隙壓力變大,有效正壓力因此減小,使滑坡速度加快,也是使得草可都滑坡速度比興培當滑坡速度大的原因。并且,草可都滑坡的前緣出現(xiàn)鼓脹裂縫及崩塌的現(xiàn)象[12],因此在實驗結(jié)果中,草可都滑坡的前緣出現(xiàn)失相干現(xiàn)象。對于滑坡A所在區(qū)域,其坡度前緣為5°~10°,中后部約為10°~25°;而對于滑坡C所在區(qū)域,其坡度前緣為15°~35°,中部為5°~15°,后部為25°~45°。雖然滑坡A的地形坡度小于滑坡C,然而滑坡A的地表滑動速度卻大于滑坡C,根據(jù)現(xiàn)場資料顯示[12],滑坡A的堆積物以礫石與少量粉土組成,而滑坡C的堆積物為碎塊石土,粉土塑性指數(shù)更大,因此滑坡A的形變更大。
為了進一步了解4個形變區(qū)域的動態(tài)變化過程,通過TCP-InSAR技術(shù)得出不同時間序列上的LOS向形變場。圖5顯示從2015-12-18—2017-02-10共16個時刻的形變。由時序結(jié)果中可看出,在滑坡A,B,C和D雖然在不同的時間形變的速度變化不一,然而主要形變?nèi)砸跃€性形變?yōu)橹?,隨著時間的推移形變量在不斷地增加。此外,在個別區(qū)域出現(xiàn)明顯的季節(jié)性形變,如圖中的P5點和P6點。
圖5 時序形變場,圖中P1~P4點為形變區(qū)域A~D中的一點,P5~P6為出現(xiàn)季節(jié)性形變區(qū)域中的一點
其中,P6點位于草可都滑坡與興培當滑坡之間,其形變在2016年7月達到約25 mm后開始下降,而后在2016年9月再次達到25 mm,且形變范圍更大;而P5點從2016年2月開始形變逐漸增大,在2016年8月達到最大35 mm的形變后,開始逐漸下降。兩點的最大形變集中于2016年的夏季,而此時屬于該地區(qū)的濕季,是降水最多的時期,推測這兩個地區(qū)的形變與降水有關(guān)。
為了驗證P1-P6點與降水之間存在的關(guān)系,從NOAA(美國國家海洋與大氣管理局)獲得該地區(qū)以南69 km的麗江市在對應時間段內(nèi)的降水與氣溫數(shù)據(jù)作進一步分析,其結(jié)果如圖6所示。氣溫為每日氣溫均值,降水量為SAR時間序列間隔內(nèi)的總降水量,形變量為LOS向累計形變量的絕對值。該地區(qū)在夏季主要受海洋性西南季風和東南季風的影響,降水充沛,而冬季則受西風帶氣流影響,天氣干燥。
圖6 不同形變點累計形變與降水/氣溫之關(guān)系
如圖6所示,位于白水河口右岸的滑坡P3點和位于興培當滑坡的P4點的形變主要以線性形變?yōu)橹鳎S時間遷移形變量不斷遞增,僅在降水較為集中的6—8月出現(xiàn)較為輕微的滑坡速度加快的現(xiàn)象,整體滑動速度則較為穩(wěn)定。而對位于白水河口與格吉河口間的不穩(wěn)定斜坡P1點和位于草可都的P2點而言,其形變主要是持續(xù)的線性形變與季節(jié)性形變共同作用的結(jié)果。兩點在大體上處于形變不斷遞增的趨勢,而在2016年6月至9月期間,即降水最集中的時間,P1點和P2點的形變量都出現(xiàn)較大增長。滑坡的發(fā)生一般與土壤滲透量有直接相關(guān),滲透量一般由土壤入滲率與降水獲得,因此此處結(jié)合相關(guān)文獻中給出的地質(zhì)土壤資料,結(jié)合降水對其進行分析。根據(jù)相關(guān)資料顯示[12],P1點和P2點地區(qū)的滑坡物質(zhì)以粘土、粉土為主,相比起P3點的碎塊石土和P4點的碎塊石砂土而言,P1點和P2點的土壤環(huán)境更為松散,由大氣降水帶來的地表水在這種土壤環(huán)境下更容易滲透,從而導致孔隙壓力增大,致使有效正壓力減小,從而使P1點和P2點比P3及P4點更易受到降水的影響。
相較于前4個點,P5點和P6點雖然在年均形變速率上較小,然而卻表現(xiàn)出十分明顯的季節(jié)性形變現(xiàn)象。P5點位于格吉河口與白水河口之間,是P1點所在山峰的另一側(cè)。該點的形變主要由降水所引起,在降水集中的6月至10月期間該點的形變值在不斷增加,隨后開始遞減,在12月出現(xiàn)少許降水時,該點又開始活躍,而后迅速遞減,與降水表現(xiàn)出較強的相關(guān)性。P6點是一個堆積體,位于含可地區(qū)。該點的形變雖然也呈現(xiàn)出明顯的季節(jié)性趨勢,然而與P5點不同的是,該點的趨勢與氣溫變化的相關(guān)性更大,雖然該點海拔只有1700 m,然而在該點附近有海拔4000 m以上的雪山,因此,推測該點的季節(jié)性主要由于氣溫的升高或降低,會導致雪山融雪或結(jié)冰,引起地下水位變化所導致的,據(jù)此,推斷該點的形變主要是氣溫與大氣降水共同作用的結(jié)果。在P1-P6點的形變時間序列中,除了P5點外,其它點在降水極少的2017-01-05—2017-01-29期間,出現(xiàn)形變速度較快的異常現(xiàn)象。由于這5個點都處于沿江地帶,初步推測有可能與庫區(qū)水位變化有關(guān),然而由于缺乏當?shù)厮臄?shù)據(jù)的緣故,此處不作進一步討論。
本文以云南省麗江市北部的金沙江中游興培當至草可都段的滑坡作為研究對象,基于新一代中高分辨率衛(wèi)星哨兵一號A星,采用TCP-InSAR技術(shù),獲取該地區(qū)2015年12月至2017年2月的LOS向形變信息。形變結(jié)果表明,該地區(qū)山體總體處于不穩(wěn)定狀態(tài),多處山體出現(xiàn)滑坡的現(xiàn)象,在沿江的草可都地區(qū),LOS向年平均形變速率最大達到了約55 mm/a,在形變大的區(qū)域更出現(xiàn)形變超出監(jiān)測梯度的現(xiàn)象。然后,根據(jù)地形及當?shù)氐刭|(zhì)情況,分析該區(qū)域內(nèi)不同滑坡形變速率不同的原因。并由NOAA提供的氣象數(shù)據(jù),分析該地區(qū)共計6個主要滑坡體的影響因素,其中,草可都滑坡與位于格吉河口與白水河口之間的不穩(wěn)定斜坡的形變速率均由降水與線性形變共同影響,較為穩(wěn)定的興培當滑坡和白水河口右岸的滑坡主要受線性形變的影響,而位于格吉河口與白水河口之間的不穩(wěn)定斜坡另一側(cè)的滑坡對于降水則較為敏感,位于含可的堆積體形變則主要與氣溫有關(guān)。
然而,該地區(qū)處于梨園水庫的上游地區(qū),滑坡極有可能與水庫水位存在一定的聯(lián)系,由于水文數(shù)據(jù)的缺乏,本文暫時不能對其進行分析,而真實地表濕度數(shù)據(jù)缺失,也會使得本文的分析可能存在一定的不足。此外,本文只獲取了該地區(qū)的一維LOS向的形變結(jié)果,僅一維的形變結(jié)果,難以對滑坡結(jié)果有充分詳實的理解[18]。因此,未來的工作將可以考慮融合多孔徑InSAR( MAI)技術(shù)或偏移量跟蹤技術(shù)[19],獲取更加全面的滑坡信息。