隐伏矿体三维可视化预测
1 引言
从国内外矿山的找矿成果和和我国大型老矿山的找矿勘探与研究程度分析，我国大多数老矿山及危机矿山仍具有巨大的资源潜力，围绕矿山基地开展找矿勘查是解决危机矿山接替资源问题的有效途径（彭省临等，2004）。但同时须看到，矿山基地潜在资源的寻找也面临着巨大的困难和风险，这是由老矿山具有勘探开发程度高、预测评价及找矿向深部三维空间发展、矿山生产对物化探异常信息产生干扰等固有特点所决定的，因此，迫切需要研究开发适应危机矿山找矿特点的矿产资源预测评价新方法。
矿产资源预测与评价经过约50年的发展，经历了起步、发展和成熟等阶段后，目前已进入科学化、定量化和数字化阶段。在起步和发展阶段，许多地质学家主要从资源总量评价、远景区评价、评价方法、计算机程序等方面进行了开创性的理论研究和应用工作（Allais et al.，1957；Slicher et al.，1960；Harris et al.，1970；Barry et al.1970；Sincair et al.,1970；Agterberg et al.，1971；Griffiths，1975；Singer，1976；McCammon 1976），在此基础上，逐步形成了较完善的矿产资源定量评价系列理论和方法，突出的代表性成果是，国际地科联IGCP 98专题推出的六种矿产资源定量预测评价方法，国际国内广泛推广应用的矿床统计预测理论及方法（Agterberg，1974；Zhao 1992；赵鹏大等，1994）、美国地质调查局倡导的“三部式”资源评价法（Singer，1993）、我国学者提出的综合信息预测（王世称等，1989,2000）和致矿地质异常预测方法（赵鹏大等，1991,1993）。自90年代开始，随着GIS空间信息技术的发展，矿产资源预测与评价进入数字化阶段，形成了以空间数据库和GIS空间分析为技术支撑、以“多元地学空间数据集成-多元成矿信息提取与融合-矿产资源潜力制图”为核心流程的矿产资源数字化预测评价方法体系（Bonham-Carter et al.，1990；Katz et al.，1991；Rencz et al.，1994；Knox-Robinson et al.，1997；Cheng et al.，1999；Harris et al.，2000，2006；Asadi et al.，2001； Chen et al.,2004；Zhou et al.,2007；Carranza et al.，2008；Cassard et al.，2008；池顺都等，1999；肖克炎等，2000；王全明等，2001；叶天竺等，2007），其理论和方法更趋完善和实用。
上述理论和方法，尤其是基于GIS技术的矿产资源定量评价方法，已成为目前开展矿产资源与评价的主流方法而广泛应用于区域矿产资源远景预测评价工作中。但上述方法在应用于危机矿山深部找矿预测时，由于需要向深部三度空间发展，会遇到矿产资源预测评价的三维空间问题，即无论是预测评价范围，还是评价模型及评价结果等，都不得不考虑真三维空间的要求。上述矿产资源定量评价方法形成和发展主要源于中小比例尺的全球性和区域性矿产预测评价,且其赖以依靠的GIS技术及软件仍然属于2维或2.5维的，因而，尚不能完全适应和满足危机矿山和老矿山可接替资源找矿向深边部三度空间发展的要求。为了将矿产资源定量评价理论和方法有效地应用于危机矿山和老矿山深部找矿预测，需要对上述理论和方法进行三维空间扩展或改造，以便支持矿山真三维空间下的隐伏矿体立体定位定量预测要求。
自90年代以来，真三维空间下的地质建模理论、方法和软件得到了飞速发展并已进入实用化阶段（Turner，1992；Houlding，1992，1994；Fisher et al.，1992；Moore et al.，2001；Sirakova et al.，2002；Wu et al.，2003；Wu，2004；Gong et al.，2004；Arens et al.，2005）。三维地质建模（3DGM）技术的实用化，为矿产资源定量评价的三维化奠定了技术基础。
本文针对危机矿山及老矿山深找矿预测的上述问题，在多年开展隐伏矿体立体定量预测探索性研究（毛先成等，1988a，1988b，1991，2009）的基础上，结合国家“十一五”科技支撑计划课题“铜陵地区危机铜矿山大比例尺定位预测技术”（编号：2006BAB01B07），以安徽铜陵铜多金属矿集区内资源紧缺的凤凰山铜矿田为示范基地，将三维地质建模（3DGM）技术引入矿产资源定量评价，通过成矿信息三维提取等关键技术的突破，成功地研究开发了一套适应于危机矿山深边部隐伏矿找矿特点的预测评价新技术—隐伏矿体三维可视化预测。
隐伏矿体三维可视化预测的核心流程：在地质数据集成和成矿规律分析的基础上，通过连续地质体（含矿化体）的三维建模与离散化，采用地质体形态分析、地质场模拟等技术进行成矿信息三维定量提取，建立控矿地质因素场模型，分析控矿地质因素和矿化分布之间的关联关系，建立反映控矿变量到矿化变量映射关系的立体定量预测模型，对预测区三维空间中的矿化分布进行定位定量预测，采用三维可视化模型表达预测成果（图14-11）。
2 矿体定位预测概念模型
为开展深边部隐伏矿体立体定量预测，需研究和总结研究区的综合地质与成矿规律，归纳矿体定位规律，建立矿体定位预测概念模型，为可视化立体定量预测提供知识驱动。
凤凰山矿田位于铜陵市东南约35km的凤凰山新屋里盆地，盆地中心为燕山晚期的新屋里岩体，主要由花岗闪长岩和石英二长闪长岩组成。盆地四周主要为二叠至三叠系的碳酸盐岩地层。围绕岩体沿接触带产出有凤凰山（又称药园山）、宝山陶、铁山头、仙人冲、清水塘和江家冲等多个铜矿床，其中，凤凰山矿床达到了中型规模，累计探明铜金属资源量约43万吨以上，其它矿床的规模较小（邵拥军等，2003；彭省临等，2004）。
根据凤凰山矿田矿床矿体的空间分布特征和成矿地质条件的分析与归纳，总结出矿体定位规律：
（1）现有矿床都分布在新屋里岩体的接触带上或靠近接触带部位。矿床基本上围绕新屋里岩体沿接触带呈等间距分布，在远离接触带的围岩中以及岩体的中心部位，未发现有成型的铜矿床，表明了新屋里岩体在铜成矿中的重要作用。
（2）从矿田的西北往东南方向，矿化深度有增加的趋势。矿田东部的仙人冲和南部的
江家冲矿床的产出深度比北部的宝山陶矿床和铁山头矿床的产出深度要大，而药园山矿床北部矿体的产出深度比南部要浅得多。
（3）北东向褶皱－冲断系统形成于印支期且在燕山晚期岩浆侵位过程继续活动，是控制岩体侵位和矿化定位的主体构造。
（4）矿体的产状基本上与所在部位接触带的产状一致。虽然矿体的形成明显晚于矽卡岩，但至今还没有发现斜切接触带的矿体，这充分表明：①接触带是矿石沉淀有利的化学环境；②岩体侵位期的韧-脆性构造系统在成矿热液聚集和定位中发挥了十分重要的作用。
（5）位于岩体西面的药园山铜矿床规模最大，其原因主要有三个方面：①在岩体西部边界附近有一条规模较大的走向北西的横向张性断层，这条断层在局部地段与接触带叠合，虽然其形成于新屋里复式岩体侵位以前，与新屋里复式褶皱是同一构造系统，但在岩体侵位期和侵位后，它一直活动，是一个重要的边界条件，对成矿流体和成矿岩浆活动起着重要的控制作用；②从整个新屋里岩体的形态和产状的变化规律来看，岩体西部的接触带较特殊，岩体超覆在围岩之上呈凹兜状；③从岩体接触带的产状、岩体和近矿围岩的变形构造特征来分析，该处是岩体侵位时的前沿，受岩体侵位冲击最大，也是成矿流体流动聚集部位。
（6）从矿床成因特征与控矿构造特征推断的成矿定位过程为：来自岩浆岩的高压成矿流体（等于或超过静岩压力），在构造变形场和温度场的双重控制下向岩体边界的扩容断裂带运移，由于扩容性断裂，特别是位于围岩中的扩容性断裂与地表连通，其间的流体主要来自大气降水（具静水压力），当岩浆来源的流体流入后，流体压力的突然降低和两种完全不同流体间的反应，导致矿石的沉淀。
从上述矿床矿体定位规律规律出发，总结出了凤凰山矿田矿体定位预测概念模型（表1）。该模型包括矿田范围内控制矿床矿体分布的主要控制因素及相关地质体，以及地质体（包括矿化）和控矿因素的建模与分析方法。
3 地质信息三维可视化建模
地质信息三维可视化建模，是指采用三维地质建模技术（3DGM），对研究的矿床对象包含的各种地质体（地层、构造、岩浆岩、矿化等）进行三维可视化模拟与表达。建立的地质体模型主要包括线框模型、块体模型。地质信息三维可视化建模是进行控矿地质因素定量分析和成矿信息定量提取的基础。
本课题采用Datamine Studio 3软件工具进行地质信息三维可视化建模，主要工作内容包括：（1）地质建模数据（Geodatabase）的构建；（2）剖面地质界线圈定；（3）地质体线框模型建立；（4）地质体块体模型建立。
3.1 地质建模数据构建
为建立地质体的三维模型，须预先准备好建模用的各种地质原始数据，如地形数据、勘
探工程数据、编录数据、化验数据、地质图件数据等。将这些专用于地质建模的数据按照一定的逻辑结构进行组织，并导入到三维地质建模软件中，统称为地质建模数据（Geodatabase）。三维地质建模软件提供地质建模数据的管理与三维可视化显示功能。
为得到凤凰山矿田的地质建模数据，进行了下述具体工作：（1）相关的原始资料的收集与整理；（2）地质图件的矢量化；（3）综合地质数据库建立；（4）地质建模数据的提取与导入。实际收集的主要原始资料包括：凤凰山矿田药园山矿床、宝山陶矿床、仙人冲矿床、铁山头矿床的各种区域地质、地质勘探与储量成果报告；钻孔、坑道等单项工程原始资料；地质综合与找矿研究资料；物化探报告与原始数据资料；药园山矿床生产探矿资料等。最终汇总导入到三维地质建模工具软件Datamine Studio中的凤凰山矿田地质数据包括：勘探工程数据、勘探线数据、DTM数据。勘探工程数据包括开孔数据表（Collar）、测斜数据表（Survey）、地质编录数据表（Geology）、样品分析数据表（Sample）。
3.2 地质体三维建模与可视化
基于地质建模数据进行三维地质建模的基本流程为：（1）导入地质建模数据（勘探工程数据、勘探线数据、DTM数据等）；（2）单项工程及勘探线剖面的三维显示；（3）按剖面人机交互圈定地质界线；（4）人机交互生成地质体线框模型；（5）生成地质体块体模型（含块体品位和金属量估计）。
三维地质建模获得的主要模型为线框模型和块体模型。在剖面地质界线圈定的基础上，通过对相同地质体的边界线依剖面顺序连接，即可建立地质体的线框模型。利用建模软件的三维栅格化功能，可对线框模型描述的地质体实体进行块体（体元）分割，即可得到地质体的块体模型。
对凤凰山矿田的全部地质体均建立了线框模型和块体模型，包括地形模型、新屋里岩体模型、地层模型、断层模型、矿体模型。图2所示凤凰山矿田地层三维模型。线框模型采用三角形文件和顶点文件描述，块体模型由原型表文件和模型表文件描述。
3.3 基于物化探技术的深部地质体推断
随着探查深度的增加，钻探工程等探查手段的花费成本也在不断增高，因此，采用成本相对偏低的物化探技术手段不失为一种合适的选择。在成熟矿集区尤其是资源危机矿山，隐伏矿体因产于地下深部，在地表或浅部开展物化探工作仅能获得微弱甚至近乎没有的指示信息，因此，物化探技术对深部隐伏矿体的探测效果较差。但深部地质体由于其体积和质量远大于隐伏矿体，在地表和浅部具有较强的物化探指示信息，所以，物化探技术可以较好地探测深部地质体，并进而根据成矿规律预测深部隐伏矿体的分布。
本研究利用凤凰山铜矿CSAMT扫面测量资料对地下-1000m深度范围内的地质体进行推断。凤凰山矿田的新屋里花岗闪长岩体（简称新屋里岩体）与灰岩地层的接触带具有低阻高极化率特征，因此，利用CSAMT可以较好地探测深部接触带与花岗闪长岩体的分布。图14-3a为根据CSAMT视电阻率剖面利用Datamine推断得到新屋里岩体的深部地质剖面；图14-3b为根据浅部地质勘探资料和深部CSAMT推断地质剖面建立的新屋里岩体三维线框模型。
图14-3 凤凰山矿田新屋里岩体的深部地质推断剖面与三维线框模型
(a)—利用CSAMT推断得到的新屋里岩体深部地质剖面；(b)—新屋里岩体三维线框模型
4 成矿信息三维定量提取
成矿信息三维定量提取，是指在综合研究成矿规律的基础上，采用各种数学建模与计算机处理手段，提取出指示矿化在三维地质空间中的质量和数量分布的各种信息指标，这些指标统称为成矿信息指标。成矿信息分为控矿因素信息（岩浆岩、地质构造、地层、岩性等地质条件）和找矿标志信息（物化探异常、遥感影像、重砂、蚀变等找矿标志）两类，相应地，成矿信息指标分为控矿因素指标（简称控矿指标）和找矿标志指标（简称标志指标）。
在危机矿山隐伏矿体预测中，由于物化遥信息对深部隐伏矿体的指示作用相对较弱，所以，本研究将物化遥手段获取的信息主要用于对深部地质现象的推断而间接地转化为深部地质信息（图14-1），因此，隐伏矿体预测成矿信息提取的关键是控矿地质因素的分析与提取。
成矿信息三维定量提取的步骤：（1）导入地质体块体模型数据作为成矿信息提取的原始数据；（2）定义地质空间和划分立体单元；（3）建立控矿地质因素场模型，实现控矿地质因素的三维栅格模型表达；（4）定量分析矿化分布与控矿地质因素的关联关系，构建控矿因素指标集。
4.1地质空间定义与立体单元划分
（1）地质空间定义
隐伏矿体预测是在一定的空间范围内进行的，当该空间范围赋予地质意义时，称为地质空间。地质空间G是地质体产出和地质作用发生的三维空间。矿化空间M是成矿作用发生的空间，为地质空间的子集。地质空间的确定，一般先定义一个巨大的立方体空间作为地质空间的包集，然后再用各种边界条件对立方体空间进行限制获得真正的地质空间。
考虑凤凰山矿田已有的地质工作程度、本次工作的目标及范围，取自地表最高点（标高600米）至地下深部-2700标高水平，为其垂直空间范围，整体垂深为3300米。定义地质空间的巨大立方体空间左下前角点（Xmin,Ymin,Zmin）的坐标为（39595200，3411000，-2700）、右上后角点（Xmax,Ymax,Zmax）的坐标为（39604000，3421800，600）。限制该立方体空间的边界条件为：（1）以地表面（由工作区DEM定义）为地质空间顶面、-2700标高水平面为地质空间底面；（2）根据工作范围定义地质空间的水平投影范围，水平投影范围多边形的角点坐标为（39597800，3411000）、（39604000，3417000）、（39603900，3417500）、（39600400，3421700）、（39595900，3418500）、（39596000，3411100）、（39592900，3411000）。考虑到预测外推的可靠性问题，将地质空间的第四系浮土层底面以下、-1000米标高水平以上的子集空间定义为矿化空间。地质信息建模和成矿信息提取在地质空间中进行，但矿化分布分析和隐伏矿体预测均限定在矿化空间范围内。
（2）立体单元划分
采用三维规则网格将地质空间划分为许多小立方体，称为地质空间的立体单元划分，小立方体称为立体单元（简称为单元，在不同场合，还称为体元、体素、块体）。
地质空间划分的精度取决于预测比例尺、建模工作精度等因素。同时考虑到预测精度的合适性要求、建模的高精度要求，本次研究选用两种精度划分地质空间：（1）建模精度（10
米精度）—立体单元尺寸为10×10×10m3
；（2）预测精度（50米精度）—立体单元尺寸为
50×50×50m3
。地质体块体模型、控矿地质因素场模型均采用建模精度，成矿信息指标提取、矿化预测均采用预测精度。
4.2 控矿地质因素场建模
为实现控矿地质因素指标的定量表达，提出控矿地质因素场的概念。控矿地质因素场是根据地质知识和地质经验建立的，反映的是控矿地质作用在地质空间中的结果与分布。控矿地质因素场与空间中某点到相关联的地质体的距离有关，即控矿地质因素场是到地质体距离的空间分布函数。在地质空间中，本研究选择欧式距离作为空间距离。地质体之间的距离或地质空间中某点到地质体的距离，用以表示和研究地质体之间的几何接近程度或地质体对空间中某点的影响程度。在研究中，我们将点到地质体的距离约定为点到地质体的最小距离，即预测空间中某单元（体元）到地质体的最近距离作为控矿地质因素场对单元的影响程度。控矿地质因素场采用三维栅格模型表达，所以，控矿地质因素场建模以地质体块体模型作为原始输入数据。为实现控矿地质因素场建模，研究开发了专门的控矿地质因素场建模算法和软件BlockModelTrans。
在控矿地质因素场模型建立过程中，为了提高分析与表达精度，地质体对象模型采用高精度的10m×10m×10m规格块体（栅格）模型，控矿地质因素场描述采用预测精度的50m×50m×50m规格栅格模型，即控矿因素指标在50m×50m×50m规格单元的赋值由位于其中的最多125个10m×10m×10m规格的有效单元的指标值进行整合计算得到。
根据凤凰山矿田矿体定位预测概念模型，控制矿体分布的主要地质体与地质因素包括新屋里岩体及其形态起伏、接触带、地层及褶皱构造、断层等，抽象出的描述这些控矿地质因素的场变量为：（1）岩体热力场因素dG；（2）岩体形态因素，即一级起伏wr1G和二级起伏wr2G；（3）接触带因素dI；（4）接触面构造因素，即原始接触面与趋势接触面的夹角aIT；（5）横向张性断层因素dF；（6）挤压远应力场因素aIP；（7）褶皱构造因素dD3。以下仅给出岩体热力场因素dG的栅格模型。
岩体热力场因素用新屋里岩体的距离场模型描述。首先从新屋里岩体块体模型获取岩体的体素模型，然后对体素模型进行三维欧式距离变换。以单元到岩体的最小距离作为距离场值来表达热力场因素，即可得到的新屋里岩体热力场因素变量dG（图14-4）。为区分岩体内外单元的距离场值的差异，将位于岩体分布范围之外的单元场值置为正（正距离场），位于岩体分布范围之外的单元场值为负（负距离场）。
4.3控矿地质因素与矿化分布关联分析
控矿地质因素与矿化分布关联分析是指，在地质空间中立体单元划分的基础上，通过计算和获取立体单元的控矿地质因素场变量值、已知单元矿化指标值，研究控矿地质因素场变量与单元矿化指标之间的关联关系，进行控矿地质因素场变量的非线性变换，提取控矿指标。
单元矿化指标包括：（1）单元铜平均品位Cu—落入某单元的全部取样样品的铜品位按样长加权平均值；（2）单元铜金属量CuOre—落入某单元的铜矿体的金属量；（3）单元含矿性IOre—落入某单元达到边界品位的取样样品的铜品位按样长加权平均值是否达到工业品位要求（1—达到，0—未达到）。控矿指标是控矿地质因素场变量经过非线性变换后得到的新变量，与矿化指标呈线性关联关系，因而，可直接用传统的统计分析方法建立控矿指标到矿化指标的关联模型。
提取控矿指标的基本方法是：（1）依据已获得的单元控矿地质因素场变量值和单元矿化指标值，生成反映单元控矿地质因素场变量与单元矿化指标关系的散点图；（2）基于散点图，构造矿化指标与控矿地质因素场变量的非线性关联数学模型，并对模型的统计效果进行检验；（3）按非线性模型对控矿地质因素场变量进行非线性变换，得到直接指示矿化分布的控矿指标。在控矿指标提取过程中，单元矿化指标计算、单元控矿地质因素场变量计算等采用自编软件（BlockModelTrans、VoxelCal等）完成，散点图的生成、非线性回归分析等采用统计分析软件SPSS 15.0完成。
控矿地质因素场变量dG与单元矿化指标Cu、CuOre的散点图如图14-5所示。
从dG-Cu、dG-CuOre散点图（图14-5）可知，单元铜平均品位Cu、单元铜金属量CuOre与因素dG之间存在着一定的关联关系，表现为：dG取值为(-250,150)的单元，其矿化指标Cu、CuOre的值最高或明显偏高，即这些单元为主要矿化富集空间（品位高、金属量大）；dG取值为(260,500) 的单元，其矿化指标Cu、CuOre的值相对较高，即这些单元为次要矿化富集空间（品位偏低、金属量偏小）；dG取值为其它值的单元，其矿化指标Cu、CuOre的值明显偏低，即这些单元为无矿或贫矿空间。
5 预测模型与可视化
根据上述成矿信息定量提取结果可知，矿化指标与控矿指标之间存在着一定的关联关系，所以，可以依据这种关系来建立预测模型。根据地质空间中的控矿地质因素分布，可以利用预测模型推断出矿化指标在空间上的分布，从而可以对分布在矿化地质空间中的隐伏矿体进行预测评价。
5.1 矿体立体定量预测模型
矿化指标在三维地质空间上的分布称为矿化分布，描述矿化指标的变量称为矿化变量。矿化变量包括：（1）Cu—单元铜平均品位；（2）CuOre—单元铜金属量。
控矿指标描述了控矿地质因素的成矿有利度，反映了地质控矿作用在三维地质空间上的分布，故称为控矿变量。不同的矿化变量对应地有不同的控矿变量：矿化变量Cu对应的控矿变量—ddG1，ddF1，daIP1，dwr1G1，aIT1，wr2G1，ddD31；矿化变量CuOre对应的控矿变量—ddG2，ddF2，daIP2，dwr1G2，aIT2，wr2G2，ddD32。
矿化指标与控矿指标的关联关系，在数学上可以表达为控矿变量（控矿指标）空间到矿化变量（矿化指标）空间的映射，而矿化变量空间中的每个矿化变量均属于某个有界的实数域，故可以将这种映射看作为泛函，其函数化表达模型为MV=f(GV)，式中MV为矿化变量空间，GV为控矿变量空间。该泛函关系可以通过多元回归等统计方法来实现函数化表达。
5.2 预测结果及可视化
利用上述预测模型即单元矿化指标回归模型和单元含矿性估计模型对未知区符合边界条件的单元的矿化指标（Cu与CuOre）和含矿性指标（Iore）均进行了预测估值，预测结果（预测结果数据表格略）存放在预测成果数据库中。
预测成果数据库是对立体单元基本属性表扩展而得到的，扩展的属性用于存储预测估值结果，包括矿化指标Cu的已知值Cu和预测值PCu、矿化指标CuOre的已知值CuOre和预测值PCuOre、含矿性指标IOre的已知值IOre和预测值PIOr。同时，为了便于对预测成果进行统一制图与可视化，还对矿化指标的已知值和预测值进行统一合并处理，矿化指标的已知值和预测值合并后用ECu、ECuOre、EIore属性表示，合并原则为：如果单元为已知单元，则ECu、ECuOre、EIore取已知值；如果单元为预测单元，则取预测值。
对预测数据库中的单元矿化指标（ECu、ECuOre）的预测估值结果，按照Datamine块体模型的格式，生成预测成果的三维栅格模型，实现了预测成果的三维可视化显示。尚未发现的隐伏矿体的单元矿化指标ECu和ECuOre预测结果栅格模型的可视化显示效果如图14-6和图14-7所示。
为了利用预测成果指导深部探矿工程布置，利用AutoCAD二次开发编程绘制了分层单元矿化指标预测成果图共27幅，利用Surfer绘制分层单元矿化指标预测等值线图共54幅，利用Datamine和AutoCAD绘制深部探矿工程设计地质剖面图5幅。为了方便三维可视化显示和查询，根据预测成果数据，编程开发了隐伏矿体立体定量预测成果三维可视化查询系统GeoPro，可实时地显示和查询每个立体单元的坐标(x,y,z)、单元地质代码（gcode）、单元矿化指标（ECu、ECuOre）和单元含矿性指标（EIore）。根据预测结果，在凤凰矿田新屋里岩体西缘即凤凰山矿区深边部圈定了4个深部找矿立体靶区（靶区具体信息略），为深部找矿工程的设计布置提供了尚未发现的隐伏矿体的位置、品位和金属量等信息的指导。
