椭球谐和球谐系数之间一个简单的转换关系

本文推导的椭球谐系数和球谐系数相互之间转换关系的核心思想是在ε~2量级下利用Legendre函数的正交性,从球谐系数求解的积分表示出发,将积分中的椭球坐标变量与球坐标变量相互转换,从而得出椭球谐系数与球谐系数之间的转换关系。本文导出的转换关系有以下优点:①对于第二类Legendre函数的计算采用Laurent级数表示,使计算第二类Legendre函数更为简单;②保留了ε~2量级下,导出的转换关系相比文献[2]的形式更简单,满足物理大地测量边值问题线性化的要求;③顾及了余纬和归化余纬的区别。     

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.     

Downward continuation of the free-air gravity anomalies to the ellipsoid using the gradient solution and terrain correction—An attempt of global numerical computations

The formulas for the determination of the coefficients of the spherical harmonic expansion of the disturbing potential of the earth are defined for data given on a sphere. In order to determine the spherical harmonic coefficients, the gravity anomalies have to be analytically downward continued from the earth's surface to a sphere—at least to the ellipsoid. The goal of this paper is to continue the gravity anomalies from the earth's surface downward to the ellipsoid using recent elevation models. The basic method for the downward continuation is the gradient solution (theg ₁ term). The terrain correction has also been computed because of the role it can play as a correction term when calculating harmonic coefficients from surface gravity data. Theg ₁ term and the terrain correction were expanded into the spherical harmonics up to180 ^th order. The corrections (theg ₁ term and the terrain correction) have the order of about 2% of theRMS value of degree variance of the disturbing potential per degree. The influences of theg ₁ term and the terrain correction on the geoid take the order of 1 meter (RMS value of corrections of the geoid undulation) and on the deflections of the vertical is of the order 0.1″ (RMS value of correction of the deflections of the vertical).     

Transformation of coordinates between two horizontal geodetic datums

The following topics are discussed in this paper: the geocentric coordinate system and its different realizations used in geodetic practice; the definition of a horizontal geodetic datum (reference ellipsoid) and its positioning and orientation with respect to the geocentric coordinate system; positions on a horizontal datum and errors inherent in the process of positioning; and distortions of geodetic networks referred to a horizontal datum. The problem of determining transformation parameters between a horizontal datum and the geocentric coordinate system from known positions is then analysed. It is often found necessary to transform positions from one horizontal datum to another. These transformations are normally accomplished through the geocentric coordinate system and they include the transformation parameters of the two datums as well as the representation of the respective network distortions. Problems encountered in putting these transformations together are pointed out.     

Space-wise approach to satellite gravity field determination in the presence of coloured noise

For many years, the gravity field of the Earth was only seen by satellite geodesy as the main factor affecting the orbit and consequently it was retrieved together with a number of other orbital perturbations. Since the advent of a new generation of accelerometers, non-gravitational perturbations can be separated from the gravity effects and a new era of gravity field estimates from space has been born. During preparatory data analysis for new missions performed by the geodetic community, three approaches have been proposed and numerically tested: the brute force method (direct approach), the semi-analytical (time-wise) method and the space-wise method. In particular, the time-wise method takes advantage of the incoming time flow of data and, after performing a Fourier transform of the observation equations, exploits the prevailing block diagonal structure of the normal equations to estimate the spherical harmonic coefficients of the gravity field. Complementary to this is the space-wise approach, which goes back to the traditional computation of the harmonic coefficients by an integration technique or by least-squares collocation. Some advantages and disadvantages are peculiar to both methods, particularly the space-wise approach, which has for a long time ignored the marked signature of the noise spectrum due to the specific measuring conditions of space-borne accelerometers. The application of a proper Wiener filter, exploiting the correlation along the orbit, embedded into an iterative scheme, seems to be the answer. The solution to this major problem of the space-wise approach is illustrated and simulation results are discussed.     

椭球界面下Neumann边值问题的积分解

应用文献 [1 ]推导出的球谐系数与椭球谐系数的转换关系 ,给出了椭球界面下Neumann边值问题的积分解     

W0: an estimate in the Finnish Height Datum N60, epoch 1993.4, from twenty-five GPS points of the Baltic Sea Level Project

Based upon a data set of 25 points of the Baltic Sea Level Project, second campaign 1993.4, which are close to mareographic stations, described by (1) GPS derived Cartesian coordinates in the World Geodetic Reference System 1984 and (2) orthometric heights in the Finnish Height Datum N60, epoch 1993.4, we have computed the primary geodetic parameter W ₀(1993.4) for the epoch 1993.4 according to the following model. The Cartesian coordinates of the GPS stations have been converted into spheroidal coordinates. The gravity potential as the additive decomposition of the gravitational potential and the centrifugal potential has been computed for any GPS station in spheroidal coordinates, namely for a global spheroidal model of the gravitational potential field. For a global set of spheroidal harmonic coefficients a transformation of spherical harmonic coefficients into spheroidal harmonic coefficients has been implemented and applied to the global spherical model OSU 91A up to degree/order 360/360. The gravity potential with respect to a global spheroidal model of degree/order 360/360 has been finally transformed by means of the orthometric heights of the GPS stations with respect to the Finnish Height Datum N60, epoch 1993.4, in terms of the spheroidal “free-air” potential reduction in order to produce the spheroidal W ₀(1993.4) value. As a mean of those 25 W ₀(1993.4) data as well as a root mean square error estimation we computed W ₀(1993.4)=(6 263 685.58 ± 0.36) kgal × m. Finally a comparison of different W ₀ data with respect to a spherical harmonic global model and spheroidal harmonic global model of Somigliana-Pizetti type (level ellipsoid as a reference, degree/order 2/0) according to The Geodesist's Handbook 1992 has been made. Received: 7 November 1996 / Accepted: 27 March 1997     

The gains of small circular,square and rectangular filters for surface waves on a sphere

Meissl has derived weighting functions for converting point gravity anomaly degree variances into mean anomaly variances over a circular cap on a sphere. If the cap is sufficiently small so that the cap on a sphere degenerates into a circle on a plane, the problem may be considered that of the gain of a circular filter for a surface wave whose wave number depends on the spherical harmonic degree. The Meissl weights then become replaced by diffraction integrals of optical physics. The expected gain for a square filter for waves coming from random directions is derived and shown to be close to the gain of a circular filter with the same area. The expected gains and cross-gains for rectangular filters are also derived. When weighted by an anomaly degree variance model, these gains and cross-gains can be used to determine rectangular anomaly variances and covariances for arbitrary bandwidths. Using the Tscherning-Rapp model, analytic gravity anomaly variances and covariances are calculated for 1°×1° blocks.     

Spherical harmonic modelling to ultra-high degree of Bouguer and isostatic anomalies

The availability of high-resolution global digital elevation data sets has raised a growing interest in the feasibility of obtaining their spherical harmonic representation at matching resolution, and from there in the modelling of induced gravity perturbations. We have therefore estimated spherical Bouguer and Airy isostatic anomalies whose spherical harmonic models are derived from the Earth’s topography harmonic expansion. These spherical anomalies differ from the classical planar ones and may be used in the context of new applications. We succeeded in meeting a number of challenges to build spherical harmonic models with no theoretical limitation on the resolution. A specific algorithm was developed to enable the computation of associated Legendre functions to any degree and order. It was successfully tested up to degree 32,400. All analyses and syntheses were performed, in 64 bits arithmetic and with semi-empirical control of the significant terms to prevent from calculus underflows and overflows, according to IEEE limitations, also in preserving the speed of a specific regular grid processing scheme. Finally, the continuation from the reference ellipsoid’s surface to the Earth’s surface was performed by high-order Taylor expansion with all grids of required partial derivatives being computed in parallel. The main application was the production of a 1′ × 1′ equiangular global Bouguer anomaly grid which was computed by spherical harmonic analysis of the Earth’s topography–bathymetry ETOPO1 data set up to degree and order 10,800, taking into account the precise boundaries and densities of major lakes and inner seas, with their own altitude, polar caps with bedrock information, and land areas below sea level. The harmonic coefficients for each entity were derived by analyzing the corresponding ETOPO1 part, and free surface data when required, at one arc minute resolution. The following approximations were made: the land, ocean and ice cap gravity spherical harmonic coefficients were computed up to the third degree of the altitude, and the harmonics of the other, smaller parts up to the second degree. Their sum constitutes what we call ETOPG1, the Earth’s TOPography derived Gravity model at 1′ resolution (half-wavelength). The EGM2008 gravity field model and ETOPG1 were then used to rigorously compute 1′ × 1′ point values of surface gravity anomalies and disturbances, respectively, worldwide, at the real Earth’s surface, i.e. at the lower limit of the atmosphere. The disturbance grid is the most interesting product of this study and can be used in various contexts. The surface gravity anomaly grid is an accurate product associated with EGM2008 and ETOPO1, but its gravity information contents are those of EGM2008. Our method was validated by comparison with a direct numerical integration approach applied to a test area in Morocco–South of Spain (Kuhn, private communication 2011) and the agreement was satisfactory. Finally isostatic corrections according to the Airy model, but in spherical geometry, with harmonic coefficients derived from the sets of the ETOPO1 different parts, were computed with a uniform depth of compensation of 30?km. The new world Bouguer and isostatic gravity maps and grids here produced will be made available through the Commission for the Geological Map of the World. Since gravity values are those of the EGM2008 model, geophysical interpretation from these products should not be done for spatial scales below 5 arc minutes (half-wavelength).     

Geoid determination using one-step integration

A residual (high-frequency) gravimetric geoid is usually computed from geographically limited ground, sea and/or airborne gravimetric data. The mathematical model for its determination from ground gravity is based on the transformation of observed discrete values of gravity into gravity potential related to either the international ellipsoid or the geoid. The two reference surfaces are used depending on height information that accompanies ground gravity data: traditionally orthometric heights determined by geodetic levelling were used while GPS positioning nowadays allows for estimation of geodetic (ellipsoidal) heights. This transformation is usually performed in two steps: (1) observed values of gravity are downward continued to the ellipsoid or the geoid, and (2) gravity at the ellipsoid or the geoid is transformed into the corresponding potential. Each of these two steps represents the solution of one geodetic boundary-value problem of potential theory, namely the first and second or third problem. Thus two different geodetic boundary-value problems must be formulated and solved, which requires numerical evaluation of two surface integrals. In this contribution, a mathematical model in the form of a single Fredholm integral equation of the first kind is presented and numerically investigated. This model combines the solution of the first and second/third boundary-value problems and transforms ground gravity disturbances or anomalies into the harmonically downward continued disturbing potential at the ellipsoid or the geoid directly. Numerical tests show that the new approach offers an efficient and stable solution for the determination of the residual geoid from ground gravity data.     

基于“浅层海水”质量法确定海洋内部层面重力值

经典物理大地测量学利用斯托克斯方法和莫洛金斯基方法解算大地测量边值问题并给出地球外部重力场表达,若忽略1~2 m量级的动力学海面地形,静止的平均海面可认为是大地水准面,后者是与平均海平面最为接近的重力等位面。经典理论无法求解海洋内部,即地球内部重力场问题,为解决这一局限,基于地表浅层法引入“浅层海水”的概念,“浅层海水”上下界面由平均海面高模型DTU21确定,利用牛顿积分和球谐展开算法确定了最优球谐分析迭代次数,分析了“浅层海水”厚度与积分区域半径大小的关系,确定了“浅层海水”厚度为100 m、500 m和1 000 m时的最优积分区域半径为1°,厚度4 000 m时为1.5°;评估了“浅层海水”质量法移去-恢复海洋表面重力值的精度,“浅层海水”厚度100 m、500 m、1 000 m和4 000 m的均方根误差分别为0.13 mGal、0.61 mGal、1.21 mGal和3.93 mGal,验证了该方法的可靠性。基于此理论,计算了不同厚度“浅层海水”下表面的层面重力值,得到了100 m、500 m、1 000 m和4 000 m深度处层面重力值与“浅层海水”上表面重力值差的均方根,分别为22.11 mGal、110.50 mGal、220.87 mGal和877.31 mGal。     

GRACE时变重力场各向异性高斯组合滤波方法

受测量误差等因素影响,直接使用GRACE时变重力场模型的地表质量变化反演结果呈现严重条带噪声,必须采用滤波消除。本文对不同滤波方法进行了试验分析,以信噪比最大为准则,确定了不同滤波方法的最优滤波参数,并在此基础上提出了一种各向异性组合滤波方法。该方法根据时变重力场模型球谐系数误差特性,结合各向异性高斯滤波和均方根滤波特点,对精度较高的低次项系数采用较大权重以保留更多有效信号,而对精度较差的高次项系数采用较小权重以压制噪声。不同于传统的两步法组合滤波,该方法仅需进行一步滤波处理。试验结果表明,本文提出的各向异性组合滤波方法计算步骤简单,能够有效消除条带噪声;与单一滤波和传统两步法组合滤波方法相比,提高了反演结果信噪比,保留了更多真实信号。     

第二大地边值问题引论EI北大核心CSCD

在空间大地测量时代,GNSS可以测定地面点的大地高,使重力扰动变成了直接观测量,以重力扰动为边界条件的第二边值问题在大地测量中得以实用化。它的解与GNSS组合正在成为一种颇有应用前景的海拔高测量方法。本文原理性地讨论了有两种不同边界面的球近似第二大地边值问题。第一种以地形面为边界面,给出了高程异常与地面垂线偏差的解析延拓解;第二种以参考椭球面为边界面,将其外部地形质量按照Helmert第二压缩法移至参考椭球面,然后将Hotine函数与从地球表面延拓至边界面的Helmert重力扰动进行卷积,并顾及地形间接影响,最后得到大地水准面高、椭球面垂线偏差、高程异常与地面垂线偏差的Helmert解。在讨论部分,进行了第二与第三大地边值问题的比较,提出了现有重力点高程从正高或正常高到大地高的改化方法,并展望了它的应用前景。     

Transforming ellipsoidal heights and geoid undulations between different geodetic reference frames

Transforming height information that refers to an ellipsoidal Earth reference model, such as the geometric heights determined from GPS measurements or the geoid undulations obtained by a gravimetric geoid solution, from one geodetic reference frame (GRF) to another is an important task whose proper implementation is crucial for many geodetic, surveying and mapping applications. This paper presents the required methodology to deal with the above problem when we are given the Helmert transformation parameters that link the underlying Cartesian coordinate systems to which an Earth reference ellipsoid is attached. The main emphasis is on the effect of GRF spatial scale differences in coordinate transformations involving reference ellipsoids, for the particular case of heights. Since every three-dimensional Cartesian coordinate system ‘gauges’ an attached ellipsoid according to its own accessible scale, there will exist a supplementary contribution from the scale variation between the involved GRFs on the relative size of their attached reference ellipsoids. Neglecting such a scale-induced indirect effect corrupts the values for the curvilinear geodetic coordinates obtained from a similarity transformation model, and meter-level apparent offsets can be introduced in the transformed heights. The paper explains the above issues in detail and presents the necessary mathematical framework for their treatment. An erratum to this article can be found at     

Transformation between surface spherical harmonic expansion of arbitrary high degree and order and double Fourier series on sphere

In order to accelerate the spherical harmonic synthesis and/or analysis of arbitrary function on the unit sphere, we developed a pair of procedures to transform between a truncated spherical harmonic expansion and the corresponding two-dimensional Fourier series. First, we obtained an analytic expression of the sine/cosine series coefficient of the \(4 \pi \) fully normalized associated Legendre function in terms of the rectangle values of the Wigner d function. Then, we elaborated the existing method to transform the coefficients of the surface spherical harmonic expansion to those of the double Fourier series so as to be capable with arbitrary high degree and order. Next, we created a new method to transform inversely a given double Fourier series to the corresponding surface spherical harmonic expansion. The key of the new method is a couple of new recurrence formulas to compute the inverse transformation coefficients: a decreasing-order, fixed-degree, and fixed-wavenumber three-term formula for general terms, and an increasing-degree-and-order and fixed-wavenumber two-term formula for diagonal terms. Meanwhile, the two seed values are analytically prepared. Both of the forward and inverse transformation procedures are confirmed to be sufficiently accurate and applicable to an extremely high degree/order/wavenumber as \(2^{30}\,{\approx }\,10^9\). The developed procedures will be useful not only in the synthesis and analysis of the spherical harmonic expansion of arbitrary high degree and order, but also in the evaluation of the derivatives and integrals of the spherical harmonic expansion.     

SGG重力场球谐系数正则解的误差估计方法

重力梯度为重力位的二阶导数,可以通过星载梯度仪进行观测。重力场球谐函数系数可以通过正则化方法由重力梯度算出。本文在对正则化方法分析的基础上提出了估计球谐函数系数正则解误差的方法,为我国今后发射重力梯度卫星提供技术准备。     

Ellipsoidal corrections to order e 2 of geopotential coefficients and Stokes’ formula

Assuming that the gravity anomaly and disturbing potential are given on a reference ellipsoid, the result of Sjöberg (1988, Bull Geod 62:93–101) is applied to derive the potential coefficients on the bounding sphere of the ellipsoid to order e ² (i.e. the square of the eccentricity of the ellipsoid). By adding the potential coefficients and continuing the potential downward to the reference ellipsoid, the spherical Stokes formula and its ellipsoidal correction are obtained. The correction is presented in terms of an integral over the unit sphere with the spherical approximation of geoidal height as the argument and only three well-known kernel functions, namely those of Stokes, Vening-Meinesz and the inverse Stokes, lending the correction to practical computations. Finally, the ellipsoidal correction is presented also in terms of spherical harmonic functions. The frequently applied and sometimes questioned approximation of the constant m, a convenient abbreviation in normal gravity field representations, by e ²/2, as introduced by Moritz, is also discussed. It is concluded that this approximation does not significantly affect the ellipsoidal corrections to potential coefficients and Stokes formula. However, whether this standard approach to correct the gravity anomaly agrees with the pure ellipsoidal solution to Stokes formula is still an open question.     

不同椭球变换的物理方法

本文采用物理方法讨论不同椭球变换的问题，论证了重力位低阶球谐系数与大地坐标系转换参数的关系，说明几何法转换模型（即大地坐标微分公式）不过是取至二阶次的谐函微分公式，为了提高转换精度，应该采用较高阶次的物理法微分公式（即球面正交多项式）作为不同系统的转换模型。     

Anholonomity when using the development method for the reduction of observations to the reference ellipsoid

For computing the geodetic coordinates ϕ and γ on the ellipsoid one needs information of the gravity field, thus making it possible to reduce the terrestrial observations to the reference surface. Neglect of gravity field data, such as deflections of the vertical and geoid heights, results in misclosure effects, which can be described using the object of anholonomity.