基于地面三维激光扫描的露天矿采剥工程量计算方法

王 森 何 群 刘善军 毛亚纯

（东北大学资源与土木工程学院，辽宁沈阳110819）

采场验收测量是露天矿山测量中的一项主要工作，而采剥工程量计算则是其中最频繁、最重要的工作，既为露天矿检查生产计划的执行情况提供依据，又为安排和调整未来生产计划提供不可或缺的资料。目前，露天采场验收测量手段主要使用全站仪对采场断面进行测量［1］，而后利用断面法计算采掘矿（岩）块体积，再与矿（岩）块体密度相乘，即可获得采掘量。由于该方法属于点式测量，因而存在采样密度低、体积计算误差大、测量周期长等不足。近年来，随着具有精度高、速度快等优势的面式测量手段（包括三维激光扫描、无人机航摄测量等）［2-3］的兴起，大大提高了矿山测量效率。

目前，三维激光扫描系统是较为先进的获取地面空间多目标三维数据长距离影像的测量技术，被称为测绘领域继GPS技术之后的又一技术革命［4］。三维激光扫描获取的数据集合——点云，是空间目标实体的真实表达，能够用于生成高精度DTM。但在大多数情况下，一次扫描范围不足以获得整个场景，需要利用部分重叠的高空间信息点云进行拼接［5］；另一方面，在进行多期点云变化监测、土石方量计算时，需要将不同坐标系统的点云归算到同一坐标系统下［6］，这就需要对点云进行精确配准。因此，点云数据精确配准便成为土石方量计算的一个关键环节。目前，Besl等［7］提出的ICP算法及改进ICP算法是应用最广、效果最佳、最为理想的点云精配准算法。但该算法对初始粗配准的输入要求较高，否则，会影响收敛速度，并且可能出现局部收敛现象［8-10］。Zheng等［11］提出了一种利用点云融合影像数据的配准方法，该方法利用相机和LiDAR扫描仪分别获取研究区域的影像和点云2种数据，将SIFT算法对影像提取出的特征点用作改进ICP算法过程中的对应点，以此进行配准迭代。该方法在配准过程中增加了一种影像数据源，并且在ICP精确配准中仅考虑了特征点配准，未顾及到全局点的配准计算。

本研究在分析上述成果的基础上，对原有SIFTICP点云配准算法进行改进，并利用露天矿三维激光扫描数据进行算法验证，实现基于地面三维激光点云数据的露天矿采剥工程量精准计算。

1 新SIFT-ICP算法原理

露天采场常用的采剥验收量计算方法有断面法、规则方格网法、等高线法、DTM（Digital terrain model）三角网法等［1］。DTM由Miller教授于1956年提出，能够反映地表真实特征，用于各种线路设计，各种工程面积、体积、坡度计算，任意两点间可视性判断以及任意断面图绘制等［12］。与其他几种方法对比，DTM能够全面真实地描述地形信息，差值计算采剥工作量的精度较高［13-14］，能够满足矿山验收测量的要求。DTM三角网法计算露天矿采掘工程量，是根据实测数据生成三角网，利用三角网计算拟定设计高程下的三棱柱体积，对指定区域内每个三棱柱体积进行求和得到整个区域在拟定高程下的体积总和，不同时间序列、同一区域、同等设计高程下的体积差值即为该区域的采掘工程量［15］。在整个过程中，点云精确配准直接影响了三棱柱体积的计算精度，进而影响了整个区域工程量的计算精度。

为构建研究区域完整的三维模型，需要将多个视点下获取的点云进行重新定位，生成一个统一坐标系下的三维数据点集，即点云配准［16］。点云配准分为粗配准和精配准，粗配准是对2个坐标系统相差较大的点云进行配准，粗配准后的点云能够达到坐标系统大体一致。ICP精配准算法是在粗配准的基础上进行全局优化配准，通过全局点不断迭代寻找对应点的关系，计算旋转矩阵R及平移系数T，使得所有对应点之间欧式距离的平方和达到最小，以完成配准。ICP算法精确配准的初始状态如果输入不佳，即会影响后续迭代解算，导致迭代收敛速度慢、迭代不收敛或无法收敛到正确结果，易导致出现局部最优解，使得配准误差偏大。

为优化配准算法，本研究首先设计了一种无靶标点的双线性内插SIFT特征点粗配准方法，然后将该方法的配准结果作为输入数据，再进行ICP配准，从而将SIFT算法与ICP算法的优势相结合，实现精确配准。由于该方法面向单一数据源（点云），而原有的SIFT-ICP算法面向2种数据源（可见光影像、点云数据）［11］，两者具有明显区别，因此该方法可称为新SIFT-ICP算法。该方法首先利用SIFT算法在实体模型上提取出特征点集，然后在特征点集上利用双线性内插方法求解点云点集对，以此进行点云粗配准。

SIFT算法由 Lowe［17］提出，该算法设计的图像特征描述子能够对于图像平移、旋转、缩放等特性保持不变。该算法的主要实现步骤为：①使用高斯函数构建不同尺度空间；②在高斯差分（Difference of Gaussian，DOG）尺度空间进行局部极值点检测；③为极值点分配128维方向参数；④生成特征点描述。

ICP算法是对SIFT粗配准结果进行基于最小二乘运算的最优匹配方法［18］，通过重复进行“确定对应点集、计算最优刚体变换过程”等步骤，直至某个表示正确匹配的收敛准则得到满足为止。该算法主要实现步骤为：①求解待配准点云在参考点云中的最近点，记为对应点集，计算最近点集距离平方和；②根据对应点集计算转换参数，由转换参数对待配准点云进行坐标转换；③转换后的点云与参考点再次求解最近点集，计算最近点集之间的距离平方和，若2次距离平方和的差值小于既定阈值，则停止迭代，否则，重复①～③步。可见，SIFT算法能够实现特征点检测，提取对应点关系，ICP算法能够完成最优匹配。

综合上述分析，新SIFT-ICP算法的主要实现步骤为：

（1）采用Delaunay三角剖分算法构建露天采场DTM，并对模型进行栅格化处理。

（2）在栅格化的模型上采用SIFT算法检测像平面直角坐标系下的特征点集。

（3）对于在像平面直角坐标系中的任意一点（x，y）（图1），采用双线性内插方法进行转换，将像平面直角坐标系x-y上的特征点（x，y）转换至点云对应的大地空间直角坐标系X-Y-Z在XY面上的投影，即XY平面直角坐标系下的坐标（X，Y），

式中，x1、y1为DTM 栅格化后的像平面长与宽，（X1，Y1）与（X2，Y2）分别为像平面坐标（图1）。根据X-Y平面直角坐标系下的坐标，在研究区域中的DEM（图1）上搜索确定出Z值，便可得到点云在大地空间直角坐标系中的特征点集。

（4）根据特征点集，对原始点云数据进行粗配准。

（5）由粗配准后的点云数据进行ICP点云数据精配准。

2 算法试验

2.1 数据获取

本研究试验区域位于辽宁省鞍山市鞍千矿矿区，该矿始建于2005年，于2006年8月正式投产。矿区由3个采区组成，分别为许东沟、哑巴岭、西大背，3个采区均已形成固定帮坡。试验区域哑巴岭采场最低露天开采标高为-132 m。试验采用Riegl VZ-4000三维激光扫描仪，最大测量距离为4 000 m，最近测量距离为5 m，垂直视场扫描范围为100°（+30°～-30°），水平视场扫描范围为0～360°，角度分辨率优于0.000 5°，测点精度为15 mm。采用该型扫描仪对露天采场进行精细扫描，分别获取2期三维激光扫描数据，数据采集时间与采场验收测量同步（图2）。

2.2 SIFT-ICP算法应用

2.2.1 SIFT算法实施

对于完成预处理后的数据进行Delaunay三角剖分建模，并对栅格化后的模型采用SIFT算法进行特征点集提取，而后进行同名点匹配，匹配结果如图3所示。

由于鞍千矿哑巴岭采场经过多年开采，外围已经形成固定边帮，模型匹配得到的点集主要分布于矿坑边帮区域。由于SIFT算法采用128维特征描述的局限性，匹配点会存在一些错误点，有必要删除。错误点识别主要考虑两个方面：

（1）局部区域是否满足平行关系，由于边帮匹配得到的特征点没有发生变形，在局部区域仅存在点云坐标系引起的系统误差，没有匹配时的离散随机误差，匹配点之间的连线应该满足平行关系，本研究将该类点判定为有效特征点。

（2）对于分布于矿坑内部的特征点，由于进行了采掘工作，匹配得到的特征点之间的连线与周围边帮难以构成平行关系，故判定为奇异点；对于图中还有部分点已经匹配到模型外围，可视其为无效特征点。本研究经过错误点剔除后保留了矿坑外围的97个有效特征点建立对应点集。

2.2.2 粗配准

依据SIFT算法匹配以及筛选后得到的特征点集，采用双线性内插方法，计算出对应点集在大地空间直角坐标系下的点，根据对应点集对原始点云数据进行粗配准。以初期点云数据为参考，对末期点云数据进行配准，结果如图4所示，其中淡色区域对应为第一期点云数据，深色区域对应为第二期点云数据。对于经过粗配准后的点云叠加结果（图4（b）），由于特征点集主要分布于右下边帮，在矿坑右下边帮中直观上能够看出深浅相间的点，粗配准效果较好。

2.2.3 ICP精配准

SIFT算法提取出的特征点数量有限，难以实现矿坑整个区域的精确配准。矿山生产过程中，由于采掘作业、道路铺设、矿车碾压等都会对矿区形态造成影响，使得不同时期点云之间的形态相差较大。为此，在SIFT算法粗配准的基础上，采用了一种基于稳定区域的全局匹配方法。该方法通过选取时间序列段内稳定不变的区域，在该区域内应用ICP精配准算法求解转换参数，再根据转换参数对全局区域进行转换配准。

考虑到矿坑内部存在非稳定点，难以使用ICP算法进行精配准，因此从点云中裁剪出边帮稳定区域，稳定区域根据粗配准结果的差值进行DEM圈定。对粗配准后的区域进行ICP算法精配准，得到边帮稳定区域的配准结果，如图5（a）所示，其中淡色区域对应为第一期点云数据，深色区域对应为第二期点云数据。在此基础上，将ICP算法得到的最终转换参数应用到整个矿区点云，得到整个区域的配准结果，如图5（b）所示。

将点云数据精配准结果与粗配准结果进行了对比分析，结果见表1。由表1可知：采用SIFT算法提取模型的特征点，用于配准的特征点有97个，该类点主要分布于边帮区域，并且在空间分布上也相对均匀，能够满足粗配准要求；利用稳定的边帮区域进行精确配准，2期配准点个数为847 874个，配准点数大大增加；2种方法的均方根误差RMSE相差1.169 105 m。

2.3 采掘边界线圈定

在同一露天矿中，一般是多个采掘区域同时进行采掘，不同的采掘区域一般是在不同的台阶上。在进行采掘工作量计算过程中，需要分别计算不同采区的工程量。为此，需要首先圈定出不同采区的范围。根据地面三维激光扫描获取的点云数据进行采区范围圈定，具体步骤为：①由配准后得到的2期三维激光点云数据，分别利用插值方法得到2期DEM，即：DEM1、DEM2；②根据 2期 DEM，计算出DEM之间的差值，ΔDEM=DEM2-DEM1；③由ΔDEM内插出指定高差的等值线，选择高差为0 m的等值线。

根据精配准后的点云数据，内插出了2期DEM，如图6所示。

对上述2期DEM计算出ΔDEM，而后圈定了高差为0 m的等值线图（图7）。由图7可知：由于受到道路铺设、矿车碾压等因素的影响，导致提取得到的等高线与实际采场范围略有偏差，故对等高线进行了修剪，进一步圈定了采场的大致范围，图7中圈定的深色区域可视为采场的主要采剥区域。

2.4 采掘量计算

采用Delaunay三角剖分算法构建露天采场DTM，图8为2期DTM叠加效果，其中圈定的不规则深色区域对应为第1期DTM，其余区域对应为第2期DTM。通过该图能够清晰看到采掘范围，与图7中等值线范围大体一致。此外，在该图上虽然能够看到有一些星星点点的变化点，并且该类点对可视化以及等值线提取会有一定的影响，但在圈定采掘边界线后，对于采场大范围工程计算来说影响较小。

根据采剥区域的不同岩矿类型，将整个区域划分为6个区域，如图9所示。图中各区域的采剥工程量计算结果见表2。由表2可知：基于点云数据计算出的结果与矿山实际记录的采矿量偏差较小，根据《有色金属矿山生产技术规程》，当矿岩验收量在20万t以上时，验收量允许误差不大于1%，该方法计算结果的相对误差为0.72%，能够满足矿山计算要求。此外，本研究采用三维激光扫描方法获取点云数据，有助于大幅提升露天矿山验收测量效率。

3 结语

针对ICP精配准对初始输入点云要求较高，以及原有SIFT-ICP算法需要使用影像数据作为SIFT算法粗配准数据源的问题，将原有SIFT-ICP算法进行了改进，提出了新的SIFT-ICP算法。该算法仅使用点云一种数据源，将SIFT算法的快速粗配准与ICP全局精配准优势相结合，实现了点云数据的快速和精确配准，能够在不需要设定靶标点、不需要影像辅助数据输入的情况下，仅使用三维激光扫描数据完成快速精确配准。基于新的SIFT-ICP算法并利用地面三维激光点云数据进行了实际露天矿山验证试验，结果表明，该算法能够快速精确完成配准，计算的采剥工作量相对误差为0.72%，能够满足矿山采场验收测量要求。