哨兵-1A合成孔徑雷達的深圳圍填海區(qū)域地面沉降監(jiān)測

周 呂，施顯健，任 超，黃遠林，朱子林

(1.桂林理工大學測繪地理信息學院，桂林 541004；2.城市空間信息工程北京市重點實驗室，北京 100038；3.廣西空間信息與測繪重點實驗室，桂林 541004；4.北部灣大學資源與環(huán)境學院，欽州 535011)

改革開放以來，深圳地區(qū)的經濟得到了快速的發(fā)展，但隨著城市的快速擴張和土地需求的迅速增長，深圳土地資源緊缺問題日益凸顯，其中圍填海成為深圳緩解土地供應的重要途徑之一[1]。資料顯示，僅在1990—2015年，深圳就通過圍填海獲得了超過80 km2的土地面積[2]。然而深圳圍填海區(qū)域的地層結構性質復雜，隨著圍填海區(qū)域地面的逐漸鞏固和地面建筑負荷的持續(xù)作用，極易引發(fā)地面沉降，進而威脅城市發(fā)展和人民生命財產安全[3]。因此，及時掌握深圳圍填海區(qū)域在填海造地工程結束后的地面沉降情況具有重大現實意義。

傳統(tǒng)的地面沉降監(jiān)測方法一般是利用全球定位系統(tǒng)(global positioning system，GPS)測量方法或水準測量方法通過布設地面監(jiān)測網的方式進行地面形變監(jiān)測[4]。利用上述方法采集到的數據的精度較高，但是需要投入較大人力和時間成本。近年來快速興起的合成孔徑雷達干涉測量技術(interferometric synthetic aperture radar, InSAR)具有大范圍[5]、高精度[6]、長時間監(jiān)測[7]等優(yōu)點，為大面積地面沉降監(jiān)測提供新思路[8]。然而InSAR技術容易受到失相關、大氣及軌道誤差的影響[9]。為了克服這些問題，中外眾多學者在合成孔徑雷達差分干涉測量技術(differential interferometric synthetic aperture radar, D-InSAR)基礎上提出了時間序列合成孔徑雷達(time-series interferometric synthetic aperture radar, TS-InSAR)技術[10-11]，通過解析在時間序列中保持穩(wěn)定散射特性的像素點來提取地面形變速率[12]，較好地消除了失相關、大氣及軌道誤差對地面形變監(jiān)測的影響[13]，實現了較高精度的地面形變時序監(jiān)測[14]。

深圳圍填海區(qū)域的SAR數據相對緊缺和商用SAR衛(wèi)星數據定制價格高昂的問題在某種程度上制約了InSAR技術用于圍填海區(qū)域的地面形變監(jiān)測。哨兵-1A(Sentinel-1A)衛(wèi)星是歐洲航天局(ESA)哥白尼計劃(global monitoring for environment and security, GMES)的新一代對地觀測衛(wèi)星，于2014年4月3日發(fā)射，具有影像分辨率高、全球大范圍監(jiān)測、重訪周期短和數據公開等優(yōu)勢[15]，可為深圳圍填海區(qū)域的地面形變時序監(jiān)測提供較好的數據支持。綜合上述分析，現利用TS-InSAR技術和深圳圍填海區(qū)域2016年12月6日—2019年10月4日的18景Sentinel-1A數據提取深圳圍填海區(qū)域的地面形變結果，以此為基礎重點分析深圳圍填海區(qū)域的地面沉降規(guī)律及沉降原因。

1 研究區(qū)和數據源

深圳位于中國華南地區(qū)、珠江口東岸，西瀕伶仃洋和珠江口，東臨大鵬灣和大亞灣。改革開放后憑借著瀕臨香港、背靠珠三角的天然地理優(yōu)勢和其獨有的經濟特區(qū)政策，深圳經濟迅猛發(fā)展，現已成為中國超一線城市中最年輕的城市以及世界一線城市。然而，隨著城市規(guī)模的急劇擴張，深圳土地資源緊缺問題日益凸顯。為了緩解土地供應緊張問題，深圳開始了較大規(guī)模的圍填海，圖1為深圳圍填海區(qū)域概況。

圖1 深圳圍填海區(qū)域概況

基于TS-InSAR技術和深圳填海區(qū)2016年12月6日—2019年10月4日的18景Sentinel-1A數據，利用ESA發(fā)放的精密定軌(precise orbit determination，POD)去除軌道誤差，利用航天飛機雷達地形測繪使命(shuttle radar topography mission，SRTM)數字高程模型(digital elevation model，DEM)Version4 90 m去除地形相位，提取了深圳圍填海區(qū)域的地面沉降結果。表1為Sentinel-1A SAR數據參數。

表1 Sentinel-1A數據參數

2 TS-InSAR原理

首先，基于研究區(qū)的18景Sentinel-1A SAR影像構造時空基線，圖2為試驗的垂直基線分布圖。

圖2 垂直基線分布

根據生成的時空基線，從研究區(qū)(t0,t1,…,tN)時刻采集的N+1景SAR影像中篩選出主影像，將主影像和剩余影像配準并生成M景干涉圖滿足

(1)

對tA時刻采集的SAR影像和tB時刻采集的SAR影像(tA

δφj(x,r)=φ(tB,x,r)-φ(tA,x,r)≈

φnoise,j(x,r)+φatm,j(x,r)+

φdef,j(x,r)+φtopo,j(x,r)

(2)

式(2)中：x、r分別為像素點的方位向坐標、距離向坐標；φ(tB,x,r)、φ(tA,x,r)分別為tB、tA時刻采集的SAR影像上某一像素點(x,r)的相位；φnoise,j(x,r)為噪聲相位；φatm,j(x,r)為大氣相位的誤差；φdef,j(x,r)為tA—tB時刻的衛(wèi)星視線向(line of sight, LOS)形變相位；φtopo,j(x,r)為外部DEM引入的地形相位誤差。

大氣相位的誤差計算公式為

φatm,j(x,r)=φatm,j(tB,x,r)-φatm,j(tA,x,r)

(3)

式(3)中：φatm,j(tB,x,r)、φatm,j(tA,x,r)分別為tB、tA時刻采集的SAR影像上的大氣相位分量。

衛(wèi)星視線向形變相位計算公式為

(4)

式(4)中：λ為波長；d(tB,x,r)、d(tA,x,r)分別為tB、tA時刻到t0時刻間衛(wèi)星LOS累積形變量。

外部DEM引入的地形相位誤差的誤差可表示為

(5)

式(5)中：B⊥j為干涉圖j的垂直基線；R為衛(wèi)星與目標點間的距離；Δz為DEM誤差；θ為視角(SAR波束和地面垂直直線形成的夾角)。

此外，式(2)干涉相位δφj(x,r)中的噪聲相位分量、大氣延遲相位誤差分量和地形相位誤差分量需要精確估算出并排除，則有

Aφ=δφ

(6)

式(6)中：A為M×N系數矩陣(M為干涉圖數量，N為SAR影像數量)，且?j=1,2,…,M，向量φT=[φ(t1),φ(t2),…,φ(tN)]由每景SAR影像上高相干點目標的相位值組成，向量δφT=[δφ(t1),δφ(t2),…,δφ(tN)]由各個差分干涉圖的解纏相位組成。

用兩景SAR影像之間的平均相位速率替換相位值以接解算影像上各高相干點目標的位移形變速率，式(6)可表示為

Bv=δφ

(7)

若B秩虧(M

3 結果與分析

3.1 形變結果與分析

利用TS-InSAR技術和深圳圍填海區(qū)域2016年12月6日—2019年10月4日的18景Sentinel-1A數據提取的地面形變結果如圖3所示，圖中正值表示地面在垂直方向上隆起，負值表示地面在垂直方向上下沉。由圖3的地面形變結果可以發(fā)現，2016—2019年，深圳西海岸的地面形變速率為-19.85～15.60 mm/a，深圳東部BA地區(qū)(寶安區(qū))和西部FT地區(qū)(福田區(qū))的地面沉降相對穩(wěn)定，于NS地區(qū)(南山區(qū))監(jiān)測到大面積的地面沉降。

圖3 深圳圍填海區(qū)域的地面形變速率圖

為檢測沿海岸線的形變情況，從圖3中沿AA′提取了地面形變剖面圖，如圖4所示。由圖4可知，AA′沿線地面自西向東呈現逐漸下沉的趨勢，AA′中最大沉降點A1、A2(分別靠近南山區(qū)與寶安區(qū)、福田區(qū)的邊界)的垂直沉降分別為-13.50、-15.03 mm。可以發(fā)現，沿海岸線的沉降漏斗主要分布在南山區(qū)，則重點對南山區(qū)進行沉降機理的分析。

圖4 沿A-A′提取的地面形變剖面圖

3.2 精度分析

由于沒有獲取到外部水準數據進行輔助驗證，因此，為了更好地評估Sentinel-1A InSAR形變反演的內部精度，通過計算速度的線性擬合偏差得到形變標準差，并對形變標準差進行統(tǒng)計分析，圖5給出本次試驗平均形變速率的標準差分布。若點目標表現出較強的非線性運動(即表現為高標準差值)，則相對于線性模型會產生較大的殘差。由圖5可以發(fā)現Sentinel-1A InSAR形變反演得出的平均形變速率的標準差為1.39 mm/a，點目標點最大標準差達到9.26 mm/a，96.75%的形變點目標的標準差小于2.25 mm/a。上述分析表明Sentinel-1A InSAR形變反演具有較高的精度和可靠性。

圖5 形變速率標準差的分布

3.3 南山區(qū)沉降機理分析

從地理上看，NS(南山區(qū))瀕臨南海，在亞熱帶海洋性季風氣候的作用下，四季溫潤、年均降雨量多，再加上該地區(qū)平地、山丘、臺地相間，地勢北高南低，地面容易產生向向南傾斜下沉的趨勢。發(fā)展上看，南山區(qū)不僅是深圳的科研、娛樂、教育以及體育中心，還是深圳新的經濟增長極，自國務院2016年5月將深圳市南山區(qū)設為中國首批雙創(chuàng)示范基地，越來越多的企業(yè)開始入駐。隨著圍填海區(qū)域經濟建設的快速推進，容易加速誘發(fā)南山區(qū)的不均勻沉降。一般的，形變點目標臨近區(qū)域具有相同的形變趨勢。為分析南山區(qū)地面沉降的時間演化及其與降水的關系，從圖3中選取形變點目標a、b、c、d與CW(赤灣)點監(jiān)測的月平均降雨量進行對比，圖6為結果分析圖。由圖6可以發(fā)現，所選形變點目標的沉降時間序列隨降雨量的季節(jié)變化而非線性遞減：平均降水量增加的時間里，可能是由于降水有效地補充了地下水，所選點目標的地面沉降速率較為明顯地減緩了。從空間位置上看，靠近海岸的形變點目標明顯受降雨的影響明顯大于其他點，這可能與深圳的沉積層厚度有關。

圖6 地面沉降的時間演化及其與降水的關系

圖7為深圳圍填海區(qū)域的沉積層厚度分布圖。從時空演變的角度來看，深圳是珠江三角洲流水的長期沉積作用而匯集成的沖積平原。深圳圍填海區(qū)域廣泛分布著第四紀的海洋沉積物，這些沉積物形成的沉積層具有孔隙率大、含水量高、壓縮性高等特點，因而沉積層厚度大的地面在施工建設和后期土體自固結過程更易發(fā)生形變。通過對比分析圖7和圖3，可以發(fā)現沉積層厚度為10～20 m的地區(qū)的地面沉降相比沉積層厚度為3～10 m的地區(qū)要明顯。

圖7 深圳圍填海區(qū)域沉積層厚度分布

4 結論

利用TS-InSAR技術和深圳填海區(qū)2016年12月6日—2019年10月4日的18景Sentinel-1A數據，提取了深圳圍填海區(qū)域的地面沉降結果，在此基礎上重點分析了圍填海區(qū)域的地面沉降規(guī)律和地面沉降原因。結果顯示：深圳西海岸的地面形變速率為-19.85～15.60 mm/a，深圳東部寶安區(qū)和西部福田區(qū)的地面沉降相對穩(wěn)定，于南山區(qū)監(jiān)測到大面積的地面沉降；與年均降雨量的對比分析顯示了南山區(qū)的地面沉降存在明顯的季節(jié)變化規(guī)律，表明該地區(qū)的地面沉降與降水有關；靠近海岸的形變點目標明顯受降雨的影響明顯大于其他點，這可能與深圳的沉積層厚度有關。為了更好地探究圍填海區(qū)域地面沉降，下一步，將采集填海區(qū)地下水數據和沉積層鉆孔測試數據進行更進一步的研究。