江培華，王金洋，唐 偉，郝延錦

(1. 華北科技學(xué)院 建筑工程學(xué)院，北京 東燕郊 065201；2. 中國礦業(yè)大學(xué)(北京) 地球科學(xué)與測繪工程學(xué)院，北京 100083)

0 引言

地面沉降是一種在自然或人為因素下導(dǎo)致地下松散巖層固結(jié)壓縮，并引起地面一定范圍內(nèi)標(biāo)高降低的地質(zhì)現(xiàn)象[1]。地面沉降的累積和空間不均勻分布會對地面產(chǎn)生破壞，對地表植被、水體、建筑產(chǎn)生威脅，影響人類生活生產(chǎn)的安全，演變成為一種地質(zhì)災(zāi)害[2-4]。天津市是中國北方最大的港口城市，地表沉降已在市區(qū)內(nèi)多個(gè)區(qū)域發(fā)生。1959-1997年，天津市的沉降情況可以分為三個(gè)階段[5]：(1)1959-1977年處于緩慢沉降階段，沉降速率為50 mm/year；(2)1978-1982年加速沉降，沉降速率為100 mm/year；(3)1983-1997年沉降減緩，沉降速率為20～30 mm/year。位于濱海新區(qū)的漢沽2000-2007年的平均沉降速率達(dá)35～45 mm/year，1959-2007年的最大累計(jì)沉降量達(dá)3.08 m，是當(dāng)前天津市四大沉降中心之一[6]。為了監(jiān)測地表沉降情況，國內(nèi)外都建立了包含眾多手段的地面沉降監(jiān)測網(wǎng)絡(luò)。地面沉降的監(jiān)測技術(shù)包括水準(zhǔn)測量、全球定位系統(tǒng)(GPS)和合成孔徑干涉雷達(dá)(InSAR)等。研究人員采用InSAR時(shí)間序列分析技術(shù)，對2007年至2009年間Envisat衛(wèi)星C波段ASAR影像與ALOS衛(wèi)星L波段PALSAR影像數(shù)據(jù)處理，來提取天津地表形變，結(jié)果顯示天津西部地區(qū)沉降主要受地下水的開采影響，而東部地區(qū)的沉降情況則受天津南、天津北和大寺斷裂控制，同時(shí)也與含水層Aquifer Ⅱ和Aquifer Ⅲ有關(guān)[7]。2018年研究人員利用多時(shí)相InSAR技術(shù)對2016年1月9日至2017年6月8日25景升軌Sentinel-1A 衛(wèi)星影像處理得出，天津市的地面形變從市區(qū)到郊區(qū)逐漸增加，不同的地質(zhì)條件也對沉降產(chǎn)生了不同的影響。此外，居住區(qū)和工業(yè)區(qū)之間明顯的沉降差異表明，即使在同一地質(zhì)區(qū)域，抽取地下水的深度差異也可能導(dǎo)致不同的沉降率[8]。2019年研究人員對2015年5月至2016年6月的25景Sentinel-1雷達(dá)影像，采用小基線集(Small baseline subset,SBAS)技術(shù)，結(jié)合城市規(guī)劃研究天津市濱海新區(qū)地表沉降特征，其中2008年存在大量的圍海造陸區(qū)域，軟土固結(jié)需要一定的周期，圍海區(qū)域大量的軟土沉積引發(fā)較大的地面沉降[9]。 針對天津市的沉降情況，目前的研究大都集中在市區(qū)范圍內(nèi)，而近年來沿海地區(qū)進(jìn)行填海造陸擴(kuò)張土地，地面沉降情況尤為突出。同時(shí)，InSAR獲取的形變信息為地表實(shí)際三維形變投影到雷達(dá)視線向(Line-of-sight,LOS)的一維形變結(jié)果。分解得到的垂直向形變值更能直觀地反映地表沉降情況，而Sentinel-1衛(wèi)星所獲得的升、降軌數(shù)據(jù)為視線向沉降量的三維分解提供了便利。

本文結(jié)合了天津市填海區(qū)的填海造陸的歷史及海岸線的變化，利用SBAS-InSAR時(shí)序分析技術(shù)對Sentinel-1衛(wèi)星SAR數(shù)據(jù)進(jìn)行分析，獲取天津沿海地面沉降時(shí)空分布信息，并結(jié)合三個(gè)GPS觀測站點(diǎn)的數(shù)據(jù)進(jìn)行對比分析。最后，從漢沽、臨港、南港三個(gè)區(qū)域及其海岸線的變化來分析研究沿海區(qū)的地面沉降時(shí)空演變特征。結(jié)論對天津市的地面沉降防治工作具有實(shí)際意義。

1 研究區(qū)與SAR數(shù)據(jù)

1.1 研究區(qū)概況

濱海新區(qū)位于中國華北平原北部，天津市中心的東邊，東臨渤海，北部與河北省唐山市相鄰，南部則與河北省黃驊市相鄰。海岸線主體長約153公里，陸地面積達(dá)2270平方公里，海域面積約3000平方公里，如圖1所示。本文將濱海新區(qū)分為漢沽區(qū)、臨港工業(yè)區(qū)和南港工業(yè)區(qū)進(jìn)行分析。漢沽區(qū)面積441.5平方公里，其中海涂面積69.85平方公里，有著“千年鹽城”之稱，航運(yùn)發(fā)達(dá)。臨港工業(yè)區(qū)是國家級重型裝備制造基地，總規(guī)劃面積約200平方公里，圍海造地面積約67平方公里。南港工業(yè)區(qū)是世界級石化產(chǎn)業(yè)和港口綜合功能區(qū)，規(guī)劃面積約200平方公里。

圖1 研究區(qū)概況

濱海新區(qū)地處華北平原的東部塌陷盆地邊緣，位于渤海盆地的西部，處在黃驊凹陷區(qū)北部區(qū)域。地形地貌類型是由沖積平原沉積而成，具有沖積平原到潮間帶之間完整的地貌分布帶規(guī)律。近年來隨著經(jīng)濟(jì)建設(shè)，填海造陸的人工地貌單元面積超過200平方公里。

1.2 InSAR數(shù)據(jù)

所采用的數(shù)據(jù)集來自Sentinel-1衛(wèi)星TOPSAR數(shù)據(jù)，其中升軌數(shù)據(jù)72期，時(shí)間范圍從2017年1月3日至2019年4月11日；降軌數(shù)據(jù)66期，時(shí)間范圍從2017年1月2日至2019年4月10日。Sentinel-1衛(wèi)星是近極地太陽同步軌道衛(wèi)星，重訪周期為12天，基于C波段采用右側(cè)視觀測，具有4種成像模式，其分辨率最高達(dá)5 m、幅寬達(dá)400 km[10]。數(shù)據(jù)范圍如圖2所示。

圖2 SAR數(shù)據(jù)范圍

2 InSAR數(shù)據(jù)處理

2.1 SBAS-InSAR解算

SBAS-InSAR是對相干目標(biāo)進(jìn)行相位分析來獲取時(shí)序形變，通過選擇合適的空間基線和時(shí)間基線閾值組成差分干涉對，并且選取相干目標(biāo)點(diǎn)利用線性相位變化模型進(jìn)行建模和解算，并通過時(shí)空濾波去除大氣延遲，在減少失相關(guān)影響及高程誤差、大氣誤差的同時(shí)獲取地表的形變時(shí)間序列。SBAS-InSAR技術(shù)可以在利用較少數(shù)據(jù)對的情況下獲取較好的結(jié)果，適用于監(jiān)測大區(qū)域、長時(shí)間序列的微小緩慢形變。本文在組成干涉對時(shí)，時(shí)間基線閾值選擇50天，空間基線閾值選擇100 m，升軌數(shù)據(jù)共形成237個(gè)干涉對，降軌數(shù)據(jù)共形成203個(gè)干涉對，干涉對組合見圖3。

SABS-InSAR算法重在應(yīng)用短基線的差分干涉圖集，集合內(nèi)的影像對基線距小，而集合間的影像對基線距大，集合構(gòu)成相位回歸分析的序列，用均值相干系數(shù)作為相干目標(biāo)識別的指標(biāo)，再根據(jù)最小范數(shù)準(zhǔn)則，利用奇異值來分解算法逐個(gè)分離相位構(gòu)成，實(shí)現(xiàn)每個(gè)相干目標(biāo)形變序列的獲取，減少空間和時(shí)間去相關(guān)影響，保證精度與可靠性[11]。

差分干涉相位組成如下：

Δφ=Δφdefo+Δφatm+Δφorbit+ΔφΔDEM+Δφnoise

(1)

式中，Δφdefo是地表形變，Δφatm是大氣延遲，Δφorbit是軌道誤差，ΔφΔDEM是DEM誤差，Δφnoise是噪聲誤差。假設(shè)不考慮大氣延遲、軌道誤差、DEM誤差和噪聲，且干涉圖中的解纏不存在誤差，可建立如下最小二乘模型

Aφdefo=Δφ

(2)

式中，A為M×N的設(shè)計(jì)矩陣，其中M為干涉圖數(shù)量，N為SAR影像數(shù)量。Δφ為已知觀測量，φdefo為未知參數(shù)向量。由于M≥N，設(shè)計(jì)矩陣A的秩為N，所以有唯一解，根據(jù)經(jīng)典最小二乘原理VTPV=min,可以得到

Δφdefo=(ATA)-1ATΔφ

(3)

圖4為利用上述SBAS-InSAR算法獲得的升降軌視線向沉降速率?？梢钥闯錾壟c降軌數(shù)據(jù)得到的沉降空間分布基本一致，天津市區(qū)范圍的沉降控制較好，沉降主要集中在濱海新區(qū)西南部與三個(gè)研究區(qū)域范圍內(nèi)，最嚴(yán)重的沉降情況出現(xiàn)在臨港工業(yè)區(qū)。

圖3 升降軌數(shù)據(jù)InSAR干涉對組合

圖4 天津市濱海新區(qū)LOS形變速率

2.2 InSAR升降軌融合提取水平和垂直向形變

對升降軌數(shù)據(jù)進(jìn)行融合分解可以獲取地表水平東西向和垂直向形變分布信息。如圖5所示，展示了SAR側(cè)視成像幾何及三維坐標(biāo)下的角度參數(shù)。

圖5 SAR側(cè)視成像幾何及角度參數(shù)

InSAR干涉圖中某個(gè)像素的視線向LOS的形變實(shí)際上是該像素點(diǎn)三維位移分量在視線向的投影之和，如下：

rLOS=uxcosφsinθ-uysinφsinθ-uzcosθ

(4)

式中，rLOS為視線向形變量；θ為雷達(dá)信號入射角，即視線向與垂直方向的夾角；φ為衛(wèi)星航向角。地面的位移矢量u=[uxuyuz]，ux為沿東西方向的形變分量，uy為沿南北方向的形變分量，uz為沿垂直方向的形變分量，可將視線向的沉降量r，由式(4)分解得

(5)

所以視線向rLOS中令uy為零，且考慮升降軌兩個(gè)視線向形S變則列出兩個(gè)方程，如下:

(6)

其中ra為升軌視線向的沉降量，rd為降軌視線向的沉降量將式(6)列為矩陣形式：

r=Pu′

(7)

u′=[PTP]-1×PT×r

(8)

得到ux與uz，即東西方向與垂直方向的形變。 圖6展示了研究區(qū)東西方向與垂直方向的形變速率。從垂直形變速率圖中能明顯的看出沉降的空間分布情況。沉降主要集中在填海造陸區(qū)域，同時(shí)此區(qū)域的東西方向形變也更為嚴(yán)重。

圖6 天津市濱海新區(qū)沉降與東西方向形變速率

3 天津沿海地面沉降特征分析

3.1 InSAR與GPS的形變結(jié)果比較

首先對InSAR結(jié)果進(jìn)行精度分析，選取試驗(yàn)區(qū)中的三個(gè)GPS監(jiān)測站點(diǎn)數(shù)據(jù)進(jìn)行對比驗(yàn)證。GPS站點(diǎn)位置如圖1中的五角星所示。

從圖7中可以看出，整體來說，二者垂直方向的結(jié)果更為吻合(RMSE(均方根誤差)小于7 mm)，略優(yōu)于東西方向的結(jié)果(RMSE小于10 mm)，印證了InSAR對垂直方向的位移敏感度最高，對東西方向位移次之。分別來說，TJBH站點(diǎn)數(shù)據(jù)吻合較好，東西方向與垂直方向RMSE分別為6.72 mm與6.60 mm；TJA1站點(diǎn)的垂直方向結(jié)果最佳，RMSE最小，為6.24 mm；KC01站點(diǎn)垂直方向RMSE為6.80 mm；而KC01與TJA1站點(diǎn)東西方向的結(jié)果并不理想，RMSE分別為9.55 mm和7.86 mm。

圖7 三個(gè)站點(diǎn)InSAR結(jié)果與GPS結(jié)果對比

3.2 濱海新區(qū)地面沉降空間分布特征

圖8展示了研究區(qū)的沉降分布情況，可以從中看出濱海新區(qū)的沉降呈現(xiàn)出不均勻的特點(diǎn)，內(nèi)陸地區(qū)的沉降較小，越靠近沿海地區(qū)，沉降量越大，而在填海造陸區(qū)，沉降則十分嚴(yán)重。原因是填海造陸區(qū)地表淺層的地質(zhì)結(jié)構(gòu)處于不穩(wěn)定的狀態(tài)，地表人類進(jìn)行各種建筑的建設(shè)活動與生產(chǎn)活動，土層固結(jié)壓縮，導(dǎo)致地表下沉較大。

圖8 濱海新區(qū)沉降速率分布圖

其中，漢沽區(qū)的情況相對較好，最大沉降速率為209.9 mm/year，南港工業(yè)區(qū)的最大沉降速率為240.7 mm/year，最嚴(yán)重的下沉情況出現(xiàn)在臨港工業(yè)區(qū)，最大下沉速率達(dá)275.4 mm/year。

根據(jù)圖8結(jié)果，在三個(gè)沉降嚴(yán)重的區(qū)域內(nèi)選取三個(gè)具有代表性的沉降點(diǎn)(圖8中的P1-P3點(diǎn))，進(jìn)行沉降時(shí)間序列分析，并進(jìn)行線性擬合(圖9)?？梢詮膱D9看到各點(diǎn)數(shù)據(jù)有明顯的下沉趨勢，且擬合情況良好。P1點(diǎn)最大下沉量為137 mm，經(jīng)過線性擬合得到下沉速率為64.54 mm/year；P2點(diǎn)最大下沉量為295 mm，經(jīng)過線性擬合得到下沉速率為155.03 mm/year；P3點(diǎn)最大下沉量為342 mm，經(jīng)過線性擬合得到下沉速率為167.62 mm/year。

圖9 三個(gè)沉降點(diǎn)垂直方向形變時(shí)間序列圖

地表沉降的原因一般分為自然因素和人為因素，自然因素有地質(zhì)構(gòu)造活動、土體變形固結(jié)、地層的壓密和地殼的升降運(yùn)動；人為因素則是地下水抽取、地下自然資源開采和地下工程施工影響[13]。資料顯示，天津市存在一條海河斷裂帶，作為一種重要的基底活動構(gòu)造，對沉降情況起了一定的作用；唐山大地震曾波及天津市及沿海地區(qū)，并且在漢沽區(qū)形成了一個(gè)下沉中心，這種地震活動，同樣也對沿海區(qū)的地表沉降情況產(chǎn)生一定影響；沿海區(qū)的土質(zhì)屬于軟弱土層，在自重壓力與地表活動影響下，地下水的壓力逐漸減小，地層固結(jié)，也導(dǎo)致了一定量的地表下沉。

3.3 地面沉降與海岸線變化分析

本文選取了1995年到2018年，每隔5年共6幅的云量低于20%的Landsat TM和Landsat ETM+影像，對天津市填海區(qū)的海岸線變化進(jìn)行提取。圖10展示了天津市填海區(qū)1995～2018年的海岸線變化和此區(qū)域垂直方向的沉降速率對比。

可以看出，圍海造陸的區(qū)域從2005年之后大規(guī)模擴(kuò)張，臨港工業(yè)區(qū)在2010形成了如今的規(guī)模，南港工業(yè)區(qū)在2015年形成如今規(guī)模，開墾土地的主要應(yīng)用為建設(shè)用地和濕地。此區(qū)域上集中了濱海新區(qū)大部分的工業(yè)園區(qū)，岸線長度占了整個(gè)海岸線長度的76.8%，并且仍在隨著圍海造陸活動與港口發(fā)展迅速向海域擴(kuò)張[14]。

由圖10可知，沉降嚴(yán)重的區(qū)域分布在近年的填海區(qū)，填海造陸活動本身對地表淺層的地質(zhì)結(jié)構(gòu)會有根本性的影響，地質(zhì)較為松軟，不斷積壓，加上海水的不斷侵蝕，加劇了地表沉降。同時(shí)大部分土地利用類型是建筑用地，地表的建設(shè)活動與載荷也會產(chǎn)生一定量的地表沉降。所以研究區(qū)的沉降分布情況與圍海造陸活動相關(guān)。根據(jù)區(qū)域監(jiān)測結(jié)果分析，臨港工業(yè)區(qū)不同填土厚度區(qū)域沉降差異較大，最初5～6年對地表沉降影響較大，之后逐年遞減[15]。

4 結(jié)論

(1) 海岸線的演變情況顯示出了填海造陸的歷史進(jìn)程，反映了沿海區(qū)地面沉降情況嚴(yán)重的原因。濱海新區(qū)的沉降呈現(xiàn)出不均勻的特點(diǎn)，沿海區(qū)比市區(qū)沉降嚴(yán)重，最嚴(yán)重的沉降情況發(fā)生在臨港工業(yè)區(qū)，最大下沉速率達(dá)275.4 mm/year。

(2) 通過填海造陸建成的臨港工業(yè)區(qū)與南港工業(yè)區(qū)，地表淺層的地質(zhì)結(jié)構(gòu)處于不穩(wěn)定的狀態(tài)，地表各種建筑的建設(shè)活動與載荷，加劇了填海區(qū)的地表沉降情況，導(dǎo)致填海造陸區(qū)比天然海岸和市區(qū)沉降嚴(yán)重，同時(shí)還展現(xiàn)出不均勻的特點(diǎn)。

(3) InSAR結(jié)果與GPS監(jiān)測數(shù)據(jù)在垂直方向上更為吻合，略優(yōu)于東西方向，印證了InSAR對垂直方向的位移敏感度最高，對東西方向位移次之。

圖10 天津市填海區(qū)1995-2018年海岸線變化與沉降速率分布圖