D-InSAR技术在矿区地表形变监测中的应用研究

马道鸣，鲁殿刚，朱文慧，张雪丽，苏永江

(河南省航空物探遥感中心，河南 郑州 453000)

0 前 言

矿区地表形变监测一直是矿山安全生产的基础保障，也是严重制约矿山升级发展的重要因素，是在人类采矿活动和自然地质条件变化的双重作用下的岩石应力变化，进而使得矿区深部采空区上部区域及邻近区域的岩矿石发生不均一的形变，甚至塌陷等，最终威胁着资源的安全开发活动，同时造成大量的矿石损耗等[1]。因此，提高矿区地表形变监测效果不仅是有效防治矿山灾害事故的主要途径，也是提高资源利用效率的基础。传统的GPS、全站仪等监测仅可获得监测点在某个时间点的沉降量，而无法获得监测区域整体的形变特征，这就限制了对监测区域整体形变规律的研究，进而造成形变特征分析出现偏差的局面，所以选择监测方法是提高监测质量的前提[2]。合成孔径雷达差分干涉测量技术(D-InSAR)是在合成孔径雷达干涉测量技术(InSAR)上发展起来的新型监测技术，该监测方法充分利用了高精度DEM成果及干涉图等基础土建，能够更精确地反映矿区地表微小的形变量，故常用于较精细及缓慢形变监测领域中，如滑坡形变监测[3-4]、冻土区公路地表形变监测[5-6]、城市地表形变监测[7-8]、铁路地表形变监测[9]、火山附近地表形变监测[10]、地震监测领域[11-13]以及矿区地表形变监测等领域，在上述形变监测中取得了良好的应用效果。鉴于此，本文以金铜矿山的多幅ENVISAT ASAR数据为基础，试图采用D-InSAR技术分析该矿区地表形变规律，综合研究矿区地表形变发展规律，为矿山资源的安全开发利用提供参考依据。

1 D-InSAR的技术流程

矿区地表形变是人类采矿活动及地质条件变化的双重作用下诱发引起的地壳浅表层高程变化以及位移的一种表现形式，D-InSAR技术就是利用了不同时段高分辨率影像数据中的微小形变特征而获取地表形变信息的。因此，D-InSAR技术在矿区地表形变监测中的应用就是利用不同时段获得的ENVISAT ASAR影像数据，通过影像数据的进一步解缠处理等，就可获得该时段内矿区地表的微小形变特征。影像数据的处理方法有二轨法、三轨法和四轨法等几种类型，其中后2种数据处理方法中引入了大气误差，容易造成监测精度降低。二轨法是最为常见的和常用的数据处理方法，本文也采用二轨法进行处理跨越形变期不同时段的ENVISAT ASAR影像数据，进而借助高精度DEM模型进行反演，将地形相位等剔除，就可获得监测区域内任何一点的形变量数据。综上所述，本文采取的二轨法在数据处理方面的技术流程见图1。

图1 D-InSAR技术的二轨法数据处理流程

2 矿区基本概况

研究区为一金铜多金属矿床，矿区面积为2.5 km2，矿区边界为一不规则的多边形。矿区地形高差变化较大，属于高山峡谷地貌与丘陵山区地貌的结合部位，导致地形地貌较为复杂。矿床成因类型为矽卡岩型金铜多金属矿床，矿山已有的开采方式为井下硐采，且局域北西区域矿体埋深较浅，南东区域矿体埋深较深。该矿山为一开采矿山，现状已有6处采空区，采空区上方区域的地表已出现不同程度的地表沉降问题，甚至在某些采空区上方已出现细微的地裂缝。因此，为了了解该矿区地表形变规律以及指导资源开发，本文以ENVISAT ASAR影像数据为基础，通过D-InSAR技术研究矿区地表形变发展规律，为矿山安全生产提供参考依据。

3 D-InSAR技术在矿区地表形变监测中的应用

3.1 配准影像数据及估算基线

前人的大量研究资料[3-9]显示，二轨法的影像数据配准方法主要包括2种类型：一是以精密星历数据以及相干系数配准影像数据；二是以地面控制点为基准点进行影像数据的配准。本次所使用的影像数据来源于COSMO-SkyMed高分辨率雷达卫星，相对而言该卫星的影像数据精度较低，故本次影像数据的配准方式选择以精密星历数据以及相干系数配准，可以通过该种配准方式有效地提高监测精度。在使用第一种配准方式过程中一般采用DORIS轨道数据对相应的ENVISAT ASAR影像数据进行基线的估算处理，在完成上述操作的基础上对矿区内地形地貌相对较为平缓的区域以条纹率进行精估计[2]，就可获得精度较高的配准后的影像数据。影像数据精配准后结合矿区高精度的DEM数据对影像数据进行相位模型处理，将影像数据中包含的地形信息、形变信息等转换为相位信息，为后期进一步展开影像数据的相位差分干涉处理以及去地形相位解缠等处理做准备。

3.2 差分处理及相位滤波

在完成ENVISAT ASAR影像数据配准和基线估算后形成了包含地形信息和形变信息的相位信息，即相位干涉图，此时的相位干涉图中包含了地形信息和形变信息。因此，需要将相位干涉图中的地形信息剔除，就可获得去地形干涉相位，该过程就是干涉影像数据的差分处理[13]。去地形干涉相位中不仅包含了大量的微小形变信息，而且包含了其他的噪声信息，如噪声相位、非线性相位等，此时需要进行噪声相位和非线性相位的剔除处理，即差分干涉相位的滤波处理。

3.3 相位解缠及地理编码

完成差分干涉相位滤波处理后可进行相位的解缠过程，该过程也是获取矿区地表形变信息的基础。本文相位的解缠采用最小费用流法，该方法能够降低影像数据中微小形变模糊的问题[5]，进而获得与线性变化地形相互一一对应的相位数据，在GAMMA软件中数据配准换算关系等解算出不同解缠点的形变量信息[12]。采集形变信息的数据并编辑地理编码等信息，并将其统一至同一坐标系统下，进而通过图件整饰处理等就可获得矿区地表形变图等成果性图件，能够将矿区地表形变直观地表达出来(见图2)。

图2 研究区地表沉降形变分布

4 矿区地表沉降形变分析

4.1 沉降形变监测结果分析

在完成相位解缠后就可在MAPGIS软件等中采集并编辑形变信息以及地理编码等，再通过图件整饰就可生成矿区地表形变分布图(见图2)，图件中包含了矿山中的形变区域、分布位置、形变范围、沉降量等信息。由图2可知：研究区地表沉降范围共识别出6个区域，地表形变范围与深部采空区具有良好的空间吻合关系；同时，矿区地表形变面积以及沉降量大小与深部资源的开发利用关系密切，其中采空区形成时间长且面积较大，地表形变面积较大且沉降量较大，如图2南东侧和北西侧的采空区所对应的地表形变区，即在垂向上地表形变分布区域深部采空区具有良好的对应关系；地表形变分布区的空间展布形态与深部采空区的展布方向和形态基本一致，形变分布区的长轴方向与采空区的长轴方向一致。由矿山沉降量形变图可知：矿区内地表形变可划分为3个等级(-25～0 mm、-50～-25 mm和-100～-50 mm)，其中矿区地表最大沉降量位于矿区南东侧的采空区上部，即矿区主矿体采空区上部，最大沉降量可达92.5 mm；矿区内最小沉降形变区分布在矿区南西侧，其沉降量一般小于25 mm；矿区内总形变区域可达0.78 km2，其中形变量大于-50 mm的面积约为0.07 km2，形变量介于-50～-25 mm的面积为0.31 km2，形变量介于-25～0 mm的面积为0.40 km2。

为了解矿区地表形变在垂向走向的变化规律，在形变最发育的范围内对垂直形变主体走向进行了剖面研究(见图3)。由图3可知：研究区矿区地表形变在剖面上具有大致为“V”形变化的规律，其最大沉降量位于采空区中心部位靠近南东一侧，与采空区的分布范围基本吻合；同时，沉降量从采空区中心部位向两侧具有大致对称分布的特征，且形变量向远离采空区范围逐渐降低。

图3 研究区A-B剖面沉降特征

4.2 精度对比分析

为了对比分析D-InSAR技术在矿区地表形变监测中的精度问题，本文利用矿区形变区内已有的控制点坐标信息等，使用传统的全站仪测量方法对D-InSAR技术的精度进行验证。本文累计收集了10个控制点的坐标信息，通过开采前后的三维坐标对比，获得了控制点的形变量数据，进而与D-InSAR形成的沉降形变图中的沉降量对比分析。本次为了获得更加准确的D-InSAR沉降形变数据，将控制点周边的约30个高相干点数据进行加权平均，进而作为该控制点的D-InSAR方法的数据，最终的数据统计见表1。由表1可知：本文所获取的D-InSAR加权平均形变量数据与全站仪测量所获形变量的变化趋势是一致的，二者的误差百分比均小于5%，说明本文所获矿区地表形变量是可靠的，其精度能够满足矿山形变监测应用。因此，可以将合成孔径雷达差分干涉测量技术(D-InSAR)应用至矿区地表形变监测应用中，其监测效果良好。

表1 D-InSAR技术监测结果与全站仪监测结果对比

5 结 论

①D-InSAR技术在矿区地表形变监测中具有良好的应用优势，该技术方法能够获取矿区微小的形变信息，对全面了解矿区地表形变特征以及综合研究形变规律，进而预测形变发展趋势等意义重大。

②研究区地表形变范围与矿区深部采空区具有良好的空间吻合关系，形变长轴与矿区采空区的长轴方向基本一致，沉降量从采空区中心向两侧呈大致对称状分布，且向远离采空区中心向两侧形变量逐渐减小，总体上具有“V”形的变化特征。

③矿区地表最大沉降量可达92.5 mm，最小沉降形变区分布在矿区南西侧，其沉降量一般小于25 mm；矿区内总形变区域可达0.78 km2，其中形变量大于-50 mm的面积约为0.07 km2，形变量介于-50～-25 mm的面积为0.31 km2，形变量介于-25～0 mm的面积为0.40 km2。