Joint interpretation method and system for gravity and magnetic data physical property inversion and interactive inversion

Abstract

The application discloses a gravity-magnetic data physical property inversion and interactive inversion alternately combined interpretation method and system, relates to the gravity-magnetic exploration field, and the method comprises the following steps: obtaining and preprocessing gravity-magnetic observation data; performing initial grid division on the underground space and executing large-grid-scale physical property inversion to obtain the physical property parameters of each initial grid unit; constructing an initial geology-geophysical model, taking the model as a starting point to execute interactive inversion to obtain an interactive inversion model; performing secondary division on the initial grid unit based on the interactive inversion model to obtain fine grid units; performing physical property parameter assignment according to the spatial position relationship between each fine grid unit and the boundary of the interactive inversion model to obtain fine grid model parameters, and constructing a physical property constraint condition; executing fine physical property inversion to obtain fine physical property inversion results, and combining the current interactive inversion model to output final interpretation achievements. The application realizes high-precision, high-efficiency and high-determinacy combined interpretation of the underground physical property distribution at the profile scale.

Description

Technical Field

[0001] This application relates to the field of gravity and magnetic exploration, and in particular to a method and system for the alternating joint interpretation of gravity and magnetic data property inversion and interactive inversion. Background Technology

[0002] Gravity and magnetic exploration is an important tool in geophysical exploration, widely used in metal mineral exploration, oil and gas resource development, geothermal resource evaluation, and geological structure research. Gravity and magnetic data inversion aims to infer the density and magnetic distribution of subsurface media based on observed gravity and magnetic anomalies, thereby revealing the spatial distribution and structural characteristics of geological bodies. However, due to the finiteness of observational data, the influence of noise, and the existence of source equivalence, gravity and magnetic data inversion suffers from severe ambiguity; that is, multiple different subsurface property distributions can produce almost identical observed anomalies. This inherent difficulty restricts the accuracy and reliability of quantitative interpretation of gravity and magnetic data and remains a core challenge that has long faced the field of geophysical inversion.

[0003] To address these challenges, geophysicists have developed various inversion methods, among which the two most representative and widely used are property inversion and interactive inversion. However, both of these methods have inherent limitations within their respective technical frameworks, hindering further improvements in the interpretation of inversion results.

[0004] Gravity and magnetic data property inversion discretizes subsurface space into regular grid cells. Based on the linear superposition relationship of potential field theory, it inverts the density or magnetic susceptibility value of each cell by solving a large set of linear equations. This method can characterize geological anomalies of arbitrarily complex shapes and quantitatively calculate information such as volume and physical properties. However, this method has the following significant drawbacks: (1) Severe ambiguity and poor stability of inversion results. The number of unknowns in the property inversion far exceeds the amount of observation data, the equations are highly underdetermined, and the influence of field source equivalence and observation errors leads to a huge solution space and prominent ambiguity. Although the Tikhonov regularization framework can make ill-conditioned problems tend to be well-conditioned to some extent by introducing model constraint terms, the choice of regularization function and the setting of regularization parameters have a significant impact on the inversion results. Different constraint strategies may lead to property models with huge differences.

[0005] (2) The skin effect leads to a severe lack of vertical resolution. The potential field kernel function decays exponentially with increasing depth, making the contribution of shallow cells to observed anomalies much greater than that of deep cells. The inversion results often over-concentrate physical property anomalies in the shallow layer, suppressing information from deeper layers. Although a depth-weighted function that monotonically increases with depth is usually introduced to compensate for the depth decay of the kernel function, the weighted function itself changes the physical nature of the inversion problem. Moreover, as the inversion depth increases, the resolution inevitably decreases gradually, resulting in poor physical property recovery of deep field sources.

[0006] (3) Regularization constraints lead to blurred property boundaries. Classical L2 norm-based regularization inversion methods tend to produce fuzzy models with smooth and continuous changes in property. They can only estimate the approximate location of the source and cannot accurately restore the sharp boundaries, attitude and morphology of the source, which seriously affects the accuracy of geological interpretation.

[0007] (4) Large-scale computing resources are required. Especially in three-dimensional physical property inversion, a large number of grid cells are often required to accurately depict the morphology of underground anomalies. Solving ultra-large-scale linear equations places extremely high demands on memory usage and computation time, which restricts the efficiency of physical property inversion in actual large-scale exploration areas.

[0008] Gravity and magnetic data interactive inversion (also known as the human-computer interaction selection method) uses visualization technology to display observed anomalies and geological models on the same screen. Interpreters adjust the model's geometry, physical properties, and spatial location using a mouse and keyboard based on geological experience. The system performs real-time forward modeling of theoretical anomalies and displays the fitting results, gradually approximating a reasonable geological interpretation through trial-and-error, feedback-iteration. This method has advantages such as flexibility, adaptability, clear physical meaning, and easy control of the inversion process. However, this method has the following significant drawbacks: (1) It heavily relies on the subjective experience and prior information of the interpreters. The establishment of the initial model, the strategy for parameter adjustment, and the judgment of iterative convergence are highly dependent on the personal experience and skill level of the interpreters. Different interpreters may obtain significantly different model results. In areas without sufficient prior information (such as borehole data, geological maps, seismic profiles, etc.), it is difficult for interpreters to quickly establish a reasonable initial model, resulting in low interpretation efficiency.

[0009] (2) The workflow is cumbersome and the interpretation cycle is long. Interactive inversion is essentially a manual trial and error process. Interpreters need to repeatedly adjust the model parameters and observe the fitting effect. When the morphology of the landmass is complex or there are many landmasses, this process is extremely time-consuming, and the interpretation efficiency is limited by the rhythm and accuracy of manual operation.

[0010] (3) It is prone to getting trapped in local optima and lacks global optimization capabilities. Interactive inversion is a local search strategy. Interpreters tend to make local adjustments near the initial model, making it difficult to effectively explore regions in the solution space far from the initial model and easily missing the globally optimal model solution.

[0011] (4) Difficult to operate in complex three-dimensional scenes. When faced with complex three-dimensional geological models, interactive operations performed on two-dimensional profiles are difficult to fully control the geometric shape and physical property distribution in three-dimensional space, which limits the operability and interpretation quality of interactive inversion.

[0012] For the inversion of 2.5D geological profile gravity and magnetic data, existing methods still mainly rely on single property inversion or interactive inversion modeling, failing to fully leverage the complementary advantages of the two in terms of resolution, boundary characterization, and interpretation efficiency.

[0013] In summary, there is an urgent need for a joint interpretation method for 2.5D gravity and magnetic data that can combine property inversion with interactive inversion to address the shortcomings of existing single methods in terms of resolution, boundary characterization, and interpretation efficiency. Summary of the Invention

[0014] The purpose of this application is to provide a method and system for the joint interpretation of gravity and magnetic data property inversion and interactive inversion, so as to achieve high-precision, high-efficiency, and high-deterministic inversion modeling of subsurface property distribution at the profile scale.

[0015] To achieve the above objectives, this application provides the following solution: Firstly, this application provides a method for the joint interpretation of gravity and magnetic data property inversion and interactive inversion, including: Acquire gravity and magnetic observation data and preprocess them to obtain preprocessed gravity and magnetic observation data; The underground space is initially meshed to obtain multiple initial mesh elements; Based on the multiple initial grid cells and the preprocessed gravity and magnetic observation data, large-scale physical property inversion is performed to obtain the physical property parameters of each initial grid cell; Based on the physical property parameters of each initial grid cell, multiple initial geological-geophysical models are constructed. Starting with the aforementioned initial geological-geophysical models, an interactive inversion is performed to obtain multiple interactive inversion models containing geological body boundary information; Based on the multiple interactive inversion models, the initial grid unit is subdivided a second time to obtain multiple fine grid units; Based on the boundary information of each fine mesh cell and multiple interactive inversion models, the physical property parameters of each fine mesh cell are assigned to obtain the fine mesh model parameters; Based on the multiple fine mesh units and the multiple interactive inversion models, physical property constraints are constructed. Based on the aforementioned physical property constraints, a fine physical property inversion is performed to obtain the fine physical property inversion results; Calculate the relative error between the parameters of the fine mesh model and the results of the fine property inversion; Determine whether the relative error meets the preset convergence condition; If so, the current interactive inversion model and the refined physical property inversion results will be used together as the final interpretation output; If not, the current interactive inversion model is corrected based on the refined property inversion results until the relative error meets the preset convergence condition.

[0016] Secondly, this application provides a joint interpretation system for gravity and magnetic data property inversion and interactive inversion, including: The acquisition and preprocessing module is used to acquire gravity and magnetic observation data and preprocess it to obtain preprocessed gravity and magnetic observation data. The initial mesh generation module is used to perform initial mesh generation on the underground space to obtain multiple initial mesh elements; The large-scale physical property inversion module is used to perform large-scale physical property inversion based on the multiple initial grid cells and the preprocessed gravity and magnetic observation data to obtain the physical property parameters of each initial grid cell. The initial geological-geophysical model construction module is used to construct multiple initial geological-geophysical models based on the physical property parameters of each initial grid cell. The interactive inversion module is used to perform interactive inversion starting from the multiple initial geological-geophysical models to obtain multiple interactive inversion models containing geological body boundary information. The secondary subdivision module is used to perform secondary subdivision on the initial mesh unit based on the multiple interactive inversion models to obtain multiple fine mesh units; The physical property parameter assignment module is used to assign physical property parameters to each of the fine mesh cells based on the boundary information of each fine mesh cell and multiple interactive inversion models, so as to obtain the fine mesh model parameters. The physical property constraint construction module is used to construct physical property constraints based on the multiple fine mesh cells and the multiple interactive inversion models; The fine property inversion module is used to perform fine property inversion based on the property constraints to obtain fine property inversion results; The relative error calculation module is used to calculate the relative error between the parameters of the fine mesh model and the fine property inversion result; The judgment module is used to determine whether the relative error meets the preset convergence condition; if yes, the current interactive inversion model and the fine property inversion result are output as the final interpretation result; if no, the current interactive inversion model is corrected based on the fine property inversion result until the relative error meets the preset convergence condition.

[0017] According to the specific embodiments provided in this application, this application has the following technical effects: This application provides a method and system for the joint interpretation of gravity and magnetic data property inversion and interactive inversion. Through large-grid-scale property inversion, reliable property distribution results are provided for interactive modeling, avoiding the blind trial-and-error process of "finding a needle in a haystack" for interpreters in blank or complex anomalies, effectively reducing modeling time. By constructing property constraints based on the fine-grid model parameters and the interactive inversion model, and performing fine-grained property inversion, fine-grained property inversion results are obtained, effectively overcoming the boundary ambiguity problem caused by regularization smoothing in traditional property inversion. Through an alternating iterative process, the property model continuously converges towards a geological model with clear boundaries, and the final generated property distribution map clearly depicts the boundaries of geological bodies, improving the ability to finely distinguish subsurface structures. Simultaneously, the results of fine-grained property inversion serve as new evidence for verifying and correcting the geological model, reducing the subjective bias of human experience and making the interactive inversion process more objective and directional. By combining the current interactive inversion model with the refined physical property inversion results as the final interpretation output, the same underground structure can be characterized from different angles and using different methods. This inherent mutual verification mechanism provides a more solid and reliable geophysical basis for subsequent geological and mineral evaluation work. Attached Figure Description

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

[0019] Figure 1 This is a flowchart illustrating the alternating joint interpretation method of gravity and magnetic data property inversion and interactive inversion in the embodiments of this application. Figure 2 This is a schematic diagram of the initial grid partitioning for the alternating joint interpretation method of gravity and magnetic data property inversion and interactive inversion in the embodiments of this application. Figure 3 This is a schematic diagram illustrating the construction of multiple initial geological-geophysical models for the alternating joint interpretation method of gravity and magnetic data property inversion and interactive inversion in the embodiments of this application. Detailed Implementation

[0020] The technical solutions of the embodiments of this application will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of this application, and not all embodiments. Based on the embodiments of this application, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of this application.

[0021] To make the above-mentioned objectives, features and advantages of this application more apparent and understandable, the application will be further described in detail below with reference to the accompanying drawings and specific embodiments.

[0022] like Figure 1 As shown, a method for joint interpretation of gravity and magnetic data property inversion and interactive inversion is provided, including the following steps: S1: Acquire gravity and magnetic observation data and preprocess them to obtain preprocessed gravity and magnetic observation data.

[0023] In a specific embodiment, the underground grid range is set according to the size of the study area, and this range should include all observation points. Observations are conducted at the observation points using gravimeters or magnetometers to obtain gravity and magnetic anomaly data, i.e., gravity and magnetic observation data, for interpretation of the observation points. The acquired gravity and magnetic observation data is further preprocessed, including data correction and noise reduction.

[0024] S2: Initially mesh the underground space to obtain multiple initial mesh elements.

[0025] In a specific embodiment, the underground space is divided into large-scale grids. The data point spacing can be set to 2-3 times the minimum target size, thus forming a series of relatively small initial grid cells with relatively large sizes. The initial grid cell regions are then activated, such as... Figure 2 As shown, the solid lines in the upper coordinate system represent the observed data curves, with the horizontal axis X representing the profile distance and the vertical axis F representing the gravity and magnetic observation data. The grid in the lower coordinate system represents the grid division of the inversion model, with the horizontal axis X representing the profile distance, the vertical axis Z representing the profile depth, and the dashed lines representing the topographic lines. Initial grid cells with a "+" sign whose grid center is above the topographic line are inactive, while initial grid cells whose grid center is below the topographic line are active.

[0026] Initial mesh cells in an active state will participate in the large-scale physical property inversion in S3, while those in an inactive state will not. In this application, large-scale physical property inversion refers to using a relatively large initial mesh for inversion. Based on the characteristics of the gravity and magnetic anomaly plane or the anomaly variation characteristics of adjacent parallel profiles, the extension length of the interpretation profile in the direction perpendicular to the profile's orientation (Y-axis) is determined (including both positive and negative directions perpendicular to the profile), and this value is set as the initial value for the model's extension length. The extension length refers to the extension length on both sides of the 2.5D polygon. This parameter will be adjusted in subsequent step S11 when the current interactive inversion model is corrected based on the refined physical property inversion results.

[0027] S3: Based on multiple initial grid cells and preprocessed gravity and magnetic observation data, perform large-scale physical property inversion to obtain the physical property parameters of each initial grid cell.

[0028] In a specific embodiment, it is assumed that the density or magnetic susceptibility is uniform within each initial grid cell. The physical property parameters of each cell are obtained by performing large-scale property inversion using preprocessed gravity and magnetic observation data. The core principle is a linear superposition relationship based on potential field theory: the anomaly generated by each initial grid cell at the observation point is equal to its physical property value multiplied by a kernel function that is only related to the cell shape, position, and observation point coordinates; the total anomaly is the sum of the contributions of all cells. This relationship can be concisely represented in matrix form: ,in It is a data vector, composed of outlier data at the observation points. These are physical property parameters, which are vectors, i.e. For the vector of physical property parameters, It is the kernel matrix.

[0029] In this case, the physical property parameter vector is the only parameter that describes the geological body and is also the only parameter that needs to be inverted. In this way, the problem of inverting the magnetic anomaly source becomes the problem of inverting the underground distribution of physical properties.

[0030] Furthermore, due to the ill-posed nature of gravity and magnetic inversion, regularization constraints must be introduced. The Tikhonov framework is typically used, constructing an objective function that includes the data fit difference and model constraint terms. The model constraints can be selected based on prior information, such as minimum model constraints (favoring simple structures), smoothing constraints (favoring continuous changes in physical properties), or sparsity constraints (favoring blocky boundaries). In addition, to overcome the "skin effect" caused by the rapid decay of the kernel function with depth, a depth weighting factor is introduced to ensure that the inversion results are reasonably distributed across different depths.

[0031] Therefore, based on the above principles, the calculation process for the physical property parameters of each initial mesh element in this application is as follows: ; in, These are physical property parameters. For kernel matrix, superscript Indicates transpose. For regularization parameters, These are the preprocessed gravity and magnetic observations. It is a weighted matrix that includes depth weighting, minimum constraint matrix and model gradient operator.

[0032] S4: Based on the physical property parameters of each initial grid cell, construct multiple initial geological-geophysical models. This includes the following steps: S401: Determine the macroscopic property contour map based on the physical property parameters of each initial grid cell.

[0033] In a specific embodiment, due to the small number of initial grid cells, the calculation speed for large-scale physical property inversion is very fast, facilitating the rapid acquisition of the macroscopic distribution of subsurface physical properties (density or magnetic susceptibility). Thus, based on the physical property parameters of each initial grid cell in step S3, a low-resolution macroscopic physical property contour map that reflects the approximate range, depth, and orientation of the main anomalies is obtained.

[0034] S402: Based on the macroscopic contour map of physical properties, determine the maximum absolute value of the physical property parameter, denoted as... and based on Determine the first preset threshold.

[0035] In a specific embodiment, the first preset threshold is: .

[0036] S403: Use the contour lines in the macroscopic physical property contour map where the absolute value of the physical property parameter is the first preset threshold as the boundary of the first major anomaly.

[0037] In a specific embodiment, the absolute value of the physical property parameter is traced on the interpretive cross-section. The contour lines are used as the boundaries of the first major anomaly.

[0038] S404: The average value of the physical property parameters of all initial mesh elements whose absolute values ​​of physical property parameters within the boundary of the first primary anomaly are greater than a first preset threshold is taken as the physical property value of the first primary anomaly. For example... Figure 3 The boundary defined by all the bold black lines in the middle is shown, with the x-axis... Represents the cross-sectional distance, ordinate Indicates depth.

[0039] S405: Based on the boundary and physical property values ​​of the first major anomaly, construct the first geological-geophysical model; the outline of the first geological-geophysical model is the boundary of the first major anomaly.

[0040] S406: The maximum absolute value of the physical property parameter of the initial mesh element that is not located within the boundary of the first major anomaly and has local physical property anomalies is denoted as... and based on Determine the second preset threshold.

[0041] In a specific embodiment, the second preset threshold is: .

[0042] S407: Use the contour lines in the macroscopic physical property contour map where the absolute value of the physical property parameter is the second preset threshold as the boundary of the next major anomaly.

[0043] S408: Take the average value of the physical parameters of all initial grid cells whose absolute values ​​of physical parameters within the boundary of the next major anomaly are greater than the second preset threshold as the physical value of the next major anomaly.

[0044] S409: Construct the next geological-geophysical model based on the boundary and physical property values ​​of the next major anomaly.

[0045] S4010: Determine whether there are still local anomaly regions on the macroscopic physical property contour map; local anomaly regions such as Figure 3 The boundary is defined by all the bold black dotted lines in the middle.

[0046] S4011: If yes, return to step S406 and repeat steps S406 to S4010 for the local physical property anomaly area until there is no longer a local physical property anomaly area on the macroscopic physical property contour map. S4012: If not, then the remaining area on the macroscopic property contour map is determined as the background area; the average property parameter of all initial grid cells located in the remaining area is taken as the property value of the background area; S4013: Construct a background geological-geophysical model based on the background region and its physical properties; S4014: The first geological-geophysical model, all subsequent geological-geophysical models, and the background geological-geophysical model are used together as the initial geological-geophysical model.

[0047] S5: Starting with multiple initial geological-geophysical models, perform interactive inversion to obtain multiple interactive inversion models containing geological body boundary information. Specifically, this includes the following steps: S501: Based on the boundary and physical property values ​​of multiple initial geological-geophysical models, construct multiple initial interactive inversion models.

[0048] In a specific embodiment, each interactive inversion model includes multiple 2.5D polygons, each corresponding to a geological body. The corner coordinates of these 2.5D polygons are obtained by selecting corner points at a certain distance along the boundaries of each initial geological-geophysical model. Corner points can be densified at locations with greater boundary curvature. The physical properties and lateral extension lengths of the 2.5D polygons are consistent with the settings of the corresponding initial geological-geophysical models.

[0049] S502: Calculate the theoretical outliers for each initial interactive inversion model.

[0050] In a specific embodiment, based on the forward modeling formula of 2.5D polygon gravity and magnetic anomaly, the gravity and magnetic anomaly generated by each 2.5D polygon at each gravity and magnetic observation point is calculated, and the forward modeling results of all polygons are summed to obtain the theoretical anomaly value of each initial interactive inversion model.

[0051] S503: Calculate the fitting error based on the theoretical outliers and the preprocessed gravity and magnetic observations. The calculation formula is as follows: ; in, For fitting error, These are the preprocessed gravity and magnetic observations. These are the theoretical outliers of the initial interactive inversion model. The smaller the value, the better the fit.

[0052] S504: Determine whether the fitting error is less than or equal to the third preset threshold.

[0053] In a specific embodiment, the third preset threshold is preset based on empirical values.

[0054] S505: If so, then the current initial interaction inversion model shall be used as the final interaction inversion model.

[0055] S506: If not, adjust the initial interactive inversion model and recalculate the fitting error until the characteristics of each initial interactive inversion model match the descriptive characteristics of the geological data and the fitting error is less than or equal to the third preset threshold.

[0056] In a specific embodiment, the interpreter adjusts the corner positions, physical property parameters, spatial positions, and lateral extension lengths of the 2.5D polygon in the initial interactive inversion model based on the fitting relationship between the theoretical and observed anomaly curves of the initial interactive inversion model, as well as the stratigraphic attitude and morphology information described in the geological data. After adjusting each parameter, the theoretical anomaly value of each initial interactive inversion model is calculated again using the 2.5D polygon gravity and magnetic anomaly forward modeling formula. If the fitting improves, the parameter is adjusted in that direction; otherwise, it is adjusted in the opposite direction. Through multiple trials, the optimal values ​​for each parameter are determined. The parameters of the 2.5D polygon are modified repeatedly until the 2.5D polygon features match the descriptive features of the geological data, and the fitting error is minimized. Less than the third preset threshold.

[0057] Compared with the interactive modeling that starts from scratch in existing technologies, step S5 above has significantly improved modeling efficiency and accuracy because it has multiple initial geological-geophysical models obtained from large-scale physical property inversion as a starting point.

[0058] S6: Based on multiple interactive inversion models, the initial grid cells are subdivided a second time to obtain multiple fine grid cells.

[0059] In a specific embodiment, the region defined by the boundary of the interactive inversion model obtained in step S5 is identified as the key research area requiring high-resolution fine characterization. The initial grid cells set in the previous property inversion are then subdivided to form a large number of smaller grid cells to achieve high-resolution imaging. If the original grid size is already less than 0.5 times the data point spacing or the minimum target size, the grid size remains unchanged.

[0060] S7: Based on the boundary information of each fine mesh element and multiple interactive inversion models, assign physical property parameters to each fine mesh element to obtain the fine mesh model parameters; specifically including the following steps: S701: Determine whether the center point of each fine mesh cell is located inside a polygon in a certain interactive inversion model.

[0061] In a specific embodiment, the method for determining which interactive inversion model each fine mesh cell belongs to in step S5 is as follows: iteratively determine whether the center point of each fine mesh cell is located inside the polygon of a certain interactive inversion model.

[0062] S702: If so, then based on the interactive inversion model corresponding to the 2.5D polygon, assign the corresponding physical property parameters to each fine mesh element to obtain the fine mesh model parameters, denoted as... .

[0063] S703: If not, then keep the original physical properties of the fine mesh cells unchanged.

[0064] S8: Construct physical property constraints based on multiple fine-grained mesh elements and multiple interactive inversion models. This includes the following steps: S801: Calculate the minimum distance from the center point of each fine mesh cell to the boundary of each polygon in the interactive inversion model; S802: Determine whether the minimum distance is less than or equal to the fourth preset threshold; In a specific embodiment, the fourth preset threshold is preset based on empirical values.

[0065] S803: If so, then define the fine mesh element as the boundary mesh element and set the gradient constraint weight coefficient of the boundary mesh element to 0; S804: If not, then define the fine mesh element as a non-boundary mesh element and set the gradient constraint weight coefficient of the non-boundary mesh element to 1; S805: Construct a constraint matrix based on the gradient constraint weight coefficients of each fine mesh element; denoted as: ;in, For the constraint matrix, These are the gradient constraint weight coefficients. This diagonal constraint matrix constitutes the constructed physical property constraints, highlighting boundary characteristics.

[0066] S9: Based on the material property constraints, perform a refined material property inversion to obtain the refined material property inversion results. This includes the following steps: S901: Use fine mesh model parameters as the initial model and reference model ; S902: Based on the initial model, reference model, constraint matrix, and preprocessed gravity and magnetic observation data, calculate the model correction amount using the following formula: ; in, This is the model correction amount. For kernel matrix, superscript Indicates transpose. These are the preprocessed gravity and magnetic observations. For the first The model parameters for the next iteration, when hour, That is, the initial model ; For regularization parameters, This is the minimum constraint matrix for the model. As a reference model, For the model gradient operator, Represents the model gradient operator Three directions, This is the constraint matrix.

[0067] S903: Update the current model parameters based on the model correction amount. The calculation formula is as follows: ; in, For the first Model parameters for the next iteration; The step size.

[0068] S904: Determine whether the norm of the model correction is less than or equal to the fifth preset threshold or whether the maximum number of iterations has been reached; In a specific embodiment, the fifth preset threshold is preset based on empirical values.

[0069] S905: If so, stop the iteration and use the current model parameters at the current iteration number as the result of the fine-grained property inversion. ; S906: If not, continue iterating until the norm of the model correction is less than the fifth preset threshold or the maximum number of iterations is reached.

[0070] After the above step S9, the fine property inversion has strong prior constraints, the multiple solutions of the inversion are greatly suppressed, and the obtained fine property inversion results have clear boundaries and reasonable property distribution.

[0071] S10: Calculate the relative error between the parameters of the fine mesh model and the results of the fine property inversion. The calculation formula is as follows: ; in, This is a relative error. For the results of fine-grained property inversion, These are the parameters for the fine mesh model.

[0072] S11: Determine whether the relative error meets the preset convergence condition. ;where the threshold It is a relatively small positive value.

[0073] If so, the current interactive inversion model (2.5D polygon) and the refined property inversion results will be used together as the final interpretation output.

[0074] The final interpretation results at this point encompass both the clear geological significance and well-defined boundaries embodied in the interactive inversion model, and the continuous property distribution information presented by the detailed property inversion results. The results of the interactive inversion model can be used for geological profile drawing, borehole layout guidance, and resource estimation; while the detailed property inversion results can identify property changes within concealed anomalies, compensate for the overly "polygonal" boundaries of the geological model, and reveal gradual transition zones of properties. The two models corroborate each other, significantly improving the reliability of the final interpretation results.

[0075] If not, the current interactive inversion model is corrected based on the refined property inversion results until the relative error meets the preset convergence condition.

[0076] In a specific embodiment, the refined property inversion results obtained in step S9 are used as a reference base map for a new round of interactive inversion, guiding interpreters to modify and optimize the current interactive inversion model a second, third, or even more times, generating an iteratively optimized interactive inversion model. Subsequently, based on the reference base map again, the geometry, property parameters, extension length, and spatial position of the current interactive inversion model are adjusted to obtain a corrected interactive inversion model, realizing an alternating iterative and feedback correction closed-loop process of "property inversion → interactive inversion → property constraint inversion".

[0077] This application also provides a joint interpretation system for alternating inversion and interactive inversion of gravity and magnetic data properties, including: The acquisition and preprocessing module is used to acquire gravity and magnetic observation data and preprocess it to obtain preprocessed gravity and magnetic observation data. The initial mesh generation module is used to perform initial mesh generation on the underground space to obtain multiple initial mesh elements; The large-scale physical property inversion module is used to perform large-scale physical property inversion based on the multiple initial grid cells and the preprocessed gravity and magnetic observation data to obtain the physical property parameters of each initial grid cell. The initial geological-geophysical model construction module is used to construct multiple initial geological-geophysical models based on the physical property parameters of each initial grid cell. The interactive inversion module is used to perform interactive inversion starting from the multiple initial geological-geophysical models to obtain multiple interactive inversion models containing geological body boundary information. The secondary subdivision module is used to perform secondary subdivision on the initial mesh unit based on the multiple interactive inversion models to obtain multiple fine mesh units; The physical property parameter assignment module is used to assign physical property parameters to each of the fine mesh cells based on the boundary information of each fine mesh cell and multiple interactive inversion models, so as to obtain the fine mesh model parameters. The physical property constraint construction module is used to construct physical property constraints based on the multiple fine mesh cells and the multiple interactive inversion models; The fine property inversion module is used to perform fine property inversion based on the property constraints to obtain fine property inversion results; The relative error calculation module is used to calculate the relative error between the parameters of the fine mesh model and the fine property inversion result; The judgment module is used to determine whether the relative error meets the preset convergence condition; if yes, the current interactive inversion model and the fine property inversion result are output as the final interpretation result; if no, the current interactive inversion model is corrected based on the fine property inversion result until the relative error meets the preset convergence condition.

[0078] This document uses specific examples to illustrate the principles and implementation methods of this application. The descriptions of the above embodiments are only for the purpose of helping to understand the methods and core ideas of this application. At the same time, for those skilled in the art, there will be changes in the specific implementation methods and application scope based on the ideas of this application. In summary, the content of this specification should not be construed as a limitation of this application.

Claims

1. A method for alternating and joint interpretation of gravity and magnetic data property inversion and interactive inversion, characterized in that, The method for alternating joint interpretation of gravity and magnetic data property inversion and interactive inversion includes: Acquire gravity and magnetic observation data and preprocess them to obtain preprocessed gravity and magnetic observation data; The underground space is initially meshed to obtain multiple initial mesh elements; Based on the multiple initial grid cells and the preprocessed gravity and magnetic observation data, large-scale physical property inversion is performed to obtain the physical property parameters of each initial grid cell; Based on the physical property parameters of each initial grid cell, multiple initial geological-geophysical models are constructed. Starting with the aforementioned initial geological-geophysical models, an interactive inversion is performed to obtain multiple interactive inversion models containing geological body boundary information; Based on the multiple interactive inversion models, the initial grid unit is subdivided a second time to obtain multiple fine grid units; Based on the boundary information of each fine mesh cell and multiple interactive inversion models, the physical property parameters of each fine mesh cell are assigned to obtain the fine mesh model parameters; Based on the multiple fine mesh units and the multiple interactive inversion models, physical property constraints are constructed. Based on the aforementioned physical property constraints, a fine physical property inversion is performed to obtain the fine physical property inversion results; Calculate the relative error between the parameters of the fine mesh model and the results of the fine property inversion; Determine whether the relative error meets the preset convergence condition; If so, the current interactive inversion model and the refined physical property inversion results will be used together as the final interpretation output; If not, the current interactive inversion model is corrected based on the refined property inversion results until the relative error meets the preset convergence condition.

2. The method for alternating joint interpretation of gravity and magnetic data property inversion and interactive inversion according to claim 1, characterized in that, Based on the multiple initial grid cells and the preprocessed gravity and magnetic observation data, a large-scale physical property inversion is performed to obtain the physical property parameters of each initial grid cell. The calculation process is as follows: ； in, These are physical property parameters. For kernel matrix, superscript Indicates transpose. For regularization parameters, These are the preprocessed gravity and magnetic observations. It is a weighted matrix that includes depth weighting, minimum constraint matrix and model gradient operator.

3. The method for alternating joint interpretation of gravity and magnetic data property inversion and interactive inversion according to claim 1, characterized in that, Based on the physical property parameters of each initial grid cell, multiple initial geological-geophysical models are constructed, specifically including the following steps: S401: Based on the physical property parameters of each of the initial grid cells, determine the macroscopic physical property contour map; S402: Based on the aforementioned macroscopic physical property contour map, determine the maximum absolute value of the physical property parameter, denoted as... and based on the Determine the first preset threshold; S403: Use the contour lines in the macroscopic physical property contour map whose absolute values ​​of physical property parameters are first preset thresholds as the boundaries of the first major anomaly. S404: The average value of the physical property parameters of all initial grid cells whose absolute values ​​of physical property parameters within the boundary of the first main anomaly are greater than the first preset threshold is taken as the physical property value of the first main anomaly. S405: Based on the boundary of the first major anomaly and the physical property values ​​of the first major anomaly, construct the first geological-geophysical model; S406: The maximum absolute value of the physical property parameter of the initial mesh element that is not located within the boundary of the first major anomaly and has local physical property anomalies is denoted as... and based on the Determine the second preset threshold; S407: Use the contour lines in the macroscopic physical property contour map whose absolute values ​​of physical property parameters are second preset thresholds as the boundaries of the next major anomalies; S408: The average value of the physical property parameters of all initial grid cells whose absolute values ​​of physical property parameters within the boundary of the next major anomaly are greater than the second preset threshold is taken as the physical property value of the next major anomaly. S409: Based on the boundary of the next major anomaly and the physical property values ​​of the next major anomaly, construct the next geological-geophysical model; S4010: Determine whether there are still local abnormal physical property areas on the macroscopic physical property contour map; S4011: If so, return to step S406 and repeat steps S406 to S4010 for the local physical property anomaly area until there is no longer a local physical property anomaly area on the macroscopic physical property contour map. S4012: If not, then the remaining area on the macroscopic physical property contour map is determined as the background area; the average value of the physical property parameters of all initial grid cells located in the remaining area is taken as the physical property value of the background area; S4013: Based on the background area and its physical properties, construct a background geological-geophysical model; S4014: The first geological-geophysical model, all subsequent geological-geophysical models, and the background geological-geophysical model are used together as the initial geological-geophysical model.

4. The method for alternating joint interpretation of gravity and magnetic data property inversion and interactive inversion according to claim 1, characterized in that, Starting with the aforementioned initial geological-geophysical models, an interactive inversion is performed to obtain multiple interactive inversion models containing geological body boundary information; specifically, the following steps are included: Based on the boundary and physical property values ​​of the multiple initial geological-geophysical models, multiple initial interactive inversion models are constructed; each of the interactive inversion models includes multiple 2.5D polygons, and each 2.5D polygon corresponds to a geological body; Calculate the theoretical outliers for each of the initial interactive inversion models; The fitting error is calculated based on the theoretical outliers and the preprocessed gravity and magnetic observations. Determine whether the fitting error is less than or equal to a third preset threshold; If so, then the current initial interaction inversion model shall be used as the interaction inversion model; If not, the initial interactive inversion model is adjusted and the fitting error is recalculated until the features of each initial interactive inversion model match the descriptive features of the geological data, and the fitting error is less than or equal to the third preset threshold.

5. The method for alternating joint interpretation of gravity and magnetic data property inversion and interactive inversion according to claim 4, characterized in that, The formula for calculating the fitting error is as follows: ； in, For fitting error, These are the preprocessed gravity and magnetic observations. These are the theoretical outliers of the initial interactive inversion model.

6. The method for alternating joint interpretation of gravity and magnetic data property inversion and interactive inversion according to claim 4, characterized in that, Based on the boundary information of each fine mesh element and multiple interactive inversion models, the physical property parameters of each fine mesh element are assigned to obtain the fine mesh model parameters. This specifically includes the following steps: Determine whether the center point of each fine mesh unit is located inside a polygon in a certain interactive inversion model; If so, then according to the interactive inversion model corresponding to the 2.5D polygon, the corresponding physical property parameters are assigned to each fine mesh unit to obtain the fine mesh model parameters; If not, then the original physical properties of the fine mesh unit remain unchanged.

7. The method for alternating joint interpretation of gravity and magnetic data property inversion and interactive inversion according to claim 6, characterized in that, Based on the multiple fine mesh elements and the multiple interactive inversion models, physical property constraints are constructed, specifically including the following steps: Calculate the minimum distance from the center point of each fine mesh cell to the boundary of each 2.5D polygon in the interactive inversion model; Determine whether the minimum distance is less than or equal to a fourth preset threshold; If so, the fine mesh cell is determined as the boundary mesh cell, and the gradient constraint weight coefficient of the boundary mesh cell is set to 0; If not, the fine mesh cell is determined as a non-boundary mesh cell, and the gradient constraint weight coefficient of the non-boundary mesh cell is set to 1; Based on the gradient constraint weight coefficients of each fine grid cell, a constraint matrix is ​​constructed; the constraint matrix is ​​the physical property constraint condition.

8. The method for alternating joint interpretation of gravity and magnetic data property inversion and interactive inversion according to claim 7, characterized in that, Based on the aforementioned property constraints, a refined property inversion is performed to obtain the refined property inversion result, specifically including the following steps: The parameters of the fine mesh model are used as the initial model and the reference model; Based on the initial model, the reference model, the constraint matrix, and the preprocessed gravity and magnetic observation data, the model correction is calculated using the following formula: ； in, This is the model correction amount. For kernel matrix, superscript Indicates transpose. These are the preprocessed gravity and magnetic observations. For the first The model parameters for the next iteration, when hour, That is, the initial model ; For regularization parameters, This is the minimum constraint matrix for the model. As a reference model, For the model gradient operator, Represents the model gradient operator Three directions, This is the constraint matrix; The current model parameters are updated based on the model correction amount, calculated as follows: ； in, For the first Model parameters for the next iteration Step size; Determine whether the norm of the model correction is less than or equal to the fifth preset threshold or whether the maximum number of iterations has been reached; If so, stop the iteration and use the current model parameters at the current iteration number as the result of the fine-grained property inversion; If not, continue iterating until the norm of the model correction is less than the fifth preset threshold or the maximum number of iterations is reached.

9. The method for alternating joint interpretation of gravity and magnetic data property inversion and interactive inversion according to claim 8, characterized in that, The formula for calculating the relative error between the parameters of the fine mesh model and the fine property inversion result is as follows: ； in, This is a relative error. For the results of fine-grained property inversion, These are the parameters for the fine mesh model.

10. A joint interpretation system for alternating inversion and interactive inversion of gravity and magnetic data properties, characterized in that, The alternating joint interpretation system for gravity and magnetic data property inversion and interactive inversion includes: The acquisition and preprocessing module is used to acquire gravity and magnetic observation data and preprocess it to obtain preprocessed gravity and magnetic observation data. The initial mesh generation module is used to perform initial mesh generation on the underground space to obtain multiple initial mesh elements; The large-scale physical property inversion module is used to perform large-scale physical property inversion based on the multiple initial grid cells and the preprocessed gravity and magnetic observation data to obtain the physical property parameters of each initial grid cell. The initial geological-geophysical model construction module is used to construct multiple initial geological-geophysical models based on the physical property parameters of each initial grid cell. The interactive inversion module is used to perform interactive inversion starting from the multiple initial geological-geophysical models to obtain multiple interactive inversion models containing geological body boundary information. The secondary subdivision module is used to perform secondary subdivision on the initial mesh unit based on the multiple interactive inversion models to obtain multiple fine mesh units; The physical property parameter assignment module is used to assign physical property parameters to each of the fine mesh cells based on the boundary information of each fine mesh cell and multiple interactive inversion models, so as to obtain the fine mesh model parameters. The physical property constraint construction module is used to construct physical property constraints based on the multiple fine mesh cells and the multiple interactive inversion models; The fine property inversion module is used to perform fine property inversion based on the property constraints to obtain fine property inversion results; The relative error calculation module is used to calculate the relative error between the parameters of the fine mesh model and the fine property inversion result; The judgment module is used to determine whether the relative error meets the preset convergence condition; if yes, the current interactive inversion model and the fine property inversion result are output as the final interpretation result; if no, the current interactive inversion model is corrected based on the fine property inversion result until the relative error meets the preset convergence condition.

Images