利用航测数据进行矿山生态修复土石方量计算模型的研究

张 鸽,垢元培

(河北省地质矿产勘查开发局第二地质大队(河北省矿山环境修复治理技术中心),河北 唐山 063000)

我国是矿产资源大国,矿山资源促进了经济的高速发展,同时矿产资源的开发利用也造成了严重的生态破坏和环境污染,成为制约我国经济、社会长远发展的重要因素。因此,开展矿山生态修复显得尤为重要。在对矿山特别是露天开采矿山进行生态修复时,要对其形成的坡面进行不同程度的处理,在保证边坡稳定性的前提下,进行生态环境修复措施[1]。矿山地质环境治理和矿山生态修复工程过程中涉及了大量的土石方量计算工作。土石方量计算结果的精确度不仅影响工程的投资预算和设计方案的优化,同时也影响到业主单位与施工单位的利益[2],以及工程的结算和施工方案的实施程度,是工程量计量的重要环节。因此,如何得到准确的土石方量计算成果显得格外重要。

对土石方量计算的研究大多为计算方法、计算原理或计算步骤,而少有从计算模型的角度进行研究。文献[3—4]研究了计算方法适用的条件、适用范围和精度要求等。文献[5—8]研究了土石方量计算原理、方法与步骤,介绍了几种主要应用软件、TIN模型模拟地面特征及DEM数据进行土方量计算的基本步骤和注意事项。

随着无人机低空摄影测量技术的快速发展,在矿山生态修复中无人机的利用也得到了空前发展。无人机摄影测量的优点是作业灵活,可快速获取高分辨率数字影像,同时受现场地质地理环境影响较小,工作中既提高了作业效率,又节约了大量的人力物力,大幅度提高了测量人员的安全保障[9]。

本文以无人机航测数据中的点云数据为基础,利用散点法的优势,以地形性质线为特征线,进行不规则三角网模型的建立,实现矿山生态修复中土石方量计算模型的构建;并通过与DEM数据计算的结果进行对比分析完成精度验证,从而快速、准确地得到土石方量计算成果。

1 土石方量计算数据类型

土石方量计算常用软件主要有南方CASS软件、飞时达软件、GIS(ArcGIS)软件和土石方大师软件4种。其中,南方CASS软件和飞时达软件均以CAD软件二次开发为基准的平台,均可采用三角网模型计算;GIS(ArcGIS)软件以DEM数据进行计算;土石方大师软件为国内创新团队自主研发软件,采用三角网模型计算,可实现三维化。这几款软件均可适用于矿山生态修复中土石方量计算。本文主要研究以无人机航测点云数据为基础,通过三角网法建立模型。

CASS软件系统是较常用的测绘软件,在土石方量计算中提供了方格网法、DTM法、等高线法和断面法等丰富的土方计算方法[10],对不同的工程和地形条件可灵活采用合适的土石方计算方法。而在矿山生态环境修复项目中主要采用DTM法建立土石方量计算模型。FastTFT(飞时达)是一款基于CAD平台二次开发的专业土石方计算绘图软件[11],针对各种复杂地形情况和工程实际要求,提供了方格网法、三角网法、断面法、道路断面法、田块法、整体估算法6种土石方量计算方法,同时对于土方挖填量的结果可进行分区域调配优化。GIS软件在土石方量计算中使用最多的为ArcGIS软件,ArcGIS 软件具有土石方量计算的工具。其计算原理与三角网法原理类似[12],采用克里金插值法,利用TIN和DEM数据通过ArcGIS软件多样化的功能模块,优化了传统方法的空间关系,从而使ArcGIS软件在土石方计算上具有优越性。土石方大师软件具有多图层、实景影像结合、三维立体化展示、模型数据种类多样、成果快速生成、智能化操作计算等优点,可计算多种三维模型的方量与面积等成果,为土石方量计算提供了丰富的三维模型展示[13]。此软件为国内创新团队自主研发土石方量计算专用软件,模型全面、计算精度高。

土石方量计算数据的类型主要有两大类,分别为散点法(点云法)和数字高程模型(DEM法),两种方法的特征见表1。

表1 土石方量计算数据类型特征

无人机点云数据是通过无人机航测手段,获取的具有位置信息和颜色信息的三维坐标数据集合。日常计算中通常对点云数据进行稀释后使用,在测量工作中主要采用三角网法。

数字高程模型(DEM)是以数字形式按一定的结构组织在一起的,能够表示实际地形特征、空间分布的一种实体地面模型,也是地形地理位置和地形起伏的数字描述[14],可用于土石方量的计算。DEM数据与无人机点云数据如图1所示。

图1 同一区域DEM数据与点云数据对照

2 计算模型建立方法的研究

无人机航测得到的点云一般用颜色表示高程,以反映地形。但在日常工作中,均需要根据精度进行点云抽稀,而抽稀后的点云数据为“规则无规律”的三维数据,无法真实反映复杂地形,具有局限性、容易造成数据处理错误、费时费力、受人为因素影响大、计算误差大等弊端。

应用点云数据计算土石方量时需要使用地性线。地性线也称地形性质线[15],是地貌形态的骨架线,是描述地貌形态时的控制线,这条线在实际地形中是不存在的,是在内业数据制图中得到的派生数据。

2.1 土石方计算模型的建立

采用Delaunay三角剖分与地性线边界虚拟点进行高程点云的有限元分析(finite element analysis,FEA)。应用此方法经三角剖分后含有最大化最小角特征及最重要的空圆特性,可生成最接近于规则化的三角形的三角网,同时具有点云中任意的4个点不共圆的唯一特性,这种方法明显优于常用的四叉树和八叉树方法及网格前沿算法。而地形性质线的使用有效弥补了航测点云的“规则无规律”性。将Delaunay方法与地形性质线相结合,弥补了因复杂地形造成的计算误差、计算错误等问题,尤其是在矿山生态修复项目等高山地地形。

空圆特性是指在一个三角形的外接圆范围内,包括外接圆范围线上,除了三角形的3个顶点,不再包含点云数据中的任何其他点。其中,四点共圆属于一个特例,空圆特性如图2所示。

图2 空圆特性示意

最大化最小角特性在本文中主要指,经抽稀后的无人机点云数据所组成的三角形中两个相邻的所组成凸四边形的对角线。最大化最小角特性如图3所示。

图3 最大化最小角特性示意

由于必须使用到地形性质线,应先将地形性质线按抽稀后的点云密度进行同间距或小于点云密度的间距插点,也可以称为地形性质线上的分割点。由于此分割点没有高程信息,需先根据地形线分割开的周围点云数据进行高程赋值。而三角网的建立采用“插值法”,即由外向内的方法。首先建立一个覆盖全部点云的三角形;然后将点云数据进行逐一插入计算,从而逐一找出符合空圆特性和最大化最小角特性的三角形;最后建立起不规则三角网数据模型。假设插入一点后的外切圆上的3个点坐标为A(x1,y1)、B(x2,y2)、C(x3,y3),则其外切圆的圆心(x0,y0)和半径(r)计算公式分别为

(1)

(2)

(3)

2.2 模型体积计算

在三角网数据中逐一选出三角形,其顶点坐标为A(x1,y1,h3)、B(x2,y2,h3)、C(x3,y3,h3),体积计算公式为

(4)

(5)

V=Sh

(6)

式中,S表示三角形面积;h表示高程;h0表示参考高程;V表示体积。

根据h的正负可以计算出挖方量和填方量,挖方量与填方量的和为总体积,最后将所有体积相加得到总方量。

3 矿山生态修复中土石方工程类型

根据工作实际和矿山生态修复工程的需要,可以将矿山生态修复中土石方工程按照矿山类型、矿山地形、破坏程度、修复方法、工作方法、施工流程等分为山体平行挖方或填方工程、依坡挖方或填方工程及山体挖削台阶工程。以上3种工程类型基本涵盖了矿山生态修复的不全类型,可真实、直观反映矿山修复后的地形地貌。

3.1 山体平行挖方或填方工程

山体平行挖方或填方工程指对独立山体进行挖方或其上平台进行填方,其工程部位不依附于周围山体边坡或地形,施工后对周围工程不产生影响,如图4所示。

图4 山体平行挖方前、后对照

3.2 依坡挖削方或填方工程

依坡挖削方或填方工程指进行工程部位一侧或四周依附于边坡或山体,当进行土石方量计算模型设计时需考虑边坡的影响。利用地形性质线划分出边坡与平台范围,形成依坡挖方或填方工程,施工后会增加或减少边坡面积,如图5所示。

图5 依坡挖削方前、后对照

3.3 山体挖削台阶工程

山体挖削台阶工程指对开采过后的矿山破碎边坡进行削坡降段,在模型建立过程中需使用地形性质线将挖削方形成的台阶与边坡进行分割,最终形成逐层台阶与边坡,如图6所示。

图6 山体挖削台阶前、后对照

4 矿山生态修复中土石方量计算模型

利用航测数据进行矿山生态修复土石方量计算模型的研究,点云数据是基础,而地形性质线(如平台边界线、坎顶线、坎底线、边坡边界线、台阶与边坡处分界线等)对土石方计算起到了非常关键的作用。

依据矿山生态修复中土石方工程类型和点云数据计算的特点,可将土石方量计算模型分成3类,分别为山体平行挖方或填方模型、依坡挖方或填方模型以及山体挖削台阶模型。模型建立流程如图7所示。

图7 模型建立流程

以唐山市某矿山生态修复治理项目中土石方量计算数据为例。本文项目航飞任务采用大疆Phantom 4 RTK无人机进行航测数据采集,采用ContextCapture软件进行空三数据处理,采用Cloud Compare或Global Mapper软件将点云数据处理成5 m间距,采用土石方大师软件进行点云数据土石方三维建模与方量计算,采用ArcGIS软件进行DEM数据的土石方量计算,最终通过对比分析得到点云数据计算的精度分析。

4.1 山体平行挖方或填方模型

山体平行挖方或填方模型主要适用于区域性整体性挖削方工程或区域性平台进行填方工程。此模型适用性最为普遍、通用性最强,由于此模型不依附任何地形,因此外部范围线即可作为地形性质线。

山体平行挖方或填方模型的建立方法是以挖填区域外部范围为地形性质线,以外圈地形性质线和点云为基础建立三角网模型,如图8所示,再根据三角网法计算两期间土石方量。

图8 山体平行挖方或填方模型

通过对实例数据的计算,使用点云法,采用山体平行挖方或填方模型所计算的土石方挖方量为25 288.1 m3。

使用ArcGIS软件对相同区域DEM数据进行土石方量计算作为“山体平行挖方或填方模型”计算方法的精度验证,如图9所示,计算结果为25 383.1 m3。

图9 同区域DEM数据模型

4.2 依坡挖方或填方模型

依坡挖方或填方模型主要适用于依附大型山体挖削方工程、基坑填方或依附山体的平台堆坡工程。此模型选择绘制地形性质线最为重要,点云间距越大地形性质线对计算的方量精度的影响越大。

依坡挖方或填方模型的建立方法是以挖填区域外部范围和坡底(顶)为地形性质线,以外圈地形性质线、坡底(顶)地形性质线和点云为基础建立三角网模型,如图10所示,再根据三角网法计算两期间土石方。

图10 依坡挖方或填方模型

通过对实例数据的计算,使用点云法,采用依坡挖方或填方模型所计算的土石方挖方量为31 978.7 m3。

使用ArcGIS软件对相同区域DEM数据进行土石方量计算作为“依坡挖方或填方模型”计算方法的精度验证,如图11所示,计算结果为32 026.4 m3。

图11 同区域DEM数据模型

4.3 山体挖削台阶模型

山体挖削台阶模型主要适用于依附山体的挖削方后形成逐层台阶的生态修复工程,也是土石方量计算模型中最复杂的一种。同时在矿山生态修复中的作用也较重要,通过点云数据与地形性质线的有效绘制,可以实现对台阶模型的精确建模。此模型与“山体平行挖方或填方模型”“依坡挖方或填方模型”的相同点为外部范围线均为地形性质线。

山体挖削台阶模型的建立方法是以挖填区域外部范围和台阶内侧、台阶外侧为地形性质线,以外圈地形性质线、台阶内测地形性质线、台阶外侧地形性质线和点云为基础建立三角网模型,如图12所示,再根据三角网法计算两期间土石方量。

图12 山体挖削台阶模型

通过对实例数据的计算,使用点云法,采用山体挖削台阶模型所计算的土石方挖方量为28 290.4 m3。

使用ArcGIS软件对相同区域DEM数据进行土石方量计算作为“山体挖削台阶模型”计算方法的精度验证,如图13所示,计算结果为28 065.5 m3。

图13 同区域DEM数据模型

5 模型精度分析

通过建立山体平行挖方或填方、依坡挖方或填方和山体挖削台阶3种模型进行土石方量计算。建立模型的点云数据为采用Cloud Compare或Global Mapper软件抽稀后所形成的5 m点云数据。无人机航测无法直接得到DEM数据,因此当地表植被和建构筑物较少时,可使用DSM代替DEM。精度分析使用的DEM数据为航测后ContextCapture软件形成的DSM数据,而非使用抽稀后的点云数据通过克里金法或三维分析法建立TIN数据后形成数字高程模型。通过实例数据分析,对比点云法与数字高程模型法计算的成果及误差,见表2。

表2 矿山生态修复中土石方量计算模型精度分析 m3

表中两种方法方量差占比计算公式为(数字高程模型法计算方量-点云法计算方量)÷点云法计算方量×100%。通过对计算结果进行分析可知,使用点云法通过对不同地形条件下模型的准确建立,与数字高程模型法进行比对分析后,山体平行挖方或填方模型的精度为0.38%,依坡挖方或填方模型的精度为0.15%,山体挖削台阶模型的精度为0.79%。对比精度可得,点云法计算方量与数字高程模型法计算方量之差的占比均小于1%,点云法计算的土石方量与数字高程模型法计算的结果精度一致,而非精度低。而点云法相比于数字高程模型法具有操作简单、适用软件多、要求专业技术水平较低等优势。因此,利用航测数据进行矿山生态修复土石方量计算模型的研究,更利于该方法的推广。

6 结 语

矿山生态修复过程中,土石方量的计算是一项基础性工作,也为生态修复后续工作提供了依据[16]。本文通过研究无人机航测点云数据特点,结合地形性质线,利用三角网建立适用于不同类型的土石方量计算模型。依据矿山生态修复中土石方工程类型和点云数据计算的特点,建立山体平行挖方或填方、依坡挖方或填方和山体挖削台阶3种模型。通过项目实例证明,在矿山生态修复中,针对不同地形条件建立合理的土石方量计算模型,采用点云数据计算的精度与采用DEM数据计算的精度相近,可快速准确地应用到土石方量计算中,有较大的推广与应用价值。