張永滎,張生鵬,2,貢寶扎西,李水平,李亞剛,省天琛
( 1. 青海省基礎測繪院, 西寧 810001;2. 青海師范大學 地理科學學院, 西寧 810008;3. 合肥工業(yè)大學 土木與水利工程學院, 合肥 230009 )
水電站是通過調控水資源時空分布成為蓄水防旱的重要樞紐工程,其運行直接關系到庫區(qū)及下游居民的重大生命財產(chǎn)安全[1]. 中國已建成各類水庫近10萬座,但是多數(shù)水庫由于建造年代久遠,受限于當時的條件,病險情況普遍存在,最嚴重時其數(shù)量占比達到30%. 形變監(jiān)測是水庫安全運營評估的基礎,目前水庫形變監(jiān)測主要采用傳統(tǒng)全站儀、水準測量等手段在重點部位進行監(jiān)測,傳統(tǒng)大地測量的工作特點決定其監(jiān)測結果空間采樣率不足及觀測成果具有周期性,導致在大范圍監(jiān)測結果中無法及時有效地發(fā)現(xiàn)新的形變區(qū)域,造成重大安全隱患的遺漏[2-3].
全球衛(wèi)星導航系統(tǒng)(GNSS)和合成孔徑雷達干涉測量(InSAR)技術作為大地測量新技術,為水庫形變監(jiān)測帶來了新的方案. 王騰等[4]采用InSAR永久散射體(PS)技術主要對三峽水庫大壩進行了穩(wěn)定性分析,得出利用此種技術可以獲得遠高于傳統(tǒng)手段監(jiān)測點的結論. 裴媛媛等[5]采用PS技術對上海部分堤壩進行了穩(wěn)定性分析,得出此方法可以得到毫米級精度,更好地用在人工線狀地物的結論. 張明晶[6]對于橫泉水庫GNSS自動化監(jiān)測進行研究,得出此類系統(tǒng)可以提高監(jiān)測的實時性和可靠性的結論. MILILLO等[7]利用InSAR時序方法研究了伊拉克摩蘇爾大壩2004—2016失穩(wěn)過程,找出了各個時段詳細的形變趨勢. 上述研究都是單獨利用GNSS或InSAR進行監(jiān)測,其中InSAR基本采用PS技術,聚焦點在水庫人工建筑,GNSS監(jiān)測未單獨驗證北斗衛(wèi)星導航系統(tǒng)(BDS)的可用性,也未融合現(xiàn)有連續(xù)運行參考站系統(tǒng)(CORS)技術的優(yōu)勢,其實時定位精度停留在厘米級. 龍羊峽水電站作為黃河上游第一座大型梯級電站,對下游水電站的水量調節(jié)起主導作用,由于其特殊的地形、地貌、巖土體類型及氣象因素的誘發(fā),成為大型滑坡的高發(fā)育區(qū),并且迄今為止還完全采用傳統(tǒng)測量方式進行形變監(jiān)測[8]. 本文綜合考慮GNSS的高時間分辨率和InSAR的高空間分辨率,結合青海省連續(xù)運行參考站系統(tǒng)(QH-CORS)數(shù)據(jù)成果采用InSAR小基線集技術對龍羊峽庫區(qū)進行了地表形變監(jiān)測分析,驗證并采用InSAR技術進行大面積形變區(qū)排查,基于CORS成果提升BDS實時動態(tài)(RTK)差分精度,進行重點區(qū)域形變RTK差分監(jiān)測方法的可行性及精度.
基于CORS成果進行BDS實時動態(tài)差分計算的步驟如下:
步驟一:將基準站的觀測值和實測氣象元素按特定頻率匯集到數(shù)據(jù)中心,建立新的對流層延遲改正模型進行RTK定位;
步驟二:選擇經(jīng)驗模型,計算所有站前一段時間每天的全球經(jīng)驗模型和實測氣象元素的天頂對流層延遲(ZTD);
步驟三:建立高程、基于全球經(jīng)驗模型和基于實測氣象元素對流層延遲的模型之間的關系,進行新的修正模型;
步驟四:利用改正后的基準站和監(jiān)測站ZTD模型得到基準站觀測值的改正量;
步驟五:處理后的基準站觀測值播發(fā)給監(jiān)測站,監(jiān)測站進行實時差分解算出三維(3D)坐標.
短基線子集(SBAS)技術于2002被提出,是當前較主流的InSAR時間序列分析方法[9]. 該方法依據(jù)時空基線的條件,自由組合基線較短的影像對,組合成若干個短基線集合,再利用奇異值分解方法,將組合成的所有短基線集聯(lián)合起來求解,從而得到研究區(qū)域內整個觀測時間的形變時間序列和平均沉降速率.SBAS的基本步驟如下.
采集同一區(qū)域按照時間順序t0, ···,tN排列的N+1幅雷達影像,選取其中一幅影像作為主影像進行配準,則形成M幅差分干涉圖.
對于從影像tA和主影像tB(tB>tA) 時刻獲取的合成孔徑雷達(SAR)影像生成的第j幅差分干涉圖,在不考慮大氣延遲相位、殘余地形相位和噪聲位的情況下,其中任意像元 (x,y) 的干涉相位為
將式(1)中相位采用兩個獲取時間之間平均相位速度和時間的乘積的方式表示,則有
式中: Φ 為殘余地形相位,得第j幅干涉圖的相位值為
利用矩陣形式表示為
B是M×N矩陣,對B進行奇異值求解,獲取最小范圍下的視線向的相位平均速率,對各時間段的沉降速率進行時間域上的積分,得到視線向上時間序列的形變量,通過三角函數(shù),將它們轉為垂直向上監(jiān)測時間段內的累計沉降量和形變速率[10-12].
龍羊峽水電站位于青海省共和縣與貴德縣之間的黃河干流上,黃河穿越其中,右岸是險峻的茶納山,左岸是莽原. 屬高原干旱-半干旱大陸性氣候,庫區(qū)植被稀少,水位基本每年7月開始上升,11月份達到最高值,然后開始逐漸回落. 水庫右岸發(fā)育北東東向活動斷裂和北西向斷裂,巖土體類型為雙層結構巖土體,下部軟弱的層狀半固結巖土體,上部為松散沖洪積一般粘性土[8,13].
研究區(qū)及GNSS監(jiān)測點布設位置如圖1所示,其中紅色框內為龍羊峽水庫庫區(qū),紅色三角點為QHCORS龍羊峽站及GNSS監(jiān)測點.
圖1 研究區(qū)及BDS站點布設圖
在龍羊峽大壩下方布設了5個GNSS監(jiān)測點,將數(shù)據(jù)通過無線網(wǎng)絡傳實時輸?shù)絈H-CORS核心服務器上. 選取QH-CORS龍羊峽站為基準站,按照1.1節(jié)中的研究方法對采集到的BDS數(shù)據(jù)進行RTK差分定位解算. 其步驟一中特定頻率為1 Hz,步驟二中全球對流層模型為GPT2w,利用薩斯塔莫寧(Saastamoinen)模型計算龍羊峽站每天各時刻的ZTD. 圖2為流程及具體方法解算出2019年BDS數(shù)據(jù)RTK差分定位形變結果.
圖2 BDS監(jiān)測點實時動態(tài)解算流程圖
利用GNSS定位GAMIT軟件解算收集的與SAR和BDS數(shù)據(jù)時間相對應的13天的靜態(tài)數(shù)據(jù)[14],為了最大程度的削弱水域對解算成果的影響,利用商業(yè)軟件預算數(shù)據(jù)得到高精度的先驗坐標,進行潮汐改正.采用目前模擬對流層狀況精度最高的模型進行大氣改正,其主要設置為:
衛(wèi)星軌道:國際GNSS服務(IGS)快速精密星歷;
衛(wèi)星截止高度角:15°;
數(shù)據(jù)采樣間隔:10 s;
坐標約束:CORS水平方向給予5 cm,垂直方向給予10 cm的約束,監(jiān)測點給予10 m的約束;
潮汐改正:otl_FES2004.grid;
氣象模型:VMF1[15].
基線重復性可以反映出基線解內部精度,是評定數(shù)據(jù)處理質量的重要指標. 表1顯示基線重復性在各方向的固定誤差均小于2 mm,相對精度都達到10-7水平,滿足形變監(jiān)測毫米級精度的要求[14].
表1 基線向量重復性統(tǒng)計表
選擇覆蓋研究區(qū)域的2017-10-15—2019-10-17獲取的25景Sentinel-1SAR影像,C波段,分辨率5 m×20 m,精密軌道文件,30 m分辨率的數(shù)字高程模型數(shù)據(jù).
采用SBAS技術解譯獲取的全部影像,構建了一個所有時相配對、連接均勻的SBAS數(shù)據(jù)集合,如圖3所示,主要參數(shù):
圖3 解纏像對連接圖
極化方式:VV;
臨界基線最大百分比:3%;
最大時間距離:60;
多視視數(shù):1∶5;
解纏方式:3D,Delaunay MCF;
相關系數(shù)閾值:0.35;
大氣低通:1 200;
大氣高通:365;
形變投影:垂直.
采用的解纏方法可以很好地處理兩個相干性高的孤立區(qū)域,更好地抑制水汽對解譯的影響.
圖4通過解譯獲得了龍羊峽庫區(qū)2017—2019年在垂直方向上的形變速率,區(qū)域內大面積沉降區(qū)主要在入庫區(qū)、查西、旅游公司區(qū)域. 最嚴重沉降區(qū)域位于入庫區(qū),最大沉降速率52.48 mm/a;最大抬升點在查西,最大抬升速率43.60 mm/a. 其全部分布于庫區(qū)的南岸,與其他學者研究龍羊峽滑坡集中發(fā)育在南岸的論述相符合[16-17].
3.3.1 地表形變原因推斷
入庫區(qū)和旅游公司所在區(qū)域巖土主要以粘土和湖積粉砂為主. 黃河在龍羊峽內直流而下進入入庫區(qū),在其沖擊作用下將部分泥沙沖入庫區(qū),并在黃河下切影響下,沖刷南岸造成沉降現(xiàn)象. 旅游公司地勢平緩,在庫區(qū)水的浸潤、浪蝕作用下形成一個沉降量較小的沉降形變. 查西在黃河下切的影響下形成臨空面,馬蹄形地勢造成沉降和抬升的同時存在水位下降,山體內地下水位下降幅度相對滯后,地下水向庫區(qū)滲流,產(chǎn)生動水壓力,增大山體的下滑力[8,13].
3.3.2 成果精度驗證
所得成果的精度采用外符合精度評定. 外符合精度是以外部提供的參考值為比對基準,主要反映觀測值與參考值之間的偏差程度. 外符合精度計算公式如下:
式中:m外是 外符合精度; Δ 是參考值減去觀測值后的殘差;n為監(jiān)測點參與計算數(shù)據(jù)期數(shù)[15].
將靜態(tài)解算成果與BDS RTK監(jiān)測成果相減得到的殘差值代入式(5),得到各個點的外符合精度值.然后對各點外符合精度值求平均值得到監(jiān)測精度. 同理,得到InSAR監(jiān)測成果的精度.
由表2~3可知,BDS RTK差分監(jiān)測成果的精度為7.1 mm,InSAR監(jiān)測成果的精度為4.4 mm. 以GNSS靜態(tài)解算成果為參考值,兩種技術的監(jiān)測精度都達到毫米級.
表3 InSAR監(jiān)測成果精度計算表 mm
通過InSAR提取高空間分辨率的地表形變信息,GNSS獲取高時間分辨率的空間3D信息. 近年來國家對BDS系統(tǒng)的建設投入越來越大,而CORS的迅猛發(fā)展推動了其廣闊的應用前景. 本文以龍羊峽水庫為研究區(qū)域,通過兩種技術獲取了整個庫區(qū)垂直方向的地表形變信息,分析結果表明:在黃河影響下,2017—2019年龍羊峽庫區(qū)存在大面積的地表形變區(qū),最嚴重沉降區(qū)域位于入庫區(qū),最大沉降速率為52.48 mm/a;最大抬升點在查西,最大抬升速率為43.60 mm/a. 以GNSS靜態(tài)解算成果為參考值進行驗證,BDS RTK差分監(jiān)測成果的精度為7.1 mm,In-SAR監(jiān)測成果的精度為4.4 mm,都滿足形變監(jiān)測的精度要求. 采用InSAR進行大面積形變區(qū)篩選,再結合CORS利用BDS進行重點區(qū)域實時監(jiān)測的模式對于水庫全方位形變監(jiān)測是可行的. 由于InSAR對平面信息不敏感,本文只對垂直方向進行了驗證分析,但是庫區(qū)全方位形變監(jiān)測在進行垂直形變監(jiān)測的同時進行平面形變監(jiān)測也十分重要[18],今后結合GNSS空間3D信息及水汽反演成果,對深入研究InSAR水汽干擾消除及3D形變場的建立很有意義.