An airborne gravity and magnetic imaging method suitable for multi-scale crustal structure analysis

Abstract

The application provides an airborne gravity and magnetic imaging method suitable for multi-scale crust structure analysis, which comprises the following steps: obtaining preprocessed research area reduction-to-pole total field data, reduction-to-pole vector component data and unified airborne gravity Bouguer gravity data; using a method based on the second-order spectrum moment of a potential field surface to extract a fracture structure, thereby obtaining a two-dimensional fracture structure based on gravity and magnetic scalar data; using a method based on the eigenvalue of a vector component covariance matrix to extract a fracture structure, thereby obtaining a two-dimensional fracture structure based on airborne magnetic vector data. A three-dimensional inversion method based on remanent magnetization intelligent discrimination is used to obtain a three-dimensional magnetic susceptibility result; a three-dimensional fracture structure and a comprehensive output of multi-scale structure imaging are obtained. The multi-scale three-dimensional imaging of the application completely describes the transition from a two-dimensional plane to a three-dimensional solid and the transition from a magnetic structure to a fracture structure, and can be applied to deep mineral resource exploration and ore prospecting target delineation, mineral resource quantity evaluation and other deep ore exploration and metallogenic geological applications.

Description

Technical Field

[0001] This invention relates to the fields of geophysics and geochemistry, and in particular to an airborne gravity and magnetic imaging method suitable for multi-scale crustal structure analysis. Background Technology

[0002] Airborne gravity and magnetic surveys, as an important tool in geophysical exploration, have shown unique advantages in deep mineralization prediction and target area selection by capturing differences in the physical properties of underground geological bodies and interpreting tectonic constraints.

[0003] In recent years, aeromagnetic vector measurement technology has gradually matured, and three-component data can provide richer structural anisotropy information. However, existing aeromagnetic and gravity data processing technologies still have shortcomings in system integration and adaptive remanent magnetization. First, methods such as two-dimensional boundary identification, three-dimensional property inversion, and three-dimensional fracture structure extraction exist mostly as individual tools, lacking a system integration scheme that organically connects them into a complete workflow. Second, under remanent magnetization conditions, traditional three-dimensional inversion methods often assume that induction magnetization is dominant, ignoring the influence of remanent magnetization, leading to distorted inversion results; and existing remanent magnetization application methods lack a data-driven objective judgment mechanism, with strategy switching relying entirely on human experience, making it difficult to guarantee the stability of the results. Finally, the potential of vector data is not fully explored. Aeromagnetic vector data contains information on the total magnetization direction and fine structural anisotropy, but existing interpretation methods mostly remain at the level of scalar data analysis, failing to fully utilize the multi-dimensional geological information brought by high-precision measurements. Summary of the Invention

[0004] This invention provides an airborne gravity and magnetic imaging method suitable for multi-scale crustal structure analysis. It comprehensively describes the study area from two-dimensional plane to three-dimensional solid, and from magnetic structure to fracture structure through multi-scale three-dimensional imaging. It can serve the deep mineral resource exploration and prospecting target area delineation, spatial distribution inversion of concealed ore-forming rock masses, interpretation of ore-controlling and ore-guiding structures, and mineral resource quantity evaluation, among other deep mineral exploration and metallogenic geology applications.

[0005] In a first aspect, embodiments of the present invention provide an airborne gravity and magnetic imaging method suitable for multi-scale crustal structure analysis, comprising:

[0006] S1, acquire the preprocessed aeromagnetic total field data, aeromagnetic vector component data, and unified aeromagnetic Bouguer gravity data of the study area after polarization;

[0007] S2, For the aeromagnetic total field data after polarization and the unified aeromagnetic-Bouguer gravity data, the fracture structure is extracted using the method based on the second-order spectral moment of the potential field surface to obtain a two-dimensional fracture structure based on gravity and magnetic scalar data.

[0008] For the aeromagnetic vector component data after polarization, the fracture structure is extracted using a method based on the eigenvalues ​​of the vector component covariance matrix, resulting in a two-dimensional fracture structure based on the aeromagnetic vector data.

[0009] Further, step S2 includes:

[0010] S211, calculate the second-order spectral moment tensor of the unified airborne weight Bouguer gravity data, and calculate two surface statistical invariants:

[0011] ;

[0012] ;

[0013] ;

[0014] in, The unified air weight Bouguer gravity data are collectively referred to as the potential field surface; For potential field surface The second-order spectral moment tensor; Represents the potential field gradient vector; Represents the tensor cross product; Represents a sliding window The mean operator; , , Represents tensor elements; This represents the first surface statistical invariant, reflecting the total intensity of the surface slope variance; This represents the statistical invariant of the second surface, reflecting the coupling and anisotropy of the slope in different directions;

[0015] S212, calculate the two-dimensional fracture structure based on gravity and magnetic scalar data according to the first surface statistical invariant and the second surface statistical invariant;

[0016] ;

[0017] in, To reflect the first result parameters of the two-dimensional fracture structure based on gravity and magnetic scalar data, Indicates the sign function, extracting the positive or negative sign of a numerical value;

[0018] Furthermore, the linearly extending regions correspond to faults or stratigraphic boundaries. Furthermore, the closed, ring-shaped area at both ends corresponds to the boundary of the rock mass or geological body.

[0019] Furthermore, step S2 further includes:

[0020] S221, based on the aeromagnetic vector component data after pole conversion, construct the second-order covariance matrix:

[0021] ;

[0022] in, express The second-order covariance matrix is ​​a real symmetric positive semi-definite matrix, which centrally represents the sliding window. Second-order statistical characteristics of aeromagnetic vector component data after internalization; This represents a sliding window, with a size of [size missing]. ; Represents a sliding window Summing over all grid points in the inner region; and All are positive integers; , and The X, Y, and Z components of the aeromagnetic vector component data after polarization;

[0023] S222, Perform eigenvalue decomposition on the second-order covariance matrix to obtain its diagonalized form and eigenvectors:

[0024] ;

[0025] in, Let be the three eigenvalues ​​of the second-order covariance matrix arranged in descending order; , , For each corresponding , , Eigenvectors with three eigenvalues Indicates transpose;

[0026] S223, Construction , , First-order and second-order elementary symmetric polynomials with three eigenvalues ​​are used as two fundamental invariants reflecting field intensity and anisotropy. Based on these two fundamental invariants, two-dimensional fracture structures based on aeromagnetic vector data are calculated.

[0027] Construct parameters that reflect the fracture structure:

[0028] ;

[0029] ;

[0030] ;

[0031] in, To reflect the second result parameter of the two-dimensional fracture structure based on aeromagnetic vector data, The field strength is an invariant, representing the sliding window. The total energy of the internal vector field; It is a second-order symmetric polynomial, representing a vector field;

[0032] Furthermore, areas extending in a linear pattern correspond to fractures or stratigraphic boundaries; Furthermore, the closed, ring-shaped area at both ends corresponds to the boundary of the rock mass or geological body.

[0033] Further, step S1 includes:

[0034] Initial aerial survey data were collected for each survey area within the study area. The aerial survey data for each survey area were then fused, spliced, and normalized to the target altitude to obtain the integrated aerial survey data for each survey area.

[0035] The integrated aerial survey data of each survey area is automatically stitched together using a stitching or hybrid method to obtain unified aerial survey data corresponding to the study area. The unified aerial survey data includes unified aeromagnetic total field data, unified aeromagnetic vector component data, and unified aeromagnetic Bouguer gravity data.

[0036] Using the unified aeromagnetic vector component data, the total magnetization direction is calculated to obtain the total magnetization tilt angle and total magnetization deflection angle of the target magnetic source in the study area;

[0037] Based on the total magnetization tilt angle and the total magnetization deflection angle, a true magnetization direction polarization operator in the frequency domain is constructed;

[0038] Based on the true magnetization direction polarization operator, the unified aeromagnetic total field data and the unified aeromagnetic vector component data are subjected to synchronous polarization processing to obtain the polarized aeromagnetic total field data and the polarized aeromagnetic vector component data.

[0039] Furthermore, step S2 is followed by:

[0040] S3. Based on the total magnetization direction and the regional geomagnetic direction, the remanent magnetization intensity of the study area is obtained, wherein the total magnetization direction includes the total magnetization tilt angle and the total magnetization deflection angle;

[0041] If the remanence intensity is less than the preset remanence threshold, then based on the total aeromagnetic field data after pole conversion, a three-dimensional focusing inversion based on Tikhonov regularization is used to obtain the three-dimensional magnetic susceptibility structure under weak magnetic conditions.

[0042] If the remanent magnetization intensity is greater than the preset remanent magnetization threshold, then based on the aeromagnetic total field data after pole conversion, the magnetic anomaly modulus comprehensive inversion with the normalized magnetic source intensity inversion result as the a priori model is adopted to obtain the three-dimensional magnetic susceptibility structure under strong magnetic conditions.

[0043] Furthermore, the aeromagnetic total field data after pole conversion is used to obtain the three-dimensional magnetic susceptibility structure under weak magnetic conditions through three-dimensional focusing inversion based on Tikhonov regularization. The steps include:

[0044] The objective function under the inversion of weak magnetic field conditions is:

[0045] ;

[0046] in, Represents the objective function under weak magnetic conditions; Represents a data weighting matrix; Forward operands; The three-dimensional magnetic susceptibility structure under the weak magnetic condition to be solved; , representing the total aeromagnetic field data after polarization; Represents the model weighting matrix; It is an adaptive regularization factor;

[0047] The objective function under weak magnetic conditions is solved iteratively by using the conjugate gradient method, resulting in the three-dimensional magnetic susceptibility structure under weak magnetic conditions.

[0048] Furthermore, based on the aeromagnetic total field data after pole conversion, a comprehensive inversion of the magnetic anomaly modulus using the normalized magnetic source intensity inversion result as a priori model is employed to obtain the three-dimensional magnetic susceptibility structure under strong magnetic conditions. The steps include:

[0049] For the three components of the aeromagnetic total field data after pole conversion, the spatial partial derivatives are calculated to obtain the magnetic gradient tensor matrix. Then, the magnetic gradient tensor matrix is ​​decomposed into three eigenvalues ​​to obtain three eigenvalues. Finally, the normalized magnetic source intensity is calculated using the three eigenvalues.

[0050] ;

[0051] in, This represents the normalized magnetic source strength; , and These represent the three eigenvalues ​​of the magnetic gradient tensor matrix;

[0052] Using the normalized magnetic source intensity as the observed value, a depth-weighted focusing three-dimensional inversion target functional is constructed. By introducing a depth weighting matrix, an adaptive regularization factor is adopted and iteratively solved using the conjugate gradient method to obtain the shallow three-dimensional magnetic structure.

[0053] Based on the total aeromagnetic field data after pole conversion, the magnetic anomaly modulus was calculated, and the shallow three-dimensional magnetic structure was used as a priori model to construct the objective function under strong magnetic conditions:

[0054] ;

[0055] ;

[0056] in, Represents the objective function under strong magnetic conditions; The three-dimensional magnetic susceptibility structure under strong magnetic conditions is to be solved. Represents a data weighting matrix; express Nonlinear forward modeling operator; To observe the magnetic anomaly modulus; It is an adaptive regularization factor; This represents the depth-weighted matrix; This refers to the shallow three-dimensional magnetic structure;

[0057] This represents the magnetic anomaly modulus. The northward component of the aeromagnetic total field data after polarization is represented in the frequency domain. The eastward component of the aeromagnetic total field data after polarization is represented in the frequency domain. The vertical component of the aeromagnetic total field data after pole conversion is represented by the frequency domain conversion.

[0058] The objective function under strong magnetic conditions is solved iteratively by using the conjugate gradient method, resulting in the three-dimensional magnetic susceptibility structure under strong magnetic conditions.

[0059] Furthermore, step S3 is followed by:

[0060] S4. Based on the three-dimensional magnetic susceptibility structure, the three-dimensional fracture structure of the study area is extracted using a method based on three-dimensional spectral moments and statistical invariants.

[0061] Further, step S4 includes:

[0062] S41, in the three-dimensional magnetic susceptibility structure, spherical sliding volume elements slide point by point, and at each volume element position, the three-dimensional second-order spectral moment tensor of the three-dimensional material volume is defined using tensor notation:

[0063] , ;

[0064] in, This represents the second-order spectral moment tensor. Indicated in spherical element The arithmetic mean within, Represents the physical field along The variance of the directional partial derivative, Spatial coupling reflecting partial derivatives Let be the gradient vector of the three-dimensional magnetic body. For tensor operators;

[0065] S42, perform eigenvalue decomposition on the second-order spectral moment tensor to obtain , and Three eigenvalues ​​are used to construct a three-dimensional local intensity. With anisotropy parameter :

[0066] ;

[0067] ;

[0068] in, It reflects the total intensity of the change in the three-dimensional physical properties of the body along three directions within the spherical element; It reflects the degree of unevenness in changes in each direction, i.e., the degree of anisotropy;

[0069] S43, using the three-dimensional local intensity and the anisotropy parameter as new input data, perform three-dimensional spectral moment and statistical invariant analysis again within the spherical volume element to obtain second-order statistical invariants. with anisotropy Based on this, the three-dimensional boundary coefficients that ultimately reflect the three-dimensional fracture structure are constructed:

[0070] ;

[0071] in, The value range is from -1 to +1, selected The values ​​are extracted using isosurface extraction, and the continuously extending undulating surfaces correspond to the three-dimensional extension of the fracture.

[0072] Secondly, embodiments of the present invention provide an airborne gravity and magnetic imaging system suitable for multi-scale crustal structure analysis, comprising:

[0073] The data acquisition module is used to acquire the preprocessed aeromagnetic total field data, aeromagnetic vector component data, and unified aeromagnetic Bouguer gravity data of the study area after polarization.

[0074] The two-dimensional fracture module is used to extract fracture structures from the aeromagnetic total field data after polarization and the unified aeromagnetic-Bougu gravity data using a method based on the second-order spectral moment of the potential field surface, and to obtain two-dimensional fracture structures based on gravity and magnetic scalar data.

[0075] For the aeromagnetic vector component data after polarization, the fracture structure is extracted using the method based on the eigenvalue of the vector component covariance matrix, resulting in a two-dimensional fracture structure based on the aeromagnetic vector data.

[0076] A three-dimensional magnetic susceptibility module is used to obtain the remanent magnetization of the study area based on the total magnetization direction and the regional geomagnetic direction, wherein the total magnetization direction includes the total magnetization tilt angle and the total magnetization deflection angle;

[0077] If the remanence intensity is less than the preset remanence threshold, then based on the total aeromagnetic field data after pole conversion, a three-dimensional focusing inversion based on Tikhonov regularization is used to obtain the three-dimensional magnetic susceptibility structure under weak magnetic conditions.

[0078] If the remanent magnetization intensity is greater than the preset remanent magnetization threshold, then based on the aeromagnetic total field data after pole conversion, the magnetic anomaly modulus comprehensive inversion with the normalized magnetic source intensity inversion result as the a priori model is adopted to obtain the three-dimensional magnetic susceptibility structure under strong magnetic conditions.

[0079] The three-dimensional fracture module is used to extract the three-dimensional fracture structure of the study area based on the three-dimensional magnetic susceptibility structure and a method based on three-dimensional spectral moments and statistical invariants.

[0080] Thirdly, embodiments of the present invention provide a computer device, including a memory, a processor, and a computer program stored in the memory and executable on the processor, wherein the processor executes the computer program to implement the steps of the above-described airborne gravity and magnetic imaging method applicable to multi-scale crustal structure analysis.

[0081] Fourthly, embodiments of the present invention provide a computer storage medium storing a computer program, which, when executed by a processor, implements the steps of the above-described airborne gravity and magnetic imaging method applicable to multi-scale crustal structure analysis.

[0082] This invention proposes an airborne gravity and magnetic imaging method suitable for multi-scale crustal structure analysis. It introduces a data-driven intelligent remanent magnetization discrimination mechanism, which automatically switches between conventional inversion under weak remanent magnetization conditions and comprehensive inversion of normalized magnetic source intensity and magnetic anomaly modulus under strong remanent magnetization conditions. By processing scalar and vector data in parallel, it provides two-dimensional fracture structures. Finally, it uses three-dimensional spectral moment tensor technology to achieve a three-dimensional and detailed characterization of underground fractures.

[0083] This invention provides a comprehensive multi-scale three-dimensional imaging description of the research area, from two-dimensional plane to three-dimensional solid, and from magnetic structure to fracture structure. It can serve the deep mineral resource exploration and prospecting target area delineation, spatial distribution inversion of concealed ore-forming rock masses, interpretation of ore-controlling and ore-guiding structures, and mineral resource quantity evaluation, among other deep mineral exploration and metallogenic geology applications. Attached Figure Description

[0084] Figure 1 A flowchart of an airborne gravity and magnetic imaging method for multi-scale crustal structure analysis provided in this embodiment of the invention;

[0085] Figure 2 A flowchart illustrating an airborne gravity and magnetic imaging method suitable for multi-scale crustal structure analysis, provided as an embodiment of the present invention;

[0086] Figure 3 This is a schematic diagram of an airborne gravity and magnetic imaging system for multi-scale crustal structural analysis, provided as an embodiment of the present invention.

[0087] The realization of the objective, functional features and advantages of the present invention will be further explained in conjunction with the embodiments and with reference to the accompanying drawings. Detailed Implementation

[0088] The embodiments of this application are described in detail below. Examples of the embodiments are shown in the accompanying drawings, wherein the same or similar reference numerals denote the same or similar elements or elements having the same or similar functions throughout. The embodiments described below with reference to the accompanying drawings are exemplary and are only used to explain this application, and should not be construed as limiting this application.

[0089] To enable those skilled in the art to better understand the solutions of this application, 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 a part of the embodiments of this application, and not all of them. All other embodiments obtained by those skilled in the art based on the embodiments of this application without creative effort are within the scope of protection of this application.

[0090] In the embodiments of this application, "at least one" refers to one or more; "multiple" refers to two or more. In the description of this application, terms such as "first," "second," and "third" are used only for descriptive purposes and should not be construed as indicating or implying relative importance or order. Furthermore, the terms "first" and "second" are used for descriptive purposes only and should not be construed as indicating or implying relative importance or implicitly specifying the number of indicated technical features. Therefore, a feature defined as "first" or "second" may explicitly or implicitly include one or more of that feature. In the description of this application, "multiple" means two or more, unless otherwise explicitly specified.

[0091] References such as “one embodiment” or “some embodiments” as described in this specification mean that one or more embodiments of this application include a specific feature, structure, or characteristic described in connection with that embodiment. Therefore, the terms “comprising,” “including,” “having,” and variations thereof, as used in this specification, mean “including, but not limited to,” unless otherwise specifically emphasized.

[0092] Figure 1A flowchart of an airborne gravity and magnetic imaging method for multi-scale crustal structure analysis provided in this embodiment of the invention is shown below. Figure 1 As shown, the method includes:

[0093] S1, acquire the preprocessed aeromagnetic total field data, aeromagnetic vector component data, and unified aeromagnetic Bouguer gravity data of the study area after polarization;

[0094] Specifically, step S1 includes:

[0095] S11: Collect initial aerial survey data for each survey area within the study area, and then fuse, stitch together, and normalize the aerial survey data of each survey area to the target altitude to obtain integrated aerial survey data for each survey area.

[0096] In this embodiment of the invention, initial aerial survey data for each survey area within the study area are first collected. The study area is the overall target area for this exploration and geological analysis, the complete geographical range to be fully imaged by the project, and the boundary of the entire work area. Dividing the entire study area into different survey areas is to facilitate field aerial surveying.

[0097] Initial aerial survey data includes initial total aeromagnetic field data. Initial aeromagnetic vector component data and initial airborne weight Bouguer gravity data The initial aeromagnetic vector component data includes three components. , and The initial aeromagnetic vector component data can be directly obtained by measuring with a three-component magnetometer; the initial total aeromagnetic field data... Hehang Heavy Bouger Gravity Data Data can be obtained through national open data platforms, project field measurements, or aerial surveys, and is in the format of ASCII grid files or raster data.

[0098] There are three major problems in each survey area: (1) there are obvious boundary steps and amplitude jumps between survey areas; (2) there are slight differences in flight altitude and correction parameters between survey areas; (3) a single grid cannot completely cover the entire study area, and it is impossible to carry out full-area three-dimensional inversion and full-area fault extraction when used alone.

[0099] Therefore, it is necessary to first perform fusion and splicing processing within each survey area, and for data at different measurement altitudes, use frequency domain extension to normalize them to the target altitude, thus obtaining the integrated aerial survey data for each survey area.

[0100] Specifically, the flight altitudes for different aerial surveys vary considerably: some fly at 500m, others at 800m. Gravity and magnetic anomalies decrease with observation altitude; at different altitudes, the anomaly morphology and amplitude are completely different, making direct comparison and joint calculation impossible. The operation method is as follows:

[0101] The initial airborne weight Bouguer gravity data were separately extended in the frequency domain and uniformly corrected to the set standard altitude; the initial aeromagnetic total field data were also separately extended in the frequency domain and similarly corrected to the same standard altitude; the extension was performed independently for the gravity and magnetic data sets, without mixing them.

[0102] S12, using a stitching or hybrid method, the integrated aerial survey data of each survey area are automatically stitched together to obtain unified aerial survey data corresponding to the study area. The unified aerial survey data includes unified aeromagnetic total field data, unified aeromagnetic vector component data, and unified aeromagnetic Bouguer gravity data.

[0103] After integrating the data within each test area, it is found that, due to the different test areas belonging to different sorties and different construction cycles, in addition to flight altitude, there are multiple layers of systematic errors that cannot be eliminated by extension:

[0104] Instrument zero drift: The baseline values ​​of the gravimeter and magnetometer show a fixed offset in different expeditions;

[0105] Differences in external correction parameters: Inconsistencies in geomagnetic field diurnal variation correction and gravity baseline connection benchmarks at different times;

[0106] The selection criteria for terrain correction and Bouguer correction parameters differ: there are subtle differences in parameter settings among personnel processing different survey areas; upward extension of the frequency domain only corrects abnormal attenuation caused by observation height, and cannot correct the overall amplitude shift caused by the instrument and calibration benchmark. Even if both survey areas are standardized to the 800m standard height, the overall amplitude of survey area A is 5nT higher and that of survey area B is 3nT lower, with a clear step appearing at the junction.

[0107] Therefore, automatic stitching between different survey areas is also required. In this embodiment of the invention, the integrated aerial survey data of each survey area are automatically stitched together using the Grid Knitting module of the Oasis Montaj geophysical data processing software, employing either a stitching method or a hybrid method to obtain unified aerial survey data covering the study area.

[0108] It should be noted that the unified aerial survey data includes unified aeromagnetic total field data, unified aeromagnetic vector component data, and unified airborne boogüter gravity data. During calculations, these three data sets are integrated both within and outside the survey area.

[0109] S13, using the unified aeromagnetic vector component data, the total magnetization direction is calculated to obtain the total magnetization tilt angle and total magnetization deflection angle of the target magnetic source in the study area;

[0110] In this embodiment of the invention, in order to accurately locate the magnetic source and eliminate the influence of oblique magnetization, it is necessary to estimate the total magnetization direction using unified aeromagnetic vector component data:

[0111] (1)

[0112] in, As a correlation indicator; The current sliding window; For the first in the sliding window The mean-removed observations of the vector magnitude anomalies of the points; , These represent the magnetization tilt angle and the magnetization deflection angle, respectively. For the first in the sliding window The point is the calculated value after the direction transformation and demeaning, under the current assumed total magnetization direction.

[0113] when The closer it gets to 1, the more it corresponds to That is, the total magnetization tilt angle and the total magnetization deflection angle corresponding to the magnetic source.

[0114] S14, Based on the total magnetization tilt angle and the total magnetization deflection angle, construct the true magnetization direction polarization operator in the frequency domain;

[0115] Based on the optimal total magnetization tilt angle and deflection angle, a true magnetization direction polarization reversal operator can be constructed in the frequency domain:

[0116] (2)

[0117] in, , Wavenumber domain coordinates; ; The direction of the regional geomagnetic field can be obtained from the IGRF international geomagnetic reference field or geomagnetic reference stations. This operator converts the magnetic field affected by oblique magnetization into a vertical magnetization field by redistributing the energy and phase in the frequency domain.

[0118] S15, based on the true magnetization direction polarization operator, synchronous polarization processing is performed on the unified aeromagnetic total field data and the unified aeromagnetic vector component data respectively to obtain polarized aeromagnetic total field data and polarized aeromagnetic vector component data.

[0119] Finally, the unified total aeromagnetic field number and unified aeromagnetic vector component data are synchronized and polarized as follows:

[0120] (3)

[0121] in, Fourier transform; This is the inverse Fourier transform.

[0122] This step synchronously corrects the morphology of anomalies in each aeromagnetic component, eliminating positional distortion caused by magnetization direction deviation. (Aeromagnetic total field data after pole shifting) And the aeromagnetic vector component data after polarization ( , , Initial flight weight Bouguer gravity data As the data input basis for subsequent two-dimensional fracture structure extraction and physical property structure inversion steps, the total magnetization tilt angle and total magnetization deflection angle of the magnetic source are used. This forms the basis for intelligent discrimination in the physical property structure inversion step.

[0123] S2, For the aeromagnetic total field data after polarization and the unified aeromagnetic-Bouguer gravity data, the fracture structure is extracted using the method based on the second-order spectral moment of the potential field surface to obtain a two-dimensional fracture structure based on gravity and magnetic scalar data.

[0124] For the aeromagnetic vector component data after polarization, the fracture structure is extracted using a method based on the eigenvalues ​​of the vector component covariance matrix, resulting in a two-dimensional fracture structure based on the aeromagnetic vector data.

[0125] In this embodiment of the invention, two-dimensional fracture structures are extracted from the aeromagnetic total field data and the aeromagnetic vector component data after pole polarization using different methods. Specifically, the fracture structure extraction from the aeromagnetic total field data after pole polarization is performed using a method based on the second-order spectral moments of the potential field surface; the fracture structure extraction from the aeromagnetic vector component data after pole polarization is performed using a method based on the eigenvalues ​​of the vector component covariance matrix. Each type of data outputs a planar distribution of the two-dimensional fracture structure in the study area.

[0126] In one implementation, step S2 includes:

[0127] S211, calculate the second-order spectral moment tensor of the unified airborne weight Bouguer gravity data, and calculate two surface statistical invariants:

[0128] (4)

[0129] (5)

[0130] in, The unified air weight Bouguer gravity data are collectively referred to as the potential field surface; For potential field surface The second-order spectral moment tensor; Represents the potential field gradient vector; Represents the tensor cross product; Represents a sliding window The mean operator; , , Represents tensor elements; This represents the first surface statistical invariant, reflecting the total intensity of the surface slope variance; This represents the statistical invariant of the second surface, reflecting the coupling and anisotropy of the slope in different directions;

[0131] The aeromagnetic total field data after pole-forming in this embodiment of the invention With Unified Airborne Weight Bouguer Gravity Data Fracture structure extraction is performed. First, the potential field surface is used... As input, in a size of sliding window Inside, the potential field surface is defined according to formula (4). The second-order spectral moment tensor is used to calculate two surface statistical invariants according to formula (5).

[0132] S212, calculate the two-dimensional fracture structure based on gravity and magnetic scalar data according to the first surface statistical invariant and the second surface statistical invariant;

[0133] (6)

[0134] in, To reflect the first result parameters of the two-dimensional fracture structure based on gravity and magnetic scalar data, Indicates the sign function, extracting the positive or negative sign of a numerical value;

[0135] Furthermore, the linearly extending regions correspond to faults or stratigraphic boundaries. Furthermore, the closed, ring-shaped area at both ends corresponds to the boundary of the rock mass or geological body.

[0136] It should be noted that, The value ranges from -1 to +1, and the entire study area is obtained during the point-by-point sliding process of the sliding window. Distribution. This step can be based on... Output a two-dimensional fracture structure distribution map of gravity and magnetic scalar data.

[0137] In another implementation, step S2 further includes:

[0138] S221, based on the aeromagnetic vector component data after pole conversion, construct the second-order covariance matrix:

[0139] (7)

[0140] in, express The second-order covariance matrix is ​​a real symmetric positive semi-definite matrix, which centrally represents the sliding window. Second-order statistical characteristics of aeromagnetic vector component data after internalization; This represents a sliding window, with a size of [size missing]. ; Represents a sliding window Summing over all grid points in the inner region; and All are positive integers; , and The X, Y, and Z components of the aeromagnetic vector component data after polarization;

[0141] S222, Perform eigenvalue decomposition on the second-order covariance matrix to obtain its diagonalized form and eigenvectors:

[0142] (8)

[0143] in, Let be the three eigenvalues ​​of the second-order covariance matrix arranged in descending order; , , For each corresponding , , Eigenvectors with three eigenvalues Indicates transpose;

[0144] S223, Construction , , First-order and second-order elementary symmetric polynomials with three eigenvalues ​​are used as two fundamental invariants reflecting field intensity and anisotropy. Based on these two fundamental invariants, two-dimensional fracture structures based on aeromagnetic vector data are calculated.

[0145] Construct parameters that reflect the fracture structure:

[0146] (9)

[0147] (10)

[0148] in, To reflect the second result parameter of the two-dimensional fracture structure based on aeromagnetic vector data, The electric field strength is an invariant, representing the sliding window. The total energy of the internal vector field; It is a second-order symmetric polynomial, representing a vector field;

[0149] Furthermore, areas extending in a linear pattern correspond to fractures or stratigraphic boundaries; Furthermore, the closed, ring-shaped area at both ends corresponds to the boundary of the rock mass or geological body.

[0150] It should be noted that, The value ranges from -1 to +1, essentially measuring the relative proportion of the principal direction energy at that location in the total three-dimensional energy. The entire study area is obtained during the point-by-point sliding of the window. distributed. Furthermore, areas extending in a linear pattern correspond to fractures or stratigraphic boundaries; Furthermore, the closed, ring-shaped area at both ends corresponds to the boundary of the rock mass or geological body.

[0151] In mathematics, when hour, Degenerate into , Degenerate into Formula (7) strictly degenerates into a two-dimensional fracture parameter based on two eigenvalues, ensuring complete compatibility of the high-dimensional form with low-dimensional data; while when hour, All the vertical anisotropy information carried is included in the normalized denominator. This method demonstrates a significant enhancement in identifying targets such as gently dipping fractures and detachment structures that exhibit clear signals only in the vertical component, a feature not found in previous methods. This stage can output a fracture location identification map of the study area. .

[0152] In summary, based on the output and The distribution map shows the planar distribution of two-dimensional fracture structures in the study area.

[0153] It should also be noted that two-dimensional fracture structures refer to planar fracture information, which only expresses the planar features of the subsurface structure projected horizontally onto the surface and does not include depth information. Extraction is based on gridded planar gravity and magnetic anomalies (polarity anomalies, Bouguer gravity planar grids), and the output is planar fracture lines and fracture structure planar maps, recording only:

[0154] Fracture plane distribution: strike, length, planar segmentation, bifurcation, and intersection.

[0155] Fault plane boundary: tectonic boundary line, stratigraphic boundary line.

[0156] Fracture plane strength: The strength of the gravity and magnetic gradient reflects the degree of fracture.

[0157] For example, the entire fault does not distinguish whether it is shallow (100m) or deep (5km), and only has a horizontal planar shape, so it is called a two-dimensional fault structure.

[0158] Figure 2 A flowchart of an airborne gravity and magnetic imaging method for multi-scale crustal structure analysis provided in an embodiment of the present invention is included. In some embodiments, step S2 is followed by:

[0159] S3. Based on the total magnetization direction and the regional geomagnetic direction, the remanent magnetization intensity of the study area is obtained, wherein the total magnetization direction includes the total magnetization tilt angle and the total magnetization deflection angle;

[0160] If the remanence intensity is less than the preset remanence threshold, then based on the total aeromagnetic field data after pole conversion, a three-dimensional focusing inversion based on Tikhonov regularization is used to obtain the three-dimensional magnetic susceptibility structure under weak magnetic conditions.

[0161] If the remanent magnetization intensity is greater than the preset remanent magnetization threshold, then based on the aeromagnetic total field data after pole conversion, the magnetic anomaly modulus comprehensive inversion with the normalized magnetic source intensity inversion result as the prior model is adopted to obtain the three-dimensional magnetic susceptibility structure under strong magnetic conditions.

[0162] In this embodiment of the invention, the angle deviation is calculated based on the estimated total magnetization direction and regional geomagnetic direction, and the remanent magnetization intensity of the study area is automatically determined.

[0163] Specifically, for each sliding window center grid point Calculate the total magnetization direction relative to the direction of the regional geomagnetic field (Angle deviation given by the IGRF International Geomagnetic Reference Field) :

[0164] (11)

[0165] Physically, This characterizes the combined deflection degree of the remanent magnetization component and the induced magnetization component of the rock. If... This indicates that the magnetization in the study area is mainly contributed by induced magnetization; with As the remanence increases, its contribution to the total magnetization increases significantly.

[0166] Calculation of the entire study area Then, the average of the entire space was used. As a quantitative proxy indicator, a preset residual magnetism threshold is set. (Recommended initial value is) ),when When the condition is determined to be weak remanence, based on the total aeromagnetic field data after pole transformation, a three-dimensional focusing inversion based on Tikhonov regularization is used to obtain the three-dimensional magnetic susceptibility structure under weak magnetic conditions; when When the condition is determined to be a strong remanent magnetization condition, the magnetic anomaly modulus comprehensive inversion with the normalized magnetic source intensity inversion result as the a priori model is adopted to obtain the three-dimensional magnetic susceptibility structure under the strong magnetic condition.

[0167] As one implementation method, the aeromagnetic total field data after pole conversion is used to obtain the three-dimensional magnetic susceptibility structure under weak magnetic conditions through three-dimensional focusing inversion based on Tikhonov regularization. The steps include:

[0168] The objective function under the inversion of weak magnetic field conditions is:

[0169] (12)

[0170] in, Represents the objective function under weak magnetic conditions; Represents a data weighting matrix; Forward operands; The three-dimensional magnetic susceptibility structure under the weak magnetic condition to be solved; , representing the total aeromagnetic field data after polarization; Represents the model weighting matrix; It is an adaptive regularization factor;

[0171] The objective function under weak magnetic conditions is solved iteratively by using the conjugate gradient method, resulting in the three-dimensional magnetic susceptibility structure under weak magnetic conditions.

[0172] As another implementation method, the aeromagnetic total field data after pole conversion is used to obtain the three-dimensional magnetic susceptibility structure under strong magnetic conditions by comprehensive inversion of the magnetic anomaly modulus with the normalized magnetic source intensity inversion result as the a priori model. The steps include:

[0173] For the three components of the aeromagnetic total field data after polarization, the spatial partial derivatives are calculated to obtain the magnetic gradient tensor matrix. Then, the magnetic gradient tensor matrix is ​​decomposed into three eigenvalues, which are used to calculate the normalized magnetic source intensity.

[0174] (13)

[0175] in, This represents the normalized magnetic source strength; , and These represent the three eigenvalues ​​of the magnetic gradient tensor matrix;

[0176] It should be noted that, in the embodiments of the present invention, ,and It is inversely proportional to the fourth power of the distance, weakly sensitive to the magnetization direction, and completely unaffected by the magnetization direction in the case of a dipole field. It is the ideal conversion quantity for inverting shallow magnetic bodies.

[0177] Using the normalized magnetic source intensity as the observed value, a depth-weighted focusing three-dimensional inversion target functional is constructed. By introducing a depth weighting matrix, an adaptive regularization factor is adopted and iteratively solved using the conjugate gradient method to obtain the shallow three-dimensional magnetic structure.

[0178] Based on the total aeromagnetic field data after pole conversion, the magnetic anomaly modulus was calculated, and the shallow three-dimensional magnetic structure was used as a priori model to construct the objective function under strong magnetic conditions:

[0179] (14)

[0180] (15)

[0181] in, Represents the objective function under strong magnetic conditions; The three-dimensional magnetic susceptibility structure under strong magnetic conditions is to be solved. Represents a data weighting matrix; express Nonlinear forward modeling operator; To observe the magnetic anomaly modulus; It is an adaptive regularization factor; This represents the depth-weighted matrix; This refers to the shallow three-dimensional magnetic structure;

[0182] This represents the magnetic anomaly modulus. The northward component of the aeromagnetic total field data after polarization is represented in the frequency domain. The eastward component of the aeromagnetic total field data after polarization is represented in the frequency domain. The vertical component of the aeromagnetic total field data after pole conversion is represented by the frequency domain conversion.

[0183] The objective function under strong magnetic conditions is solved iteratively by using the conjugate gradient method, resulting in the three-dimensional magnetic susceptibility structure under strong magnetic conditions.

[0184] in, .

[0185] Furthermore, the magnetic anomaly modulus It is another type of magnetic conversion quantity that is weakly sensitive to the magnetization direction. It is inversely proportional to the cube of the distance and of the same order of magnitude as the magnetic anomaly. It does not significantly amplify noise during the conversion process and can effectively preserve long-wavelength low-frequency information, thus providing a reliable data foundation for the inversion of deep magnetic structures.

[0186] It should be noted that the three-dimensional magnetic susceptibility structure is a three-dimensional distribution model of the magnetic susceptibility of each geological body unit in the three-dimensional space (X, Y, Z) of the entire underground area. At the same time, it combines the results of remanent magnetization discrimination to distinguish between induced magnetism and remanent magnetism, and to three-dimensionally depict the spatial occurrence state of underground magnetic materials.

[0187] A simple breakdown of the three-dimensional magnetic susceptibility structure into its three dimensions:

[0188] X, Y: Lateral range in the plane (the entire study area, without division into survey areas);

[0189] Z: Vertical depth, dividing the earth's surface into multiple depth slices; The model discretizes the underground space into a large number of three-dimensional cubic meshes (three-dimensional voxels), each voxel is assigned a magnetic susceptibility value, and all voxels are combined to form a three-dimensional magnetic structure.

[0190] See Figure 2 In some embodiments, step S3 is followed by:

[0191] S4. Based on the three-dimensional magnetic susceptibility structure, the three-dimensional fracture structure of the study area is extracted using a method based on three-dimensional spectral moments and statistical invariants.

[0192] Specifically, step S4 includes:

[0193] S41, in the three-dimensional magnetic susceptibility structure, spherical sliding volume elements slide point by point, and at each volume element position, the three-dimensional second-order spectral moment tensor of the three-dimensional material volume is defined using tensor notation:

[0194] , (16)

[0195] in, This represents the second-order spectral moment tensor. Indicated in spherical element The arithmetic mean within, Represents the physical field along The variance of the directional partial derivative, Spatial coupling reflecting partial derivatives Let be the gradient vector of the three-dimensional magnetic body. For tensor operators;

[0196] S42, perform eigenvalue decomposition on the second-order spectral moment tensor to obtain , and Three eigenvalues ​​are used to construct a three-dimensional local intensity. With anisotropy parameter :

[0197] (17)

[0198] in, It reflects the total intensity of the change in the three-dimensional physical properties of the body along three directions within the spherical element; It reflects the degree of unevenness in changes in each direction, i.e., the degree of anisotropy;

[0199] S43, using the three-dimensional local intensity and the anisotropy parameter as new input data, perform three-dimensional spectral moment and statistical invariant analysis again within the spherical volume element to obtain second-order statistical invariants. with anisotropy Based on this, the three-dimensional boundary coefficients that ultimately reflect the three-dimensional fracture structure are constructed:

[0200] (18)

[0201] in, The value range is from -1 to +1, selected The values ​​are extracted using isosurface extraction, and the continuously extending undulating surfaces correspond to the three-dimensional extension of the fracture.

[0202] It should be noted that, compared to two-dimensional fracture structures, three-dimensional fracture structures possess depth information.

[0203] This invention also integrates the products of the aforementioned steps as the final output, specifically including the following three types of output:

[0204] (1) Multidimensional fracture structure planar distribution. Two-dimensional fracture parameters extracted from the multi-source data output in step S2. and The distribution map reflects the distribution characteristics of the fault and geological body boundary on the shallow surface.

[0205] The three-dimensional boundary coefficients output from step S4 The study area is presented with three-dimensional fracture extension, and the two-dimensional planar fracture and the three-dimensional fracture together constitute the multi-dimensional fracture structure of the study area.

[0206] (2) Multi-scale three-dimensional magnetic structure. The three-dimensional magnetic susceptibility structure of the study area output from step S3. The weak remanent magnetization branch output Obtained by focused inversion; strong remanent magnetization branch output or Obtained by comprehensive inversion of normalized magnetic source intensity (NSS) and magnetic anomaly modulus (Ta), the shallow part is obtained by... The inversion results are anchored with high spatial resolution, and the deep layers are... Inversion preserves long-wavelength information, enabling multi-scale magnetic structure characterization that integrates shallow and deep layers.

[0207] (3) Comprehensive output format. The comprehensive output of this invention is presented in the form of three types of datasets: "two-dimensional fracture plan view + three-dimensional magnetic structure + three-dimensional fracture structure". These datasets can be displayed in geophysical data visualization software in various forms such as plan view, depth slice, profile, and three-dimensional volume rendering, for use in subsequent geological interpretation and exploration deployment.

[0208] This invention provides a comprehensive multi-scale three-dimensional imaging description of the research area, from two-dimensional plane to three-dimensional solid, and from magnetic structure to fracture structure. It can serve the deep mineral resource exploration and prospecting target area delineation, spatial distribution inversion of concealed ore-forming rock masses, interpretation of ore-controlling and ore-guiding structures, and mineral resource quantity evaluation, among other deep mineral exploration and metallogenic geology applications.

[0209] This invention proposes a multi-scale structural imaging method for airborne gravity and magnetic data, the innovation of which lies in:

[0210] (1) Innovation Point 1: System integration of multi-module integrated system framework for airborne gravity and magnetic multi-scale imaging. Multiple independent method modules, such as two-dimensional fracture extraction from gravity and magnetic scalar data, two-dimensional fracture extraction from vector aeromagnetic data, intelligent discrimination and inversion branch decision of remanent magnetization, adaptive three-dimensional physical property inversion under remanent magnetization conditions, and three-dimensional fracture structure extraction based on inverted physical property volume, are organically integrated into a complete workflow in a six-step serial structure of "data input → common preprocessing → two-dimensional fracture extraction → physical property inversion (remanent magnetization discrimination + dual branch) → three-dimensional fracture extraction → comprehensive output", forming an airborne gravity and magnetic multi-scale, multi-dimensional imaging system, overcoming the shortcomings of fragmented links in traditional methods.

[0211] (2) Innovation Point 2: Quantitative discrimination criteria and automatic decision-making mechanism for inversion branch of remanence intensity. This is achieved by calculating the deviation angle between the total magnetization direction and the geomagnetic field direction of the vector data. This enables the objective quantification of the impact of remanence. The criterion drives the automatic switching of the inversion strategy: under weak remanence conditions, three-dimensional focusing inversion is invoked; under strong remanence conditions, NSS-Ta integrated inversion is automatically triggered.

[0212] This mechanism effectively avoids the drawbacks of traditional reliance on human experience, significantly improves the automation level, cross-regional adaptability, and consistency of results of the inversion, and is especially suitable for accurate inversion of deep structures in intermediate-acidic magmatic rock development areas.

[0213] (3) Innovation Point 3: Two-dimensional fracture identification using aeromagnetic vector data. Breaking through the limitations of traditional methods that rely on directional derivative calculations and easily amplify high-frequency noise, this method directly utilizes the three-component data after polarization to construct a three-dimensional second-order covariance matrix. This matrix innovatively introduces vertical component information, filling the mathematical gap that traditional two-dimensional frameworks cannot express anisotropic structures in three-dimensional space. The fracture characterization parameters constructed based on this not only perfectly match traditional low-dimensional methods but also significantly enhance the identification accuracy and noise resistance of gently dipping fractures, detachment structures, and geological boundaries with weak vertical gradients.

[0214] (4) Innovation Point 4: Multi-scale three-dimensional magnetic structure inversion method of NSS-Ta integrated inversion. The depth-weighted focusing inversion results of NSS are innovatively used as a priori model and integrated into the Ta inversion objective function in the form of "reference model deviation term" to realize multi-scale inversion and integrated reconstruction of shallow high-frequency information (anchored by NSS) and deep low-frequency information (driven by Ta data).

[0215] The integrated inversion method fully leverages the complementary advantages of NSS (which is weakly sensitive to magnetization direction and can distinguish shallow parts) and Ta (which preserves deep low-frequency signals), overcoming the bottleneck that single-conversion inversion cannot characterize multi-scale magnetic bodies. It can directly support mineral exploration drilling deployment and target area selection.

[0216] Figure 3 This is a schematic diagram of an airborne gravity and magnetic imaging system suitable for multi-scale crustal structure analysis, provided as an embodiment of the present invention. Figure 3 As shown, the system includes:

[0217] Data acquisition module 310 is used to acquire preprocessed aeromagnetic total field data, aeromagnetic vector component data, and unified aeromagnetic Bouguer gravity data of the study area after polarization.

[0218] The two-dimensional fracture module 320 is used to extract the fracture structure from the aeromagnetic total field data after polarization and the unified aeromagnetic-Bougu gravity data using a method based on the second-order spectral moment of the potential field surface, and obtain the two-dimensional fracture structure based on the gravity and magnetic scalar data.

[0219] For the aeromagnetic vector component data after polarization, the fracture structure is extracted using the method based on the eigenvalue of the vector component covariance matrix, resulting in a two-dimensional fracture structure based on the aeromagnetic vector data.

[0220] The three-dimensional magnetic susceptibility module 330 is used to obtain the remanent magnetization of the study area based on the total magnetization direction and the regional geomagnetic direction, wherein the total magnetization direction includes the total magnetization tilt angle and the total magnetization deflection angle;

[0221] If the remanence intensity is less than the preset remanence threshold, then based on the total aeromagnetic field data after pole conversion, a three-dimensional focusing inversion based on Tikhonov regularization is used to obtain the three-dimensional magnetic susceptibility structure under weak magnetic conditions.

[0222] If the remanent magnetization intensity is greater than the preset remanent magnetization threshold, then based on the aeromagnetic total field data after pole conversion, the magnetic anomaly modulus comprehensive inversion with the normalized magnetic source intensity inversion result as the prior model is adopted to obtain the three-dimensional magnetic susceptibility structure under strong magnetic conditions.

[0223] The three-dimensional fracture module 340 is used to extract the three-dimensional fracture structure of the study area based on the three-dimensional magnetic susceptibility structure and a method based on three-dimensional spectral moments and statistical invariants.

[0224] This embodiment is a system embodiment corresponding to the above method embodiment. Its specific implementation process is the same as that of the above method embodiment, and this system embodiment will not be described again here.

[0225] The modules in the aforementioned airborne gravity and magnetic imaging system applicable to multi-scale crustal structural analysis can be implemented entirely or partially through software, hardware, or a combination thereof. These modules can be embedded in or independent of the processor in a computer device, or stored in the computer device's memory as software, so that the processor can call and execute the corresponding operations of each module.

[0226] This invention provides a computer device, which may be a server. The computer device includes a processor, memory, network interface, and database connected via a system bus. The processor provides computational and control capabilities. The memory includes a computer storage medium and internal memory. The computer storage medium stores an operating system, computer programs, and a database. The internal memory provides an environment for the operation of the operating system and computer programs in the computer storage medium. The database stores data generated or acquired during the execution of an airborne gravity and magnetic imaging method suitable for multi-scale crustal structure analysis, such as post-polarization aeromagnetic total field data, post-polarization aeromagnetic vector component data, and unified airborne gravity Bouguer gravity data. The network interface communicates with external terminals via a network connection. When the computer program is executed by the processor, it implements an airborne gravity and magnetic imaging method suitable for multi-scale crustal structure analysis.

[0227] In one embodiment, a computer device is provided, including a memory, a processor, and a computer program stored in the memory and executable on the processor. When the processor executes the computer program, it implements the steps of an airborne gravity and magnetic imaging method suitable for multi-scale crustal structure analysis as described in the above embodiment. Alternatively, when the processor executes the computer program, it implements the functions of various modules / units in this embodiment of an airborne gravity and magnetic imaging system suitable for multi-scale crustal structure analysis.

[0228] In one embodiment, a computer storage medium is provided, on which a computer program is stored. When executed by a processor, the computer program implements the steps of an airborne gravity and magnetic imaging method suitable for multi-scale crustal structure analysis as described in the above embodiment. Alternatively, when executed by a processor, the computer program implements the functions of each module / unit in this embodiment of an airborne gravity and magnetic imaging system suitable for multi-scale crustal structure analysis.

[0229] Those skilled in the art will understand that all or part of the processes in the methods of the above embodiments can be implemented by a computer program instructing related hardware. The computer program can be stored in a non-volatile computer-readable storage medium. When executed, the computer program can include the processes of the embodiments of the above methods. Any references to memory, storage, databases, or other media used in the embodiments provided in this application can include non-volatile and / or volatile memory. Non-volatile memory may include read-only memory (ROM), programmable ROM (PROM), electrically programmable ROM (EPROM), electrically erasable programmable ROM (EEPROM), or flash memory. Volatile memory may include random access memory (RAM) or external cache memory. By way of illustration and not limitation, RAM is available in a variety of forms, such as static RAM (SRAM), dynamic RAM (DRAM), synchronous DRAM (SDRAM), dual data rate SDRAM (DDRSDRAM), enhanced SDRAM (ESDRAM), synchronous link DRAM (SLDRAM), RAMbus direct RAM (RDRAM), direct memory bus dynamic RAM (DRDRAM), and memory bus dynamic RAM (RDRAM), etc.

[0230] Those skilled in the art will clearly understand that, for the sake of convenience and brevity, the above-described division of functional units and modules is used as an example. In practical applications, the above functions can be assigned to different functional units and modules as needed, that is, the internal structure of the device can be divided into different functional units or modules to complete all or part of the functions described above.

[0231] The above-described embodiments are only used to illustrate the technical solutions of the present invention, and are not intended to limit it. Although the present invention has been described in detail with reference to the foregoing embodiments, those skilled in the art should understand that modifications can still be made to the technical solutions described in the foregoing embodiments, or equivalent substitutions can be made to some of the technical features. Such modifications or substitutions do not cause the essence of the corresponding technical solutions to deviate from the spirit and scope of the technical solutions of the embodiments of the present invention, and should all be included within the protection scope of the present invention.

Claims

1. An airborne gravity and magnetic imaging method suitable for multi-scale crustal structure analysis, characterized in that, include: S1, acquire the preprocessed aeromagnetic total field data, aeromagnetic vector component data, and unified aeromagnetic Bouguer gravity data of the study area after polarization; S2, For the aeromagnetic total field data after polarization and the unified aeromagnetic-Bouguer gravity data, the fracture structure is extracted using the method based on the second-order spectral moment of the potential field surface to obtain a two-dimensional fracture structure based on gravity and magnetic scalar data. For the aeromagnetic vector component data after polarization, the fracture structure is extracted using a method based on the eigenvalues ​​of the vector component covariance matrix, resulting in a two-dimensional fracture structure based on the aeromagnetic vector data.

2. The airborne gravity and magnetic imaging method for multi-scale crustal structure analysis according to claim 1, characterized in that, Step S2 includes: S211, calculate the second-order spectral moment tensor of the unified airborne weight Bouguer gravity data, and calculate two surface statistical invariants: ； ； ； in, The unified air weight Bouguer gravity data are collectively referred to as the potential field surface; For potential field surface The second-order spectral moment tensor; Represents the potential field gradient vector; Represents the tensor cross product; Represents a sliding window The mean operator; , , Represents tensor elements; This represents the first surface statistical invariant, reflecting the total intensity of the surface slope variance; This represents the statistical invariant of the second surface, reflecting the coupling and anisotropy of the slope in different directions; S212, calculate the two-dimensional fracture structure based on gravity and magnetic scalar data according to the first surface statistical invariant and the second surface statistical invariant; ； in, To reflect the first result parameters of the two-dimensional fracture structure based on gravity and magnetic scalar data, Indicates the sign function, extracting the positive or negative sign of a numerical value; Furthermore, the linearly extending regions correspond to faults or stratigraphic boundaries. Furthermore, the closed, ring-shaped area at both ends corresponds to the boundary of the rock mass or geological body.

3. The airborne gravity and magnetic imaging method for multi-scale crustal structure analysis according to claim 1, characterized in that, Step S2 further includes: S221, based on the aeromagnetic vector component data after pole conversion, construct the second-order covariance matrix: ； in, express The second-order covariance matrix is ​​a real symmetric positive semi-definite matrix, which centrally represents the sliding window. Second-order statistical characteristics of aeromagnetic vector component data after internalization; This represents a sliding window, with a size of [size missing]. ; Represents a sliding window Summing over all grid points in the inner region; and All are positive integers; , and The X, Y, and Z components of the aeromagnetic vector component data after polarization; S222, Perform eigenvalue decomposition on the second-order covariance matrix to obtain its diagonalized form and eigenvectors: ； in, Let be the three eigenvalues ​​of the second-order covariance matrix arranged in descending order; , , For each corresponding , , Eigenvectors with three eigenvalues Indicates transpose; S223, Construction , , First-order and second-order elementary symmetric polynomials with three eigenvalues ​​are used as two fundamental invariants reflecting field intensity and anisotropy. Based on these two fundamental invariants, two-dimensional fracture structures based on aeromagnetic vector data are calculated. Construct parameters that reflect the fracture structure: ； ； ； in, To reflect the second result parameter of the two-dimensional fracture structure based on aeromagnetic vector data, The electric field strength is an invariant, representing the sliding window. The total energy of the internal vector field; It is a second-order symmetric polynomial, representing a vector field; Furthermore, areas extending in a linear pattern correspond to fractures or stratigraphic boundaries; Furthermore, the closed, ring-shaped area at both ends corresponds to the boundary of the rock mass or geological body.

4. The airborne gravity and magnetic imaging method for multi-scale crustal structure analysis according to claim 1, characterized in that, Step S1 includes: Initial aerial survey data were collected for each survey area within the study area. The initial aerial survey data for each survey area were then fused, spliced, and normalized to the target altitude to obtain the integrated aerial survey data for each survey area. The integrated aerial survey data of each survey area is automatically stitched together using a stitching or hybrid method to obtain unified aerial survey data corresponding to the study area. The unified aerial survey data includes unified aeromagnetic total field data, unified aeromagnetic vector component data, and unified aeromagnetic Bouguer gravity data. Using the unified aeromagnetic vector component data, the total magnetization direction is calculated to obtain the total magnetization tilt angle and total magnetization deflection angle of the target magnetic source in the study area; Based on the total magnetization tilt angle and the total magnetization deflection angle, a true magnetization direction polarization operator in the frequency domain is constructed; Based on the true magnetization direction polarization operator, the unified aeromagnetic total field data and the unified aeromagnetic vector component data are subjected to synchronous polarization processing to obtain the polarized aeromagnetic total field data and the polarized aeromagnetic vector component data.

5. The airborne gravity and magnetic imaging method for multi-scale crustal structure analysis according to claim 4, characterized in that, Step S2 is followed by: S3. Based on the total magnetization direction and the regional geomagnetic direction, the remanent magnetization intensity of the study area is obtained, wherein the total magnetization direction includes the total magnetization tilt angle and the total magnetization deflection angle; If the remanence intensity is less than the preset remanence threshold, then based on the total aeromagnetic field data after pole conversion, a three-dimensional focusing inversion based on Tikhonov regularization is used to obtain the three-dimensional magnetic susceptibility structure under weak magnetic conditions. If the remanent magnetization intensity is greater than the preset remanent magnetization threshold, then based on the aeromagnetic total field data after pole conversion, the magnetic anomaly modulus comprehensive inversion with the normalized magnetic source intensity inversion result as the a priori model is adopted to obtain the three-dimensional magnetic susceptibility structure under strong magnetic conditions.

6. The airborne gravity and magnetic imaging method for multi-scale crustal structure analysis according to claim 5, characterized in that, The aeromagnetic total field data after pole conversion is used to obtain the three-dimensional magnetic susceptibility structure under weak magnetic conditions through three-dimensional focusing inversion based on Tikhonov regularization. The steps include: The objective function under the inversion of weak magnetic field conditions is: ； in, Represents the objective function under weak magnetic conditions; Represents a data weighting matrix; Forward operands; The three-dimensional magnetic susceptibility structure under the weak magnetic condition to be solved; , representing the total aeromagnetic field data after polarization; Represents the model weighting matrix; This is an adaptive regularization factor; The objective function under weak magnetic conditions is solved iteratively by using the conjugate gradient method, resulting in the three-dimensional magnetic susceptibility structure under weak magnetic conditions.

7. The airborne gravity and magnetic imaging method for multi-scale crustal structure analysis according to claim 5, characterized in that, The aeromagnetic total field data after pole conversion is used to obtain the three-dimensional magnetic susceptibility structure under strong magnetic conditions by comprehensive inversion of the magnetic anomaly modulus with the normalized magnetic source intensity inversion result as the a priori model. The steps include: For the three components of the aeromagnetic total field data after polarization, the spatial partial derivatives are calculated to obtain the magnetic gradient tensor matrix. Then, the magnetic gradient tensor matrix is ​​decomposed into three eigenvalues, which are used to calculate the normalized magnetic source intensity. ； in, This represents the normalized magnetic source strength; , and These represent the three eigenvalues ​​of the magnetic gradient tensor matrix; Using the normalized magnetic source intensity as the observed value, a depth-weighted focusing three-dimensional inversion target functional is constructed. By introducing a depth weighting matrix, an adaptive regularization factor is adopted and iteratively solved using the conjugate gradient method to obtain the shallow three-dimensional magnetic structure. Based on the total aeromagnetic field data after pole conversion, the magnetic anomaly modulus was calculated, and the shallow three-dimensional magnetic structure was used as a priori model to construct the objective function under strong magnetic conditions: ； ； in, Represents the objective function under strong magnetic conditions; The three-dimensional magnetic susceptibility structure under strong magnetic conditions is to be solved. Represents a data weighting matrix; express Nonlinear forward modeling operator; To observe the magnetic anomaly modulus; This is an adaptive regularization factor; This represents the depth-weighted matrix; This refers to the shallow three-dimensional magnetic structure; This represents the magnetic anomaly modulus. The northward component of the aeromagnetic total field data after polarization is represented in the frequency domain. The eastward component of the aeromagnetic total field data after polarization is represented in the frequency domain. The vertical component of the aeromagnetic total field data after pole conversion is represented by the frequency domain conversion. The objective function under strong magnetic conditions is solved iteratively by using the conjugate gradient method, resulting in the three-dimensional magnetic susceptibility structure under strong magnetic conditions.

8. The airborne gravity and magnetic imaging method for multi-scale crustal structure analysis according to claim 5, characterized in that, Step S3 is followed by: S4. Based on the three-dimensional magnetic susceptibility structure, the three-dimensional fracture structure of the study area is extracted using a method based on three-dimensional spectral moments and statistical invariants.

9. The airborne gravity and magnetic imaging method for multi-scale crustal structure analysis according to claim 8, characterized in that, Step S4 includes: S41, in the three-dimensional magnetic susceptibility structure, spherical sliding volume elements slide point by point, and at each volume element position, the three-dimensional second-order spectral moment tensor of the three-dimensional material volume is defined using tensor notation: ， ； in, This represents the second-order spectral moment tensor. Indicated in spherical element The arithmetic mean within, Represents the physical field along The variance of the directional partial derivative, Spatial coupling reflecting partial derivatives Let be the gradient vector of the three-dimensional magnetic body. For tensor operators; S42, perform eigenvalue decomposition on the second-order spectral moment tensor to obtain , and Three eigenvalues ​​are used to construct a three-dimensional local intensity. With anisotropy parameter : ； ； in, It reflects the total intensity of the change in the three-dimensional physical properties of the body along three directions within the spherical element; It reflects the degree of unevenness in changes in each direction, i.e., the degree of anisotropy; S43, using the three-dimensional local intensity and the anisotropy parameter as new input data, perform three-dimensional spectral moment and statistical invariant analysis again within the spherical volume element to obtain second-order statistical invariants. with anisotropy Based on this, the three-dimensional boundary coefficients that ultimately reflect the three-dimensional fracture structure are constructed: ； in, The value range is from -1 to +1, selected The values ​​are extracted using isosurface extraction, and the continuously extending undulating surfaces correspond to the three-dimensional extension of the fracture.

Images