王 勇 楊 軍 閆 勇 李 鎖 陳 巖

1 天津城建大學(xué)地質(zhì)與測(cè)繪學(xué)院，天津市津靜路26號(hào)，300384 2 中國科學(xué)院精密測(cè)量科學(xué)與技術(shù)創(chuàng)新研究院大地測(cè)量與地球動(dòng)力學(xué)國家重點(diǎn)實(shí)驗(yàn)室，武漢市徐東大街340號(hào)，430077 3 天津市地質(zhì)工程勘測(cè)設(shè)計(jì)院有限公司，天津市紅旗南路261號(hào)，300191

由于多年露天開采礦石，天津薊州五名山人工堆積斜坡及臨空面曾多次發(fā)生滑坡、崩塌等地質(zhì)災(zāi)害，但未進(jìn)行有效監(jiān)測(cè)和預(yù)警，造成一定經(jīng)濟(jì)損失，嚴(yán)重威脅人民群眾安全。滑坡災(zāi)害具有破壞性大、隱蔽性強(qiáng)等特點(diǎn)，對(duì)滑坡進(jìn)行早期識(shí)別和長(zhǎng)期監(jiān)測(cè)值得重點(diǎn)研究[1]。目前，滑坡監(jiān)測(cè)主要依賴傳統(tǒng)的地表形變監(jiān)測(cè)技術(shù)，如GNSS監(jiān)測(cè)和定期水準(zhǔn)測(cè)量等，雖然精度較高，但只能針對(duì)滑坡體邊緣布設(shè)，無法有效地大面積覆蓋整個(gè)滑坡體，且對(duì)于潛在未知隱患區(qū)監(jiān)測(cè)非常有限[2]。逐漸成熟的時(shí)序InSAR技術(shù)具有精度高、不受天氣限制、監(jiān)測(cè)范圍廣等特點(diǎn)[3]，可以克服該限制。近年來，國內(nèi)外學(xué)者利用InSAR技術(shù)獲取已知或未知的潛在滑坡位置、邊界、運(yùn)動(dòng)趨勢(shì)，回溯評(píng)估滑坡體的運(yùn)動(dòng)學(xué)特征等[4-5]。Goorabi[6]利用InSAR技術(shù)探測(cè)到伊朗北部區(qū)域2017年地震后誘發(fā)的滑坡仍在繼續(xù)移動(dòng)；Bayer等[7]發(fā)現(xiàn)人類工程活動(dòng)已激活意大利某山脈潛在的深層滑坡，并用GPS等數(shù)據(jù)與InSAR數(shù)據(jù)進(jìn)行對(duì)比，結(jié)果具有良好的一致性；王之棟等[8]采用多種InSAR監(jiān)測(cè)技術(shù)進(jìn)行研究，結(jié)果表明可正確識(shí)別70%的形變區(qū)域，證明了不同InSAR技術(shù)可為地質(zhì)災(zāi)害探測(cè)提供幫助。然而不同滑坡體朝向會(huì)使InSAR監(jiān)測(cè)精度受限，可利用滑坡特殊性，認(rèn)為滑坡主要發(fā)生在垂向與坡向，結(jié)合升降軌技術(shù)獲取滑坡體垂向和坡向二維形變，可全面探查天津薊州五名山滑坡災(zāi)害分布情況，為后續(xù)集中監(jiān)測(cè)提供良好基礎(chǔ)。

本文利用2017-03～2021-04 Sentinel-1A升軌和Sentinel-1B降軌數(shù)據(jù)，使用SBAS-InSAR技術(shù)對(duì)天津薊州五名山滑坡隱患進(jìn)行識(shí)別和探測(cè)。同時(shí)聯(lián)合解算滑坡隱患區(qū)垂向和坡向形變場(chǎng)，與GNSS形變數(shù)據(jù)進(jìn)行對(duì)比分析，評(píng)價(jià)SBAS-InSAR技術(shù)的適用性和精度，并結(jié)合區(qū)域內(nèi)降水?dāng)?shù)據(jù)回溯滑坡的運(yùn)動(dòng)規(guī)律和特征，對(duì)天津薊州五名山滑坡災(zāi)害的穩(wěn)定性和范圍進(jìn)行初步評(píng)估，為有關(guān)部門提供精確的形變數(shù)據(jù)資料。

1 研究區(qū)概況及實(shí)驗(yàn)數(shù)據(jù)

1.1 區(qū)域概況

五名山位于天津市薊州區(qū)漁陽鎮(zhèn)與官莊鎮(zhèn)交界處，山腳周圍分布小刀剪營村、大星峪村等多個(gè)村莊。整體地貌屬于低山丘陵，由于受長(zhǎng)年采石影響，部分區(qū)域地貌已發(fā)生較大變化，形成多個(gè)臨空面。其中已知的部分滑坡類型為中深層牽引式巖質(zhì)順層滑坡，大部分滑坡走向與坡向基本一致，曾多次發(fā)生局部巖體順層滑動(dòng)[9]。薊州區(qū)降水主要集中在6～9月，占全年降水量的80%以上，其中汛期主要集中在7、8月，此時(shí)正是強(qiáng)降水引發(fā)地質(zhì)災(zāi)害的多發(fā)期。

圖1為初步圈定的五名山可能存在的滑坡隱患區(qū)域，其中A、D為主要的采石區(qū)，山頂已被削平深挖，呈不規(guī)則漏斗狀，B、C多是在山腳處存在部分開采區(qū)域。已有資料顯示，災(zāi)害類型主要為淺表層風(fēng)化殘坡積土質(zhì)滑坡、尾礦堆碎屑泥石流，以及廢棄礦山采面滑坡和崩塌、山坡崩塌等。

圖1 五名山主要采石區(qū)Fig.1 Main quarrying areas of Wuming mountain

1.2 實(shí)驗(yàn)數(shù)據(jù)

本文選用歐空局的Sentinel-1A/B衛(wèi)星數(shù)據(jù)，成像模式為IW(interferometric wide-swath)模式、SLC(single look complex)、C波段SAR影像數(shù)據(jù)，分辨率為5 m×20 m，采用VV極化方式。Sentinel-1A/B成像時(shí)間為2017-03～2021-04，Sentinel-1A升軌數(shù)據(jù)量為51幅，入射角為35.82°，軌道方位角為-13.38°；Sentinel-1B降軌數(shù)據(jù)量為48幅，入射角為41.76°，軌道方位角為-166.52°。數(shù)字高程模型(DEM)采用30 m空間分辨率的ALOS World 3D(AW3D30)作為外部DEM數(shù)據(jù)，也用于提取坡度坡向信息。

2 SBAS基本原理與二維形變獲取

SBAS技術(shù)主要是通過設(shè)置時(shí)空基線閾值，將參與計(jì)算的所有SAR影像分成多個(gè)子集，選取短基線像對(duì)集合進(jìn)行干涉處理，可有效克服時(shí)空失相干，提高InSAR精度和探測(cè)能力。其核心是利用最小二乘原理獲得各子集的形變時(shí)間序列，通過奇異值分解(SVD)方法估算形變參數(shù)，具體原理及處理流程可參見文獻(xiàn)[10]。SAR衛(wèi)星多數(shù)為側(cè)視成像，無論是升軌還是降軌，對(duì)于單一軌道獲得的形變結(jié)果并不是地表真實(shí)形變，而是三維形變?cè)谝暰€向(LOS)的投影，因此可將升降軌視線向形變分解成3個(gè)維度分量dU、dN、dE(圖2)，即

dlos=dUcosθinc-dEsinθincsin(αazi-3π/2)-

dNsinθinccos(αazi-3π/2)

(1)

式中，dlos為視線向形變；dU、dN、dE分別為地表某點(diǎn)的垂向、水平南北、水平東西方向形變；θinc為SAR衛(wèi)星入射角；aazi為SAR軌道方位角。

圖2 坡面形變分解示意圖Fig.2 Diagram of slope deformation decomposition

基于研究區(qū)先驗(yàn)信息，滑坡大部分屬于順層坡向運(yùn)動(dòng)特性，假設(shè)滑坡體只受重力影響沿坡面均勻滑動(dòng)，而沿坡面垂直坡向方向的運(yùn)動(dòng)可忽略。由圖2可知，坡面坐標(biāo)系由dasp坡向形變、dT垂直坡向形變和dI垂直坡面形變構(gòu)成，即認(rèn)為SAR衛(wèi)星在滑坡面視線向形變主要由dU垂向和dasp坡向形變組成[11]：

(2)

式中，β為坡向。因此，可利用升降軌2個(gè)軌道依據(jù)式(2)聯(lián)合解算獲得研究區(qū)滑坡體垂向和坡向二維形變。其中垂向形變以向上為正、向下為負(fù)；坡向形變則順坡方向?yàn)檎?，反之則為負(fù)。

3 形變結(jié)果分析

3.1 InSAR結(jié)果精度分析

為評(píng)估SBAS-InSAR監(jiān)測(cè)結(jié)果的精度，聯(lián)合研究區(qū)周圍2021年最新建立的GNSS連續(xù)監(jiān)測(cè)站數(shù)據(jù)，經(jīng)處理轉(zhuǎn)化為對(duì)應(yīng)的垂向、坡向相對(duì)形變結(jié)果，并與SBAS-InSAR數(shù)據(jù)進(jìn)行對(duì)比。經(jīng)計(jì)算可知，GNSS數(shù)據(jù)與SBAS-InSAR數(shù)據(jù)的垂向相關(guān)系數(shù)為0.65，均方根誤差(RMSE)為1.44 mm；坡向相關(guān)系數(shù)為0.54，均方根誤差為1.93 mm，SBAS-InSAR結(jié)果與GNSS監(jiān)測(cè)結(jié)果的基本趨勢(shì)相同，具有較好的一致性。

3.2 滑坡隱患分析

通過SBAS-InSAR技術(shù)初步獲得2017-03～2021-04五名山的二維形變結(jié)果。研究結(jié)果表明，五名山滑坡隱患區(qū)最大垂向累積形變達(dá)到-38.28 mm，臨空面坡體最大坡向累積形變?yōu)?2.10 mm。五名山形變區(qū)域主要分布在A、B、C區(qū)域，以下將對(duì)A、B、C區(qū)域的形變特征展開詳細(xì)分析。A區(qū)域以前為五名山主要的采石區(qū)域，其形態(tài)大體呈凹陷漏斗狀，因?qū)Νh(huán)境破壞較大，近幾年逐漸關(guān)停。由圖3可知，A區(qū)內(nèi)不穩(wěn)定區(qū)域集中分布在采石區(qū)域，因A區(qū)大部分地區(qū)受采石影響嚴(yán)重，高程變化巨大，轉(zhuǎn)化成坡向形變時(shí)可能存在錯(cuò)誤，因此只討論垂向形變結(jié)果。

圖3 五名山滑坡A區(qū)垂向形變Fig.3 Vertical deformation of area A of Wuming mountain landslide

A區(qū)采石場(chǎng)西側(cè)形變沿著礦區(qū)邊緣長(zhǎng)度約為760 m，其中最大累積沉降為-23.70 mm；南側(cè)區(qū)域形變集中分布在開采區(qū)兩側(cè)，最大累積沉降為-25.06 mm；北部臨空面上部累積垂向形變?yōu)?17.28 mm。其中不穩(wěn)定區(qū)域面積共約31.5萬 m2。東側(cè)坡體暫未發(fā)現(xiàn)明顯形變。

B、C區(qū)為五名山坡腳區(qū)域，其二維形變?nèi)鐖D4所示。B區(qū)沿公路走向的采石區(qū)具有明顯的滑坡運(yùn)動(dòng)趨勢(shì)，最大垂向形變速率為-4.48 mm/a，最大垂向累積形變?yōu)?23.42 mm，最大坡向形變速率為3.68 mm/a，最大坡向累積形變?yōu)?8.11 mm，形變區(qū)域面積約3.2萬 m2。中部存在2個(gè)臨空面，未出現(xiàn)滑動(dòng)跡象，最南端采石區(qū)坡腳處的最大垂向形變速率為-4.64 mm/a，最大垂向累積形變?yōu)?20.87 mm，面積約0.5萬m2。臨空面山腳處采挖區(qū)累積垂向形變達(dá)到30 mm以上，坡向形變具有沿周邊向中心移動(dòng)的趨勢(shì)。原因?yàn)锽區(qū)山腳處曾在采石區(qū)向下深挖后，為恢復(fù)原有的地質(zhì)環(huán)境，又進(jìn)行簡(jiǎn)單的回填處理，由于回填土沉降，形成現(xiàn)在的漏斗狀沉降分布，因此在InSAR結(jié)果中主要表現(xiàn)為垂向沉降形式。圖4(b)為五名山C區(qū)滑坡，該滑坡位于小刀剪營村北側(cè)。由于坡體前緣臨空，坡體特征為前緣切割和后緣裂縫，存在發(fā)生滑坡的可能。通過本次監(jiān)測(cè)發(fā)現(xiàn)，最大垂向形變速率大于-4 mm/a，最大垂向累積形變?yōu)?25.08 mm，坡向形變速率超過6 mm/a，最大坡向累積形變?yōu)?5.67 mm，形變區(qū)沿坡向最大長(zhǎng)度約227 m，寬度約158 m，面積約2.7萬 m2?；码[患的產(chǎn)生是受坡腳處持續(xù)采石影響，導(dǎo)致上部滑坡體前緣臨空逐漸加大，具有向下蠕滑趨勢(shì)，穩(wěn)定性較差。形變最大區(qū)距山頂約150 m，遠(yuǎn)離坡腳處較穩(wěn)定。如果發(fā)生滑坡會(huì)對(duì)下方居民區(qū)形成巨大的安全隱患，因此后續(xù)應(yīng)重點(diǎn)關(guān)注該區(qū)域。

圖4 五名山B、C滑坡區(qū)Fig.4 Areas B and C of Wuming mountain landslide

3.3 滑坡形變與降水影響因素分析

研究表明，季節(jié)性強(qiáng)降水容易引發(fā)滑坡。分別在B、C坡體上提取形變區(qū)垂向、坡向形變，分析其與日均降水量的關(guān)系，對(duì)比結(jié)果如圖5所示。從圖中可以看出，每年6～9月研究區(qū)月均降水量在200 mm以上，較大降水后滑坡體均表現(xiàn)出不穩(wěn)定狀態(tài)，強(qiáng)降水后坡體形變存在加速趨勢(shì)。

圖5 累積形變與日均降水對(duì)比Fig.5 Comparison of cumulative deformation and daily precipitation

前期降水量較少時(shí)(2017-11、2018-05)，強(qiáng)降水會(huì)導(dǎo)致坡體形變位移快速增加，這可能與土體含水率極速增多導(dǎo)致坡體的穩(wěn)定性變差有關(guān)；而連續(xù)降水(2019-06～2020-06)對(duì)滑坡形變影響較小，并未明顯促使滑坡位移出現(xiàn)加速；滑坡形變峰值對(duì)于降水峰值具有滯后性，但滯后時(shí)間存在差異，需要一定的時(shí)間延遲才能對(duì)降水作出響應(yīng)。從圖5累積形變量可以看出，2018年雨季過后，五名山主要滑坡隱患區(qū)的垂向和坡向累積形變均明顯快速增大。以C區(qū)滑坡體為例，繪制2018-07-24強(qiáng)降水前后垂向、坡向形變結(jié)果(圖6)。

圖6 2018年強(qiáng)降水前后C區(qū)垂向、坡向形變對(duì)比Fig.6 Comparison of vertical and slope deformation in area C before and after heavy rainfall in 2018

雨季過后坡體中下部位移逐漸加大，且潛在滑坡體面積增大，說明降水會(huì)加劇潛在滑坡形變。隨著形變的逐漸積累，坡體可能在某次強(qiáng)降水過后，從蠕滑狀態(tài)達(dá)到臨界值而形成滑坡，對(duì)坡下居民區(qū)造成嚴(yán)重影響。這表明大量降水會(huì)導(dǎo)致潛在滑坡體被激活，位移逐漸增大，因此降水是影響五名山坡體滑坡的主要外部因素之一。

4 結(jié) 語

本文利用SBAS-InSAR技術(shù)處理2017～2021年共99景覆蓋天津市薊州區(qū)五名山區(qū)域的Sentinel-1A/B影像，獲取垂向和坡向形變，識(shí)別滑坡災(zāi)害以及隱患點(diǎn)，并開展滑坡原因分析，獲得以下結(jié)論：

1)在2017-03～2021-04期間，五名山滑坡隱患區(qū)域中最大垂向累積形變達(dá)到-38.28 mm，臨空面坡體最大坡向累積形變?yōu)?2.10 mm。隱患區(qū)域中A區(qū)形變主要集中在露天采石區(qū)；B區(qū)坡向形變表現(xiàn)為沿周邊向中心移動(dòng)趨勢(shì)；C區(qū)由于滑坡體前緣臨空，處于緩慢蠕滑階段；D區(qū)較穩(wěn)定，未出現(xiàn)明顯形變。

2)強(qiáng)降水是影響滑坡形變的重要自然因素，降水量大小及降水持續(xù)時(shí)間對(duì)坡體的穩(wěn)定性具有不同影響，突發(fā)大量降水會(huì)加速滑坡位移，而持續(xù)少量降水對(duì)滑坡影響較小，且降水峰值與滑坡形變峰值并不重合，存在一定的時(shí)間延遲。