基于时序InSAR的新疆阿希矿区地表形变监测与分析
虎小强1, 杨树文1,2,3, 闫恒1, 薛庆1, 张乃心1
1.兰州交通大学测绘与地理信息学院,兰州 730070
2.地理国情监测技术应用国家地方联合工程研究中心,兰州 730070
3.甘肃省地理国情监测工程实验室,兰州 730070
摘要:
新疆阿希矿区地理环境复杂,长期的矿产资源开采致使矿区的地面发生了较为严重的沉降和变形,造成了采矿生产安全隐患和周边生态环境的破坏。为了进一步研究分析阿希矿区地面沉降时空变化特征和地表形变规律。首先,利用小基线集合成孔径雷达干涉测量(small baseline subset-interferometric synthetic aperture Radar,SBAS-InSAR)技术对2017年2月9日—2021年4月25日获取的127景升轨Sentinel-1A影像数据进行沉降计算; 其次,将InSAR沉降监测结果与水准测量结果进行对比验证; 最后,分析了阿希矿区近5 a地面沉降时空变化特征,并研究了影响矿区地面沉降的驱动因素。结果表明,在整个监测期间内,阿希矿区地表形变呈现出整体基本趋于稳定状态、局部沉降较为明显的特征。造成矿区沉降的主要因素为矿产开采、地质构造、降水及露天矿井蓄水等。本研究可为矿区地面沉降监测以及今后地下矿物的合理开采提供科学依据。
0 引言
地面沉降是在自然因素和人为因素共同影响下,导致一定范围内地面高程降低的地质现象,是一种无法自我恢复的环境和资源损失[1⇓-3]。新疆阿希地区的高强度矿产开发引发了不同程度的地表形变,破坏了当地生态环境,造成潜在的采矿安全隐患。由此,探索和研究阿希矿区的地面形变特征以及时空变化规律具有非常重要的现实意义。
近年来,应用于地表形变监测的方法主要有基于点的水准测量和全球定位系统等[4],然而,这些传统监测手段通常耗时耗力且空间分辨率较低[3-4]。基于合成孔径雷达数据的干涉测量(interferometric synthetic aperture Radar,InSAR)技术发展迅速,与传统方法相比,该方法具有效率高、精度高等优点,已经广泛应用于区域性地表形变监测、滑坡灾害预测及地震形变[3]。差分雷达干涉测量技术(differential interferometric synthetic aperture Radar,D-InSAR)技术可获得精度达毫米级的地表形变信息,但其精度容易受时空失相干和大气延迟等因素的影响,且时间分辨率也比较低[4⇓⇓-7]。时间序列InSAR技术有效地克服了时空相干、大气延迟异常以及地形对差分的影响等问题[8⇓⇓⇓⇓⇓⇓-15]。由此,许强等[3]利用时序InSAR技术对延安新区地面沉降时空演化进行了研究与分析; Dong等[7]使用时间序列InSAR技术结合开采沉陷对称特征对煤矿区进行三维位移反演; Yang等[8]运用多时相InSAR技术对东胜矿区地表形变进行了研究与分析。总结上述研究并结合其他文献可知,时序InSAR技术在矿区地表形变监测中具有广泛的应用前景。因此,本文利用时序InSAR技术对新疆阿希矿区地面沉降进行时空变化监测。
地表形变是新疆阿希矿区严重的次生灾害,对采矿生产造成严重安全隐患,致使周边生态环境也发生破坏。但针对新疆阿希矿区地表形变相关研究较少,未能有效、科学地反映矿区地表形变规律。因此,为了更好地研究阿希矿区的形变特征以及时空变化规律,本文选取Sentinel-1A升轨雷达影像数据,采用SBAS(small baseline subset)-InSAR技术获得2017—2021年5 a间新疆阿希矿区的形变信息,分析了地面沉降与开山取金的时空相关性,并对其沉降驱动因素进行分析与探讨,进而为今后阿希金矿有效开采提供科学依据。
1 研究区概况及数据源
1.1 研究区概况
阿希矿区位于新疆维吾尔自治区伊宁县境内,地处西天山西段的吐拉苏火山断陷盆地内[16],其范围大致分布在E81°36' ~ 81°38',N44°13'~44°13',是大型的浅层低温黄金矿山[17]。矿区地质构造发育,以断裂为主,区域地层及其构造线均呈西北方向展布,地形切割强烈[17]。研究区概况如图1所示。
1.2 数据源
本研究应用的数据源有Sentienl-1A影像数据、SRTM DEM(digital elevation model)数据、降水数据、地质数据和Google Earth影像等。其中,Sentienl-1A升轨影像数据用于对矿区地面沉降特征进行研究,其空间分辨率为5 m×20 m,重访周期为12 d,入射角为39°,极化方式为VV[3]。为了更好地研究与分析地面沉降的时空特性以及沉降机理,选取了2017年2月9日—2021年4月25日的127景升轨Sentinel-1A单视复数(single look complex,SLC)影像。并采用空间分辨率为30 m的SRTM DEM数据去除地形相位分量。此外,采用2017年和2020年2期的Google Earth影像对矿区工程建设的空间变化进行解译,用于沉降原因分析。
2 研究方法
2.1 SBAS-InSAR沉降监测原理
SBAS-InSAR是一种利用多幅主影像与基于短时空基线准则的InSAR时间序列分析方法,该方法在一定程度上可以有效地克服传统InSAR技术空间基线过长,导致的时间和空间失相干现象,并能有效地削弱大气延迟的影响,同时也增加了形变监测的时间分辨率[18]。
沉降监测原理为[19⇓⇓⇓⇓⇓-25]: 利用覆盖同一区域的N+1幅SAR影像,选一幅影像作为超级主影像进行配准,通过设置恰当的时空基线阈值组成不同的干涉对,生成差分干涉图。
2.2 数据处理及形变计算
利用ENVI/SARscape软件对所获取的127景Sentinel-1A数据进行时间序列处理和形变计算,具体处理流程如图2所示。主要包括以下5个步骤: 1)生成连接图和干涉对组合。根据阿希矿区的地形和气候等因素分析,将空间基线阈值设为60%,以防止完全空间失相干; 同时,为避免无相干数据对,将时间基线阈值设为300,便于估算低相干性区域,建立新的连接; 然后,选择超级主影像并进行处理,生成时空基线(图3)。
图3 SBAS-InSAR时空基线
2)干涉工作流处理。首先利用复数像对共轭相乘得到干涉图,再利用参考DEM模拟并剔除其地形相位,采用Goldstein算法平滑相干斑噪声,提高干涉条纹清晰度。选择最小费用流(minimum cost flow,MCF)算法获得差分干涉序列的解缠相位图,并剔除不理想的干涉像对。
3)SBAS-InSAR轨道精炼和重去平。选取相位稳定且远离形变区的35个地面控制点(ground control points,GCP),采用多项式拟合优化方法。达到消除斜坡效应,修正相位偏移,估算以及去除残余恒定相位的目的。
4)SBAS-InSAR反演。先通过建立线性模型计算出所有像对的残余地形和地表形变速率,利用SVD法解算相关参数,得出干涉像对的平均形变速率和高差。然后再通过大气滤波估算、去除大气相位和地形残余相位,最终获得时间序列位移图。
5)地理编码。将基于SBAS坐标系的反演结果转化为地理坐标,并显示地理坐标系下视线向形变和垂直向形变。
3 地面沉降时空特征及成因分析
3.1 地面沉降空间分布特征
利用SBAS-InSAR技术对2017年2月9日—2021年4月25日期间获取的覆盖新疆阿希矿区的Sentinel-1A升轨影像数据进行数据处理,得到了新疆阿希矿区长时间序列的地表形变信息(图4),其中,负速率值和正速率值分别表示地面沉降和抬升。通过对图4分析发现,在监测时间内,阿希矿区地面平均形变速率范围主要分布在-15~-5 mm/a之间,最大平均形变速率可达-21 mm/a; 阿希矿区地面沉降变形主要分布在重点监测点N1~N5(图4中小圆点)附近。各监测点的平均形变速率如表1所示,平均形变速率最大区域分布在N3点附近,平均形变速率最小区域分布N5点附近。监测区边缘地带在监测时间内基本处于稳定状态,平均形变速率在-2~1 mm/a之间。
表1 2017—2021年监测点平均形变速率 (mm·a-1)
3.2 InSAR监测精度验证
为验证SBAS-InSAR技术监测结果的可靠性与一致性,选用了研究区内8个地面水准测量基准点(图4中小五星),以每个基准点为中心特定半径(选取半径为10 m)内像元的平均形变速率作为对应InSAR测量结果。
研究中将2019年3月和2020年3月实测水准数据与Sentinel-1A数据获取的沉降速率进行对比分析,总体上呈现了良好的一致性,其对比结果如表2所示。对比发现,二者之间存在一定的差异。为此,计算了每个水准点的InSAR和水准测量之间差异的均方根误差RMSE,并计算了2个测量结果之间的最大和最小误差。RMSE为0.75 mm/a,最大和最小误差分别为3.3 mm/a和1.0 mm/a。
表2 SBAS-InSAR监测结果精度验证(mm·a-1)
进一步分析表明,二者之间的差异主要归因于水准测量和InSAR观测之间的位置和采集时间点并非完全吻合,及水准测量时外界条件等因素对其精度均有一定的影响,这种误差对InSAR监测结果的使用影响相对较小。因此,可以利用InSAR监测结果来分析矿区地面沉降时空演化特征及地面沉降与各驱动因素之间的关系。
3.3 地面沉降时空特征分析
3.3.1 地面沉降空间演化特征分析
以2017年2月9日为起始日期,2021年4月25日为截止日期,以6个月为时间段,选取了研究区的9幅Sentinel-1A升轨影像数据来研究矿区地面沉降的空间演化过程(图5)。通过图5分析发现,2017年2月—2018年2月,矿区地面沉降分布面积明显增加,累积形变量范围在-25~-8 mm的区域增加尤为显著; 2018年8月—2019年8月,矿区地面沉降分布总面积持续增加,且局部区域(红色图斑)沉降面积变化更加迅速,累积沉降量增加明显,其中重点监测点N1~N4面积变化最为明显; 2020年2月—2021年2月,研究区内累积形变量为-15~-8 mm(黄色)之间的区域有所减小; 在局部区域,累积形变量为-20~-9 mm(红色)之间的区域沉降面积分布扩大,且沉降速率变化更加明显。
从整个空间演化过程可知,2017—2021年研究区内地面沉降分布总面积呈现出先增加后减少现象; 但在局部沉降分布面积持续增加,沉降速率也大幅上升,最大累积形变量至-90 mm左右。据此表明,沉降严重的每个区域都位于中度沉降区域内,说明严重沉降只是一种局部现象,人为采矿是致使局部地面严重沉降的主要因素。
3.3.2 地面沉降时间演化特征
为了更好地研究新疆阿希矿区地面沉降时间演化特征,将对图4中的5个重点监测点进行时序形变信息分析,得到了2017—2021年间矿区累积形变信息随时间变化折线图,如图6所示。分析可知,在监测时间内,各监测点均以不同的沉降速率呈现出逐年增加的趋势,其中,N3在监测时段内形变速率变化最大,沉降最为严重,累积形变量约为-87 mm; N5在监测时段内沉降变化速率最小,其累积形变量为约为-37 mm,此区域沉降速率缓慢,沉降范围基本趋于稳定。
此外,利用线性拟合模型对各个重点监测点的时序形变信息进行拟合,结果如图7所示。由此分析可知,在2017—2021年间,地面累积形变量随时间的增加而增加,且基本呈线性关系。因此推测,在未来较长的时间段内,地表形变尚未达到稳定状态,沉降将依然继续发生,此后可能会继续加剧。
3.4 地面沉降驱动因素分析
近年来,随着阿希矿区矿产资源储量进一步探明以及黄金价格不断攀升,阿希矿区生产规模逐年扩大,致使矿区地表发生了严重的形变。此外,地质构造、降水补给以及露天矿井蓄水等因素都不同程度影响着矿区的地表形变,因此仍需进一步研究矿区地面沉降与各种驱动因素之间的关系。
3.4.1 采矿工程
在监测时间内,选取覆盖矿区的2期历史 (2017年7月15日和2020年7月3日)Google Earth光学影像对矿区工程建设过程中的空间变化进行解译,其可以基本反映出矿区空间演化过程,并将其与2017—2021年累积形变量图进行叠加(图8)。由图8(a)分析可知,从2017—2020年矿区工程区域面积明显增加,各重点监测区域也伴随着不同程度的增加。根据图8(b)发现,该地区形变最大的区域在空间分布上与重点监测区域吻合,且随着采矿工程的进行,矿区沉降面积不断扩大,沉降速率持续加快。因此,在2017—2021年,采矿工程可能是影响并加剧新疆阿希矿区地面沉降速率的主要因素。
3.4.2 地质构造
相关研究表明,断层构造、地层岩性特征对地面沉降的产生和发展有一定的影响[26]。因此,本研究将阿希矿区地质图与年平均形变速率图进行叠置分析(图9)。阿希金矿床南北向展布的由3个矿化带构成,构造以断裂为主,褶皱不发育[27⇓⇓-30]。矿带总长为1 280 m,西缘矿化带是金矿床最主要的成矿带,是矿山主要开发对象。主矿体受控于弧形断裂F2,总体呈向南西凸出的弧形带状分布; 矿体东南部被F4断裂带强烈切割,发生明显断层现象。矿区发育的断裂构造和复杂岩性组合。在成矿前期F2断裂为长期活动的张性断裂,火山活动末期岩浆房收缩导致火山口下陷,主要是来自火山活动晚期的次火山岩浆和期后的火山气液依次沿该断裂裂隙上侵瀑裂围岩、温度转高的热液使围岩发生硅化、黄铁绢英岩化等蚀变作用,然后在区域构造应力作用下,该断裂再度张开,深部含矿热液再次沿蚀变围岩裂隙贯入最终结束矿化作用。随着矿物资源的开发,受F2断层控制,上下盘未发现明显形变。F2断裂既是导矿构造,又是控矿构造。此外,由于长期的大气降水渗透,成矿物质的火山岩对基底火山岩进行改造,火山下部受热上升。上述过程持续进行形成了地下热液的对流循环,上升的溶液侵蚀了导矿构造岩,使得断裂节理裂隙向外移动,进而使得矿体中部区域的地表形变加剧。由此,通过上述分析表明,地质构造是新疆阿希矿区地面沉降的控制性因素。
3.4.3 降水
据相关文献分析表明区域降水是诱发沉降的重要因素之一[26]。因此,为了分析降水量与地面形变之间的关系,分别选取伊宁县2019年的降水数据与2019年阿希矿区4个季节的平均形变速率数据进行研究。2019年的地面形变速率可通过年末和年初相邻日期的地面形变累积结果进行相减得到,每个季节的形变量除以时间得到其平均形变速率。通过研究分析2019年降水可知,伊宁县降水主要集中现在夏季(6—8月),春季(3—5月)、秋季(9—11月)、冬季(12—2月)降水量相对较少; 将矿区夏季(6—8月)的地面形变速率与春季、秋季、冬季进行对比。由表3可知,冬季的地面形变速率最大,夏季和秋季次之,春季的平均地面形变速率最小。
表3 2019年各季节形变速率(mm·a-1)
此外,为了进一步分析地面沉降与月降水量的关系,本文针对5个重点监测点的地面累积形变量与月降水量进行比较分析,结果如图10所示。在降水量较多的范围内,重点监测点的形变量增加缓慢,反之在降水量减小的时间范围内,重点监测点的形变量增加显著,说明降水在一定程度上减缓了矿区地面沉降速率。
由上述地面沉降与降水的关系可知,降水在一定程度上减缓了矿区地面沉降速率。但在某些特殊区域(露天矿坑与其周围地区),降水却对地面沉降速率有着不同的影响,如图11所示。一是采矿过程中产生的碎石、沙砾、泥土形成的沙土推,在多雨季节,当有大量的雨水注入时,软化致使其地面下沉,造成地面沉降。另外,露天矿井蓄水也是引起矿区地表形变的原因之一。如图11(a)和(c)可知,2017—2020年工程区域面积、体积基本没有发生明显变化。但由图11(b)和(d)可知,在监测时间范围内工程区域发生了明显的地面沉降现象。由此,造成上述现象的过程可能是大量积水浸泡使矿井底层下陷,周围岩壁滑落,致使矿井周围地面发生塌陷,导致地面沉降。说明在某些特殊区域降水加速了矿区地面沉降速率。
4 结论
采用SABS-InSAR技术对2017—2021年覆盖研究区的127景Sentinel-1A升轨影像数据进行地面沉降计算,获得了新疆阿希矿区地表平均形变速率分布,分析了矿区地面沉降的时空演化特性以及沉降驱动因素,得到的主要结论如下:
1)将InSAR监测结果与水准测量结果进行交叉验证,其结果基本一致,表明监测结果可靠性。
2)在观测期间,整个研究区基本处于稳定状态,但在局部区域沉降显著,其年均最大形变速率大约可达-21 mm/a。
3)由时序形变信息拟合结果可知,各监测点沉降量随时间的增加而增加,且基本呈线性关系。可以预测,研究区局部地面沉降将在未来一段时间内将持续增加,沉降漏斗将不断变大。
4)通过对阿希矿区地面沉降驱动因素分析可知,采矿工程是引起矿区地面沉降的主要原因。此外,地质构造是阿希矿区地面沉降的控制性因素; 降水在一定程度上减缓了矿区地面沉降速率,但在某些特殊区降水却加速了矿区地面沉降速率。
（原文有删减）
【作者简介】虎小强(1995-),男,硕士研究生,研究方向为InSAR形变监测。Email: 1938204120@qq.com。
【通讯作者】杨树文(1975-),男,博士,教授,研究方向为遥感数字图像处理及信息自动提取、灾害遥感。Email: 825198827@qq.com。
【】国家自然科学基金项目“西北重点城市彩钢板建筑群与产业园区时空关联关系”(42161069);国家自然科学基金项目“基于高分辨率卫星影像的彩钢板建筑与城市空间结构演变关系研究”(41761082)
【】虎小强, 杨树文, 闫恒, 薛庆, 张乃心. 基于时序InSAR的新疆阿希矿区地表形变监测与分析[J]. 自然资源遥感, 2023, 35(1): 171-179.
转自：“测绘学术资讯”微信公众号
如有侵权，请联系本站删除！
