InSAR监测连云港地面沉降与因素分析

摘要： 针对江苏连云港地区在多因素共同作用下形成的复杂地面沉降问题，该文提出一套以合成孔径雷达干涉(InSAR)为主、有机结合地面常规测量、地下水位及分层位移监测的“星—地—内”立体综合观测与驱动因素协同分析框架，并成功应用于获取研究区2012—2020年的地面沉降时空变化特征。结果表明，研究区沉降受不良地质条件、独特沉积特征及人为活动共同作用，监测时段内最大累计沉降量达到571 mm;综合时空变化特征来看，连云区沉降主要由软土固结和工程建设活动引起，呈持续发育状态；灌南和灌云区沉降由地下水开采引起，但不同阶段采水层位和需求存在明显差异：2012—2016年以开采深层地下水满足工业和生活需求为主；2017—2020年则以水产养殖过程中的浅层地下水开采为主。

引言

沿海地区处于海陆交互的特殊地带，地质条件复杂，由于成陆时间晚，本身地质环境较为脆弱[1]。近些年，随着沿海地区开发步伐的加快，人类活动对其影响日益加剧，进而引发了一系列较为严重的地质环境问题[2-3]。地面沉降作为最具代表性的一类，其危害主要体现在两个方面：①增加相对海平面上升风险。地面沉降会造成高程资源的损失，导致相对海平面上升；降低海堤防护能力，加大海岸带城市内涝、海水倒灌等灾害风险的发生[4]。《2020年中国海平面公报》指出，在全球变暖的大背景下，1980—2020年，中国沿海海平面上升速率达到3.4 mm/a, 且呈加速上升趋势，随着社会经济发展进程的加快，沿海地区所面临的海平面上升风险进一步提升。②增加人工目标失稳风险。地面沉降会削弱建、构筑物基础稳定性，缩短使用寿命。尤其对空间跨度较大的高等级公路、铁路等线性工程的影响最为突出，不仅降低运营期道路的安全性，还为待建项目的规划设计、施工等环节增加难度[5-6]。因此，不断优化沿海地区地面沉降监测网布局、提升地面沉降监测水平，全面掌握其空间分布、动态追踪其演化过程，查清其驱动因素和贡献程度对规避形变病害引发的灾害风险、实现地面沉降精准防控具有重要作用。

合成孔径雷达干涉(interferometric synthetic aperture radar, InSAR)测量、水准/全球定位系统(global positioning system, GPS)、分层标组等不同测量方法凭借各自的技术优势在地面沉降多源综合监测方面已有较多的研究和应用成果，文献[7]采用永久散射体InSAR(PS-InSAR)技术对2017—2019年的珠海市地表沉降进行了获取，并利用连续运行参考站(continuous operation of reference stations, CORS)设置的连续观测水准点数据对结果进行了验证；文献[8]采用PS-InSAR技术获取了常州市2018年1月—2019年12月的地面沉降时空演变规律，并结合水准测量、地下水水位变化、单组分层标数据进行了成因分析。而在“星—地—内”立体监测网构建方面，文献[9]对上海市地面沉降监测网络的探索和实践过程进行了总结；文献[10]对江苏沿海地区地面沉降多手段融合立体监测网络的发展过程及应用现状进行了介绍；文献[11]对郑州市地面沉降监测网的建设提出了建议。前人在地面沉降多源综合观测网络建设方面的研究已取得了长足进展，但在构建“星—地—内”地面沉降立体监测与驱动因素协同分析的一体化框架方面探索较少，缺乏从时空变化特征监测到驱动因素分解的整体性、系统性研究，难以实现地面沉降的标本兼治，影响防控工作的实效性。

针对目前存在的问题，本文以时序InSAR技术为主，辅以水准/GPS常规地面测量手段、分层标组及地下水监测方法共同构建了“三位一体、一主两辅”的“星—地—内”地面沉降综合立体观测网络与驱动力协同分析框架；以此为基础，首先获取了连云港地区2012—2020年的长时序地面沉降时空变化特征；然后对各主要沉降区的驱动因素进行了详细分析。该研究工作的开展对进一步优化区域地面沉降立体监测网络的构建、持续保障地面沉降防控工作的精准、高效实施具有重要指导意义。

1.1 地面沉降立体综合观测网

“星—地—内”立体综合观测网在构成上，分为星上、地表和地下3部分。其中，星上部分由满足InSAR监测条件的SAR卫星与高分辨率光学卫星构成；地表部分由水准和GPS等常规观测构成；地下部分则包含分层标组和水位自动化监测井。3部分有机结合表现在地表与地下观测成果为星上InSAR监测提供校正与协同对比分析数据；而InSAR则指导地表和地下监测点的合理布设。在方法侧重上，采用以星上InSAR全域常态化监测为主，地表与地下重点区域稀疏点位监测为辅的组合模式。这是因为InSAR具有经济成本低、监测效率高等突出优势，高密度空间覆盖的监测点位能够实现地面沉降时空变化特征的全面获取，实现长时序历史地面沉降的追溯及未来监测的延续。

1.2 地面沉降驱动力协同分析

立体综合监测网采集的时空观测数据除了能够定量获取地面沉降的时空发育特征外，还能够对诱发地面沉降的驱动因素进行协同分析，并定性评价各影响因素的贡献程度，见图1。

1)卫星部分。

首先，高分辨率光学卫星影像可以发现由地表变化而引起的沉降现象。此外，还能够依据地物类型辅助判断可能产生沉降的原因。其次，InSAR能够通过时空形变特征分析地面沉降的潜在诱发因素。空间上，不同驱动因素导致的沉降区在大小和形态存在一定差异；时间上，形变演化过程具有不同的发育规律。

2)地面部分。

水准或GPS观测的数据在分析地面沉降驱动因素方面的作用有限，主要受二者的时空采样能力限制。

3)地下部分。

首先，不同埋深的分层标能够精确获取土层的压缩量，判断沉降发生的位置，并通过浅部和深部的变形特征区分土体自然压实和抽采地下水引起的土体压缩。其次，依据地下水位多层自动化监测数据判断因不同层位开采地下水形成的沉降。

2.1 数据集

2.1.1 SAR影像

本研究对象为江苏省连云港市。该区域拥有独特的沉积特征，为地面沉降的发育提供了先天条件。此外，在近些年《江苏沿海地区发展规划》大力推进的背景下，人类活动强度持续加大，进一步驱动该区域地面沉降的发生。为此，收集RADARSAT-2卫星2012—2016年37景30 m中分辨率以及2017—2020年51景5 m高分辨率时序SAR数据。由于两阶段数据均为同一颗卫星拍摄，因此除空间分辨率外，参数基本保持一致，这为后续地面沉降监测成果的长时序拼接奠定了基础。分别选取2014年9月18日和2018年8月4日的SAR影像作为主影像，并计算相邻影像拍摄的时间间隔以及相对于主影像的空间垂直基线距，如图2所示，从中看出SAR数据集时空基线分布比较合理，满足时序InSAR处理需求。

2.1.2 光学卫星影像

研究采用的光学卫星影像由行星(Planet)小卫星群采集获取，其空间分辨率为3 m, 能够完全覆盖连云港市全域。采集频率为每月一次，可以实现对研究区内地物及其变化过程进行动态跟踪，及时发现由地表变化引起的形变现象。

2.1.3 地表水准、GPS数据

依据研究区历史地面沉降空间发育情况布设水准和GPS点，主要集中在连云港东南区域，且以水准测量为主。

1)水准数据。

共收集到研究区2013年12月—2016年12月和2018年12月—2019年12月两个时段的水准测量成果，有效观测标志点数量分别为36和86个，用于开展对2014—2016年以及2019年InSAR形变速率结果的对比。

2)GPS数据。

该成果资料获取时段为2015年12月—2016年12月，相比水准点而言，其空间分布更为分散，监测点共有9处，主要用于评价2016年中分辨率InSAR形变速率成果的可靠性。

2.1.4 地下分层标组、水位监测数据

分层标组和地下水位监测井具有多层自动化观测能力，能够提供地下土层和水位的动态变化数据。

1)分层标组。

收集到徐圩和燕尾港两处分层标组中的地面标2016年12月1日—2020年12月1的形变发育曲线，观测频率为半年一次，用于与高分InSAR形变过程的一致性分析。除地面标外，这两座分层标组还包含1个基岩标，5个分层标，以徐圩为例，5个分层标的埋深分别为18.98 m(F1)、46.78 m(F2)、62.58 m(F3)、66.58 m(F4)和94.18 m(F5)。

2)地下水监测井。

收集到燕尾港和灌南地区两处地下水位监测井观测数据。其中，堆沟港监测井采集到2017年9月1日—2020年12月1日共40次Ⅱ承压水水位数据；灌南监测井自2017年3月1日—2020年12月1日共采集43次Ⅱ承压水水位数据。该资料主要用于与InSAR形变演化趋势的耦合分析。

2.2 SAR数据处理

2.2.1 主影像选取

主影像以时间和空间基线作为选取标准，同时考虑季节对相干性的影响。分别以影像集中的一景为参考计算其与其他影像之间的空间基线与时间基线和，并设定时空基线等权重，将两者相加值最小(无量纲)的参考影像确定为主影像，见式(1)。

ΤB ⊥ i = ∑ j≠i Ν B ⊥ i,j + ∑ j≠i Ν Τ i,j ⋅w i (1)

式中：TBi⊥代表第i期影像的时空基线质量评分，分数越低说明质量越好；N代表影像数量；Bi,j⊥和Ti,j分别代表第i期影像与第j期影像的空间垂直基线距与时间间隔；wi代表参考影像的季节权重，且 0

表1 SAR影像基本参数Tab.1 Basic Parameters of SAR Images

2.2.2 多主影像干涉组合

多主影像干涉组合的选择需要通过合理的时空基线阈值进行控制。阈值越小，一方面会造成干涉网络无法充分连接，形成孤立子网不便于后续的处理，另一方面会降低形变信息在干涉相位中的占比成分，损失监测精度；阈值越大，将会引入质量较差的干涉对，同样对结果精度造成不利影响。考虑到夏季影像受地物、大气等失相干现象影响更为明显，为此，根据季节设置时间基线阈值，见式(2)。

(1-wi)·(1-wj)·Tthre>Ti,j (2)

式中：Tthre代表时间基线固定阈值。处理过程中将空间基线阈值设为300 m; 时间基线固定阈值设为150 d。

2.2.3 相干目标提取

多主影像干涉构型决定了除相干点目标外，质量较好的分布式目标同样可以作为地表形变信息获取的载体。本文采用幅度差分离差+相位一致性分析+多主影像干涉网络闭合环残差检测的多层级组合方法提取相干目标。

与幅度离差指数相同[12],幅度差分离差指数同样通过SAR幅度影像初选相干目标[13],不同之处在于后者能够兼顾非城市区域分布式目标的提取，因而更加适合多主影像干涉处理。DΔ,A计算过程见式(3)。

D Δ,A = σ ΔA μ A = ∑ i=1 Κ (| Μ i |-| S i |) 2 Κ 1 2Κ ∑ i=1 Κ | S i | (3)

式中：σΔ,A和μA分别为主、辅影像幅度差的标准偏差和均值；Mi和Si分别为第i个干涉对主、辅影像的幅度；K代表干涉组合个数。DΔ,A的阈值设置与干涉对时间间隔以及影像总体时间跨度有关。结合使用的数据情况，分别设置中、高分辨率影像处理的DΔ,A为0.75和0.6。

相位空间一致性分析方法用于对上步提取的初始相干目标作进一步精化[14]。它依赖空间一定范围内相干目标相位存在相关性的假设，组合低通和带通滤波器实现空间相关相位的估计，通过对残余相位建立模型达到评价每个目标质量的目的，相位一致性因子γx见式(4)。

γ x = 1 Κ | ∑ i=1 Κ e -1 (φ x,i - φ ¯ x,i -Δφ θ,x,i u ) | (4)

式中：φx,i为相干目标x第i幅的差分干涉相位；φ ¯ x,i为由滤波器估计得到的空间相关相位；Δφ θ,x,i u 为空间非相关视角误差，通过其与空间垂直基线的关系进行计算[15]。

多主影像干涉网络闭合环残差通常用来检测相位解缠错误[16-17],而相干目标质量不佳是造成解缠处理出现问题的原因之一，因此，可以通过对闭合环残差的分析再次精化相干目标，见式(5)。

ε= ∑ i=1 Η U W{φ x,i,ε } (5)

式中：UW为相位解缠操作算子；H为单条闭合环上的干涉图总数；φx,i,ε 为闭合环的每个干涉相位。由于受噪声及多视处理等非线性操作的影响，即使无解缠错误发生的情况下，ε也不会为0。本文处理中将两个时段分别设为 | ε |<1 和 | ε |<0.5 。

2.2.4 地表形变信息获取

先分阶段各自解算得到形变量和各年度形变速率，而后利用同名点将两时段形变量进行拼接，得到研究区长时序地表形变演化过程。

1)分阶段形变量解算。

通过三维解缠算法对差分干涉相位进行周期还原[13],为降低解缠难度，处理前先将相位一致性分析中得到的空间非相关视角误差项从中减去，如式(6)所示。

UW{φx,i-Δφ θ,x,i u }=φdef, x,i+φatm, x,i+φorb, x,i+Δφ θ,x,i c +φnoise, x,i (6)

式中：φdef, x,i为雷达视向形变相位；φatm, x,i为大气传播延迟相位；φorb, x,i为残余轨道相位；Δφ θ,x,i c 为空间相关视角误差相位；φnoise, x,i为噪声相位。在得到解缠相位后，采用加权最小二乘法将多主解缠相位转为单主形式；其次通过二次多项式趋势面拟合得到并去除残余轨道信息[18];最后在时间域使用低通滤波得到地表形变序列[19-20]。其中，阶段一的低通滤波时间窗口设置为150 d; 阶段二设置为90 d。

2)两阶段形变量拼接。

该过程主要分为空间同名点查找和时序形变曲线连接两部分。首先，以阶段一提取的各相干目标位置为圆心，设置50 m半径查找两阶段同名点，采用距离加权平均对阶段二落入阶段一同一点位附近的相干目标地表形变信息进行处理。其次，利用式(7)计算两阶段中间空白时刻每日的形变量Def1,2。

式中：下标beg、end分别表示第一期和最后一期影像；上标1和2分别表示监测阶段。将每日形变量与两阶段空白时段相乘并加至2017—2020年地表形变序列上，便得到2012—2020年长时序相干目标形变信息。

3.1 地面沉降时空变化特征

连云港市2012—2020年逐年度的地表形变速率显示，监测时段内研究区地面沉降主要位于连云区、灌云以及灌南县的部分区域。其中，连云区集中在猴嘴街道、连云新城以及徐圩新区等地，地面沉降的主要特征表现为速率值大、普遍超过20 mm/a, 且呈零星分布、处于持续发育状态。灌云县地面沉降共有两处，分别表现为“一增一减”,改善区域集中在沿海的燕尾港化工园区附近，2012—2013年处于快速沉降状态，该时期漏斗中心处速率超过100 mm/a; 2014—2017年沉降逐渐趋缓；2018—2020年随着地下水位的逐渐回升，出现了轻微反弹迹象。新增区域位于圩丰镇-五图河农场一带，2018年出现沉降；2019年沉降范围和量值迅速增加，漏斗中心速率超过100 mm/a; 2020年略有减弱。灌南县地面沉降分为两块，一块位于田楼镇-堆沟港镇的狭长区域，且形成了以堆沟港镇工业园区为中心的剧烈沉降漏斗，2012—2016年处于快速沉降阶段，该时期中心处速率超过65 mm/a; 2017—2020年为沉降减缓及地表轻微反弹阶段。另一块位于灌南县城区的东北位置，2012—2016年沉降速率逐渐加大；2017—2020年快速减缓并逐步趋稳。

3.2 监测网数据对比分析

利用上述人工标志点监测数据分别与两阶段InSAR成果进行对比分析，验证前已将InSAR成果由雷达视向转至垂向并完成了基准偏差纠正。

1)2012—2016年。

受SAR影像分辨率的影响，该阶段InSAR相干点数量较少，因此在查找对比点对时，以人工标志点为中心，将搜索半径距离设置为80 m, 分别提取到9组水准-InSAR和7组GPS-InSAR点对，它们在地表形变趋势上的对比情况如图3所示，可以看到吻合程度较好。表2列出了该阶段点对数值的差异统计情况，两类对比的均方根误差分别为4.8和4.2 mm/a, 除个别点位偏差绝对值超过5 mm/a外，总体数值较为接近。

2)2017—2020年。

高分辨率InSAR结果在相干点数量上有明显提升，对比时以人工标志点为中心，设置空间搜索半径50 m, 同时顾及该范围内建筑物形变可能对验证结果造成的影响，尽量剔除自身形变特征明显的相干目标，最终提取到63组水准-InSAR点对，趋势对比如图4(a)所示。可以看到，两者不仅在形变速率较小的情况下具有较好的一致性，在一些局部快速沉降漏斗附近区域的监测数值同样较为接近，水准-InSAR对比的均方根误差为4.7 mm/a。此外，将该阶段InSAR形变曲线与徐圩和灌云两处分层标组中的地面标数据进行对比，见图4(b)。可以看到，除徐圩地面标在个别时间点上的数值与InSAR偏差超过20 mm外，两者的总体形变趋势相关程度较高。

表2 两阶段点对数值差异统计Tab.2 Statistics of Numerical Difference of Point Pairs in Two Stages

3.3 地面沉降驱动因素分析

3.3.1 自然因素

研究区表层软土分布广泛，且厚度自西北向东南沿海递增，如图5所示。其中，连云新城大部分区域软土厚度超过9 m; 而南部的徐圩新区软土平均厚度大于12 m, 部分区域超过20 m, 呈零星状分布。此外，研究区软土多为具高压缩、大孔隙比的海相沉积类软弱淤泥质土，不仅含水量高、强度低，而且灵敏度高、固结程度低，压缩系数大，因此具有显著的流变特性和弹塑性变形能力。将研究区2012—2020年的累计沉降量与软土厚度超过10 m的区域进行叠加显示，可以看到沿海地区呈零星状发育的地面沉降主要位于软土较厚区域，两者在空间分布上吻合度较高。图5(a)和图5(b)的分层标监测数据反映沿海地区地面沉降主要集中在地面标-F1标之间，即在20 m以浅的软土区发育。

3.3.2 人为因素

3.3.2.1 工程活动

工程活动进一步加速了研究区浅表层软土的固结压实，通过基础结构埋深较浅的建、构筑物发生明显沉降得以体现。如图6(a)和图6(b)所示，连云区猴嘴街道的地面沉降不均匀特征明显且量值较大，主要发生在基础结构较浅的道路上，除软土引起的自然沉降和车辆动荷载外，还包含道路周边建筑施工引起的附加沉降，图6(c)和图6(d)分别为2016年1月和2020年12月获取的两期光学卫星影像，通过对比可以发现，该区域在监测时段内正处于高强度开发阶段，工程建设活动较为频繁。

3.3.2.2 深层地下水开采

连云港地区深层地下水的开采分为满足居民日常生活和化工园区工业生产两方面用途，前者引发的地面沉降在空间上形成灌南县城—灌南堆沟港镇的区域性地面沉降带，后者形成以燕尾港和堆沟港两大化工园区为中心的沉降漏斗，随着地表水城乡区域供水一体化工程的建成使用，以及深层地下水井封填压采工作的实施，由深层地下水开采引发的地面沉降现象得到明显改善。自2015—2020年底，连云港累计完成封填(存)深层地下水井408眼，深层地下水压采量累计9.5×106 m3,超采区域大部分水位埋深回升至25 m限采红线以上。连云港市2016—2020年的Ⅱ、Ⅲ承压地下水开采量如表3所示。

表3 连云港市Ⅱ、Ⅲ承压地下水开采量Tab.3 Mining Volume of Ⅱ and Ⅲ Confined Groundwater in Lianyungang City

在燕尾港和堆沟港化工园区地下水超采漏斗区完成永久填埋水井84口、封存备用水井203口。此外，连云港南部于2016年底实现城乡区域供水一体化。图7给出了连云港南部2012—2020年的地面累计沉降量及两类深层地下水开采区的形变发育过程。

3.3.2.3 浅层地下水开采

连云港市浅层地下水开采以淡水养殖为主，与深层地下水开采相比，由其引发的地面沉降在空间上具有明显的局部漏斗特征，时间上表现出短期快速沉降的发育过程。以连云港南部具有代表性的堆沟港镇养殖区为例，详细分析浅层地下水开采引发的地面沉降发育规律。图8(a)给出了该养殖区2012—2020年逐年度地表形变速率以及沉降漏斗中心区域点位的时序形变曲线，可以清晰看到其地表变化分为3个阶段：前两个阶段分别为2012—2016年由深层地下水开采引起的区域地面沉降以及2016—2018年深层地下水井封填后出现的地表稳定及轻微反弹迹象；第三个阶段的沉降由水产养殖基地投产后大量抽采地下水引起，漏斗中心处3个阶段总体沉降速率超过130 mm/a。光学卫星影像和InSAR时序曲线对比分析显示，2019年初已有少量水产养殖大棚建成，并随即引发沉降现象，且随着养殖规模的不断扩大，大棚周边房屋进入快速沉降期。除确定该区域沉因外，时序光学影像还反映了养殖区的建设过程(图8(b)):2018年10月—2019年4月为初步建设期；2019年10月—2020年4月为快速建设期，养殖范围迅速扩大；2020年4月以后为建成期，总体规模保持不变。

经野外调研发现，连云港水产养殖场均以饲养南美白对虾为主，企业规模一般在133 333 m2左右，并以开采60 m以上的浅层地下水为主，部分区域伴有微承压+Ⅰ承压水，且每个养殖场内开采井数量在100～200口。此外，养殖基地周边人工地面及房屋受快速不均匀沉降影响，已出现不同程度的损毁，建筑物失稳风险明显增加，如图9所示。