Regional geological three-dimensional model construction method based on multi-source data fusion

Abstract

The invention relates to the technical field of geological modeling, in particular to a regional geological three-dimensional model construction method based on multi-source data fusion, and solves the technical problems that in the prior art, in the geological model construction process, the form of a model is distorted in a drilling-free region, a topological structure is not connected, and boundary extraction is inaccurate. The method comprises the following steps: constructing a multi-source feature field fusing regional geophysical data and geological structure information; the multi-source characteristic field is used for representing geological attribute impedance and a geological structure direction; according to the multi-source characteristic field and the multiple discrete mineralization centers, determining an ore body communication framework connecting the multiple discrete mineralization centers; and according to the ore body connected skeleton and the geological attribute impedance, reconstructing a three-dimensional geological entity model through interface evolution.

Description

Technical Field

[0001] This invention relates to the field of geological modeling technology, specifically to a method for constructing a regional geological 3D model based on multi-source data fusion. Background Technology

[0002] In the field of deep metal mineral exploration and geological modeling, constructing high-precision three-dimensional geological models is the core foundation for resource assessment and engineering decision-making. However, the data in the current exploration area generally exhibits a polarized characteristic: borehole data has precise location and clear lithology, but can only provide discrete point or line constraints, and cannot reflect the spatial continuity of geological bodies between boreholes; geophysical inversion data can achieve full coverage, but is affected by the multiple solutions and smoothing effect of inversion algorithms, resulting in unclear boundaries and difficulty in directly indicating the fine connectivity structure of ore bodies.

[0003] In view of the characteristics of the above-mentioned multi-source data, existing geological modeling methods are mainly divided into two categories: one is explicit modeling based on geometric interpolation, which directly connects sparse borehole control points through Euclidean distance; the other is implicit modeling based on voxels, which relies on a single attribute threshold to extract geological bodies. Both methods attempt to achieve three-dimensional reconstruction of geological bodies through their respective technical paths, but neither has effectively solved the core contradiction in the process of multi-source data fusion.

[0004] Existing explicit modeling methods, by ignoring the anisotropic characteristics of the subsurface medium, are prone to producing erroneous connections that cross layers or violate geological laws in borehole-free areas with large borehole spacing, resulting in distorted model morphology. Implicit modeling methods, on the other hand, are prone to problems such as model fragmentation, disconnection, or inclusion of a large amount of surrounding rock background due to the inability of a single threshold to take into account the gradient differences between the mineralization center and the edge. Overall, it is difficult to effectively utilize fuzzy geophysical inversion data to construct a 3D model that conforms to geological laws under the hard constraints of sparse boreholes. Summary of the Invention

[0005] To address the technical problems in existing geological model construction techniques, such as morphological distortion, disconnected topology, and inaccurate boundary extraction in areas without boreholes, the present invention aims to provide a method for constructing regional geological 3D models based on multi-source data fusion. The specific technical solution adopted is as follows: Firstly, a method for constructing a regional geological 3D model based on multi-source data fusion is provided, including: constructing a multi-source feature field that integrates regional geophysical data and geological structural information; using the multi-source feature field to characterize geological attribute impedance and geological structural direction; determining the ore body connectivity framework connecting multiple discrete mineralization centers based on the multi-source feature field and multiple discrete mineralization centers; and reconstructing a 3D geological entity model through interface evolution based on the ore body connectivity framework and geological attribute impedance.

[0006] Based on the above technical solution, in the method for constructing a regional geological 3D model based on multi-source data fusion provided by the present invention, a multi-source feature field that integrates geological attribute impedance and geological structural direction is constructed. The ore body connectivity skeleton connecting discrete mineralization centers is determined by relying on the constraint information implicit in the multi-source data. Then, the 3D geological entity is reconstructed through interface evolution by combining geological attribute impedance. This effectively solves the problems of morphological distortion of poreless areas, poor model topological connectivity, and inaccurate boundary definition in traditional modeling. It realizes the deep fusion and utilization of multi-source heterogeneous data and finally constructs a high-precision regional geological 3D model with topological connectivity, morphology conforming to geological laws, and accurate boundaries.

[0007] In conjunction with the first aspect mentioned above, in one possible implementation, the method for constructing a multi-source feature field that integrates regional geophysical data and geological structural information specifically includes: establishing a three-dimensional grid covering the target area; constructing a geological attribute impedance field to characterize the ease of medium connectivity based on regional geophysical data; constructing a macroscopic structural tensor field to characterize macroscopic structural directions based on regional geological structural information; and constructing a microscopic texture tensor field to characterize local texture directions based on the geological attribute impedance field.

[0008] In conjunction with the first aspect above, in one possible implementation, the method for constructing a geological attribute impedance field to characterize the ease of medium connectivity specifically includes: mapping regional geophysical data to a three-dimensional grid; normalizing the mapped geophysical data; performing a nonlinear transformation on the normalized geophysical field to enhance the difference between the mineralized area and the surrounding rock area, thereby obtaining the geological attribute impedance field; the magnitude of the geological attribute impedance field is positively correlated with the difficulty of medium connectivity.

[0009] In conjunction with the first aspect above, in one possible implementation, the method for determining the orebody connectivity framework connecting multiple discrete mineralization centers specifically includes: determining the path connectivity cost based on the geological attribute impedance field, macroscopic structural tensor field, and microscopic texture tensor field; the path connectivity cost is used to characterize the comprehensive cost of connecting two points in a three-dimensional grid; constructing a spatial adjacency graph based on multiple discrete mineralization centers; constructing minimum cost connectivity paths corresponding to each connecting edge in the spatial adjacency graph based on the path connectivity cost, obtaining multiple initial connectivity paths; filtering the multiple initial connectivity paths based on the path connectivity cost, and marking the grid positions traversed by the retained connectivity paths after filtering, thereby generating the orebody connectivity framework.

[0010] In conjunction with the first aspect above, in one possible implementation, the method for determining the path connectivity cost specifically includes: weighted fusion of the macroscopic structural tensor field and the microscopic texture tensor field to obtain a hybrid guiding tensor field; and determining the single-step path connectivity cost from the current grid position to the adjacent grid position based on the geological attribute impedance value of the current grid position and the projection value of the hybrid guiding tensor field on the current connection direction.

[0011] In conjunction with the first aspect above, in one possible implementation, the method for filtering multiple initial connected paths based on path cost value specifically includes: determining the cost value per unit length of each initial connected path based on the total cost value and geometric length of each initial connected path; determining the cost value threshold based on the distribution of the cost value per unit length of multiple initial connected paths using a preset segmentation threshold determination method; and removing initial connected paths whose cost value per unit length is greater than the cost value threshold.

[0012] In conjunction with the first aspect mentioned above, in one possible implementation, the method for reconstructing a three-dimensional geological entity model based on the ore body connectivity framework and geological property impedance through interface evolution specifically includes: generating an initial geological interface distance field based on the ore body connectivity framework; the initial geological interface distance field is used to implicitly characterize the initial geological body boundary; generating a boundary evolution rate field based on the geological property impedance field; the boundary evolution rate field is used to control the rate at which the geological interface expands outward or contracts inward; driving the geological interface evolution through a level set evolution method based on the initial geological interface distance field and the boundary evolution rate field until a stable state is reached; and extracting the geological interface after the evolution has stabilized to obtain a three-dimensional geological entity model.

[0013] In conjunction with the first aspect above, in one possible implementation, the method for generating the boundary evolution rate field specifically includes: determining an impedance reference threshold based on the numerical distribution of the geological attribute impedance field using a preset segmentation threshold determination method; determining the boundary evolution rate value for each grid location based on the impedance reference threshold and the geological attribute impedance value for each grid location, thereby constituting the boundary evolution rate field; when the geological attribute impedance value is lower than the impedance reference threshold, the corresponding boundary evolution rate value is positive, and when the geological attribute impedance value is higher than the impedance reference threshold, the corresponding boundary evolution rate value is negative.

[0014] In conjunction with the first aspect above, in one possible implementation, the method further includes, during the process of driving the geological interface evolution through the level set evolution method, forcibly constraining the geological interface distance value of the grid position corresponding to the ore body's connected skeleton to a preset internal distance value.

[0015] In conjunction with the first aspect above, in one possible implementation, the method for determining multiple discrete mineralization centers specifically includes: obtaining the mineralization location in borehole exploration data and mapping the mineralization location to a three-dimensional grid; in regional geophysical data, searching for geophysical attribute extreme points within a preset neighborhood of the mapped mineralization location, and establishing the grid location corresponding to the extreme point as the mineralization center.

[0016] The present invention has the following beneficial effects: By constructing a multi-source feature field that integrates geological attribute impedance and geological structural direction, and relying on the constraint information implicit in the multi-source data, the ore body connectivity framework connecting discrete mineralization centers is determined. Then, by combining geological attribute impedance and reconstructing three-dimensional geological entities through interface evolution, the problems of morphological distortion in pore-free areas, poor model topological connectivity, and inaccurate boundary definition in traditional modeling are effectively solved. This achieves deep integration and utilization of multi-source heterogeneous data, and finally constructs a high-precision regional geological three-dimensional model with topological connectivity, morphology conforming to geological laws, and accurate boundaries. Attached Figure Description

[0017] To more clearly illustrate the technical solutions and advantages in the embodiments of the present invention or the prior art, the drawings used in the description of the embodiments or the prior art will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0018] Figure 1 A flowchart illustrating a method for constructing a regional geological 3D model based on multi-source data fusion, as provided in one embodiment of the present invention; Figure 2 This is a schematic diagram of the hardware structure of a regional geological three-dimensional model construction device based on multi-source data fusion, provided as an embodiment of the present invention. Detailed Implementation

[0019] To further illustrate the technical means and effects adopted by the present invention to achieve its intended purpose, the following, in conjunction with the accompanying drawings and preferred embodiments, details the specific implementation, structure, features, and effects of a regional geological three-dimensional model construction method based on multi-source data fusion proposed according to the present invention. In the following description, different "one embodiment" or "another embodiment" do not necessarily refer to the same embodiment. Furthermore, specific features, structures, or characteristics in one or more embodiments can be combined in any suitable form.

[0020] Unless otherwise defined, all technical and scientific terms used herein have the same meaning as commonly understood by one of ordinary skill in the art to which this invention pertains.

[0021] The following description, in conjunction with the accompanying drawings, details a specific scheme for constructing a regional geological three-dimensional model based on multi-source data fusion, provided by the present invention.

[0022] Please see Figure 1 The diagram illustrates a flowchart of a method for constructing a regional geological 3D model based on multi-source data fusion, according to an embodiment of the present invention. This method includes: S1. Construct a multi-source feature field that integrates regional geophysical data and geological structural information.

[0023] Multi-source characteristic fields are used to characterize geological property impedance and geological structural orientation.

[0024] In some implementation methods, the method of constructing a multi-source feature field that integrates regional geophysical data and geological structural information can be specifically implemented through the following steps S11 to S14, which are explained in detail below: S11. Create a three-dimensional mesh covering the target area.

[0025] In some implementation methods, the spatial range of the three-dimensional grid must first be accurately defined by combining the topography of the target exploration area, the scope of geological exploration work, and the coverage boundary of geophysical data. This ensures that the grid can completely cover the target exploration area and the necessary surrounding buffer zone, avoiding boundary distortion caused by data truncation.

[0026] Next, the core parameters of the grid are set, including the origin coordinates, voxel resolution, and number of grids: the origin coordinates can be selected from the surface projection point of the southwest corner of the target area as the starting reference datum of the grid in the three-dimensional Cartesian coordinate system; the voxel resolution needs to be determined according to the exploration accuracy requirements and the balance of computing resources. For example, if the target area needs to balance detail and efficiency, the voxel resolution along the x, y, and z axes can be set to 5m; the number of grids is derived from the area range and voxel resolution. Assuming that the target area has a length of 1000m in the x direction, a length of 1000m in the y direction, and a depth of 750m in the z direction, the corresponding number of grids in the x direction is 200, the number of grids in the y direction is 200, and the number of grids in the z direction is 150.

[0027] Based on the above parameters, any voxel cell in the mesh can be represented by a unique integer index vector. The identifier, where i, j, and k are integer indices in the x, y, and z axes respectively, and T is the transpose symbol, converts a row vector into a column vector, changing only the vector arrangement without changing the numerical values, in order to adapt to the standard specifications of linear algebra operations such as tensors and matrices.

[0028] This step constructs a unified spatial computing container, eliminating the spatial reference difference between the continuous coordinates of borehole data and the voxel structure of geophysical data. It provides a standardized spatial framework for the subsequent alignment, storage, and computation of multi-source data, ensuring that all data can be processed interactively under the same coordinate system.

[0029] S12. Based on regional geophysical data, construct a geological property impedance field to characterize the ease or difficulty of medium connectivity.

[0030] In some implementations, regional geophysical data is mapped to a three-dimensional grid; the mapped geophysical data is normalized, and the normalized geophysical field is nonlinearly transformed to enhance the difference between the mineralized area and the surrounding rock area, thus obtaining a geological property impedance field; the magnitude of the geological property impedance field is positively correlated with the difficulty of medium connectivity.

[0031] Specifically, the original 3D geophysical inversion data (such as resistivity and polarizability volume data) covering the target area are first mapped to the established 3D grid. The mapping process can select an appropriate resampling technique based on the data characteristics: if the original data is in a regular grid format, trilinear interpolation can be used to obtain the attribute values ​​of the target voxels by linear interpolation of the attribute values ​​of adjacent voxels, ensuring a smooth data transition; if the original data has obvious spatial correlation, kriging interpolation can be used to achieve an interpolation mapping that better fits the geological laws by utilizing the spatial variogram model of the data, and finally obtain a discretized attribute field that perfectly matches the 3D grid.

[0032] Next, data cleaning and truncation are performed on the discretized attribute field. First, the histogram distribution of the attribute field is statistically analyzed to determine the distribution characteristics of the data. Then, a truncation threshold is set. For example, for resistivity data, if the statistical analysis shows that the value corresponding to the 1% quantile is 20 Ω·m and the value corresponding to the 99% quantile is 500 Ω·m, then all data less than 20 Ω·m are truncated to 20 Ω·m, and data greater than 500 Ω·m are truncated to 500 Ω·m. This effectively eliminates extreme outliers in the inversion data caused by instrument errors and geological noise, avoiding such data interference with the accuracy of subsequent calculations and ensuring the rationality of the data distribution.

[0033] Then, the truncated discretized attribute field was subjected to max-min normalization to unify geophysical attributes with different dimensions and numerical ranges into the same numerical range, thus solving the problem that the original data had large differences in magnitude and could not be directly fused.

[0034] It should be noted that, unless otherwise specified, the normalization method and normalization function norm mentioned in the embodiments of this invention all adopt maximum and minimum value normalization. The maximum and minimum values ​​are preset empirical extreme values ​​derived from a large amount of historical experimental data. If the calculation result exceeds the interval [0, 1], it is restricted to the range [0, 1] by a truncation function (i.e., if the result is less than 0, it is taken as 0; if it is greater than 1, it is taken as 1) to eliminate the influence of outliers on the result.

[0035] It should be understood that during the normalization process, polarity unification needs to be performed according to the geophysical data type: if the data value is positively correlated with the mineralization probability (such as polarizability), the reverse mapping needs to be performed first, that is, the difference between 1 and the maximum and minimum normalization results is taken as the final normalization result.

[0036] After normalization, a nonlinear transformation is performed on the normalized attribute field to enhance the difference between the mineralized region and the surrounding rock region, thus obtaining the geological attribute impedance field. : ; In the formula, Non-zero infinitesimal quantity (e.g., 10). -4 This is used to avoid calculation errors caused by the impedance field value being zero; For normalized attribute fields; The compactness factor of the medium (taken as an example of 2.0) is used to adjust the intensity of nonlinear enhancement; for the normalized property field conduct The exponentiation amplifies the numerical difference between the mineralized zone and the surrounding rock zone, and then... Control the weight of this enhancement, and finally... Adding them together ensures the result is not zero, thus yielding... A property field that indirectly characterizes the ease or difficulty of medium connectivity, with values ​​approaching... A value close to 1 indicates low connectivity resistance and good connectivity; a value close to 1 indicates high connectivity resistance and poor connectivity.

[0037] S13. Based on regional geological structure information, construct a macroscopic structural tensor field to characterize the direction of macroscopic structures.

[0038] In some implementation methods, the dominant structural features of the target area are first identified based on regional geological survey reports, geological mapping data, and field measurement data. The normal unit vector of the dominant structural surface (such as the main stratigraphic bedding plane and the main fault bedding plane) is determined. This vector is the core benchmark representing the macroscopic structural direction. For example, if the dominant structural surface of the target area is trending east-west, dips north, and has a dip angle of 60°, the corresponding normal unit vector can be calculated as (0, 0.866, 0.5). Its physical meaning is that it points to the direction with the greatest resistance to the extension of the structural surface, that is, the direction that cuts through the structural surface.

[0039] Next, the outer product operator is used to construct the primitive construction tensor for each voxel in the 3D mesh. : ; In the formula, It is a unit tensor used to maintain the basic structure of tensors; This is the macroscopic anisotropy strength coefficient (taken as 10.0 for example), used to adjust the constraint strength in the macroscopic structural direction; The normal unit vector of the dominant structural surface; for The transpose; through the unit tensor To ensure the positive definiteness of tensors, and at the same time through The effect of magnifying the outer product term of the normal unit vector, where Essentially, a term is a tensor component constructed along the normal direction. The larger the value of , the more significant the influence of this component, making the original construction tensor... The larger the value along the normal direction, the stronger the resistance to connectivity in that direction, and the more likely it is to result in... It is a second-order tensor used to initially quantify the anisotropic characteristics of macroscopic structures and clarify the differences in connectivity resistance in different directions.

[0040] To ensure that the scale of the macroscopic structural tensor is consistent with that of the subsequent microscopic texture tensor, spectral norm normalization is performed on the original structural tensor. Each element of the original structural tensor is divided by its largest eigenvalue, so that the largest eigenvalue of the normalized tensor is 1, thus standardizing the tensor magnitude and obtaining the macroscopic structural tensor field. This field accurately represents the directional constraints of the macroscopic structure. The normalized tensor has a larger eigenvalue along the normal to the dominant structural surface and a smaller eigenvalue along the stratigraphic strike surface, which can guide the subsequent path to extend along the structural strike. At the same time, because it is consistent with the scale of the microscopic tensor, it can be directly fused for calculation.

[0041] S14. Based on the geological property impedance field, construct a micro-texture tensor field to characterize the local texture direction.

[0042] In some implementations, a local neighborhood window is first defined for each voxel in the 3D mesh. The window size needs to be determined based on the texture features of the geological body and the mesh resolution. For example, if the geological texture of the target area is relatively fine, a 3×3×3 voxel window can be selected; if the texture is relatively coarse and contains large-scale continuous structures, a 5×5×5 voxel window can be selected. The function of this window is to cover the local spatial region around the current voxel in order to extract the texture direction information within this range and ensure the local continuity of texture features.

[0043] Next, the gradient vector of each voxel within the neighborhood window is calculated using the central difference method. The calculation of this gradient vector is based on the spatial distribution of the geological property impedance field. Specifically, for any voxel in the grid, the variation components in the x, y, and z coordinate axes are obtained by the central difference method: the x-axis component is obtained by dividing the impedance difference between the voxels to the right and left by twice the x-axis voxel resolution; the y-axis component is obtained by dividing the impedance difference between the voxels to the front and back by twice the y-axis voxel resolution; and the z-axis component is obtained by dividing the impedance difference between the voxels to the top and bottom by twice the z-axis voxel resolution. These three components are combined and transposed into a column vector (to adapt to subsequent tensor operations), which gives the gradient vector corresponding to the voxel. Its physical meaning is to characterize the spatial rate of change and direction of change of the geological property impedance field at the corresponding voxel position. The larger the magnitude of the gradient vector, the more drastic the impedance change at that position, and the more likely it is to be close to the boundary of the geological body.

[0044] Then construct the structure tensor matrix : ; In the formula, For voxel x, it represents the local neighborhood window; This is a Gaussian smoothing kernel function (specifically a three-dimensional Gaussian filter function), used to perform weighted smoothing on the calculation results in the neighborhood to suppress noise interference; The parameters of the Gaussian smoothing kernel correspond to the characteristic scale of the geological body (for example, take twice the voxel width, i.e., 10m). is the spatial offset between voxel x and its neighboring voxel y; Let y be the gradient vector of the geological property impedance field corresponding to voxel y; This is the transpose of the gradient vector. For each voxel within the neighborhood window, first calculate the outer product of its gradient vectors to obtain a matrix representing the gradient direction of that voxel. Then multiply this matrix by the Gaussian smoothing kernel weights determined by the spatial offset (the smaller the offset, the larger the weights, and the more significant the influence of neighboring voxels). Finally, sum the calculation results for all neighboring voxels to integrate and smooth the local gradient information, resulting in a positive semi-definite symmetric matrix. , is a dimensionless parameter, and its physical meaning is to characterize the directional distribution characteristics of the impedance field gradient in a local area, reflecting the trend of the micro-texture. Due to the weighting effect of the Gaussian smoothing kernel, high-frequency noise is effectively suppressed and the continuity of texture features is enhanced.

[0045] Subsequently, the structure tensor matrix was analyzed. Perform eigenvalue decomposition to extract the unit eigenvector corresponding to the largest eigenvalue. The physical meaning of this vector is that it represents the normal direction of the local texture, that is, the direction in which the impedance field changes most drastically, and the corresponding texture extension direction is perpendicular to this vector.

[0046] Next, construct the original gradient tensor. : ; In the formula, It is a unit tensor used to maintain the basic structure of tensors; This is the micro-guidance strength coefficient (taken as an example of 5.0), used to adjust the constraint strength of local texture direction; The unit eigenvector corresponding to the largest eigenvalue; for The transpose; through the unit tensor To ensure the positive definiteness of tensors, and at the same time through The effect of the outer product term of the amplified texture normal vector, where The term constructs tensor components along the texture normal direction. The larger the value of , the more significant the influence of this component, resulting in a larger value of the original gradient tensor along the texture normal, reflecting a stronger connectivity resistance in that direction. It is a second-order tensor, which preliminarily quantifies the anisotropic characteristics of micro-texture and clarifies the differences in connectivity resistance in different directions within a local region.

[0047] Finally, spectral norm normalization is performed on the original gradient tensor to obtain the normalized microtexture tensor field, which standardizes the magnitude of the microtexture tensor and accurately represents the directional constraints of the local texture. The normalized tensor has larger eigenvalues ​​along the texture normal and smaller eigenvalues ​​along the texture extension direction, which can guide the subsequent path to extend along the local texture direction. At the same time, because it is consistent with the scale of the macroscopic construction tensor, it can be directly fused and calculated.

[0048] S2. Based on the multi-source feature field and multiple discrete mineralization centers, determine the ore body connectivity framework connecting the multiple discrete mineralization centers.

[0049] In some implementations, methods for determining multiple discrete mineralization centers include: obtaining mineralization locations from borehole exploration data and mapping these locations to a three-dimensional grid; and in regional geophysical data, searching for geophysical attribute extreme points within a preset neighborhood of the mapped mineralization location and establishing the grid location corresponding to the extreme point as the mineralization center.

[0050] Specifically, the process begins by acquiring complete mineralization-related data from the borehole exploration database. This includes the three-dimensional geometric center coordinates of the mineralized section, lithological information, and corresponding geophysical attribute data. The three-dimensional geometric center coordinates of the mineralized section are calculated by interpolating the start and end depths of the mineralized section with the borehole trajectory coordinates, ensuring that the coordinates accurately reflect the core location of the mineralized body. Subsequently, these continuous three-dimensional coordinates of the mineralized locations are mapped onto an established three-dimensional grid. Using the grid's coordinate datum and voxel indexing rules, each mineralization center coordinate is converted into a corresponding initial grid index.

[0051] Next, a preset neighborhood range is set for the ore-bearing mapping location. The size of the neighborhood needs to be determined based on the spatial scale of the geological body and the grid resolution. For example, if the voxel resolution is 5m, a 3×3×3 voxel neighborhood (i.e., a cubic area with a spatial range of ±5m) can be selected. This neighborhood range can cover the potential anomaly area around the ore-bearing location while avoiding the introduction of irrelevant background noise over an excessively large area. Then, within this preset neighborhood, based on the mapped regional geophysical data (which has been converted into a discretized attribute field), a preset extreme point detection algorithm is used to search for geophysical attribute extreme points. If the target mineralization body corresponds to a low resistivity anomaly (such as a sulfide ore body), the minimum attribute value point in the neighborhood is searched; if it corresponds to a high polarization anomaly, the maximum attribute value point in the neighborhood is searched. For example, a local extreme value scanning algorithm can be used to compare the attribute values ​​of each voxel with all other voxels in the neighborhood, and select the voxels with the largest or smallest attribute values ​​as candidate extreme points.

[0052] Meanwhile, to avoid extreme point drift caused by local noise, the validity of candidate extreme points needs to be verified: calculate the variance of geophysical attribute values ​​in the preset neighborhood, set a preset variance threshold (for example, 10% of the variance of geophysical attributes in the whole area), if the variance in the neighborhood is greater than the threshold, it indicates that there is a significant attribute change in the neighborhood, and the candidate extreme point is the real mineralization anomaly center, and its corresponding grid position is established as the mineralization center; if the variance in the neighborhood is less than or equal to the threshold, it indicates that the neighborhood is an attribute flat area or a background noise area, and at this time, the grid position corresponding to the initial grid index is directly established as the mineralization center to avoid noise interference causing the mineralization center to deviate from the real mineralization location.

[0053] The above steps, through neighborhood extremum search and variance verification, corrected the slight spatial misalignment between borehole data and geophysical data, ensuring that the mineralization center is accurately located in the core area of ​​the geophysical anomaly.

[0054] In some implementations, the method for determining the orebody connectivity framework connecting multiple discrete mineralization centers can be specifically implemented through the following steps S21 to S24, which are explained in detail below: S21. Determine the path connectivity value based on the geological attribute impedance field, macroscopic structural tensor field, and microscopic texture tensor field.

[0055] Path connectivity cost is used to characterize the overall cost of connecting two points in a 3D mesh.

[0056] In some implementations, the macroscopic structural tensor field and the microscopic texture tensor field are weighted and fused to obtain a hybrid directional tensor field; based on the geological attribute impedance value of the current grid position and the projection value of the hybrid directional tensor field on the current connection direction, the single-step path connectivity cost of moving from the current grid position to the adjacent grid position is determined.

[0057] Specifically, since the macroscopic tensor field constrains the macroscopic structural orientation of the region and the microscopic texture tensor field constrains the local texture direction, it is necessary to balance the constraint weights of the two through weighted fusion. A scale balance coefficient β (range [0, 1], for example, 0.5, so that the macroscopic and microscopic constraint weights each account for 50%) is set, and the hybrid guided tensor field is calculated based on this coefficient: ; In the formula, This is the hybrid directional tensor field corresponding to the current grid position u, used to simultaneously characterize macroscopic and microscopic directional constraints; To construct a tensor field for macroscopic purposes, representing the directional constraints of the macroscopic construction; It is a micro-texture tensor field, representing the constraint of local texture direction; The scale balance coefficients are used to balance the constraint weights of the macroscopic construction tensor field and the microscopic texture tensor field; after weighted fusion, a tensor containing both macroscopic construction and microscopic texture direction information is obtained. As a positive definite second-order tensor, its physical meaning is to integrate macroscopic and microscopic anisotropic directional constraints. Its function is to allow subsequent paths to conform to both regional tectonics and local texture directions, thereby improving the geological rationality of the path.

[0058] After constructing the hybrid guided tensor field, determine the single-step path connectivity cost from the current grid position u to the adjacent grid position v: First, calculate the unit displacement vector d between u and v (obtained by dividing the coordinate difference between v and u by the Euclidean distance (i.e., the L2 norm, such as the Euclidean norm, representing the displacement direction); then, read the geological attribute impedance value corresponding to u, combine it with the projection value of the hybrid guided tensor field in the d direction, and the geometric distance between u and v, to calculate the single-step cost. : ; In the formula, The cost of single-step path connectivity; The geological property impedance value of u represents the material constraint. It is a unit displacement vector. Transpose it; Let u be the geometric distance between u and v (L2 norm). By integrating macroscopic and microscopic directional constraints and combining them with material property constraints, a path substitution value is obtained that conforms to both the material characteristics of the geological body and the spatial extension law. The final result is a length-weighted value with a range of []. ,+∞)( (This is the minimum non-zero value of the impedance field), its physical meaning is the comprehensive cost of single-step connection, and its function is to guide the path to extend preferentially along areas with low impedance and in accordance with geological direction.

[0059] S22. Construct a spatial adjacency graph based on multiple discrete mineralization centers.

[0060] In some implementations, the three-dimensional mesh coordinates of all discrete mineralization centers are first obtained to form a set of mineralization centers. If pairwise path calculations are performed on this set, the computational complexity increases quadratically with the number of mineralization centers, leading to computational explosion on large-scale grids. Therefore, it is necessary to construct a spatial adjacency graph to reduce the computational load.

[0061] Specifically, a K-nearest neighbor strategy can be used: Set the number of nearest neighbors K (e.g., 6), and for each mineralization center, search for the 6 nearest neighbors in Euclidean space and establish connecting edges. To ensure graph connectivity, supplement the edges to construct the minimum spanning tree, ensuring all mineralization centers are in the same connected component. Alternatively, a 3D Delaunay triangulation strategy can be used: Perform 3D Delaunay triangulation on the set of mineralization centers, using the edges of the generated tetrahedral mesh as connecting edges. This strategy guarantees the spatial proximity of adjacent edges.

[0062] S23. Based on the path connectivity cost, construct the minimum cost connectivity path corresponding to each connecting edge in the spatial adjacency graph to obtain multiple initial connectivity paths.

[0063] In some implementations, the 3D mesh is first treated as a weighted directed graph, with voxels as nodes and the weight of the edges connecting adjacent voxels as the cost of a single-step path connection. For each edge in the adjacency graph (corresponding to the starting point)... and the finish line ), Execution path search: First, initialize: Let The cumulative cost of one voxel is 0, while the cumulative cost of the other voxels is infinite; create a priority queue to... Insert into the queue (the queue is sorted in ascending order of cumulative cost).

[0064] Then perform iterative search: pop the node u with the minimum cumulative cost from the queue, if u is... If the search ends, then the search continues; otherwise, iterate through the 26 neighboring volees of u (adjacent volees within ±1 range of the x, y, and z axes), calculate the single-step generation value from u to v; if the cumulative generation value of u plus the single-step generation value is less than the current cumulative generation value of v, then update the cumulative generation value of v, record its parent node as u, and insert v into the queue. Finally, backtrack the path: from... Tracing back to the parent node Extract discrete paths and record the total value of the agent.

[0065] Repeat the above process for all connected edges to obtain an initial set of paths. Guided by cost, the minimum cost path connecting the mineralization center is obtained. The path prioritizes extending along areas with low impedance and in line with geological direction, which conforms to the connectivity law of geological bodies.

[0066] S24. Filter multiple initial connectivity paths according to their path value, and mark the grid positions traversed by the retained connectivity paths after filtering to generate the ore body connectivity skeleton.

[0067] In some implementations, the method for filtering multiple initial connected paths based on path cost value includes: determining the cost value per unit length of each initial connected path based on the total cost value and geometric length of each initial connected path; determining the cost value threshold by using a preset segmentation threshold determination method based on the distribution of cost values ​​per unit length of multiple initial connected paths; and removing initial connected paths whose cost value per unit length is greater than the cost value threshold.

[0068] Specifically, first calculate the cost per unit length for each initial path: traverse the set of initial paths, extract the total cost and geometric length (number of voxels) of each path, divide the total cost by the geometric length to obtain the cost per unit length (for example, a total cost of 200, a length of 40 voxels, and a cost per unit length of 5).

[0069] Subsequently, the cost threshold is determined: a frequency distribution histogram is constructed based on the cost of all units of length, and the segmentation threshold is automatically calculated using the Otsu method (e.g., a threshold of 6 is obtained). This method distinguishes between connected and disconnected paths by maximizing the ratio of inter-class to intra-class variance.

[0070] Next, the paths are filtered: paths with a unit length cost exceeding a threshold are removed, and valid paths are retained. It should be understood that before removing paths with a unit length cost exceeding the threshold, it is checked whether the path is a minimum spanning tree (MST) edge of the spatial adjacency graph; if it is an MST edge, the path is forcibly retained and not removed to ensure the overall connectivity of the ore body's connecting framework.

[0071] Finally, the skeleton is generated: an all-zero field of the same size as the mesh is initialized, and the voxels traversed by the effective path are marked as 1 to obtain the ore body connected skeleton.

[0072] Through the above steps, erroneous paths that cross high-resistivity surrounding rocks are eliminated, and valid paths that conform to geological genesis are retained. The generated skeleton topology conforms to geological laws, and its geometry is controlled by multi-source data, providing a reliable geometric centerline for subsequent volume reconstruction.

[0073] S3. Based on the ore body's connectivity framework and geological property impedance, reconstruct a three-dimensional geological entity model through interface evolution.

[0074] In some implementations, the method for reconstructing the three-dimensional geological entity model can be specifically implemented through the following steps S31 to S34, which are explained in detail below: S31. Generate the initial geological interface distance field based on the ore body connectivity framework.

[0075] The initial geological interface distance field is used to implicitly characterize the boundary of the initial geological body.

[0076] In some implementations, the orebody connectivity framework is first read. This is a binary field of the same size as the 3D geological discrete grid, where voxels with a value of 1 correspond to the central region of the orebody connectivity framework and serve as the core benchmark for subsequent geological body growth. Next, for each voxel in the 3D geological discrete grid, its Euclidean distance to the nearest framework voxel (i.e., the voxel marked as 1) is calculated. The grid is traversed, and the Euclidean distance from each voxel to the nearest target voxel is output, which represents the spatial distance between each voxel and the orebody framework. The smaller the value, the closer it is to the framework region.

[0077] Then, an initial channel radius is set, typically 1.5 to 2.0 times the voxel width (e.g., when the voxel width is 5m, the initial channel radius is 7.5m), to define the basic channel width of the initial geological body. The initial geological interface distance field is obtained by subtracting the Euclidean distance from the voxel to the framework, defining the spatial extent of the initial geological body: regions with a result less than 0 are the interior of the geological body, zero isosurfaces equal to 0 are the boundaries of the initial geological body, and regions greater than 0 are the surrounding rock.

[0078] S32. Generate the boundary evolution rate field based on the geological property impedance field.

[0079] Boundary evolution rate fields are used to control the rate at which geological interfaces expand outward or contract inward.

[0080] In some implementations, an impedance reference threshold is determined by a preset segmentation threshold determination method based on the numerical distribution of the geological attribute impedance field. The boundary evolution rate value of each grid location is determined based on the impedance reference threshold and the geological attribute impedance value of each grid location to form a boundary evolution rate field. When the geological attribute impedance value is lower than the impedance reference threshold, the corresponding boundary evolution rate value is positive, and when the geological attribute impedance value is higher than the impedance reference threshold, the corresponding boundary evolution rate value is negative.

[0081] Specifically, the impedance reference threshold is first determined based on the numerical distribution of the geological attribute impedance field: the geological attribute impedance field is read, and the data of the geological attribute impedance field of the whole area is analyzed by the Ottoman method, and the optimal segmentation threshold for distinguishing between mineralized areas and surrounding rock areas is calculated (for example, for the geological attribute impedance field, the optimal segmentation threshold is calculated to be 0.6). This threshold is the impedance boundary value between mineralization and surrounding rock.

[0082] Based on the optimal segmentation threshold Calculate the adaptive reference threshold : ; In the formula, The impedance sensitivity coefficient (a constant greater than 0, e.g., 10.0) controls the sensitivity of the velocity field to impedance changes; the final result... is a dimensionless parameter with a range of (0, 1). Its physical meaning is the reference threshold of the velocity field, used to distinguish the positive and negative intervals of the velocity. Its function is to allow the velocity field to adaptively match the distribution characteristics of the impedance field.

[0083] Based on adaptive reference threshold With geological property impedance field Determine the boundary evolution rate value for each grid location: ; In the formula, Let be the boundary evolution rate value corresponding to voxel x, be a dimensionless parameter with a range of (-1, 1), and physically represent the evolution trend of the geological interface: when < (Mineralized area), the first item is greater than , Greater than 0, drives interface expansion; when Greater than (Surrounding rock area), the first item is less than , <0, the driver interface shrinks; when ≈ (Boundary area) ≈0, interface evolution stagnates.

[0084] Finally, for the voxels at the edges of the 3D mesh, force settings are applied. =-1.0, construct an artificial absorbing boundary to prevent the model evolution from overflowing the computational domain.

[0085] The above steps generate a boundary evolution rate field that precisely matches geological properties, enabling adaptive control of the interface as it expands in the mineralized zone, contracts in the surrounding rock zone, and stagnates in the boundary zone. At the same time, it avoids the deadlock problem of extreme evolution across the entire field and provides rate constraints for the reasonable growth of the interface.

[0086] S33. Based on the initial geological interface distance field and boundary evolution rate field, the geological interface is driven to evolve through the level set evolution method until a stable state is reached.

[0087] In some implementations, during the evolution process, the geological interface distance value of the grid position corresponding to the ore body's connecting skeleton can be forcibly constrained to a preset internal distance value.

[0088] Specifically, the initial geological interface distance field and boundary evolution rate field are first read, and the geological interface evolution is driven by the level set evolution method: the geological interface is represented as the zero isosurface of the distance field, and the interface deformation is achieved by iteratively updating the distance field. The numerical discretization form of the evolution equation is as follows: ; In the formula, This represents the geological interface distance field after the (t+1)th iteration; The distance field after the t-th iteration; The time step must satisfy the Courant-Friedrich-Levy (CFL) stability condition, for example, taking... ≤0.5); The gradient magnitude of the range field (calculated using an upwind difference scheme, such as a first-order upwind difference, to ensure evolution stability); final results The updated distance field has the dimension of length and the physical meaning of the updated geological body spatial extent. Its function is to iteratively allow the interface to evolve along the direction constrained by the rate field, gradually approaching the boundary of the real geological body.

[0089] After each iteration of the distance field update, a forced constraint operation is performed: the geological interface distance value of the grid position (i.e., the voxel marked as 1) corresponding to the ore body's connected skeleton is forcibly set to a preset internal distance value. ; In the formula, A binary field representing the connected framework of the ore body; The preset internal distance value (corresponding to the internal region of the geological body) is used; the updated distance field value is compared with the preset value: if the current voxel is a skeleton voxel, then the preset value is forcibly taken. If not, the updated value is retained. The final result ensures that the skeleton voxels are always inside the geological body (distance value less than 0). The physical meaning is to forcibly constrain the connectivity of the ore body skeleton. Its function is to prevent the skeleton from breaking due to local high impedance during the evolution process and to maintain the topological rationality of the geological body.

[0090] Finally, the evolution termination criterion is set: calculate the relative change rate of the volume enclosed by the zero isosurface between two adjacent iterations (the volume is obtained by multiplying the number of voxels with a statistical distance field less than 0 by the voxel volume). When the change rate is less than a preset threshold (e.g., 0.1%), it is determined that the evolution has reached a physical equilibrium state.

[0091] S34. Extract the geological interface after it has stabilized to obtain a three-dimensional geological entity model.

[0092] In some implementations, once the evolution reaches a stable state, the final geological interface distance field is read, and the zero isosurface of the distance field is extracted using the Marching Cubes algorithm: This algorithm traverses each cube cell in the three-dimensional mesh, and interpolates the triangular facets of the zero isosurface within the cell based on the distance field sign of the cell vertex (less than 0 for the interior of the geological body, greater than 0 for the surrounding rock), and finally splices all the triangular facets into a closed three-dimensional surface.

[0093] This closed triangular mesh surface is the regional geological three-dimensional solid model, which integrates the hard constraints of borehole data (ore body skeleton), the soft constraints of geophysical data (impedance field velocity), and geological topological features (interface connectivity), and can accurately characterize the spatial morphology and boundaries of geological bodies.

[0094] Based on the above technical solution, by constructing a multi-source feature field that integrates geological attribute impedance and geological structural direction, and relying on the constraint information implicit in the multi-source data to determine the ore body connectivity framework connecting discrete mineralization centers, and then combining geological attribute impedance to reconstruct three-dimensional geological entities through interface evolution, the problems of morphological distortion in pore-free areas, poor model topological connectivity, and inaccurate boundary definition in traditional modeling are effectively solved. This achieves deep integration and utilization of multi-source heterogeneous data, and finally constructs a high-precision regional geological three-dimensional model with topological connectivity, morphology conforming to geological laws, and accurate boundaries.

[0095] It should be noted that the order of the above embodiments of the present invention is merely for descriptive purposes and does not represent the superiority or inferiority of the embodiments. The processes depicted in the accompanying drawings do not necessarily require a specific or sequential order to achieve the desired result. In some embodiments, multitasking and parallel processing are also possible or may be advantageous.

[0096] The various embodiments in this specification are described in a progressive manner. The same or similar parts between the various embodiments can be referred to each other. Each embodiment focuses on describing the differences from other embodiments.

[0097] In this embodiment of the invention, the functional units of the regional geological 3D model construction device based on multi-source data fusion can be divided according to the above method example. For example, each function can be divided into separate functional units, or two or more functions can be integrated into one processing unit. The integrated unit can be implemented in hardware or as a software functional unit. It should be noted that the unit division in this embodiment is illustrative and only represents one logical functional division; other division methods may be used in actual implementation.

[0098] This invention also provides a schematic diagram of the hardware structure of a regional geological 3D model construction device based on multi-source data fusion, see [link / reference]. Figure 2The regional geological 3D model construction device 200 based on multi-source data fusion includes a processor 201, and optionally, a memory 202 connected to the processor 201.

[0099] In the first possible implementation, see Figure 2 The regional geological 3D model construction device 200 based on multi-source data fusion also includes a transceiver 203. The processor 201, memory 202, and transceiver 203 are connected via a bus. The transceiver 203 is used to communicate with other devices or communication networks. Optionally, the transceiver 203 may include a transmitter and a receiver. The device in the transceiver 203 that implements the receiving function can be considered as a receiver, which is used to perform the receiving steps in the embodiments of the present invention. The device in the transceiver 203 that implements the transmitting function can be considered as a transmitter, which is used to perform the transmitting steps in the embodiments of the present invention.

[0100] Based on the first possible implementation method Figure 2 The structural diagram shown can be used to illustrate the structure of the regional geological three-dimensional model construction device based on multi-source data fusion involved in the above embodiments.

[0101] in, Figure 2 The diagram can also illustrate the system chip in a regional geological 3D model construction device based on multi-source data fusion. In this case, the actions performed by the aforementioned regional geological 3D model construction device based on multi-source data fusion can be implemented by this system chip. The specific actions performed can be found above and will not be repeated here.

[0102] Although the invention has been described herein in conjunction with various embodiments, those skilled in the art, by reviewing the accompanying drawings and the disclosure, will understand and implement other variations of the disclosed embodiments in carrying out the claimed invention. In this invention, the word "comprising" does not exclude other components or steps, and "a" or "an" does not exclude a plurality. A single processor or other unit can implement several of the functions listed in this invention.

[0103] Although the invention has been described in conjunction with specific features and embodiments, it is obvious that various modifications and combinations can be made therein without departing from the spirit and scope of the invention. Accordingly, this specification and drawings are merely illustrative of the invention and are to be considered as covering any and all modifications, variations, combinations, or equivalents within the scope of the invention. Clearly, those skilled in the art can make various alterations and modifications to the invention without departing from its spirit and scope. Thus, if such modifications and modifications of the invention fall within the scope of the invention and its equivalents, the invention is also intended to include such modifications and modifications.

Claims

1. A method for constructing a regional geological 3D model based on multi-source data fusion, characterized in that, include: Construct a multi-source feature field that integrates regional geophysical data and geological structural information; The multi-source feature field is used to characterize geological property impedance and geological structural direction; Based on the multi-source feature field and multiple discrete mineralization centers, a ore body connectivity framework connecting the multiple discrete mineralization centers is determined. Based on the ore body's connectivity framework and geological property impedance, a three-dimensional geological entity model is reconstructed through interface evolution.

2. The method for constructing a regional geological three-dimensional model based on multi-source data fusion according to claim 1, characterized in that, Constructing a multi-source feature field that integrates regional geophysical data and geological structural information, including: Establish a three-dimensional mesh covering the target area; Based on regional geophysical data, a geological property impedance field is constructed to characterize the ease or difficulty of media connectivity. Based on regional geological structure information, a macroscopic structural tensor field is constructed to characterize the direction of macroscopic structures; Based on the geological property impedance field, a micro-texture tensor field is constructed to characterize the local texture direction.

3. The method for constructing a regional geological three-dimensional model based on multi-source data fusion according to claim 2, characterized in that, Constructing a geological property impedance field to characterize the ease of connectivity of a medium, including: Map the geophysical data of the region onto the three-dimensional grid; The mapped geophysical data is normalized, and the normalized geophysical field is nonlinearly transformed to enhance the difference between the mineralized area and the surrounding rock area, thus obtaining the geological attribute impedance field; the magnitude of the geological attribute impedance field is positively correlated with the difficulty of medium connectivity.

4. The method for constructing a regional geological three-dimensional model based on multi-source data fusion according to claim 2, characterized in that, Determining the orebody connectivity framework connecting the multiple discrete mineralization centers includes: The path connectivity cost is determined based on the geological property impedance field, the macroscopic structural tensor field, and the microscopic texture tensor field; the path connectivity cost is used to characterize the comprehensive cost of connecting two points in the three-dimensional mesh. A spatial adjacency graph is constructed based on the multiple discrete mineralization centers; Based on the path connectivity cost, construct the minimum cost connectivity path corresponding to each connecting edge in the spatial adjacency graph to obtain multiple initial connectivity paths; The multiple initial connectivity paths are filtered based on their path value, and the grid positions traversed by the retained connectivity paths are marked to generate the ore body connectivity skeleton.

5. The method for constructing a regional geological three-dimensional model based on multi-source data fusion according to claim 4, characterized in that, Determining the path connectivity cost includes: The macroscopic construction tensor field and the microscopic texture tensor field are weighted and fused to obtain a hybrid guided tensor field. Based on the geological attribute impedance value of the current grid location and the projection value of the hybrid guiding tensor field in the current connection direction, the single-step path connectivity cost of moving from the current grid location to the adjacent grid location is determined.

6. The method for constructing a regional geological three-dimensional model based on multi-source data fusion according to claim 4, characterized in that, The multiple initial connected paths are filtered based on their path cost, including: The unit length cost of each initial connected path is determined based on the total cost and geometric length of each initial connected path. Based on the unit length cost distribution of the multiple initial connected paths, the cost threshold is determined by a preset segmentation threshold determination method; Initial connected paths whose unit length cost is greater than the cost threshold are removed.

7. The method for constructing a regional geological three-dimensional model based on multi-source data fusion according to claim 2, characterized in that, Based on the connectivity framework and geological property impedance of the ore body, a three-dimensional geological entity model is reconstructed through interface evolution, including: Based on the ore body connectivity framework, an initial geological interface distance field is generated; the initial geological interface distance field is used to implicitly characterize the boundary of the initial geological body. Based on the geological property impedance field, a boundary evolution rate field is generated; the boundary evolution rate field is used to control the rate at which the geological interface expands outward or contracts inward. Based on the initial geological interface distance field and the boundary evolution rate field, the geological interface is driven to evolve using the level set evolution method until a stable state is reached; The geological interface after evolution and stabilization is extracted to obtain the three-dimensional geological entity model.

8. The method for constructing a regional geological three-dimensional model based on multi-source data fusion according to claim 7, characterized in that, Generate the boundary evolution rate field, including: Based on the numerical distribution of the geological property impedance field, an impedance reference threshold is determined using a preset segmentation threshold determination method. Based on the impedance reference threshold and the geological attribute impedance value of each grid location, the boundary evolution rate value of each grid location is determined to form the boundary evolution rate field; when the geological attribute impedance value is lower than the impedance reference threshold, the corresponding boundary evolution rate value is positive, and when the geological attribute impedance value is higher than the impedance reference threshold, the corresponding boundary evolution rate value is negative.

9. The method for constructing a regional geological three-dimensional model based on multi-source data fusion according to claim 7, characterized in that, The process of driving geological interface evolution through level set evolution also includes: The geological interface distance value of the grid position corresponding to the ore body connecting skeleton is forcibly constrained to a preset internal distance value.

10. The method for constructing a regional geological three-dimensional model based on multi-source data fusion according to claim 2, characterized in that, Identify multiple discrete mineralization centers, including: Obtain the mineralization locations from the borehole exploration data and map the mineralization locations onto the three-dimensional grid; In the regional geophysical data, geophysical attribute extreme points are searched within a preset neighborhood of the mineralization mapping location, and the grid location corresponding to the extreme point is established as the mineralization center.

Images