神东矿区三维地形定量变化检测研究

刘威辉，余代俊，李鸿洲，刘 凯

（1.成都理工大学，四川 成都 610059；2.国家测绘地理信息局卫星测绘应用中心，北京 100048）

神东矿区三维地形定量变化检测研究

刘威辉1，余代俊1，李鸿洲2，刘 凯1

（1.成都理工大学，四川 成都 610059；2.国家测绘地理信息局卫星测绘应用中心，北京 100048）

针对二维变化检测方法不再适用于当前三维变化检测的问题，提出一种三维地形定量变化检测方法。方法通过构建三维场景确定变化区域，采用离散化积分方法计算出变化区域的面积和土方量。实验结果表明该方法不仅可以直观地显示出变化区域的三维场景，而且能够定量计算出面积和土方量，可应用于矿区的遥感动态监测。

三维地形定量变化检测；离散化积分；遥感动态监测

我国是一个矿产资源丰富的国家，为了保证社会经济的可持续发展，必须对矿区地表变化进行有效检测和管理[1]。实验结果表明，利用遥感技术对矿区进行变化检测更为有效可行[2,3]。目前，二维变化检测研究已经非常成熟[2-5]，而三维变化检测的相关理论和方法还存在一定缺陷[6-13]。本文提出一种三维地形定量变化检测方法，来实现研究区的三维地形定量变化检测。

1 三维地形定量变化检测方法

1.1 三维地形定量变化检测流程

本文三维地形定量变化检测方法包括定性研究和定量计算两个过程[14]。首先，通过构建新旧两个时相的模拟三维场景，确定发生变化区域的位置，即定性研究过程；然后计算得到变化区域的面积和土方量，即定量计算过程。流程如图1所示。

图1 三维地形定量变化检测流程图

1.2 高程变化阈值的确定

数字表面模型（Digital Surface Model，DSM）相减后得到的差值即为实际高程差值。由于在DSM匹配过程中存在一定的高程误差，为得到关注的变化量，需要设定一个阈值ε作为高程变化的认定指标。当变化值超出这一阈值时，则认定该对象区域为变化区域，否则为未变化区域[15,16]。由研究区高程差值影像的直方图可以知道，研究区高程差值服从正态分布，再根据最大似然估计法，可以得出总体的均值和标准差的最大似然估计值是样本的均值和标准差。因为发生高程变化的区域面积占研究区总面积比例较小，未发生高程变化区域的一个样本的均值和标准差可以作为总体的均值和标准差的最大似然估计。又根据正态分布中的 3σ法则，研究区内的点落在（μ-3σ ,μ+3σ）范围内的概率高达99.74%，未落在此置信区间的点可视为发生高程变化的点。因此，本文采用未发生高程变化区域的一个样本的均值μ和标准差σ作为高程变化阈值的确定参数，阈值的区间为（μ-3σ, μ+3σ）。

1.3 变化区域面积和土方量计算

DSM存储的是离散的栅格点数据，采用离散化积分方法，选取研究区的DSM为数据样本，根据高差提取影像计算变化区域的面积和土方量[17,18]。面积的计算方法如下：

变化区域的面积=变化区域的像元数×影像分辨率×影像分辨率。

土方量的计算方法如下：假设某一研究区域为D，其中某点高程的变化函数为：

式中，ΔH表示某点高程的变化值，表示高程变化函数。由积分公式可以求得体积，函数表示为：

式中，V表示体积，是一个计算体积的三重积分表达式。具体计算过程是：假设用网格把闭区域D分割成n个小的闭区域σ1，σ2，…，σn，将σi的面积也记作σi，乘以相应小闭区域的高程变化ΔH，可得体积变化ΔVi。将研究区域的每一个细小圆柱体的行列号记作k和l，体积表示为：

1.4 三维地形定量变化检测方法的优势

1）采用DSM而不是DEM进行高程计算。DEM只包含了地形的高程信息，而DSM是在DEM的基础上，进一步涵盖其他地表信息如地表建筑物、桥梁等的高程，用DSM进行高程计算更准确。

2）使用的资源三号卫星数据高程精度比较高。资源三号卫星无控制点高程定位精度优于5 m，带控制点高程定位精度优于3 m，完全满足1︰50 000测图精度要求。

2 实验与分析

2.1 研究区及数据概况

本文研究区神东矿区地处内蒙古自治区鄂尔多斯市的伊金霍洛旗及东胜区南部和准格尔旗的西南部，在38°50'～39°50'N、110°00'～111°00' E 之间。由于矿区活动导致的地表形变非常明显，因此本文选取该地区作为研究区。本文所用实验数据是该地区2012-08和2013-05两个时期的资源三号卫星影像数据。

2.2 资源三号卫星

资源三号卫星是我国第一颗高分辨率光学立体测绘卫星，正视影像分辨率为2.1 m，多光谱影像分辨率为5.8 m，前后视影像分辨率为3.5 m，前后视相机与正视的夹角为22°，基高比为0.89。采用推扫式拍摄模式，可稳定获取同轨立体影像。

2.3 利用PixelGrid软件制作DSM、DOM

利用PixelGrid软件制作DSM、DOM的流程如下：

1）对原始影像进行连接点的匹配，再进行区域网平差，平差结果的单位权中误差在0.5像素以内，最后生成DSM。

2）对全色和多光谱影像分别进行数字正射影像（Digital Orthophoto Map，DOM）纠正，控制点的点位误差控制在0.5像素以内，然后进行影像融合处理，得到DOM。

2.4 三维地形变化检测

将生成的2012年和2013年的DOM、DSM 在ArcScene中构建模拟三维场景，叠加显示可以清楚地看到高程增加、高程减少和高程未发生变化的各个区域。例如某一区域的三维地形变化如图2所示，可以清楚地看到该区域的高程增加了。根据三维地形变化检测的结果，可以确定变化区域的位置、类型和大小。

图2 某一区域的三维地形变化图

2.5 三维地形定量变化检测

根据2.4中确定变化区域的位置在影像上建立感兴趣区。研究区2013年的感兴趣区影像如图3所示。

图3 感兴趣区影像图

将感兴趣区2012年和2013年的DSM相减，得到高程差值影像。在高程差值影像上均匀地选取了74个控制点，剔除错误点和发生高程变化的点，剩余61个点。这61个控制点是未发生高程变化区域的一个样本。根据最大似然估计法，将这61个控制点的高差值的均值和标准差作为总体的均值和标准差。实验得到，这61个控制点的高差值的均值为0.28 m，标准差为1.24 m。根据 3σ法则，得到高程变化阈值区间为（-3.44，4.01）。

图4 高差提取影像

图4 显示了高差提取的结果，黄色（例如区域A）表示高程未发生变化的区域；浅绿色（例如区域B）表示高程减少量在3.44 m以上的轻度沉降区域；绿色（例如区域C）表示高程减少量在10 m以上的严重沉降区域，这种高程变化极有可能是采煤区的挖方作业导致，是关注的重点；橙色（例如区域D）表示高程增加量在4.01 m以上的轻度沉降区域；红色（例如区域E）表示高程增加量在10 m以上的严重沉降区域，这种高程变化极有可能是煤渣堆积区的填方作业导致，是关注的重点。

2.6 变化区域面积和土方量计算

通过式 （3） 计算得到各高程变化区域的面积和土方量。高程减少区域有5个，面积和土方量如表1所示。高程增加区域有6个，面积和土方量如表2所示。为了对高程变化区域有更直观的了解，本文制作了高程变化区域面积和土方量三维柱状图（如图5所示），其中，黄色表示高程未发生变化的区域，绿色表示高程减少区域，红色表示高程增加区域。从图5不仅可以直观地看出哪些地方发生了变化，发生变化的地方高程是增加了还是减少了，而且可以清楚地看到各个变化区域的面积和土方量的大小以及各个区域之间的对比。

表1 高程减少区域的面积和土方量

表2 高程增加区域的面积和土方量

图5 高程变化区域面积和土方量三维柱状图

3 结 语

本文利用神东矿区的资源三号卫星影像数据，通过生成该地区的DSM和DOM，构建出研究区的模拟三维地形场景，确定了变化区域的位置、类型和大小；对感兴趣区变化前后两个时期的DSM作差值运算，运用统计学方法确定了高程变化的阈值，提取出变化区域；并利用离散化积分方法，较为准确地计算出变化区域的面积和土方量，实现了研究区三维地形定量变化检测，证明所提方法的正确性与可行性。

[1]孔改红,李富平.采煤沉陷区环境保护与治理研究[J].河北煤炭,2006(1):4-5

[2]侯鹏.基于遥感和地理信息系统的矿区环境监测与评价[D].青岛:山东科技大学,2004

[3]李德仁.利用遥感影像进行变化检测[J].武汉大学学报(信息科学版) , 2003, 28(特刊) :7-12

[4]赵英时.遥感应用分析原理与方法[M].北京:科学出版社,2003

[5]Peter A R. Change Detection Thresholds for Remotely Sensed Images[J]. Journal of Geographical Systems, 2002(4) : 85-97

[6]卜丽静,寇可心.基于资源三号数据的三维变化检测可行性研究[J].测绘通报,2013,1(12):31-35

[7]Tian J,Chaabouni-Chouayakh H, Reinartz P, et al. Automatic 3D Change Detection Based on Optical Satellite Stereo Imagery[C]. ISPRS TC VII Symposium 100 Year ISPRS,Vienna, Austria, 2010

[9]冯甜甜,王密,潘俊.基于物方影像匹配的DEM聚类变化检测[J].测绘信息与工程,2008,33(1):35-37

[10]曹建军,骆莹,李德仁,等. 利用多源空间数据进行的DEM变化检测方法研究[J].矿山测量,2011(5):17-20

[11]杜培军,柳思聪.融合多特征的遥感影像变化检测[J].遥感学报,2012,16(4) :663-677

[12]Lu D, Mausel P, Brondizio, et al. Change Detection Techniques[J]. International Journal of Remote Sensing, 2004,25(12) :2 365-2 407

[13]Radke R J, Andra S, Alkofahi O, et al. Image Change Detection Algorithms:a Systematic Survey[J]. IEEE Transactions on Image Processing,2005,14(3) :294-307

[14]夏松,李德仁,巫兆聪.利用多源空间数据进行地形的三维变化检测[J].测绘科学,2007,32(1):49-50

[15]马国锐,李平湘,秦前清.基于融合和广义高斯模型的遥感影像变化检测[J].遥感学报,2006,10(6):847-853

[16]Li D R. Remotely Sensed Images and GIS Data Fusion for Automatic Change Detection[J]. International Journal of Image and Data Fusion. 2010, 1(1):99-108

[17]Malpica J, Alonso M. A Method for Change Detection with Multi-temporal Satellite Images Using the RX Algorithm[J].International Archives of Photogrammetry, Remote Sensing and Spatial Information Sciences, 2008(37):1 631-1 635

[18]李畅, 熊昊, 陶顺勇,等.地质灾区光学立体影像植被与三维地形定量变化检测[J].遥感学报,2014,18(6):1 258-1 267

P237.4

B

1672-4623（2017）12-0022-03

10.3969/j.issn.1672-4623.2017.12.007

2016-08-03。

国家自然科学基金资助项目(41271394) 。

刘威辉，硕士研究生，研究方向为遥感原理与应用。