基于地质统计学的采空区储量估算
报告比较少，地质统计学资源储量估算方法还没有在地质勘查和矿山行业得到普遍的应用。但是随着矿业全球化的发展，应用地质统计学方法进行地质矿产资源储量的评估已成为一种趋势[2-4]。采空区的三维模型的构建，可以为矿体开采定位提供详实可靠的依据，对整个矿山的开采进度进行合理控制。应用地质统计学计算采空区储量有助于企业用户准确估算矿山的开采量、损失量、保有资源储量、拟动用的资源储量等，提高估算精度和工作效率，减少开采风险。文中从模型切割出发，通过OBB树碰撞检测，找到和采空区相交的块体，然后进行布尔求交，通过体积加权进行储量计算。
1 采空区建模方式
1.1 三维地学建模技术三维地学模拟技术（3D GMS，3D Geosciences Modeling System）是固体矿产勘探三维可视化储量计算的基础和核心，也是三维GIS领域研究的热点问题。综合所有的方法，3D地学建模可以分为面元模型（Facial Model）、体元模型（Volumetric Model）和混合元模型（Mixed Model）3类[5]。1）面元模型的构建方法主要是侧重在空间实体的表面、地质层面、地下工程的轮廓，其中不规则网TIN模型和规则Grid模型，它们通常用于不封闭的表面模型，如DEM，DTM等;线框模型（Wire Frame）和边界（B-Rep）模型主要用于封闭表面或者外部轮廓模拟。采用面元表示的三维地学模型，便于显示和数据更新[6];
2）体元模型是基于3D空间的体元分割和实体表达，体元的属性可以独立描述和存储，可以分为规则和不规则2类。规则体元，如矢量建模方法结构实体几何（CSG）建模、八叉树（Octree）模型;不规则体元，如四面体格网（TEN）建模，用三维Delaunay法则将空间中的散乱点剖分成四面体;实体（Solid）建模，采用多边形网格描述地质和开采过程形成的形体边界，并用传统的块段模型描述地质体内部的品位或质量的分布;
3）混合模型：TEN+Octree建模，其中Octree作整体描述，TEN作局部描述;TIN+CSG建模，TIN模型表示地形表面，CSG模型表示建筑物实体，2种模型的数据是分开存储的[7]。
1.2 采空区建模采空区模型的建立主要是利用矿山提供的各类剖面图、中段平面图、钻孔及坑道测量编录资料而进行的。建立采空区模型，将采空区模型和任意采区以及勘探工程进行复合，可以清楚地看到矿体的的已采区域和未采区域，为矿山设计开采计划提供依据，还可以计算采空区模型的体积和储量，并与相应的矿体储量进行比较核实。这有利于矿山每年进行新增储量和消耗储量核算[8]。近年来，基于实测采空模型的CMS空区监测系统也是一种有效的探测手段，尤其是在危险和人员无法进入采空区，主要适用于井下采场及空区的探测和精密测量，将测量好的模型可以导入到其他系统中，生成模型，采用CMS系统探测空区效率高，探测结果可视化效果好，但是相对于价格较为昂贵，数据量大，操作也相对复杂[9]。采空区是通过线框模型构造，主要是基于在矿山开采过程中，采空区的设计都是较为规则，拟下一时间段的开采区域，不仅要分析近区开采的岩体稳定性，更要确定待开采区域的各种金属量是否符合要求，这需要不断调整区域的范围，对模型更新度要求比较高。根据采空模型的规划，可以将建模方式分为剖面法和中段法，线框法就是通过在二维视图勾划出一系列的采空区截面数据，然后通过轮廓线拼接算法，将相邻属于同一采空区的边界进行拼接，最终得到采空区的模型[10-12]。不同之处在于，剖面法是基于勘探线剖面，一般应用于矿山开采前期阶段;而中段法，主要应用在矿山深部开采的过程中。如图1和图2分别为2种采空区边界和线框实体模型。
2 地质统计学储量估算方法IExploration-EM是基于MAPGIS，综合了传统矿产资源储量估算方法、地质统计学的克里格法与3D建模技术，面向全国危机矿山接替资源以及其他固体矿产勘查项目的资源储量估算系统，实现储量计算的方法主要有：平行剖面法、地质块段法和地质统计学方法。
1）平行剖面和地质块段法，都是把形状复杂的矿体简化为与该矿体体积大致相等的简单几何体，并将矿化复杂状态变为在影响范围内的均匀化状态，以便采用简单的数学公式计算其体积和储量[13]。其优点在于简便和易于掌握，特别当工程数量少，对矿产储量进行概略估算。缺点是可靠性差，特别是当矿体形态和矿化复杂，工程控制不是很密集时[14-15];
2）地质统计学是以区域化变量理论作为基础，以变差函数作为主要工具，对既具有随机性又具有结构性的变量（如品位值）进行统计学研究。地质统计学核心为克里格估值，它是一种无偏的最小误差的储量计算方法。IExploration-EM系统在地质统计学资源储量方面，实现了普通克里格、泛克里格、指示克里格等方法。具体处理步骤为：首先对原始样品进行组合样划分，然后对组合后的样品进行数据分析确定样品分布形态，然后通过结构分析理解空间样品分布的相关性，最后通过合适的估值算法对属性模型进行块体赋值[16]。地质统计学建模流程如图3所示。
3 采空区储量计算参考矿体的储量计算方法，采空区可以应用平行剖面法和地质块段法，在二维剖面和中段视图中计算储量。但是由于采空区在矿山开采过程中，对于矿山储量的平衡和管理，矿山科学开采有重要的指导意义，而基于传统方法计算的矿山储量准确性较低。因此本论文，采用基于地质统计学储量计算得到的三維矿体晶胞模型，和采空区模型进行布尔求交，计算求交后晶胞模型的体积，通过体积加权方法的计算储量，累加统计整个采空区的储量信息[17-19]。计算流程如图4所示。
根据地质统计学的储量统计结果，采空区储量计算可以采用以下方法。1）采空区单元素储量统计，每次只计算采空区中一种元素的储量;2）基于表达式的多元素储量统计，将属于同一空块模型的多种元素，按照规则表达式分区间显示地质统计学统计结果模型，实现一次切割，统计采空区多种元素的储量。关键技术
1）型切割：地质统计学主要采用规则晶胞块体建模，数据量较大，而基于线框建模的采空区模型相对较为复杂，数据较少，进行模型切割之前，首先通过OBB树碰撞检测，找到和采空区相交的块体，然后进行布尔求交，求交的过程实际就是根据2个模型碰撞信息，计算边界，重新组成新的模型（切割后模型）;
2）体积计算：三维采空区模型的储量主要是通过体积加权的方法计算，因此求交后模型的体积计算方法，直接决定了采空区的储量计算精度。这里将采空区和晶胞块体求交后的每个模型进行四面体划分，然后基于四面体计算切割后模型的体积。这种方法主要适合于凸多面体模型，而对于凹多面体等极不规则的采空区模型，可以通过增加晶胞模型的精度，数据量虽然增长，切割过程较为缓慢，但是产生切割后模型都可以为凸多面体的，体积计算误差小，采空区的储量估算叶更为准确;
3）储量统计：根据地质统计学储量估算的结果，通过对切割后每个晶胞的体积加权，计算金属量。采空区某元素金属量计算见式（1）。
4 实际应用以新疆阿舍勒矿区为例，新疆阿舍勒矿床为铜锌硫化物矿产，伴生有Au，Ag，Pb，Ga，Se，S等6种有益的矿产。文中主要对该区I号铜矿体进行矿体品位分布与赋存储量的研究，该矿体主要产于细碧岩与石英角斑火山碎屑岩间的接触面上，在水平断面中为月牙状，横断面中为鱼钩状。矿体与上下地层整合接触，同步褶皱，形态受向斜构造的控制，呈向北倾伏、向南扬起、矿体东翼向西倒转的紧闭形态[20]。采空区储量统计，采用基于表达式的多元储量计算方法。表达式设置见表1.
部分储量统计输入见表2.h40012西采场、h45012#东采场、h50011 3个采空区的计算结果和实际勘探结果对比信息见表3，采空区切割剖面法矿体模型是指基于平行剖面法形成矿体，然后进行采空区和三维矿体模型布尔求交，并计算切割后采空区模型的储量。通过表2，表3可以计算采空区储量估算的2种方法的误差，这里为了直观体现，下面以铜（Cu）的储量建立直方图，如图7所示，其中A方法为采空区切割矿体模型计算方法，B为采空区切割地质统计学模型计算方法，C为勘探储量，通过对比可以发现，通过基于地质统计学方法计算的采空区储量结果，更为接近勘探储量的结果，这为矿山的开采设计可以提供有力的依据。
5 结 论
1）利用布尔求交方法得到采空区的晶胞模型，用插值法计算晶包的矿产元素含量，采用四面体方法对切割模型进行划分，减少了体积计算误差，通过体积加权方法实现地质统计学采空区估算。
2）应用于新疆某矿区验证，基于地质统计学模型估算结果与勘探结果接近率约为90%，与传统模型比较，地质统计学模型估算结果更精确，对随机性和机构性的变量适应性更强，使用该模型有助于精确估算资源量，提高估算精度，减少开采风险，为资源开发设计及精准开采提供指导。
参考文献（References）：
[1] 王 议，孙 伟，马冬梅，等.基于3D管理系统的煤炭矿山三维建模及实践——以陕北某煤矿为例[J].中国矿业，2016，25（7）：148-152.WANG Yi，SUN Wei，MA Dong-mei，et al.Three-dimensional modeling and practice of coal mine based on 3D management system：taking a coal mine in northern Shaanxi as an example[J].China Mining Industry，2016，25（7）：148-152.
[2]孫玉建.非参数地质统计学进展[J].地质与勘探，2007，43（4）：79-82.SUN Yu-jian.Advances in nonparametric geostatistics[J].Geology and Exploration，2007，43（4）：79-82.
[3]罗周全，王中民，刘晓明，等.基于地质统计学与Surpac的某铅锌矿床储量计算[J].矿业研究与开发，2010，30（2）：4-6，69.LUO Zhou-quan，WANG Zhong-min，LIU Xiao-ming，et al.Reserve calculation of a lead-zinc deposit based on geostatistics and Surpac[J].Mining Research and Development，2010，30（2）：4-6，69.
[4]姚改委，张静娜，李增辉，等.矿山采空区调查中地质手段方法应用浅析[J].科技创业，2013，26（3）：170-171，174.YAO Gai-wei，ZHANG Jing-na，LI Zeng-hui，et al.Brief analysis of the application of geological means and methods in mine goaf survey[J].Scientific and Technological Entrepreneurship，2013，26（3）：170-171，174.
[5]吴立新，朱旺喜，张瑞新.数字矿山与我国矿山未来发展[J].科技导报，2004（7）：29-31，28.WU Li-xin，ZHU Wang-xi，ZHANG Rui-xin.Digital mines and future development of mines in China[J].Science and Technology Report，2004（7）：29-31，28.
[6]刘修国，陈国良，候卫生，等.基于线框架模型的三维复杂地质体建模方法[J].地球科学，2006（5）：668-672.LIU Xiu-guo，CHEN Guo-liang，HOU Wei-sheng，et al.Three-dimensional complex geological body modeling method based on line frame model[J].Geosciences，2006（5）：668-672.
[7]査达剑，李乐林，江万寿，等.LiDAR三维重建中基于CSG方法的扩展研究[J].国土资源遥感，2016，28（4）：35-42.CHA Da-jian，LI Le-lin，JIANG Wan-shou，et al.Extension of CSG based on LiDAR three-dimensional reconstruction[J].Remote Sensing for Land and Resources，2016，28（4）：35-42.
[8]王钰东，余传涛，刘欣玥.煤矿双层采空区高密度电阻率法正演模拟[J].煤矿安全，2017，48（10）：212-215.WANG Yu-dong，YU Chuan-tao，LIU Xin-yue.The forward modeling of high density resistivity method in double-layer goaf of coal mine[J].Safety in Coal Mines，2017，48（10）：212-215.
[9]劉晓明，罗周全，杨承祥，等.基于实测的采空区稳定性数值模拟分析[J].岩土力学，2007，28（S1）：521-526.LIU Xiao-ming，LUO Zhou-quan，YANG Cheng-xiang，et al.[J].Geotechnical Mechanics，2007，28（S1）：521-526.
[10]Yavuz H Yavuz.An estimation method for cover pressure re-establishment distance and pressure distribution in the goaf of longwall coal mines[J].International Journal of Rock Mechanics and Mining Sciences，2004（2）：193-205.
[11]Singh A K Singh，Rajendra Singh，Maiti J.Assessment of mining induced stress development over coal pillars during depillaring[J].International Journal of Rock Mechanics and Mining Sciences，2011，（5）：805-818.
[12]J.N.van der Merwe.Fundamental analysis of the interaction between overburden behaviour and snook stability in coalmines[J].Journal of the South African Institute of Mining & Metallurgy，2005（1）：63-73.
[13]龚国清，刘修国，倪平泽.基于iExploration-EM的数字矿产勘查成果编制方法研究[J].金属矿山，2009，23（1）：34-38. GONG Guo-qing，LIU Xiu-guo，NI Ping-ze.Research on the compilation method of digital mineral exploration results based on iExploration-EM[J].Metal Mines，2009，23（1）：34-38.
[14]高帮飞.资源量估算的边界分析与矿化体圈定[J].地质与勘探，2018，54（2）：415-425.GAO Bang-fei.Boundary analysis of resource estimation and mineralization delineation[J].Geology and Exploration，2018，54（2）：415-425.
[15]苏 昆，周钟生，易 飞.基于地质统计学方法的某金矿采空区储量计算[J].现代矿业，2016，32（11）：144-146.SU Kun，ZHOU Zhong-sheng，YI Fei.Calculation of reserves in a gold mine based on geostatistics method[J].Modern Mining，2016，32（11）：144-146.
[16]刘海英，刘修国，李超岭.基于地质统计学法的三维储量估算系统研究与应用[J].吉林大学学报（地球科学版），2009，39（3）：541-546.LIU Hai-ying，LIU Xiu-guo，LI Chao-ling.Research and application of three-dimensional reserves estimation system based on Geostatistics[J].Journal of Jilin University（Geoscience Edition），2009，39（3）：541-546.
[17]田玉川.基于3DMine软件下地质统计学在矿山储量计算中的应用[J].中国矿业，2018，27（S1）：353-357.TIAN Yu-chuan.Application of geostatistics in mine reserves calculation based on 3DMine software[J].China Mining，2018，27（S1）：353-357.
[18]李 华，宁化展，丁希旺.基于Delaunay三角网矿体和采空区储量计算算法与实现[J].全球定位系统，2012，37（3）：82-85.LI Hua，NING Hua-zhan，DING Xi-wang.Algorithm for calculating reserves of ore bodies and mined areas based on Delaunay triangulation and implementation of[J].Global Positioning System，2012，37（3）：82-85.
[19]许 磊，刘永团.地质统计学在矿山储量估算中的应用——以内蒙古大苏计钼矿为例[J].资源与产业，2017，19（5）：65-70.XU Lei，LIU Yong-tuan.Application of geostatistics in estimation of mine reserves：taking Inner Mongolia Da Su Ji molybdenum mine as an example[J].Resources and Industry，2017，19（5）：65-70.
[21]张新宇，肖克炎，刘光胜，等.阿舍勒铜矿可视化储量计算的指示克里格法应用研究[J].吉林大学学报（地球科学版），2006，36（2）：305-308.ZHANG Xin-yu，XIAO Ke-yan，LIU Guang-sheng，et al.The applied research of indicator kriging method in visualized reserve computing of ashele copper mine[J].Journal of Jilin University（Earth Science Edition），2006，36（2）：305-308.