百年煤城地表沉降融合PS/SBAS InSAR监测——以徐州市为例

姜德才，张继贤，张永红，吴宏安，康永辉

(中国测绘科学研究院，北京 100830)

百年煤城地表沉降融合PS/SBASInSAR监测
——以徐州市为例

姜德才，张继贤，张永红，吴宏安，康永辉

(中国测绘科学研究院，北京 100830)

融合PS/SBAS算法成为当前InSAR技术的研究热点，本文通过融合PS/SBAS的优势，利用128景ALOS PALSAR影像和24景Radarsat-2影像获取了百年煤城徐州地区2007—2015年的地表形变场及变化趋势。试验结果表明：①2007—2011年，徐州主要有4个沉降区，分布在沛县、丰县、铜山区和贾汪区。②2012—2015年，沛县、丰县和铜山区地表沉降范围有所扩大；贾汪区地表沉降范围和速率明显减小，说明近年来贾汪区重点治理采矿塌陷工作已初步见效；另外，睢宁县县城出现明显沉降。③已有文献对沛县大屯中心区的沉降监测结果验证了试验结果是可靠的。④首次提取了徐州近8年间地表沉降的信息，可为该地区开展地表沉降调查监测与防治工作提供参考；试验方法可为全国开展地表沉降调查监测工作提供参考。

雷达差分干涉测量;煤城;地表沉降;徐州

随着我国人口的持续增长和社会经济的高速发展，地下水、煤、石油或天然气等自然资源的大量开采已步入常态化，超高层建筑、地铁已成为大城市不可或缺的组成部分。这些人类活动诱发并加剧了地表沉降。截至2009年，我国有50多个城市出现了较严重的地表沉降[1]。徐州市位于华北平原的东南部，作为中国重要的煤矿产地，由于近百年以来煤炭资源的开采以及近30年以来地下水的超采，已造成了不同程度的地表沉降。

传统的地表沉降监测主要是采用水准测量，缺点是成本高、周期长、区域小。对于较大区域地表沉降监测可采用全球定位系统(GPS)进行测量,但是GPS监测站空间密度稀疏，并且成本高，不适宜推广。经过数十年的发展，合成孔径雷达差分干涉测量技术(DInSAR)，成为国际上公认的开展地表沉降监测的最有效、最先进手段。但是在实际应用中，传统DInSAR还存在着受时间、空间失相干和大气效应的制约[2-4]，为了解决这些问题，国内外学者探索出一些新型的技术,总体上可以概括为两类： 以PS-InSAR[5-6]为代表的单一主影像时间序列InSAR技术和以SBAS-InSAR[7-9]为代表的多主影像时间序列InSAR技术。前者多采用像元幅度选取点目标，优点是在提取非线性形变和消除大气相位影响方面算法简单[10]，缺点是要求足够多的SAR影像数据(通常25景以上)；后者多采用相干系数作为判定条件提取点目标，优点是可以显著增加有效干涉像对的数目，缺点是在反演线性形变时没有考虑大气影响。

PS/SBAS InSAR方法的比较和融合的研究成为当前热点问题。针对两种时间序列InSAR方法的特点，开展了融合研究。基于相干系数法、幅度差分离差法[11]和平均幅度法3种方法选取高相干点目标，线性形变提取通过构建局部Delaunay三角网连接高相干点目标，选择可靠连接边，从某一参考点开始对形变速率增量做积分，得到各点的线性形变速率。而非线性形变提取采用类似PS-InSAR的思路。

以往的InSAR研究都是针对局部地区、较短时段的监测，本试验依托地理国情监测项目，开展大区域、较长时段的地表沉降监测研究，利用128景ALOS PALSAR影像和24景Radarsat-2影像，反演了徐州大区域较长时段(2007—2015年)的地表形变场,分析了地表沉降的时空变化。使用试验区内泰山(毗邻云龙山)上保存良好的国家二等水准点D026对线性反演结果进行整体定标。通过文献[12—13]中沛县大屯中心区沉降结果对试验结果进行了验证。经国务院批准，国土资源部、水利部在2012年3月2日印发了《全国地面沉降防治规划(2011—2020年)》，试验方法可为全国开展地面沉降调查监测与防治工作提供参考。

1 研究方法与数据准备

1.1 研究方法

融合PS/SBAS InSAR核心算法简述如下：

1.1.1 高相干点目标提取

由于对高相干点目标的认识不同，采用的点目标提取方法也不尽相同。点目标提取的质量影响反演结果的精度。目前，高相干点目标提取主要分为两类：基于像元幅度提取和基于相干系数提取。前者多用于PS-InSAR中，后者多用于SBAS InSAR中。

幅度差分离差法是依据SAR影像同一像元在时间序列上强度值的离散特性进行点目标的识别，计算公式[11]如下

DΔA=σΔA/μA

(1)

式中，σΔA与μA分别是干涉图幅度的标准差和均值。此方法计算快捷，但要求有足够多的SAR影像(一般大于25景)以满足像元幅度的统计信息。另外，该方法只考虑点目标的幅度稳定性，忽略了点目标的强反射性及相干性，易造成点的误选，如水体、冰面。

平均幅度法依据SAR影像强度值的均值进行点目标的筛选，确保所提取的点目标具有强后向散射特征，剔除后向散射较弱的目标被误选。此方法没有考虑点目标的相干性，易误选长时间序列失相干的点目标。

为提高点目标提取的正确率，针对平均相干系数、幅度差分离差和平均幅度3种方法的特点，提取同时满足初始阈值(经验值，相干系数大于0.3、平均幅度大于1.2、幅度差分离差小于0.5)的点作为PS候选点，将点叠加至平均幅度上，根据点的空间分布情况，调节各参数阈值，反复试验求取理想阈值。

1.1.2 线性形变反演

从M幅干涉相位图相位中减去平地相位和地形相位，即可获得M幅差分干涉图相位，第k幅差分干涉图中的第i个点目标的相位可表示如下

(2)

在此基础上，通过局部Delaunay三角网连接所有的点目标，构建整体较密、反演效果更佳的三角网[14]。当三角网两点间距在大气影响相关距离(一般约1～3 km)内时，那么这两点的大气影响相位近似相等。因此对于三角网上的任一条边，其两顶点i、j之间的相位差为线性模型

(3)

式中，Δv和Δξ分别是形变速率和高程误差增量，采用吴宏安等人提供的整体相位相干系数法求取Δv和Δξ。

在所有的边中，选择同时满足整体相位相干系数大于0.7，残余相位的标准差和均值分别小于1.2和1.0的边(均为经验值)作为可靠的连接边，从某一参考点开始对Δv和Δξ做积分，得到各点的线性形变速率和高程误差绝对量。

1.1.3 非线性形变反演

从原始差分干涉图相位减去模型相位式(4)后，得到点目标上的残余相位φres，它包括大气影响相位φatm，非线性形变相位φnon-linear及噪声相位φnoise。为了得到完整的形变信息，需要对解缠后的残余相位进行时空频谱特征分析[5-6,10]，以分离出非线性形变相位，再将得到的非线性形变相位与线性形变叠加可得到全部形变信息。

1.2 数据准备

选择时间跨度从2007年7月8日—2011年3月3日的128景ALOS PALSAR影像，2012年2月11日—2015年8月6日的24景Radarsat-2影像。数据格式都是单视复数据(SLC数据)。PALSAR影像工作频率采用L波段，成像中心入射角是38.7°，每景影像覆盖范围是64 km×64 km，研究区范围是徐州，8景影像能够全覆盖；Radarsat-2影像工作频率采用C波段，成像中心入射角是26.2°，影像覆盖面积约为2.9万km2，基本全覆盖研究区。不同SAR影像的覆盖范围如图1所示。为消除干涉相位中的地形影响，需要研究区的数字高程模型数据(DEM)模拟地形相位[15]。为此，下载了研究区约90 m分辨率的SRTM DEM数据。

图1 徐州地区地表沉降监测SAR影像覆盖范围

2 徐州市地表沉降分析

对PALSAR影像和Radarsat-2影像分别作10×2、5×1多视处理，地面分辨率分别为31 m×29 m、26 m×27 m。8景PALSAR影像可将徐州全覆盖，在每景覆盖区获取了16期时间序列影像，对每景覆盖区单独反演线性形变，8景覆盖区的线性反演结果按照大地坐标拼接。利用上一节中介绍的方法，选取高相干点目标时，最终确定的平均相干系数阈值、平均幅度阈值、幅度差分离差阈值和提取的PS点数见表1。

表1 3种方法的参数阈值与提取的PS点数

2.1 试验区地表沉降时空变化分析

基于融合PS/SBAS InSAR方法提取了徐州大区域较长时段(2007—2015年)的地表沉降信息。考虑山顶稳定性较好，利用试验区泰山(毗邻云龙山)山顶保存良好的国家二等水准点D026对监测结果进行整体定标，得到最终形变反演结果。图2、图3分别为徐州市2007年7月—2011年3月、2012年2月—2015年8月年平均沉降速率分布图。

图2 徐州市2007年7月—2011年3月年平均沉降速率分布(A区域为沛县—丰县，B区域为部分铜山区；背景为Radarsat-2的平均幅度)

2007—2011年，徐州主要有4个沉降区，分布在沛县、丰县、铜山区和贾汪区。2012—2015年，沛县—丰县(A区)、铜山区(B区)沉降范围有扩大趋势，下文将具体分析A区和B区的地表沉降空间分布；贾汪区沉降范围明显减小，沉降速率明显减缓，这说明近年来贾汪区对塌陷地的综合治理工作已初步见效；睢宁县城出现了明显的区域性沉降，为此，睢宁县自2013年8月启动了为期7个月的城区地表沉降专项调查。

图3 徐州市2012年2月—2015年8月年平均沉降速率分布(A区域为沛县—丰县，B区域为部分铜山区；背景为Radarsat-2的平均幅度)

下面将具体分析沛县—丰县、铜山区地表沉降的空间分布特征。

2.1.1 沛县—丰县地表沉降空间分布

图4为沛县—丰县地区两个时间段的地表沉降对比图。2007—2011年，沛县地表沉降分布在县城及以北的地区，在县城以南的地区地面略有抬升；主要沉降漏斗有4个，从北往南依次为：龙固镇龙固矿区、大屯镇张双楼矿区、刘河涯、县城西环路东；大屯中心区也出现地表沉降。丰县地表沉降分布在县城及城周边，沿县城中心往外，年平均沉降速率逐渐增大。

2012—2015年，两县沉降范围明显扩大，主要分布在沛县县城、大屯镇、鹿楼镇、丰县县城；三河尖矿区出现了明显沉降，形成了局部沉降漏斗；大屯中心区沉降加剧；刘河涯、张双楼矿区沉降速率略有减缓。

从地理位置上看，沛县虽然东临微山湖，但是淡水资源供需缺口较大。因全县工业(如县城西环路东)、农业生产(如鹿楼镇的大片耕地)及生活(如县城及大屯中心区)过度利用地下水[16]，沛县一般年份缺水1.0多亿m3，特枯年份缺水达1.8亿m3。沛县煤炭资源非常丰富，主要分布于县城以北，由图4可知，主要矿区都有沉降发生。丰县地表沉降主要是因过度开采地下水[16]。至2012年年底，丰县地区已形成了205 km2的漏斗区，水位逐年下降，漏斗区逐年扩大，被江苏省列为限采区，局部地区禁采。

图4 沛县、丰县地区两个时间段的地表沉降对比图(背景为Radarsat-2的平均幅度)

2.1.2 铜山区地表沉降空间分布

图5为铜山区两个时间段的地表沉降对比图。2007—2011年，该地区形成了若干个沉降漏斗：柳新镇冯庄村、庞庄—张小楼矿区、刘集镇车村矿区、刘集镇丁场村、大彭镇土楼新村。2011—2015年，该地区地表沉降范围略有扩大，主要分布在九里区；柳新镇的柳新矿区有了明显沉降；大彭镇土楼新村沉降中心向北偏移；车村矿区近年来开展了大规模的拆房填坑工作，导致该地区在本时段提取的PS点较少。

铜山区地表沉降主要是因开采岩溶水、煤矿引起的[17-18]，沉降已造成当地许多村庄房屋开裂，村民被迫迁移。且京台高速(G3)的部分路段(共计超过2 km)处在沉降区，应加强对该路段的保护。

2.2 试验区2012—2015年地表沉降时间变化

因为多景影像覆盖区的每景SAR影像的时间序列通常不一致，累积沉降量自然不同。多景覆盖区的形变结果镶嵌通常针对平均沉降速率[19]，不适用于累计沉降量。ALOS PALSAR共128景影像，平均每期8景，时间序列不完全一致。因此，基于24景Radarsat-2影像，在反演线性形变的基础上，继续反演非线性形变。将线性形变和非线性形变结果叠加，得到徐州2012—2015年的累计沉降量(如图6所示)。从图中可以看出，随时间变化，徐州5个主要沉降区逐渐显现出来，分别为：沛县大屯镇，-274 mm；丰县县城，-262 mm；铜山区刘集镇，-317 mm；贾汪区潘安湖，-115 mm；睢宁县县城，-113 mm。

图5 铜山区两个时间段的地表沉降对比图(背景为Radar sat-2的平均幅度)

图6 2012—2015年徐州市各时间段的累计沉降量

2.3 试验区地表沉降结果可靠性验证

分析沛县大屯中心区(如图4所示)的地表沉降空间分布特征，2007—2011年，大屯中心区大部分地区为沉降趋势，75%以上的区域沉降速率达到-5 mm/a，40%以上的区域的沉降速率在-10～-37 mm/a，最大沉降速率达到-37 mm/a。2011—2015年，85%以上的区域沉降速率达到-5 mm/a，65%以上的区域的沉降速率在-10～-51 mm/a，最大沉降速率为-51 mm/a。可见，大屯中心区地表沉降范围、沉降速率随时间有所增加。

张茂永等利用1988—2005年的徐州大屯中心水准观测数据，根据灰色理论建立模型预测2005—2010年最大沉降速率达到-31 mm/a[12]。何敏等利用两景ENVISAT ASAR数据(20070327，20090120)进行分析，监测结果表明，2007—2009年大屯中心区70%以上的区域的沉降速率达到5 mm/a，最大沉降速率达到-34 mm/a。

试验第一时段(2007-07-08—2011-03-03)的监测结果对比张茂永、何敏等的研究结果，地表沉降范围、年平均最大沉降速率基本一致，且随时间有所增加，验证了试验结果是可靠的。

3 结论与探讨

本融合PS/SBAS InSAR算法提取了百年煤城徐州大区域地表形变场，全面分析了2007—2015年徐州地表沉降的时空变化。监测结果表明：①2007—2011年，徐州主要有4个沉降区，分布在沛县、丰县、铜山区和贾汪区；②2012—2015年，沛县、丰县和铜山区地表沉降范围有所扩大，贾汪区地表沉降范围和速率明显减小，睢宁县城出现明显沉降；③已有文献对沛县大屯中心区的沉降监测结果验证了试验结果是可靠的；④试验结果可服务于资源型城市转型，试验方法为全国开展地理国情普查——地面沉降监测提供参考。

综合考虑徐州地区地表沉降的人为因素，主要是地下水、煤过度开采。徐州地区碳酸岩地层分布较广，多为岩溶发育[20]。岩溶水是徐州市区、铜山区、沛县和丰县的主要供水源，在干旱年份、季节，易造成地下水源的过量开采，诱发地面塌陷。软土零星分布于市区、刘集等地[20]，软土具有高压缩性、触变性等特点，易发生地面形变、地基滑移等灾害。徐州是全国主要煤矿能源基地，但煤炭资源的大规模开发利用给徐州带来了分布广泛的采煤塌陷地。

由试验结果可知，徐州地区现有5个沉降区，除贾汪区地表沉降范围和速率有所减缓外，其他4个沉降区的沉降范围都有扩大的趋势。贾汪区政府采取治理复垦与生态修复相结合治理采煤塌陷地，并重点治理了潘安湖采煤塌陷地。通过引运河水形成湖泊、湿地来进行采煤塌陷地整治和生态修复，将占地1160 hm2的潘安湖采煤塌陷地建成了综合型国家4A级湿地景区。借鉴贾汪区对塌陷地的综合治理经验，徐州应全面加强对已有地表沉降地区的治理和复垦工作；同时，加强水资源管理，合理利用地下水，引进其他补寄水源。

致谢： 日本宇宙航空研究开发机构(JAXA)提供了ALOS PALSAR数据，加拿大MDA公司提供了RADARSAT-2数据，国家科学数据服务平台提供了SRTM数据，在此一并表示感谢。

[1] 王寒梅.上海市地面沉降风险评价体系及风险管理研究[D].上海：上海大学,2013：1-10.

[2] ZEBKER H A，VILLASENOR J.Decorrelation in Interferometric Radar Echoes[J].IEEE Transaction on Geoscience and Remote Sensing，1992,30(S)：950-959.

[3] ZEBKER H A,ROSEN P A and HENSLEY S.Atmospheric Effects in Interferometric Synthetic Aperture Radar Surface Deformation and Topographic Maps[J].Journal of Geophysical Research,1997,102(B4)：7547-7563.

[4] CHANG Z Q，GONG H L，ZHANG J F，et al.A Feasible Approach for Improving Accuracy of Ground Deformation Measured by D-InSAR[J].Journal of China University of Mining amp; Technology,2007,17(2):262-266.

[5] 张永红,张继贤,龚文瑜，等.基于SAR干涉点目标分析技术的城市地表形变监测[J].测绘学报，2009，38(6)：482-487.

[6] 张永红,吴宏安,孙广通.时间序列InSAR技术中的形变模型研究[J].测绘学报，2012，41(6)：864-869.

[7] MORA O，MALLORQUI J J，BROQUETAS A.Linear and Nonlinear Terrain Deformation Maps from a Reduced Set of Interferometric SAR Images[J].IEEE Transactions on Geoscience and Remote Sensing，2003，41(10):2243-2253.

[8] CASU F，MANZO M，PEPE A，et al.SBAS-DInSAR Analysis of Very Extended Areas：First Results on a 60000 km2Test Site[J].IEEE Transactions on Geosciences and Remote Sensing Letters，2008，5(3)：438-442.

[9] LIU G X，JIA H G，NIE Y J，et al.Detecting Subsidence in Coastal Areas by Ultrashort-Baseline TCPInSAR on the Time Series of High-Resolution TerraSAR-X Images[J].IEEE Transactions on Geoscience and Remote Sensing，2014，52(4)：1911- 1923.

[10] 吴宏安，张永红，陈晓勇,等.基于小基线DInSAR技术监测太原市2003～2009年地表形变场[J].地球物理学报，2011，54(3)：673-680.

[11] HOOPER A.A Multi-temporal InSAR Method Incorporating Both Persistent Scatterer and Small Baseline Approaches[J].Geophysical Research Letters，2008，35(16)：96-106.

[12] 张茂永,孟磊,王传林,等.徐州大屯中心区地面沉降趋势预测[J].地质灾害与环境保护，2007，18(13)：11-13.

[13] 何敏,陆晓燕,何秀凤.利用D-InSAR二轨法监测徐州大屯中心区地表形变[J].地理空间信息，2011，9(5):3-6.

[14] ZHANG L，DING X L，LU Z.Modeling PSInSAR Time Series without Phase Unwrapping[J].IEEE Transactions on Geoscience and Remote Sensing，2011,49(1)：547-556.

[15] 吴宏安，张红，王超，等.基于SRTM DEM的InSAR高分辨率山区地表高程重建算法(英文)[J].遥感学报，2009，13(1)：145-151.

[16] 蔺宗宗，吴海波，田玉新,等.徐州市地面沉降成因分析及防治对策[J].中国科技博览，2010(26):256-256.

[17] 花修权，黄敬军，缪世贤，等.徐州市地质灾害分布特征及成因分析[J].地质灾害与环境保护,2015，26(2):74-80.

[18] 花修权，黄敬军，魏永耀，等.江苏徐州城市地质调查需求分析及关键问题[J].地质学刊，2015，39(2):322-328.

[19] ZHANG Y H，WU H A，KANG Y H，et al.Ground Subsidence in the Beijing-Tianjin-Hebei Region from 1992 to 2014 Revealed by Multiple SAR Stacks[J].Remote Sensing，2016,8(8).

[20] 邢雪,周丹,罗跃.江苏徐州岩溶塌陷及其防治对策[J].中国地质灾害与防治学报，2014，25(4)：51-58.

GroundDeformationoverCentury-longCoalMiningCityMonitoredthroughIncorporatingBothPSandSBASInSAR:ACaseStudyofXuzhou

JIANG Decai,ZHANG Jixian,ZHANG Yonghong,WU Hongan,KANG Yonghui

(Chinese Academy of Surveying and Mapping, Beijing 100830, China)

The algorithm fusion of permanent scatterers interferometric SAR (PS InSAR) and small baseline subset interferometric SAR (SBAS InSAR) has become a hot spot of InSAR technology. By fusing the advantages of PS/SBAS, 128 ALOS PALSAR images (from 2007 to 2011) and 24 Radarsat-2 images (from 2012 to 2015) are used to extract the ground deformation field and trend in recent 8 years over Xuzhou region with century-long coal mining history. The experiment results demonstrate that:① From 2007 to 2011, four main subsidence areas have been found in Xuzhou, locating in Peixian County, Fengxian County, Tongshan District and Jiawang District. ② From 2012 to 2015, the subsidence range of Peixian County, Fengxian County and Tongshan District is expanding; however, the subsidence range and rate of Jiawang District decrease significantly. The phenomenon shows that the governance of mining subsidence in Jiawang District has achieved good results in recent years. In addition, there comes into apparent subsidence in Suining County. ③ The results of this experiment has been validated by the published research results on ground subsidence in Datun Center District, Peixian County. ④ The results represents the first-ever subsidence information in the recent 8 years in Xuzhou, which can be applied to monitor and control the subsidence in this area.The monitoring result is of high significance for monitoring and control, prevention and mitigation of ground subsidence disaster over the region.

differential interferometric SAR(DInSAR); coal mining city; ground deformation; Xuzhou city

P237

A

0494-0911(2017)01-0058-07

姜德才，张继贤，张永红，等.百年煤城地表沉降融合PS/SBAS InSAR监测——以徐州市为例[J].测绘通报，2017(1)：58-64.

10.13474/j.cnki.11-2246.2017.0013.

2016-04-20

国家973计划(2012CB719905);自然科学基金(41271430；41304010)

姜德才(1990—),男,硕士生，研究方向为合成孔径雷达干涉测量。E-mail:jiangdc_rs@sina.com

张永红。E-mail：yhzhang@casm.ac.cn