Geoid calculation method, system, equipment and medium
By combining the Earth's gravity field model with measured gravity anomaly data, the removal-restoration method was used to calculate the marine geoid, solving the problem of high-precision measurement in marine areas and achieving high-precision and high-resolution calculation of the geoid.
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- NAT MARINE DATA & INFORMATION SERVICE
- Filing Date
- 2026-01-27
- Publication Date
- 2026-05-15
AI Technical Summary
Existing technologies cannot achieve high-precision geoid measurement in marine areas. The resolution and accuracy of global gravity field models are insufficient, and existing methods lack a systematic data processing framework, resulting in unreliable marine geoid calculation results.
By combining the Earth's gravity field model with measured gravity anomaly data, and using the removal-recovery method, the model perturbation potential and model gravity anomaly are first calculated. Then, the measured gravity anomaly is subtracted from the model gravity anomaly to obtain the residual gravity anomaly. Finally, the residual geoid height is calculated by Stokes integration and summed with the model geoid height.
It has achieved the acquisition of a high-precision geoid model in the ocean region, combining the overall advantages of the global gravity field model with the ability to depict the details of local gravity data, thereby improving the reliability and operability of the model.
Smart Images
Figure CN122045562A_ABST
Abstract
Description
Technical Field
[0001] This invention relates to the field of marine surveying technology, and more specifically to a method, system, equipment, and medium for calculating the geoid. Background Technology
[0002] The accurate determination of the marine geoid is a fundamental task in fields such as marine surveying, geophysics, and space navigation. Its results directly serve engineering and scientific research such as seabed topography inversion, marine resource exploration, and satellite altimetry calibration.
[0003] In existing technologies, the determination of the geoid mainly relies on the mature GPS / leveling method for land or nearshore areas. This method obtains the geoid anomaly by connecting the geometric leveling height of the measurement points with the GPS geodetic height. However, in the vast ocean areas, it is almost impossible to carry out high-precision geometric leveling measurements in engineering, which makes this method unsuitable for direct application to the ocean, leaving a large data gap. Although global gravity field models (such as EGM2008) can provide global gravity potential information, their own resolution and accuracy are limited, especially in reflecting the fine structure of the local gravity field in the ocean. If they are directly used to calculate the ocean geoid, it will lead to the loss of shortwave signals, making it difficult to meet the requirements of high-precision applications.
[0004] To address the issue of missing ocean data, existing technologies often employ single-source gravity data combined with potential models for estimation. However, these methods often lack a systematic data processing framework. They fail to establish a standardized and rigorously derived process for effectively integrating models and observational data of different scales and accuracies. For example, from observational data preprocessing and coordinate transformation to the separation and fusion of model and observational values, the error propagation and control between each stage lack clear definition, making it difficult to guarantee the reliability and consistency of the final results.
[0005] Furthermore, in the core step of using gravity anomalies to invert the geoid, namely the application of Stokes integration, existing methods often face specific technical challenges such as how to reasonably handle the integration kernel function, determine the integration radius, and optimize algorithm efficiency. Many implementation schemes either oversimplify the integration process and introduce errors, or have excessively high computational complexity and are difficult to use practically. They have failed to form a rigorous and efficient computational system, which restricts their practical application in large-scale, high-resolution marine geoid modeling.
[0006] Therefore, how to design a geoid calculation method, system, equipment and medium that integrates the Earth's gravity field model with multi-source ocean gravity observation data, and make up for the shortcomings of existing technologies in terms of processing flow, accuracy and practicality, is a problem that urgently needs to be solved by those skilled in the art. Summary of the Invention
[0007] In view of this, the present invention provides a geoid calculation method, system, equipment and medium, which aims to solve the problem that the traditional geometric leveling method cannot be implemented in the ocean. By integrating the Earth's gravity field model and measured gravity anomaly data, a logically clear and step-by-step calculation system is constructed, and finally a high-precision geoid model that can reflect the gravity field characteristics of the ocean area is obtained.
[0008] To achieve the above objectives, the present invention adopts the following technical solution:
[0009] In a first aspect, the present invention provides a method for calculating a geoid, comprising the following steps: S1. Convert the geographic coordinates of the observation point to geocentric coordinates to obtain the geocentric sphere coordinate parameters; S2. Based on the geocentric coordinate parameters, calculate the model perturbation potential of the observation point using the Earth's gravity field model; S3. Based on the model perturbation position, calculate the model gravity anomaly at the observation point; S4. Subtract the measured gravity anomaly from the model gravity anomaly to obtain the residual gravity anomaly at the observation point; S5. Substitute the residual gravity anomaly into the Stokes integral model to calculate the residual geoid height; S6. Based on the Earth's gravity field model, the height of the model's geoid is calculated; S7. Sum the residual geoid height with the model geoid height to obtain the geoid height.
[0010] Preferably, S1 includes: Convert geographic coordinates (longitude L, latitude B, elevation H) to geocentric coordinates (geocentric distance r, geocentric latitude). and longitude :
[0011] Where N is the radius of curvature of the ramusoidal circle, and e is the first eccentricity of the reference ellipsoid.
[0012] Preferably, in S2, the Earth's gravity field model is represented as follows:
[0013] in, For model perturbation bits, For the remaining latitude, GM is the gravitational constant. With reference to the semi-major axis of the ellipsoid, and Here are the fully normalized potential coefficients of the Earth's gravity field model, where n is the order and m is the degree. To fully normalize the association of Legendre functions, This represents the highest order of the model expansion.
[0014] Preferably, in S3, the model exhibits a gravity anomaly. Represented as:
[0015] in, For model perturbation bits, The distance from the center of the Earth is the coordinate system of the Earth's center.
[0016] Preferably, in S5, the Stokes integral model is expressed as:
[0017] in, Let be the residual geoid height, R be the Earth's mean radius, γ be the normal gravity value, and σ be the integration region. For Stokes functions, To calculate the spherical angular distance between the point and the flow surface element, Let be a infinitesimal element of the spherical surface area.
[0018] Preferably, the Stokes function Represented as:
[0019] in, It is an nth-order Legendre polynomial. This is the cutoff coefficient.
[0020] Preferably, in step S6, the model geoid height Represented as:
[0021] in, This is the perturbation bit for the model.
[0022] Secondly, the present invention provides a geoid calculation system, comprising: Coordinate transformation module: used to convert the geographic coordinates of the observation point to geocentric coordinates to obtain geocentric sphere coordinate parameters; Model perturbation potential calculation module: used to calculate the model perturbation potential of the observation point based on the geocentric sphere coordinate parameters and the Earth's gravity field model; Model gravity anomaly calculation module: used to calculate the model gravity anomaly at the observation point based on the model perturbation position; Residual gravity anomaly calculation module: used to subtract the measured gravity anomaly from the model gravity anomaly to obtain the residual gravity anomaly at the observation point; Stokes integral calculation module: used to substitute the residual gravity anomaly into the Stokes integral model to calculate the residual geoid height; Model geoid calculation module: used to calculate the height of the model geoid based on the Earth's gravity field model; Geoid height synthesis module: used to sum the residual geoid height with the model geoid height to obtain the geoid height.
[0023] Thirdly, the present invention provides an electronic device, including a memory, a processor, and a computer program stored in the memory and executable on the processor, wherein the processor executes the computer program to implement the above-described geoid calculation method.
[0024] Fourthly, the present invention provides a computer-readable storage medium storing a computer program that, when executed by a processor, implements the above-described geoid calculation method.
[0025] The descriptions of the second to fourth aspects of this invention can be referred to the detailed description of the first aspect; and the beneficial effects described in the second to fourth aspects can be referred to the analysis of the beneficial effects of the first aspect, which will not be repeated here.
[0026] As can be seen from the above technical solution, compared with the prior art, the present invention has the following beneficial effects: 1. This method addresses the challenge of directly conducting geometric leveling measurements at sea by providing a computational path based entirely on gravity field data. By systematically processing all aspects from observation data to the final geoid height, it enables the acquisition of high-precision geoid models even in vast sea areas lacking traditional leveling data, thus solving a key technical bottleneck in marine surveying.
[0027] 2. The core idea of the removal-restoration method is to first use a high-order Earth gravity field model to calculate and remove the model part containing medium and long wave information, and then perform Stokes integration calculation on the remaining residual gravity anomaly that reflects local short wave details. This frequency division processing strategy effectively combines the overall advantages of the global gravity field model and the ability of local gravity data to characterize details, thereby obtaining geoid results that have both wide-area consistency and local high resolution.
[0028] 3. In the specific steps, from the coordinate reduction of the original observation data to the theoretical calculation of the disturbance potential and gravity anomaly, and then to the integration processing and final synthesis of the residual signal, the input, output and processing methods of each step are clearly defined, which enhances the repeatability of the method and the convenience of engineering implementation, and provides clear technical guidance for practical applications. Attached Figure Description
[0029] To more clearly illustrate the technical solutions in the embodiments of the present invention or the prior art, the drawings used in the description of the embodiments or the prior art will be briefly introduced below. Obviously, the drawings described below are only embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on the provided drawings without creative effort.
[0030] Figure 1 A flowchart of a geoid calculation method provided in an embodiment of the present invention; Figure 2 This is a spatial distribution map of the original observed gravity anomaly and the model gravity anomaly provided in an embodiment of the present invention; Figure 3 This is a Stokes integral residual geoid map based on residual gravity anomaly provided in an embodiment of the present invention. Figure 4 This is a final result map of the synthetic gravity geoid height provided in an embodiment of the present invention; Figure 5 A framework diagram of a geoid calculation system provided in an embodiment of the present invention; Figure 6 This is a schematic diagram of the electronic device structure provided in an embodiment of the present invention. Detailed Implementation
[0031] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.
[0032] The homomorphic comparison method provided in this application can be applied to a homomorphic comparison server. The homomorphic comparison server can be hardware or software. When the homomorphic comparison server is hardware, it can be implemented as a distributed server cluster providing homomorphic comparison services, or it can be implemented as a single server. When the homomorphic comparison server is software, it can be installed on the servers listed above. It can be implemented as multiple software programs or software modules, or it can be implemented as a single software program or software module; no specific limitations are made here.
[0033] Example 1; like Figure 1 As shown, this embodiment provides a method for calculating the geoid, including the following steps: S1. Convert the geographic coordinates of the observation point to geocentric coordinates to obtain the geocentric sphere coordinate parameters; S2. Based on the geocentric coordinate parameters, calculate the model perturbation potential of the observation point using the Earth's gravity field model; S3. Based on the model perturbation position, calculate the model gravity anomaly at the observation point; S4. Subtract the measured gravity anomaly from the model gravity anomaly to obtain the residual gravity anomaly at the observation point; S5. Substitute the residual gravity anomaly into the Stokes integral model to calculate the residual geoid height; S6. Based on the Earth's gravity field model, the height of the model's geoid is calculated; S7. Sum the residual geoid height with the model geoid height to obtain the geoid height.
[0034] This method addresses the challenge of directly conducting geometric leveling measurements in marine areas. Through a systematic removal-restore process, it effectively determines the marine geoid by combining a global gravity field model with local gravity observation data. This not only improves the overall accuracy and local detail resolution of the geoid model, but also enhances the operability and reliability of engineering implementation through standardized calculation steps.
[0035] The following provides a further explanation of each step and related features in the above method; In this embodiment S1, the geographic coordinates of the observation point are converted into geocentric coordinates to obtain geocentric sphere coordinate parameters; including: Convert geographic coordinates (longitude L, latitude B, elevation H) to geocentric coordinates (geocentric distance r, geocentric latitude). and longitude :
[0036] Where N is the radius of curvature of the ramusoidal circle, and e is the first eccentricity of the reference ellipsoid.
[0037] This step converts the geographic coordinates of the observation points into geocentric coordinates, providing a unified coordinate reference for subsequent physical calculations based on spherical harmonics. In practice, to ensure the accuracy of the conversion and the efficiency of large-scale data processing, the radius of curvature N of the primordial circle needs to be accurately calculated. Parallel computing processes and coordinate validity verification mechanisms can be adopted to ensure data quality and control error propagation.
[0038] In this embodiment S2, based on the geocentric sphere coordinate parameters, the model perturbation potential of the observation point is calculated using the Earth's gravity field model; The Earth's gravitational field model is represented as follows:
[0039] in, For model perturbation bits, For the remaining latitude, GM is the gravitational constant. With reference to the semi-major axis of the ellipsoid, and Here are the fully normalized potential coefficients of the Earth's gravity field model, where n is the order and m is the degree. To fully normalize the association of Legendre functions, This represents the highest order of the model expansion.
[0040] Using a global gravity field model, such as the fully normalized potential coefficients of EGM2008, the model perturbation potential T at the observation point is calculated. For higher-order associated Legendre functions, a stable standard recursive algorithm is usually used for calculation. Furthermore, by fusing regional refined potential coefficients or weighting multiple global models, the medium- and long-wave gravity field background of a specific sea area can be optimized, thereby enhancing the adaptability of this method to different marine geological structures.
[0041] In this embodiment, S3, the model gravity anomaly at the observation point is calculated based on the model perturbation position. Among them, the model gravity anomaly Represented as:
[0042] Furthermore, , For model perturbation bits, The distance from the center of the Earth is the coordinate system of the Earth's center.
[0043] This step converts the model perturbation T into the corresponding model gravity anomaly. In actual programming implementation, it can be directly calculated using a closed analytical formula derived from the potential coefficients and which is of the same origin as the perturbation potential expression in S2, thereby avoiding the additional errors that may be introduced by numerical differentiation. The gravity anomaly of this model characterizes the long-wave variation features of the gravity field on a global scale.
[0044] In this embodiment, S4, the difference between the measured gravity anomaly and the model gravity anomaly is calculated to obtain the residual gravity anomaly at the observation point; Among them, residual gravity anomaly value Represented as:
[0045] in, This represents the measured gravity anomaly value. These are the gravity anomalies in the model.
[0046] Furthermore, the measured gravity anomaly Δg at this observation point is derived from shipboard gravity measurement data, airborne gravity measurement data, or gravity anomaly data obtained by inversion from satellite altimetry data.
[0047] This step is central to the removal process, achieved by analyzing measured gravity anomalies. Subtracting the model gravity anomaly The residual gravity anomaly was obtained. It can effectively separate the shortwave components mainly caused by local topography and geological structure from the observation signal, so that the subsequent Stokes integration mainly processes local high-frequency components, which not only improves the integration efficiency, but also reduces the influence of far-field data, and enhances the stability and resolution of local inversion.
[0048] like Figure 2 As shown, the horizontal and vertical axes represent longitude and latitude, respectively, defining the planar location of the target sea area. The figure uses contour lines and color scales to display the spatial distribution of two types of gravity anomalies side by side: one is the original observed gravity anomaly obtained from actual measurements. Another type is the model gravity anomaly calculated based on the Earth's gravity field model. The visual comparison of the two directly presents the overall spatial agreement and differences between the observed values and the background field of the global model, providing an intuitive basis for subsequent calculation of residual anomalies.
[0049] In this embodiment, S5, the residual gravity anomaly is substituted into the Stokes integral model to calculate the residual geoid height. The Stokes integral model is expressed as follows:
[0050] in, Let be the residual geoid height, R be the Earth's mean radius, γ be the normal gravity value, and σ be the integration region. For Stokes functions, To calculate the spherical angular distance between the point and the flow surface element, Let be a infinitesimal element of the spherical surface area.
[0051] Furthermore, the Stokes function mentioned above Represented as:
[0052] in, It is an nth-order Legendre polynomial. This is the cutoff coefficient.
[0053] Substituting the residual gravity anomaly into the Stokes integral formula, the residual geoid height is calculated using surface integral. The output of this step directly characterizes the fine geoid undulations determined by the local gravity shortwave signal, and its spatial distribution is as follows: Figure 3 As shown, the method's ability to recover local details is intuitively demonstrated.
[0054] In this embodiment, S6, the height of the geoid is calculated based on the Earth's gravity field model. Among them, the model geoid height Represented as:
[0055] Furthermore, , This is the perturbation bit for the model.
[0056] Model geoid height It represents the long-wavelength portion of the smooth geoid as determined by the global gravity field model. The step-by-step calculation strategy decouples the contributions of long waves from those of medium and short waves, which is beneficial for separate analysis and quality control, and also provides flexibility for users with different accuracy requirements.
[0057] In this embodiment, S7, the residual geoid height is summed with the model geoid height to obtain the geoid height. The geoid height N is represented as:
[0058] in, The model's geoid height, This represents the residual geoid height.
[0059] This is the final restoration step, where the model geoid height is added to the residual geoid height to obtain the geoid height N, which integrates global trends and local details; the spatial distribution of the synthesized result is as follows. Figure 4 As shown, it retains all local details while overlaying smooth macroscopic trends, forming a more realistic and complete geoid shape.
[0060] The geoid calculation method in this embodiment overcomes the limitations of single methods in terms of spectral coverage. The final result integrates the overall consistency of the global model with the high-resolution details of local integrals, and can directly serve the establishment of marine mapping vertical benchmarks, satellite altimetry data correction, and marine geophysical research.
[0061] Example 2; like Figure 5 As shown, this embodiment provides a geoid calculation system, including: Coordinate transformation module: used to convert the geographic coordinates of the observation point to geocentric coordinates to obtain geocentric sphere coordinate parameters; Model perturbation potential calculation module: used to calculate the model perturbation potential of the observation point based on the geocentric sphere coordinate parameters and the Earth's gravity field model; Model gravity anomaly calculation module: used to calculate the model gravity anomaly at the observation point based on the model perturbation position; Residual gravity anomaly calculation module: used to subtract the measured gravity anomaly from the model gravity anomaly to obtain the residual gravity anomaly at the observation point; Stokes integral calculation module: used to substitute the residual gravity anomaly into the Stokes integral model to calculate the residual geoid height; Model geoid calculation module: used to calculate the height of the model geoid based on the Earth's gravity field model; Geoid height synthesis module: used to sum the residual geoid height with the model geoid height to obtain the geoid height.
[0062] Through the collaborative work of its various modules, the system constructs a complete and logically rigorous calculation process, realizing the systematic solution of the geoid based on gravity anomalies and the Earth's gravity field model. It is particularly suitable for establishing vertical benchmarks in areas such as oceans where geometric leveling is difficult to implement, and has good engineering practicality and computational reliability.
[0063] Example 3; like Figure 6 As shown, this embodiment provides an electronic device, including a memory, a processor, and a computer program stored in the memory and executable on the processor. When the processor executes the computer program, it implements the above-described geoid calculation method.
[0064] Example 4; This embodiment provides a computer-readable storage medium storing a computer program that, when executed by a processor, implements the above-described geoid calculation method.
[0065] In the embodiments provided in this application, it should be understood that the disclosed methods, systems, devices, and media can be implemented in other ways. The embodiments of methods, systems, devices, and media described above are merely illustrative. For example, the division of modules or units is only a logical functional division, and there may be other division methods in actual implementation. Each functional unit can be integrated into one processing unit, or each unit can exist physically separately, or two or more units can be integrated into one unit.
[0066] The units and algorithm steps of the various examples described in conjunction with the embodiments disclosed herein can be implemented in electronic hardware, or a combination of computer software and electronic hardware. Whether these functions are implemented in hardware or software depends on the specific application and design constraints of the technical solution. Those skilled in the art can use different methods to implement the described functions for each specific application, but such implementation should not be considered beyond the scope of this application.
[0067] Computer programs include computer program code, which can be in the form of source code, object code, executable files, or certain intermediate forms. Computer-readable media can include: any entity or device capable of carrying computer program code, recording media, USB flash drives, portable hard drives, magnetic disks, optical disks, computer memory, read-only memory, random access memory, electrical carrier signals, telecommunication signals, and software distribution media, etc.
[0068] The above embodiments are only used to illustrate the technical solutions of this application, and are not intended to limit them. Although this application has been described in detail with reference to the foregoing embodiments, those skilled in the art should understand that modifications can still be made to the technical solutions described in the foregoing embodiments, or equivalent substitutions can be made to some of the technical features. Such modifications or substitutions do not cause the essence of the corresponding technical solutions to deviate from the spirit and scope of the technical solutions of the embodiments of this application, and should all be included within the protection scope of this application.
Claims
1. A method for calculating the geoid, characterized in that, Includes the following steps: S1. Convert the geographic coordinates of the observation point to geocentric coordinates to obtain the geocentric sphere coordinate parameters; S2. Based on the geocentric coordinate parameters, calculate the model perturbation potential of the observation point using the Earth's gravity field model; S3. Based on the model perturbation position, calculate the model gravity anomaly at the observation point; S4. Subtract the measured gravity anomaly from the model gravity anomaly to obtain the residual gravity anomaly at the observation point; S5. Substitute the residual gravity anomaly into the Stokes integral model to calculate the residual geoid height; S6. Based on the Earth's gravity field model, the height of the model's geoid is calculated; S7. Sum the residual geoid height with the model geoid height to obtain the geoid height.
2. The geoid calculation method according to claim 1, characterized in that, S1 includes: Convert geographic coordinates (longitude L, latitude B, elevation H) to geocentric coordinates (geocentric distance r, geocentric latitude). and longitude : Where N is the radius of curvature of the ramusoidal circle, and e is the first eccentricity of the reference ellipsoid.
3. The geoid calculation method according to claim 1, characterized in that, In S2, the Earth's gravity field model is expressed as: in, For model perturbation bits, For the remaining latitude, GM is the gravitational constant. With reference to the semi-major axis of the ellipsoid, and Here are the fully normalized potential coefficients of the Earth's gravity field model, where n is the order and m is the degree. To fully normalize the association of Legendre functions, This represents the highest order of the model expansion.
4. The geoid calculation method according to claim 1, characterized in that, In S3, the model experiences gravity anomalies. Represented as: in, For model perturbation bits, The distance from the center of the Earth is the coordinate system of the Earth's center.
5. The method for calculating the geoid according to claim 1, characterized in that, In S5, the Stokes integral model is expressed as: in, Let be the residual geoid height, R be the Earth's mean radius, γ be the normal gravity value, and σ be the integration region. For Stokes functions, To calculate the spherical angular distance between the point and the flow surface element, Let be a infinitesimal element of the spherical surface area.
6. The geoid calculation method according to claim 1, characterized in that, The Stokes function Represented as: in, It is an nth-order Legendre polynomial. This is the cutoff coefficient.
7. The geoid calculation method according to claim 1, characterized in that, In S6, the model geoid height Represented as: in, This is the perturbation bit for the model.
8. A geoid calculation system, characterized in that, include: Coordinate transformation module: used to convert the geographic coordinates of the observation point to geocentric coordinates to obtain geocentric sphere coordinate parameters; Model perturbation potential calculation module: used to calculate the model perturbation potential of the observation point based on the geocentric sphere coordinate parameters and the Earth's gravity field model; Model gravity anomaly calculation module: used to calculate the model gravity anomaly at the observation point based on the model perturbation position; Residual gravity anomaly calculation module: used to subtract the measured gravity anomaly from the model gravity anomaly to obtain the residual gravity anomaly at the observation point; Stokes integral calculation module: used to substitute the residual gravity anomaly into the Stokes integral model to calculate the residual geoid height; Model geoid calculation module: used to calculate the height of the model geoid based on the Earth's gravity field model; Geoid height synthesis module: used to sum the residual geoid height with the model geoid height to obtain the geoid height.
9. An electronic device comprising a memory, a processor, and a computer program stored in the memory and executable on the processor, characterized in that, When the processor executes the computer program, it implements the geoid calculation method as described in any one of claims 1 to 7.
10. A computer-readable storage medium storing a computer program, characterized in that, When the computer program is executed by the processor, it implements the geoid calculation method as described in any one of claims 1 to 7.