Global geopotential models from Satellite Laser Ranging data with geophysical applications : A review

How to cite this article: Botai CM, Combrinck L. Global geopotential models from Satellite Laser Ranging data with geophysical applications: A review. S Afr J Sci. 2012;108(3/4), Art. #662, 10 pages. http:// dx.doi.org/10.4102/sajs. v108i3/4.662 The launch of artificial satellites (as early as in 1957), specifically the launch of the first laser tracked satellite, Beacon-B, in 1964, has provided data sets which have allowed researchers to probe the long to medium components of the gravitational field of the Earth. In particular, observational data recorded at satellite laser ranging tracking stations have since been used to develop models that quantify the global long-wavelength and medium-wavelength gravity field of the Earth. Currently, literature reviewing gravity field models with geophysical applications is scarce and not up to date. The most recent review paper was published more than a decade ago. In the interim, there has been an unprecedented increase in gravity field modelling, which can be attributed to the deployment of new and dedicated satellite missions. As a result, a number of existing geopotential models have been improved and new models have been developed. Each of these models differs in accuracy and spatial-temporal scale. This review extends the earlier review of gravity field models, by incorporating up-to-date research efforts in geopotential modelling with geophysical applications in oceanography, hydrology, geodesy and solid Earth science.

The launch of artificial satellites (as early as in 1957), specifically the launch of the first laser tracked satellite, Beacon-B, in 1964, has provided data sets which have allowed researchers to probe the long to medium components of the gravitational field of the Earth.In particular, observational data recorded at satellite laser ranging tracking stations have since been used to develop models that quantify the global long-wavelength and medium-wavelength gravity field of the Earth.Currently, literature reviewing gravity field models with geophysical applications is scarce and not up to date.The most recent review paper was published more than a decade ago.In the interim, there has been an unprecedented increase in gravity field modelling, which can be attributed to the deployment of new and dedicated satellite missions.As a result, a number of existing geopotential models have been improved and new models have been developed.Each of these models differs in accuracy and spatial-temporal scale.This review extends the earlier review of gravity field models, by incorporating up-to-date research efforts in geopotential modelling with geophysical applications in oceanography, hydrology, geodesy and solid Earth science.

Introduction
The concept of the Earth's gravity field is often described through gravity or gravitational potential.Alternatively, the definition of gravity can be viewed in terms of the cause and effect of gravity.To this end, gravity could be described as an attraction that causes acceleration amongst objects on or near the Earth's surface.The gravitational potential of the Earth is the quantity of energy that is associated with the position of a unit mass in the gravitational field of the Earth. 1 Earth itself is a complex dynamic system driven by many geophysical processes.These include the coupled atmosphere-ocean system, varying mass distribution of ice and the isostatic correction from the glacial loading of the last Ice Age and mobile tectonic plates. 2 In addition, internal mass distribution is often controlled by thermal convection of the core mantle. 2Some of the geophysical processes taking place within Earth's system act to redistribute the Earth's mass thereby changing the motion of the solid Earth relative to the centre, as well as causing spatial and time-dependent variations of the gravitational field of the Earth.Observations of this variability of the Earth's gravity field using artificial satellites via the Satellite Laser Ranging (SLR) technique and other geodetic techniques, Global Navigation Satellite Systems (GNSS), Very Long Baseline Interferometry and Doppler Orbitography and Radiolocation Integrated by Satellite (DORIS) can be used to study a wide variety of geophysical processes that involve changes in mass. 3e Earth's gravity field plays an important role in understanding dynamic processes taking place within the Earth system.These processes include interactions between the cryosphere, hydrosphere, atmosphere and ocean at spatial scales ranging from a few metres to continental and global scales.Temporal scales of these dynamic processes range from an hour to geological time. 4The gravity field of the Earth can also be used to determine global ocean circulation which relates to global climate change.Gravitational field changes may be used in detecting mass shifts in the Earth's interior, which might be associated with movements on the Earth's surface. 4haracteristics of the gravitational field are often defined from artificial satellite tracking data.In particular, artificial satellites are used to detect long-wavelength components of the gravitational potential.A gravitational potential is often expressed as a series expansion of spherical harmonics known as a global geopotential model (GGM).In this paper, a general review of GGMs, released from 1990 to date, is presented.

The concept of Satellite Laser Ranging
Gravity field models are often derived by use of data collected from the SLR observations technique.This is a technique that measures the two-way travel time of a short laser pulse which is reflected by an orbiting satellite.SLR as an observational method is possible through satellites equipped with retroreflectors made from glass prisms.An example of retroreflectors is illustrated in Figure 1.In a typical SLR system, a transmitting telescope emits short laser pulses with energy between 10 mJ and 100 mJ at a repetition frequency ranging between 5 Hz and 20 Hz.Some modern systems have lower power levels and higher firing rates (of up to 2 kHz).Laser pulses which illuminate any of the retroreflectors are reflected back to the ground station where they are collected via the receiving telescope and detected by a photomultiplier or a solid state photo diode.The measurement of laser ranges from laser tracking stations to a retroreflector on an orbiting satellite, for example, the time-of-flight (TOF) is often measured by either a time interval counter or an epoch timer.
A basic equation representing the approximate TOF is given by: where c is the speed of light, d is the round-trip distance from the SLR station to the target satellite retroreflector and t is the TOF.In order to obtain the best possible range precision from the ground station to the satellite, additional parameters and numerous corrections corresponding to internal delays in the transmission and detection systems need to be considered.Taking into account such corrections, the basic range equation given by [Eqn 1] can be expanded as in [Eqn 2] as reported in Seeber 5 : In [Eqn 2], ∆t is the measured TOF and is mostly affected by uncertainties in the signal identification.The preferred resolution for the measured TOF is often a few picoseconds.In addition, the measured TOF needs to be tied to universal time (because of the satellite's motion relative to the Earth).The ∆d 0 term corresponds to the eccentric correction on the ground, which is the intersection of the vertical axis and horizontal axis and is used as a reference point in the laser system.Similarly, ∆d S corresponds to the eccentric correction at the satellite and gives a geometrical relationship between the centre of the corner cube and the centre of mass of the satellite; the accuracy of this parameter is very difficult to obtain on satellites with irregular shapes (e.g.satellites equipped with solar panels and antennas).The ∆d b term in [Eqn 2] corresponds to the signal delay in the ground system -the geometric reference point 0 to the electrical 0 point and is often not exactly at the same point; this correctional parameter is often determined through calibration with older systems that were calibrated with respect to a defined terrestrial target.Furthermore, ∆d r is the refraction correction as a result of atmospheric conditions which affect the propagation velocity of laser pulses.Laser pulses experience a delay in the lower part of the atmosphere, which makes measurement of these parameters along the total path difficult.Therefore atmospheric models are used that incorporate variables such as SLR site pressure and temperature and are supported by measured data at the laser site. 5Lastly, η are random systematic and observation errors related to unmodelled residual effects.
The first SLR experiment campaign began in the 1960s with the development of Ruby-based SLR stations tracking satellites such as the Beacon Explorer-B.Since then numerous satellite missions have been launched for different applications, such as geodetic, Earth sensing and radio navigation, and a global network of SLR stations has been established, replacing the old Baker-Nunn optical camera tracking network. 6A historical overview of such missions is summarised in Table 1.The current global network of SLR stations involved in satellite tracking consists of over 40 stations and their global distribution is depicted in Figure 2. Most of the stations are located in the Northern Hemisphere leaving the Southern Hemisphere with weak coverage.In Africa, there are two stations: Helwan in Egypt and MOBLAS-6 (see Figure 3) located at Hartebeesthoek Radio Astronomy Observatory (HartRAO) in South Africa.The space geodetic fundamental station HartRAO is involved with the International Laser Ranging Service activities as  where is the position vector of the centre of mass of the satellite, ā g is the sum of the gravitational forces acting on the satellite, ā ng is the sum of the non-gravitational forces acting on the surface of the satellite and ā emp represents the unmodelled forces which act on the satellite because of either a functionally incorrect or an incomplete description of the various forces acting on the satellite. 5The gravitational forces (ā g ) acting on an orbiting satellite are composed of a series of perturbations, expressed as: where Ρ̅ geo is the geopotential force as a result of the gravitational attraction of the Earth; Ρ̅ set and Ρ̅ ot are perturbations as a result of solid Earth tides and ocean tides, respectively; Ρ̅ rd is the perturbation caused by the rotational deformation of the Earth; Ρ̅ smp is the perturbation caused by the Sun, Moon and other planets; and Ρ̅ rel is the perturbation as a result of general relativity. 5e non-gravitational forces acting on an orbiting satellite are given by: where Ρ̅ drag is the atmospheric drag acting on a satellite; Ρ̅ solar is the perturbation as a result of solar radiation pressure; Ρ̅ earth is the perturbation caused by Earth radiation pressure; and Ρ̅ thermal is the perturbation as a result of thermal radiation imbalances resulting from non-uniform temperature distribution on different satellite surfaces.The forces described in [Eqn 4] and [Eqn 5], together with unmodelled forces, are solved for during geopotential modelling and therefore are central, in particular, to the derivation of precise geopotential models.

Precise satellite orbit determination
Precise satellite orbit determination (POD) is one of the most essential applications of geopotential modelling.The POD process involves the estimation of position and velocity of an orbiting satellite at a specific time epoch. 7POD is used for geolocation of the satellite sensors and to measure the gravity field and its variations in time.There are currently three ways in which a satellite orbit can be calculated: dynamic, kinematic and reduced dynamic.

Dynamic orbit determination
Dynamic orbit determination 7 requires a set of tracking observations and mathematical models acting on an orbiting satellite.Here the force and satellite models are used to compute a model of satellite acceleration over a given time.
A nominal trajectory is generated analytically or numerically by integrating the acceleration model.The orbit solution is compared with the one predicted by the observations.Selected parameters of the force models acting on the satellite may be adjusted along with an initial satellite position and velocity in order to minimise the difference between the actual observations and the predicted ranges (this difference is called O-C residuals) in a least-squares sense. 7The accuracy of the dynamic orbit determination approach is highly dependent on the satellite force models.

Kinematic orbit determination
Kinematic orbit determination is a purely geometric technique that depends only on GNSS (e.g.GPS) measurements and cannot be used by SLR. 8 It does not take into account the dynamic properties (e.g.gravity field or air drag) of an orbiting satellite.Here the errors emanating from the satellite force models do not affect the accuracy of the kinematic orbit determination, but its accuracy does depend on the availability and accuracy of GNSS data.

Reduced-dynamic orbit determination
In the dynamic and kinematic methods, the accuracy of the solution may be reduced as a result of modelling errors and GNSS measurement noise, respectively.The reduceddynamic technique proposed by Yunck et al. 9 may be defined as a method that exhibits a combination of dynamic and kinematic components and that minimises the errors caused by each method.In the reduced-dynamic orbit determination approach, the kinematic components of the dynamic force models are introduced in the form of a process noise model containing two parameters -the correlation time constant T (which defines the correlation in the dynamic model error over one update interval) and the dynamic model steadystate variance V. Accuracy in this method depends on proper adjustment of the two parameters.

Global geopotential models
The perturbing potential function for the solid-body mass distribution of the Earth is often expressed in terms of a spherical harmonic expansion, obtained when solving a Laplace equation in spherical coordinates described in Tapley et al. 10 by: Here (r,φ,λ) represent the magnitude of the radius vector, the latitude and the longitude, respectively; {C̅ nm , S nm } are fully normalised spherical harmonic coefficients of degree n and order m; and Ρ̅ nm is a fully normalised associated Legendre function.The adopted gravity mass constant GM is set to 398600.4415km 3 /s 2 in most recent geopotential models. 11A typical geopotential model is often described by {C̅ nm , S nm } spherical harmonic coefficients.The values of {C̅ nm , S nm } coefficients decrease as the degree increases.For satellite-based global gravity field models the accuracy of the lower degree coefficients is typically higher than the higher degree coefficients.
A number of spherical harmonic models have been developed over the years by different research groups, for example, Ohio State University (OSU), GeoForschungs Zentrum (GFZ) Potsdam, Goddard Earth Models (GEM), Joint Gravity Models (JGM), Texas Earth Gravity models (TEG) and European Improved Gravity Model of the Earth by New Techniques (EIGEN).The existing GGMs representing the Earth's gravitational field can be classified into three groups: satellite-only, combined and tailored gravity field models.The satellite-only GGMs are derived from the analysis of the orbits of artificial Earth satellites.Numerous factors have been attributed to the inherent inaccuracies of the satellite-only models.These factors include weakening of the gravitational field with altitude, precession of the Earth-based range measurements to the satellites, the lack of continuous tracking data from the existing stations and difficulties in modelling non-gravitational and third body perturbations. 12he satellite-only models are often combined with terrestrial gravity data and marine gravity anomalies to yield highdegree combined GGMs.The combined GGMs are subjected to the same deficiencies as in satellite-only GGMs and also to other errors emanating from terrestrial gravity anomalies. 13n tailored GGMs, the spherical harmonic coefficients of the satellite-only or the combined models are often adjusted and extended to higher degrees by using previously used or unused gravity data. 14rrently a number of GGMs have been derived and made available freely to the scientific community. 15A review of gravity field models derived between 1970 and 1997 can be found in Rapp 16 .This paper reviews developments undergone in the gravity field modelling for the last two decades (i.e.1990-2010).Characteristics of these models are summarised in Table 2.The first considered model is a combined gravity field model, GRIM4C1, reported by Schwintzer et al 17 .This model was computed as a joint collaboration between the German Geodetic Research Institute (DGFI) and Groupe de Recherches de Geodesie Spatiale (GRGS).The GRIM4C1 model was derived up to degree and order 50 in terms of spherical harmonics.It incorporated GRIM4S1 satellite solution, mean gravity anomalies and Seasat altimeter derived mean geoid undulations.The OSU91A geopotential model was reported by Rapp et al 18 .This model was an upgraded version of OSU89a and OSU89b.It was computed complete to degree and order 360 in terms of spherical harmonics in a blended form.In the computation of the OSU91A, coefficients to degree 50 were based on a combined solution from the GEM-T2 model, surface gravity data and GEOSAT altimeter data.The remaining coefficients (51-360) were derived from a combined solution computed from terrestrial data, altimeter-derived anomalies and topographic or isostatic anomalies.
The Joint Gravity Model 3 (JGM3) released in 1994 was reported by Tapley et al. 19 This model was developed by NASA Goddard Space Flight Center (GSFC) and the University of Texas at Austin as part of the Topex Poseidon (T/P) project.This combined model was derived by adding the geopotential coefficients from the prelaunch model, JGM1 and their associated error covariance with GPS, SLR, DORIS tracking of T/P, laser ranging tracking of LAGEOS 2 and Stella and DORIS tracking of SPOT 2. The model was derived complete to degree and order 70.The GRIM4S4 and GRIM4C4 reported by Schwintzer et al. 20 were developed jointly by GFZ Potsdam and GRGS Toulouse/Grasse for requirements of geodetic and altimeter satellite missions.The GRIM4S4 model was derived solely from satellite tracking data complete to degree and order 70.On the other hand, the GRIM4C4 model was derived based on a least squares adjustment involving a combined solution from the GRIM4S4 model and surface gravity data from gravimetric and altimeter measurements.This model was computed complete to degree and order 72, corresponding to a spatial resolution of 555 km at the surface of the Earth. 20The GRIM4S4 and GRIM4C4 models were thought to be efficient for satellite orbit computations especially with orbit altitudes exceeding about 800 km. 20The GFZ96 geopotential model which was an upgrade of the GFZ93 and GFZ95 models was reported to provide high resolution in the history of GFZderived models. 21This combined model was computed from the then improved terrestrial data derived from a 3-year ERS-1 mean sea surface and PMG055 solution.The solution was also combined with altimeter-derived gravity anomalies and normal equations and potential coefficients of the GRIM4S4 model as the a priori model.The GFZ96 model was derived to degree and order 359.
Lemoine et al. 22 described the combined spherical harmonic model, EGM96, which is complete to degree and order 360 and corresponds to a global resolution of about 55 km.The EGM96 model was developed based on a joint collaboration between NASA-GSFC, the National Imagery and Mapping Page 5 of 10 Agency (NIMA) and the Ohio State University.EGM96 is a blend model in which three computational procedures were used.The spherical harmonic coefficients from 2-70 were derived based on a least squares adjustment involving satellite tracking data, terrestrial data and altimeter data of the ocean surface from the T/P, ERS-1 and GEOSAT missions and fill-in gravity anomalies in areas lacking data. 12rom degree 71-359, the coefficients were computed from a combined solution based on normal equations derived from the satellite tracking data which were used as a priori values.The remaining coefficients at degree 360 were taken from a quadrature combined solution derived from the a priori satellite model and ERS-1/GEOSAT altimeterderived anomalies.The EGM96 geopotential model was believed to provide a more accurate reference surface for the topography, as well as an improved orbit determination for low orbiting satellites. 22The GRIM5C1 gravity field model reported by Gruber et al. 23 was derived in a German-French joint collaboration between GFZ Potsdam and GRGS Toulouse.The model was computed up to degree and order 120.It incorporated terrestrial and airborne mean gravity anomalies, altimetric gravity anomalies from NIMA and mean gravity anomalies derived from the GRIM5S1 model.
Most of the geopotential models released from 2000 onwards were derived mostly from CHAMP, GRACE and GOCE missions plus other satellites, terrestrial and altimeter data.Geopotential models generated from the inclusion of the three satellite missions data are believed to be more accurate than prior models (e.g. they allow, with an unprecedented accuracy and resolution, the recovery of the mean sea surface topography from the difference between an altimetry-based mean sea surface height model and the gravity model's derived geoid). 24The first CHAMP geopotential model, EIGEN-1, reported by Reigber et al. 25 , was derived in a German-French collaboration complete to degree and order 119.This model was derived by use of GPS tracking and 3 months of on-board accelerometer data from CHAMP.
The EIGEN-1 geopotential model was reported to resolve the geoid and gravity with an accuracy of about 20 cm and 1 mGal, respectively, at a half-wavelength resolution of 550 km. 25 The EIGEN-2 model reported by Reigber et al. 26 was also derived in collaboration between Germany and France.This satellite-only model was derived complete to degree and order 140.The model incorporated gravity orbit perturbations exploiting GPS CHAMP satellite-to-satellite tracking and 6 months of on-board accelerometer data.The accuracy in terms of geoid and gravity for the EIGEN-2 model was reported to be about 10 cm and 0.5 mGal, respectively.
Like the CHAMP mission, the GRACE mission data set has enabled a homogeneous determination of the geopotential gravity field modelling.The first model that was derived from the GRACE data was a 'satellite-only model' called the GGM01S reported by Tapley et al. 10 This model, derived to complete degree and order 120, incorporated GRACE tracking data spanning from April to November 2002 (summing to a total of 111 selected days) and used least squares adjustment.The authors reported error estimates to an accuracy of about 2 cm over the land and ocean regions.An improvement on the geopotential model GGM01, called GGM02, was released in 2005.GGM02 exists both in the GRACE-based satellite-only GGM02S and the combined model GGM02C. 27The combined geopotential model incorporated the GRACE-only model GGM02S with EGM96 plus 14 months of GRACE data spanning from April 2002 to December 2003.The GGM02C model was computed to maximum degree and order 200 in terms of spherical harmonics.Improvements by a factor of two were reported with error estimates of less than 1 cm geoid height to spherical harmonic at degree 70.
The satellite-only model EIGEN-GL04S1, described by Foerste et al. 28 , has a maximum degree and order of 150.It incorporated GRACE-only (EIGEN-GRACE04S) and GRACE-LAGEOS (EIGEN-GL04S) solutions.EIGEN-GL04S1 was later combined with surface gravity data from altimetry over the oceans and gravimetry over the continents to derive a high-resolution gravity model, EIGEN-GL04C, released in 2006. 28This combined gravity field model is an outcome of the joint gravity field processing between GRGS Toulouse and GFZ Potsdam.The satellite part of EIGEN-GL04C is based on GRACE and LAGEOS data and the maximum degree and order of this model is 360 in terms of spherical harmonics.
EIGEN-5C 29 was also a joint collaboration between GFZ Potsdam and GRGS Toulouse.EIGEN-5C is an upgrade of EIGEN-GL04C and has a maximum degree and order of 360.
The model is again a combination of GRACE and LAGEOS tracking data combined with gravimetry and altimeter surface data.The combination of the satellite and surface data has been done by combining normal equations obtained from observation equations for the spherical harmonic coefficients.The National Geospatial-Intelligence Agency released the first ever global model capable of resolving the Earth's gravity field beyond a spherical harmonic degree of 2000; this model is called EGM2008. 30A description of this model can be found in Pavlis et al. 30 The EGM2008 gravity field model has a maximum degree and order of 2159.It incorporates improved gravity anomaly data, altimetryderived gravity anomalies and GRACE-based satellite solutions.It allows proper computation of quasi-geoid heights, gravity anomalies and vertical deflections and has a spatial resolution of ~5 arc minutes or ~9 km in the latitudinal direction. 30

Improvements in gravity field models
Improvement in gravity field modelling in terms of accuracy and spatial resolution is necessary in order to understand the physics of the interior of the Earth; the dynamics of the ocean and the interaction of continents; to study the sea levels of ice and oceans; and to better determine satellite orbits and height systems in science and engineering. 31uch improvements are expected owing to the availability of qualitative data, especially from the low Earth-orbiting satellites.Satellite missions such as CHAMP, GRACE and GOCE (launched in 2000, 2002 and 2009, respectively), are believed to have improved the spatial resolution sensitivity and accuracy of the newly developed GGMs. 32he CHAMP, GRACE and GOCE missions were designed to resolve the long-wavelength part of the gravity field and hence provide unprecedented accuracy. 32In contrast to the sporadic tracking by the SLR network of the ILRS, the three satellite missions carry GPS receivers on board that allow continuous orbit tracking.Furthermore, these satellites are equipped with accelerometers which provide direct measurements of the non-conservative forces (e.g.air drag).
In the case of GOCE, six accelerometers are installed in a gradiometer arrangement which additionally allows for direct measurement of the Earth's gravity gradients, which gives an improvement in the medium-wavelength part of the gravity.
Improvement in gravity field modelling has already been noticed in several models as the resolution in such models is increased to reach higher degree and order in terms of spherical harmonics.Such improvements can be measured by studying characteristics of the GGMs based on several factors.For example, the behaviour of GGMs can be analysed by performing orbit adjustment tests on artificial satellites and GPS or levelling tests, and by comparing the spectral behaviour of the models or ocean geoid. 33Whilst old geopotential models derived up to degree and order 70 can resolve spatial features (geoid computation) at a half wavelength of about 290 km, models (particularly the most recent) computed up to degree and order 360 can resolve spatial features down to 55 km. 34Early evaluations of gravity field models by Zhang and Featherstone 35 reported that the OSU91A geopotential model provided better fits to the gravity field over the Australian region than did prior models.
In contributions by Pearse and KearsIey 36 and Kirby et al. 37 , the OSU91A model has been ousted by the EGM96 model, which reportedly gives better solutions to the computation of geoid heights.
Evaluations of GGMs released between 1996 and 2002 by Amos and Featherstone 12 , based on comparisons of gravity anomalies, free-air gravity anomalies, geoid heights and GPS or levelling tests, found that EIGEN-1S was the best satelliteonly GGM when applied in the Australian and New Zealand region, whilst the best combined GGM over the same region was reported to be PGM2000A.The quality of the GGM01 model was assessed by Ellmann 38 by comparing it with the combined gravity field model EGM96.It was reported that the GGM01 model gives better solutions of gravity anomaly and geoidal heights over Fennoscandia (e.g.Finland, Germany, Norway and Sweden) and the Baltic Sea region.
In a comparison study of 10 geopotential models (EGM96, GGM02C, GGM03S, ITG-GRACE03, JEM01-RL03B, EIGEN-GL04C, EIGEN-5C/5S and EGM2008) using geoid heights and GPS or levelling data points, the EGM2008 model was reported to provide the best solution compared to the other models at degree 360. 33A much improved solution was also reported for EGM2008 when its coefficients were increased to degree 2190.A similar study evaluating the GGMs EGM96, EIGEN-5C and EGM2008, based on the comparison of geoid heights to the GPS or levelling over Afyonkarahisar in Western Turkey, has also confirmed the improvements in the EGM2008 model in the computation of geoid heights. 39Botai and Combrinck 40 investigated the general improvement of 13 gravity field models released between 1990 and 2008 based on LAGEOS data and they found that gravity field modelling had improved during this period by a factor of two (Figure 4).In this review, the GRACE-only gravity field model, AIUB-GRACE01S, has been shown to provide lower root mean square orbital errors compared to the other tested models.

Geophysical applications of gravity field modelling
According to Newton's law, changes in the gravity field are evidence of changes in mass and density distribution in the Earth system.Any movement of masses in, on or above the Earth will therefore introduce variations in the gravity field of the Earth.often in the order of 10 -6 N/kg (variation from the mean) and occur on a scale ranging from hours to thousands of years. 42uch temporal variations are caused by several phenomena that redistribute mass, which include tides caused by the Sun and Moon, and postglacial rebound caused by isostatic correction.Surface mass changes in the atmosphere, oceans, hydrosphere and cryosphere are dominated by seasonal and interannual variations, whilst processes such as isostatic glacial recovery and sea-level change give rise to long-term secular or quasi-secular signatures. 43veral studies have investigated the long-term and the seasonal variations of the Earth's gravity field using data collected from different satellite missions. 44,45,46In particular, the lower order harmonic component of the gravity field with n = 2 and m = 0 (hereafter J 2 ), which characterises the gravitational oblateness of the Earth has attracted a lot of interest from the scientific community.Early studies of J 2 , for example, that by Yoder et al. 47 , showed a secular increase in J 2 that was consistent with a steady migration of mass from low latitudes towards high latitudes, resulting in a linearly decreasing trend.Such a trend was thought to be related to postglacial rebound, the Earth's ongoing response to the removal of the ice loads at the end of the last Ice Age.However, long-term studies by Cox and Chao 48 have discovered that J 2 started to increase from about 1997.
The detected increasing trend reported in Cox and Chao 48 is believed to have turned again, with J 2 once more decreasing from late 1997.Several mechanisms have been suggested to be the cause in the sudden changes of the J 2 coefficient.These mechanisms include processes involved in the surge in subpolar glacial melting and mass shifts in the Southern, Pacific and Indian Oceans. 1 In addition to the increasing trend of the J 2 coefficient, Nerem et al. 49 reported that the J 2 coefficient might be exhibiting seasonal variability as a result of a combination of atmospheric pressure variations and variations in the distribution of water in the oceans and on land.Furthermore, Dickey et al. 2 detected interannual variability in the J 2 which they attributed to climatically driven oscillations in the ocean, storage of water, snow and ice on land, and also partly to the consequence of the effects of anelasticity on the 18.6-year solid Earth tide, as suggested by Benjamin et al. 50rth orientation changes, often represented by polar motion, X-equatorial and Y-equatorial polar components in a geographical reference frame and variations in the length of day (LOD), are often explained by studying variations of atmospheric and/or oceanic angular momentum.Such variations are caused by the exchange of angular momentum between the solid Earth and its geophysical fluid envelope.Eubanks 51 found that variations in the Earth's rotation rate corresponded to changes in LOD and amounted to a few parts in 10 8 .Studies by Ponsar et al. 52 suggested that variations in LOD are caused by the interaction between the Earth's core and mantle.Similar studies by Gross et al. 53 related the LOD variations with tidal variations, exhibiting periods between 12 h and 18.6 years.Such variations were believed to be as a result of the deformation of solid Earth and changes in the strength and direction of the winds.
Geopotential models can be used to derive the geoid, the equipotential surface of the Earth's gravity field that corresponds closely with mean sea level in the open oceans, ignoring oceanographic effects as well as the geoidal height which is the separation between the geoid and the ellipsoid. 54etermination of the geoid has been one of the main research areas in Geodesy for decades, because of its various applications, which include vertical data for orthometric heights, understanding of ocean circulation patterns and dynamics, refinement of satellite orbits and the modelling of geodynamical phenomena.To this end, geoid heights at any point on the Earth's surface can be determined with an accuracy ranging from 30 cm to a few metres. 55A number of researchers have addressed the precise determination of geoid height on local and regional scales for oceanographic and geophysical applications. 56,57At a local scale, the geoid can be determined by a combination of GPS-derived heights and levelled heights, through gravimetric and geometric approaches.For instance, the quasi-geoid for southern Africa based on SLR-derived geopotential gravity models has been reported by Merry 58 .In Merry's 58 study, gravity data for South Africa was combined with different geopotential models (based on the remove-restore procedure) to derive a quasigeoid model, UCT2006.In addition, quasi-geoids produced in South Africa by use of different GGMs were compared with GPS or levelling data points to assess the suitability and reliability of the considered models.Merry 58 concluded that the UCT2006 model gives a better solution (with a root mean square fit of 15 cm) compared to the EGM96, EIGEN-CG03C and GGM02C models.A slight improvement of 4 cm was also reported when the UCT2006 geoid model was allowed to tilt in two directions (north-south and east-west).
Global gravity change has also attracted particular attention in the scientific community as it is often related to global sea-level changes. 59,60The sources of global sea-level rise often involve the redistribution of mass from the continents to the ocean.The use of gravity field measurements allows the discrimination of several sources through the continual monitoring of geoid changes on both global and regional scales as well as on basin scales.Gravity field solutions can be used to numerically estimate components such as thermal expansion (eustatic) and freshwater influx influencing global sea level.Measurements of temporal gravity variations can also be used to determine water storage change in the hydrological system.Since the launch of the GRACE mission in 2002, numerous articles assessing the potential of GRACE recovering hydrological signals have been published.For example, Andersen and Hinderer 61 have investigated the potential of inferring interannual gravity field changes caused by continental water storage change, as determined from GRACE observations between 2002 and 2003.
Contributions from continental water storage change were compared to the output from global hydrological models.Andersen et al. 62 and Neumeyer et al. 63 correlated large scale hydrological events with the estimated change in the gravity field for certain areas of the world to an accuracy of 0.4 µGa corresponding to 9 mm of water.On a regional scale, Winsemius et al. 64 compared hydrological model outputs for the Zambezi river basin with estimates derived from GRACE.Monthly storage depths produced by the hydrological model displayed larger amplitudes and were partly out of phase compared to the estimates based on GRACE data.

Summary
The continuous design and deployment of satellite missions dedicated to gravity field measurements and the availability of high-precision data have led to the availability of gravity information with unprecedented spatial-temporal resolution and accuracy.In particular, the advent of satellite data has made it possible to determine the gravity field of the Earth via modelling.To this end, these data sets are the basis of robust gravity field modelling with more than 100 gravity models being released to the scientific community since the early 1960s.Research dedicated to assessing the accuracy of the existing gravity field models has been reported in the literature.Different gravity field models could be characterised by various degrees spatial-temporal resolution.Such improvements are as a result of the availability of quantitative and qualitative SLR and terrestrial gravity data.With the development of gravity field models, dedicated review papers that report on the chronology of gravity field modelling for geophysical applications have been lacking.Review papers known thus far are more than a decade old and therefore cannot provide a complete account of up-to-date gravity field modelling activities.This review has explored the various gravity field modelling efforts with specific geophysical applications.It is concluded that gravity field modelling algorithms have improved over time partly due to the availability of specialised gravity mission satellite data and partly because of the advancement of technology and Earth and orbit system modelling techniques or approaches.

FIGURE 1 :
FIGURE 1:The satellite LAGEOS 1 contains an array of retroreflective mirrors covering its surface.

FIGURE 2 :
FIGURE 2: The global distribution of the International Laser Ranging Service (ILRS) tracking network.

FIGURE 3 :
FIGURE 3: The Satellite Laser Ranging station MOBLAS-6 is located at the Hartebeesthoek Radio Astronomy Observatory in South Africa.

TABLE 1 :
Timeline of artificial satellites that were tracked by global Satellite Laser Ranging stations.

well as the other services of the International Association of Geodesy. This SLR tracking station is relatively isolated in Africa and more active than Helwan; hence HartRAO plays a very important role as far as data coverage is concerned. Forces acting on an orbiting satellite
inertial reference frame perturbed by gravitational, non-gravitational and unmodelled forces is often expressed by [Eqn 3], as reported in Seeber 5 : [Eqn 3]

TABLE 2 :
Summary of some of the global geopotential models released between 1990 and 2008.Geophysical applications of these models include gravity field, satellite orbit determination, station coordinates, reduction of altimeter data, Earth's rotation and computation of geoid undulations.
2,41Temporal variations in Earth's gravity field are