A method and system for detecting gravity linear anomalies with all-around coverage
By employing a comprehensive horizontal derivative calculation method combined with two-dimensional fast Fourier transform, the problem of full coverage in gravity linear anomaly detection was solved, enabling detailed analysis of fracture structures and effective identification of hidden fractures, thus improving the application effect of gravity exploration.
Patent Information
- Application Number
- CN202511027946.1
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2025-07-24
- Publication Date
- 2026-04-24
- Estimated Expiration
- 2045-07-24
AI Technical Summary
Existing technologies are insufficient to comprehensively detect linear gravity anomalies, especially in deep fractures where the attenuation effect of the gravity field leads to weak signals, making it difficult to effectively identify and analyze the spatial distribution and activity characteristics of fracture structures.
A comprehensive horizontal derivative calculation method is adopted. By setting the angle change and the upper extension height, and combining two-dimensional fast Fourier transform and inverse transform, the horizontal derivative in each direction is calculated to highlight the linear anomaly characteristics of gravity. The fracture structure is displayed through two result plot methods.
It enables detailed characterization of linear gravity anomalies at different depths and directions, improves the overall understanding and analytical capabilities of fracture structures, and particularly enhances the ability to identify concealed fractures, thus enriching the theory of gravity data processing.
Smart Images

Figure CN120892666B_ABST
Abstract
Description
Technical Field
[0001] This invention relates to the field of gravity exploration technology, and in particular to a calculation method and system for omnidirectional detection of linear gravity anomalies. Background Technology
[0002] Gravity exploration is a classic geophysical exploration method with the most mature theoretical approach, the widest range of applications, and the greatest importance. It determines the spatial location, size, and shape of geological bodies by observing gravity anomalies caused by density differences between the geological body and the surrounding rock, thereby making judgments about geological structures and mineral distribution. In the application of gravity to the study of related geological problems, using regional gravity to study fault structures is one of the most effective geological applications of gravity.
[0003] Faults are one of the most widespread and common structural forms and modes in the Earth's crust. To a certain extent, the history of Earth's tectonic evolution is the history of the occurrence and development of faults. Most large and medium-sized fault structures directly control the distribution and metallogenic types of mineral resources. The study of faults is of great significance for understanding Earth's tectonic evolution, regional geotectonics, and the distribution patterns of oil and gas resources and mineral resources.
[0004] Faults cause geological bodies to shift and dislocate in three-dimensional space, creating abrupt transition zones at stratigraphic density interfaces. The larger the fault, the more significant the density differences between interfaces, the larger the abrupt transition zone, and the more pronounced the gravity anomaly gradient. Linear gravity anomalies or gradient zones are the most basic gravity anomaly characteristics of faults. Years of resource exploration practice have shown that gravity anomalies are a very effective means of studying regional faults.
[0005] Linear gravity anomalies often manifest as gravity gradient zones, i.e., areas of dense gravity isopleths, reflecting abrupt changes in subsurface density boundaries. Lithological differences on either side of a fault zone (such as contact between high-density bedrock and low-density sedimentary layers) or the presence of fracture zones can lead to significant changes in the gravity field gradient, forming linear anomaly zones. Linear gravity anomalies can reveal fault characteristics at different depths. Shallow faults correspond to high-frequency anomalous signals, while deep, large faults exhibit low-frequency gradient zones.
[0006] Gravity linear anomalies, as surface responses to fault structures, provide important evidence for fault identification, activity assessment, and resource exploration through their gradient changes, strike extensions, and anomaly combinations. Gravity linear anomalies are a crucial indicator of subsurface structures in geophysical exploration, and their correlation with fault structures is primarily revealed through changes in gravity field gradients, which demonstrate the spatial distribution, density differences, and activity characteristics of faults.
[0007] Previous researchers have utilized the characteristic that fractures in gravity anomalies mostly manifest as linear anomalies or gradient zones to study various methods for delineating fractures and determining their locations using gravity. Fault zones typically exhibit density difference boundaries, and the directional derivative of gravity represents the gradient change of the gravity potential function along a specific direction. By capturing the gradient extrema in these density abrupt change regions, the fracture location perpendicular to the derivative direction can be highlighted. However, while the directional derivative emphasizes linear and gradient zones, it also causes anomalies without a definite orientation to resemble fractures as linear or gradient zone anomalies, resulting in anomaly distortion. In vertical or wide fracture models, the location of the maximum value of the directional derivative of gravity is basically consistent with the fracture surface projection; for dipping fractures, its maximum value is closer to the fracture top. The directional derivative anomaly is more pronounced in shallow fractures, while deep fractures, due to the attenuation effect of the gravity field, require continuation or filtering processing to enhance the signal. Summary of the Invention
[0008] This invention provides a calculation method and system for omnidirectional detection of linear gravity anomalies. It effectively utilizes multi-directional horizontal derivatives to characterize the linear gravity anomaly features in various directions, providing an effective means for overall fracture classification and analysis of fracture structural features and the relationship between fractures in various directions.
[0009] To achieve the above objectives, the present invention adopts the following technical solution:
[0010] A computational method for detecting linear gravity anomalies with comprehensive coverage includes:
[0011] Step 1: Set the angle change and upper extension height used to calculate the horizontal derivative;
[0012] Step 2: Perform a two-dimensional fast Fourier transform on the original gravity grid data to obtain the original gravity spectrum;
[0013] Step 3: Based on the original gravity spectrum and the upper extension frequency response, calculate the upper extension gravity spectrum corresponding to the upper extension height, and then perform a two-dimensional fast Fourier inverse transform on the upper extension gravity spectrum to obtain the upper extension gravity anomaly.
[0014] Step 4: Based on the angle change, determine the direction angles within the range of 0 to 180 degrees from which the horizontal derivative needs to be calculated, and based on the upward gravity anomaly, calculate the horizontal derivatives corresponding to each direction angle;
[0015] Step 5: Based on the horizontal derivatives corresponding to the angles in each direction, calculate the maximum value among the absolute values of the horizontal derivatives in each direction at the grid point to obtain the calculation result for the first grid point;
[0016] Step 6: Based on the horizontal derivatives corresponding to the angles in each direction, calculate the horizontal derivative value corresponding to the maximum absolute value of the horizontal derivatives in each direction at the grid point, and obtain the calculation result for the second grid point.
[0017] In this manual, step 7: Based on the upper gravity anomaly, the horizontal derivatives corresponding to the angles in each direction, the calculation results of the first grid point and the calculation results of the second grid point, draw and output the corresponding planar diagram.
[0018] In this specification, in step 1, when the set upper extension height is 0, the upper extension gravity anomaly is the original gravity grid data.
[0019] In this specification, in step 3, the extended gravity spectrum is the product of the original gravity spectrum and the extended frequency response.
[0020] In this specification, in step 4, the horizontal derivative in the range of 0 to 180 degrees covers the detection requirements of gravity linear anomalies in all directions, and the horizontal derivative in the range of 180 to 360 degrees is the opposite value of the horizontal derivative corresponding to 0 to 180 degrees.
[0021] In this specification, the calculation results of the first grid point are used to reflect the apparent linear gravity anomaly and the geological body boundary anomaly; the calculation results of the second grid point are used to reflect the hidden linear-like gravity anomaly.
[0022] In this specification, the horizontal derivatives corresponding to each directional angle are calculated as follows: first, the corresponding horizontal derivative frequency response is generated based on the directional angle; then, the horizontal derivative frequency response is multiplied by the upper gravity spectrum to obtain the horizontal derivative spectrum for each direction; finally, a two-dimensional fast Fourier inverse transform is performed on the horizontal derivative spectrum for each direction to obtain the horizontal derivatives corresponding to each directional angle.
[0023] In this manual, the magnitude of the angle change is determined through experiments; the smaller the angle change, the more directional angles whose horizontal derivatives need to be calculated within the range of 0 to 180 degrees, and the higher the directional resolution accuracy.
[0024] In this specification, the maximum angle of each direction is 180 degrees and the difference between ((number of direction angles - 1) multiplied by the angle change) does not exceed the angle value of the angle change.
[0025] A computational system for omnidirectional coverage detection of linear gravity anomalies, employing the computational method for omnidirectional coverage detection of linear gravity anomalies described in any one of the above-mentioned methods; the computational system for omnidirectional coverage detection of linear gravity anomalies includes:
[0026] The setting module is used to set the angle change and upper extension height used to calculate the horizontal derivative;
[0027] The forward transform module is used to perform a two-dimensional fast Fourier forward transform on the original gravity grid data to obtain the original gravity spectrum;
[0028] The inverse transform module is used to calculate the upper gravity spectrum corresponding to the upper extension height based on the original gravity spectrum and the upper extension frequency response, and then perform a two-dimensional fast Fourier inverse transform on the upper extension gravity spectrum to obtain the upper extension gravity anomaly.
[0029] The horizontal derivative calculation module is used to determine the direction angle for calculating the horizontal derivative within the range of 0 to 180 degrees based on the angle change, and to calculate the horizontal derivative corresponding to each direction angle based on the upward gravity anomaly.
[0030] The first grid point calculation module is used to calculate the maximum value of the absolute value of the horizontal derivative in each direction at the grid point based on the horizontal derivative corresponding to the angle in each direction, and to obtain the calculation result of the first grid point.
[0031] The second grid point calculation module is used to calculate the horizontal derivative value corresponding to the maximum absolute value of the horizontal derivative in each direction at the grid point based on the horizontal derivative corresponding to the angle in each direction, and obtain the calculation result of the second grid point.
[0032] In summary, the present invention has at least the following beneficial effects:
[0033] This invention meticulously characterizes linear gravity anomalies at different depths and in different directions through the upward extension of frequency-domain gravity anomalies, the calculation of horizontal derivatives covering all directions, and two unique methods of representing the resulting maps. This provides an innovative processing method for the precise division and analysis of fault structures in different directions. The two innovative mapping methods of this invention integrate linear gravity anomalies (including geological body boundary anomalies) in various directions into a single map, greatly facilitating the overall understanding and analysis of fault structures in different directions. In particular, the second information extraction mapping method further expands the potential for identifying fault structures using gravity, effectively highlighting quasi-linear anomalies hidden within gravity anomalies, providing an effective method and technology for delineating hidden faults using gravity. This invention enriches and develops the theory of gravity data processing. Attached Figure Description
[0034] To more clearly illustrate the technical solutions of the embodiments of the present invention, the drawings used in the following description of the embodiments will be briefly introduced. Obviously, the drawings described below are only some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.
[0035] Figure 1 This is a schematic diagram of the calculation method for omnidirectional coverage detection of linear gravity anomalies involved in this invention.
[0036] Figure 2 This is a schematic diagram of the Bouguer gravity anomaly plan view of the SLPD area involved in this invention.
[0037] Figure 3 This is a schematic diagram of the 0° direction derivative plane of gravity in the SLPD region involved in this invention.
[0038] Figure 4 This is a schematic diagram of the 30° directional derivative plane of gravity in the SLPD region involved in this invention.
[0039] Figure 5 This is a schematic diagram of the 120° directional derivative plane of gravity in the SLPD region involved in this invention.
[0040] Figure 6 This is a schematic diagram of the 150° directional derivative plane of gravity in the SLPD region involved in this invention.
[0041] Figure 7 This is a schematic diagram of the plan view of the processing result (calculation result of the first grid point) of the first method involved in this invention.
[0042] Figure 8 This is a schematic diagram of the plan view of the processing result (second grid point calculation result) of the second method involved in this invention.
[0043] Figure 9 This is a schematic diagram of the SLPD region gravity extension 2km involved in this invention.
[0044] Figure 10 This is a schematic diagram of the 0° directional derivative plane after the gravity of the SLPD region extends 2km upwards, as described in this invention.
[0045] Figure 11 This is a schematic diagram of the 30° directional derivative plane of gravity in the SLPD region, which is extended upwards by 2km as described in this invention.
[0046] Figure 12 This is a schematic diagram of the 120° directional derivative plane of gravity in the SLPD region, which is extended upwards by 2km as described in this invention.
[0047] Figure 13 This is a schematic diagram of the 150° directional derivative plane of gravity in the SLPD region, which is extended upwards by 2km as described in this invention.
[0048] Figure 14 This is a schematic diagram of the plan view of the first method's processing results after the gravity extension of 2km in the SLPD area involved in this invention.
[0049] Figure 15 This is a schematic diagram of the plan view of the second method processing result after the gravity extension of 2km in the SLPD area involved in this invention. Detailed Implementation
[0050] In the following description, only certain exemplary embodiments are briefly described. As those skilled in the art will recognize, the described embodiments can be modified in various ways without departing from the spirit or scope of the embodiments of the invention. Therefore, the drawings and description are considered to be exemplary in nature and not restrictive.
[0051] The following disclosure provides many different implementations or examples for carrying out different structures of the embodiments of the present invention. To simplify the disclosure of the embodiments of the present invention, specific examples of components and arrangements are described below. Of course, these are merely examples and are not intended to limit the embodiments of the present invention. Furthermore, reference numerals and / or reference letters may be repeated in different examples of the embodiments of the present invention; such repetition is for simplification and clarity and does not in itself indicate a relationship between the various implementations and / or arrangements discussed.
[0052] The embodiments of the present invention will now be described in detail with reference to the accompanying drawings.
[0053] like Figure 1 As shown, this embodiment provides a calculation method for omnidirectional coverage detection of linear gravity anomalies, including:
[0054] Step 1: Set the angle change and upper extension height used to calculate the horizontal derivative;
[0055] Step 2: Perform a two-dimensional fast Fourier transform on the original gravity grid data to obtain the original gravity spectrum;
[0056] Step 3: Based on the original gravity spectrum and the upper extension frequency response, calculate the upper extension gravity spectrum corresponding to the upper extension height, and then perform a two-dimensional fast Fourier inverse transform on the upper extension gravity spectrum to obtain the upper extension gravity anomaly.
[0057] Step 4: Based on the angle change, determine the direction angles within the range of 0 to 180 degrees from which the horizontal derivative needs to be calculated, and based on the upward gravity anomaly, calculate the horizontal derivatives corresponding to each direction angle;
[0058] Step 5: Based on the horizontal derivatives corresponding to the angles in each direction, calculate the maximum value among the absolute values of the horizontal derivatives in each direction at the grid point to obtain the calculation result for the first grid point (the result processed by the first method);
[0059] Step 6: Based on the horizontal derivatives corresponding to the angles in each direction, calculate the horizontal derivative value corresponding to the maximum absolute value of the horizontal derivatives in each direction at the grid point, and obtain the calculation result of the second grid point (the result processed by the second method).
[0060] In some embodiments, step 7: Based on the upward gravity anomaly, the horizontal derivatives corresponding to the angles in each direction, the calculation results of the first grid point and the calculation results of the second grid point, draw and output the corresponding planar diagram.
[0061] In some embodiments, in step 1, when the set upper extension height is 0, the upper extension gravity anomaly is the original gravity grid data.
[0062] In some embodiments, in step 3, the extended gravity spectrum is the product of the original gravity spectrum and the extended frequency response.
[0063] In some embodiments, in step 4, the horizontal derivative in the range of 0 to 180 degrees covers the detection requirements of gravity linear anomalies in all directions, and the horizontal derivative in the range of 180 to 360 degrees is the opposite value of the horizontal derivative corresponding to 0 to 180 degrees.
[0064] In some embodiments, the calculation results of the first grid point are used to reflect the apparent linear gravity anomaly and the geological body boundary anomaly; the calculation results of the second grid point are used to reflect the hidden linear-like gravity anomaly.
[0065] In some embodiments, the horizontal derivatives corresponding to each directional angle are calculated as follows: first, the corresponding horizontal derivative frequency response is generated based on the directional angle; then, the horizontal derivative frequency response is multiplied by the upper-extended gravity spectrum to obtain the horizontal derivative spectrum for each direction; finally, a two-dimensional fast Fourier inverse transform is performed on the horizontal derivative spectrum for each direction to obtain the horizontal derivatives corresponding to each directional angle.
[0066] In some embodiments, the magnitude of the angle change is determined experimentally; the smaller the angle change, the more directional angles whose horizontal derivatives need to be calculated within the range of 0 to 180 degrees, and the higher the directional resolution accuracy.
[0067] In some embodiments, the difference between the maximum angle of each direction and ((number of direction angles - 1) multiplied by the angle change) is no more than the angle value of the angle change.
[0068] In some embodiments, the number of directional angles is calculated by dividing 180 by the angle change amount and adding 1. Each directional angle is a series of angles starting from 0 and increasing sequentially according to the angle change amount to the maximum angle.
[0069] In some embodiments, the upper-extended frequency response is an attenuation function generated based on the upper-extended height and frequency parameters.
[0070] In some embodiments, the set upper extension height is used to reflect geological structural features at different depths. The greater the upper extension height, the deeper the underground geological structure it reflects.
[0071] In some embodiments, the generation of the horizontal derivative frequency response is calculated based on the sine and cosine values of the corresponding directional angle to capture the gravitational field gradient change in that direction.
[0072] In some embodiments, the two grid point results are stored in SURFER format so that subsequent mapping software can directly call and plot them.
[0073] In some embodiments, the magnitude of the angle change needs to be determined in conjunction with the scale of the gravity anomaly and the upper extension height. The optimal value can be selected after experimentally calculating the processing effect of different angle changes.
[0074] In some embodiments, the drawn planar diagrams include an upward gravity anomaly planar diagram, horizontal derivative planar diagrams in each direction, a planar diagram corresponding to the first type of grid point results, and a planar diagram corresponding to the second type of grid point results. All types of planar diagrams are used together for the comprehensive analysis of fracture structures.
[0075] In some embodiments, after the horizontal derivatives in each direction are calculated, they need to be stored on disk for subsequent generalization and plotting.
[0076] In some embodiments, the second grid point results can be used to identify concealed fractures, providing a basis for analyzing the spatial distribution and extension characteristics of concealed fractures by highlighting hidden linear gravity anomalies.
[0077] In some embodiments, when performing up-extension processing on gravity data, the up-extension frequency response is calculated by taking the square root of the sum of the squares of the up-extension height and frequency parameters, so as to achieve depth filtering of gravity anomalies.
[0078] In some embodiments, the horizontal derivative is calculated in the frequency domain. The spatial domain gravity data is converted into a spectrum by a fast Fourier transform, processed by frequency response, and then converted back to the spatial domain by an inverse fast Fourier transform to improve computational efficiency.
[0079] In some embodiments, the series values of the angles in each direction are continuous angles starting from 0 and increasing sequentially according to the amount of angle change, with the maximum angle not exceeding 180 degrees, to ensure coverage of all possible fracture directions.
[0080] In some embodiments, after the upward gravity anomaly is stored on the disk, it can be directly used to plot the planar distribution of gravity anomalies at different depths, providing basic data for analyzing deep geological structures.
[0081] In some embodiments, the results of the first grid point can be combined with geological data for the preliminary division of fault structures and the delineation of geological body boundaries, and its prominent apparent anomalies provide intuitive evidence for geological interpretation.
[0082] In some embodiments, the calculation of the horizontal derivative can capture gradient changes in areas of abrupt changes in underground density, thereby enhancing linear anomaly signals and highlighting gravity field variations associated with faults.
[0083] In some embodiments, the mapping software used is SURFER planar mapping software, which can visualize the processing results in the form of contour lines, clearly showing the trend, extension and gradient changes of anomalies.
[0084] In some embodiments, the hidden linear anomalies highlighted by the second type of grid point results correspond to zones of drastic changes in gravity anomalies, and can be verified by comparison with known fractures, thereby improving the reliability of hidden fracture identification.
[0085] In some embodiments, the experimental selection process for the change in angle requires comparing the overall results at different angles to determine the fineness of the fracture characterization, and prioritizing the selection of angle values that can clearly distinguish fractures in different directions.
[0086] The technical concept of this invention is as follows:
[0087] In order to better obtain deep gravity anomalies, it is often necessary to perform upward extension processing. By selecting an appropriate upward extension height, the gravity field at different depths can be reflected. This gravity field at different depths, combined with the horizontal derivative, is used to delineate fractures and determine their locations.
[0088] In gravity and magnetic exploration theory, the directional derivative of the potential field (gravity and magnetic field) is defined as follows:
[0089]
[0090] Where: V represents the potential field (gravity, magnetism) anomaly; ∠V is the directional derivative of the potential field (gravity, magnetic) anomaly; ∠V is the gradient of the potential field (gravity, magnetic) anomaly; Let be the unit vector in the direction t; Gradient and unit vector of potential field (gravity, magnetic) anomalies The inner product of.
[0091] because:
[0092]
[0093] so:
[0094]
[0095] Equation (2) is the formula for calculating the directional derivative using the first-order horizontal derivatives in the x and y directions.
[0096] θ is the azimuth angle for the direction of differentiation. It is 0° for due north, 90° for due east, 180° for due south, 270° for due west, and 0°≤θ≤360° for other directions.
[0097] This invention employs a fast and efficient algorithm in the frequency domain for calculating the horizontal derivative of gravity. The specific implementation process of this fast and efficient method will be described in detail in the technical solution.
[0098] As discussed earlier, the directional derivative can only highlight linear gravity anomalies caused by fracture structures perpendicular to its derivative direction, and it also produces some false distorted linear anomalies. To overcome this drawback, this invention uses a 360-degree omnidirectional coverage method to obtain the horizontal derivative, which can provide enhanced features of fracture gravity linear anomalies in any direction. The fineness of the direction characterization depends on the size (Δθ) of the given segmentation azimuth angle.
[0099] In fact, for a 360-degree omnidirectional scan, the horizontal derivative only needs to be calculated for the 0-180 degree range, and the horizontal derivative for the 180-360 degree range is simply the inverse of the horizontal derivative for the 0-180 degree range. For linear anomaly zones with prominent gravity anomalies, full coverage is already achieved from 0-180 degrees, so there is no need to calculate the horizontal derivative for the 180-360 degree range.
[0100] For a given azimuth angle (Δθ), the number of directions for calculating the directional derivative within the range of 0-180 degrees is: n θ =(180 / Δθ) 取整 +1 The angle for calculating the horizontal derivative of the scan is θ=[0,Δθ,2Δθ,------,(n θ -1)Δθ]. The level of detail in describing fractures in each direction depends on the magnitude of the given Δθ. The smaller Δθ is, the more directional derivatives can be calculated, resulting in higher directional resolution and more detailed fracture characterization.
[0101] To more comprehensively demonstrate the characteristics of the overall linear gravity anomaly zone based on omnidirectional directional derivatives, this invention employs a method of determining which direction's directional derivative and how it is represented as the grid node result after maximizing the absolute value of the horizontal derivative at each grid point (is it the absolute value of the directional derivative or the original directional derivative value?). These two different results reflect different types of linear gravity anomalies: one reflects the apparent linear gravity anomaly and the boundary anomaly of the geological body; the other reflects a quasi-linear gravity anomaly hidden within the gravity anomaly itself. The method for characterizing the second type of hidden quasi-linear gravity anomaly is unprecedented in previous literature and plays a crucial role in discovering concealed faults.
[0102] Furthermore, the method of this invention also introduces frequency domain gravity upward continuation technology, and further combines the method of determining the full coverage directional derivative and grid node results, which can use gravity anomalies at different depths to characterize fracture structures at different depths.
[0103] The method of this invention effectively utilizes multi-directional horizontal derivatives to characterize the linear anomaly features of gravity in various directions, providing an effective means for overall fracture classification and analysis of fracture structural features and the relationship between fractures in various directions.
[0104] The purpose of this invention is to comprehensively characterize the fracture-related gravity anomalies contained within them by extracting multi-directional linear gravity anomalies from multi-directional gravity directional derivatives. It cleverly utilizes two mapping techniques—one using the maximum or the other using the absolute maximum value of the horizontal derivative in each direction—to integrate all linear gravity anomalies into a single map. This provides an important methodological means for comprehensively classifying fractures and studying their structural characteristics. The processing method employed has unique and innovative ideas and has achieved highly effective application results, showing broad application prospects. This invention plays a crucial role in deeply exploring the linear anomaly information of potential gravity fractures.
[0105] Before describing the technical solution, to facilitate the writing of formulas and the description of the method of this invention, several necessary symbols are defined and the meanings represented by these symbols are briefly described.
[0106] i and j: represent the i-th row and j-th column of the gravity grid point, respectively.
[0107] M and N represent the total number of rows in the gravity grid data and the total number of points in each row, respectively.
[0108] θ: represents the angle used to calculate the horizontal derivative, which is the angle between the horizontal derivative and the coordinate axis.
[0109] Δθ: The angular change between a series of horizontal derivatives.
[0110] FFT: Represents Fast Fourier Transform; FFT -1 This represents the inverse of the Fast Fourier Transform.
[0111] G(i,j): represents gravity grid data; G(u,v) represents the spectral function of the two-dimensional Fourier transform of G(i,j), i.e., G(u,v) = FFT(G). Here, u is the frequency in the direction of the measurement point, and v is the frequency along the direction of the measurement line.
[0112] Gravity anomaly after extending gravity G(i,j) up to a height h.
[0113] Upward Gravity The spectral function of the two-dimensional Fourier transform.
[0114] H(u, v): represents the filter response in the frequency domain.
[0115] This represents the horizontal derivative of θ in the direction of kΔθ.
[0116] This represents the spectrum of the horizontal derivative of θ in the direction of kΔθ.
[0117] A single horizontal derivative can highlight and enhance linear gravity anomalies associated with faults or geological body boundaries that are perpendicular to the derivative direction. Fault-related linear gravity anomalies are characterized by significant elongation, while quasi-linear anomalies associated with geological body boundaries often exhibit a closed loop shape. However, linear gravity anomalies associated with geological body boundaries are distorted in single-directional horizontal derivatives and cannot accurately represent the boundaries of geological bodies. The multi-directional horizontal derivative and its integrated display method of this invention overcome these shortcomings. It can not only highlight and enhance all linear gravity anomalies associated with faults in all directions, but also accurately obtain pseudo-linear gravity anomalies associated with geological body boundaries. These pseudo-linear gravity anomalies are actually the boundaries of geological bodies. Therefore, the multi-directional horizontal derivative of this invention can identify both fault-related linear gravity anomalies and pseudo-linear gravity anomalies associated with geological body boundaries.
[0118] Calculating the directional derivative is not the subject of this invention, but it is a method required for extracting fracture-line gravity anomalies or determining pseudo-linear gravity anomalies at the boundaries of geological bodies.
[0119] The present invention proposes a method for comprehensively discovering and extracting linear gravity anomalies and pseudo-linear gravity anomalies reflecting the boundaries of geological bodies, and how to comprehensively, accurately, and reliably display these linear gravity anomalies.
[0120] The most fundamental aspect of this invention is how to calculate the directional derivative in any direction. Based on equation (2) and the filtering characteristics in the frequency domain, it is not difficult to derive the following process for obtaining the horizontal first derivative in the frequency domain. Since the object of the process is gravity (G), v in equation (2) is g. Therefore, the directional derivative of gravity G is written in the following form:
[0121]
[0122] Applying a two-dimensional FFT to both sides of equation (3) yields:
[0123]
[0124] For gravity extended upwards by h, G in the formula becomes... As can be seen from equation (4), the directional derivative of gravity GDH θ The spectrum GDH of (i, j) θ(u, v) is equal to the product of the frequency response H(u, v) = [iusin(θ) + ivcos(θ)] and the spectrum G(u, v) of the gravity anomaly [G(i, j)]. Given the angle (θ) between the differentiation direction and the measurement direction, the spectrum of the differentiation direction can be obtained using equation (4). Based on the properties of FFT, [GDH] θ Perform inverse Fourier transform (FFT) on [u, v] -1 GDH θ [u, v] gives the spatial domain directional derivative [GDH] in the θ direction. θ (i, j)].
[0125] Given an interval Δθ for calculating the directional derivative of gravity within the range of 0-180 degrees, calculate n θ =(180 / Δθ) 取整 +1, then the angle for calculating the horizontal derivative is θ=[0,Δθ,2Δθ,---,kΔθ,---,(n θ -1)Δθ], the directional derivatives at each angle are obtained by calculating the directional derivative in the frequency domain. The magnitude of Δθ depends on the precision with which the direction of the linear gravity anomaly is resolved, and is also related to the scale and upper extension height of the gravity anomaly. It can be determined through experimental calculations.
[0126] This invention uses two methods to synthesize and map the horizontal derivatives that provide comprehensive coverage, and the results of these two mapping methods have special geological significance.
[0127] 1. Use the maximum absolute value of the horizontal derivative in each direction of the grid point as the calculation result for the grid point.
[0128] 2. Use the directional derivative corresponding to the maximum absolute value of the horizontal derivative in each direction at the grid point as the calculation result for the grid point.
[0129] The processing results can effectively highlight the fracture-related linear anomalies hidden in gravity anomalies (actually, the zones of drastic changes in gravity anomalies), which is quite useful for discovering hidden fractures using gravity.
[0130] Furthermore, this invention introduces a gravity frequency domain extension method, the purpose of which is to obtain gravity anomalies that reflect geological structural characteristics, and then use the core content of this invention to conduct research on the identification of fault structures and geological body boundaries at different depths.
[0131] Specific steps of the present invention:
[0132] Step 1: Based on the required level of detail in the analysis of linear gravity anomalies, determine the angular variation Δθ between a series of horizontal derivatives and the given height h of the gravity extension. When h = 0, the extension operation has no effect, and the result remains the surface gravity anomaly G(i, j).
[0133] Step 2: Perform a two-dimensional forward FFT on G(i,j) to obtain the two-dimensional spectrum G(u,v) of G(i,j), i.e., G(u,v) = FFT[G(i,j)].
[0134] Step 3: Utilize the extended frequency response The spectrum obtained by applying gravity to G(u, v) and extending it upwards by a height h is as follows: right Perform an inverse Fourier transform to obtain the gravity anomaly after the upper extension height h. Right now The gravity anomaly after the extension is stored on the disk for use in drawing.
[0135] Step 4: Based on the given Δθ, use n θ = (180 / Δθ) rounded up + 1 to find n θ The angle θ between each direction and the direction of the measuring point is θ = [0, Δθ, 2Δθ, ..., kΔθ, ..., (n θ -1)Δθ], and calculate n respectively θ The frequency response of the directional derivative in each direction is H(u, v) = [iusin(θ) + ivcos(θ)], and according to H(u, v) (u, v) obtain n θ The directional derivative spectrum of gravity in each direction extended upwards to a height h right Perform inverse Fourier transforms to obtain n. θ Directional derivatives in each direction Right now n θ The directional derivatives in each direction are stored on disk for use in fracture structure analysis by plotting the data.
[0136] Step 5: Use the maximum absolute value of the horizontal derivative in each direction of the grid point as the calculation result for the first grid point;
[0137] After calculation Then, the maximum absolute value of the horizontal derivative in each direction at the grid point is taken as the result at grid point (i, j), i.e.
[0138] The processing results effectively highlight linear gravity anomalies associated with faults and their geological body boundaries. These results can be used to further delineate faults and delineate geological bodies in conjunction with other geological data. The results should be saved as a SURFER format file for drawing or analysis and research on other geological issues.
[0139] Step 6: Use the directional derivative corresponding to the largest absolute value of the horizontal derivative in each direction at the grid point as the calculation result for the second grid point;
[0140] After calculation Then, the maximum absolute value of the horizontal derivative in each direction at each grid point corresponds to... (i, j) represents the result at grid point (i, j), and L is the index of the horizontal derivative corresponding to the maximum absolute value of the directional derivative. That is, in calculating... During the process, mark the k that achieves the maximum value and use this k value as L.
[0141] The processing results effectively highlight fault-related linear anomalies (actually zones of dramatic gravity anomaly changes) hidden within gravity anomalies, proving highly useful for discovering concealed faults using gravity. The results should be saved in SURFER format for plotting or analysis and research on other geological issues.
[0142] Step 7: Using SURFER plane mapping software, draw the gravity anomaly plane diagram extending h upwards, the horizontal derivative plane diagrams in each direction after extending h upwards, the plane diagram of the maximum absolute value of the horizontal derivative, and the plane diagram of the horizontal derivative calculated at the maximum absolute value of the horizontal derivative, etc., to complete the processing calculation and result output of the method of the present invention.
[0143] A computational system for omnidirectional coverage detection of linear gravity anomalies, employing the computational method for omnidirectional coverage detection of linear gravity anomalies described in any one of the above-mentioned methods; the computational system for omnidirectional coverage detection of linear gravity anomalies includes:
[0144] The setting module is used to set the angle change and upper extension height used to calculate the horizontal derivative;
[0145] The forward transform module is used to perform a two-dimensional fast Fourier forward transform on the original gravity grid data to obtain the original gravity spectrum;
[0146] The inverse transform module is used to calculate the upper gravity spectrum corresponding to the upper extension height based on the original gravity spectrum and the upper extension frequency response, and then perform a two-dimensional fast Fourier inverse transform on the upper extension gravity spectrum to obtain the upper extension gravity anomaly.
[0147] The horizontal derivative calculation module is used to determine the direction angle for calculating the horizontal derivative within the range of 0 to 180 degrees based on the angle change, and to calculate the horizontal derivative corresponding to each direction angle based on the upward gravity anomaly.
[0148] The first grid point calculation module is used to calculate the maximum value of the absolute value of the horizontal derivative in each direction at the grid point based on the horizontal derivative corresponding to the angle in each direction, and to obtain the calculation result of the first grid point.
[0149] The second grid point calculation module is used to calculate the horizontal derivative value corresponding to the maximum absolute value of the horizontal derivative in each direction at the grid point based on the horizontal derivative corresponding to the angle in each direction, and obtain the calculation result of the second grid point.
[0150] To illustrate the role of the method of the present invention in highlighting the linear anomaly information of fracture gravity, gravity data from the SLPD region were processed using the method described in the present invention, and the processing effect was analyzed.
[0151] The gravity data of the SLPD region were processed using the method of this invention, including ground gravity and gravity data extending 3 km upwards. With Δθ = 30°, simple directional derivative processing was performed in six directions, resulting in two innovative maps from this invention.
[0152] 1. Ground gravity processing results and corresponding descriptions
[0153] right Figure 2 The gravity anomaly shown yielded directional derivatives in six directions. The results for the directional derivatives in four directions (0°, 30°, 120°, and 150°) are shown below.
[0154] Figures 3 to 6 To illustrate the directional derivatives of gravity in the four directions, these four figures are compared with... Figure 2 The comparison reveals that they effectively highlight and characterize the linear gravity anomalies perpendicular to their respective directional derivatives and the boundary anomalies of local gravity.
[0155] The first processing result (calculation result of the first grid point) comprehensively reflects the set of boundaries between linear anomalies reflected by all directions of the guide numbers and the boundaries reflecting local gravity anomalies of the geological body, such as... Figure 7 The integrated anomaly information is complete, detailed, and without redundancy, making it an important map for comprehensively interpreting faults and delineating the boundaries of geological bodies. It also fully demonstrates the effectiveness of the method of this invention, which has a strong ability to characterize gravity linear anomalies and quasi-linear anomalies. However, it has not been able to effectively characterize incomplete linear anomalies that reflect faults hidden in gravity anomalies.
[0156] A comparison of the second processing results (calculation results of the second grid points) with the first processing results and the ground gravity anomaly clearly shows that the implicit fracture information not displayed in the first processing method is well highlighted. Figure 8 The clearly reflected near-east-west, northeast, and northwest-trending banded anomalies are all hidden fault gravity anomalies related to faults. These anomalies often reflect abrupt change zones or zonal distributions in surface gravity anomalies. This method is highly effective in identifying hidden faults reflected by such gravity anomalies. This processing method has not been found or reported in any literature. It adds an innovative processing method for dividing and extracting fault information using gravity anomalies, making more effective use of gravity for geological interpretation.
[0157] 2. Results of ground gravity treatment extending 2km upwards
[0158] Figures 9 to 15 This also demonstrates that the two methods of this invention have the same effect on studying deep fractures using extended gravity. The corresponding description is the same. Figures 2 to 8 The description.
[0159] The embodiments described above are for illustrative purposes only and are not intended to limit the invention. Therefore, any changes in numerical values or substitutions of equivalent elements should still fall within the scope of this invention.
Claims
1. A calculation method for omnidirectional detection of linear gravity anomalies, characterized in that, include: Step 1: Based on the required level of detail in the analysis of linear gravity anomalies, determine the angular variation and upward extension height between a series of horizontal derivatives; Step 2: Perform a two-dimensional fast Fourier transform on the original gravity grid data to obtain the original gravity spectrum; Step 3: Based on the original gravity spectrum and the upper extension frequency response, calculate the upper extension gravity spectrum corresponding to the upper extension height, and then perform a two-dimensional fast Fourier inverse transform on the upper extension gravity spectrum to obtain the upper extension gravity anomaly. The upper-extended frequency response is an attenuation function generated based on the upper-extended height and frequency parameters; Step 4: Based on the angle change, determine the direction angles within the range of 0–180 degrees from which the horizontal derivative needs to be calculated, and calculate the horizontal derivatives corresponding to each direction angle based on the upward gravity anomaly. The calculation method for the horizontal derivatives corresponding to each direction angle is as follows: first, generate the corresponding horizontal derivative frequency response based on the direction angle; then, multiply the horizontal derivative frequency response with the upward gravity spectrum to obtain the horizontal derivative spectrum for each direction; finally, perform a two-dimensional inverse fast Fourier transform on the horizontal derivative spectrum for each direction to obtain the horizontal derivatives corresponding to each direction angle. The calculation of the horizontal derivatives is completed in the frequency domain. The spatial domain gravity data is converted into a spectrum through fast Fourier transform, and after frequency response processing, it is converted back to the spatial domain through inverse fast Fourier transform. Step 5: Based on the horizontal derivatives corresponding to the angles in each direction, calculate the maximum value among the absolute values of the horizontal derivatives in each direction at the grid point to obtain the calculation result for the first grid point; Step 6: Based on the horizontal derivatives corresponding to the angles in each direction, calculate the horizontal derivative value corresponding to the maximum absolute value of the horizontal derivatives in each direction at the grid point, and obtain the calculation result for the second grid point; Step 7: Based on the upper gravity anomaly, the horizontal derivatives corresponding to the angles in each direction, the calculation results of the first grid point and the calculation results of the second grid point, draw and output the corresponding planar diagram.
2. The calculation method for omnidirectional coverage detection of linear gravity anomalies according to claim 1, characterized in that, In step 1, when the set upper extension height is 0, the upper extension gravity anomaly is the original gravity grid data.
3. The calculation method for omnidirectional coverage detection of linear gravity anomalies according to claim 1, characterized in that, In step 3, the extended gravity spectrum is the product of the original gravity spectrum and the extended frequency response.
4. The calculation method for omnidirectional coverage detection of linear gravity anomalies according to claim 1, characterized in that, In step 4, the horizontal derivative in the range of 0 to 180 degrees covers the detection requirements of gravity linear anomalies in all directions, and the horizontal derivative in the range of 180 to 360 degrees is the opposite value of the horizontal derivative corresponding to 0 to 180 degrees.
5. The calculation method for omnidirectional coverage detection of linear gravity anomalies according to claim 1, characterized in that, The calculation results of the first grid point are used to reflect the apparent linear gravity anomaly and the boundary anomaly of the geological body; the calculation results of the second grid point are used to reflect the hidden linear-like gravity anomaly.
6. The calculation method for omnidirectional coverage detection of linear gravity anomalies according to claim 1, characterized in that, The magnitude of the angle change is determined through experiments; the smaller the angle change, the more directional angles whose horizontal derivatives need to be calculated within the range of 0 to 180 degrees, and the higher the directional resolution accuracy.
7. The calculation method for omnidirectional coverage detection of linear gravity anomalies according to claim 1, characterized in that, The maximum angle in each direction is 180 degrees and ( -1) The difference does not exceed Angle value, For the number of directions and angles, This represents the change in angle.
8. A computational system for omnidirectional detection of linear gravity anomalies, characterized in that, The calculation method for detecting linear gravity anomalies with full coverage according to any one of claims 1 to 7; The computational system for detecting linear gravity anomalies with full coverage includes: The setting module is used to determine the angular variation and upper extension height between a series of horizontal derivatives based on the level of detail in the analysis of linear gravity anomalies. The forward transform module is used to perform a two-dimensional fast Fourier forward transform on the original gravity grid data to obtain the original gravity spectrum; The inverse transform module is used to calculate the upper gravity spectrum corresponding to the upper extension height based on the original gravity spectrum and the upper extension frequency response, and then perform a two-dimensional fast Fourier inverse transform on the upper extension gravity spectrum to obtain the upper extension gravity anomaly; the upper extension frequency response is an attenuation function generated based on the upper extension height and frequency parameters; The horizontal derivative calculation module is used to determine the direction angles from 0 to 180 degrees where the horizontal derivative needs to be calculated based on the angle change, and to calculate the horizontal derivatives corresponding to each direction angle based on the upward gravity anomaly. The calculation method for the horizontal derivatives corresponding to each direction angle is as follows: first, the corresponding horizontal derivative frequency response is generated based on the direction angle; then, the horizontal derivative frequency response is multiplied by the upward gravity spectrum to obtain the horizontal derivative spectrum for each direction; finally, a two-dimensional inverse fast Fourier transform is performed on the horizontal derivative spectrum for each direction to obtain the horizontal derivatives corresponding to each direction angle. The calculation of the horizontal derivatives is completed in the frequency domain. The spatial domain gravity data is converted into a spectrum through fast Fourier transform, and after frequency response processing, it is converted back to the spatial domain through inverse fast Fourier transform. The first grid point calculation module is used to calculate the maximum value of the absolute value of the horizontal derivative in each direction at the grid point based on the horizontal derivative corresponding to the angle in each direction, and to obtain the calculation result of the first grid point. The second grid point calculation module is used to calculate the horizontal derivative value corresponding to the maximum absolute value of the horizontal derivative at each grid point based on the horizontal derivative corresponding to the angle in each direction, and obtain the calculation result of the second grid point; based on the upward gravity anomaly, the horizontal derivative corresponding to the angle in each direction, the calculation result of the first grid point and the calculation result of the second grid point, the corresponding planar map is drawn and output.
Citation Information
Patent Citations
Physical geography exploration gravity and magnetic data processing method
CN101285896A
Methods and systems for data collection, learning, and streaming of machine signals for part identification and operating characteristics determination using the industrial internet of things
US20200133254A1