摘 要 本文先阐述了Kriging插值方法的基本概念；再基于离散的地质勘探取样数据，在实现三维地质体空间网格化的基础上，采用Kriging方法进行网格样品插值，并对插值结果进行了比较分析。最后结合MC构建等值面技术和OpenGL可视化技术开发了三维地质体空间网格插值计算与可视化模拟程序。
关键词 Kriging法；网格插值；三维可视化
0 引言
三维地质建模是三维地学研究的一项基础性工作，它对地质体形态构造研究、矿体品位储量估算、地下矿产勘查与采矿等都具有非常重要的意义[1-4]。由于地质条件以及勘查技术的局限性，三维地质建模与可视化研究所使用的往往是不完备的地质勘探数据，如何通过空间数据插值来生成较为完备的地质体数据，快速建立接近真实的地质体三维空间可视化模型，已经成为目前地学研究中的研究难题。
本文将结合实例地质勘探数据，基于Kriging空间插值方法[5-6]，在三维地质体空间网格化的基础上进行插值计算，并对插值结果进行比较与分析。再以此为基础，结合MC构建等值面技术和OpenGL可视化技术，实现三维地质体的三维形态模拟和可视化。
1 Kriging插值方法
克里金法（Kriging）最初由南非矿业工程师D. G. Krige针对金矿品位特点结合地质推估问题创立，后由法国数学家Georges Matheron完善，形成一套完整的理论体系。Kriging法基于地质统计学的思想，将空间上连续性变化的属性称为“区域性变量”，用协方差函数和半变异函数描述区域性变量的变化规律。Kriging法插值的数学模型如下：
2 三维地质体空间网格插值计算
2.1 实例数据处理与分析
本文所使用的实例地质勘探数据共有原始样品铜品位化验数据708条，以样品平均长度1.6m进行样品组合，以去除异常值并保证样品的均一性，共得到组合样品744个。分别对原始样品和组合样品进行统计分析，统计结果如表1所示。
2.2 空间网格插值过程
先对待插值网格进行遍历，以当前网格为中心，变程a为半径进行样品网格搜索。对落入搜索范围内的样品，代入（3）式计算各矩阵，再代入（2）式，进行线性方程组求解，得到权重向量。最后由（1）式求出当前待插值网格的估算品位。如果落入搜索范围的样品数为零，则当前待插值网格的品位等于零。
2.3 结果验证与分析
对Kriging插值方法进行交叉验证，即每次去掉一个已知网格的铜品位值，用剩余的已知网格样品对其品位进行估算，直到对所有已知网格样品都进行了一次估算，最终得到一组和实际值对应的估算值，对两组数值进行分析，进而对不同的插值方法进行分析与评价。
为了准确地比较两者的结果，采用了以下统计模型：均方根（Root Mean Square，RMS）和决定系数R2（相关系数的平方）。分别在50m×50m×50m规格网格和25m×25m×25m规格网格下进行交叉验证，结果如下：50m网格，RMS=0.283，R2=0.056；25m网格，RMS=0.3556，R2=0.0292。
统计结果较好地说明了Kriging空间插值方法在实例数据上的表现，在不同的网格间距下，Kriging仍能提供较优的线性无偏估计，较低的均方根和较高的决定系数表明具有较好的预测性能。
图2为25m网格间距下的Kriging网格插值交叉验证结果散点图，发现验证点均比较均匀地分布在理想直线两边，两者均有相当数量的验证点靠近理想直线。
3 三维可视化模拟
基于实例数据的Kriging空间网格插值结果，本文使用移动立方体（Marching Cubes，MC）算法[2]，结合OpenGL三维图形工具，采用VC++编程，使用Access数据库为数据源，开发了三维可视化模拟程序。Kriging法空间网格插值结果如图3所示，其中网格大小为50m×50m×50m，MC算法显示（铜品位阈值设为工业品位0.3%）。
4 结论
1）Kriging方法充分考虑了研究区域内变量的空间结构性和随机性，它比较适用于已知数据较多且相对平稳的插值场合。
2）基于三维空间插值算法（如Kriging法），并结合三维模拟算法（如MC算法）和三维可视化技术（如OpenGL技术），可以实现对地质体的快速三维形态模拟与可视化。工作周期短，而且能快速判断地质体的大致分布情况。
参考文献
[2]邹艳红，何建春.移动立方体算法的地质体三维空间形态模拟[J].测绘学报，2012，41（6）：910-917.
[3]郑通科，陈庆，袁峰，等.三维矿体模型快速构建方法研究[J].合肥工业大学学报：自然科学版，2015，（1）.
[5]郭艳军，潘懋，燕飞，等.自然邻点插值方法在三维地质建模中的应用[J].解放军理工大学学报（自然科学版），2009，10（6）：650-655.
[6]陈靓，蒋春曦，杨欣，等.基于Kriging方法的三维地质可视化系统构建[J].华南地震，2014，34（4）：61-66.