基于高分辨率PSP-InSAR技術(shù)的深圳市長嶺皮水庫壩體形變監(jiān)測應(yīng)用研究

黃 松，倪紹文，習(xí)樹峰，董 凱

（深圳市城市公共安全技術(shù)研究院有限公司，廣東深圳，518000）

0 前言

為了解決徑流在時(shí)間上和空間上的重新分配問題，充分開發(fā)利用水資源，深圳市內(nèi)修建了大量水庫工程。但是，水庫工程建設(shè)、水位調(diào)節(jié)或其他相關(guān)事件會對水庫壩體邊坡及周圍區(qū)域的穩(wěn)定性產(chǎn)生影響，進(jìn)而導(dǎo)致安全隱患。一般而言，在重大安全事故發(fā)生前，目標(biāo)都會存在一定程度的微小形變，對水庫壩體邊坡及其周圍區(qū)域進(jìn)行合理有效的形變監(jiān)測已經(jīng)成為水利行業(yè)重點(diǎn)關(guān)注的問題。

相比于傳統(tǒng)的接觸式測量技術(shù)（水準(zhǔn)儀測量和GPS測量等），星載合成孔徑雷達(dá)干涉測量（Spaceborne-InSAR）技術(shù)屬于非接觸測量范疇，存在其他技術(shù)不具備的優(yōu)勢。它利用兩幅或多幅SAR圖像中的相位數(shù)據(jù)獲取地表形變信息，每景SAR圖像的覆蓋范圍可達(dá)數(shù)百甚至上千平方公里，形變測量精度能達(dá)到mm量級，且不受光照和天氣條件限制，能實(shí)現(xiàn)全天時(shí)、全天候?qū)Φ赜^測。因此，星載InSAR技術(shù)特別適合對水庫進(jìn)行定期大面積形變監(jiān)測[1]。

目前，最具應(yīng)用潛力的InSAR形變測量技術(shù)是永久散射體（PSInSAR）技術(shù)，它先在多幅SAR圖像中選出能長時(shí)間保持高相干性的散射點(diǎn)，并將其定義為PS（Persistent Scatterer）點(diǎn)，然后基于PS點(diǎn)的相位數(shù)據(jù)，獲取目標(biāo)區(qū)域在監(jiān)測時(shí)間范圍內(nèi)的形變速率信息和形變演化歷史信息[2-3]。2009年，意大利Constantini等人在此基礎(chǔ)上提出了PSP-InSAR（PS PairsInSAR）技術(shù)，能在SAR圖像序列中選出更多的PS點(diǎn)，且形變反演精度更高。PSP-InSAR技術(shù)的提出進(jìn)一步拓展了InSAR技術(shù)的實(shí)用能力[4-5]。

COSMO-SkyMed星座是由意大利政府投資，目前世界上唯一有4顆衛(wèi)星同軌的軍民兩用SAR衛(wèi)星星座[6]。它發(fā)射X波段的微波信號，最高分辨率能達(dá)到1 m，且具有多種數(shù)據(jù)獲取模式，最短能在6 h內(nèi)提供目標(biāo)區(qū)域的信息。COSMO-SkyMed系統(tǒng)在災(zāi)害監(jiān)測、科研、商業(yè)和國防等方面的應(yīng)用具有巨大的社會和經(jīng)濟(jì)效益以及發(fā)展?jié)摿Α?/p>

以COSMO-SkyMed數(shù)據(jù)為輸入，利用PSP-In-SAR技術(shù)對深圳市全市大部分水庫邊坡及其周邊區(qū)域開展了形變監(jiān)測。在長嶺皮水庫壩體邊坡形變的重點(diǎn)分析過程中，通過結(jié)合水位變化數(shù)據(jù)，建立了形變隨水位變化和時(shí)間變化的綜合模型。在此基礎(chǔ)上，根據(jù)參數(shù)估計(jì)理論，成功分解出水庫壩體不同坡面由水位變化導(dǎo)致的形變和由自身沉降引起的形變。通過對數(shù)據(jù)的綜合分析，證實(shí)了In-SAR技術(shù)在水庫穩(wěn)定性評估方面的優(yōu)勢。

1 PSP-InSAR形變監(jiān)測算法

圖1顯示了InSAR形變測量技術(shù)的幾何模型。M和S分別為形變發(fā)生前后SAR衛(wèi)星的位置，其距離為空間基線B。在形變發(fā)生前，目標(biāo)點(diǎn)位于位置A，形變發(fā)生后，目標(biāo)點(diǎn)移動到位置A'。當(dāng)從干涉相位中剔除模擬的地形相位后，目標(biāo)點(diǎn)在視線（LOS）方向的形變量Δr與形變相位φdef的關(guān)系可以表示為：式（1）中，λ表示雷達(dá)信號的波長。

Ferretti等在1999年提出了PSInSAR方法。PSInSAR方法首先在SAR圖像中選出那些在長時(shí)間范圍內(nèi)能保持高相關(guān)性的目標(biāo)點(diǎn)，并定義為PS點(diǎn)。通過數(shù)據(jù)處理，最終獲取所有PS點(diǎn)的形變信息和三維位置信息。

Persistent Scatterer Pair（PSP）方法是傳統(tǒng)PSIn-SAR算法的進(jìn)一步升級，它的的核心思想是定義和分析PS點(diǎn)對。PS點(diǎn)對的聯(lián)合分析能夠降低空間相關(guān)性誤差（如大氣相位誤差）對形變反演結(jié)果的影響。因此，這種方法能夠克服傳統(tǒng)PSInSAR算法的限制，并且獲取密度更高的PS點(diǎn)。

圖1 InSAR形變測量技術(shù)的幾何模型Fig.1 Geometric model of InSAR deformation measurement

PSP算法的核心步驟是迭代建立PSP網(wǎng)格，并估計(jì)PS對之間的相對形變速率和相對高程差。由于形成PS對的兩個(gè)PS點(diǎn)距離較近，在干涉相位中的大氣相位項(xiàng)Δφatmo和軌道誤差相位項(xiàng)Δφtrack基本相同，那么這兩個(gè)PS點(diǎn)的干涉相位差可以建模為：

式中，Ti表示第i幅SAR圖像相對于參考SAR圖像的時(shí)間；δva表示第a個(gè)PS對中兩個(gè)PS點(diǎn)的相對形變速率；Bi表示第i幅干涉圖的有效基線；δha表示第a個(gè)PS對中兩個(gè)PS點(diǎn)的相對高程差；r表示目標(biāo)點(diǎn)的斜距；θ表示下視角；εa,i表示噪聲和未建模的誤差。一般來說，如果εa,i較小，那么這個(gè)像素點(diǎn)對就可以被認(rèn)為是PS對，并稱之為PSP。

為了從數(shù)學(xué)模型上判斷PSP，可以先定義第a個(gè)PSP的時(shí)間相關(guān)系數(shù)：

式中，ωa,i表示第a個(gè)像素點(diǎn)在第i幅干涉圖像中所對應(yīng)的權(quán)重值，簡單情況下可全設(shè)置為1，在處理過程中，先對未知參數(shù)δva和δha進(jìn)行最佳估計(jì)，然后計(jì)算γa。同時(shí)，設(shè)定時(shí)間相關(guān)系數(shù)門限γThreshold，當(dāng) γa＞γThreshold時(shí)，這個(gè)像素點(diǎn)對可以被判定為PS對。

2 基于InSAR技術(shù)的水庫形變分析基本原理

基于水利工程的基礎(chǔ)知識，水庫水位的變化會導(dǎo)致水庫邊坡壓力發(fā)生改變，進(jìn)而會引起邊坡及周邊區(qū)域發(fā)生一定程度的形變[7-8]。根據(jù)經(jīng)驗(yàn)?zāi)Ｐ停伤蛔兓鸬男巫兛梢员硎緸椋?/p>

式中，Δhwater表示水位的相對變化值；k表示水位變化到形變的轉(zhuǎn)換因子。

在同時(shí)考慮水庫存在由沉降和由水位變化引起形變的情況下，水庫邊坡上測得的PS點(diǎn)形變量可以建模為：

式中，PS點(diǎn)的沉降采用二階運(yùn)動模型，其中d0表示初始形變值；v0表示初始形變速率；a表示形變加速度；t表示相對于初始時(shí)刻的時(shí)間；e則表示未建模的形變和形變測量誤差。在式中，未知的參數(shù)包括d0、v0、a和k，那么，最終水庫形變分析的問題就可以轉(zhuǎn)化為這四個(gè)未知參數(shù)的估計(jì)問題。同時(shí)，為了簡化分析，在參數(shù)估計(jì)后，PS點(diǎn)的平均形變速率可以近似為：

式中，T表示InSAR技術(shù)整體的監(jiān)測時(shí)間長度。

3 深圳市水庫基本概況和實(shí)測數(shù)據(jù)介紹

截至2015年底，深圳市共有蓄水工程168座，其中大型水庫1座、中型水庫13座、小（1）型水庫65座、小（2）型水庫89座，總控制集雨面積593.95 km2，總庫容7.76億m3。水庫安全一直是深圳市重點(diǎn)關(guān)注的問題。

為了對深圳全市的水庫開展形變監(jiān)測，基于COSMO-SkyMed系統(tǒng)獲取的深圳西部51幅和深圳東部45幅3 m分辨率條帶模式下的SAR干涉圖像序列，利用PSP-InSAR技術(shù)，獲取了深圳全市2013年9月～2016年9月3年內(nèi)11 668 810個(gè)PS點(diǎn)的三維位置信息、形變速率信息和形變歷史信息。其中，數(shù)據(jù)的基本信息和衛(wèi)星的照射范圍分別如表1和圖2所示。

InSAR技術(shù)能獲取長時(shí)間范圍內(nèi)水庫及周邊區(qū)域的高精度形變信息，最終為水庫邊坡及周邊目標(biāo)的形變機(jī)理研究和安全評估提供基礎(chǔ)數(shù)據(jù)支持。

表1 深圳InSAR數(shù)據(jù)的基本參數(shù)Table 1 Basic parameters of Shenzhen InSAR data

圖2 COSMO-SkyMed系統(tǒng)照射示意圖和影像覆蓋范圍Fig.2 The irradiation range of COSMO-SkyMed system and the image coverage area

4 深圳市長嶺皮水庫壩體邊坡形變監(jiān)測分析

長嶺皮水庫位于深圳市南山區(qū)大沙河支流長嶺皮河的上游，始建于1977年，1981年竣工建成，1985年投入運(yùn)行，原為總庫容786萬m3的?。?）型水庫。2007年8月，啟動水庫加固及擴(kuò)容工程，主體工程已于2010年竣工，擴(kuò)建后總庫容為1 754萬m3，升級為中型水庫。長嶺皮水庫主要任務(wù)為供水、調(diào)蓄及防洪，供水對象為紅木山水廠及高峰水廠，防護(hù)對象為大沙河上游段的大學(xué)城校區(qū)、田寮村、福光村、大馬村等。

一般而言，在水庫工程的建設(shè)期及竣工后一段時(shí)間內(nèi)，壩體及周邊區(qū)域都會存在一定程度的形變。因此，考慮到長嶺皮水庫的擴(kuò)建工程數(shù)年前才竣工，該水庫的形變監(jiān)測就成為研究建筑工程對水庫及周邊影響的熱點(diǎn)案例。

圖3顯示了InSAR監(jiān)測時(shí)間范圍內(nèi)長嶺皮水庫的水位變化數(shù)據(jù)。可以看出，在2015年1月后，水庫水位每年內(nèi)都存在較大幅度的波動。根據(jù)PSPInSAR的計(jì)算結(jié)果，在長嶺皮水庫的邊坡上共獲取570個(gè)PS點(diǎn)。采用式（5）的模型進(jìn)行參數(shù)估計(jì)后，獲取了每個(gè)PS點(diǎn)的水位轉(zhuǎn)換因子、平均形變速率和形變加速度，其結(jié)果分別如圖4～6所示。同時(shí)，為了更好地分析壩體不同邊坡位置處的形變，將長嶺皮水庫主壩和2號壩分為4個(gè)部分。

圖3 長嶺皮水庫水位變化數(shù)據(jù)Fig.3 Graph of water level of Changlingpi reservoir

圖4 長嶺皮水庫大壩的水位轉(zhuǎn)換因子監(jiān)測結(jié)果Fig.4 Monitoring result of the water level conversion factor of Changlingpi reservoir

圖4顯示了整個(gè)水庫大壩水位轉(zhuǎn)換因子的形變監(jiān)測結(jié)果。通過對比可知，1號區(qū)域的形變受水位變化的影響非常明顯，大多數(shù)PS點(diǎn)的水位轉(zhuǎn)換因子分布在0.8～1.1 mm/m之間，即當(dāng)水位變化1 m時(shí)，這些PS點(diǎn)在LOS方向的形變可達(dá)0.8～1.1 mm。同時(shí)，3號區(qū)域和4號區(qū)域部分PS點(diǎn)的水位轉(zhuǎn)換因子分布在0.2～0.6 mm/m之間，表明這些區(qū)域也會隨水位變化而發(fā)生一定程度的形變。相反，2號區(qū)域處大部分PS點(diǎn)的水位轉(zhuǎn)換因子較小，說明這塊區(qū)域的形變受水位變化的影響較小。根據(jù)初步分析，這類現(xiàn)象可能與邊坡的朝向相關(guān)。

圖5顯示了長嶺皮水庫大壩形變加速度圖。形變加速度為負(fù)值表示形變速率沿背離雷達(dá)方向增大，形變加速度為正值表示形變速率沿朝向雷達(dá)方向增大。從結(jié)果可以看出，除了3號區(qū)域背水坡的部分PS點(diǎn)，整個(gè)水庫大壩大多數(shù)PS點(diǎn)的形變加速度都小于2 mm/yr2，說明邊坡的大部分區(qū)域不存在明顯形變加速的趨勢。

圖5 長嶺皮水庫邊坡的形變加速度監(jiān)測結(jié)果Fig.5 Monitored deformation acceleration rate of Changlingpi dam slope

圖6 長嶺皮水庫邊坡的形變速率監(jiān)測結(jié)果Fig.6 Monitored deformation rate of Changlingpi dam slope

圖6顯示了長嶺皮水庫大壩的平均形變速率。形變速率為負(fù)值表示PS點(diǎn)沿遠(yuǎn)離雷達(dá)的方向運(yùn)動，形變速率為正值表示PS點(diǎn)沿朝向雷達(dá)的方向運(yùn)動。從監(jiān)測結(jié)果可知，這四部分區(qū)域都存在一定程度的形變。通過對比分析可知，在3號區(qū)域處的形變速率和形變區(qū)域都最大，最大形變速率達(dá)到7.4 mm/yr。

為了更好地分析大壩不同位置處的形變演化歷史，在1號和3號區(qū)域處分別選擇了A和B兩個(gè)PS點(diǎn)，其位置如圖6所示。圖7和圖8則顯示了這兩個(gè)PS點(diǎn)的實(shí)測形變歷史和模型分析形變歷史。其中，模型形變是運(yùn)動形變與水位變化形變之和。通過對比兩個(gè)PS點(diǎn)的實(shí)測形變曲線和模型形變曲線，兩者在形變趨勢上非常相似，也從一定程度上驗(yàn)證了式（5）所對應(yīng)模型的準(zhǔn)確性。同時(shí)，將這兩個(gè)PS點(diǎn)的形變歷史進(jìn)行對比，可以看出A點(diǎn)由水位變化引起的形變更加明顯，而B點(diǎn)由自身運(yùn)動引起的形變速率更大，呈現(xiàn)出更明顯的運(yùn)動形變特征。

圖7 長嶺皮水庫A點(diǎn)的形變演化歷史Fig.7 The deformation history of monitoring point A of Changlingpi reservoir

圖8 長嶺皮水庫B點(diǎn)的形變演化歷史Fig.8 The deformation history of monitoring point B of Changlingpi reservoir

5 結(jié)語

以COSMO-SkyMed衛(wèi)星采集的SAR數(shù)據(jù)為輸入，利用PSP-InSAR形變測量技術(shù)，對深圳市內(nèi)大部分水庫開展了時(shí)間跨度為3年（2013年9月～2016年9月）的形變監(jiān)測。在對長嶺皮水庫大壩進(jìn)行重點(diǎn)形變分析的過程中，結(jié)合水位變化數(shù)據(jù)，建立了目標(biāo)點(diǎn)形變隨水位變化和時(shí)間變化的綜合模型。在此基礎(chǔ)上，采用信號參數(shù)估計(jì)理論，成功從實(shí)測形變演化歷史數(shù)據(jù)中分解出由水位變化引起的形變和由目標(biāo)自身沉降引起的形變。同時(shí)，通過對比分析不同壩體坡面的形變結(jié)果，初步確認(rèn)了長嶺皮水庫大壩受水位變化影響較大的坡面和自身沉降較為明顯的坡面。下一步的工作可能需要結(jié)合地面測量手段，對形變較大的區(qū)域開展重點(diǎn)監(jiān)測，并綜合分析形變原因。

基于COSMO-SkyMed數(shù)據(jù)的InSAR形變測量技術(shù)能在觀測區(qū)域內(nèi)獲取高密度的測量點(diǎn)，完全有能力對目標(biāo)進(jìn)行高精度的形變監(jiān)測。在未來的應(yīng)用過程中，它能對水庫及周邊區(qū)域進(jìn)行安全普查，并與傳統(tǒng)測量手段相結(jié)合，共同實(shí)現(xiàn)對水庫區(qū)的安全評估和災(zāi)害預(yù)警。 ■

[1]Massonnet D,Rossi M,Carmona C,Adagna F,Peltzer G,Feigl K,Rabaute T.The displacement field of the Landers earthquake mapped by radar interferometry[J].Nature,1993,64(8)：138-142.

[2]Ferretti A,Prati C,Rocca F.Nonlinear subsidence rate estimation using permanent scatterers in differential SAR interferometry[J].IEEE Trans.On Geosci.Remote Sensing,2000,38(5)：2202-2212.

[3]Ferretti A,Prati C,Rocca F.Permanent scatterers in SAR interferometry[J].IEEE Trans.Geosci.Remote Sens.,2001,39(1)：8-20.

[4]Costantini M,Falco S,Malvarosa F,Minati F.A new method for identification and analysis of persistent scatterers in series of SAR images[C].Proc.Int.Geosci.Remote Sens.Symp.(IGARSS),Boston,MA,USA,2008,2：449-452.

[5]Costantini M,Falco S,Malvarosa F,Minati F,Trillo F,Vecchioli F.Persistent Scatterer Pair Interferometry：Approach and Application to COSMO-SkyMed SAR Data[J].IEEE JOURNAL OF SELECTED TOPICS IN APPLIED EARTH OBSERVATIONS AND REMOTE SENSING,2014,7(7)：2869-2879.

[6]F.Covello,F.Battazza,A.Coletta,E.Lopinto,C.Fiorentino,L.Pietranera,G.Valentini,S.Zo_oli.COSMO-SkyMed an existing opportunity for observing the Earth[J].Journal of Geodynamics,2010,49(3-4)：171-180.

[7]Milan Lazecky,Daniele Perissin,Wang Zhiying,Lei Ling,and Qin Yuxiao.Observing Dam's Movements with Spaceborne SAR Interferometry[J].Engineering Geology for Society and Territory,2014：131-136.

[8]JJMD Sousa,M Lazecky,I Hlavacova,M Bakon,G Patricio.Satellite SAR Interferometry for Monitoring Dam Deformations in Portugal[J].Dam World,2015.