Gravity data downward continuation method, equipment, medium and product of discretization model based on inverse band-limited Poisson integral

Abstract

The invention discloses a gravity data downward continuation method, device, medium and product of a discretization model based on inverted band limit Poisson integral, and relates to the technical field of gravity data processing, gridding is carried out on discrete gravity point data, and grid limit frequency, the truncation order of grid gravity abnormal data and the type of the discretization model of the inverted band limit Poisson integral are determined; when a direct discretization model in a recovery mode is removed, a coefficient matrix, a parameter vector and a constant vector are determined based on grid gravity abnormal data, continuation height and an earth gravity field model; determining a target calculation formula based on the coefficient matrix, the parameter vector and the constant vector; when the model is an indirect discretization model in a near-region integration mode, determining a target calculation formula of the discretization model based on the grid gravity anomaly data and the continuation height; and obtaining a downward continuation result based on a target calculation formula of the discretization model. According to the invention, stable and high-precision downward continuation can be carried out on the gravity data.

Description

Technical Field

[0001] The present application relates to the field of gravity data processing, and particularly to a method, device, medium and product for downward continuation of gravity data based on a discretization model of the rewind-limited Poisson integral. Background Art

[0002] Gravity data is an important potential field data in the field of physical geodesy, and its acquisition methods include absolute gravity measurement, relative gravity measurement, satellite altimetry, shipborne gravity measurement, airborne gravity measurement, etc. When the measurement height surface of gravity data is inconsistent with the actual application height surface, it is necessary to transform the gravity data to the required height surface through an extension algorithm.

[0003] The Poisson integral is an important theoretical tool for the extension of potential field data and is widely used in practice. The Poisson integral is a continuous integral, but since the integrated gravity data is grid data, it is necessary to discretize the Poisson integral when performing downward continuation. The rewind-limited Poisson integral has a filtering property and stable numerical operations, and has received increasing attention in the field of physical geodesy. It has been found that the discretization method has an important influence on the downward continuation accuracy of the rewind-limited Poisson integral. In order to obtain a stable and high-precision downward continuation result, there is an urgent need for a discretization method of the rewind-limited Poisson integral that can perform efficient continuation according to the characteristics of the spectral information contained in the gravity data. Summary of the Invention

[0004] The purpose of the present application is to provide a method, device, medium and product for downward continuation of gravity data based on a discretization model of the rewind-limited Poisson integral, which can perform stable and high-precision downward continuation of gravity data.

[0005] To achieve the above purpose, the present application provides the following solutions:

[0006] In a first aspect, the present application provides a method for downward continuation of gravity data based on a discretization model of the rewind-limited Poisson integral, where the method for downward continuation of gravity data based on the discretization model of the rewind-limited Poisson integral includes:

[0007] Obtain discrete gravity point data to be continued, grid resolution, and continuation height;

[0008] Perform grid processing on the discrete gravity point data according to the grid resolution to obtain grid gravity anomaly data;

[0009] Based on the grid resolution, determine the grid frequency limit, and based on the grid frequency limit and the discrete gravity point data, determine the truncation order of the grid gravity anomaly data;

[0010] Determine the type of the discretization model of the rewind-limited Poisson integral based on the grid-limited frequency and the truncation order of the grid gravity anomaly data; the type of the discretization model is a direct discretization model in the remove-restore mode or an indirect discretization model in the near-zone integral mode;

[0011] When the type of the discretization model is a direct discretization model in the remove-restore mode:

[0012] Obtain the Earth gravity field model for reference in the remove-restore mode;

[0013] Based on the grid gravity anomaly data and the continuation height, determine the coefficient matrix of the discretization model; the coefficient matrix includes diagonal elements and non-diagonal elements;

[0014] Based on the grid gravity anomaly data, the continuation height, and the Earth gravity field model, determine the parameter vector and the constant vector of the discretization model;

[0015] Based on the coefficient matrix, the parameter vector, and the constant vector of the discretization model, determine the target calculation formula of the discretization model;

[0016] When the type of the discretization model is an indirect discretization model in the near-zone integral mode:

[0017] Based on the grid gravity anomaly data and the continuation height, determine the coefficient matrix and the parameter vector of the discretization model;

[0018] Based on the coefficient matrix and the parameter vector of the discretization model, determine the target calculation formula of the discretization model;

[0019] Based on the target calculation formula of the corresponding discretization model, obtain the downward continuation result.

[0020] In a second aspect, the present application provides a computer device, including: a memory, a processor, and a computer program stored on the memory and executable on the processor, where the processor executes the computer program to implement the gravity data downward continuation method based on the discretization model of the rewind-limited Poisson integral as described in any one of the above.

[0021] In a third aspect, the present application provides a computer-readable storage medium, on which a computer program is stored, and when the computer program is executed by a processor, it implements the gravity data downward continuation method based on the discretization model of the rewind-limited Poisson integral as described in any one of the above.

[0022] In a fourth aspect, the present application provides a computer program product, including a computer program, and when the computer program is executed by a processor, it implements the gravity data downward continuation method based on the discretization model of the rewind-limited Poisson integral as described in any one of the above.

[0023] According to the specific embodiments provided in this application, this application has the following technical effects:

[0024] This application discloses a gravity data downward continuation method, device, medium, and product based on a discretized model of the rewind-limited Poisson integral. First, grid the discrete gravity point data according to the grid resolution; determine the grid limited frequency based on the grid resolution, and determine the truncation order of the grid gravity anomaly data based on the grid limited frequency and the discrete gravity point data; determine the type of the discretized model of the rewind-limited Poisson integral based on the grid limited frequency and the truncation order; the type of the discretized model is a direct discretized model in the remove-and-recover mode or an indirect discretized model in the near-zone integral mode; secondly, for different discretized models, determine the target calculation formula of the corresponding discretized model; finally, obtain the downward continuation result based on the target calculation formula of the discretized model. This application uses different discretized models to perform rewind-limited Poisson integral processing on the gravity data to be continued according to different situations, significantly improving the accuracy of the obtained downward continued gravity data. Description of the Drawings

[0025] In order to more clearly illustrate the technical solutions in the embodiments of this application or the prior art, the following will briefly introduce the drawings required for use in the embodiments. Obviously, the drawings in the following description are only some embodiments of this application. For those of ordinary skill in the art, other drawings can be obtained based on these drawings without creative efforts.

[0026] Figure 1 It is an application environment diagram of the gravity data downward continuation method based on the discretized model of the rewind-limited Poisson integral in an embodiment of this application;

[0027] Figure 2 It is a schematic flowchart of the gravity data downward continuation method based on the discretized model of the rewind-limited Poisson integral provided in an embodiment of this application;

[0028] Figure 3 It is a schematic structural diagram of a computer device provided in an embodiment of this application.

[0029] Reference Signs:

[0030] Terminal - 102, Server - 104. Detailed Embodiments

[0031] The following will clearly and completely describe the technical solutions in the embodiments of the present application with reference to the accompanying drawings in the embodiments of the present application. Obviously, the described embodiments are only a part of the embodiments of the present application, rather than all the embodiments. All other embodiments obtained by those of ordinary skill in the art based on the embodiments in the present application without creative efforts belong to the scope of protection of the present application.

[0032] To make the above objects, features, and advantages of the present application more obvious and understandable, the present application will be further described in detail below with reference to the accompanying drawings and specific embodiments.

[0033] The downward continuation method of gravity data based on the discretized model of the rewind-limited Poisson integral provided by the embodiments of the present application can be applied to an application environment as Figure 1 shown. Among them, the terminal 102 communicates with the server 104 through a network. The data storage system can store the data that the server 104 needs to process. The data storage system can be set separately, integrated on the server 104, or placed on the cloud or other servers. The terminal 102 can send the discrete gravity point data to be continued, the grid resolution, and the continuation height to the server 104. After receiving the discrete gravity point data to be continued, the grid resolution, and the continuation height, for the discrete gravity point data to be continued, the grid resolution, and the continuation height, the server 104 performs the downward continuation of gravity data based on the discretized model of the rewind-limited Poisson integral. The server 104 can feedback the obtained downward continuation result to the terminal 102. In addition, in some embodiments, the downward continuation method of gravity data based on the discretized model of the rewind-limited Poisson integral can also be implemented by the server 104 or the terminal 102 alone. For example, the terminal 102 can directly perform the downward continuation of gravity data based on the discretized model of the rewind-limited Poisson integral on the discrete gravity point data to be continued, the grid resolution, and the continuation height. The server 104 can also obtain the discrete gravity point data to be continued, the grid resolution, and the continuation height from the data storage system and perform the downward continuation of gravity data based on the discretized model of the rewind-limited Poisson integral for the discrete gravity point data to be continued, the grid resolution, and the continuation height.

[0034] Among them, the terminal 102 can be, but is not limited to, various desktop computers, laptop computers, smartphones, and tablet computers. The server 104 can be implemented by an independent server or a server cluster composed of multiple servers, and can also be a cloud server.

[0035] In an exemplary embodiment, as Figure 2As shown, a downward continuation method of gravity data based on a discretization model of the rewind-limited Poisson integral is provided. This method is executed by a computer device, which can be specifically executed by a computer device such as a terminal or a server alone, or jointly executed by a terminal and a server. In the embodiments of the present application, taking this method applied to Figure 1 server 104 in

[0036] Step S1, obtain the discrete gravity point data to be continued, the grid resolution, and the continuation height.

[0037] Step S2, perform grid processing on the discrete gravity point data according to the grid resolution to obtain grid gravity anomaly data.

[0038] Step S3, determine the grid limited frequency based on the grid resolution, and determine the truncation order of the grid gravity anomaly data based on the grid limited frequency and the discrete gravity point data.

[0039] Step S4, determine the type of the discretization model of the rewind-limited Poisson integral based on the grid limited frequency and the truncation order of the grid gravity anomaly data; the type of the discretization model is a direct discretization model in the remove-and-restore mode or an indirect discretization model in the near-zone integral mode.

[0040] When the type of the discretization model is a direct discretization model in the remove-and-restore mode, execute Steps S5 to S8:

[0041] Step S5, obtain the Earth's gravity field model referred to in the remove-and-restore mode.

[0042] Step S6, determine the coefficient matrix of the discretization model based on the grid gravity anomaly data and the continuation height; the coefficient matrix includes diagonal elements and non-diagonal elements.

[0043] Step S7, determine the parameter vector and the constant vector of the discretization model based on the grid gravity anomaly data, the continuation height, and the Earth's gravity field model.

[0044] Step S8, determine the target calculation formula of the discretization model based on the coefficient matrix, the parameter vector, and the constant vector of the discretization model.

[0045] When the type of the discretization model is an indirect discretization model in the near-zone integral mode, execute Steps S9 to S10:

[0046] Step S9, determine the coefficient matrix and the parameter vector of the discretization model based on the grid gravity anomaly data and the continuation height.

[0047] Step S10, determine the target calculation formula of the discretization model based on the coefficient matrix and the parameter vector of the discretization model.

[0048] Step S11: Obtain the downward continuation result based on the target calculation formula of the corresponding discretization model.

[0049] As an optional implementation manner, step S3 specifically includes:

[0050] Step S31: Calculate the grid aliasing limit based on the grid resolution using the formula M' = π / Δ; where M' represents the grid aliasing limit; Δ represents the grid resolution in radians.

[0051] Step S32: Calculate the truncation order of the grid gravity data based on the grid aliasing limit and the discrete gravity point data using the formula M = M'·min{n 1 / n 2 , 1}; where M represents the truncation order of the grid gravity data; n 1 represents the number of discrete gravity points; n 2 represents the number of grid gravity data.

[0052] As an optional implementation manner, step S4 specifically includes:

[0053] Step S41: Determine whether the difference between the grid aliasing limit and the truncation order of the grid gravity data is greater than or equal to a preset threshold to obtain a judgment result. Specifically, if the grid resolution is 3′, the difference threshold is 900 orders; if the grid resolution is 2′, the difference threshold is 2200 orders, and the higher the resolution, the larger the difference threshold.

[0054] Step S42: If the judgment result is yes, determine that the type of the discretization model of the rewind-limited Poisson integral is the direct discretization model in the remove-restore mode.

[0055] Step S43: If the judgment result is no, determine that the type of the discretization model of the rewind-limited Poisson integral is the indirect discretization model in the near-zone integral mode.

[0056] As an optional implementation manner, when the discretization model is the direct discretization model in the remove-restore mode, in step S6, the calculation formula for the non-diagonal elements in the coefficient matrix of the discretization model is

[0057]

[0058] where represents the element in the i-th row and j-th column of the coefficient matrix of the direct discretization model in the remove-restore mode; s j represents the area of the j-th integration grid, s j =Δ 2 sinθ j , θ jis the co-latitude of the center point of the j-th integration grid, and Δ represents the grid resolution in radians; N 1，1 represents the low-order truncation order of the rewind-limited Poisson kernel function when the discretization model is the direct discretization model in the remove-restore mode, N 1，1 = L, where L represents the truncation order of the reference Earth gravity field model; N 1,2 represents the high-order truncation order of the rewind-limited Poisson kernel function when the discretization model is the direct discretization model in the remove-restore mode, N 2,1 = M; R represents the average radius of the Earth; r represents the geocentric radius vector of the calculation point; P n represents the Legendre polynomial of degree n; ψ ij represents the angular distance between the center point of the i-th calculation grid and the center point of the j-th integration grid; ψ 1 represents the radius of the near-region integration domain.

[0059] In step S6, the calculation formula for the diagonal element in the coefficient matrix of the discretization model is

[0060]

[0061] where represents the diagonal element of the i-th row in the coefficient matrix of the direct discretization model in the remove-restore mode; I n (ψ 0 ) represents the function value of degree n with parameter ψ 0 ; ψ 0 represents the radius of the circular approximation region of the calculation grid. Among them, ψ 0 can be determined according to the principle that the area of the calculation grid is equal to that of its circular approximation region, and the specific calculation formula is

[0062]

[0063] where θ represents the co-latitude of the calculation point.

[0064] I n (ψ 0 ) has the calculation formula of

[0065]

[0066] In step S7, the calculation formula for the parameter vector of the discretization model is

[0067]

[0068] where represents the i-th element in the parameter vector; Δg i represents the gravity anomaly data of the i-th grid; Denote the reference gravity anomaly of the i-th grid calculated using the Earth's gravity field model.

[0069] In step S7, the calculation formula for the constant vector of the discretized model is

[0070]

[0071] where, b i denotes the i-th element in the constant vector; N denotes the highest order affected by the high-frequency far zone; Ω i denotes the spherical coordinates of the i-th calculation point; denotes the n-th order gravity anomaly calculated using the Earth's gravity field model; denotes the n-th order far zone truncation coefficient of the rewind-limited Poisson kernel function.

[0072] where, The calculation formula for

[0073]

[0074] where, K IBL denotes the rewind-limited Poisson kernel function. When the discretized model is the direct discretized model in the remove-and-restore mode, its calculation formula is

[0075]

[0076] In step S8, the target calculation formula for the discretized model is

[0077] Δg High (R) = A IBL1 Δg High (r) + b(9)

[0078] where, Δg H igh(R) denotes the output vector of the direct discretized model in the remove-and-restore mode; A IBL1 denotes the coefficient matrix of the direct discretized model in the remove-and-restore mode; Δg High (r) denotes the parameter vector; b denotes the constant vector.

[0079] As an optional implementation, when the type of the discretized model is the indirect discretized model in the near zone integration mode:

[0080] In step S9, the calculation formula for the non-diagonal elements of the coefficient matrix of the discretized model is

[0081]

[0082] where, Denote the element in the \(i\)-th row and \(j\)-th column of the coefficient matrix of the indirect discretization model in the near-zone integration mode; \(N\) 2，1 Denote the low-order truncation order of the rewind-limited Poisson kernel function when the type of the discretization model is the indirect discretization model in the near-zone integration mode; \(N\) 2，1 = -1; \(N\) 2,2 Denote the high-order truncation order of the rewind-limited Poisson kernel function when the type of the discretization model is the indirect discretization model in the near-zone integration mode; \(N\) 2,2 = \(M\), that is \(N\) 2,2 The value of is equal to the truncation order of the grid gravity anomaly data.

[0083] In step S9, the calculation formula for the diagonal element in the coefficient matrix of the discretization model is

[0084]

[0085] where Denote the diagonal element in the \(i\)-th row of the coefficient matrix of the indirect discretization model in the near-zone integration mode; \(J\) denotes the number of near-zone grids.

[0086] At this time, use the grid gravity anomaly data to construct the parameter vector \(\Delta g(r)\) of the discretization model.

[0087] In step S10, the objective calculation formula of the discretization model is

[0088] \(\Delta g\) IBL2 (R) = A IBL2 \(\Delta g(r)(12)\)

[0089] where \(\Delta g\) IBL2 (R) represents the output vector of the indirect discretization model in the near-zone integration mode, \(A\) IBL2 represents the coefficient matrix of the indirect discretization model in the near-zone integration mode; \(\Delta g(r)\) represents the parameter vector constructed from the grid gravity anomaly data.

[0090] As an optional implementation manner, when the type of the discretization model is the direct discretization model in the remove-restore mode, step S11 specifically includes:

[0091] Step S1111, perform rewind-limited Poisson integral processing based on the objective calculation formula of the direct discretization model in the remove-restore mode to obtain the output vector of the discretization model.

[0092] Step S1112, based on the output vector of the discretization model and the Earth's gravity field model, use the formula \(\Delta g\) IBL1 (R) = \(\Delta g\) High (R) + \(\Delta g\) G(R), obtaining the output vector after restoring the reference model value; where, Δg IBL1 (R) represents the output vector after restoring the reference model value; Δg G (R) is the reference gravity anomaly vector on the boundary surface with radius R calculated using the Earth's gravity field model.

[0093] Step S1113, based on the output vector after restoring the reference model value, eliminate the data at the edge of the calculation area to obtain the downward continuation result.

[0094] As an optional implementation manner, when the type of the discretized model is an indirect discretized model in the near-zone integration mode, step S11 specifically includes:

[0095] Step S1121, based on the target calculation formula of the direct discretized model in the remove-restore mode, perform the rewind-limited Poisson integral processing to obtain the output vector of the discretized model;

[0096] Step S1122, based on the output vector of the discretized model, eliminate the data at the edge of the calculation area to obtain the downward continuation result.

[0097] Rewind-limited Poisson integral

[0098] The rewind-limited Poisson integral uses the gravity data of the measurement surface at a certain height as the integrand physical quantity, and its expression is

[0099]

[0100] where, Ω 0 represents the global integration domain, dΩ' is the spherical differential unit, Δg(R,Ω) is the gravity anomaly at the calculation point on the boundary surface, Δg(r,Ω') is the gravity anomaly at the flowing integration point on the measurement surface, K IBL represents the rewind-limited Poisson kernel function, and its calculation formula is

[0101]

[0102] where, N 1 is the low-order truncation order, N 2 is the high-order truncation order, ψ is the geocentric angular distance between the calculation point and the integration point, and n is the spherical harmonic order.

[0103] Since the downward continuation using the rewind-limited Poisson integral can obtain the gravity data on the boundary surface without inverse operation, the downward continuation process based on the rewind-limited Poisson integral has numerical stability. When applying the rewind-limited Poisson integral, it needs to be discretized. If the kernel function value at the central grid center point is used as the approximate average value of the kernel function within the grid during discretization, it will cause relatively large discretization errors. To avoid large discretization errors in the central area, direct or indirect discretization schemes can be adopted in practical applications. The direct discretization method is a calculation method that directly calculates the integral value of the central grid using the analytical formula of the continuous integral in the central grid, while the indirect discretization method is a calculation method that uses a modified integral formula that makes the integral value of the central grid zero to avoid directly calculating the integral value in the central area.

[0104] In the past, the rewind-limited Poisson integral was usually applied to the remove-and-replace mode. In fact, the rewind-limited Poisson integral can also be applied to the non-remove-and-replace mode, and the non-remove-and-replace mode can be further divided into the near-zone integral mode and the integral mode considering the influence of the far zone. The discretization method of the rewind-limited Poisson integral has an important impact on the continuation accuracy. The traditional discretization algorithm of the rewind-limited Poisson integral is essentially an indirect discretization algorithm that includes the influence of high-frequency far zones in the remove-and-replace mode. To further improve the downward continuation accuracy of the rewind-limited Poisson integral, this application first derives the direct and indirect discretization formulas of the rewind-limited Poisson integral in the near-zone integral mode, the integral mode considering the influence of the far zone, and the remove-and-replace mode, and then verifies the effectiveness of the discretization algorithms in different calculation modes through numerical experiments.

[0105] Near-zone integral mode

[0106] Direct discretization method:

[0107] Theoretically, the Poisson integral is a global integral. However, due to the large amount of calculation for global integration, the Poisson integral is usually calculated using the near-zone integral in practical applications.

[0108] To derive the analytical formula of the rewind-limited Poisson integral for the central grid, the central grid is approximated as a circular region C with the grid center point as the origin. 0 C 0 The radius ψ 0 of C 0 can be determined according to the principle that the area of C is equal to the area of the central grid (i.e., ). Through derivation, it can be obtained that the diagonal elements of the discretization coefficient matrix of the near-zone rewind-limited Poisson integral when using the direct discretization method are

[0109]

[0110] where Denotes the diagonal element of the \(i\)-th row of the discretization coefficient matrix of the rewind-limited Poisson integral, \(I\) n (\(\psi\) 0 ) The calculation formula is

[0111]

[0112] The off-diagonal elements of the discretization coefficient matrix of the rewind-limited Poisson integral are

[0113]

[0114] Wherein Denotes the element of the \(i\)-th row and \(j\)-th column in the discretization coefficient matrix of the rewind-limited Poisson integral, \(\psi\) 1 Denotes the radius of the near-region integration domain.

[0115] Indirect discretization method:

[0116] The low-order truncation order \(N\) in the rewind-limited Poisson kernel function 1 Is an integer greater than or equal to -1. Using the orthogonality of Legendre polynomials:

[0117]

[0118] The spherical surface integral of the rewind-limited Poisson kernel function can be obtained as

[0119]

[0120] Multiply both sides of Equation (19) by \(r / R\cdot\Delta g(r,\Omega)\), and then subtract it from the rewind-limited Poisson integral formula, the modified form of the rewind-limited Poisson integral can be obtained as follows:

[0121]

[0122] Changing the integration domain in Equation (20) to the near-region integration domain is the near-region rewind-limited Poisson integral in the modified form. According to Equation (20), the diagonal elements of the indirect discretization coefficient matrix of the near-region rewind-limited Poisson integral are

[0123]

[0124] Where \(J\) denotes the number of near-region grids.

[0125] It can be seen by comparing equations (20) and (13) that the calculation formulas for the non - diagonal elements of the indirect discretization coefficient matrix of the near - zone rewind - limited Poisson integral are the same as those of the direct discretization coefficient matrix. Additionally, since the near - zone integration mode ignores the far - zone integration value, neither the direct nor the indirect discretization formulas in the near - zone integration mode contain a constant vector.

[0126] Integration mode considering far - zone influence

[0127] Direct discretization method:

[0128] Since directly performing the far - zone Poisson integral requires far - zone measured gravity data and involves a large amount of calculation, in practical applications, the far - zone Poisson integral value is usually calculated using the Earth's gravity field model. The expression for the rewind - limited Poisson integral considering far - zone influence is

[0129]

[0130] where \(C\) 1 represents the near - zone integration domain, \(L'\) is the cut - off order of the far - zone influence, and \(\Delta g_{n}\) is the \(n\) - th order model gravity anomaly at the calculation point.

[0131] It can be seen from this that the difference between the direct discretization formula of the rewind - limited Poisson integral in the integration mode considering far - zone influence and that in the near - zone integration mode lies only in the addition of a constant term composed of the far - zone influence.

[0132] Indirect discretization method:

[0133] Through derivation, the modified formula for the rewind - limited Poisson integral considering far - zone influence is

[0134]

[0135] It can be known from this that the diagonal elements of the indirect discretization coefficient matrix of the rewind - limited Poisson integral considering far - zone influence are

[0136]

[0137] By comparing equations (22) and (23), it can be seen that the constant vectors of the indirect discretization formula of the Poisson integral considering far - zone influence are the same as those of the direct discretization formula of the Poisson integral considering far - zone influence. In addition, the calculation formulas for the non - diagonal elements of the direct and indirect discretization coefficient matrices of the rewind - limited Poisson integral considering far - zone influence are the same as those of the non - diagonal elements of the near - zone rewind - limited Poisson integral.

[0138] Remove - restore mode

[0139] Direct discretization method:

[0140] When applying the remove-restore technique, it is first necessary to remove the low-frequency information of the gravity data to be extended, so that the gravity data to be extended in the Poisson integral only contains high-frequency information. In this way, the gravity data obtained by the Poisson integral only contains high-frequency information, and restoring its low-frequency information can obtain the final extension result. Among them, the low-frequency information of the gravity data is calculated by the Earth gravity field model. The expression of the rewind-limited Poisson integral in the remove-restore mode is

[0141]

[0142] In the formula, Δg High (R,Ω) is the high-frequency gravity anomaly at the projection position of the calculation point on the boundary sphere, and Δg High (r,Ω') is the high-frequency gravity anomaly at the flowing integration point on the measurement sphere.

[0143] If the low-order truncation order of the high-frequency gravity data in the remove-restore mode is lower than the highest order of the gravity field model, the influence of the high-frequency far zone can be further increased. The expression of the rewind-limited Poisson integral considering the influence of the high-frequency far zone in the remove-restore mode is

[0144]

[0145] In the formula, N represents the highest order of the influence of the high-frequency far zone, and L represents the truncation order of the reference Earth gravity field model in the remove-restore mode. According to formulas (25) and (26), it can be seen that the calculation formulas of the diagonal and non-diagonal elements of the direct discretization coefficient matrix of the rewind-limited Poisson integral in the remove-restore mode are the same as the corresponding formulas of the near-zone rewind-limited Poisson integral.

[0146] Indirect discretization method:

[0147] The modified form of the rewind-limited Poisson integral in the remove-restore mode is

[0148]

[0149] The modified rewind-limited Poisson integral including the influence of the high-frequency far zone in the remove-restore mode is

[0150]

[0151] According to Equations (27) and (28), it can be seen that the calculation formula for the diagonal elements of the indirect discretization coefficient matrix of the rewind-limited Poisson integral without the influence of the high-frequency far zone in the removal and restoration mode is the same as the corresponding formula for the near-zone rewind-limited Poisson integral; while the calculation formula for the diagonal elements of the indirect discretization coefficient matrix of the rewind-limited Poisson integral with the influence of the high-frequency far zone in the removal and restoration mode is the same as the corresponding formula for the rewind-limited Poisson integral considering the far-zone influence. For the non-diagonal elements of the discretization coefficient matrix of the rewind-limited Poisson integral, regardless of the discretization algorithm and calculation mode used, the calculation formula is the same.

[0152] Numerical experiment

[0153] The discrete gravity point values obtained from actual measurements can contain gravity information in the full frequency band, but the spectral content that can be expressed by the gridded gravity data is limited by the spatial distribution of the discrete gravity points and the highest order of the spectral domain that the grid can accommodate. Among them, the highest order that the grid can accommodate is related to the grid resolution. In order to conduct the downward continuation experiment of the rewind-limited Poisson integral on gravity data containing different spectral domain information under different resolution conditions, this application first generated the model gravity anomaly data with different resolutions and truncated to different orders on the 0m and 4000m altitude planes as statistically shown in Table 1 using the XGM2019e model. The test area is a mountainous area in the range of 106°E to 109°E and 33°N to 36°N. The terrain of the test area with a resolution of 2′×2′ has a minimum of 365.24m, a maximum of 3423.56m, an average value of 1300.84m, and a root mean square value of 1383.16m.

[0154] Table 1 Statistics of model gravity anomalies truncated to different orders at different resolutions (unit: mGal)

[0155]

[0156]

[0157] As can be seen from Table 1, the signal strength of gravity data mainly depends on the truncation order of the gravity data. The higher the truncation order, the stronger the signal of the gravity data on the same height plane. When conducting the downward continuation test, the model gravity anomaly data without error and with a 2 mGal random error at a height of 4000 m are used as the gravity data to be continued respectively, while the model gravity anomaly data at a height of 0 m are used as the verification data. In the test, the Poisson integral radius is taken as 1°, and to avoid the influence of edge effects, the data within a range of 1° from the edge of the test area are removed after continuation. The error in the downward continuation using the discretized Poisson integral includes two parts: the error of the algorithm itself and the error caused by measurement noise. The error in the downward continuation of noiseless gravity data actually reflects the error of the algorithm itself. To compare and illustrate the downward continuation effect of the discretized backward-limited Poisson integral, this application first conducts a downward continuation test using the standard Poisson integral considering the influence of the far zone. Table 2 lists the standard deviation statistics of the continuation error, where M represents the truncation order of the gravity anomaly data (the same hereinafter).

[0158] Table 2 Downward Continuation Error of the Standard Poisson Integral (unit: mGal)

[0159]

[0160] The downward continuation process has an amplifying effect on the noise of gravity data. As can be seen from Table 2, when using the standard Poisson integral to conduct downward continuation on gravity data with the same resolution and different truncation orders, the amplification degree of measurement noise is very close. This is because when using the standard Poisson integral without filtering attributes for downward continuation, the noise of the gravity data is amplified to the limiting frequency corresponding to its grid resolution. According to the Nyquist sampling theorem, the highest order of the frequency domain (limiting frequency) that can be accommodated by grid gravity data with resolutions of 5′, 4′, 3′, and 2′ are 2160 orders, 2700 orders, 3600 orders, and 5400 orders respectively. The higher the resolution of the gravity data, the more serious the amplifying effect of measurement noise. The test in Table 2 shows that it is not advisable to use the standard Poisson integral for downward continuation of high-resolution gravity data containing noise. Table 3 lists the standard deviation of the error in the downward continuation of gravity data with different resolutions using the discretized backward-limited Poisson integral in the non-remove-and-recover mode. Among them, N 1 is taken as -1, and N 2 is taken as M, and the cut-off order of the far zone influence in the integral mode considering the far zone influence is taken as 360.

[0161] Table 3 Downward Continuation Error of the Backward-Limited Poisson Integral in the Non-Remove-and-Recover Mode (unit: mGal)

[0162]

[0163] Comparing Table 3 and Table 2, it can be seen that when using N 2 The downward continuation accuracy of the rewind-limited Poisson integral of M for non-full-order gravity data with measurement noise is significantly higher than that of the standard Poisson integral. This is because when using the standard Poisson integral for downward continuation, the measurement noise of the gravity data is amplified to the grid limit frequency, while the noise of the gravity data in the downward continuation value obtained by the rewind-limited Poisson integral with filtering properties is at most amplified to N 2 order. The downward continuation accuracy of the rewind-limited Poisson integral using the direct discretization method for non-full-order gravity data is significantly higher than its continuation accuracy for full-order gravity data of the same resolution, which reflects that N 2 When approaching the grid limit frequency, the diagonal element formula of the direct discretization coefficient matrix is 0.5∑(2n + 1)(r / R) n+2 I n (ψ 0 ) has a very large discretization error. When using the direct discretization method, the difference in the downward continuation error between the rewind-limited Poisson integral in the near zone and considering the influence of the far zone is very small, indicating that this term has a very small influence on the downward continuation solution. However, when using the indirect discretization method, there is a significant difference in the downward continuation accuracy of the rewind-limited Poisson integral in these two calculation modes. From this, it can be inferred that in the diagonal element formula of the indirect discretization coefficient matrix, 0.5∑(2n + 1)(r / R) n+2 I n (ψ 1 ) has a very large influence on the downward continuation accuracy. Since the magnitude of this term is independent of the grid resolution, the indirect discretization algorithm of the rewind-limited Poisson integral considering the influence of the far zone is not very sensitive to the resolution change of the gravity data truncated to the same order. 0.5∑(2n + 1)(r / R) n+2 I n (ψ 1 ) The order error increases with the increase of the order, and its order error is the largest at the highest order (i.e., N 2 order). Therefore, when using the indirect discretization method for the rewind-limited Poisson integral considering the influence of the far zone, the difference in the downward continuation error calculated by the kernel functions with the same N 2 value under different resolution conditions is relatively small. As shown in Table 3, when using the rewind-limited Poisson integral with indirect discretization for downward continuation, the continuation accuracy of the near zone integral mode is always significantly higher than that of the integral mode considering the influence of the far zone. This phenomenon indicates that when using the indirect discretization method for the rewind-limited Poisson integral, in the integral mode considering the influence of the far zone, 0.5∑(2n + 1)(r / R) n+2 I n (ψ 1) The error of this term far exceeds the far - zone truncation error caused by neglecting this term in the near - zone integration mode. According to the test results in Table 3, for non - full - order gravity data with a cut - off order significantly lower than the grid aliasing frequency, in the non - remove - restore mode, a high - precision downward continuation solution can be obtained by using the direct discretization formula of the rewind - limited Poisson integral considering the far - zone influence; while for full - order or nearly full - order gravity data, the extension effect of applying the indirect discretization formula of the near - zone rewind - limited Poisson integral is the most ideal. Finally, Table 4 statistics the standard deviation of the extension error of the rewind - limited Poisson integral for discretizing gravity data with 2 mGal random noise under different resolution conditions in the remove - restore mode. Among them, N takes 720, L takes 360, N 1 equals L, N 2 takes the value of M.

[0164] Table 4 Extension error of the rewind - limited Poisson integral in the remove - restore mode (unit: mGal)

[0165]

[0166] It can be seen from Table 4 that, like the non - remove - restore mode, when performing the rewind - limited Poisson integral on non - full - order gravity data with a cut - off order much smaller than the grid aliasing frequency in the remove - restore mode, a high - precision downward continuation solution can be obtained by using the direct discretization method. Considering that the kernel function with N 1 equals L has low - order spectral leakage, at this time, the rewind - limited Poisson integral including the influence of the high - frequency far - zone should be carried out. When using the indirect discretization method to perform the rewind - limited Poisson integral in the remove - restore mode, for gravity data with a truncation order M equal to 2160, the downward continuation accuracy including the influence of the high - frequency far - zone is slightly higher than that without including the influence of the high - frequency far - zone; for gravity data with M exceeding 2700, the downward continuation accuracy without including the influence of the high - frequency far - zone is significantly higher than that including the influence of the high - frequency far - zone. This shows that when N 2 is large, 0.5∑(2n + 1)(r / R) n+2 I n (ψ 1) The influence of the error of item [] on the downward continuation accuracy exceeds that of the spectral leakage phenomenon. Based on the test results in Table 3 and Table 4, it can be concluded that when downward continuing full-order or nearly full-order gravity data, the best continuation effect is achieved by using the indirect discretization formula of the near-zone rewind-limited Poisson integral; while for non-full-order gravity data with a cut-off order significantly lower than the grid aliasing limit frequency, it is recommended to use the direct discretization formula of the rewind-limited Poisson integral considering the influence of high-frequency far zone in the remove-and-recover mode for downward continuation. According to the test results in Table 3 and Table 4, it can be inferred that for a 3′ resolution grid, when the difference between the grid aliasing limit frequency and the cut-off order of the grid gravity data is greater than or equal to 900 orders, the direct discretization method considering the influence of high-frequency far zone in the remove-and-recover mode should be used; for a 2′ resolution grid, the difference threshold is approximately 2200 orders. The higher the grid resolution, the larger the difference threshold.

[0167] The traditional discretization algorithm of the rewind-limited Poisson integral is an indirect discretization formula that includes the influence of high-frequency far zone in the remove-and-recover mode. According to the tests of this application, when downward continuing non-full-order gravity data truncated to 2700 orders and containing 2 mGal random noise under the 2′ resolution condition, the continuation accuracy of the rewind-limited Poisson integral considering the influence of high-frequency far zone using the direct discretization method in the remove-and-recover mode is 68.27% higher than that of the traditional algorithm. For non-full-order gravity data truncated to 3600 orders and containing 2 mGal random noise under the 2′ resolution condition, the continuation accuracy of the near-zone rewind-limited Poisson integral using the indirect discretization method is 59.86% higher than that of the traditional algorithm. When downward continuing full-order gravity data containing 2 mGal random noise under the 5′, 4′, 3′, 2′ resolution conditions, the continuation accuracy of the near-zone rewind-limited Poisson integral using the indirect discretization method is 47.52%, 45.48%, 52.04%, 51.95% higher than that of the traditional algorithm respectively.

[0168] This application also provides an application scenario, which applies the above-mentioned downward continuation method of gravity data based on the discretization model of the rewind-limited Poisson integral. Specifically: The downward continuation method of gravity data based on the discretization model of the rewind-limited Poisson integral provided in this embodiment can be applied to the processing of airborne gravity data. Obtain airborne gravity data, and input the airborne gravity data into the discretization model based on the rewind-limited Poisson integral to obtain ground gravity data.

[0169] In an exemplary embodiment, a computer device is provided, including: a memory, a processor, and a computer program stored on the memory and executable on the processor, and the processor executes the computer program to implement the downward continuation method of gravity data based on the discretization model of the rewind-limited Poisson integral.

[0170] In an exemplary embodiment, a computer device is provided. The computer device can be a server or a terminal, and its internal structure diagram can be as shown in Figure 3 . The computer device includes a processor, a memory, an input / output interface (Input / Output, abbreviated as I / O), and a communication interface. Among them, the processor, the memory, and the input / output interface are connected through a system bus, and the communication interface is connected to the system bus through the input / output interface. Among them, the processor of the computer device is used to provide computing and control capabilities. The memory of the computer device includes a non-volatile storage medium and an internal memory. The non-volatile storage medium stores an operating system, a computer program, and a database. The internal memory provides an environment for the operation of the operating system and the computer program in the non-volatile storage medium. The input / output interface of the computer device is used to exchange information between the processor and external devices. The communication interface of the computer device is used to communicate with an external terminal through a network connection. When the computer program is executed by the processor, it implements a gravity data downward continuation method based on a discretized model of the rewind-limited Poisson integral.

[0171] Those skilled in the art can understand that Figure 3 the structure shown in is only a block diagram of some structures related to the solution of this application, and does not constitute a limitation on the computer device to which the solution of this application is applied. The specific computer device may include more or fewer components than those shown in the figure, or combine some components, or have different component arrangements.

[0172] It should be noted that the user information (including but not limited to user device information, user personal information, etc.) and data (including but not limited to data for analysis, stored data, displayed data, etc.) involved in this application are all information and data authorized by the user or fully authorized by all parties, and the collection, use, and processing of relevant data need to comply with relevant regulations.

[0173] Those of ordinary skill in the art can understand that all or part of the processes of implementing the methods in the above embodiments can be completed by instructing relevant hardware through a computer program. The computer program can be stored in a non-volatile computer-readable storage medium. When the computer program is executed, it can include the processes of the embodiments of the above methods. Among them, any reference to a memory, database, or other medium used in the embodiments provided in the present application can include at least one of non-volatile and volatile memories. Non-volatile memories can include read-only memory (ROM), magnetic tapes, floppy disks, flash memories, optical memories, high-density embedded non-volatile memories, resistive random access memories (ReRAM), magnetoresistive random access memories (MRAM), ferroelectric random access memories (FRAM), phase change memories (PCM), graphene memories, etc. Volatile memories can include random access memory (RAM) or external cache memories, etc. By way of illustration and not limitation, RAM can be in various forms, such as static random access memory (SRAM) or dynamic random access memory (DRAM), etc.

[0174] The databases involved in the embodiments provided in the present application can include at least one of relational databases and non-relational databases. Non-relational databases can include distributed databases based on blockchain, etc., without limitation. The processors involved in the embodiments provided in the present application can be general-purpose processors, central processors, graphics processors, digital signal processors, programmable logics, data processing logics based on quantum computing, etc., without limitation.

[0175] The technical features of the above embodiments can be combined arbitrarily. For the sake of brevity of description, not all possible combinations of the technical features in the above embodiments are described. However, as long as there is no contradiction in the combination of these technical features, it should be considered as the scope described in this specification.

[0176] In the present application, specific examples are used to elaborate on the principles and implementation manners of the present application. The descriptions of the above embodiments are only used to help understand the method and its core idea of the present application; at the same time, for those of ordinary skill in the art, according to the idea of the present application, there will be changes in the specific implementation manners and application scopes. In summary, the content of this specification should not be construed as a limitation to the present application.

Claims

1. A method for downward extension of gravity data based on a discretized model of rewind-limited Poisson integral, characterized in that: The downward extension method of gravity data based on the discretized model of the rewind-limited Poisson integral includes: Obtain the discrete gravity point data, grid resolution and extension height to be extended; The discrete gravity point data are gridded according to the grid resolution to obtain grid gravity anomaly data; Based on the grid resolution, the grid frequency limit is determined, and based on the grid frequency limit and the discrete gravity point data, the truncation order of the grid gravity anomaly data is determined; Determining the type of discretization model of the rewind limited Poisson integral based on the grid frequency limit and the truncation order of the grid gravity anomaly data; the type of the discretization model is a direct discretization model in the removal recovery mode or an indirect discretization model in the near zone integration mode; When the discretization model type is a direct discretization model in remove recovery mode: Obtaining the Earth's gravity field model without the reference of the recovery mode; Based on the grid gravity anomaly data and the extended height, a coefficient matrix of a discretization model is determined; the coefficient matrix includes diagonal elements and off-diagonal elements; Based on grid gravity anomaly data, extended height and earth gravity field model, the parameter vector and constant vector of the discretization model are determined; Based on the coefficient matrix, parameter vector and constant vector of the discretization model, determine the target calculation formula of the discretization model; When the discretization model type is an indirect discretization model in near-zone integration mode: Based on the grid gravity anomaly data and extended height, the coefficient matrix and parameter vector of the discretization model are determined; Based on the coefficient matrix and parameter vector of the discretization model, determine the target calculation formula of the discretization model; Based on the target calculation formula of the corresponding discretization model, the downward extension result is obtained.

2. The method for downward extension of gravity data based on the discretization model of the rewind-limited Poisson integral according to claim 1 is characterized in that: Based on the grid resolution, the grid frequency limit is determined, and based on the grid frequency limit and discrete gravity point data, the truncation order of the grid gravity anomaly data is determined, including: Based on the grid resolution, the grid frequency limit is calculated using the formula M'=π / Δ; where M' represents the grid frequency limit; Δ represents the grid resolution in radians; Based on the grid frequency-limited and discrete gravity point data, the truncation order of the grid gravity data is calculated using the formula M=M'·min{n1 / n2,1}; where M represents the truncation order of the grid gravity data; n1 represents the number of discrete gravity points; and n2 represents the number of grid gravity data.

3. The method for downward extension of gravity data based on the discretization model of the rewind-limited Poisson integral according to claim 2 is characterized in that: The type of discretization model of the rewind limited Poisson integral is determined based on the grid frequency limit and the truncation order of the grid gravity anomaly data, including: Determine whether the difference between the grid frequency limit and the truncation order of the grid gravity data is greater than or equal to a preset threshold, and obtain a determination result; If the judgment result is yes, then the type of the discretization model of the rewind-limited Poisson integral is determined to be a direct discretization model under the removal recovery mode; If the judgment result is no, then it is determined that the type of the discretization model of the rewind-limited Poisson integral is an indirect discretization model in the near-zone integration mode.

4. The method for downward extension of gravity data based on the discretization model of the rewind-limited Poisson integral according to claim 3 is characterized in that: When the discretization model is a direct discretization model in the removal recovery mode: The calculation formula for the off-diagonal elements in the coefficient matrix of the discretized model is in, represents the element in the i-th row and j-th column of the coefficient matrix of the direct discretization model without the recovery mode; s j represents the area of ​​the jth integration grid; N 1，1 represents the low-order cutoff of the rewind-limited Poisson kernel function when the discretization model is a direct discretization model under the removal recovery mode, N 1，1 =L, L represents the truncation order of the reference earth gravity field model; N 1,2 represents the high-order truncation order of the rewind-limited Poisson kernel function when the discretization model is a direct discretization model under the removal recovery mode, N 1,2 =M; R represents the average radius of the earth; r represents the diameter of the calculation point from the center of the earth under the spherical approximation; P n represents the n-th order Legendre polynomial; ψ ij represents the angular distance between the center point of the i-th computational grid and the center point of the j-th integral grid; ψ1 represents the radius of the near-zone integral domain; The calculation formula for the diagonal elements in the coefficient matrix of the discretized model is in, I represents the diagonal element of the i-th row in the coefficient matrix of the direct discretization model with the recovery mode removed; n (ψ0) represents the value of the n-order function with parameter ψ0; ψ0 represents the radius of the circular approximation area of ​​the computational grid; The calculation formula for the parameter vector of the discretized model is in, represents the i-th element in the parameter vector; Δg i Represents the gravity anomaly data of the i-th grid; represents the reference gravity anomaly of the ith grid calculated using the Earth gravity field model; The calculation formula of the constant vector of the discretized model is Among them, b i represents the i-th element in the constant vector; N represents the highest order of high-frequency far-zone influence; Ω i represents the spherical coordinates of the i-th calculation point; represents the nth-order gravity anomaly calculated using the Earth gravity field model; represents the nth-order far-zone cutoff coefficient of the reversal-limited Poisson kernel function; The objective calculation formula of the discretization model is: Δg High (R)=A IBL1 Δg High (r)+b; Where, Δg H igh(R) represents the output vector of the direct discretization model without the recovery mode; A IBL1 represents the coefficient matrix of the direct discretization model with the recovery mode removed; Δg High (r) represents a parameter vector; b represents a constant vector.

5. The method for downward extension of gravity data based on the discretization model of the rewind-limited Poisson integral according to claim 3 is characterized in that: When the discretization model type is an indirect discretization model in near-zone integration mode: The calculation formula for the off-diagonal elements in the coefficient matrix of the discretized model is in, represents the element in the i-th row and j-th column of the coefficient matrix of the indirect discretization model in the near-zone integral mode; R represents the average radius of the earth; r represents the geocentric diameter of the calculation point under the spherical approximation; s j represents the area of ​​the jth integration grid; N 2，1 It represents the low-order cutoff number of the reversal-limited Poisson kernel function when the discretization model type is an indirect discretization model in the near-zone integral mode, N 2，1 =-1; N 2,2 Indicates the high-order truncation order of the reversal-limited Poisson kernel function when the discretization model type is an indirect discretization model in the near-zone integral mode, N 2,2 =M;P n represents the n-th order Legendre polynomial; ψ ij represents the angular distance between the center point of the i-th computational grid and the center point of the j-th integral grid; ψ1 represents the radius of the near-zone integral domain; The calculation formula for the diagonal elements in the coefficient matrix of the discretized model is in, represents the diagonal element of the i-th row in the coefficient matrix of the indirect discretization model in the near-zone integration mode; J represents the number of near-zone grids; The objective calculation formula of the discretization model is: Δg IBL2 (R)=A IBL2 Δg(r); Where, Δg IBL2 (R) represents the output vector of the indirect discretization model in the near-zone integration mode, A IBL2 represents the coefficient matrix of the indirect discretization model in the near-zone integration mode; Δg(r) represents the parameter vector constructed by the grid gravity anomaly data.

6. The method for downward extension of gravity data based on the discretization model of the rewind-limited Poisson integral according to claim 1 is characterized in that: When the type of the discretization model is a direct discretization model in the removal recovery mode, the downward continuation result is obtained based on the target calculation formula of the corresponding discretization model, including: Based on the target calculation formula of the direct discretization model under the removal recovery mode, the reversal limited Poisson integral processing is performed to obtain the output vector of the discretization model; Based on the output vector of the discretized model and the Earth gravity field model, the formula Δg IBL1 (R) = Δg High (R)+Δg G (R), and the output vector after restoring the reference model value is obtained; where Δg IBL1 (R) represents the output vector after restoring the reference model value; Δg High (R) represents the output vector of the direct discretization model without the recovery mode; Δg G (R) a reference gravity anomaly vector on a boundary surface with a radius of R calculated using the Earth gravity field model; Based on the output vector after restoring the reference model value, the data at the edge of the calculation area is eliminated to obtain the downward extension result.

7. The method for downward extension of gravity data based on the discretization model of the rewind-limited Poisson integral according to claim 1 is characterized in that: When the type of the discretization model is an indirect discretization model in the near-zone integral mode, the downward extension results are obtained based on the target calculation formula of the corresponding discretization model, including: Based on the target calculation formula of the direct discretization model under the removal recovery mode, the reversal limited Poisson integral processing is performed to obtain the output vector of the discretization model; Based on the output vector of the discretized model, the data at the edge of the calculation area is eliminated to obtain the downward extension result.

8. A computer device comprising: 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 downward extension method of gravity data based on a discretized model of a rewind-limited Poisson integral as described in any one of claims 1 to 7.

9. A computer-readable storage medium having a computer program stored thereon, characterized in that: When the computer program is executed by a processor, the method for downward extension of gravity data based on a discretized model of rewind-limited Poisson integral described in any one of claims 1 to 7 is implemented.

10. A computer program product, comprising a computer program, characterized in that When the computer program is executed by a processor, the method for downward extension of gravity data based on a discretized model of rewind-limited Poisson integral described in any one of claims 1 to 7 is implemented.

Images