On estimating gravity anomalies— A comparison of least squares collocation with conventional least squares techniques

The least squares collocation algorithm for estimating gravity anomalies from geodetic data is shown to be an application of the well known regression equations which provide the mean and covariance of a random vector (gravity anomalies) given a realization of a correlated random vector (geodetic data). It is also shown that the collocation solution for gravity anomalies is equivalent to the conventional least-squares-Stokes' function solution when the conventional solution utilizes properly weighted zero a priori estimates. The mathematical and physical assumptions underlying the least squares collocation estimator are described.     

A modified Wiener-type filter for geodetic estimation problems with non-stationary noise

 One of the most basic and important tools in optimal spectral gravity field modelling is the method of Wiener filtering. Originally developed for applications in analogue signal analysis and communication engineering, Wiener filtering has become a standard linear estimation technique of modern operational geodesy, either as an independent practical tool for data de-noising in the frequency domain or as an integral component of a more general signal estimation methodology (input–output systems theory). Its theoretical framework is based on the Wiener–Kolmogorov linear prediction theory for stationary random fields in the presence of additive external noise, and thus it is closely related to the (more familiar to geodesists) method of least-squares collocation with random observation errors. The main drawback of Wiener filtering that makes its use in many geodetic applications problematic stems from the stationarity assumption for both the signal and the noise involved in the approximation problem. A modified Wiener-type linear estimation filter is introduced that can be used with noisy data obtained from an arbitrary deterministic field under the masking of non-stationary random observation errors. In addition, the sampling resolution of the input data is explicitly taken into account within the estimation algorithm, resulting in a resolution-dependent optimal noise filter. This provides a more insightful approach to spectral filtering techniques for noise reduction, since the data resolution parameter has not been directly incorporated in previous formulations of frequency-domain estimation problems for gravity field signals with discrete noisy data. Received: 1 November 2000 / Accepted: 19 June 2001     

Integration of the Monte Carlo Covariance Estimation Strategy into Tailored Solution Procedures for Large-Scale Least Squares Problems

Large-scale least squares problems require tailored numerical techniques to overcome the computational burden. For these types of problems, iterative strategies are suitable because of their flexibility and effectiveness. The only shortcoming of iterative strategies in least squares estimation is that the inverse of the normal equation matrix as the carrier of the covariance information is either unavailable or very expensive to compute. This paper presents algorithms based on Monte Carlo integration, which can be incorporated very efficiently into iterative solvers and which are demonstrated to close the aforementioned gap. Tailored strategies for different types of solution techniques with respect to normal equations, observation equations, and combined models are treated. Finally, the paper presents new criteria to define confidence regions for the estimated covariance matrix of the parameters, as well as for all additional derived quantities. In a case study these techniques are applied to simulated GOCE data, where satellite gravity gradiometry and satellite-to-satellite tracking information are combined for reconstructing the gravity field. The problem of deriving the covariance matrix of gravity fields with high spatial resolution by combined iterative estimation processes, unsolved until now, is treated.     

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.     

最小二乘法求解三类卫星重力梯度边值问题的研究

研究了最小二乘法求解3类卫星重力梯度边值问题的理论和方法，给出了3类梯度观测值{Гzz}、{Гxz、Гyz}和{Гxx-Гyy,2Гxy}对应边值问题解的核函数严密表达式。模拟试算结果表明，最小二乘法求解的卫星重力梯度积分公式用于恢复地球重力场是有效而严密的。     

Strategies for spatial and temporal extrapolation and interpolation of wet delay

. The excess radio-path delay due to the atmospheric water vapor, the wet delay, can be derived from water vapor radiometer (WVR) measurements. WVR data used for external calibration of space geodetic measurements are not always acquired in the directions of the space geodetic signal sources, thus extrapolation and interpolation methods for the wet delay are needed. We evaluate three different methods using approximately 10 days of WVR measurements. Two methods, the gradient method and turbulence method, use the directional information in the data, while the third method used is linear regression in time regardless of the direction of the observations. The turbulence method yielded at least 10% less RMS estimation error than the errors from the other two methods. Received: 20 May 1997 / Accepted: 15 December 1997     

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.     

A method of error adjustment for marine gravity with application to Mean Dynamic Topography in the northern North Atlantic

International compilations of marine gravity, such as the International Gravity Bureau (BGI) contain tens of millions of point data. Lemoine et al. (The Development of the Joint NASA GSFC and the National Imagery and Mapping Agency (NIMA) Geopotential Model EGM96, NASA/TP-1998-206861) chose not to include any marine gravity in the construction of the global gravity model EGM96. Instead they used synthetic anomalies derived from altimetry, so that no independent information about Mean Dynamic Topography (MDT) can be deduced. Software has been developed not only to identify and correct those aspects of marine gravity data that are unreliable, but to do so in a way that can be applied to very large, ocean-wide data sets. First, we select only straight-line parts of ship-tracks and fit each one with a high-degree series of Chebyshev polynomials, whose misfit standard deviation is σ _line and measures the random error associated with point gravity data. Then, network adjustment determines how the gravity datum is offset for each survey. A free least squares adjustment minimises the gravity anomaly mismatch at line-crossing points, using σ _line to weight the estimate for each line. For a long, well crossed survey, the instrumental drift rate is also adjusted. For some 42,000 cross-over points in the northern Atlantic Ocean, network adjustment reduces the unweighted standard deviation of the cross-over errors from 4.03 to 1.58 mGal; when quality weighted, the statistic reduces from 1.32 to 0.39 mGal. The geodetic MDT is calculated combining the adjusted gravity anomalies and satellite altimetry, and a priori global ocean model through a new algorithm called the Iterative Combination Method. This paper reports a first demonstration that geodetic oceanography can characterise the details of basin wide ocean circulation with a resolution better than global ocean circulation models. The result matches regional models of ocean circulation from hydrography measurements (Geophys Res Lett 29:1896, 2002; J Geophys Res 108:3251, 2003).     

由星载GPS相位数据确定地球重力场模型若干问题研究

人造地球卫星在地球引力场中运动,可以探测地球重力场的长波信息。随着GPS技术的发展,星载GPS技术日趋成熟,因此由星载GPS相位数据确定地球重力场模型是当前国际地学研究的热点之一。本文给出了确定地球重力场模型中的星载GPS星地相位双差观测量,阐述了Cowell II数值轨道积分公式,导出了参数估计中星地双差观测量的偏导数,利用分块Bayes最小二乘参数估计地球引力场位系数等有关参数。     

The DNSC08GRA global marine gravity field from double retracked satellite altimetry

Satellite radar altimetry has been monitoring the earth’s oceans from space for several decades. However, only the GEOSAT and ERS-1 geodetic mission data recorded more than a decade ago provide altimetry with adequate spatial coverage to derive a high-resolution marine gravity field. The original geodetic mission data suffer from degradation in quality and coverage close to the coast and in Polar Regions as well as the occasionally wrongly retracking of these, even in the open ocean. In order to improve the quality of these geodetic mission data and to derive a new improved global marine gravity field called DNSC08GRA, a new double retracking technique for analyzing the waveform data has been developed. Multiple retracking allows the system to retrack more data to increase the spatial coverage of the data. Subsequently, a second retracking run is used to enhance the SSH determination by using information from the first fitting to inform the second set of retrackers about smoothly varying sea state parameters. The development of the new global marine gravity field DNSC08GRA is described in this paper. Besides application of new retracking techniques the radar altimetry has been processed using EGM2008 as reference and augmented with ArcGP gravity data and laser altimetry from ICESat to close the Polar gap. DNSC08GRA is seen to perform significantly better than previous global marine gravity field like KMS02. The improvement in accuracy is better than 20% in general, but in coastal regions, the improvement is in many places of the order of 40–50% compared to older global marine gravity field KMS02.     

GSPP a general surface representation module designed for geodesy

GSPP is a computer program system which has been developed for the purposes of automatically determining and representing gravity field surfaces like the geoid, the field of gravity anomalies or deviations of the vertical at prescribed altitude, etc. The system processes gravity field information given by a heterogeneous set of linear functionals of the anomalous potential superimposed by noise, and provides automatically gravity field surfaces in terms of profiles, contour maps and/or 3-dimensional representations. The solution is generally based on least-squares collocation; for a homogeneous data set, a simple weighted average interpolation is available as well. Based on the given data, surface function values at the grid points of a regular rectangular grid are predicted. The representation of the surfaces is smooth using bicubic spline functions. GSPP has a control unit which performs all necessary decision processes and such reduces the user’s decision making to a minimum. The system has been designed for geodetic purposes only; however, because of its versatility and flexibility it presents itself also for applications in other geosciences.     

Improvement of theories and methods in geodetic joint inversion

This paper first establishes the prior globe dynamical model by geophysics,which is a solid earth elastic deformation model.Then,the parameters of the globe dynamic model can be obtained by inverting the geodetic data.The inverse method can be used in seismology and geology,and to make earthquake prediction.     

Improvement of theories and methods in geodetic joint inversion

This paper first establishes the prior globe dynamical model by geophysics, which is a solid earth elastic deformation model. Then, the parameters of the globe dynamic model can be obtained by inverting the geodetic data. The inverse method can be used in seismology and geology, and to make earth-quake prediction.     

Using radial basis functions in airborne gravimetry for local geoid improvement

Radial basis functions (RBFs) have been used extensively in satellite geodetic applications. However, to the author’s knowledge, their role in processing and modeling airborne gravity data has not yet been fully advocated or extensively investigated in detail. Compared with satellite missions, the airborne data are more suitable for these kinds of localized basis functions especially considering the following facts: (1) Unlike the satellite missions that can provide global or near global data coverage, airborne gravity data are usually geographically limited. (2) It is also band limited in the frequency domain. (3) It is straightforward to formulate the RBF observation equations from an airborne gravimetric system. In this study, a set of band-limited RBF is developed to model and downward continue the airborne gravity data for local geoid improvement. First, EIGEN6c4 coefficients are used to simulate a harmonic field to test the performances of RBF on various sampling, noise, and flight height levels, in order to gain certain guidelines for processing the real data. Here, the RBF method not only successfully recovers the harmonic field but also presents filtering properties due to its particular design in the frequency domain. Next, the software was tested for the GSVS14 (Geoid Slope Validation Survey 2014) area in Iowa as well as for the area around Puerto Rico and the US Virgin Islands by use of the real airborne gravity data from the Gravity for the Redefinition of the American Vertical Datum (GRAV-D) project. By fully utilizing the three-dimensional correlation information among the flight tracks, the RBF can also be used as a data cleaning tool for airborne gravity data adjustment and cleaning. This property is further extended to surface gravity data cleaning, where conventional approaches have various limitations. All the related numerical results clearly show the importance and contribution of the use of the RBF for high- resolution local gravity field modeling.     

利用最小二乘直接法反演卫星重力场模型的MPI并行算法

针对海量卫星重力数据反演高阶次地球重力场模型的密集型计算任务与高内存耗用问题,基于MPI实现了最小二乘直接法恢复高阶次位系数的并行算法。引入并行读写、分块存储与分块计算等方式完成了设计矩阵的构建、法方程的形成与求解等密集型计算任务的并行算法,数值计算结果表明三者的并行相对效率峰值可分别达到95%、68%、63%。利用GOCE轨道跟踪和径向扰动重力梯度数据(共518 400个历元)分别反演了120、240阶次地球重力场模型,计算时间仅为40 min、7 h,内存耗用峰值仅为290 MB、1.57 GB;采用与GOCE同等噪声水平的观测数据恢复的重力场模型精度与GOCE已发布模型的解算精度相一致,联合GRACE和GOCE的解算模型能够实现二者独立信息的频谱互补,表明本文方法可高效稳定地恢复高阶次地球重力场模型。     

Height datum definition,height datum connection and the role of the geodetic boundary value problem

Vertical datum definition is identical with the choice of a potential (or height) value for the fundamental bench mark. Also the connection of two adjacent vertical datums poses no principal problem as long as the potential (or height) value of two bench marks of the two systems is known and they can be connected by levelling. Only the unification of large vertical datums and the connection of vertical datums separated by an ocean remains difficult. Two vertical datums can be connected indirectly by means of a combination of precise geocentric positions of two points, as derived by space techniques, their potential (or height) value in the respective height datum and their geoid height difference. The latter requires the solution of the linear geodetic boundary value problem under the assumption that potential and gravity anomalies refer to a variety of height datums. The unknown off-sets between the various datums appear in the solution inside and outside the Stokes integral and can be estimated in a least squares adjustment, if geocentric positions, levelled heights and adequate gravity material are available for all datum zones. The problem can in principle also be solved involving only two datums, in case a precise global gravity field becomes available purely from satellite methods.     

Least squares collocation and regularization

Certain geodetic problems such as the downward continuation of gravity information from satellite or aerial altitudes to the surface of the earth or the inverse Stokes problem are improperly posed in the sense that the best approximate solution does not continuously depend on the given observations. In order to obtain a stable solution a technique of regularization is discussed which can be shown to be identical to the method of least squares collocation. The characteristic features of regularization are analysed transferring the problem through a singular value decomposition into the spectral domain which allows an easy interpretation as a special method of filtering. In practical applications the stability depends not only on the observational errors including the computer round-off but in the same way on the number of observations and their distribution. The regularized solution should achieve a proper trade-off between sufficient smoothness and highest possible resolution with a limit defined through the internal accuracy of the computer.     

我国大地测量学的进展和展望

回顾了我国大地测量工作的进展。面向 2 1世纪前期的我国经济和国防建设及科技和社会发展 ,展望了我国新世纪的大地测量 ,提出应逐步进入精确、动态、实时的现代化体系 ,即完善国家三维空间大地网 ;建立 GPS综合服务体系 ;提供导航和定位服务 ;测定地壳运动、电离层参数、大气中可降水份等信息 ;精化中国地区重力场参数 ;建立新的国家重力基准网 ;完成分米级精度的中国似大地水准面的推算 ;积极开展海洋和空间大地测量 ,为资源、环境的管理以及防灾监测做出应有的贡献。     

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.     

Solving three types of satellite gravity gradient boundary value problems by least-squares

The principle and method for solving three types of satellite gravity gradient boundary value problems by least-squares are discussed in detail. Also, kernel function expressions of the least-squares solution of three geodetic boundary value problems with the observations {Гzz}, {Гxz, Гyz} and {Гxx ? Гyy, 2Гxy} are presented. From the results of recovering gravity field using simulated gravity gradient tensor data, we can draw a conclusion that satellite gravity gradient integral formulas derived from least-squares are valid and rigorous for recovering the gravity field.