0 引 言
随着煤炭工业的发展，煤炭开采所造成的地表沉陷问题越发严峻，破坏耕地，损毁公共基础设施，严重时还将威胁当地居民的生活安全[1-3]。 我国华北平原矿区大部分属于高潜水位矿区，在开采过程中极易形成大面积的水体，沉陷较深，威胁性更大[4-5]。 因此快速高效识别沉陷范围，对矿区土地复垦工作具有重要意义。 目前针对高潜水位矿区沉陷信息提取研究并不多。 早期的矿区沉陷水体信息提取主要是实地测量和依靠遥感影像目视解译，费时费力[6-7]。 后期随着影像分辨率的提高，各种分类方法相应提出，各种水体指数也不断改进。 徐涵秋[8]提出的改进归一化水体指数MNDWI(Modified Normalized Difference Water Index)在水体提取方面具有较好的优势，可以很容易地区分阴影和水体，解决了水体提取中难于消除阴影的难题。 李晶等[9]针对整个济宁市的情况，基于Landsat TM 影像提出了一种集成MNDWI 法和HIS 空间水体模型优点的综合性水体识别方法，在识别地表沉陷形成地表不规则的零星积水方面具有较好的优势。 肖武等[10]通过监督分类验证高分一号数据对矿区水体监测精度较高。 针对非积水区域的地表形变研究，近几年迅速发展起来的合成孔径雷达干涉测量技术(InSAR)是目前监测地表形变最为高效且精度可达毫米级的技术，具有数据处理流程化的优势[11-12]。 其中SBAS-InSAR 技术能够有效解决时空失相关和大气延迟影响的问题，使监测地表形变精度得以大幅提高，在探测地表长期累积的缓慢形变方面表现出了极大的潜力[13]。 胡乐银等[14]证明 SBAS-InSAR 技术获取到的地表形变时间序列在空间上相对更加连续。 周吕等[15] 采用 Sentinel-1 数据，在 SBASInSAR 技术的支持下得到武汉地表形变速率精度为6 mm/a。 董少春等[16]利用 SBAS-InSAR 技术获取常州监测时期内年均沉降速率和地表累积形变图，揭示了监测期内常州市各地区的地面沉降时空分布和变化特征。
目前针对高潜水位矿区的地表形变监测研究通常仅获取沉陷水体信息或者沉陷非积水区信息，并没有将两者相结合。 虽然高分辨率影像提取水体结果精度较高，但是其覆盖面较小且价格较高，适用性较低。 而Landsat 系列数据具有时间和空间上覆盖全面、数据免费的优势。 Sentinel-1A 数据在监测地表变形方面具有很好的优势，但在采煤沉陷区域的应用较少。 因此，这2 种数据和方法组合可以获得更全面的矿区沉陷信息。 笔者以济宁市鲍店煤矿为例，针对高潜水位矿区地表沉陷特点，结合沉陷积水区和非积水区的提取方法，利用多源多时相数据综合分析，获取整个矿区地表沉陷信息，以期为矿区土地复垦提供参考依据。
1 矿区概况与数据选取
1.1 矿区概况
鲍店煤矿是山东省济宁市邹城西部矿区的重要组成部分(图1)，地处兖州煤田西南端，地理坐标范围为北纬 35°23′13.2″～35°28′8.4″，东经 116°48′7.2″～116°52′26.4″，面积约 35.71 km2[17]。 该矿区地势平坦，东北高西南低，地面标高为+39.97—+44.75 m。区域地貌类型属冲积平原，大部分土地为农耕区。地下水位埋深较浅，在2～5 m，属于典型的高潜水位煤矿区。 鲍店矿区是我国自行设计施工的大型现代化矿井，于 1977 年10 月 14 日动工兴建，1986 年 6月10 日建成投产。 经过30 多年的开采，造成了矿区大范围的土地沉陷，由于潜水位较高，又汇集形成了大面积水体。
1.2 数据源及预处理
利用Landsat 系列数据提取水体信息(表1)。行列号是获取数据时的一个重要参数，可以根据行列号直接定位到研究区的位置。 为区分矿区自然水体，选取该矿区开采前后的1973 年、1979 年和1984 年的数据进行处理。 矿区沉陷水体主要选取2015 年至2017 年5 景影像质量良好，云量稀少的Landsat 8 OLI 影像进行提取。 Landsat 8 卫星是目前唯一一颗在轨运行的Landsat 系列卫星，性能各方面相对于之前都有了很大提高，并且在长期观测陆地表面状态及其变化监测方面有较好的优势。 利用ENVI 5.3 进行数据的辐射校正，大气校正等处理。
图1 鲍店矿区
Fig.1 Baodian Mining Area
表1 Landsat 系列数据主要信息
Table 1 Main information of Landsat series data
年份 行列号 获取时间 分辨率/m 卫星1973 1979 1984 122/35 122/35 122/35 1973-11-19 1979-10-28 1984-04-19 30 30 30 Landsat1 Landsat2 Landsat5 2015 122/35 2015-11-03 30 Landsat8 2016 122/35 2016-03-10 30 Landsat8 2016 122/35 2016-10-04 30 Landsat8 2017 2017 122/35 122/35 2017-04-30 2017-11-08 30 30 Landsat8 Landsat8
Sentinel-1A 卫星是由2014 年欧空局发射的地球观测卫星，不受天气的影响，可以提供连续影像，在海洋环境监视，监测地面运动风险，森林制图，水和土壤管理和测绘等方面具有较好的优势。 利用Sentinel-1A 数据进行SBAS-InSAR 处理，获取研究区非积水地表形变信息。 选取2016 年1 月至2017年12 月共24 景哨兵数据，极化方式为VV，轨道号为升轨142(表2)。 利用高级雷达图像处理软件(SARScape)进行数据预处理，主要包括裁剪和获取研究区高程数据。
2 沉陷信息提取方法
2.1 沉陷积水区信息提取
2.1.1 改进归一化水体指数
徐涵秋[8]在对Mcfeeters 提出的归一化水体指数(NDWI)的基础上，对构成该指数的波长组合进行了修改，将近红外波段(NIR)改为中红外波段(MIR)，提出了改进归一化水体指数(MNDWI)，即
表2 哨兵1 号A 星数据主要信息
Table 2 Main information of Sentinel-1A data
序号 获取时间 轨道数 时间基线/d 1234567891 0 0 11 12 13 14 15 16 17 18 19 2016-01-14 2016-03-02 2016-03-26 2016-04-19 2016-05-13 2016-06-30 2016-07-24 2016-08-17 2016-10-04 2016-10-16 2016-11-09 2016-12-15 2017-01-08 2017-02-13 2017-03-21 2017-04-14 2017-05-20 2017-06-13 2017-07-19 142 142 142 142 142 142 142 142 142 142 142 142 142 142 142 142 142 142 142 48 72 96 120 168 192 216 264 276 300 336 360 396 432 456 492 516 552 20 2017-08-12 142 576 21 2017-09-17 142 612 22 2017-10-11 142 636 23 2017-11-16 142 672 24 2017-12-22 142 708
MNDWI ＝(Green - MIR)/(Green + MIR)
式中：Green 代表绿光波段。
2.1.2 增强型改进归一化水体指数
由于高潜水位煤矿区采煤造成的沉陷区域基本为耕地，且自然形成的未经整治的沉陷水体没有明显的坑沿，与周边耕地混成一片，植被覆盖较多，矿区因采煤导致土壤中含有很多粉煤灰，而在利用传统方式提取水体信息时往往将这些土坎和粉煤灰误识为水体，仅依靠绿光波段(Green)和中红外波段(MIR)所提取的水体易模糊水体边界，最终所获取的水体面积明显偏大。 因此需要找到更为适合高潜水位矿区的提取沉陷水体的方法。
经过多次试验和各种水体指数对比，在徐涵秋提出的MNDWI 指数的基础之上，针对高潜水位采煤沉陷区沉陷水体区域的提取，提出了增强型改进归一化差异水体指数(E-MNDWI)。 对比水体、土坎、粉煤灰和植被的地物光谱曲线发现，四者在可见光波段的影像值大致相当，而在红外波段范围内其影像值趋势一致，但差距明显，水体影像值最低，其他地物均高于水体(图2)。
图2 地物光谱曲线
Fig.2 Spectral curves of ground object
在OLI 影像中，短波红外(SWIR)对水体具有强烈的吸收[18-20]，而植被、湿润土壤、矿物等非水体地物在该波段的反射均高于水体，因此在识别植被覆盖、湿润土壤和矿物方面具有较好的优势。 因此在MNDWI 的基础上，增加OLI 数据的短波红外波段，公式如下：
2.2 沉陷非积水区信息提取
2.2.1 小基线集技术
目前针对地表形变应用最广泛的方法是InSAR技术。 本研究主要利用SBAS-InSAR 技术来监测矿区非积水区域的沉陷情况。 2003 年，Berardino等[21]提出 SBAS-InSAR 方法，后经 LANARI 等[22]不断完善。 SBAS-InSAR 是将所有参加计算的SAR数据按照空间基线和时间基线分成不同的短基线子集，再对各个子集里的影像处理，其基本思想是连接由长基线造成的相互独立的SAR 影像，形成短基线SAR 影像集合，目的是为了增加数据获取的采样率，从而可以在已有的SAR 影像数据集中形成若干小集合，每个小集合内SAR 影像间的基线较小，集合间SAR 影像的基线较大。 SBAS-InSAR技术对于跨度范围大、需求精度高以及时间序列上的地表形变监测有很好的效果，其监测精度可以达到毫米级。
2.2.2 SBAS-InSAR 处理过程中的主要参数
将24 景影像导入SARscape 软件，利用已有影像获取覆盖研究区的SRTM-3 Version4 DEM 影像数据，以辅助研究区地形信息的生成。 考虑到采煤沉陷区分布及影响范围的特殊性，将24 景影像空间基线设为临界基线的0～30%，时间基线设为0 ～365天，采用德洛尼处理方法进行3D 解缠，冗余度设置为高。 干涉处理过程中，选取2017 年1 月8 日为主影像，总共生成了207 组干涉图，3D 解缠的结果为60 个像对。 在 2016 年 12 月 15 日 ～ 2016 年 10 月16 日生成的干涉图和解缠结果图上选择120 个GCP 点进行轨道精炼和重去平，其余参数默认即可。
3 结果与分析
3.1 水体分类结果精度评价
为验证E-MNDWI 指数的提取精度，利用混淆矩阵对MNDWI 和E-MNDWI 分类结果进行精度评价。 在分辨率较高的高分一号数据上选择40 个感兴趣区作为验证样本区，得到MNDWI 和E-MNDWI的分类结果的总体精度分别为87.5%和92.5%，Kappa 系数分别为75.37%和85.07%(Kappa 系数是通过混淆矩阵所计算出来的，和前面所说的总体精度类型一致)。 对比观察图3 中红色方框标记的区域，2017 年7 月的高分影像中可以明显看到类似于土坎的地物，而利用MNDWI 提取水体时，将土坎混淆为水体；E-MNDWI 则将土坎和水体很好地区分出。 相比之下，E-MNDWI 分类结果精度更好，更能准确地区分出水体和非水体，消除其他地物的影响，针对高潜水位矿区沉陷水体提取具有较好的优势。
为了定量分析这2 种指数分类结果，将高分一号数据监督分类后的结果作为真实值，采用相对误差来反映数据间的可信度，分别选取2017 年7 月的数据作为研究对象，结果见表3。
表3 2 种水体指数相对于高分一号的相对误差
Table 3 Relative error of two water indexes relative to GF-1
相对误差/%数据类型 获取时间MNDWI E-MNDWI高分一号 2017-07-08 — —Landsat 8 2017-07-19 34.67 7.44
由表3 可知，E-MNDWI 提取的水体结果相对误差小，即其分类精度更高。 而 MNDWI 的误差高达34.67%，其原因有2 点：①数据选取的时间问题，高分数据要比Landsat 8 数据早11 d，在这11 d 内可能由于降水的增加导致水体面积的大幅度变化；②MNDWI 将水边植被、土坎以及粉煤灰误识为水体。 对于高潜水位矿区而言，能够精准识别矿区水体范围，对于后期矿区沉陷治理具有重要意义，而30 m 分辨率的Landsat 8 数据通过E-MNDWI 指数提取水体结果和16 m 分辨率的高分一号数据利用监督分类提取水体结果具有一致性，精度相当。
图3 MNDWI 和E-MNDWI 水体信息提取对比
Fig.3 Comparison of MNDWI and E-MNDWI water information extraction
3.2 水体分类结果分析
根据煤矿区动工前后的Landsat 数据，发现在开采前矿区的水体只有泗河河段，未发现其他水域。因此在后期提取矿区沉陷水体时直接去除河流面积即可。
利用E-MNDWI 分别提取 2015 年 11 月，2016年 3 月，2016 年 11 月，2017 年 4 月以及 2017 年 11月矿区水体(图4)，通过目视解译以及实地数据考察，将水体阈值设置为 0 ～1。 借助 ArcGis10.2 统计分析，获得矿区沉陷水体的面积(表4)。
表4 E-MNDWI 提取矿区沉陷水体面积及其占研究区面积比例
Table 4 E-MNDWI extraction of subsidence water area in mining area and its proportion in study area
时间 沉陷水体面积/km2 所占比例/%2015-11 3.81 10.67 2016-03 3.99 11.17 2016-10 4.05 11.34 2017-04 4.12 11.54 2017-11 4.60 12.88
图4 不同时间E-MNDWI 提取矿区水体结果
Fig.4 Results of water body extraction from mining area by E-MNDWI
统计分析显示，春季或秋季，2016 年较2015 年沉陷水体面积变化较小，为0.24 km2，占整个矿区面积的 0.67%；2017 年比 2016 年增长了 0.55 km2，变化幅度较大，为2016 年变化量的2 倍。 总体来看，矿区沉陷水体面积在监测期内还在继续扩大，但是速度缓慢，年均增长面积为0.4 km2。
3.3 沉陷非积水区提取结果分析
Sentinel-1A 影像通过SBAS-InSAR 技术处理，得到矿区地表累计沉陷量图和矿区地表年均形变速率图(图5)，其中负值代表地表沉陷，正值代表地表抬升。 根据该矿区设置的水准监测桩，选取10 个监测点验证SBAS-InSAR 监测结果(图1)，得到两者最大误差为2.30 mm/a，最小误差为0.33 mm/a，误差率在 0.04 ～0.49。 SBAS-InSAR 监测值与水准数据的具有较好的相关关系，决定系数R2为0.92。
图5 SBAS-InSAR 处理结果
Fig.5 SBAS-InSAR processing results
结果表明，该矿区在监测期内沉陷速率范围为-41.69～12.83 mm/a，沉陷速率较快且较为严重；沉陷累计沉陷范围为30.16 ～-85.16 mm，该矿区非积水地表沉陷面积为5.5 km2。 由图5 看出，目前地表发生沉陷的区域大部分位于矿区的东部和南部，且靠近中部沉陷水坑的地方沉陷明显，沉陷水坑有扩张的趋势，极有可能发展为积水区域；矿区中西部积水区域较多，非积水区域地表形变不明显，表明该区域基本达到稳沉状态。 经过实地考察发现，较少区域产生的地表抬升情况与当地村庄的施工建设和地下水位的升高有关。
3.4 煤矿区整体沉陷范围
地表下沉量10 mm 即可认定为沉陷地。 通过克里金法将SBAS-InSAR 在监测时期内的地表累计沉陷量和沉陷水体结果进行插值处理[23]，填补相邻区域的空值，获得非积水区域的监测时期内下沉10 mm 及以上的等值线(图6)。 图6 中等值线围绕水体分布，越靠近水体边缘等值线越密集，说明地表沉陷剧烈，沉陷水坑有扩张的趋势，这和当地地表形变实际情况相符合。 在已有的沉陷水体之间形成了很多小的沉陷漏斗，推测后期会形成沉陷积水区域。等值线具有连续性，通过插值处理将周边沉陷量相近的范围划为一体，因而其面积较SBAS-InSAR 和E-MNDWI 直接监测的沉陷面积大，从图6 可以看出矿区的大部分区域都处于沉陷之中。
4 矿区沉陷信息提取方法讨论
图6 煤矿区沉陷信息
Fig.6 Coal mining area subsidence information
虽然利用水体指数获取研究区的水体信息的方式简单有效，但是在对高潜水位矿区水体提取时发现，各种水体指数其分类效果并不佳，可以推测各种水体指数具有一定的环境适用性。 水体在MIR 和SWIR 均有强烈的吸收，而植被、土坎和粉煤灰等非水体地物在这几个波段范围内的反射则均高于水体。 因此本研究在提取水体时，考虑到积水区域边缘大量的植被影响，以及矿区开采产生大量的粉煤灰和煤矸石，经过多种指数和阈值范围变化的试验，最终 改 进 MNDWI 指 数， 增 加 SWIR， 提 出 了E-MNDWI指数。 结果发现，E-MNDWI 能够更为清晰获取矿区沉陷水体的边缘信息，有效消除了水体周边植被的影响，从而使分类结果更为准确。 但该指数是否适用于其他积水区域还有待进一步验证。
SBAS-InSAR 监测结果可以达到毫米级，精度高，在提取地表形变时具有较好的优势，特别是针对沉陷较浅，肉眼难以辨别的区域。 能够减少人工成本，减少资金投入，提高了工作效率。 但是SBASInSAR 在监测矿区地表形变时，监测结果图斑破碎，推测可能是以下3 个方面的影响：①SBAS-InSAR受到相干性的影响，有些区域会出现空值，导致了结果图斑较为破碎；②推测研究时间较短，数据密度不够，导致有些地区发生微小形变但未监测到；③研究区有些地区可能因为地质原因并没有受到煤炭开采的影响，地表没有发生形变。 对此，图斑破碎问题还有待进一步的研究。 此外，哨兵数据量大，SBASInSAR 处理时间长，对计算机硬件设备的要求更高，如何提高处理速度也需要继续探究。
经过克里金插值弥补地表形变存在空值，因此获取的数值和真实值之间存在一定的误差。 等值线的获取是在综合SBAS-InSAR 结果和提取水体结果上处理的，因监测期较短，有些地表在短期内未发生变化，或是已经达到稳沉状态，形成的等值线会存在一定的误差。 如何消除以上2 种误差以及如何利用遥感技术快速监测矿区地表在整个采煤矿区内的沉陷情况还有待进一步探究。
5 结 论
1)针对高潜水位矿区的水体和其他地物的光谱特征，提出了E-MNDWI 沉陷水体提取指数，与MNDWI 相比，该指数能消除水体边缘植物以及矿区粉煤灰的影响，使水体范围更为清晰，获取面积更为准确。 截至2017 年，矿区的沉陷水体面积占整个煤矿区面积的13.08%。 结果表明，矿区沉陷水体面积依旧在扩大，年增长率约0.4 km2/a。
2)从空间角度看，矿区东部和南部存在面积的地表沉陷，监测期内鲍店煤矿公司所在地地表沉陷最为明显，且地表累计沉陷量由西南向东北逐渐增加，最大地表沉陷速率为41.69 mm/a；矿区西南有些区域地表表现为抬升，最大地表抬升速率为12.83 mm/a。 监测期内矿区整体非积水沉陷面积为10.1 km2，占整个矿区的28.28%。 从时间角度看，鲍店煤矿具有30 年的开采历史，在监测期内其地表依旧发生沉陷，地表沉陷问题严峻，矿区土地复垦工作迫在眉睫。
3)结合多源多时相数据，利用 E-MNDWI 和SBAS-InSAR 技术能够有效获取高潜水位矿区的地表沉陷信息，且精度较高。 这对快速高效提取采煤沉陷区域信息具有重要意义，也为后续沉陷区域的整治工作提供了参考依据。
[1]杨 璐， 胡振琪， 赵艳玲， 等.基本农田保护与煤炭资源开采的协调发展[J].中国煤炭， 2007， 33(9)：16-17.YANG Lu， HU Zhenqi， ZHAO Yanling， et al.Coordinated development of basic farmland protection and coal resource exploitation[J].China Coal， 2007， 33(9)：16-17.
[2]QU F， ZHANG Q， LU Z， et al.Land subsidence and ground fissures in Xi'an， China 2005 - 2012 revealed by multi-band InSAR time-series analysis[J].Remote Sensing of Environment，2014， 155：366-376.
[3]DUAN G Y，GONG H L，LIU H，et al.Monitoring and analysis of land subsidence along Beijing-Tianjin inter-city railway[J].Journal of the Indian Society of Remote Sensing， 2016， 44：915-931.
[4]闵祥宇， 李新举.高潜水位矿区不同复垦方式下土壤热导率及其影响因素[J].水土保持学报， 2017， 31(3)：176-181.MIN Xiangyu， LI Xinju.Soil thermal conductivity and influencing factors under different reclamation modes in a mining area[J].Journal of Soil and Water Conservation.， 2017， 31(3)：176-181.
[5]付艳华， 胡振琪， 肖 武， 等.高潜水位煤矿区采煤沉陷湿地及其生态治理[J].湿地科学， 2016， 14(5)：671-676.FU Yanhua， HU Zhenqi， XIAO Wu， et al.Subsidence wetlands in coal mining areas with high water level and their ecological restoration[J].Wetland Science， 2016， 14(5)：671-676.
[6]鹿献章， 杨义忠， 喻 根， 等.淮北市北区地表水体及采空塌陷遥感调查[J].安徽地质， 1997(1)：82-90.LU Xianzhang， YANG Yizhong， YU Gen， et al.Remote sensing survey of surface water and mining subsidence in northern area of Huaibei city[J].Geology of Anhui， 1997(1)：82-90.
[7]彭苏萍， 王 磊， 孟召平， 等.遥感技术在煤矿区积水塌陷动态监测中的应用：以淮南矿区为例[J].煤炭学报， 2002，27(4)：374-378.PENG Suping， WANG Lei， MENG Zhaoping， et al.Monitoring the seeper subside in coal district by the remote sensing：examples from Huainan coal district[J].Journal of China Coal Society，2002，27(4)：374-378.
[8]徐涵秋.利用改进的归一化差异水体指数(MNDWI)提取水体信息的研究[J].遥感学报， 2005， 9(5)：589-595.XU Hanqiu.A study on information extraction of water body with the Modified Normalized Difference Water Index(MNDW I)[J].Journal of Remote Sensing， 2005， 9(5)：589-595.
[9]李 晶， 李恩来， 苗慧玲， 等.基于 LandsatTM 数据的济宁市水体识别方法研究[J].环境科学与技术， 2013(9)：175-179.LI Jing， LI Enlai， MIAO Huiling，et al.Identification of water bodies in Jining city based on LandsatTM dataset[J].Environmental Science ＆ Technology， 2013(9)：175-179.
[10]肖 武， 李素萃， 徐建飞， 等.GF-1 与 TM 影像对高潜水位矿区地表监测的比较分析[J].中国煤炭， 2016，42(6)：119-125.XIAO Wu， LI Sucui， XU Jianfei，et al.Comparative analysis on surface monitoring high physical water level in mining area foundation up GF-1 and TM images[J].China Coal， 2016，42(6)：119-125.
[11]GABRIEL A K， GOLDSTEIN R M， ZEBKER H A.Mapping small elevation changes over large areas： Differential radar interferometry[J].Journal of Geophysical Research Solid Earth，1989， 94(B7)：9183-9191.
[12]LIU G， DING X， CHEN Y， et al.Ground settlement of Chek Lap Kok Airport，Hong Kong，detected by satellite synthetic aperture radar interferometry[J].Chinese Science Bulletin， 2001(21)： 1778-1782.
[13]CASU F，MANZO M，LANARI R.A quantitative assessment of the SBAS algorithm performance for surface deformation retrieval from DInSAR data[J].Remote Sensing of Environment， 2006，102：195-210.
[14]胡乐银，张景发，商晓青.SBAS-InSAR 技术原理及其在地壳形变监测中的应用[J].地壳构造与地壳应力文集， 2010(1)：88-95.HU Leyin， ZHANG Jingfa， SHANG Xiaoqing.SBAS-InSAR technology and its application in monitoring the crustal deformation[J].Bulletin of the Institute of Crustal Dynamics，2010(1)：82-89.
[15]ZHOU L，GUO J，HU J，et al.Wuhan Surface Subsidence Analysis in 2015-2016 Based on Sentinel-1A Data by SBAS-InSAR[J].Remote Sensing， 2017， 9(10)：982.
[16]董少春， 种亚辉， 胡 欢， 等.基于时序 InSAR 的常州市2015—2018 年地面沉降监测[J].南京大学学报：自然科学版， 2019，55(3)：370-380.DONG Shaochun， ZHONG Yahui， HU Huan，et al.Ground subsidence monitoring during 2015—2018 in Changzhou based on timeseries in SAR method[J].Journal of Nanjing University：Natural Science， 2019，55(3)：370-380.
[17]孙问娟， 李新举.煤矿区土壤有机碳含量的遥感反演与分布特征[J].水土保持学报， 2018，32(3)：328-333，339.SUN Wenjuan， LI Xinju.Remote sensing inversion and distribution characteristics of soil organic carbon content in coal mining area[J].Journal of Soil and Water Conservation， 2018，32(3)：328-333，339.
[18]丁 凤.基于新型水体指数(NWI)进行水体信息提取的实验研究[J].测绘科学， 2009，34(4)：155-157.DING Feng.Study on information extraction of water body with a new water index (NWI)[J].Science of Surveying and Mapping，2009，34(4)：155-157.
[19]肖艳芳， 赵文吉，朱 琳.利用TM 影像Band1 与Band7 提取水体信息[J].测绘科学， 2010，35(5)：226-227，216.XIAO Yanfang， ZHAO Wenji， ZHU Lin.A study on information extraction of water body using band1 and band7 of TM imagery[J].Science of Surveying and Mapping， 2010，35(5)： 226-227，216.
[20]倪 衡， 刘晓颖.一种基于Landsat8 遥感数据快速提取水体信息的新方法[J].科技资讯，2015，13(19)：7-9.NI Heng， LI Xiaoying.A new method for fast information extraction of water bodies using remotely sensed data based on Landsat8[J].Science ＆ Technology Information，2015，13(19)：7-9.
[21]BERARDINO P， FORNARO G， LANARI R， et al.A new algorithm for surface deformation monitoring based on small baseline differential SAR interferograms [J].Geoscience ＆ Remote Sensing IEEE Transactions， 2003， 40：2375-2383.
[22]LANARI R， ÓSCAR Mora Sacristan， MANUNTA M， et al.A small baseline DIFSAR approach for investigating deformations on full resolution SAR interferograms[J].IEEE Transactions on Geoscience ＆ Remote Sensing， 2012， 42(7)：1377-1386.
[23]周祺超.基于InSAR 的老采空区地表变形监测与分析[D].西安：西安科技大学，2017.ZHOU Qichao.Monitoring and analysis of surface defermation in old goaf based on InSAR[D].Xi’ an： Xi’ an University of Science and Technology， 2017.
