# 3D geological modeling for mineral resource assessment of the Tongshan Cu deposit, Heilongjiang Province, China

### Introduction

Three-dimensional geological modeling (3DGM) is a developing technology for geological studies, mineral resource exploration, and quantitative estimation of mineral resources (Houlding, 1994; Mallet, 2002; Wu and Xu, 2004; Fallara et al., 2006; Sprague et al., 2006; Pan et al., 2007; Wang et al., 2007; Calcagno et al., 2008; Kaufmann and Martin, 2008; Pouliot et al., 2008; Wang et al., 2009a; Xiao, 2009; Zanchi et al., 2009). Nowadays, 3DGM can integrate two-dimensional (2D) GIS, database, statistics methods, and 3D visual technology (e.g., Micromine, GoCAD, and Surpac software, as well as VC++ and OpenGL) developing software with three-dimensional visualization modules (Wang et al., 2009b). 3DGM allows modeling of complex and irregular geological objects in a 3D environment using geological maps, geologic survey records, structural information, geophysical and geochemical data (Mallet, 2002; Zanchi et al., 2009). Therefore, 3DGM can represent sophisticated 3D geological objects in three spatial dimensions as a continuous function of the geological coordinates: V = f(x, y, z). 3DGM can be applied to construct 3D structural models of geological objects involving stratum, structure, rock body, geophysical anomaly, geochemical anomaly, and ore body (Kaufmann and Martin, 2008; Xiao, 2009). 3DGM can also be used to illustrate the genesis of geological objects through spatial analysis or logical calculation with respect to regional metallogenic theory and typical ore deposit models (Pouliot et al., 2008).

Researches on the Tongshan Cu deposit in Heilongjiang Province, China have utilized geological, geochemical, geophysical, and remote sensing methods (Du et al., 1988; Zhao et al., 1995; Han et al., 2007; Wang et al., 2007; Cui et al., 2008) (Fig. 1). However, two unresolved geological issues are (1) the relationship between orebody III and deep magmatic-hydrothermal mineralization, and (2) the effect of the Tongshan Fault on orebody II of the deposit (Fig. 1). In this paper, 3DGM and mathematical methods including BP (Back Propagation) network, virtual borehole, virtual section data, and an improved interpolation method are used to construct 3D models of the Tongshan Cu deposit (including mineralization stratum, granodiorite mineralization, altered rock, magnetic anomaly, and orebody).

### Geological setting

The Tongshan Cu deposit is situated in the Daxinganling fold system of Inner Mongolia and is located within 4 km of the well-known Duobaoshan Cu deposit. The 3D study area is 4248400–4249000 N, 3764000–3767000 E, at 620–1000 m elevation.

The study area consists of up to 6000 m of Ordovician marine volcanic (basalt, andesite, rhyolite)-sedimentary rocks formed in an active continental margin, island arc environment (Du et al., 1988, Cui et al., 2008). The Tongshan Cu deposit occurs within the southwestern limb of an inverted anticline. The core of the inverted anticline consists of Middle Ordovician rocks of the Tongshan Formation, whereas the limbs comprise andesite and andesitic tuff of the Middle Ordovician Duobaoshan Formation. Volcanic rocks of the Duobaoshan Formation host the main stratum of Cu mineralization (106–132 ppm Cu), over a thickness of some 3000 m (Du et al., 1988). Early Ordovician granodiorite is considered to be the Cu-ore source rock. SHRIMP U-Pb dating of zircon in the granodiorite yield an age of 479.15 Ma, and an Re-Os isotopic isochron age of molybdenite is 506 ± 14 Ma (Cui et al., 2008). The granodiorite is associated with tonalite with a marginal facies of quartz diorite porphyry developed along the footwall of Tongshan Fault. The mineralized Duobaoshan Formation exhibits strong hydrothermal alteration, similar to typical porphyry-type alteration zoning with intense potassium, silica, chlorite, and propylite (locally developed illite-carbonate) zonation from the center outward. On average, the amount of Cu in the alteration rock is more than 0.1%. The ore body is hosted mainly within zones of potassium, silica and phyllic alteration of the quartz diorite porphyry. The Tongshan Cu deposit is a typical magmatic-hydrothermal deposit and has a similar geological setting, alteration sequence, and mineralization characteristics as the Duobaoshan Cu deposit.

There are three main orebodies (I–III) in study area, with a total Cu reserve of the Tongshan Cu deposit of 1,007,992 t (Tongshan detailed exploration report in 1993, the second geological survey of Heilongjiang Province, China). Orebody II is the longest orebody (2.2 km), with a width of 25–713 m. Orebody I occurs at the highest elevation (530 m), whereas orebody III is concealed and in contact with granodiorite at a depth of 862.5 m. The Tongshan Fault is a 30–40° dipping reverse fault, strikes nearly EW for more than 10 km and is >10 m wide (Fig. 1). This fault transects and displaces orebodies II and III and therefore postdates the formation of the Tongshan Cu deposit.

### 3D modeling

In this paper, 3D modeling includes 3D geological object modeling, 3D modeling of orebody grade using a 3D kriging interpolation method, and a 3D modeling of mineralization (i.e., Cu-grade is not less than 0.1%) using an improved IDW interpolation method. The methodology is based on 3DGM and mathematical methods and includes several steps of processing data depending on the type of data. Some steps require establishing a virtual section based on geological knowledge in order to control the boundary of a geological object, or establishing virtual boreholes by combining statistics with a BP network for 3D interpolation calculation (Wang et al., 2009a,b). Other steps such as data interpretation, information extraction, or 3D object model validation require user interaction.

To construct an accurate 3D geological model from geological data (e.g., geological maps at different scales, cross-sections and boreholes), it was necessary to develop a methodology that also takes into account magnetic data. Geological maps synthesize geological information but they do not give a complete representation of the subsurface geology. Cross-sections and borehole logs add the third dimension to give a more detailed interpretation of subsurface structure. However, if geophysical information is available, a better constraint on the interpretation of structures or intrusive rocks is possible. By combining contact locations and orientation of geology and geophysical information, 3D geological models can be constructed.

#### 3D geological objects modeling

The methodology that was followed for 3D geological object modeling involved the following steps:

(1)Geoscience data acquisition and compilation: The study area data include a 1:2000 scale geological map, 21 exploration sections, e.g., 1080 exploration section (Figure 1, Figure 3), 94 borehole logs, a 1:5000 scale magnetic survey data, 23,800 Cu-grade surveys of boreholes (Fig. 2), and ETM+ images. Geological, geophysical, topographical and geochemical data are quantified and standardized in the same 3D coordinate system.

(2)Geoscience information interpretation and extraction: This is efficient and effective for deriving information from the amount of geosciences data on the basis of metallogenic theory and an ore deposit model, e.g., many complex borehole logs can be interpreted on the basis of typical porphyry alteration zoning.

(3)Virtual section and virtual borehole construction based on geological information (e.g., all sections are not completely accordant in length or depth, and some boundaries of deep geological objects) are effectively controlled from sparse borehole data. Orientation data of surface and subsurface rocks can be applied to infer the shape of a deep-seated object, and contact relationships of various surfaces or subsurface geological structures in association with geophysical inversion information are usually applied to infer the boundaries of deep-seated objects. The BP network and statistic methods (e.g., nonlinear simulation with least squares) can be applied to simulate and infer vertical or equal Cu-grade variation at some point on the basis of continuous Cu-grade surveys from boreholes. In this paper, the greatest depth of Cu mineralization is 1732 m using the BP network method based on deepest borehole data (1392.5 m). The depth value of 1732 m was used to define the maximum depth limit studied in this paper.

(4)Construction and integration of regional sections used for the homogenization and simplification of regional geological contacts. Geological surfaces were constructed from contact curves and dip vectors derived from surface geological maps and cross-sections, and topographical data (digital elevation model, DEM). Orientation data were available since contact points and orientation vectors control the geological interfaces, and we polarized the 3D unit vectors of the orientation field of the interfaces along the most recent direction. Geophysical data were used to construct additional 3D models to help constrain depths, dips and boundaries of subsurface geological bodies (orebodies, strata, faults, folds and rock masses).

(5)Validation of the 3D geological model using unconstrained and constrained geophysical inversion, geological map and field survey data according to quantitative relationships between geological and geophysical information; all of the 3D geological object models are therefore related in a 3D environment with same 3D coordinate system.

### 3D orebody modeling

#### Interpolation method

Interpolation is widely used for both predictive and visualization purposes in geoscience studies (Falivene et al., 2010). A variety of algorithms have been developed to construct such interpolations (Morrison, 1974), e.g., inverse distance weighting (IDW) (Kane et al., 1982), kriging (Matheron, 1963), splines (Ahlberg et al., 1967; Mitasova and Mitas, 1993) or polynomial regression. We considered two commonly used interpolation algorithms and the statistical histogram of all Cu sample surveys in the study area (Fig. 2): IDW and ordinary kriging (OK) provide an estimate of the studied variable at every unsampled location, by means of a linear combination of N-observed values of

where is the measured value at the ith location, is an unknown weight for the measured value at the ith location, and N is the number of measured values.

In IDW, depends solely on the distance to the prediction location. However, with OK, the weights are based not only on the distance between the measured points and the prediction location but also on the overall spatial arrangement of the measured points. To use the spatial arrangement of the weights, spatial autocorrelation must be quantified. Thus, in OK, depends on a fitted model to the measured points, the distance to the prediction location, and the spatial relationships among the measured values around the prediction location.

The IDW algorithm for interpolation of 3D scattered data can be written as (Shepard, 1968):

where Zi is the value of the original function, Zx,y,z is the interpolated function value on the x, y, z-grid, represents the effective distance between the nodes, and the interpolant, , , is the distance of the grid node to . It is possible to include three weighting parameters: b is the weighting or power parameter (it is also possible to use b in a manner that a point on a grid node has a beta of 1 and points off the grid have a beta value of zero); is the smoothing parameter.

In this paper, the IDW interpolation method is improved by addition of a constraint weight (with a value of 0.8 and 1) to Zi. If Zi is located in the geological object of mineralization, the constraint weight is 1, if Zi is not located in the geological object of mineralization, the constraint weight is 0.8.

Kriging is a statistical method of interpolation for estimating unknown values from an original sample data set. It differs from other techniques such as IDW in that it uses the concept of spatial continuity between data. The aim is to minimize variances of the estimates and arrive at the best data value for each unknown point. Part of the input is a semi-variogram model that fits the input data. The validity of the kriging result depends on how well the semi-variogram models the data. It is easy to include a distance criterion in the IDW method. A well-established statistical method to estimate this value relies on the use of a semi-variogram (or experimental variogram). A semi-variogram is the graph that is most commonly used in applied Geostatistics to explore spatial interdependence (Kitanidis, 1997). The semi-variogram is calculated by averaging one half of the squared difference of the z-values over all pairs of observations with the specified separation distance and direction. According to Kitanidis (1997), the semi-variogram equation can be written as

where is the number of pairs of measurements, k is the index of the consecutive intervals and index i refers to each pair of measurements and for which

As the range is an estimate of the correlation length of the data, it is an appropriate value for reducing the search distance within the interpolation algorithm.

In a 3D environment, the kriging method regards the 3D cell (e.g., hexahedron, cube) as a point, and the 3D coordinate of the point is generally that of the focus or vertex of a 3D cell. In order to decrease the error from the orebody boundary, the 3D cell is usually divided into several sub units. One practical rule is 4n, where n is the dimension 3 (Jounrel and Huijbergts, 1978).

### 3D interpolation model

On the basis of a mining exploration grid net (100 m × 50 m) of the Tongshan Cu deposit (Fig. 3), combined 3D geological object models using Geostatistics analysis (Fig. 4), the 3D geological cell unit are determined: row × line × layer for a regular block is 20 m × 20 m × 20 m, and sub-block cells may be 20 m × 10 m × 10 m, 20 m × 20 m × 10 m, and 10 m × 10 m × 10 m. The sub-cells are not only helpful for accurate construction of a 3D grade model for the orebody, but also for identification of continuous sensitive differences of Cu-grade in a single orebody. The latter can be applied to infer a causal process of magmatic-hydrothermal fluid relations with reference to a concept of metallogenesis (Figure 4, Figure 5).

On the basis of the improved IDW method, using VC++ and an OpenGL developing environment, the 3D study area of the mineralization model (with a threshold Cu grade of 0.1%) is constructed (Fig. 6). The results show that: (1) Cu mineralization has a northwestern trend but the main orebody lies to the east of the study area; (2) the higher grade orebody is located in the northern and southern parts of the study area; (3) the concealed orebody is located southeast of orebody III at depths between 300 and 500 m.

Comparison of Fig. 1 with Figure 3, Figure 4, Figure 5, Figure 7, shows that the 3D geological model of the Tongshan Cu deposit is accurate. As indicated in Fig. 3, the study area is characterized by a typical porphyry alteration zonation. Orebodies I and II occur in the phyllitic zone above the Tongshan Fault. Orebody III occurs in the potassium, silica (biotite-K-feldspar) alteration zone below the Tongshan Fault. On basis of ore-forming intrusive rock buffer analysis within 500 m (the buffer distance is determined from geochemical data of alteration zoning of the intrusive rock) in the 3D model of the study area, orebody III is the mineralization center of the Tongshan Cu deposit. Ore body II was formed by magmatic-hydrothermal fluid activity on the basis of metallogenic theory, and orebody I forms the upper part of orebody III based on the Tongshan Cu deposit model. Orebody III was destroyed by Tongshan Fault in Jurassic period.

### Discussion and conclusions

On basis of the methodology and research results described in this study, the necessary geological constraints and feasible regular block cell provide for accurate construction of a 3D Cu-grade model for several orebodies and for identification of continuous sensitive differences of Cu-grade in a single orebody, which can be used to infer a magmatic-hydrothermal fluid origin for the Tongshan Cu deposit. The 3D improved IDW modeling associated with virtual boreholes provides for accurate estimation of potential targets of mineralization. Although 3D kriging modeling and 3D improved IDW modeling are different algorithms, they can be used for cross-validation of results to assist accurate identification of potential targets of mineralization in the study area.

From comparison of 3D orebodies, 3D geological object models, and 3D potential targets of mineralization in study area, six new insights have been determined:

(1)Orebodies I–III of the Tongshan Cu deposit form a single ore-forming system and can be explained by an integrated metallogenic model.

(2)Orebody II is displaced by the Tongshan Fault, so that part of the orebody remains along the footwall of the fault as shown by the 1080–1120 exploration line at a depth of between 800 and 1000 m.

(3)Orebody III is the mineralization center of the Tongshan Cu deposit, and occurs at a depth of 1732 m based on BP network estimation using the longest borehole data and the 3D object model.

(4)In the 3D IDW model of the Tongshan Cu deposit, potential targets of mineralization are significantly larger than those of the 3D kriging model, e.g., the deposit may also have potential targets of Au mineralization.

(5)The geophysical anomaly reflects the northwestern distribution of the orebodies, and the southernmost anomaly is a potential target of mineralization (Fig. 7). Combining information on the Tongshan Fault with orebodies I–III, indicates that the displaced part of orebody II is an important potential target for exploration.

(6)A 500 m buffer based on orebody III involves the entire Tongshan Cu deposit in one 3D environment (Fig. 7). The buffer distance of 500 m was estimated by statistical analysis of geochemical data of mineralized rock in the alteration zone and indicates that orebody III is the center of mineralization in study area. Orebodies I and II have spatio-temporal and casual relationships with orebody III in accord with the magmatic-hypothermal theory of mineralization. On basis of our results, the original 3D model of orebody II model is reconstructed in Fig. 8.

(7)The proposed 3D geological model requires project management and the acquisition, extraction, interpretation, and validation from theoretical to digital modeling of geological and geophysical data 3D geological modeling methodology can be developed on the basis of integrated 3DGM and mathematics methods. Integrated methods can make use of 2D and 3D geological, geophysical, and geochemical data to accurately derive and integrate geosciences information for mineral resources assessment. With the development and implementation of second mineral resource exploration space at a depth of 500–2000 m in China, 3DGM will be one of the most important and necessary technologies for the deep exploration of mineral resources.

Source: 3D geological modeling for mineral resource assessment of the Tongshan Cu deposit, Heilongjiang Province, China

Authors: Gongwen Wang, Lei Huang