基于时序InSAR技术的矿区形变监测与分析

陈志达，庞校光

（1.绍兴市上虞区自然资源监测中心，浙江 绍兴 312399；2.中国科学院精密测量科学与技术创新研究院 大地测量与地球动力学国家重点实验室，湖北 武汉 430077；3.中国科学院大学，北京 100049）

InSAR技术是20世纪后期发展起来的卫星遥感新方法，可以获取高精度地形和地表形变信息[1-2]。与常规地面水准测量、全球定位系统（GPS）、光学遥感等测绘手段相比，InSAR技术进行地面形变监测时优势明显[3-5]，且已在地表形变[6-8]、矿区[9-10]、滑坡[11]及地震形变[12]等研究与监测中得到了广泛的应用。

矿区开采引起的地面沉降是非常复杂的过程，其危害十分严重，研究矿区塌陷具有重要的意义[13-14]。2011年，Jiang[15]等利用时序InSAR对矿区煤火引起的地表形变进行了监测，结果显示时序InSAR技术对矿区煤火引起的地表形变具有可靠的应用潜力。2014年，Dong[16]等利用差分干涉测量短基线集时序分析技术（SBAS-InSAR）对淮南地下煤矿开采引起的地面沉降进行时空分析，证明了L波段InSAR可用于监测矿区地下开采活动引起的地面沉降。同年，Zhao[17]等通过InSAR技术对大同不同采矿沉降阶段及时间演变进行分析，发现1992—2011年间约80%煤矿发生形变。2015年，LI[3]等对徐州矿区进行实验，再次验证了In-SAR技术用于矿区形变监测结果的可靠性。2019年，朱建军[18]利用InSAR技术求取了煤矿开采沉陷的三维形变场，并进行了实验模拟。

山西省大同市浑源县煤炭资源丰富，从20世纪70年代以来，煤矿无序开采导致周边地区大范围地面沉降，严重破坏了周围的生态环境，直接威胁着当地居民的生命及财产安全。为此，本文利用覆盖大同市浑源县46景SAR影像，采用时序InSAR技术进行分析处理，获取了2007—2013年间浑源县平均年形变速率结果和时序结果，并根据实地调查资料分析了地面沉降时空变化特征及其与采空区分布的关系。

1 研究区概况及实验数据

1.1 研究区概况

浑源县（39.67°～39.83°N、113.72°～113.92°E）地处黄土高原的边缘地带，位于山西省东北部大同盆地东南边缘，东临广灵县，南邻繁峙县，西毗应县，北接云州区。浑源县地处桑干河支流浑河中上游，全县南北长57 km，东西宽约52 km，总面积1 968 km2。浑源县地形复杂，东、南、北三面环山，境内起伏较大，海拔最高的卧羊场山主峰海拔2 333 m，最低的浑河西部河床海拔1 026 m，相对高差1 307 m。境内山脉较多，海拔均在2 000 m以上，山势陡峭，沟壑纵横[19]。

浑源县矿产资源丰富，储量大，包含煤炭、花岗岩、铁矿等，已探明的煤矿储量约17亿t。自上世纪70年代以来，粗放式煤矿开采破坏了当地的生态环境，威胁着人民的生命财产安全[17,20]。本文选取整个浑源县为研究区域，图1为浑源县概况及数据覆盖情况。

图1 研究区域和数据覆盖情况

1.2 实验数据

ALOS-1于2006-01-24发射升空，是日本新一代的对地观测卫星。其携带的L波段SAR适用于地表破坏严重的煤矿开采区域，可有效提取浑源县矿区地表形变范围和形变量信息。RADARSAT-2是一颗搭载C波段传感器的高分辨率商用雷达卫星，于2007-12-14发射升空，是RADARSAT-1的后续卫星，在新的图像获取能力及性能方面，有长足进展。该卫星具有最高1 m高分辨率成像能力，多种极化方式。

本文 获取2007年6月—2011年1月36景ALOS PALSAR影 像 和2012年1月—2013年4月10景RADARSAT-2影像数据，以获取浑源县长期地表高程变化速率及时间序列。ALOS/PALSAR数据分辨率为10 m，入射角为34.3°。RADARSAT-2数据分辨率为30 m，入射角为21.25°。DEM数据采用NASA提供的SRTM，其空间分辨率为90 m，用以去除干涉对中的地形相位。此外，还搜集5大煤矿区域、99座小煤矿坐标以及煤矿的地质报告等煤矿开采资料，用于InSAR监测结果的验证与分析。通过获取浑源县地面形变范围和形变量信息，结合煤矿开采存档、野外实地调查等辅助数据，进行相互验证。

2 SBAS-InSAR数据处理

2.1 SBAS-InSAR基本原理

SBAS-InSAR方法原理是将D-InSAR形变结果作为观测值，通过SVD方法将多个小基线集联合进行最小二乘解求解，获得的形变时间序列可有效反演出地表形变范围和形变量信息。对某一区域不同时间获取N+1幅SAR影像，选择主影像后将其他影像配准至主影像，为避免时空失相干造成观测误差及影响监测结果可靠性，设置时空基线阈值以保证干涉对的有效性，后利用外部DEM去除参考椭球相位与地形相位，对形变相位进行解缠及地理编码得到地理坐标系下研究区域形变结果[1-2,21]。两景SAR影像干涉相位可用φint表示为：

式中，φflat为参考椭球相位；φtopo为地形相位；φatm为大气延迟相位；φorbit为轨道误差；φnoise为噪声相位；φdef为地表形变相位。对去除各项误差后的差分干涉对地表形变相位进行解缠后，可获取视线向的位移量研究区域各像元形变量。

2.2 数据处理

本文对ALOS/PALSAR与RADARSAT-2影像数据分别处理，干涉对主辅SAR数据精配准精度均控制在0.1个像元内，空间基线均为最大基线的45%，最大基线分别为1 984.26 m和-271.49 m，时间基线分别为365 d和180 d，最终分别生成70和24个干涉对。对各干涉对分别进行精密轨道数据优化空间基线，消除轨道误差，去除地平效应，利用外部SRTM-DEM去除地形相位，干涉图滤波，相位解缠，消除大气相位延迟等处理得到地表形变相位，对获取的形变干涉图采用多基线方法获取时间序列[22-23]，流程图如图2所示。

图2 SBAS-InSAR方法流程

3 结果与讨论

3.1 地表形变结果

利用覆盖浑源县的ALOS/PALSAR和RADARSAT-2卫星雷达遥感数据，采用SBAS-InSAR技术获取了浑源县地表形变结果，包括浑源县年均形变速率与重点区域的累积形变量。根据InSAR地表形变结果，结合99座小煤矿范围、5大煤矿位置和野外地面调查资料，圈定了49个采空区分布范围，如图3所示。针对49个采空区获取了采空区历年累计区域平均形变量统计结果，其时序形变过程如图4所示。

图4 浑源县49个形变区域中心2007-2013年间时序结果

从时序图上可知各采空区7 a间地表高程持续下降，最大累计形变达-211.21 cm。其中15、17、18、19、20等区域平均形变量在2010年前后发生了剧烈变化，从图3可知，5个区域位于黄花滩乡-大磁窑镇-青瓷窑乡交界处，属浑源县东部沉降群。上述5个区域自2010年后煤矿开采规模与速度均有所加剧，其中20号区域平均形变量为当年浑源县最为严重区域。2012年地面沉降进一步加剧，区域平均沉降量超过-50 cm的有3个区域，分别为20号的-62.01 cm、36号的-60.34 cm、41号的-56.32 cm。2013年平均沉降量超过50 cm的区域为13号和20号，其值分别为-53.33 cm、-52.45 cm。

图3 浑源县年均沉陷速率及沉陷区分布范围

为验证InSAR监测结果的可靠性，结合实地调查结果与煤矿开采记录资料等辅助数据对结果进行验证。实地调查了沙圪坨镇西坪村、大仁庄乡吴家圪坨村、蔡村镇白道村、青磁窑乡青磁窑村，如图5所示，结合形变时间序列结果可知7 a间累计形变量分别达-90.0 cm、-141.22 cm、-125.51 cm、-75.77 cm。受煤矿开采影响，煤矿区域周围居民遭遇不同程度的地面沉降危害，村庄房屋建筑出现不同程度受损情况。

图5 部分村庄房屋损毁情况

3.2 采空区形变分布特征

从图3可知，浑源县7 a间在东部、北部、西北部分别形成了三大沉降漏斗群，其中东部范围最广，包括大磁窑镇-黄花滩乡-青瓷窑乡北部-千佛岭乡北部-大仁庄乡中部等广阔区域，该区域分布着果子园、古瓷窑、白庄等众多煤矿。对各区域时序形变结果中历年不同形变量比例进行统计得到图6，从图中可知，随着时间变化，浑源县49个区域中历年形变量大的比重不断增加。2007—2013年间，沉降量＞-10 cm的区域从75.5%下降至2%，年沉降量超过-50 cm区域从无到有，其中2012年共3个，占6.1%；2013年2个，占4.1%。但从2007年看形变量较2008年大，分析其原因可能是当年数据从2007年6月起，仅半年数据，获取的形变量受季节性变化影响较大，难以准确反映出全年的形变规律；2013年最大沉降值较2012年更小，可能有两个原因，一是2013年采用的数据不足半年，由于数据量不足及季节性影响可能无法完整反映出浑源县2013年形变结果；另一个原因是2012年政府对煤矿开采加以控制，地面形变情况在2013年得到改善，针对原因不确定性，后续采用完整数据验证具体影响因素。

图6 浑源县历年不同形变量所占比例

根据时序InSAR结果可知浑源县地表形变在空间上呈散点式漏斗形状，面积小而多，均呈持续下沉状态。结合时序InSAR结果获取了浑源县采空区中心随时间变化的空间演变结果（图7）。从图上可知浑源县北部与西部在2007—2013年间中心形变发生明显形变，其中部分采空区中心形变具有明显阶段性特征，以北部41号采空区为例，该区形变以2011年为分界线；2007—2011年间呈稳步线性趋势变化，发生突然加速。浑源县东部在2007—2010年间采空区中心形变较为平稳，无明显突变，但2011—2013年该区域地表形变异常，通过与时序结果对比发现在2010年后地表发生明显变化的15、17、18、19、20号采空区均位于东部漏斗群，7 a间最大累计形变超过-2.3 m。

图7 浑源县采空区中心历年累计变化量

2007—2013年间，浑源县三大沉降漏斗群均出现不同程度形变。在时间纬度，49个采空区均呈持续的、大量级的形变，最大累计沉降超过-2.3 m；在空间纬度，采空区分布于11个乡镇，其中以东部的黄花滩乡、大磁窑镇、青瓷窑乡、千佛岭乡、大仁庄乡最为严重。7 a间浑源县地面沉降无论从沉降量级还是沉降范围均发生了巨大的变化。

3.3 讨论

为分析时序InSAR监测结果可靠性，本文搜集了浑源县采煤沉陷区村庄治理、大同市原煤产量与浑源县煤炭销量[24]等资料。首先对涉及采煤沉陷区治理的村庄分析，发现村庄主要分布在大仁庄乡、大磁窑镇、东坊城乡、黄花滩乡、青瓷窑乡、千佛岭乡、蔡村镇、沙圪坨镇等乡镇，分布于浑源县东部、北部、西北部，其空间分布上与根据时序InSAR结果圈定的49个采空区空间位置分布相一致。说明煤矿开采是浑源县出现形变漏斗的主要原因。

此外，通过对大同市原煤产量与浑源县煤炭销量进行分析发现，自2008年起，当地煤矿开采呈现爆炸式增长。大同市原煤销量从2007年的7 394万t增长至2010年的9 403万t，2015年达122 59万t；浑源县煤炭销量从2008年的15.46万t，增长至2011年1 060万t，年均增长幅度达309%。当地煤矿的产量与销量变化情况与部分采空区形变量在2010年前后发生突变趋势相吻合（图8）。经调查发现，自2008年起，煤炭价格飞速上涨，煤矿开采日益严重，尤其开采难度较低的露天煤矿已初具规模，使得多数区域年形变量开始呈指数加速增长状态。

图8 大同市原煤产量与浑源县煤炭销量及沉降中心年均形变时序图

4 结论

本文分析了2007年6月—2011年1月36景ALOS PALSAR影像和2012年1月—2013年4月10景RADARSAT-2影像数据，获取了浑源地区地表形变速率结果和时序形变结果。通过对结果分析可知：

1）根据时序InSAR地表形变结果，结合99座小煤矿范围、5大煤矿位置和野外地面调查资料，圈定49个采空区范围及各采空区历年累计区域平均形变量统计结果。7 a间浑源县地表呈持续亏损状态，累计形变量不断加剧，且在部分区域2010年后形变呈突然加剧状态。

2）2007—2013年间，浑源县出现严重的地表形变，在东部、北部、西北部分别形成了三大沉降漏斗群，其中范围最广的东部区域分布着果子园、古瓷窑、白庄等众多煤矿。随着时间推移，浑源县49个区域中历年形变量大的比重不断增加，沉降量级与沉降范围均呈恶化趋势。

3）根据时序InSAR结果圈定的形变区空间分布与采煤沉陷区开采量一致性较高，表明煤矿开采是造成地表沉陷加剧的主要原因。