Accurate computation of gravitational field of a tesseroid

We developed an accurate method to compute the gravitational field of a tesseroid. The method numerically integrates a surface integral representation of the gravitational potential of the tesseroid by conditionally splitting its line integration intervals and by using the double exponential quadrature rule. Then, it evaluates the gravitational acceleration vector and the gravity gradient tensor by numerically differentiating the numerically integrated potential. The numerical differentiation is conducted by appropriately switching the central and the single-sided second-order difference formulas with a suitable choice of the test argument displacement. If necessary, the new method is extended to the case of a general tesseroid with the variable density profile, the variable surface height functions, and/or the variable intervals in longitude or in latitude. The new method is capable of computing the gravitational field of the tesseroid independently on the location of the evaluation point, namely whether outside, near the surface of, on the surface of, or inside the tesseroid. The achievable precision is 14–15 digits for the potential, 9–11 digits for the acceleration vector, and 6–8 digits for the gradient tensor in the double precision environment. The correct digits are roughly doubled if employing the quadruple precision computation. The new method provides a reliable procedure to compute the topographic gravitational field, especially that near, on, and below the surface. Also, it could potentially serve as a sure reference to complement and elaborate the existing approaches using the Gauss–Legendre quadrature or other standard methods of numerical integration.     

Optimized formulas for the gravitational field of a tesseroid

Various tasks in geodesy, geophysics, and related geosciences require precise information on the impact of mass distributions on gravity field-related quantities, such as the gravitational potential and its partial derivatives. Using forward modeling based on Newton’s integral, mass distributions are generally decomposed into regular elementary bodies. In classical approaches, prisms or point mass approximations are mostly utilized. Considering the effect of the sphericity of the Earth, alternative mass modeling methods based on tesseroid bodies (spherical prisms) should be taken into account, particularly in regional and global applications. Expressions for the gravitational field of a point mass are relatively simple when formulated in Cartesian coordinates. In the case of integrating over a tesseroid volume bounded by geocentric spherical coordinates, it will be shown that it is also beneficial to represent the integral kernel in terms of Cartesian coordinates. This considerably simplifies the determination of the tesseroid’s potential derivatives in comparison with previously published methodologies that make use of integral kernels expressed in spherical coordinates. Based on this idea, optimized formulas for the gravitational potential of a homogeneous tesseroid and its derivatives up to second-order are elaborated in this paper. These new formulas do not suffer from the polar singularity of the spherical coordinate system and can, therefore, be evaluated for any position on the globe. Since integrals over tesseroid volumes cannot be solved analytically, the numerical evaluation is achieved by means of expanding the integral kernel in a Taylor series with fourth-order error in the spatial coordinates of the integration point. As the structure of the Cartesian integral kernel is substantially simplified, Taylor coefficients can be represented in a compact and computationally attractive form. Thus, the use of the optimized tesseroid formulas particularly benefits from a significant decrease in computation time by about 45 % compared to previously used algorithms. In order to show the computational efficiency and to validate the mathematical derivations, the new tesseroid formulas are applied to two realistic numerical experiments and are compared to previously published tesseroid methods and the conventional prism approach.     

A study on the reference frame consistency in recent Earth gravitational models

All gravity field functionals obtained from an Earth gravitational model (EGM) depend on the underlying terrestrial reference frame (TRF), with respect to which the EGM’s spherical harmonic coefficients refer to. In order to maintain a coherent framework for the comparison of current and future EGMs, it is thus important to investigate the consistency of their inherent TRFs, especially when their use is intended for high precision studies. Following the methodology described in an earlier paper by Kleusberg (1980), the similarity transformation parameters between the associated reference frames for several EGMs (including the most recent CHAMP/GRACE models at the time of writing this paper) are estimated in the present study. Specifically, the differences between the spherical harmonic coefficients for various pairs of EGMs are parameterized through a 3D-similarity spatial transformation model that relates their underlying TRFs. From the least-squares adjustment of such a parametric model, the origin, orientation and scale stability between the EGMs’ reference frames can be identified by estimating their corresponding translation, rotation and scale factor parameters. Various aspects of the estimation procedure and its results are highlighted in the paper, including data weighting schemes, the sensitivity of the results with respect to the selected harmonic spectral band, the correlation structure and precision level of the estimated transformation parameters, the effect of the estimated differences of the EGMs’ reference frames on their height anomaly signal, and the overall feasibility of Kleusberg’s formulae for the assessment of TRF inconsistencies among global geopotential models.     

Evaluation of gravitational curvatures of a tesseroid in spherical integral kernels

Proper understanding of how the Earth’s mass distributions and redistributions influence the Earth’s gravity field-related functionals is crucial for numerous applications in geodesy, geophysics and related geosciences. Calculations of the gravitational curvatures (GC) have been proposed in geodesy in recent years. In view of future satellite missions, the sixth-order developments of the gradients are becoming requisite. In this paper, a set of 3D integral GC formulas of a tesseroid mass body have been provided by spherical integral kernels in the spatial domain. Based on the Taylor series expansion approach, the numerical expressions of the 3D GC formulas are provided up to sixth order. Moreover, numerical experiments demonstrate the correctness of the 3D Taylor series approach for the GC formulas with order as high as sixth order. Analogous to other gravitational effects (e.g., gravitational potential, gravity vector, gravity gradient tensor), numerically it is found that there exist the very-near-area problem and polar singularity problem in the GC east–east–radial, north–north–radial and radial–radial–radial components in spatial domain, and compared to the other gravitational effects, the relative approximation errors of the GC components are larger due to not only the influence of the geocentric distance but also the influence of the latitude. This study shows that the magnitude of each term for the nonzero GC functionals by a grid resolution 15\(^{{\prime } }\,\times \) 15\(^{{\prime }}\) at GOCE satellite height can reach of about 10\(^{-16}\) m\(^{-1}\) s\(^{2}\) for zero order, 10\(^{-24 }\) or 10\(^{-23}\) m\(^{-1}\) s\(^{2}\) for second order, 10\(^{-29}\) m\(^{-1}\) s\(^{2}\) for fourth order and 10\(^{-35}\) or 10\(^{-34}\) m\(^{-1}\) s\(^{2}\) for sixth order, respectively.     

Harmonic analysis of the Earth's gravitational field by means of semi-continuous ephemerides of a low Earth orbiting GPS-tracked satellite. Case study: CHAMP

An algorithm for the determination of the spherical harmonic coefficients of the terrestrial gravitational field representation from the analysis of a kinematic orbit solution of a low earth orbiting GPS-tracked satellite is presented and examined. A gain in accuracy is expected since the kinematic orbit of a LEO satellite can nowadays be determined with very high precision, in the range of a few centimeters. In particular, advantage is taken of Newton's Law of Motion, which balances the acceleration vector with respect to an inertial frame of reference (IRF) and the gradient of the gravitational potential. By means of triple differences, and in particular higher-order differences (seven-point scheme, nine-point scheme), based upon Newton's interpolation formula, the local acceleration vector is estimated from relative GPS position time series. The gradient of the gravitational potential is conventionally given in a body-fixed frame of reference (BRF) where it is nearly time independent or stationary. Accordingly, the gradient of the gravitational potential has to be transformed from spherical BRF to Cartesian IRF. Such a transformation is possible by differentiating the gravitational potential, given as a spherical harmonics series expansion, with respect to Cartesian coordinates by means of the chain rule, and expressing zero- and first-order Ferrer's associated Legendre functions in terms of Cartesian coordinates. Subsequently, the BRF Cartesian coordinates are transformed into IRF Cartesian coordinates by means of the polar motion matrix, the precession–nutation matrices and the Greenwich sidereal time angle (GAST). In such a way a spherical harmonic representation of the terrestrial gravitational field intensity with respect to an IRF is achieved. Numerical tests of a resulting Gauss–Markov model document not only the quality and the high resolution of such a space gravity spectroscopy, but also the problems resulting from noise amplification in the acceleration determination process.     

Effects of topographic–isostatic masses on gravitational functionals at the Earth’s surface and at airborne and satellite altitudes

Topographic–isostatic masses represent an important source of gravity field information, especially in the high-frequency band, even if the detailed mass-density distribution inside the topographic masses is unknown. If this information is used within a remove-restore procedure, then the instability problems in downward continuation of gravity observations from aircraft or satellite altitudes can be reduced. In this article, integral formulae are derived for determination of gravitational effects of topographic–isostatic masses on the first- and second-order derivatives of the gravitational potential for three topographic–isostatic models. The application of these formulas is useful for airborne gravimetry/gradiometry and satellite gravity gradiometry. The formulas are presented in spherical approximation by separating the 3D integration in an analytical integration in the radial direction and 2D integration over the mean sphere. Therefore, spherical volume elements can be considered as being approximated by mass-lines located at the centre of the discretization compartments (the mass of the tesseroid is condensed mathematically along its vertical axis). The errors of this approximation are investigated for the second-order derivatives of the topographic–isostatic gravitational potential in the vicinity of the Earth’s surface. The formulas are then applied to various scenarios of airborne gravimetry/gradiometry and satellite gradiometry. The components of the gravitational vector at aircraft altitudes of 4 and 10 km have been determined, as well as the gravitational tensor components at a satellite altitude of 250 km envisaged for the forthcoming GOCE (gravity field and steady-state ocean-circulation explorer) mission. The numerical computations are based on digital elevation models with a 5-arc-minute resolution for satellite gravity gradiometry and 1-arc-minute resolution for airborne gravity/gradiometry.     

GOCE gravitational gradients along the orbit

GOCE is ESA’s gravity field mission and the first satellite ever that measures gravitational gradients in space, that is, the second spatial derivatives of the Earth’s gravitational potential. The goal is to determine the Earth’s mean gravitational field with unprecedented accuracy at spatial resolutions down to 100 km. GOCE carries a gravity gradiometer that allows deriving the gravitational gradients with very high precision to achieve this goal. There are two types of GOCE Level 2 gravitational gradients (GGs) along the orbit: the gravitational gradients in the gradiometer reference frame (GRF) and the gravitational gradients in the local north oriented frame (LNOF) derived from the GGs in the GRF by point-wise rotation. Because the V _XX , V _YY , V _ZZ and V _XZ are much more accurate than V _XY and V _YZ , and because the error of the accurate GGs increases for low frequencies, the rotation requires that part of the measured GG signal is replaced by model signal. However, the actual quality of the gradients in GRF and LNOF needs to be assessed. We analysed the outliers in the GGs, validated the GGs in the GRF using independent gravity field information and compared their assessed error with the requirements. In addition, we compared the GGs in the LNOF with state-of-the-art global gravity field models and determined the model contribution to the rotated GGs. We found that the percentage of detected outliers is below 0.1% for all GGs, and external gravity data confirm that the GG scale factors do not differ from one down to the 10⁻³ level. Furthermore, we found that the error of V _XX and V _YY is approximately at the level of the requirement on the gravitational gradient trace, whereas the V _ZZ error is a factor of 2–3 above the requirement for higher frequencies. We show that the model contribution in the rotated GGs is 2–35% dependent on the gravitational gradient. Finally, we found that GOCE gravitational gradients and gradients derived from EIGEN-5C and EGM2008 are consistent over the oceans, but that over the continents the consistency may be less, especially in areas with poor terrestrial gravity data. All in all, our analyses show that the quality of the GOCE gravitational gradients is good and that with this type of data valuable new gravity field information is obtained.     

A comparison of the tesseroid,prism and point-mass approaches for mass reductions in gravity field modelling

The calculation of topographic (and iso- static) reductions is one of the most time-consuming operations in gravity field modelling. For this calculation, the topographic surface of the Earth is often divided with respect to geographical or map-grid lines, and the topographic heights are averaged over the respective grid elements. The bodies bounded by surfaces of constant (ellipsoidal) heights and geographical grid lines are denoted as tesseroids. Usually these ellipsoidal (or spherical) tesseroids are replaced by “equivalent” vertical rectangular prisms of the same mass. This approximation is motivated by the fact that the volume integrals for the calculation of the potential and its derivatives can be exactly solved for rectangular prisms, but not for the tesseroids. In this paper, an approximate solution of the spherical tesseroid integrals is provided based on series expansions including third-order terms. By choosing the geometrical centre of the tesseroid as the Taylor expansion point, the number of non-vanishing series terms can be greatly reduced. The zero-order term is equivalent to the point-mass formula. Test computations show the high numerical efficiency of the tesseroid method versus the prism approach, both regarding computation time and accuracy. Since the approximation errors due to the truncation of the Taylor series decrease very quickly with increasing distance of the tesseroid from the computation point, only the elements in the direct vicinity of the computation point have to be separately evaluated, e.g. by the prism formulas. The results are also compared with the point-mass formula. Further potential refinements of the tesseroid approach, such as considering ellipsoidal tesseroids, are indicated.     

Four dimensional studies in Earth space

Summary A system of reference which is directly related to observations, is proposed for four dimensional studies in Earth space. The requisite data is used to define both global control network and also polar wandering. The determination of variations of the Earth’s gravitational field with time also forms part of such a system. Techniques are outlined for the unique definition of the motion of the geocenter, and the changes in the location of the axis of rotation of an instantaneous Earth model, in relation to values at some epoch of reference. The instantaneous system referred to is directly related to a fundamental equation in geodynamics. The reference system defined would provide an unambiguous frame for long period studies in Earth space, provided the scale of the space were specified. Presented at the I.A.G./A.G.U. Symposium on Earth Gravity Models related problems, St Louis, Missouri, U.S.A., 16–18 August 1972.     

弹性地球自转动力学理论的研究进展

地球自转的精确测定是高精度大地参考系建立的理论基础,也是天体测量学、大地测量学和地球物理学共同关注的研究领域。研究了弹性地球自转动力学的基本理论,主要内容包括:修正了固体潮和地球自转对地球惯量张量元影响的表达式;首次给出了高阶岁差章动力矩对地球自转的影响;介绍了大气对地球自转影响的数值积分法,并给出了负荷理论方法;总结了大陆水分布对地球自转影响的数值原理;研究了海洋负荷潮汐对地球自转影响的理论和方法,在推导过程中抛弃了传统的直接数值积分法,直接以负荷引力位为基础,并将海潮潮高引入表达式中,推导出了海潮对地球自转影响的有关公式。本文给出的理论公式可为动力大地测量学和天文地球动力学的研究提供理论参考和依据。     

Gravitational gradients by tensor analysis with application to spherical coordinates

This contribution deals with the derivation of explicit expressions of the gradients of first, second and third order of the gravitational potential. This is accomplished in the framework of tensor analysis which naturally allows to apply general formulae to the specific coordinate systems in use in geodesy. In particular it is recalled here that when the potential field is expressed in general coordinates on a 3D manifold, the gradient operation leads to the definition of the covariant derivative and that the covariant derivative of a tensor can be obtained by application of a simple rule. When applied to the gravitational potential or to any of its gradients, the rule straightforwardly provides the expressions of the higher-order gradients. It is also shown that the tensor approach offers a clear distinction between natural and physical components of the gradients. Two fundamental reference systems—a global, bodycentric system and a local, topocentric system, both body-fixed—are introduced and transformation rules are derived to convert quantities between the two systems. The results include explicit expressions for the gradients of the first three orders in both reference systems.     

The spherical horizontal and spherical vertical boundary value problem – vertical deflections and geoidal undulations – the completed Meissl diagram

 In a comparison of the solution of the spherical horizontal and vertical boundary value problems of physical geodesy it is aimed to construct downward continuation operators for vertical deflections (surface gradient of the incremental gravitational potential) and for gravity disturbances (vertical derivative of the incremental gravitational potential) from points on the Earth's topographic surface or of the three-dimensional (3-D) Euclidean space nearby down to the international reference sphere (IRS). First the horizontal and vertical components of the gravity vector, namely spherical vertical deflections and spherical gravity disturbances, are set up. Second, the horizontal and vertical boundary value problem in spherical gravity and geometry space is considered. The incremental gravity vector is represented in terms of vector spherical harmonics. The solution of horizontal spherical boundary problem in terms of the horizontal vector-valued Green function converts vertical deflections given on the IRS to the incremental gravitational potential external in the 3-D Euclidean space. The horizontal Green functions specialized to evaluation and source points on the IRS coincide with the Stokes kernel for vertical deflections. Third, the vertical spherical boundary value problem is solved in terms of the vertical scalar-valued Green function. Fourth, the operators for upward continuation of vertical deflections given on the IRS to vertical deflections in its external 3-D Euclidean space are constructed. Fifth, the operators for upward continuation of incremental gravity given on the IRS to incremental gravity to the external 3-D Euclidean space are generated. Finally, Meissl-type diagrams for upward continuation and regularized downward continuation of horizontal and vertical gravity data, namely vertical deflection and incremental gravity, are produced. Received: 10 May 2000 / Accepted: 26 February 2001     

Construction of spherical harmonic series for the potential derivatives of arbitrary orders in the geocentric Earth-fixed reference frame

The derivatives of the Earth gravitational potential are considered in the global Cartesian Earth-fixed reference frame. Spherical harmonic series are constructed for the potential derivatives of the first and second orders on the basis of a general expression of Cunningham (Celest Mech 2:207–216, 1970) for arbitrary order derivatives of a spherical harmonic. A common structure of the series for the potential and its first- and second-order derivatives allows to develop a general procedure for constructing similar series for the potential derivatives of arbitrary orders. The coefficients of the derivatives are defined by means of recurrence relations in which a coefficient of a certain order derivative is a linear function of two coefficients of a preceding order derivative. The coefficients of the second-order derivatives are also presented as explicit functions of three coefficients of the potential. On the basis of the geopotential model EGM2008, the spherical harmonic coefficients are calculated for the first-, second-, and some third-order derivatives of the disturbing potential T, representing the full potential V, after eliminating from it the zero- and first-degree harmonics. The coefficients of two lowest degrees in the series for the derivatives of T are presented. The corresponding degree variances are estimated. The obtained results can be applied for solving various problems of satellite geodesy and celestial mechanics.     

On the spherical and spheroidal harmonic expansion of the gravitational potential of the topographic masses

This paper is devoted to the spherical and spheroidal harmonic expansion of the gravitational potential of the topographic masses in the most rigorous way. Such an expansion can be used to compute gravimetric topographic effects for geodetic and geophysical applications. It can also be used to augment a global gravity model to a much higher resolution of the gravitational potential of the topography. A formulation for a spherical harmonic expansion is developed without the spherical approximation. Then, formulas for the spheroidal harmonic expansion are derived. For the latter, Legendre’s functions of the first and second kinds with imaginary variable are expanded in Laurent series. They are then scaled into two real power series of the second eccentricity of the reference ellipsoid. Using these series, formulas for computing the spheroidal harmonic coefficients are reduced to surface harmonic analysis. Two numerical examples are presented. The first is a spherical harmonic expansion to degree and order 2700 by taking advantage of existing software. It demonstrates that rigorous spherical harmonic expansion is possible, but the computed potential on the geoid shows noticeable error pattern at Polar Regions due to the downward continuation from the bounding sphere to the geoid. The second numerical example is the spheroidal expansion to degree and order 180 for the exterior space. The power series of the second eccentricity of the reference ellipsoid is truncated at the eighth order leading to omission errors of 25 nm (RMS) for land areas, with extreme values around 0.5 mm to geoid height. The results show that the ellipsoidal correction is 1.65 m (RMS) over land areas, with maximum value of 13.19 m in the Andes. It shows also that the correction resembles the topography closely, implying that the ellipsoidal correction is rich in all frequencies of the gravity field and not only long wavelength as it is commonly assumed.     

GNSS科学发展与前景

GNSS早已成为多学科研究的强大工具。着重介绍了近几年来GNSS系统及其观测技术,特别是数据处理方法和技术的重大发展,阐述了数据格式标准化的进展、全球与区域参考框架现状、参考框架对地壳运动研究的影响和IGS的产品等。同时对利用GNSS观测研究成果,主要是地壳运动的研究成果也作了较为详细的介绍。     

A comparison of different mass elements for use in gravity gradiometry

Topographic and isostatic mass anomalies affect the external gravity field of the Earth. Therefore, these effects also exist in the gravity gradients observed, e.g., by the satellite gravity gradiometry mission GOCE (Gravity and Steady-State Ocean Circulation Experiment). The downward continuation of the gravitational signals is rather difficult because of the high-frequency behaviour of the combined topographic and isostatic effects. Thus, it is preferable to smooth the gravity field by some topographic-isostatic reduction. In this paper the focus is on the modelling of masses in the space domain, which can be subdivided into different mass elements and evaluated with analytical, semi-analytical and numerical methods. Five alternative mass elements are reviewed and discussed: the tesseroid, the point mass, the prism, the mass layer and the mass line. The formulae for the potential, the attraction components and the Marussi tensor of second-order potential derivatives are provided. The formulae for different mass elements and computation methods are checked by assuming a synthetic topography of constant height over a spherical cap and the position of the computation point on the polar axis. For this special situation an exact analytical solution for the tesseroid exists and a comparison between the analytical solution of a spherical cap and the modelling of different mass elements is possible. A comparison of the computation times shows that modelling by tesseroids with different methods produces the most accurate results in an acceptable computation time. As a numerical example, the Marussi tensor of the topographic effect is computed globally using tesseroids calculated by Gauss–Legendre cubature (3D) on the basis of a digital height model. The order of magnitude in the radial-radial component is about  ± 8 E.U. Electronic supplementary material  The online version of this article (doi:) contains supplementary material, which is available to authorized users.     

The reference frames of Mercury after the MESSENGER mission

We report on recent refinements and the current status for the rotational state models and the reference frames of the planet Mercury. We summarize the performed measurements of Mercury rotation based on terrestrial radar observations as well as data from the Mariner 10 and the MESSENGER missions. Further, we describe the different available definitions of reference systems for Mercury and obtain the corresponding reference frame using data provided by instruments on board MESSENGER. In particular, we discuss the dynamical frame, the principal-axes frame, the ellipsoid frame, as well as the cartographic frame. We also describe the reference frame adopted by the MESSENGER science team for the release of their cartographic products, and we provide expressions for transformations from this frame to the other reference frames.     

Integral transformations of gradiometric data onto a GRACE type of observable

Integral transformations of gravitational gradients onto a Gravity Recovery And Climate Experiment (GRACE) type of observable are derived in this article. The gravitational gradients represent components of the gravitational tensor in the local north-oriented frame. The GRACE type of observable corresponds to a difference between two gravitational vectors as projected onto the line of sight between the two GRACE satellites. In total, three integral transformations relating vertical–vertical, vertical–horizontal and horizontal–horizontal gravitational gradients with the GRACE type of observable are provided. Spectral and closed forms of corresponding isotropic kernels are derived for each transformation. Special cases show that the integral transformations are general and relate gravitational gradients to many other quantities of the gravitational field, such as the gravitational vector, and its radial and tangential components. Correctness of the mathematical derivations is validated in a closed-loop simulation using synthetic data.     

同时顾及章动和极移的地球自转方程

通过引进章动坐标系相对惯性参照系的转动角速度随时间的变化 ,导出了一个可同时解出章动和极移的地球自转方程 ,用这个方程可同时研究地球的强迫和自由转动。与现行研究地球自转的惯用方法相比 ,该方法综合性强 ,易于理解。     

On the spherical harmonics series expansion of the geopotential

Summary The spherical harmonics development of the gravitational potential at the earth's surface and on the geoid is discussed. First, the effect of Kelvin transformations is studied; secondly, numerical information as obtained from satellite and other data is investigated. Institut für astronomische und physikalische Geod?sie, Publ. No. 66, Techn. Hochschule, Munich.