基于时序InSAR的珠江口大面积地面沉降监测

高 辉，罗孝文*，吴自银，阳凡林

(1.山东科技大学 测绘科学与工程学院，山东 青岛 266590；2.自然资源部 第二海洋研究所，浙江 杭州 310012；3.自然资源部 海底科学重点实验室, 浙江 杭州 310012)

0 引言

珠江三角洲位于中国广东省中南部，毗邻港澳，区域内城市密集，经济发达，有“南海明珠”之称。河口等沿海地区同时受海洋、陆地及河流的相互作用，生态环境非常脆弱[1]。近30 a来珠江口发展迅速，大量基础设施建设、地下资源开采、软土固结压实等人类活动和自然因素的影响，使该区域发生大面积的地面沉降，带来地面塌陷、基础设施破坏等地质灾害，全球海平面上升[2]加剧了咸潮上溯的影响及风暴潮灾害的风险[3-4]，愈发威胁该区域的生态环境及人类的生产生活。因此，及时准确地揭示珠江河口地区地面沉降的时空分布和形变趋势具有重要的意义。

合成孔径雷达干涉测量(Interferometric Synthetic Aperture Radar，InSAR)技术具有全天时、全天候、覆盖范围广等特点，能够获取高分辨率、高精度的地形及地表形变信息[5-6]。近年来，以永久散射体技术(Persistent Scatterers，PS)[7]和小基线集技术(Small BAseline Subset，SBAS)[8]为代表的时序InSAR技术已广泛应用于河口海岸带地区的地面沉降监测研究。刘勇[9]、张金芝 等[10-11]及刘一霖 等[12]分别用不同的方法研究了黄河口地区的地面沉降；JIANG et al[13-14]利用SAR干涉点目标分析法(IPTA)监测香港国际机场及澳门填海区的地面沉降，并揭示了填海与沉降的关系；WANG et al[15]利用多轨SAR数据采用π-rate方法监测了珠江三角洲深圳、东莞地区的地面沉降，得到沿岸500 m内年沉降速率为5 mm/a；ERBAN et al[16]探测了湄公河三角洲的地面沉降，并认为地下水抽取是引起大范围沉降的原因。

这些时序InSAR方法对监测点目标的选取十分严格，要求在所涉及的干涉纹图中都保持较高的相干性[17]。而河口海岸带地区植被覆盖密，人工目标稀疏，大气扰动相对复杂，难以准确反映地面沉降的时空变化[18-19]，限制了这些方法在河口海岸带地面沉降检测中的应用。

FERRETTI et al[20]提出SqueeSAR技术，将同分布散射体(Distributed Scatterer，DS)融入到传统时序InSAR技术中，大大提高了有效监测点的密度。但该方法求解非线性最优相位需要迭代，影响计算效率。针对上述问题，HAN et al[21]提出一种基于特征分解估算DS候选点最优相位的方法，该方法大大提高了计算效率，并通过杭州湾沿岸上虞区的应用验证了在沿海地区的有效性。鉴于目前关于珠江口地区的地面沉降研究都在小区域范围内[15,22-23]，随着Sentinel-1A数据的免费公开，本文采用HAN et al[21]提出的方法，利用2015—2018年间67景Sentinel-1A SAR影像对珠江口地区进行大范围的地面沉降监测，获得了该区域的地面沉降信息及时空分布，并对几处沉降严重地区进行了分析。

1 技术方法

采用基于Kolmogorov-Smirnov(KS)检验和特征分解的DS-InSAR时序分析方法，其基本数据处理如下[21]：首先将Single Look Complex(SLC)影像配准至主影像，采用KS检验算法选择稳定的DS候选点，利用DS候选点目标的相位信息进行空间自适应滤波并构建协方差矩阵；然后通过特征分解方法估算最优相位并代替DS目标原始相位；最后根据拟合优度检测确定最终的DS点目标，联合PS点进行时序InSAR处理，经过地形误差改正、大气误差改正、相位解缠等步骤得到最终的形变速率和形变时间序列。下面重点介绍DS点选取过程中基于特征分解方法估算最优相位的关键步骤。

基于特征分解的最优相位估计方法最初由FORNARO et al[24]提出，通过从协方差矩阵或相干矩阵中提取与最大特征值相关的特征向量的相位估计最优相位。对通过KS检验算法选取并进行空间自适应滤波后构建的相干矩阵进行特征分解，得：

(1)

其中：T为候选点的最优相位，N为干涉相位，λi为特征值且按降序排序，Ui为特征值对应的特征向量。假设U1为最大特征分量，可表示为最大值约束问题，将最大特征值对应的特征向量的相位作为DS候选点的最优相位[25]：

(2)

约束条件：

(3)

2 研究区概况与实验数据

2.1 研究区概况

珠江三角洲位于珠江入海口处，由西江、北江和东江入海时沉积而成，地理坐标为112°30′～114°15′E， 21°17′～23°30′N。珠江三角洲呈倒置三角形，地貌类型以三角洲平原为主体，三角洲前缘及西北部有低山、丘陵分布，中部则缀有零星低丘及台地，平原上河道汊港纵横交错，水塘泽地星罗棋布。地层以第四纪地层分布最广，大面积展布于三角洲平原，岩性以松软土层沉积为主要特征。南部为海岸带，是区内海陆相互作用最强烈的地带，具有海岸线长、海岸曲折、岛屿众多等特点，但软土厚度较大，容易发生地面沉降地质灾害。珠江三角洲地处亚热带季风气候区，温暖潮湿，雨量充沛。近年来人类活动的影响和自然条件的变化使该区域发生大面积地面沉降，严重危害该地区的经济社会可持续发展。本文选择珠江口地面沉降较为严重的地区作为监测区，面积约3 000 km2，有效覆盖了珠江口地区。

2.2 实验数据

使用67景Sentinel-1A SAR影像进行数据处理，Sentinel-1A卫星由欧空局发射，具有短重复周期、大幅宽、多极化等特点。数据获取时间跨度为2015年6月—2018年3月，其采用干涉宽幅(IW)数据成像模式，极化方式为VV极化，所处波段为C波段，具有250 km的扫描宽度和5 m×20 m的地面分辨率。同时采用美国航空航天局(NASA)提供的SRTM 30 m DEM数据，去除地形相位与DEM误差。

3 监测结果

3.1 InSAR地面沉降结果

利用基于KS检验和特征分解的DS-InSAR时序分析方法，处理SAR影像数据，设定时间基线100 d和空间基线200 m的阈值作为小基线集干涉图选取限制，生成了178对干涉组合，其组合如图1所示。采用最小费用流法(minimum cost flow, MCF)进行相位解缠。通过数据处理，获得了珠江口地区2015年6月—2018年3月地面沉降的年均形变速率(图2)和累计形变分布(图3)。其中，负值代表下降，正值代表抬升。从图中可以看出，监测区时序InSAR获取的地面形变结果总体呈现地面沉降趋势，且沉降在空间上分布不均匀，部分区域存在明显沉降现象，抬升区范围较小。珠江口西岸沉降相对东岸分布更广，沉降量更大。珠江三角洲西北部(广州、佛山)和东南部(东莞南部)为主要沉降区，其中监测区西北部佛山禅城区南庄镇及顺德区西部龙江镇、东莞市南部长安镇局部存在严重的沉降现象，最大沉降速率可达25 mm/a。监测稳定区主要集中在东莞市西部虎门镇附近，该区域没有明显的下沉现象。

图1 短基线集干涉组合Fig.1 Combination of small-baseline interferograms

图2 珠江口地区InSAR地面年均形变速率Fig.2 Annual mean deformation rate in the Pearl River Estuary from InSAR time series

图3 珠江口地区InSAR累积形变量Fig.3 Cumulative deformation in the Pearl River Estuary from InSAR time series

通过DS-InSAR方法在实验区获得的DS目标分布，可以看出在非城市区域有着较高密度的高相干点，主要集中在沥青路面、裸土及其他有DS点特性的地物上，有效增强了高质量InSAR监测点的空间密度，进而提高了时序InSAR形变反演的精度。

3.2 沉降严重地区分析

3.2.1 东莞市长安镇

长安镇位于东莞市南端，东邻深圳市，南临珠江口，西连虎门港，是地面沉降信号明显的重点区域。由图4可以看出，长安镇整体都有明显沉降信号，沉降速率在8～13 mm/a范围内。长安镇沉降由西北(靠近虎门)向东南表现为沉降加剧，尤其是靠近珠江口区域，包括南部深圳市宝安区沙井街道，其最大沉降速率达24.1 mm/a，3 a累积最大沉降量可达45.3 mm。近年来东莞市经济飞速发展，随之而进行的大量填海活动、基础设施建设、工业科技园区建设是该地区发生沉降的主要原因。

3.2.2 广州市白云区

广州市白云区位于广州市城区北部，东邻增城区，西界佛山市南海区，南连荔湾区、越秀区、天河区、黄埔区4城区，北接花都区和从化区。由图5可以看出，白云区沙贝立交附近存在严重不均匀分布的地面沉降现象，年均沉降速率为8～15 mm/a；尤其是沙贝立交南部以及东部沙涌村下沙附近，地面沉降中心最大沉降速率可达15.2 mm/a，3 a最大累计沉降量达40 mm；沙贝立交北部为相对稳定区域，年均沉降速率大部分为3～8 mm/a，仅有零星沉降分布。广州白云区位于广花盆地，其地质环境脆弱，相关调查结果表明，基坑工程、地铁工程、隧道工程等地下空间的开发引起的岩溶地面塌陷是该地区地面沉降的主要原因。

3.2.3 佛山市南海区

佛山市南海区地处珠江三角洲腹地，东连广州市芳村区、番禺区，南接顺德区、鹤山市、新会区，西邻三水区、高明区，北邻花都区，环抱佛山禅城区。由图6可知对江沙北部整体较为稳定，部分地区存在地面沉降现象，平均沉降速率为3～8 mm/a，最大沉降中心在海南立交南部，沉降速率为10.2 mm/a；对江沙南部存在3个沉降速率大于15 mm/a的沉降中心，从西至东分别位于平胜大桥东部，新五斗大桥南部，三江郊野公园东南部。其中三江郊野公园东南部最为严重，沉降范围最大，最大沉降速率达23.9 mm/a，3 a累积沉降量最大可达65 mm。近年来，南海区的经济发展迅速，大规模城市、地铁建设和工业园开发以及区内广佛江珠、广珠西线高速公路等荷载量大增，为该地区地面沉降的主要原因。

图4 东莞市长安镇地面年均形变速率Fig.4 Annual mean deformation rate in Chang’an Town, Dongguan

图5 广州市白云区地面年均形变速率Fig.5 Annual mean deformation rate in Baiyun District, Guangzhou

图6 佛山市南海区地面年均形变速率 Fig. 6 Annual mean deformation rate in Nanhai District, Foshan

4 小结

利用基于KS检验和特征分解的DS-InSAR方法，对67景Sentinel-1A影像数据进行时序分析，获得了2015—2018年间珠江口地区地面沉降场的时空演变过程，重点分析了3处沉降较明显区域的沉降时空变化特征。采用的DSI方法在沥青路面、裸土等非城市区域提取了较大密度的高相干点。研究结果表明，珠江口地区总体呈地面沉降趋势，珠江三角洲西北部(广州、佛山)和东南部(东莞南部)为主要沉降区，最大沉降速率可达25 mm/a；珠江口西岸沉降相对东岸分布更广，沉降量更大；东莞市西部虎门镇附近比较稳定，该区域没有明显的下沉现象。

本文监测了珠江口地区的地面沉降，后续将结合海岸线、河岸线的变化，基础设施建设，地下水抽取及填海造陆情况，对该地区地面沉降的因素进行详细的分析，并研究其和海岸线变化的关系，进一步为地面沉降的预防和治理提供参考。