Surface Gravity Disturbances Research Articles

Regional gravity field recovery of the void areas using SGG-derived surface residual gravity disturbances based on least-squares collocation: a case study in Iran

Objective of this paper is regional gravity field modelling of the Earth combining surface gravity disturbances and satellite gravity gradiometry (SGG) data based on the least-squares collocation (LSC) over the void areas in Iran, where the existing gravity data are not intensified or available. The gravity disturbances are obtained at minimum geocentric sphere, which is inside the topographic masses by joint inversion of the four high-accuracy components (Vxx, Vyy, Vzz, Vxz) of the GOCE SGG data using second derivatives of the Hotine integral. The GOCE-only global geopotential model (GGM) GO_CONS_GCF_2_TIM_R5 up to degree and order 60 is applied to reduce the effect of the omitted distant-zone, which is restored after calculation. The residual terrain model (RTM) is removed to estimate the effect of omission error of the GGM and provide mass-free space outside the boundary surface. The GOCE-based gravity disturbances at minimum geocentric sphere are continued upward on the Earth's surface based on planar approximation of the Poisson integral and subtracted from terrestrial gravity measurements to compute the residual gravity disturbances. The surface residual gravity disturbances are gridded in different cell-sizes according to the distribution of input gravity data and the locations of the empty cells that represent the void areas are determined. The height-uncorrelated residual values are applied to predict the residual gravity disturbance for every empty cell based on the LSC method using analytical local covariance function. The calculated results yielded an improvement of about 18.1% to 24.6% in gravity field recovery compared to the only RTM-corrected GGM-based solutions for two selected test areas in Iran.

IRG2016: RBF-based regional geoid model of Iran

The aim of this research is the optimal determination of the regional geoid model of Iran based on radial basis functions (RBFs). In this case, the type and number of RBFs, their horizontal positions, depths, and unknown coefficients must be properly determined. The quality of calculations strongly depends on the correct choice of these unknown parameters. Given the precise geocentric position of any point on the Earth’s surface with the beginning of the global navigation satellite system (GNSS), the surface gravity disturbances were used to calculate the height anomaly according to Molodensky’s theory. The residual surface gravity disturbances derived by subtracting the global gravitational model EIGEN-6C4 up to degree and order 360 were applied to determine the unknown RBF parameters using the stabilized orthogonal matching pursuit (SOMP) algorithm. Based on this iterative sparse approach, non-zero components of unknown RBF parameters having the maximum recoverable energy for the desired signal were found at each iteration. The SOMP algorithm was applied for optimal determination of the proper basis functions since each unknown RBF coefficient is related to a specific basis function. Only the RBFs representing the best solution to the problem were selected at each iteration, then several new RBFs were added at suitable positions to enhance the calculation result. The new RBF-based regional geoid model entitled IRG2016 was calculated by applying the geoid-to-quasigeoid corrections to the height anomaly. The IRG2016 was fitted to 1288 GNSS/levelling control points over Iran, by applying the polynomial corrector surface. Relying on this new strategy, the calculated height reference surface shows an RMS value of approximately 0.23 m for the difference in geoidal height at the independent control points, which is comparable with the last Stokes-based geoid model.

A 10 km-resolution synthetic Venus gravity field model based on topography

A high resolution gravity field model is extremely important in the exploration of Venus. In this paper, we present a 3-dimensional Venus gravity field VGM2014 constructed by using the latest gravity and topography models, residual terrain model (RTM) and the Airy–Heiskanen isostatic compensation model. The VGM2014 is the first 10km scale Venus gravity field model; the final results are representations of the 3-dimensional surface gravity accelerations and gravity disturbances for Venus. We found that the optimal global compensation depth of Venus is about 60km, and the crustal density is potentially less than the commonly accepted value of 2700–2900kgm−3. This model will be potentially beneficial for the precise orbit determination and landing navigation of spacecraft around Venus, and may be utilized as a priori model for Venus gravity field simulation and inversion studies. The VGM2014 does not incorporate direct gravity information beyond degree 70 and it is not recommended for small-scale geophysical interpretation.

A novel scheme for solving the oblique derivative boundary-value problem

The paper is devoted to a novel scheme for solving the boundary-value problem (BVP) with the oblique derivative boundary condition (BC). In this approach, the oblique derivative in the BC is decomposed into its normal and two tangential components which are approximated by means of numerical solution values. Then the numerical scheme by the finite volume method is developed and testing numerical experiments are done. The obtained numerical solutions are compared to the exact one to show that the proposed method is second order accurate. Afterwards, the algorithm is applied to solving the fixed gravimetric BVP, namely, the numerical solution is sought in a domain bounded by a part of the Earth’s surface (Himalaya region or Slovakia), four side boundaries and a corresponding upper boundary at the satellite level. On the Earth’s surface, the oblique derivative BC in the form of surface gravity disturbances from the EGM2008, DTU10-GRAV or the detailed gravity mapping is taken into account. On the upper and side boundaries, the Dirichlet BC from the EGM2008 or GOCO03s is applied. The disturbing potential as a direct numerical result is compared with the solution to the more common BVP with the Neumann BC considered on the Earth’s surface. All numerical experiments show better agreement of the solution to the BVP with the oblique derivative BC than solution to the BVP with the Neumann BC in comparison with the disturbing potential obtained by a different mathematical approach. In area of Slovakia, when applying the GPS/levelling test at 61 points, we have gained 1.7 cm improvement in favour of the standard deviation of residuals of quasigeoidal heights obtained by solving the BVP with the oblique derivative BC in comparison with the solution to BVP with the Neumann BC.

Ellipsoidal Neumann geodetic boundary-value problem based on surface gravity disturbances: case study of Iran

An ellipsoidal Neumann type geodetic boundary-value problem (GBVP) for the computation of disturbing potential on the surface of the Earth based on the surface gravity disturbance as the boundary data is formulated. The solution methodology of the GBVP can be algorithmically summarized as follows: (i) using global navigation satellite systems (GNSS) coordinates of the gravity stations, the surface gravity disturbances are generated as the boundary data. (ii) Applying the deflection correction to the gravity disturbances to arrive at the derivative of the surface disturbing potential along the ellipsoidal normal. (iii) Removing the low frequencies part of the gravity field using harmonic expansion to degree and order 110. (iv) Using the short wavelength part of the corrected gravity disturbances derived in the previous section as the boundary data within the constructed GBVP to derive the short wavelength disturbing potential over the Earth surface. (v) The computed shortwave length signals of disturbing potentials are converted to disturbing potential values by restoring the removed effects.

Application of the BEM approach for a determination of the regional marine geoid model and the mean dynamic topography in the Southwest Pacific Ocean and Tasman Sea

Application of the BEM approach for a determination of the regional marine geoid model and the mean dynamic topography in the Southwest Pacific Ocean and Tasman SeaWe apply a novel approach for the gravimetric marine geoid modelling which utilise the boundary element method (BEM). The direct BEM formulation for the Laplace equation is applied to obtain a numerical solution to the linearised fixed gravimetric boundary-value problem in points at the Earth's surface. The numerical scheme uses the collocation method with linear basis functions. It involves a discretisation of the Earth's surface which is considered as a fixed boundary. The surface gravity disturbances represent the oblique derivative boundary condition. The BEM approach is applied to determine the marine geoid model over the study area of the Southwest Pacific Ocean and Tasman Sea using DNSC08 marine gravity data. The comparison of the BEM-derived and EGM2008 geoid models reveals that the geoid height differences vary within -25 and 18 cm with the standard deviation of 6 cm. The DNSC08 sea surface topography data and the new marine geoid are then used for modelling of the mean dynamic topography (MDT) over the study area. The local vertical datum (LVD) offsets estimated at 15 tide-gauge stations in New Zealand are finally used for testing the coastal MDT. The average value of differences between the MDT and LVD offsets is 1 cm.

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).

Direct BEM for high-resolution global gravity field modelling

The paper presents a high-resolution global gravity field modelling by the boundary element method (BEM). A direct BEM formulation for the Laplace equation is applied to get a numerical solution of the linearized fixed gravimetric boundary-value problem. The numerical scheme uses the collocation method with linear basis functions. It involves a discretization of the complicated Earth’s surface, which is considered as a fixed boundary. Here 3D positions of collocation points are simulated from the DNSC08 mean sea surface at oceans and from the SRTM30PLUS_V5.0 global topography model added to EGM96 on lands. High-performance computations together with an elimination of the far zones’ interactions allow a very refined integration over the all Earth’s surface with a resolution up to 0.1 deg. Inaccuracy of the approximate coarse solutions used for the elimination of the far zones’ interactions leads to a long-wavelength error surface included in the obtained numerical solution. This paper introduces an iterative procedure how to reduce such long-wavelength error surface. Surface gravity disturbances as oblique derivative boundary conditions are generated from the EGM2008 geopotential model. Numerical experiments demonstrate how the iterative procedure tends to the final numerical solutions that are converging to EGM2008. Finally the input surface gravity disturbances at oceans are replaced by real data obtained from the DNSC08 altimetryderived gravity data. The ITG-GRACE03S satellite geopotential model up to degree 180 is used to eliminate far zones’ interactions. The final high-resolution global gravity field model with the resolution 0.1 deg is compared with EGM2008.

The quasigeoid modelling in New Zealand using the boundary element method

We compile a quasigeoid model at the study area of New Zealand using theboundary element method (BEM). The direct BEM formulation for the Laplace equationis applied to obtain a numerical solution to the linearized fixed gravimetric boundaryvalueproblem in points at the Earth’s surface. The numerical scheme uses the collocationmethod with linear basis functions. It involves a discretisation of the Earth’s surfacewhich is considered as a fixed boundary. The surface gravity disturbances represent theoblique derivative boundary condition. The geocentric positions of the collocation pointsare determined combining the digital elevation data and the a priori quasigeoid model(onshore) and the mean sea surface topography (offshore). In our numerical realization,we use the global elevation data from SRTM30PLUS V5.0, the detailed DTM of NewZealand, the EGM2008 quasigeoid heights, and the mean sea surface topography fromthe DNSC08 marine database. The gravity disturbances are computed using two heterogeneousgravity data sets: the altimetry-derived gravity anomalies from the DNSC08gravity database (offshore) and the observed ground gravity anomalies from the GNSScience gravity database (onshore). The transformation of gravity anomalies to gravitydisturbances is realized using the quasigeoid heights calculated from the EGM2008 globalgeopotential model. The new experimental quasigeoid model NZQM2010 is compiled atthe study area of New Zealand bounded by the parallels of 34 and 47.5 arc-deg southernlatitude and the meridians of 166 and 179 arc-deg eastern longitude. The least-squaresanalysis is applied to combine the gravimetric solution with GPS-levelling data using a7-parameter model. NZQM2010 is validated using GPS-levelling data and compared withthe existing regional and global quasigeoid models NZGeoid2009 and EGM2008.&nbsp

Finite element method for solving geodetic boundary value problems

The goal of this paper is to present the finite element scheme for solving the Earth potential problems in 3D domains above the Earth surface. To that goal we formulate the boundary-value problem (BVP) consisting of the Laplace equation outside the Earth accompanied by the Neumann as well as the Dirichlet boundary conditions (BC). The 3D computational domain consists of the bottom boundary in the form of a spherical approximation or real triangulation of the Earth’s surface on which surface gravity disturbances are given. We introduce additional upper (spherical) and side (planar and conical) boundaries where the Dirichlet BC is given. Solution of such elliptic BVP is understood in a weak sense, it always exists and is unique and can be efficiently found by the finite element method (FEM). We briefly present derivation of FEM for such type of problems including main discretization ideas. This method leads to a solution of the sparse symmetric linear systems which give the Earth’s potential solution in every discrete node of the 3D computational domain. In this point our method differs from other numerical approaches, e.g. boundary element method (BEM) where the potential is sought on a hypersurface only. We apply and test FEM in various situations. First, we compare the FEM solution with the known exact solution in case of homogeneous sphere. Then, we solve the geodetic BVP in continental scale using the DNSC08 data. We compare the results with the EGM2008 geopotential model. Finally, we study the precision of our solution by the GPS/levelling test in Slovakia where we use terrestrial gravimetric measurements as input data. All tests show qualitative and quantitative agreement with the given solutions.

Numerical solution of the linearized fixed gravimetric boundary-value problem

The fixed gravimetric boundary-value problem (FGBVP) represents an exterior oblique derivative problem for the Laplace equation. Terrestrial gravimetric measurements located by precise satellite positioning yield oblique derivative boundary conditions in the form of surface gravity disturbances. In this paper, we discuss the boundary element method (BEM) applied to the linearized FGBVP. In spite of previous BEM approaches in geodesy, we use the so-called direct BEM formulation, where a weak formulation is derived through the method of weighted residuals. The collocation technique with linear basis functions is applied for deriving the linear system of equations from the arising boundary integral equations. The nonstationary iterative biconjugate gradient stabilized method is used to solve the large-scale linear system of equations. The standard MPI (message passing interface) subroutines are implemented in order to perform parallel computations. The proposed approach gives a numerical solution at collocation points directly on the Earth’s surface (on a fixed boundary). Numerical experiments deal with (i) global gravity field modelling using synthetic data (surface gravity disturbances generated from a global geopotential model (GGM)) (ii) local gravity field modelling in Slovakia using observed gravity data. In order to extend computations, the memory requirements are reduced using elimination of the far-zone effects by incorporating GGM or a coarse global numerical solution obtained by BEM. Statistical characteristics of residuals between numerical solutions and GGM confirm the reliability of the approach and indicate accuracy of numerical solutions for the global models. A local refinement in Slovakia results in a local (national) quasigeoid model, which when compared with GPS-levelling data, does not make a large improvement on existing remove-restore-based models.

Obtaining Earth surface gravity disturbances from a GPS-based ‘high-low’ satellite-to-satellite tracking experiment

A ‘high-low’ satellite-to-satellite tracking measurement process, involving a single low orbiter (300 km altitude) and the GPS constellation (20 189 km altitude), has been proposed by Jekeli & Upadhyay (1990) to estimate Earth surface gravity disturbances. This paper examines the potentiality of such a mission by initially conducting an error analysis based on signal-to-noise ratios of the spectra of intersatellite line-of-sight (los) accelerations possessing an assumed white noise process. Using a 60 day mission duration, los acceleration data rates ranging from 1 to 60s and data noise levels from 0.1 to 0.3 mGal, corresponding errors of commission and omission propagating into the spectrum of a desired set of 1° and/or 2° mean Earth surface vertical gravity disturbances (residual to harmonic degree and order n=m=6) were obtained. The resulting total estimated errors (in a global average sense) of the 1° mean disturbances ranged from 6 to 12 mGal, and for 2° values the range was 1.5 to 6 mGal, the lower error bounds corresponding to a 1 s data rate and a 0.1 mGal data noise level. An optimal least-squares collocation simulation was next carried out. The simulated los accelerations reflected all gravimetric perturbations due to an adopted true gravity field extending to n=m= 180. A downward continuation approach was used in conjunction with singular value decomposition techniques to predict 1° mean disturbances over two separate regions on the Earth's surface. For a 10° by 10° region possessing a fairly tranquil local gravity field, the overall RMS of 100 ‘true-predicted’ mean differences, using the aforementioned los data rates and noise levels, ranged from 5 to 9 mGal and was only 1 mGal when the data were noise-free. For a more turbulent 15° by 20° region, the RMS of 300 differences was 1.5 mGal using perfect data and was in the 6 to 11 mGal range using noisy data. Again, the lower error bounds corresponded to 1 s and 0.1 mGal data rate and noise levels. If the low satellite takes on a polar orbit one could therefore obtain respectable 1° mean disturbances and highly reliable 2° mean values over all previously inaccessible Earth surface areas.