北京来广营地区地面沉降时空演化特征

杨翠玉，王彦兵，赵亚丽，郭亚航，黄世强

(1.首都师范大学 资源环境与旅游学院，北京 100048；2.首都师范大学 三维信息获取与应用教育部重点实验室，北京100048；3.首都师范大学 城市环境过程与数字模拟国家重点实验室培育基地，北京 100048)

0 引言

由于长期超量开采地下水、城市建筑群容积率增加等原因，北京地区地面沉降加剧。从20世纪60年代，北京平原区地面沉降呈快速发展，已形成昌平沙河—八仙庄、朝阳八里庄—大郊亭、朝阳—来广营、顺义—平各庄、大兴榆垡—礼贤5大沉降区[1]。2004—2009年间，各沉降区的累积沉降量分别达到昌平沙河—八仙庄为272.8 mm，朝阳八里庄—大郊亭为460.5 mm，朝阳—来广营为434.0 mm，顺义—平各庄为281.6 mm，大兴榆垡—礼贤为435.8 mm[2]。2012—2015年间，朝阳—来广营沉降区、朝阳八里庄—大郊亭沉降区的沉降速率分别达到119 mm/a、115 mm/a[3]。

针对地面沉降的研究，目前主要围绕遥感监测技术、演化特征及成因分析等方面进行。遥感监测技术目前以差分干涉测量(differential interferometry synthetic aperture radar， D-InSAR)、永久散射体干涉测量(presistent scatterers for sar interferometry，PS-InSAR)、小基线干涉测量(small baseline subset，SBAS)为主要的地面沉降监测方法[4-6]。针对北京地区的地面沉降演化与成因，大量学者进行了研究[2，7-13]。张雯等[2]利用统计方法对北京5大沉降区进行了时序演化特征分析，发现沉降区的沉降量均呈平稳上升趋势，沉降演化的空间特征不明显。周超凡等[12]利用全局空间自相关指数和局部空间自相关指数对北京地面沉降格局特征进行了分析，结果表明,严重沉降区和轻微沉降区在空间上均存在集聚现象，时间特性不明显。Zuo等[13]利用标准差椭圆方法，研究发现北京东部地区的地面沉降由西北向东南方向发展的空间演化特征，发现其时间演化特征较不明显。Li等[14]应用区域重心法对比分析了北京平原区地面沉降演化特征与地下水、可压缩层厚度之间的关系。

本文为进一步挖掘北京地面沉降的演化规律，提出利用区域重心法定量描述沉降演化特征。区域重心法是一种描述地理现象空间集聚程度的一种方法，它能对地理现象的时空演化特征进行定量描述，被广泛应用于人口、经济等领域的空间分布演化规律研究[15-17]。从地理学角度，地面沉降是一种地理现象，用区域重心法可定量描述其时空演化特征。

1 区域重心法

重心是物体在重力场中保持平衡的点。在地学领域，重心的位置可以通过地理坐标和地理空间现象的组合计算获取，可以用来定量描述地理现象的空间分布情况及演化规律[18-19]。本文采用区域重心法定量描述区域地面沉降演化特征。

假设研究区的沉降监测点分布大致均匀，每个沉降监测点的属性描述包括沉降量S和地理坐标(Xi，Yi)，则该区域的地面沉降的重心地理坐标如式(1)所示。

(1)

式中：X、Y代表该区域地面沉降重心的地理坐标；Xi、Yi表示第i个沉降监测点的地理坐标；Si是第i个沉降监测点的沉降量。

根据公式(1)，当研究区沉降监测点空间分布均匀，且均匀沉降时，重心与几何中心位置相同；当沉降监测点空间分布均匀，但沉降不均匀时，则重心偏离几何中心。因此在后续的计算中，首先要对沉降监测点进行处理，使得其空间分布均匀。

为描述区域地面沉降重心的偏移情况，本文以重心移动的距离Lt2-t1和方位角θt2-t1定量化表达。

(2)

(3)

式中：t1、t2表示2个不同的时间；Lt2-t1表示由t1到t2时段的区域地面沉降重心迁移的距离；θt2-t1表示由t1到t2时段的区域地面沉降重心迁移的方位角[15]。

2 地面沉降数据处理及精度验证

本文以北京来广营沉降区为研究区域(图1)，利用2004年1月至2010年8月的49景ENVISAT ASAR数据和2010年11月至2015年12月的40景Radarsat-2数据，基于PS-InSAR技术获取地面沉降信息。

首先，在PS-InSAR处理过程中，由于SAR影像范围受到轨道位置、视角、大气状况等条件的影响，需要对研究区的SAR影像进行统一裁剪和配准等处理；其次以辅影像和主影像构成的干涉影像对为基础，利用原始数据的共轭相乘获取干涉相位，根据轨道数据和外部数字高程模型(digital elevation model， DEM)数据计算并去除地形和地平相位；然后通过振幅离差指数大于0.65且时间相干系数大于0.70筛选沉降监测点(来广营沉降区(128 km2)共获取3 002个沉降监测点)，并利用大气延迟相位(atmospheric delay phase， APS)在空间上低频和时间上高频，而非线性形变在空间高频时间低频的特征，通过滤波的方法将大气延迟相位去除，从而获取来广营地区2004年1月至2015年12月的地面沉降数据(图1)；接着对获取的沉降监测点进行均匀化处理；最后用区域重心法定量计算地面沉降重心的移动量与移动方向，分析该地区的地面沉降时空演化特征，并结合水文地质剖面对成因进行了分析。具体流程如图2所示。

注：该图基于自然资源部标准地图服务网站下载的审图号为GS(2019)3333号的标准地图制作，底图无修改。 图1 2004—2015年北京地区年均沉降速率图

图2 数据处理流程图

根据获取的地面沉降数据初步发现，北京地区存在明显的不均匀沉降，沉降严重的区域主要集中在东部的朝阳、通州以及顺义区，最大年均沉降速率为109 mm/a；主城区、丰台区和石景山区相对比较稳定，年均沉降速率小于10 mm/a。

为验证PS-InSAR方法获取的地面沉降数据的精度，本文用2005—2013年间同一时段、相同点位的北京市水文地质大队以水准测量方法监测的沉降数据进行对比分析。

首先，以水准监测点为圆心，选取半径200 m范围内的所有PS点的年均沉降速率作为该点的PS-InSAR监测值；然后计算研究区内8个水准监测点与对应PS-InSAR监测值的相对误差，并绘制相关关系图(图3)。结果表明，PS-InSAR监测值与水准监测值相对误差在10 mm/a以内，相关系数达到0.997，2种监测结果具有较好的一致性，反映了PS-InSAR监测结果的可靠性。

图3 2005—2013年InSAR监测结果与水准结果对比

3 地面沉降演化特征分析

3.1 地面沉降演化趋势分析

在PS-InSAR监测的来广营地区地面沉降数据基础上，本文分析了2004年1月至2015年12月该区域地面沉降的逐年演化趋势(图4)。通过分析可知，来广营地区最大累积沉降量达1 278 mm，年均沉降速率为-98～-25 mm/a。来广营地区地面沉降漏斗向东部发展，从2008年开始，在研究区的西北、东南2处形成2个沉降漏斗中心，以较高的沉降速率累积，并且东南的沉降漏斗中心沉降速率比西北的沉降漏斗中心高，存在明显的不均匀沉降。

注：该图基于自然资源部标准地图服务网站下载的审图号为GS(2019)3333号的标准地图制作，底图无修改。图4 2004—2015年来广营逐年地面沉降演化过程

3.2 区域重心法分析

根据图4的地面沉降演化过程，可以发现地面沉降的发展趋势，但难以对其进行定量表达。本文采用区域重心法，以累积沉降量为数据源，进一步定量描述来广营地面沉降的时空演化特征。

首先，对沉降监测数据进行均匀化处理。根据式(1)可知，区域重心的位置受研究区沉降监测点的空间分布和沉降监测点的沉降量二者的共同影响，为保证计算的准确性，研究的沉降监测点要均匀分布。但是，由于PS-InSAR算法获取的沉降监测点通常分布在人工建筑物、路灯、桥梁、人工角反射器等相干性较高的硬目标上，所以沉降监测点的空间不均匀分布。本文通过对研究区的沉降监测点空间插值，以100 m为间隔均匀提取样本点，使得点位分布均匀。

其次，根据式(2)、式(3)计算2004—2015年来广营地面沉降重心的迁移距离与方位角(表1)，绘制区域地面沉降重心迁移轨迹图(图5)。根据图5和表1的结果发现来广营沉降区具有如下的时空演化特征。

1)2004—2015年，来广营地面沉降重心由西北向东南方向迁移，累积迁移1 516.43 m，迁移的平均方位角为110.56°。

表1 2004—2015年来广营沉降重心迁移量统计表

2)来广营地区地面沉降重心迁移速度呈现逐年递减趋势，2014年地面沉降重心迁移速度达极小值，2015年沉降重心迁移速度有迅速回升的趋势，可能由于东南和西北地区地下水开采量发生变化导致。

3)来广营地面沉降重心以东西向迁移为主，南北向迁移为辅，2004年1月至2015年12月期间地面沉降重心累积向东迁移1 388.17 m，东西向的迁移速度在59～200 mm/a间波动，2007年达极大值，2014年达极小值；累积向南迁移586.73 m，南北向迁移速度在1～138 mm/a间波动，2005年达极大值，2014年达极小值。

4)来广营地面沉降重心东西向迁移速度的变化趋势和整体迁移速度的变化趋势保持一致，逐年递减，南北向迁移速度先缓慢递增，从2011年开始迅速减缓，到2014年南北向迁移速度仅为1 mm/a，2015年东西向和南北向的迁移速度都有回升趋势。

进一步基于区域重心法对来广营地区不同年份季度之间地面沉降重心的迁移距离和方位角进行了计算，发现其迁移量较小，且规律性不明显。

注：该图基于自然资源部标准地图服务网站下载的审图号为GS(2019)3333号的标准地图制作，底图无修改。图5 2004—2015年来广营地面沉降重心迁移轨迹图

由于地面沉降受地质环境条件、地下水开采强度、城市开发建设等因素的影响[20]。本文从中国地质科学院水文地质环境地质研究所编著的《华北平原地下水可持续利用图集》中获取了昌平区-通州区的水文地质剖面图(图6)，进一步分析来广营的地面沉降成因。

注：该图基于自然资源部标准地图服务网站下载的审图号为GS(2019)3333号的标准地图制作，底图无修改。图6 水文地质剖面图[21]

水文地质剖面图以线段AA′为剖面线，横跨来广营地区(图6)。以AA′为剖面线同时绘制累积沉降量的剖面图(图7)。通过图6 和图7的对比分析可知：来广营地区水文地质条件基本一致，岩性以含水砂层和粉土为主，水文地质剖面沿线累积沉降量与深层地下水水头具有较好的对应关系，与浅层地下水水位的关系不明显。结合现有研究成果，北京地面沉降漏斗形成的主要原因是地下水超采[2，7，22]，且开采深度逐渐增加[23]。此外，结合地面沉降监测数据和天地图影像，发现来广营地区地面建筑均以低矮村庄为主，沉降漏斗重心的变化与地物的相关性不明显。因此，来广营地区不均匀沉降可能是由于深层地下水开采差异导致。

图7 水文地质剖面沿线累积沉降量剖面图

4 结束语

本文利用永久散射体干涉测量技术(PS-InSAR)获取了来广营地面沉降数据，通过区域重心法对来广营地面沉降时空演化特征进行了定量分析，得到如下结论。

1)2004年1月到2015年12月，北京来广营地面沉降年均沉降速率为-98～-25 mm/a，空间上呈现出明显的不均匀沉降，从2008年开始，来广营地区西北、东南2处形成2个沉降漏斗，且东南沉降漏斗的沉降速率明显高于西北沉降漏斗。

2)来广营地面沉降重心逐年由西北向东南方向移动，累积迁移1 516.43 m，迁移的平均方位角为110.56°；其中以东西向迁移为主，南北向迁移为辅，期间累积向东迁移1 388.17 m，迁移速度在59～200 mm/a之间波动，累积向南迁移586.73 m，迁移速度在1～138 mm/a间波动。

3)来广营地区不均匀沉降可能由于深层地下水开采差异导致。

来广营地区覆盖北京地铁首都机场线、首都机场高速公路等多个重要轨道交通线，随着不均匀沉降的加剧，需要对该沉降区内轨道交通的安全运营进行重点关注。另外，本次研究侧重于从地理学的角度分析来广营地区地面沉降年际的时空演化特征，年内季度之间地面沉降重心迁移变化规律不明显。此外，受限于水文地质数据，对于地面沉降时空演化的成因后续须进一步展开研究。