基于PS-InSAR与SBAS-InSAR的地面沉降监测分析

陈媛媛 朱晨玮 郑加柱 王诗涵

摘要：

为了全面掌握长江岸线区域的地面沉降情况，为维护长江防洪安全以及河势稳定提供参考依据，以长江南京段沿岸5 km缓冲区域为研究区域，采用PS-InSAR和SBAS-InSAR两种时序分析技术对2017年4月至2021年4月期间30景Sentinel-1A影像进行处理与分析比较，并分析了造成沉降的原因。结果表明：两种方法所得的特征点的形变曲线变化趋势一致，通过对沉降点的年平均地表形变速率进行线性拟合以及对其差值进行正态分布分析，验证了两种方法所得结果具有一致性和可靠性；两种方法所得沉降区域分布情况大致吻合，该区域有5个沉降漏斗，分别位于浦口区的长江沿岸区域、雨花台区东北部与建邺区西南部的交界处、建邺区东北部与鼓楼区西南部交界区域、八卦洲区域以及六合区张营水库区域；地质构造、城市改造建设及承压水头下降均对该区域地面沉降产生一定的影响。

关 键 词：

地面沉降； PS-InSAR； SBAS-InSAR； 時序分析； 长江岸线； 长江南京段

中图法分类号： P237；P642.26

文献标志码： A

DOI：10.16232/j.cnki.1001-4179.2023.08.022

0 引 言

江苏地处长江下游，是南水北调东线工程的水源地，也是长江黄金水道水上运输大动脉的重要组成部分。江苏省防洪保安任务重、河势控制难度大、水资源保护要求高、岸线利用强度大。确保长江江苏段防洪安全、维护河势稳定、促进长江生态向好发展等问题是长江治理的关键问题之一［1］，其中，确保沿江岸线及其附近区域地面的稳定对于长江治理与保护至关重要，有必要开展长江沿岸区域的长周期持续地面监测工作。

目前，地面沉降监测的手段和技术主要有精密水准测量和全球定位系统（GPS）技术。这些传统的监测技术存在工作量大、效率低、易受测区范围影响等问题，且只能获取离散监测点的下沉量，无法监测地面下沉趋势和具体范围。合成孔径雷达差分干涉测量技术（InSAR）是一种全新的对地观测技术，它使用两幅或多幅合成孔径雷达影像图，根据卫星或飞机接收到的回波相位差来生成地表形变图，以探取地表的微小形变信息［2］。InSAR技术具有全天时全天候探测、探测范围广、探测精度高、探测成本低等优点，但常规InSAR技术常面临时间失相关、空间失相关和大气延迟的问题。为了解决这些问题，就产生了永久性散射体（Persistent Scatterer，PS）和小基线集（Small Baseline Subsets，SBAS）时序分析等技术［3-5］，并且在城市地面沉降、矿山沉降、地震形变、火山活动等领域得到广泛的应用。目前，已有不少专家基于不同的InSAR技术对南京市的地面沉降情况进行了研究与分析，如：黄其欢等采用13景ERS2 SAR影像，基于短基线D-InSAR方法对南京河西漫滩区进行监测，结果表明该区域最大累计沉降量达12 cm［6］。高二涛等采用2015～2017年覆盖南京地区的23幅Sentinel-1A影像，基于PS-InSAR和SBAS-InSAR对南京主城区的地面沉降进行监测，结果发现该区域的最大沉降速率达-30 mm/a［7］。杨振等利用SBAS-InSAR技术对2015～2018年的19景南京河西的Sentinel-1A数据进行处理，发现该区域存在4个沉降漏斗［8］。但是，鲜有针对长江岸线南京段区域的时间序列地面监测研究。

本文以长江岸线南京段附近区域为典型研究区域，获取了30景欧洲航天局（以下简称“欧空局”）（Sentinel-1A）卫星雷达影像，基于PS-InSAR和SBAS-InSAR两种时序分析技术对南京段长江岸线区域进行地面沉降监测，并通过对特征点在两种方法结果中的形变曲线进行趋势分析，以及对8 659个经纬度完全相同的特征点的年平均沉降速率及其差值进行相关性分析和正态分布拟合，验证两种方法所得结果的一致性，从而保证研究区域沉降监测结果的可靠性与真实性。此外，还从地质构造、城市改造建设等方面对造成长江南京段地面沉降的原因进行了分析。

1 研究区域及数据简介

1.1 研究区域概况

南京市位于长江下游中部，江苏省西南部。市域地理坐标为北纬31°14′～32°37′、东经118°22′～119°14′。南京市跨江而居，北边是广阔的江淮平原，东边靠近长江三角洲。《南京市长江岸线保护办法》对长江岸线做出了解释：长江岸线指本市行政区域内长江（含洲岛）水陆边界一定范围内的带状区域。本文将长江沿岸5 km范围内的区域视为长江岸线区域（见图1红色边界内区域）。

南京的地质在全国大地构造单元上属扬子古陆的北部边缘，基底主要是轻变质的片岩和变质的火成岩。对于南京长江沿岸地区而言，其地质条件主要特点是地质松散、承载力差、地下水位埋藏浅等。近年来地铁建设施工及地下空间开发等一系列工程项目的开展，造成南京市地表以及长江岸线区域的沉降日益明显。

1.2 实验数据

Sentinel-1A卫星由欧空局于2014年4月在法属圭亚那库鲁发射，轨道高度693 km，搭载的是C波段合成孔径雷达。它具有12 d的重访周期，并且有着多种极化方式和工作模式。本文选用2017年4月至2021年4月南京市30景Sentinel-1A干涉宽幅（IW）模式VV极化的影像进行实验，影像空间分辨率为5 m×20 m（方位向×距离向）。由于IW模式的影像覆盖范围比较大，本文在实验前对数据进行了裁剪。此外，实验还获取了Sentinel-1A卫星对应成像时间的精密轨道数据，地形参考数据使用的是南京地区12.5 m高分辨率DEM产品数据。

2 数据处理

2.1 PS-InSAR技术处理流程

永久散射体的干涉测量技术（PS-InSAR）于2000年由Ferretti等人首次提出［3］，指通过对研究区域内长时间保持稳定散射特征的物体进行研究，即研究一组影像数据中具有高相干性与稳定反射特性的地物点作为PS-InSAR点以提取相位变化，从而避免时间失相关与空间失相关的问题，并且更好地削弱了大气效应对形变信息的影响，以得到更高精度的地面形变信息与DEM［9-10］。PS-InSAR技术数据处理流程如图2所示。

本实验选择2019年9月2日的SAR影像作为主影像，其余的SAR影像与之组成29对干涉影像对。空间基线主要分布于-150～110 m之间，最长空间基线为-140.39 m，最短空间基线为-4.48 m。对每一对像对进行干涉流数据处理，分别为影像配准、干涉图生成、去平地效应以及振幅离差指数计算。

采用振幅离散指数方法［11］，显示单个像元值的相干性，辅助PS-InSAR点的筛选，完成高相干性目标点的初步选择。生成的干涉图需要剔除偏移相位，偏移相位根据振幅离差指数筛选多个参考点计算得到，参考点的个数可根据研究区域的面积大小设置，本次实验设置为每25 km2的区域选择一个参考点，重叠比例参数设置为30%。使用ENVI提供的线性反演模型，基于选择出的参考点，对地表形变速率与地形残差进行估算。

根据大气影响在空间上具有高相关性、在时间上具有低相干性的特性，將反演得到的地表形变速率与地形残差在空间上采用高通滤波，在时间上采用低通滤波，完成对大气延迟相位的估算，消除干涉相位中的大气延迟相位，克服D-InSAR的大气相干与时空相干等问题。再次使用ENVI提供的线性反演模型重新对地表形变速率与地形残差进行二次估算，获得PS-InSAR点最终的形变速率与残余高程［12］。

为便于在地理空间中对PS-InSAR点进行几何信息与形变信息的特征分析，需要将上述过程中获取的雷达坐标系下的PS点进行地理编码。本次实验使用12.5 m高分辨率DEM产品对PS点进行地理编码，相干系数阈值设置为0.75，经过数据处理得到具有高相干性的PS点，在ArcGIS中进行时序分析与空间分析。

2.2 SBAS-InSAR技术处理流程

小基线集算法（SBAS-InSAR）于2002年由Berardino等人提出［4］，该技术仅需要少量的SAR影像数据组合成若干个短基线的干涉相对，从而得到多个差分干涉集合以增加干涉图的数量，再利用最小二乘法与奇异值分解法在解决时空失相关问题的同时增加采样率，以此获取连续且高精度的地表形变信息［13-14］。SBAS-InSAR的主要技术流程如图3所示。

该部分实验选择2019年11月1日的SAR影像作为超级主影像，将其他的SAR影像与其进行配准。配准完成之后，需要对所有的干涉对进行时间基线与空间基线的筛选，由于本次研究获取的SAR影像数据集时间跨度较长（48个月），故设置时间基线阈值为550 d。此外，因为Sentinel-1A卫星轨道较为稳定，产生的偏差较小，所以设置空间基线为临界基线的45%。本研究通过引入12.5 m的DEM数据来消除平地效应，采用最小费用流（MCF）方法进行相位解缠［15］。对解缠后的结果进行轨道精炼与重去平处理，选取高质量的地面控制点，通过多项式法来消除轨道参数误差与去除平地效应。随后，选择线性模型来初步估算形变速率与地形误差，再通过时间域高通滤波与空间域的低通滤波去除大气延迟误差。

将SBAS-InSAR反演结果进行地理编码，将斜距坐标系转换为地理坐标系。此外，需要将结果的栅格数据转换为矢量数据，以便在ArcGIS中进行空间分析。

3 结果分析与比较

3.1 PS-InSAR结果分析

通过PS-InSAR数据处理流程得到研究区域的年平均形变速率图，以长江为对象构建一个5 km的缓冲区对研究区域的地面沉降监测结果进行裁剪，最终得到长江南京段沿线区域的年平均形变速率图（见图4），其中红色矩形表示主要沉降区域，紫色十字符表示特征点位置。

从图4可以看出，南京市长江流域的年平均形变速率分布范围为-25～15 mm/a，采用红色区域表示地面沉降，蓝色区域表示地面抬升。研究区域中的主要沉降区域分布于浦口区的长江沿岸区域（A）、雨花台区东北部与建邺区西南部的交界处（B）、建邺区东北部与鼓楼区西南部交界区域（C）、八卦洲区域（D）以及六合区张营水库区域（E），沉降严重的区域大多分布于长江沿岸，最大沉降速率达到-22 mm/a，最大沉降量达到-33 mm，最大累计沉降量为-101 mm。

根据主要沉降区域的分布位置，发现沉降区域A位于浦口区的江浦街道与顶山街道的连接地区，该区域范围内有地铁十号线通过且设有3个站点，其最大沉降速率为-19.3 mm/a。包含地铁二号线与地铁十号线的建邺区地面沉降程度最严重，其西南地区与东北地区分别有两个沉降漏斗B与C，且在空间上具有相连成片的趋势，其中B区域最大沉降速率为-18.2 mm/a，C区域最大沉降速率为-19.8 mm/a。由于长江冲积作用形成的八卦洲地区的南部地区以及位于八卦洲北侧六合区的张营水库附近分别存在大面积沉降区域D与E，其中D区域最大沉降速率为-13.3 mm/a，E区域最大沉降速率为-12.8 mm/a。

3.2 SBAS-InSAR结果分析

图5为采用SBAS-InSAR技术处理所得的长江南京段缓冲区的形变速率图，其中红色矩形表示与图4对应的沉降区域，紫色十字符表示对应的特征点位置。由该图可以看出，南京市长江岸线年平均形变速率范围为-19～10mm/a，在整个长江岸线，有5个主要沉降区，其位置分别在浦口区长江岸线（A）、建邺区东北部与鼓楼区西南部交界区域（C）、八卦洲区域（D）以及六合区张营水库区域（E）。最大沉降速率达到-19 mm/a，最大沉降量达到-44 mm，最大累计沉降量为-88 mm，其中，A区（浦口区长江岸线）沉降范围最广，C区是建邺区东北部与鼓楼区西南部交界处的秦淮区莫愁湖区域，其沉降最严重。

根据SBAS-InSAR结果中主要沉降区域的分布位置，可以看出A区域是由江浦街道、顶山街道、泰山街道以及沿江街道共4个街道组成的沉降区，该范围内共有3条地铁经过且设有7个站点，A区域最大沉降速率为-19.2mm/a；建邺区存在两个沉降漏斗B和C，B区域最大沉降速率为-11.9 mm/a，C区域最大沉降速率为-14.4 mm/a；沉降区域D和E的最大沉降速率分别为-12.4 mm/a和-11.8 mm/a。

3.3 两种技术手段结果对比分析

由于本次实验缺乏南京市长江岸线同期水准数据，故通过对比分析两种处理方法的结果来验证结果的可靠性。选取沉降速率较大的5个区域（A，B，C，D，E区域，见图4～5），并在每个区域内选取可以代表该区域沉降特征的特征点（图4～5中的紫色点位），通过绘制5个特征点在整个时间跨度中的累计沉降量变化图（见图6）来对比分析两种技术方法获取的形变结果。

从图6可以直观地看出：特征点的形变曲线整体呈非线性变化，都具有形变缓慢、沉降速率变化不大的特点。PS-InSAR处理方法所得特征点的沉降曲线较为粗糙，形变曲线抖动较大。而SBAS-InSAR处理所得形变曲线更光滑，趋势也更为稳定。从形变曲线的趋势和形变结果来看，特征点处的PS-InSAR的最大沉降量为-60 mm，稍小于SBAS-InSAR的-65 mm。两种方法获取的结果在沉降量级上虽然存在微小差异，但总体而言具有较强的一致性，沉降区域大致重合，因此实验结果具有较高的可靠性。

此外，本文还对PS-InSAR与SBAS-InSAR两种方法提取的8 659个经纬度完全相同的特征点的年平均沉降速率及其差值进行正态分布拟合和相关性分析，来进一步验证两种方法获取结果的可靠性和一致性。如图7所示，PS-InSAR和SBAS-InSAR两种方法获取的年平均沉降速率差值的平均值和标准差分别为0.42 mm/a和0.9 mm/a，满足CH/T 6006-2018《时间序列InSRA地表形变监测数据处理规范》中的精度要求［16］。图8中以PS-InSAR监测结果为纵坐标，以SBAS-InSAR监测结果为横坐标，8 659个沉降点均匀分布在拟合直线两侧，决定系数R=0.789 2，均方根误差为0.456 1 mm/a，由此可得两种处理结果具有一致性。

3.4 成因分析

（1） 地质构造。长江岸线南京段两侧漫滩宽广，据其地形地貌单元与岩土组合，可将其归为河漫滩软土地区，主要由长江冲击、淤积作用形成的第四系沉积软弱土层构成，故其具有地质松散、承载力差等特点。实验结果中研究区域有5个主要沉降区，均分布在长江岸线附近，岸线的沉降与地质构造具有密切联系。其中A区域（浦口区岸线）沉降范围最大，最大沉降速率达-16 mm，累计沉降量达-40 mm。B区与C区均为三角洲区域，长期受到河水冲击作用，使得岸底累积大量沉积物，是造成区域沉降的重要因素之一。

（2） 城市改造建设。近年来，南京市地铁建设及地下工程施工十分频繁，使得部分区域的沉降日益明显，其中C区周边地下有多条地铁线路运营，该区域内的最大沉降速率达-19 mm/a，累计沉降量达-60 mm。A区域为南京江北新区，作为南京着力打造的长三角区域金融中心核心功能区及国家级发展新区，近几年建筑改造施工频繁，地下开发严重，而开发区域都是以河滩为主的软基区域，这是造成该区域大面积沉降的重要原因之一。

（3） 其他因素。相关研究表明，长江岸线南京段地面沉降与该区域附近承压水头下降存在较高的相关性［17］。此外，由于其他人为与自然因素，在某种程度上对区域沉降也有着重要的影响。E区为六合区张营水库，该区域内的地面沉降与坝体的自身结构有关。

4 结 论

本文利用从欧空局获取的30景Sentinel-1A卫星雷达影像，采用PS-InSAR和SBAS-InSAR时序分析技术对南京长江岸线区域进行地面沉降监测。实验结果发现：该区域有5个沉降漏斗，分别位于浦口区的长江沿岸区域、雨花台区东北部与建邺区西南部的交界处、建邺区东北部与鼓楼区西南部交界区域、八卦洲区域以及六合区张营水库区域。对两种方法所获结果的特征点沉降速率进行比较，分析了年平均地表形变速率的线性拟合曲线及差值正态分布曲线，发现两种方法获取的地面形变结果具有较高的相关性，沉降趋势与沉降区域的空间分布具有一致性，侧面验证了这两种方法的可靠性。本文研究可为政府对长江岸线环境保护制定相关决策提供强有力的数据支持，为确保长江岸线区域的防洪安全、维护河势稳定提供一定的理论依据。

参考文献：

［1］ 朱宇驰.近50年长江南京段河道变迁研究［D］.南京：南京大学，2016.

［2］ 朱建军，李志伟，胡俊.InSAR变形监测方法与研究进展［J］.测绘学报，2017，46（10）：1717-1733.

［3］ FERRETTI A，PRATI C，ROCCA，et a1.Nonlinear subsidence rate estimation using Permanent Scatters in differential SAR Interferometry［J］.IEEE Transactions on Geoscience and Remote Sensing，2000，38（5）：2202-2212.

［4］ BERARDINO P，FORNARO G，LANARI R，et a1.A new algorithm for surface deformation monitoring based on small baseline Differential SAR interferograms［J］.IEEE Transactions on Geoscience and Remote Sensing，2002，40（11）：2375-2383.

［5］ 侯安業，张景发，刘斌，等.PS-InSAR与SBAS-InSAR监测地表沉降的比较研究［J］.大地测量与地球动力学，2012，32（4）：125-128，134.

［6］ 黄其欢，徐佳.短基线DInSAR法长江漫滩区大面积沉降监测研究［J］.山东科技大学学报（自然科学版），2011，30（1）：7-10，15.

［7］ 高二濤，范冬林，付波霖，等.基于PS-InSAR和SBAS技术监测南京市地面沉降［J］.大地测量与地球动力学，2019，39（2）：158-163.

［8］ 杨振，徐佳.利用Sentinel-1A数据监测南京河西地区地面沉降［J］.测绘通报，2020（1）：61-65，75.

［9］ REHAN K.Monitoring subsidence in urban area by PSInSAR：A case study of Abbottabad City，Northern Pakistan［J］.Remote Sensing，2021，13（9）：1651-1651.

［10］ 朱邦彦，唐超，任志忠，等.基于PS-InSAR技术的珠海市地表形变监测与驱动力分析［J］.测绘通报，2022（6）：108-113.

［11］ 刘昱彤.PS-InSAR技术的地表沉降监测研究［D］.唐山：华北理工大学，2020.

［12］ 罗三明，杜凯夫，畅柳，等.基于PS-InSAR方法反演北京地区地表沉降速率［J］.大地测量与地球动力学，2014，34（1）：43-46.

［13］ 娄连惠，刘强，谭玉敏.基于高分辨率SAR的三峡库区地表形变分析：以巴东城区为例［J］.人民长江，2021，52（5）：114-119.

［14］ 高文锦.基于SBAS-InSAR技术的深圳西部地区地表形变研究［D］.西安：西安科技大学，2021.

［15］ 王志勇，孟庆颖.InSAR相位解缠算法在矿区沉降监测中的对比与分析［J］.遥感信息，2014，29（6）：47-51，57.

［16］ 国家测绘地理信息局.时间序列InSAR地表形变监测数据处理规范：CH/T 6006-2018［S］.北京：测绘出版社，2018.

［17］ 朱邦彦，姚冯宇，孙静雯，等.利用InSAR与地质数据综合分析南京河西地面沉降的演化特征和成因［J］.武汉大学学报（信息科学版），2020，45（3）：442-450.

（编辑：高小雲）

Abstract：

In order to fully understand the land subsidence in the riparian area along the Changjiang River，and provide references for guaranteeing the flood control safety and the river regime stability of Changjiang River，a 5km buffer area along the Nanjing section of Changjiang River was taken as the study area，two time series analysis technologies，PS InSAR and SBAS InSAR，were used to process，analyze and compare with 30 Sentinel-1A images from April 2017 to April 2021，and the causes of subsidence were analyzed.The results showed that the deformation curves of feature points obtained by the two methods were basically consistent，and the consistency and reliability of the two results were also verified by linear fitting for annual mean ground deformation rate and the normal distribution analysis for the difference of the ground deformation rate.The distribution of the subsidence areas obtained by the two methods coincided roughly.There were 5 subsidence funnels in this area，which were located along the Changjiang River in Pukou District，the junction of the northeast of Yuhuatai District and the southwest of Jianye District，the junction of the northeast of Jianye District and the southwest of Gulou District，Baguazhou and Zhangying reservoir area of Liuhe District.It was found that geologic structures，urban renewal construction and drop of confined water all had an effect on land subsidence in the region.

Key words：

land subsidence；PS-InSAR；SBAS-InSAR；time series analysis；Changjiang River bank；Nanjing section of the Changjiang River