摘要
克立金估值技术,源于采矿设计阶段对矿块品位进行准确估计的问题,最早由南非采矿工程师D G Krige提出,即根据待估矿块内外不同矿石样品对矿块的平均品位进行最佳拟合与估计。法国数学家Matheron (1963)对Krige等人的研究成果进一步的系统化和理论化,形成地质统计学理论。经过数十年的研究与发展,地质统计学已逐渐形成一套较为完整的理论和工作方法应用于矿产、地质、地球物理、地球化学、水文地质、环境地质、气象和农业等各领域。目前国际矿业界已将地质统计学克立金估值方法作为固体矿产资源储量估算与评价的标准方法。
地质统计学引入我国矿业与研究领域已近30年,但据报道,目前以地质统计学估值方法提交的固体矿产地质调查/勘查资源储量报告才只有十几份，可以说地质统计学资源储量估算与评价方法还没有在我国地质矿产勘查和矿山行业得到普遍的应用。但是随着矿业全球化的发展,应用地质统计学方法进行地质矿产资源储量的评估已成为一种趋势,目前国外矿业市场只承认利用克立金法估算证实的矿产资源储量,所以若不使用地质统计学法,则地质勘查成果和储量资产就无法得到国际同行的认同,这对于我国开展国际矿业合作是非常不利的。然而地质统计学资源储量估算方法在地勘和矿山行业一直无法得到广泛的应用,究其原因主要有以下几点：
(1)传统方法仍然是我国地质矿产勘查资源储量估算的主要方法,地质统计学法要被广大地质矿产勘查研究人员接受还需要一段过程；自建国以来,我国一直沿用原苏联提出的地质块段法、剖面法等(简述“传统方法”)资源储量估算方法进行固体矿产资源储量估算。传统方法由于操作性较强、应用简单、计算方便,经过了近50年的应用已形成一套比较完整的实践方法和理论体系,而且传统方法在预查、普查等勘探程度较低的阶段有着其不可替代的优势。因此地质统计学作为一种新的资源储量估算与评价方法,要被地质矿产研究人员普遍的接受还需要进一步的过程。
(2)没有形成一套完整的地质统计学矿体储量估算的应用方案。地质统计学方法在我国地质矿产储量评价应用中难以推广一个很重要的原因是由于地质研究人员对克立金估值方法的选择及模型构建中参数的确定存在疑问,缺少相关的参考准则,参数的选取无经验可依,不同参数设置、估值方法和约束边界设置所获得估计结果往往相差很大。虽然目前已有很多学者从多个角度对地质统计学应用中的影响估值精度的因素进行过探讨：如Matheron(1984)从理论上讨论了数据分布和样品支撑对估值精度的影响；Journel(1977)和Pitard(1994)详细研究了块金值在操作中的生成原因及理论依据。Andre (1996)探讨了块金值在对偶克立金法估计精度的影响。但是目前任然缺少一套完整的地质统计学矿产资源储量估算的应用方案,特别是应用地质统计学法资源储量估算前如何对原始数据进行预处理保证后期估值精度；如何通过品位的分布信息、矿床的类型选择合适的克立金估值方法；怎样通过结构分析获得满意的理论变异函数模型等问题都有待进一步地讨论。
(3)缺少克立金估算方法与矿床类型及成矿因素之间的适用关系模型；虽然地质统计学在方法与理论上都有了很大的突破,但是目前还没有真正在地质矿产勘查中形成一套通用的地质统计学矿床建模指导方法,缺乏针对具体某一类矿床的地质统计学模型构建分析与研究是主要原因之一；针对具体的某一类型矿床,如何选择合适的克立金方法进行变异函数分析,获得矿体空间品位分布模型与资源储量、不同的估值方法与模型精度之间的关系等问题都还有待进一步研究。
综上所述,地质统计学资源储量评估方法在我国一直无法得到普遍的应用,其主要原因是目前还没有形成一套完整的应用方案供地质矿产勘查人员在进行地质统计学法资源储量估算时进行参考。因此本文以中国地质调查局发展研究中心“资源量估算子系统模块开发2006-2008”、全国危机矿山资源接替找矿项目管理办公室“新技术新方法”项目：“2008年危机矿山勘查项目成果报告编制GIS系统开发课题”、“2009年危机矿山勘查项目成果报告编制GIS系统完善课题”等研究课题为依托,结合地质统计学理论、矿床建模理论和三维地质建模技术,以实际铜矿床地质勘查数据为例,研究斑岩型(文中代号BYCu)、矽卡岩型(文中代号XKCu)、海相火山型(文中代号CJCu)等三类不同成因的铜矿床的地质统计学矿床模型构建过程。从地质勘查数据建库和预处理开始,详细讨论地质统计学法矿产资源储量估算中,如何结合矿床成因和数据分布选择合适的区域化变量,如何对数据进行处理获得稳健的实验变异函数、如何根据变异函数构建合适的拟合模型及如何根据数据统计特征选择合适的克立金估值方法对矿体进行赋值,最后研究总结归纳一套适用于斑岩型、矽卡岩型、海相火山岩型铜矿床的地质统计学模型构建方法及其应用模式。获得的阶段性成果和认识如下：
(1)论文结合地质矿产勘查实际业务流程,在现有矿床建模技术基础上,提出了一套更为完善的、面向数字矿产调查业务处理的矿床建模技术方案。应用于三实验铜矿区原始勘探数据的数据建库及后期的业务处理。具体过程如下：依据勘探编录规范构建地质矿产勘查信息数据库→根据工业指标自动圈定单工程矿体→在勘探剖面图上构建地质体及矿体的边界线→构建矿体、地质体的三维实体模型→矿体及地质体的空间域内进行规则块体的细分,构建出待估值的空块模型→基于地质统计学及其它统计学理论完成区域化变量分析、模型套合及块体赋值；
(2)在进行地质统计学矿床模型构建之前,对矿区数据进行预处理与组合划分。在数据预处理阶段,区域化变量的选择、特高品位处理、不同性质样品处理等三方面的问题需要引起足够的重视。根据对三矿区不同矿化带品位直方统计,获得了三实验数据基本特征：(a)BYCu矿区：主矿体位于矿化中心相,赋存于花岗闪长斑岩与右英闪长斑岩中；(b) CJCu矿区：块状硫化物Ⅰ号铜锌矿体的品位直方分布存在明显的混合分布特征,矿体在Cu品位2%-6%区间内有二次富集的现象；(c) XKCu矿区：三个主要矿化带均呈现相似的品位直方分布,因此在地质统计学模型构建时对其进行整体处理。
(3)取样间距、特高品位、比例效应和混合分布是影响铜矿床数据变异函数稳定性的主要因素。随着取样间距的增大,变异函数的平稳性逐渐提高,但同时也掩盖了矿体中品位的微观变异信息。可以通过统计变异函数曲线值的离散程度来确定曲线的平稳性。而特高品位是使得实验变异函数平稳性下降的最主要因素之一。因此建议在变异函数计算之前,将特高品位数据按品位上限值进行统一替换或剔除处理。比例效应的存在会使实验变异函数波动性变大,抬高块金值、基台值及估计方差,其中以正比例效应最为常见。对XKCu矿区矿化数据按内接触带、接触带、围岩带、似层状含铜黄铁矿等四类数据分别计算方位角10。,倾角35。方向上的实验变异函数,对比发现内接触带、围岩带、接触带存在较为明显的正比例效应。因此对XKCu矿区数据进行对数变换,获得较为稳定的实验变异函数曲线,有利于模型后期的拟合。
(4)对确立的三矿区实验变异函数进行了地质解释与分析,并在此基础上利用理论模型对其进行拟合。获得几下主要认识：(a) BYCu矿区矿化中心相矿体的品位变化呈现几何各向异性的特征,块金值为0.032，局部基台值为0.246。主轴方向,即矿化相变接触带方向的矿体,在滞后距为180m之内的矿化连续性较好,之后则呈现一定程度的空穴效应,推断是由于矿体与夹石层交替出现造成。由于成矿岩体受一背形穹窿构造控制,导致矿体在近水平方向(次轴方向)的变异性增大,稳定性较其它两方向稍差。而垂直轴方向代表了含矿岩体侵入方向上矿体品位变异性,变异函数在滞后距为280m处开始突变,因此推断矿体在该方向上的平均厚度在280m左右；(b) CJCu矿区Id矿化带的品位变化呈现带状各向异性的特征,垂直轴方向上的基台值高于主轴与垂直轴,且变程各不同。矿体在倾向方向上有较为明显的矿物相分带从而导致变异函数在该方向上自滞后距50m之后呈现明显空穴效应。同时由于混合分布的存在,导致变异函数自滞后距100m之后曲线出现上翘的情况。Id矿化带在矿体倾向与走向方向具有各向同性的特征。而在矿体的厚度方向。拟合模型的局部基台值为3.442，高于主轴与次轴的基台值。根据拟合曲线的变程推断矿体平均厚度在85m左右；(c)XKCu矿区矿化数据在三方向上呈现几何各向异性。碳酸盐岩体与中酸性斑岩体的接触带上,沿该方向主要发育侵染状、脉状铜矿体或透镜体。根据主轴变程推断,矿体的平均长度在232m左右。同时沿碳酸岩体倾向(次轴)与走向(垂直轴)方向矿体的矿化连续性较好,具有相同的变异程度；
(5)进行克立金估值时,搜索椭球的设置对估值精度有一定的影响。搜索椭球半径的确定也是估值样品选择的主要参数,搜索椭球半径的设置一般与变异函数三方向上的变程一致。但是当拟合模型与实验变异函数拟合程度不高时,建议尽量使理论模型与实验变异曲线的前几个点吻合或相近,同时缩小搜索椭球半径以保证搜索范围内的样品品位变异性与拟合模型一致。另外在样品搜索时需要考虑矿体空间各方向上的变异性,尤其是前几个滞后距内数据的变异长度。如果数据在前几个滞后距内的变异性很大。呈带状各向异性则建议通过设置容差范围只搜索理论模型三方向上的样品进行估值；反之,则可以使用八分圆方式进行样品搜索。
(6)利用XKCu矿区数据对不同数据方差、品位上限值及估值方法下获得的估值精度进行对比实验,建立了待估数据的方差与克立金估值方法选择的对应关系,应用于铜矿床地质统计学品位模型构建：当数据即不服从正态分布也不服从对数正态分布时,可以以数据方差和品位上限值为依据选择估值方法。当待估数据的方差大于1.0时,如果数据品位上限值大于数据均值+9倍标准差,则建议应用指示克立金法；如果数据品位上限值小于均值+9倍标准差,则建议应用普通克立金；当待估数据的方差介于0.5-1.0时,建议使用普通克立金；当待估数据的方差小于0.5且品位上限值小于均值+5倍标准差时,建议不剔除原始数据中的特高品位直接进行计算,普通克立金与对数克立金获得的估值精度基本一致；如果品位上限值大于均值+5倍标准差,则建议进行特高品位处理；若要使用对数正态克立金进行估值,如果待估数据方差大于0.5,则建议将待估数据的品位上限制设置为均值+2倍标准差左右进行估值；将上述结论应用到BYCu、CJCu矿区进行验证实验,初步证明了上述结论在铜矿床地质统计学建模与资源储量估算中的普遍性。
Kriging estimation techniques, from the design phase for the mining of ore block grade estimate accurately. The problem was first proposed by the South African mining engineer DG Krige. It made that according to inside and outside of the ore block to be estimated in different blocks of ore samples, the average grade of ore optimal fitting with the estimated. French mathematician Matheron (1963) gave Krige's theory a further systematic research, so that it found the formation of geostatistics theory. After several decades of research and development, geostatistics has been gradually formed a relatively complete set of theories and methods of work used in mining, geology, geophysics, geochemistry, hydrogeology, environmental geology, meteorology and agriculture and other fields. At present the international mining industry has build Geostatistical Kriging estimation method as the standard method of solid mineral resources estimation.
Geostatistics, the introduction of China's mining industry and research for almost 30 years, but it was reported that currently the estimation method submitted by geostatistics solid mineral geological survey/exploration resource reserves are only a dozen copies of, it can be said that geostatistical resource/reserves estimation and evaluation methods have not yet widely application in China geology and mineral exploration and mining industry. But with the global mining development, applying geostatistical evaluation methods in geological mineral resources estimation has become a trend. The current foreign mining market only confirm the mineral resource/reserve reports that evaluated by kriging estimation method estimates, so that if do not use geostatistics method, geological exploration results and reserves of assets can not be recognized by international peers, which for our international mining cooperation is very negative. However, geostatistics reserves estimation method in geological exploration and mining industry has been unable to obtain a wide range of applications, the main reasons are as following:
(1) The traditional method is still the primary means of estimating reserves in China geology mineral resources, geostatistics should be general accepted by geology mineral exploration researchers also need a process; The geological block segment, profiling method, etc. (brief "traditional method") reserves estimation method of solid mineral resources are still primary methords for reserves estimation. Traditional methods due to strong operational, the application simple and easy calculation, after nearly 50 years of application has formed a relatively completive system of practices and theories. Furthermore, traditional methods in pre-investigation, survey and other exploration and other extent of the lower stage has its unique advantages. Therefore, geostatistics as a new resource/reserves estimation and evaluation methods be generally accepted by geology and mineral resources need to be further process.
(2) Have not form a complete set of geostatistical mineral reserve estimation applications. A very important reason why geostatistics of geology mineral reserves is so difficult to popularized in China's geology mineral reserves estimation is that:①geology researchers are very unsure about the Kriging method choice;②The confirms of parameters in model construction are really doubt. Have no standard refence criteria;③The selection of parameters are inexperienced to follow; the estimation results made by a different parameter setting, valuation methods or constraints boundary set are often have a big difference. Although there are many scholars from various angles to discuss the factors that affect accuracy in the application of geostatistics. Discussed from both a theoretical support for data distribution and sample the impact on the valuation accuracy. Journel and Pitard detailed study the reasons that creates the nugget values in the operation and its the theoretical basis. Andre (1996) discuss the nugget on dual kriging estimation accuracy.However, there's still has no a completly set of application preject for geostatistical mineral reserve estimation. Particularly in original data pretreatment to makesure the estimation precision, how to choose a appropriative kriging evaluation method by grade distribution and deposit type, how to obtain theoretics variance by construct analysis are still need to further discuss.
(3) Have not a detailed research on application mode of kriging estimation methods in a specify types of deposits; Although the methods and theory of geostatistics has a great breakthrough, there is not a common set of geological deposits statistical modeling guidance methods in the geology mineral exploration and absent a certain category-specific statistical model of the geological deposits to analysis and research; Due to certain specific types of deposits, how to choose a suitable kriging method to analysis variogram, spatial grade distribution model to obtain ore reserves and resources, different valuation methods and model accuracy of the relationship between the formation of deposits for a certain type of deposit types and mineralization Geostatistical between the factors relational model for such problems are yet to be further studied.
The geostatistical resource reserves assessment method has been unable to be widely applied in China, mainly because there is no formation of a complete set of applications for the geology mneral resources researchers to reference when estimating reserves. Therefore, this article supported by China Geological Survey Development Research Center "Solid mineral resource/ resverse estimation subsystem module development from 2006 to 2008," national crisis prospecting mining resources project management office:"New technologies and new methods" subject:" GIS system development on crisis in mine exploration project report on the outcome of the preparation 2008 ""the improve subject on crisis of mining exploration project in reporting the outcome by GIS system 2009 ", etc combined with geostatistics theory, deposit modeling theory and three-dimensional geological modeling techniques used the actual copper deposit geology survey data, to research geostatistical model construction process on porphyry-type (code-named BYCu), skarn-type (code-named XKCu), marine volcanic-type (code-named CJCu) three different copper deposit model bed. From the data base construction process to discuss the geostatistical method mineral resources estimates in detail. Mainly foucs on how to select the appropriate regional variable by data distribution and ore genesis, how to construct adpetable theory variance fuction by experimental variogram fuction and how to choose the right Kriging valuation method of ore-body estimation by data characteristics, and finally summarizes an applicable geological modeling construction method on porphyry-type, skarn type, marine volcanic copper deposits and their application mode. Obtained results and knowledge are as follows:
(1) This research conbined with the actual business processes of mineral exploration, based on existing deposit modeling techniques, gave a more sophisticated mineral deposit modeling technology programs for the digital geology survey processing. And Apply it to three exploration data of copper deposit the data to build the original database and the latter part of business processing. The specific process is as follows: Based on the exploration and catalog standard bilding geology mineral resources information database→Automatic delineate orebodies base on Industry index→profile the boundary line of geological body and ore body in the exploration section maps→build three-dimensional solid body of orebody and geological model→to create block model with the orebody solid model restrict→based on Geostatistics and other statistics to complete the regionalized variable theory, analysis, model fitting and block assignment;
(2) Before geostatistics deposit modeling construct, the data pre-processing and sample conbined are very important. In the data preprocessing stage, the choice of regionalized variables, special high-grade processing and different types of sample handling need to pay enough attention. Based on the three different mineralized mining belt grade straight statistics, obtained basic characteristics of three experimental data:(a) BYCu mining areas:the main ore mineralization in the central phase, occur in the granodiorite porphyry and quartz diorite porphyry; (b) CJCu mine: The grade straight-square of I massive sulphide copper-zinc orebody distribution exists an obvious mixture distribution characteristics. The Cu ore grade within the range from 2% to 6% has secondary enrichment of the phenomenon; (c) XKCu mining: the three main mineralized zone showed a similar grade of straight-square distribution, and therefore modeling, geostatistics when dealing with them as a whole set.
(3) sample distance, ultra-high grade, the proportion effect, and data with mixed distribution are the main factors that impact the stability of variogram in copper deposit data main. As the sampling distance grow and the smooth variation of the function has a gradual increase, but it also conceals the micro-grade ore body variation information. The the value of geostatistical variation can be a refence to judge the smoothness of the discrete function curve. The ultra-high grade is the most important factor to make experimental variogram steady declining. It's important to handle the ultra-grade before experimental variogram calculation, usually the ultra-grade can be replaced or delete by an setting-max grade. the proportion effect can also impact the experimental variogram to a fluctuate, raising the nugget, sill and estimate variance of which the most common effects of direct proportion.XKCu mining data are classified within four types of data as inner-contact zone mineralization, contact mineralization, wall rock mineralization and similar-layered copper pyrite, And than calculated the experimental variogram by azimuth angle 10°and dip 35°.It denotes that inner-contact zone mineralization, contact mineralization, wall rock mineralization have an obvious effects of direct proportion. Therefore process the logarithmic transformation to mining data to obtain a more stable experimental variogram curve so that it's conducive to the later model fitting.
(4) Established three experimental variogram mining carried out geological interpretation and analysis, and on this basis, the use of theoretical models fit them. Access a few key understanding: (a) BYCu mine orebody grade mineralization centers relative to a characteristics of geometric anisotropy with nugget value of 0.032, sill of 0.246. Principal axis direction, which is the direction of mineralized phase change in contact zone of ore body, within a lag distance of 180m manifest a continuative mineralization, and then become a degree of cavitation effect is due to conclude with alternation between stone layers and the ore body. As the ore-forming rock back by a dome-shaped structure control, resulting in nearly horizontal ore body (sub-axis direction) of the variability increased somewhat less stable than the other two directions. The direction of the vertical axis represents the direction of the ore-bearing intrusive rock ore grade variability, variation function in the lag distance of 280m at the beginning mutations, so inferred that the average thickness of ore body in this direction is 280m; (b) CJCu mine Id grade mineralized zone showed changes in the characteristics of band anisotropy, the vertical axis of the base station is higher than the spindle and vertical axis, and also have the different variable-range. Ore body in the direction of incline has a more obvious tendency of the mineral phase sub-zone and thus leading to mutations function in the direction away from the lagged effect of holes became apparent, after 50m. At the same time due to the presence of mixed distribution, leading to mutations function of distance from the hysteresis curve after 100m upturned situation. Id mineralized zone in the direction of orientation and the direction it has isotropic characteristics. In the ore body thickness direction. Fitting the model of local base station values 3.442, higher than the spindle and sub-axis value of the base station. Fitting curve according to variable-range inferred average thickness of ore bodies in the 85m or so; (c) XKCu three mineralized mining data is presented to geometric anisotropy. The infection along the direction of the main shape development vein copper ore body or lens on carbonate rock and the acidic porphyry contact zone, According to variable-range axis inferred ore body, the average length of about 232m. At the same time along the carbonate rock tendencies (sub-axis) and direction (vertical axis) direction of the mineralized ore body continuity is better, have the same degree of variation;
(5) when it valuate by kriging, the search ellipsoid on the settings have a certain impact to the accuracy of valuation. Search ellipsoid is also a valuation to determine the radius of the main parameters of the sample selection, the search ellipsoid with the variation of the radius of the general function of the tripartite setting of the variable-range line up. However, when fitting model with the experimental variogram fitting degree is not high, it is recommended as far as possible theoretical models and experimental variation of the curve fit the first few points, or similar, while narrow your search to ensure the search ellipsoid radius within the scope of variability of the sample grade consistent with the fitting model. At the same time need to be considered when searching in the sample space for parties up ore variability, particularly within the first few lags from the variation of the length of the data. If the data from the last few lag variability within a large band anisotropy was suggested by setting the tolerance extends only to search theory model of the tripartite upward valuation of the samples; the other hand, you can use the eight round way sample search.
(6) Taking XKCu data as exsample, the author gave a series of comparative experiments to compare the estimation precision under different data variance, upper limit of grade and estimation methods. When the data that is not subject to the normal distribution does not obey log-normal distribution in variance, grade variance can be the upper limit of the method of valuation is based on choice. When the variance of the data to be estimated is greater than 1.0, if the data is greater than the upper limit of quality data on the mean +9 times the standard deviation, it is recommended that application of the Indicator Kriging; if the data quality +9 times the upper limit of less than the mean standard deviation, then the the proposal of using Ordinary Kriging; when the variance of the data to be estimated between 0.5-1.0, it is recommended the use of Ordinary Kriging; when the variance of the data to be estimated is less than 0.5, if the grade is less than the mean+5 times the upper limit of the standard deviation, it is recommended not remove the original data on the ultra-high grade with directing calculation, with the number of Ordinary Kriging was basically the same as the valuation of accuracy; if grade is greater than the mean+5 times the upper limit of the standard deviation, it is recommended to conduct special high-grade treatment; to the use of lognormal kriging valuation, if the data variance to be estimated is greater than 0.5, it is recommended to be assessed on the quality of data limitations is set to +2 times the standard deviation about the mean valuation; Applied above conclusions to BYCu and CJCu data, generally verified the universality of those conclusions in the copper deposit modeling and geostatistics reserves estimates.
引文
[1]Journel,A.G. Kriging in terms of projections[J]. Mathematical Geology,1977,9(6):563-586.
[2]Krige,D.G. A statistical approach to some mine valuation and allied problems on the watersrand, MSc(Eng)thesis, University of the Witwatersrand, Johannesburg, South Africa,1951.
[3]Matheron G. Principles of geostatistics [J]. Economic Geology,1963,58:1246-1266.
[4]孙玉建.地质统计学在固体矿产资源评价中的若干问题研究[D].中国地质大学(北京)博士学位论文,2008.
[5]孙玉建.以地质统计学为基础的矿业软件在中国的历史和现状[J].中国矿业,2007,16(123)：79-82.
[6]国土资源部矿产资源储量司.矿产资源储量计算方法汇编[M]匕京：地质出版社,1991.
[7]赵鹏大,矿产勘查理论与方法[M].北京：中国地质大学出版社,2001.
[8]Matheron,G.,The selectivity of the distributions and"the second principle of Geostatistics",in GVerly et al.(editors),Geostatistics for natural resources characterization,by D.Reidel Publishing Company,part Ⅰ,1984,p.421-433.
[9]Pitard,F.F. Exploration of the"nugget effect",in Geostatistics for the next century, edited by Roussos Dimirakopoulos,Kluwer Academic Publishers,the Netherland,1994.
[10]Andre, Trochu B.,F., Dansereau J. Approach for the smoothing of three dimensional reconstructions of the human spine using dual Kriging interpolation[J]. Medical & Biological Engineering & Computing,1996,34:185-191.
[11]王瑞江,王义天,王高尚,等.世界矿产勘查态势分析[J].地质通报,2008,27(1)：154-162
[12]屠雯.矿业复苏可待中国势不可挡[J].中国有色金属,2009，23：32-35.
[13]刘树臣,王淑玲,崔荣国.全球矿产资源供需形势及勘查动向[J].地质通报,2009，28(2-3)：160-165.
[14]周进生,罗晓玲.近十年全球矿产勘查趋势分析[J].中国国土资源经济,2008,9：13-15.
[15]2009年前三季度我国固体矿产勘查情况[J].现代矿业,2009，487(11)：5.
[16]陈喜峰,彭润民.加强隐伏矿勘查的必要性和意义[J].矿业研究与开发,2007，27(4)：87-89.
[17]龚国清,刘修国,倪平泽等.基于iExploration-EM的数字矿产勘查成果编制方法研究[J].金属矿山,2009,23(1)：34-38.
[18]李超岭,杨东来,李丰丹,等.中国数字地质调查系统的基本构架及其核心技术的实现[J]. 地质通报,2008,27(7)：923-944.
[19]王莉,武法东.GIS在矿产资源勘查开发的应用现状和前景展望[J].中国矿业,2008，17(7)：18-29.
[20]叶天竺,薛建玲.金属矿床深部找矿中的地质研究[J].中国地质,2007，34(5)：855-569.
[21]黄竞先,侯景儒.泛克立格法和指示克立格法在地球化学探矿中的应用[J].地球科学中国地质大学学报,1994，19(3)：321-328.
[22]肖斌,赵鹏大,侯景儒等.山东归来庄金矿床金异常分布及其时空演化的地质统计学研究[J].现代地质,1991,13(4)：419-424.
[23]肖斌,潘懋,赵鹏大等.山东归来庄金矿区g(Au/Ag)的空间信息统计学特征[J].地质科学,1991,36(4)：500-508.
[24]燕永锋.“矿床数学-经济模型”软件包的开发及其在元江金矿的应用[D].昆明理工大学博士学位论文,2000.
[25]黎应书,秦德先,蔡燕等.云南省大红山铜矿床三维矿床数学模型[J].有色金属(矿山部分),2005,57(5)：20-22.
[26]邓明国,李文昌,秦德先等.克立格方法在个旧矿区芦塘坝10-9号矿体储量计算中的应用[J].地质与勘探,2006,42(6)：67-70.
[27]余海军,李文昌,尹光侯等.普朗铜矿床铜品位分布地质统计学研究[J].地质与勘探,2009,45(4)：437-443
[28]潘国成,李秀峰,译,地质统计学储量估算方法在矿业中应用的若干实际问题(一)[J].国外金属矿山,1996,7：7-12.
[29]潘国成,李秀峰,译,地质统计学储量估算方法在矿业中应用的若干实际问题(二)[J].国外金属矿山,1996,8：5-9.
[30]黄竞先,李春霞,矿业地质统计学在我国的应用现状及对影响克立格估计的若干因素的研究[J].有色金属矿产与勘查,1997，6：38-43(增刊)
[31]刘永社,印兴耀,贺维胜,2004,空间相关分析因素对储层建模中克立金估计结果的影响[J].石油大学学报(自然科学版),282(24)：27-37.
[32]新华网,http://www.lm.cn/zjtg/societyDiscussion/200801/t20080124_192277.htm
[33]芮宗瑶,王龙生.中国铜矿床分类新方案[J].有色金属矿产与勘查,1994，3(2)：96-97.
[34]王全明.我国铜矿勘查程度及资源潜力预测[D].中国地质大学(北京)：博士学位论文,2005.
[35]袁见齐,朱上庆,翟裕生.矿床学[M].北京：地质出版社.1985.
[36]夏斌,涂光炽,陈根文,喻享祥.超大型斑岩铜矿床形成的全球地质背景[J].矿物岩石地球化学通报,2000,19(4)：406-408.
[37]芮宗瑶,张洪涛,陈仁义等.斑岩铜矿研究中若干问题探讨[J].矿床地质,2006，25(4)：491-500.
[38]黄崇轲等,中国铜矿床,北京：地质出版社,2001.
[39]陈毓川、朱裕生等,中国矿床成矿模式[M].北京：地质出版社,1993。
[40]Journel A G, Huijbregts C J, Mining Geostatistics[M]. Academic Press, New York,1978.
[41]Isaaks,E.H.,and R.M.Srivastava. An introduction to applied geostatistics[M]. Oxford Univ. Press,New York,1989,561.
[42]Cressie,N. Statistics for Spatial Data[M]. Revised Edition.Wiley,New York,1993.
[43]Goovaerts,P. Geostatistics for Natural Resources Evaluation[M]. Oxford University Press, 1997.
[44]侯景儒,黄竞先.地质统计学及其在矿产储量计算中的应用[M].北京：地质出版社,,1981.
[45]侯景儒,黄竞先,等译.矿业地质统计学[M].北京：冶金工业出版社,1982.
[46]王仁铎、胡光道.线性地质统计学[M].北京：地质出版社,1988.
[47]侯景儒,尹镇南,李维明.实用地质统计学[M].北京：地质出版社,1998.
[48]侯景儒,黄竞先.地质统计学在固体矿产资源/储量分类中的应用[J].地质与勘探,2001,37(6)：61-66.
[49]侯景儒,古梅.矿石储量分级的地质统计学研究[J].地质与勘探,1987，23(3)
[50]侯景儒,黄竞先.对数正态克立格法理论及其应用[J].北京科技大学学报,1989,11(5)：391--398.
[51]王志民,侯景儒.协同克里格法及其在矿产储量计算中的应用[J].地质与勘探,1994，30(3)：39-48.
[52]肖斌,赵鹏大,侯景儒.时空域中的指示克立格理论研究[J].地质与勘探,1999，35(4)：25-28.
[53]肖斌,赵鹏大,侯景儒.时空泛克立格法的理论研究[J].江汉石油学院学报,2001,23(1)：81-87.
[54]贾明涛,潘长良,王李管.客服地质统计学矿床中的主观因素影响技术研究[J].地质与勘探,2003,39(4)：73-77.
[55]孙玉建,孟伟,万会.矿产资源储量估算中工程控制程度划分的探索[J].地质与勘探,2006,42(6)：81-84.
[56]http://www.microminechina.com/.
[57]http://www.surpac.com.cn/.
[58]侯景儒.中国地质统计学(空间信息统计学)发展的回顾及前景[J].地质与勘探,1996，32(1)：20-25.
[59]陈爱兵,秦德先，张学书,等.基于MICROMINE矿床三维立体模型的应用.地质与勘探,2004,40(5)：77-80.
[60]陈国旭.支持传统方法的矿产资源储量估算软件研究进展[J].国土资源信息化,2009，2：8-12.
[61]中国地质调查局数字地质调查技术支持网站WWW.DGST.CGS.GOV.CN
[62]李超岭,杨东来,李丰丹,等,2008.中国数字地质调查系统的基本构架及其核心技术的实现.地质通报,27(7)：923-944.
[63]吴健生,朱谷昌,曾新平,等,2004.三维GIS技术在固体矿产勘探和开发中的应用.地质与勘探,40(1)：68-72.
[64]陈爱兵,秦德先,张学书,等,2004.基于MICROMINE矿床三维立体模型的应用.地质与勘探,40(5)：77-80.
[65]黄艳丽,秦德先,李连伟,2007.Surpac Vision及其在数字矿山中的应用.矿业快报,457：27-30.
[66]刘修国,陈国良,侯卫生,等,2006.基于线框架模型的三维复杂地质体建模方法.地球科学--中国地质大学学报,31(5)：668-672.
[67]Ekoule A B, Peyrin F C, Odet C L,1991.A triangulation algotrithm from arbitrary shaped multiple planar contours. ACM Transactions on Graphics,10 (2):182-199.
[68]Meyers D, Skinner S, Sloan K.1992. Surfaces from Contours[J]. ACM Transactions on Graphics,11(3):228-258.
[69]国家地质矿产行业标准煤、泥炭地质勘查规范(DZ/T0215-2002)
[70]固体矿产勘查原始地质编录规程(DD 2006-01)
[71]固体矿产勘查地质资料综合整理研究规定(DZ/T0079-1993)
[72]地质行业标准DZ/T0124-2002《铜、铅、锌、银、镍、钼矿地质勘查规范》
[73]徐兵,马小刚,田宜平.递归算法在单工程矿体边界圈定中的应用[J].有色金属(矿山部分),2006,58(6)：21-24.
[74]陈国旭.垂直剖面法矿产资源储量估算软件实现技术及自动分类方法研究[C].中国地质大学(武汉)硕士学位论文,2008.
[75]中华人民共和国国家技术监督局,GBl3687-92,固体矿产普查总则(第一分册)[M].北京：中国标准出版社。
[76]中华人民共和国国家技术监督局,GBl3688-92,固体矿产详查总则(第二分册)[M].北京：中国标准出版社。
[77]罗周全,刘晓明,苏家红等.基于Surpac的矿床三维模型构建[J].金属矿山,2006.358(4)：33-36.
[78]徐静,胡乃联.Surpac软件在某金铜矿山的应用研究[J].黄金科学技：术,2007.15(1)：54-58.
[79]Boisvert J B, Leuangthong O, Oritiz J M, Decutsch C V. A methodology to construct training images for vein-type deposits[J]. Computers & Geosciences,2008,34: 491-502.
[80]Bastante F G, Taboada J, Alejano C, etc. Evaluation of the resources of a slate deposit using indicator kriging[J].Engineering Geology,2005,81:407-418.
[81]史文中,吴立新,李清泉等,2007.三维空间信息系统模型与算法.北京：电子工业出版社.
[82]张宝一,尚建嘎,吴鸿敏,等,2007.三维地质建模及可视化技术在固体矿产储量估算中的应用.地质与勘探,43(2)：76-81.
[83]Meyers D, Skinner S, Sloan K. Surfaces from Contours[J]. ACM Transactions on Graphics, 1992, Vol.11(3):228-258.
[84]Keppel,1975. Approximations of complex surfaces by triangulation of contour lines. IBM Journal of Research and Development,19:2-11.
[85]Christiansen H N, Sederberg TW,1978.Conversion of complex contour line definitions into polygonal element mosaics. Computer Graphics,12:187-192.
[86]Ganapathy S, Dennephy T G,1982. A new general triangulation method for planar contours. Computer Graphics,16(3):69-75.
[87]李清泉,杨必胜,史文中,等.三维空间数据模型的实时获取、建模与可视化[M].武汉：武汉大学出版社,2003.
[88]陈慧新.变异函数的稳健性及特异值的处理方法[J].内蒙古农业大学学报,2000,21(2)：84-90.
[89]Clow,G. Why gold mines fail,northen miner[J]. February,1991,31-34.
[90]Rendu,J.M. How the geologist can prevent a geostatistical study from running out of control some suggestions,in Applied mining geology:ore reserve estimation,edited by D.E. Ranta[J]. Society of Mining Engineers of AIME,,1986,11-17.
[91]Diamond,P.and M.Armstrong, Robustness of variograms and conditioning of kriging rices[J]. mathematical geology,1984,16(8):809-822.
[92]Armstrong M. Improving the estimation and modeling of the variogram[M]. Geostatistics for natural resources characterization:D.Reidel, Amsterdam,1984,1-19.
[93]Bardossy A. Notes on the robustness of the kriging system[J]. mathematical geology,1988, 20(3):189-203.
[94]Brooker P I. A parametric study of robustness of kriging variance as a function of range and relative nugget effect for a spherical semivariogram [J]. mathematical geology,1986,18(5): 477-488.
[95]赵鹏大,矿产勘查的若干重要思路及途径[J].矿产储量管理,1995,Vol.5：4-10.
[96]赵鹏大,“三联式”资源定量预测与评价[J].地球科学,2002。
[97]赵鹏大,2002.“三联式”资源定量预测与评价-数字找矿理论与实践探讨.地球科学-中国地质大学学报,27(5)：482-489.
[98]Krige,D.G. A statistical approach to some mine valuation and allied problems on the watersrand, MSc(Eng)thesis, University of the Witwatersrand, Johannesburg, South Africa,1951.
[99]Krige,D.G. Lognormal-de Wijsian geostatistics for ore evaluation, South African Institute of Mining and Metallurgy, Monograph Series,1978,40.
[100]Krige,D.G.,1981,Lognormal-de Wijsian geostatistics for ore evaluation,South African Institute of Mining and Metallurgy,Monograph Series,51.
[101]Krige,D.G, Magri,E.J., Geostatistical case studies of the advantages of lognormal-De Wijsian kriging with mean for a base metal mine and a gold mine, Mathematical geology, 1982,14(6):547-555.
[102]Matheron G.Traitede geostatistique applique,tome 1,Memoires de Bureau de Recherches Geologiques et Minijeres, No.14,Editions Technip, Paris, English translation:Treatise, on applied geolstatistics,Kennecott Copper Corporation,Salt Lake City,Utah,1960.
[103]Matheron G. Traitede Geostatistique Appliquee, Editions Technip, Paris,1962, vols(1,2): 333.