Fast calculation method of near-zone terrain correction value based on gauss-legendre integral

Abstract

The application belongs to the field of geophysical exploration, and particularly relates to a near-zone terrain correction value fast calculation method based on Gauss-Legendre integral, comprising the following steps: taking a survey point as a center, determining a rectangular calculation range of near-zone terrain correction, and obtaining elevation data of terrain nodes in the rectangular calculation range; uniformly dividing the rectangular calculation range into four rectangular sub-regions with equal areas; in each rectangular sub-region, applying two-dimensional Gauss-Legendre integral of a current Gauss integral order to calculate terrain correction values of the rectangular sub-regions respectively; and adding the terrain correction values of the four rectangular sub-regions to obtain a total near-zone terrain correction value. The method significantly reduces the calculation complexity and shortens the calculation time of terrain correction values.

Description

Technical Field

[0001] This application belongs to the field of geophysical exploration, and in particular relates to a rapid calculation method for near-field topographic corrections based on Gauss-Legend integral. Background Technology

[0002] Gravity exploration, as a fundamental and efficient geophysical method, is widely used in mineral resource exploration, oil and gas detection, and hydrological engineering surveys. In surface gravity observations, topographic effects caused by terrain undulations significantly interfere with gravity values ​​at measuring points. To obtain gravity field data that reflects the true density anomalies underground, rigorous topographic corrections must be applied to the observed values. Among these, the gravitational effect produced by near-field topography (typically referring to the area within tens to hundreds of meters of the measuring point) is the most significant, and the accuracy of its correction directly affects the reliability of the final geological interpretation.

[0003] When calculating near-field topographic corrections, the topography is typically discretized into multiple regular prisms to approximate the true continuous topography. For each prism, its gravity effect needs to be calculated using triple volume integrals, which is computationally intensive, inefficient, and fails to meet the needs of rapid processing of large-scale gravity measurement network data.

[0004] Therefore, researching a fast algorithm that can significantly improve computational efficiency while ensuring the accuracy of near-area topographic correction calculations is of great practical significance and application value for promoting the efficient processing of gravity exploration data. Summary of the Invention

[0005] This application provides a fast calculation method for near-field terrain correction based on Gauss-Legendal integral, which solves the problem of large computational load and low efficiency caused by traditional near-field terrain correction due to its reliance on grid subdivision.

[0006] This application provides a method for fast calculation of near-field terrain correction values ​​based on Gauss-Legends integral, including: Centered on the measuring point, determine the rectangular calculation range for the near-area terrain correction, and obtain the elevation data of each terrain node within the rectangular calculation range; The rectangular calculation range is evenly divided into four rectangular sub-regions with equal areas; Initialize the Gaussian integral order to 2 and set the convergence threshold; Within each rectangular sub-region, apply the two-dimensional Gaussian-Legendal integral of the current Gaussian integral order to calculate the terrain correction value for each rectangular sub-region; The terrain correction values ​​of the four rectangular sub-regions are added together to obtain the total near-field terrain correction value; Determine whether the absolute error of the total near-field terrain correction value calculated in two adjacent iterations is less than a preset threshold. If so, stop the iteration and use the Gaussian integral order of the last iteration as the final Gaussian integral order. If not, increase the Gaussian integral order by 1 and return to the step of calculating the terrain correction value of each rectangular sub-region to continue the calculation.

[0007] Further, obtaining the elevation data of each terrain node within the rectangular calculation range specifically includes: Load the raw terrain data of the near area, which includes a terrain elevation matrix and a corresponding coordinate grid matrix, the coordinate grid matrix including a horizontal coordinate grid matrix and a vertical coordinate grid matrix; Based on the planar coordinates of the measuring points and the preset extension distance, the boundary of the rectangular calculation range is determined. Based on the boundary of the rectangular calculation range, the coverage area of ​​the rectangular calculation range in the coordinate grid matrix is ​​determined. The boundary includes the starting coordinates of the horizontal axis, the ending coordinates of the horizontal axis, the starting coordinates of the vertical axis, and the ending coordinates of the vertical axis. The terrain nodes are then mapped to the network nodes in the coordinate grid matrix. The column index range is determined by comparing the coordinate values ​​in the horizontal coordinate grid matrix with the boundaries in the horizontal direction; The row index range is determined by comparing the coordinate values ​​in the vertical coordinate grid matrix with the boundaries of the vertical axis direction; Based on the row index range and the column index range, a corresponding submatrix is ​​extracted from the terrain elevation matrix. The submatrix is ​​the elevation data of the area covered by the rectangular calculation range in the coordinate grid matrix.

[0008] Furthermore, by applying the two-dimensional Gauss-Legendal integral of the current Gaussian integral order, the terrain correction values ​​for each rectangular sub-region are calculated, including: Obtain the coordinates of the Gaussian nodes and their corresponding node weights on the standard one-dimensional interval for the current Gaussian integral order; By using a linear transformation, the standard one-dimensional interval is mapped to the current rectangular sub-region, resulting in two-dimensional Gaussian nodes distributed on the rectangular sub-region. The combined weights corresponding to each two-dimensional Gaussian node are calculated based on the product of the node weights in the horizontal and vertical directions. The contribution value of each two-dimensional Gaussian node to the terrain correction of the measuring point is calculated based on the combined weights. The terrain correction contribution values ​​of all two-dimensional Gaussian nodes within the rectangular sub-region are summed to obtain the terrain correction value of the current rectangular sub-region.

[0009] Furthermore, for the actual planar coordinates of each two-dimensional Gaussian node, the Gaussian elevation of the two-dimensional Gaussian node is obtained through cubic interpolation, including: Determine the position of the two-dimensional Gaussian node in the coordinate grid matrix corresponding to the rectangular calculation range; The Gaussian elevation at the two-dimensional Gaussian node is calculated using the bicubic interpolation method and the original elevation data of several adjacent terrain nodes around the location.

[0010] Furthermore, the contribution of each two-dimensional Gaussian node to the terrain correction of the measurement point is calculated, including: Calculate the Jacobian determinant of the linear coordinate transformation from the standard interval to the current rectangular subregion; Calculate the integral kernel function for terrain correction; The terrain correction contribution value is obtained by multiplying the integral kernel function, Jacobian determinant, gravitational constant, medium density, and combined weight.

[0011] Furthermore, the integral kernel function is: , in, Here are the plane coordinates of the measuring point, and the reference point for terrain correction calculation. The elevation of the measuring point, These are the actual coordinates of the two-dimensional Gaussian nodes. Gaussian elevations for two-dimensional Gaussian nodes. The result of the integral kernel function, Represents the first [unit] within the rectangular sub-region One integral sampling point.

[0012] Furthermore, the Jacobian determinant is: ,in, The x-coordinate range of the rectangular sub-region. The range of the ordinates of the rectangular sub-region. This is the result of the Jacobian determinant.

[0013] Compared with the prior art, the advantages of this application are as follows: The proposed method significantly reduces computational complexity and shortens the computation time for terrain correction values. It alleviates the efficiency bottleneck caused by the reliance on grid partitioning in traditional terrain correction methods. Its computational advantage lies not only in the effective control of computational accuracy but also in transforming the terrain correction integration process into efficient standard interval calculations through the combination of Gaussian nodes and weights. This reduces the data processing intensity and time cost of field gravity exploration, providing a reliable technical path for the theoretical refinement and engineering application of high-resolution terrain impact correction. Attached Figure Description

[0014] Figure 1This is a schematic diagram of the algorithm flow of a method for fast calculation of near-field terrain correction values ​​based on the Gauss-Legend integral formula provided in an embodiment of this application; Figure 2 This is a schematic diagram of the Gaussian nodes and weights of a standard one-dimensional interval provided in this application, mapped to the actual region. Detailed Implementation

[0015] To make the objectives, technical solutions, and advantages of this application clearer, the following detailed description is provided in conjunction with embodiments. It should be understood that the specific embodiments described herein are merely illustrative and not intended to limit the scope of this application.

[0016] Near-field topographic correction refers to the correction process that eliminates the influence of topographic undulations within a certain range around the measuring point on gravity observation values. With the continuous improvement of the accuracy of gravity instruments and geodetic instruments, the accuracy of topographic correction has become a decisive factor affecting the accuracy of gravity measurement data, especially in mountainous areas with complex terrain and large elevation differences.

[0017] Traditional near-field terrain correction methods rely heavily on grid density for computational efficiency. To ensure accuracy, fine grids are often required, leading to an exponential increase in computational load and extremely long computation time over large computational areas.

[0018] This application addresses the aforementioned problems by using the Gauss-Legendal integral formula to calculate near-field terrain correction values, thereby significantly improving computational efficiency while maintaining accuracy, and thus enhancing the computer's computational efficiency and running speed.

[0019] See Figure 1 As shown, a fast method for calculating near-field terrain correction values ​​based on Gauss-Legends integrals is presented. This method includes the following steps: S101, with the measuring point as the center, determine the rectangular calculation range of the near-area terrain correction, and obtain the elevation data of each terrain node within the rectangular calculation range.

[0020] Gravimeters measure the Earth's total gravitational pull, which includes the gravitational pull from the normal strata below the measuring point, as well as the additional gravitational pull from the topographic mass above or below the measuring point's plane. When the measuring point is located on a topographic protrusion, the excess mass above generates an upward additional gravitational pull, causing the measured gravity value to be lower. Conversely, when the measuring point is located on a topographic depression, the missing mass below generates a downward missing gravitational pull, causing the measured gravity value to be higher. The essence of near-field topographic correction is to eliminate the additional influence of topographic undulations around the measuring point on gravity observations.

[0021] The planar coordinates of the measuring point serve as the spatial reference point for near-area shape correction calculations. The positional descriptions of all terrain nodes are relative to this measuring point. For any terrain node, its planar position relative to the measuring point can be represented by horizontal distance and azimuth. The closer the terrain node, the greater its contribution to the near-area shape correction value. The elevation of the measuring point serves as a reference surface for measuring terrain undulation.

[0022] The rectangular calculation range for near-field topographic correction is determined based on the scale of the gravity survey and the complexity of the terrain surrounding the measuring point. Specifically, it is determined using the plane coordinates of the measuring point (…). Centered on the horizontal axis, extend a preset distance in both the positive and negative directions. This yields the starting coordinates in the horizontal direction. End coordinates in the horizontal direction Extend a preset distance along the vertical axis in both the positive and negative directions. This yields the initial coordinates along the vertical axis. End coordinates along the vertical axis .Depend on , , , The four boundary coordinates together form a rectangular calculation range.

[0023] For example, in a gravity measurement at a certain work area, based on a 1:2000 scale topographic map, the near-area topographic correction range is set to 200 meters × 200 meters, that is, extending 100 meters outwards from the measuring point. rice, Meters, thus determining the boundary of the rectangular calculation range.

[0024] In this embodiment, the original terrain data of the near area is loaded, including the terrain elevation matrix. and the corresponding coordinate grid matrix ( .

[0025] Topographic Elevation Matrix Used to store the terrain elevation values ​​at each terrain node. The terrain elevation size is M×N, where M is the number of rows (corresponding to the number of grid nodes in the vertical axis direction) and N is the number of columns (corresponding to the number of grid nodes in the horizontal axis direction).

[0026] x-axis grid matrix : Same dimension as the terrain elevation matrix, used to store the x-coordinate values ​​of each terrain node; y-coordinate grid matrix. : Same dimension as the terrain elevation matrix, used to store the ordinate values ​​of each terrain node.

[0027] During the actual loading process, the terrain elevation matrix, horizontal coordinate grid matrix, and vertical coordinate grid matrix are loaded into the computer memory by reading the pre-stored data file for subsequent calculations.

[0028] Based on the determined rectangular computational range boundary, the corresponding boundary coordinates are extracted from the coordinate grid matrix to determine the coverage area of ​​the rectangular computational range in the coordinate grid.

[0029] The specific extraction method is as follows: First, determine the column index range along the horizontal axis, and then traverse the horizontal axis grid matrix. The first row retrieves the x-coordinate values ​​for each column. Then, the x-coordinate values ​​of each column are compared with the starting coordinates of the x-axis. End coordinates in the horizontal direction Compare and filter out the starting coordinates whose x-coordinate values ​​are greater than or equal to those in the x-axis direction. And less than or equal to the end coordinate in the horizontal direction The columns record the indexes of these columns, forming a column index range.

[0030] Next, determine the row index range along the vertical axis. Traverse the vertical axis grid matrix. In the first column, obtain the y-coordinate value for each row. Then, compare the y-coordinate value of each row with the starting coordinate of the y-axis. End coordinates along the vertical axis Compare and filter out the starting coordinates whose ordinate values ​​are greater than or equal to those in the ordinate direction. And less than or equal to the end coordinate in the vertical direction The rows are recorded, and the indexes of these rows are recorded to form a row index range.

[0031] Based on the determined row index range and column index range, from the terrain elevation matrix Extract the corresponding submatrix, which contains the elevation data of the area covered by the rectangular calculation range in the coordinate grid matrix. This submatrix corresponds one-to-one with the rectangular calculation range, where the value of each element is the elevation at the corresponding terrain node. The terrain nodes are mapped to network nodes in the coordinate grid matrix, and there is a one-to-one correspondence between them.

[0032] S102, the rectangular calculation range is evenly divided into four rectangular sub-regions with equal areas; Specifically, the division of rectangular sub-regions includes the following steps: Calculate the horizontal center point of the rectangle's calculation range. , , , The x-coordinate of the horizontal center point The ordinate is the horizontal center point.

[0033] Based on the horizontal center point, the rectangular calculation range is divided into four rectangular sub-regions: First sub-region: to , to Second sub-region: to , to Third sub-region: to , to Fourth sub-region: to , to .

[0034] S103, Initialize the Gaussian integral order n=2, set the convergence threshold. ; Before performing iterative calculations for near-field terrain correction, the Gaussian integral order is first initialized. The Gaussian integral order n determines the number of Gaussian nodes used for numerical integration within each rectangular sub-region, specifically employing n×n two-dimensional Gaussian nodes for integration calculations.

[0035] In this embodiment, the Gaussian integral order is initialized to n=2. The reasons for choosing n=2 as the starting order are as follows: Prioritizing computational efficiency: Starting the calculation from the lowest order avoids wasting computational resources on high-order integrals in flat terrain areas with low accuracy requirements; Iterative increment mechanism: For complex terrain areas, the order is automatically increased through subsequent iterations until the accuracy requirements are met, achieving adaptive adjustment of the order; Avoiding overcomputation: For most flat terrain areas, the second-order Gaussian integral can obtain sufficiently accurate terrain correction values, eliminating the need for high-order calculations.

[0036] The convergence threshold is used to determine whether the absolute error of the total near-field terrain correction value calculated in two adjacent iterations meets the accuracy requirements.

[0037] S104. Within each rectangular sub-region, apply the two-dimensional Gauss-Legend integral of the current Gaussian integral order to calculate the terrain correction value for each rectangular sub-region. In this embodiment, for each rectangular sub-region, a two-dimensional Gaussian-Legendal integral is used to calculate the terrain correction value, including: Obtain the coordinates of the Gaussian nodes and their corresponding node weights on the standard one-dimensional interval for the current Gaussian integral order; By mapping the standard one-dimensional interval to the current rectangular sub-region through a linear transformation, two-dimensional Gaussian nodes distributed on the rectangular sub-region are obtained. See [link to documentation]. Figure 2As shown, multiple two-dimensional Gaussian nodes are distributed in a rectangular sub-region formed by the horizontal and vertical directions. The combined weights corresponding to each two-dimensional Gaussian node are calculated based on the product of the node weights in the horizontal and vertical directions. The contribution value of each two-dimensional Gaussian node to the terrain correction of the measuring point is calculated based on the combined weights. The terrain correction contribution values ​​of all two-dimensional Gaussian nodes within the rectangular sub-region are summed to obtain the terrain correction value of the current rectangular sub-region.

[0038] The specific process is as follows: Obtaining the order of Gaussian integrals Obtain the coordinates of Gaussian nodes on the standard one-dimensional interval [-1,1]. and node weight Where i = 1, 2, 3, ..., n. Gaussian node coordinates. It is the zero of the nth Legendre polynomial on the standard one-dimensional interval [-1, 1].

[0039] By using a linear transformation to map the standard one-dimensional interval [-1, 1] to the current rectangular sub-region, we obtain the distribution within the rectangular sub-region. Two-dimensional Gaussian nodes ( ), where i=1,2,3,…,n,j=1,2,3,…,n, and the actual planar coordinates of each two-dimensional Gaussian node ( At point ), the Gaussian elevation of the two-dimensional Gaussian node is obtained through cubic interpolation. .

[0040] For each 2D Gaussian node, the Gaussian elevation is obtained through cubic interpolation at its actual planar coordinates, including: Determine the position of the two-dimensional Gaussian node in the coordinate grid matrix corresponding to the rectangular calculation range; The Gaussian elevation at the two-dimensional Gaussian node is calculated using the bicubic interpolation method and the original elevation data of several adjacent terrain nodes around the location.

[0041] Specifically, firstly, an index mapping relationship for the coordinate grid matrix is ​​established based on the coordinate boundaries and grid spacing of the rectangular calculation area. For a given two-dimensional Gaussian node's actual planar coordinates, its floating-point row and column index values ​​in the coordinate grid matrix are calculated, thereby determining the location of the enclosing cell of the two-dimensional Gaussian node in the coordinate grid matrix. Based on this enclosing cell, four adjacent 4×4 terrain nodes are selected to form the interpolation support region, ensuring that the 16 elevation data points required for bicubic interpolation completely cover the target node.

[0042] Calculate the terrain correction contribution of each 2D Gaussian node to the measuring point, including: Calculate the Jacobian determinant of the linear coordinate transformation from the standard interval to the current rectangular subregion; Calculate the integral kernel function for terrain correction; The terrain correction contribution value is obtained by multiplying the integral kernel function, Jacobian determinant, gravitational constant, medium density, and combined weight.

[0043] Among them, the integral kernel function of terrain correction is calculated based on the elevation relationship between the measuring point and the two-dimensional Gaussian node; , in, Here are the plane coordinates of the measuring point, and the reference point for terrain correction calculation. The elevation of the measuring point, These are the actual coordinates of the two-dimensional Gaussian nodes. Gaussian elevations for two-dimensional Gaussian nodes. The result of the integral kernel function, Represents the first [unit] within the rectangular sub-region One integral sampling point.

[0044] Calculate the Jacobian determinant ,in, The x-coordinate range of the rectangular sub-region. The range of the ordinates of the rectangular sub-region. This is the result of the Jacobian determinant.

[0045] Calculate the contribution of the 2D Gaussian node to the terrain correction value of the measuring point. ,in The gravitational constant is... = For combined weights, The node weights are in the horizontal direction. The weights are the node weights in the vertical direction. The combined weights for each two-dimensional Gaussian node are calculated by multiplying the node weights in the horizontal direction with the node weights in the vertical direction.

[0046] The terrain correction value for the rectangular sub-region is obtained by summing the terrain correction contributions of all two-dimensional Gaussian nodes within that sub-region. .

[0047] S105, add the terrain correction values ​​of the four rectangular sub-regions to obtain the total near-field terrain correction value. ; in, , , , These correspond to the terrain correction values ​​for the four rectangular sub-regions. This summation process achieves a complete integral calculation of the terrain effect in the near-field rectangular region centered on the measurement point. The division into four rectangular sub-regions ensures the completeness of the integration area coverage while maintaining the independence and parallelism of the integration calculations in each rectangular sub-region.

[0048] S106: Determine whether the absolute error of the total near-field terrain correction value calculated in two consecutive iterations is less than a preset threshold. If so, stop the iteration and take the Gaussian integral order of the last iteration as the final Gaussian integral order. If not, increase the Gaussian integral order by 1 and return to the step of calculating the terrain correction value for each rectangular sub-region to continue the calculation.

[0049] When the calculation results tend to stabilize (i.e., the difference between adjacent order results is less than a preset threshold), the iteration stops, and the current Gaussian integral order is used as the final Gaussian integral order.

[0050] In summary, the method presented in this application, based on the Gauss-Legendal integral formula, achieves rapid and high-precision calculation of near-field terrain correction values, alleviating the efficiency bottleneck caused by the reliance on grid partitioning in traditional terrain correction methods. Its computational advantages lie not only in the effective control of computational accuracy but also in transforming the terrain correction integration process into efficient standard interval calculations through the combination of Gaussian integral nodes and weights. This reduces the data processing intensity and time cost of field gravity exploration, providing a reliable technical path for the theoretical improvement and engineering application of high-resolution terrain impact correction.

[0051] The above description is only a preferred embodiment of this application and is not intended to limit this application. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of this application should be included within the protection scope of this application.

Claims

1. A fast calculation method of near-zone terrain correction value based on Gauss-Legendre integral, characterized in that, include: Centered on the measuring point, determine the rectangular calculation range for the near-area terrain correction, and obtain the elevation data of each terrain node within the rectangular calculation range; The rectangular calculation range is evenly divided into four rectangular sub-regions with equal areas; Initialize the Gaussian integral order to 2 and set the convergence threshold; Within each rectangular sub-region, apply the two-dimensional Gaussian-Legendal integral of the current Gaussian integral order to calculate the terrain correction for each rectangular sub-region, including: Obtain the coordinates of the Gaussian nodes and their corresponding node weights on the standard one-dimensional interval for the current Gaussian integral order; By using a linear transformation, the standard one-dimensional interval is mapped to the current rectangular sub-region, resulting in two-dimensional Gaussian nodes distributed on the rectangular sub-region. The combined weights corresponding to each two-dimensional Gaussian node are calculated based on the product of the node weights in the horizontal and vertical directions. The contribution value of each two-dimensional Gaussian node to the terrain correction of the measuring point is calculated based on the combined weights. The terrain correction contribution values ​​of all two-dimensional Gaussian nodes within the rectangular sub-region are summed to obtain the terrain correction value of the current rectangular sub-region. The terrain correction values ​​of the four rectangular sub-regions are added together to obtain the total near-field terrain correction value. Determine whether the absolute error of the total near-field terrain correction value calculated in two adjacent iterations is less than a preset threshold. If so, stop the iteration and use the Gaussian integral order of the last iteration as the final Gaussian integral order. If not, increase the Gaussian integral order by 1 and return to the step of calculating the terrain correction value of each rectangular sub-region to continue the calculation.

2. The method for fast calculating near-zone terrain correction value based on Gauss-Legendre integral according to claim 1, characterized in that, Obtaining the elevation data of each terrain node within the rectangular calculation range specifically includes: Load the raw terrain data of the near area, which includes a terrain elevation matrix and a corresponding coordinate grid matrix, the coordinate grid matrix including a horizontal coordinate grid matrix and a vertical coordinate grid matrix; Based on the planar coordinates of the measuring points and the preset extension distance, the boundary of the rectangular calculation range is determined. Based on the boundary of the rectangular calculation range, the coverage area of ​​the rectangular calculation range in the coordinate grid matrix is ​​determined. The boundary includes the starting coordinates of the horizontal axis, the ending coordinates of the horizontal axis, the starting coordinates of the vertical axis, and the ending coordinates of the vertical axis. The terrain nodes are then mapped to the network nodes in the coordinate grid matrix. The column index range is determined by comparing the coordinate values ​​in the horizontal coordinate grid matrix with the boundaries in the horizontal direction; The row index range is determined by comparing the coordinate values ​​in the vertical coordinate grid matrix with the boundaries of the vertical axis direction; Based on the row index range and the column index range, a corresponding submatrix is ​​extracted from the terrain elevation matrix. The submatrix is ​​the elevation data of the area covered by the rectangular calculation range in the coordinate grid matrix.

3. The method for fast calculating near-zone terrain correction value based on Gauss-Legendre integral according to claim 2, characterized in that, For each 2D Gaussian node, the Gaussian elevation is obtained through cubic interpolation at its actual planar coordinates, including: Determine the position of the two-dimensional Gaussian node in the coordinate grid matrix corresponding to the rectangular calculation range; The Gaussian elevation at the two-dimensional Gaussian node is calculated using the bicubic interpolation method and the original elevation data of several adjacent terrain nodes around the location.

4. The method for rapid calculation of near-field terrain correction values ​​based on Gauss-Legends integral according to claim 3, characterized in that, Calculate the terrain correction contribution of each 2D Gaussian node to the measuring point, including: Calculate the Jacobian determinant of the linear coordinate transformation from the standard interval to the current rectangular subregion; Calculate the integral kernel function for terrain correction; The terrain correction contribution value is obtained by multiplying the integral kernel function, Jacobian determinant, gravitational constant, medium density, and combined weight.

5. The method for rapid calculation of near-field terrain correction values ​​based on Gauss-Legendal integrals according to claim 4, characterized in that, The integral kernel function is: ， in, Here are the plane coordinates of the measuring point, and the reference point for terrain correction calculation. The elevation of the measuring point, These are the actual coordinates of the two-dimensional Gaussian nodes. Gaussian elevations for two-dimensional Gaussian nodes. The result of the integral kernel function, Represents the first [unit] within the rectangular sub-region One integral sampling point.

6. The method for rapid calculation of near-field terrain correction values ​​based on Gauss-Legends integral according to claim 4, characterized in that, The Jacobian determinant is: ,in, The x-coordinate range of the rectangular sub-region. The range of the ordinates of the rectangular sub-region. This is the result of the Jacobian determinant.

Images