Flechtner F.M., Gruber Th., G?ntner A., Mandea M., Rothacher M., Sch?ne T., Wickert J. (Eds.) System Earth via Geodetic-Geophysical Space Techniques

Подождите немного. Документ загружается.

GOCE Data Analysis: From Calibrated Measurements to the Global Earth Gravity Field 229

Kaula W (1966) Theory of Satellite Geodesy, Blaisdell, Waltham, MA.

Klees R, Ditmar P, Broersen P (2003) How t o handle colored noise in large least-squares problems.

J. Geod. 76, 629–640.

Koch KR, Kusche J (2002) Regularization of the geopotential determination from satellite data by

variance components. J. Geod. 76, 259–268.

Kusche J, Ilk KH (2000) The polar gap problem. In: Sünkel H (ed.), From Eötvös to Milligal, Final

Report. ESA/ESTEC Contract No. 13392/98/NL/GD.

Lemoine FG, Kenyon SC, Factor JK, Trimmer RG, Pavlis NK, Chinn DS, Cox CM, Klosko

SM, Luthcke SB, Torrence MH, et al. (1998) The Development of the Joint NASA GSFC

and the National Imagery and Mapping Agency (NIMA) Geopotential Model EGM96.

NASA/TP-1998-206861, Greenbelt.

Metzler B, Pail R (2005) GOCE data processing: The spherical cap regularization approach. Stud.

Geophys. Geod. 49, 441–462.

Miggliacchio F, Reguzzoni M, Sansó F, Tscherning CC, Veicherts M (2007) The latest tests of the

space-wise approach for GOCE data analysis. Proceedings of the 3rd International GOCE User

Workshop in Frascati, Italy.

Plank G (2004) Numerical solution strategies for the GOCE mission by using cluster technologies.

PhD thesis. Institut für Navigation und Satellitengeodäsie, TU Graz.

Reigber C (1969) Zur Bestimmung des Gravitationsfeldes der Erde aus Satellitenbeobachtungen.

Deutsche Geodätische Kommission, Reihe C, No. 137, Munich.

Schuh WD (1995) SST/SGG tailored numerical solution strategies. CIGAR III/Phase 2, Final

Report, Part 2.

Schuh WD (1996) Tailored numerical solution strategies for the global determination of the earth’s

gravity ﬁeld. Mitteilungen der Geodätischen Institute der TU Graz, No. 81, Graz.

Schuh WD (2002) Improved modeling of sgg-data sets by advanced ﬁlter strategies. In: Sünkel H

(ed.), From Eötvös to mGal+, Final Report. ESA/ESTEC Contract No. 14287/00/NL/DC.

Schuh WD (2003) The processing of band-limited measurements – Filtering techniques in the least

squares context and in the presence of data gaps. Space Sci. Rev. 108(1–2), 67–78.

Schuh WD, Kargoll B (2004) The numerical treatment of the downward continuation problem

for the gravity potential. In: Sansó F (ed.), V Hotine-Marussi Symposium on Mathematical

Geodesy, IAG Symposia 127, Springer, Berlin, pp. 22–31.

Schwarz HR (1970) Die Methode der konjugierten Gradienten in der Ausgleichsrechnung. ZfV

95, 130–140.

Siemes C (2008) Digital ﬁltering algorithms for decorrelation within large least squares problems.

PhD thesis. http://hss.ulb.uni-bonn.de/diss_online/landw_fak/2008/siemes_christian/

Sneeuw N, van Gelderen M (1997) The polar gap. In: Sansó F, Rummel R (eds.), Geodetic

Boundary Value Problems in View of the One Centimeter Geoid, Lecture Notes in Earth

Sciences, Vol. 65, Springer, Berlin, pp. 559–568.

Sneeuw N (2000) A semi-analytical approach to gravity ﬁeld analysis from satellite observations.

Deutsche Geodätische Kommission, Reihe C, No. 527, Munich.

Svehla D, Rothacher M (2002) Kinematic orbit determination of LEOs based on zero or dou-

ble difference algorithms using simulated and real SST GPS data. In: Adam J, Schwarz KP

(eds.), Vistas for Geodesy in the New Millennium, IAG Symposia 125, Springer, Heidelberg,

pp. 322–328.

Visser P (2008) GOCE gradiometer: Estimation of biases and scale factors of all six individual

accelerometers by precise orbit determination. J. Geod. 83, 69–85.

Wermuth M (2008) Gravity ﬁeld analysis from the satellite missions CHAMP and GOCE. PhD

thesis. TU Munich.

Yunck TP, Wu SC, Wu JT, Thornton CL (1990) Precise tracking of remote sensing satellites with

the global positioning system. IEEE Trans. Geosci. Remote Sens. 28(1), 108–116.

GOCE and Its Use for a High-Resolution

Global Gravity Combination Model

Richard Shako, Christoph Förste, Oleg Abrikosov, and Jürgen Kusche

1 Pre-GOCE Satellite-only Models

The determination of t he global Earth gravity ﬁeld is one of the major tasks of

geodesy. Until the end of the last century, a satisfactory solution for the Earth gravity

ﬁeld could only be obtained by combining data from a large number of satellites at

different altitudes and inclinations. The processing of these different data resulted

in the so-called “satellite-only” normal equation matrix (for instance cf. Biancale

et al., 2000).

Important milestones on the way to a more precise, satellite-based global

Earth gravity model were then set with the CHAMP and especially the GRACE

missions. By using data of these satellites alone, the Earth gravity ﬁeld could be

resolved with an accuracy of 1 cm for half-wavelengths down to approximately

270 km (e.g. EIGEN-GRACE02S; Reigber et al., 2005 or GGM02S; Tapley et al.,

2005). Furthermore, solutions had been derived which made use of combined

CHAMP/GRACE-data or combined GRACE/LAGEOS-data (e.g. EIGEN-GL04S1;

Förste et al., 2008b) which resulted in an even more accurate Earth gravity model.

Compared to the past, where a conglomerate of satellites contributed to the solution,

it is of great advantage that the long to medium wavelength features of the Earth

gravity ﬁeld are now homogeneously derived from only few satellites’ data. This is

one part of the milestone. The other part is the improved accuracy of the GRACE-

models of more than two orders of magnitude compared to the latest pre-CHAMP

satellite-only model GRIM5-S1 (Biancale et al., 2000).

Figure 1 illustrates the progress by means of some GFZ and GFZ/GRGS derived

satellite-only models plotted in their internal unregularized versions. In this ﬁgure,

the degree amplitudes between degrees 90 and 160 are displayed together with the

combination solution EIGEN-5C (Förste et al., 2008a), serving as a reference for a

realistic spectral behaviour of the satellite-only solutions. The degree variances of

R. Shako (B)

Helmholtz Centre Potsdam, GFZ German Research Centre for Geosciences,

Department 1: Geodesy and Remote Sensing, Telegrafenberg, 14473 Potsdam, Germany

e-mail: rst@gfz-potsdam.de

231

F. Flechtner et al. (eds.), System Earth via Geodetic-Geophysical Space Techniques,

Advanced Technologies in Earth Sciences, DOI 10.1007/978-3-642-10228-8_18,



Springer-Verlag Berlin Heidelberg 2010

232 R. Shako et al.

Fig. 1 Signal amplitudes per degree in terms of geoid heights for various satellite-only grav-

ity ﬁeld models derived at GFZ and GFZ/GRGS with the combination model EIGEN-5C as

reference

the satellite-only models show a run-off at higher degrees. This run-off and its unre-

alistic increase in power is an indication for the limited sensitivity of the satellite

data, and therefore a measure for the stability at higher degrees of a satellite-only

model. For the earliest version, EIGEN-GRACE02S, the run-off starts at about 125,

while for EIGEN-GL04S1 it starts at about 140. The version EIGEN-5S does not

show a run-off at all, but only a spike at the end of the spectrum, probably a trun-

cation error effect. The fact that from version to version the run-off is gradually

suppressed can be simply explained by the amount of included GRACE data, which

increased over the course of the time as well. The latest version, EIGEN-5S, is

the best and most accurate performing model, since its degree variances are most

close to the combination solution. While for EIGEN-GL04S1 the sensitivity limit

is around degree/order (d/o) 140, Fig. 1 gives no information about the sensitivity

limit of EIGEN-5S. It can be concluded that GRACE-based satellite-only gravity

models are stable beyond degree 150, if enough observations are accumulated.

2 GOCE and Satellite-only Models

The recently launched GOCE-mission with its gradiometry observations is to be

expected the next milestone in geopotential modeling.

For geosciences, the satellite-gravity-gradiometry (SGG) measurements are a

very new type of observations. Already in early studies (for instance CIGAR IV,

GOCE and Its Use for a High-Resolution Global Gravity Combination Model 233

1996, or Abrikosov et al., 2004) it was shown that the gravity gradients observed by

GOCE will allow for the determination of the geopotential with a resolution of 80–

100 km half-wavelength (corresponding to about d/o 200–250 in terms of spherical

harmonics). Thus GOCE is expected to increase the quality of Earth gravity ﬁeld

models signiﬁcantly in the short and medium wavelength ranges. Unfortunately this

is not true for the long wavelengths due to the characteristic gradiometer bandwidth

and the speciﬁc GOCE orbit conﬁguration. This problem, however, can be simply

tackled by combining GOCE with GRACE data (which are much better for the long

wavelengths).

Another issue is caused by the GOCE mission proﬁle. Since GOCE will have

two separate observation periods of only six months each, its data are not perfect

for the determination of a static gravity ﬁeld. Indeed, the data will be contaminated

by seasonal and high-frequency signals caused by mass variations of the oceans,

in the atmosphere and by hydrological effects. Neglecting these seasonal and high-

frequency effects will result in fundamental systematic distortions. Therefore, data

corrections by de-aliasing and reduction techniques are necessary to use the GOCE-

SGG data at their best (Stubenvoll, 2006). During the last years, GFZ worked

intensely on the development of software and strategies in order to de-alias and

reduce these annual and high-frequency signals in the GOCE-SGG data. It has

turned out that, again, GRACE will be of big beneﬁt: by using GRACE time-variable

models for the reduction of seasonal hydrology mass variations and by using proper

models for the high-frequency signals in the ocean and the atmosphere, the GOCE-

data can be cleaned in a pre-processing step such that they contain only information

of the static gravity ﬁeld. The GRACE-models seem to be more reliable than, for

instance, hydrology derived models (Petrovic et al., 2004; Schmidt et al., 2008).

Note that the use of GRACE models is possible, because the two missions overlap

in time.

In order to predict the quality of a future GOCE gravity model and to develop

software and processing strategies, a simulation study has been performed on behalf

of ESA by the GOCE High Performance Facility HPF (Rummel et al., 2004). Within

the GOCE-HPF processing system three different strategies for GOCE gravity ﬁeld

determination have been adapted, these are the timewise approach (Pail et al.,

2005), the spacewise approach (Migliaccio et al., 2004), and the direct approach

(Abrikosov and Schwintzer, 2004). In the following the focus is set on the direct

approach which is used at GFZ. The simulated GOCE observations used for the

study were the E2E data as provided by ESA within the GOCE-HPF community

(Catastini et al., 2007). They had been processed by the HPF teams in the con-

text of the Acceptance Review 3 (AR3) in the framework of the HPF-development,

and cover 58 days, representative for a nominal mission scenario. The satellite-to-

satellite-tracking (SST) and SGG data were generated with the EGM96 model to

d/o 200 and 360, respectively. In fact, the GOCE SST normal equations (up to d/o

150) and the SGG normal equations (up to d/o 200) both computed for AR3 were

used in this study. In order to close the polar gaps, these data were combined with

simulated surface gravity anomalies (using EGM96 (Lemoine et al., 1998) gravity

model at 30



×30



grid covering polar areas bounded by latitudes ±80.25

◦

). Random

234 R. Shako et al.

Fig. 2 Results of a study using simulated GOCE observations: (a) Spectrum of geoid-errors in

mm w.r.t. EGM96 (left); (b) Errors of harmonic coefﬁcients w.r.t. EGM96 (right)

noise with an rms of 2.5 mGal was added to these simulated values. Then gravity

anomalies computed from EGM96 coefﬁcients up to d/o 45 were subtracted from

these data, and normal equations for harmonics up to d/o 200 have been computed

using the latitude-dependent weighting and standard deviation of 2 mGal.

In this study, the GOCE mission objectives (1–2 cm accumulated geoid error at

degree 200) cannot be attained, because only two months of simulated data were

available. Figure 2a, b s how the study results in the spectral domain. The spec-

trum of the resulting geoid-errors (Fig. 2a) demonstrates that the individual errors

coincide almost perfectly with the formal errors. The mentioned problems at low

degrees (up to d/o 15) are visible (due to the fact that the study did not use any

GRACE data). The typical run-off (as known from Fig. 1) does not show up, which

emphasizes that GOCE will be able to give accurate information for wavelengths

beyond d/o 200. The accumulated error passes the 2 cm threshold at about d/o 150.

Figure 2b shows the errors of the individual harmonic coefﬁcients. The low order

coefﬁcients have about the same error as the high order ones. Also, the zonal coef-

ﬁcients are less accurate. Figure 2b backs up the conclusion that GOCE will give

accurate information for wavelengths beyond d/o 200.

Figure 3 shows the results in the spatial domain. The polar regions exhibit much

higher errors. The global rms is 5 cm, while the rms of the non-polar regions is

4 cm. This is a very good result, because only 58 days of data had been used. It can

be estimated that – once the complete period of GOCE observations is available –

the mission goal is reachable. Finally, it can be concluded that by proper treatment

of the GOCE-SGG data and by optimal combination with GRACE-SST and GOCE-

SST observations, GOCE will increase the resolution and quality of satellite-only

Earth gravity models signiﬁcantly. While the current pre-GOCE models have an

accuracy of about 4 cm at 200 km half-wavelength, the expected GOCE accuracy i s

about 1–2 cm at 100 km half-wavelength. Thus, GOCE will open the door to a new

age of geopotential modeling.

GOCE and Its Use for a High-Resolution Global Gravity Combination Model 235

Fig. 3 Results of a study with simulated GOCE observations: geoid-errors w.r.t. EGM96 in cm

3 GOCE and Global Gravity Field Combination Models

3.1 Surface Data

The focus of the preceding section was set on satellite-only models. In order to

enhance the resolution of the gravity model to degree and order 360 (and higher) the

shortest wavelengths must be determined from surface gravity data. In the following,

the GFZ strategy for the combination of GOCE with surface data is discussed.

The currently existing surface data sets (compiled from satellite altimetry, ship

and airborne gravimetry over the oceans, and airborne and terrestrial gravimetry

over land) provide – except for Antarctica – an almost complete global cover-

age. Due to inconsistencies between the various data sets and regionally varying

accuracies, the surface data do not contain precise long to medium wavelengths

gravity information. However, if properly preprocessed and thoroughly combined

with satellite-only normal equations, the resolution of the global model can be

extended down to at least 55 km half-wavelength.

For preprocessing the surface data, it is of great importance to: (a) unify the

different data sets in order to eliminate errors in the geodetic datum; (b) apply the

necessary correction terms (for instance free air correction, linear and quadratic

terms of the normal gravity ﬁeld gradient, ellipsoidal correction terms, atmospheric

correction); (c) choose a high quality sea surface topography model, if sea surface

heights derived from satellite altimetry are used.

After these steps the different surface data sets are condensed to regular gridded

block mean values (30



×30



grids). Then, these gridded data can be used for the

combination with satellite-only normal equations. Basically two approaches exist:

236 R. Shako et al.

1. Combine the satellite data with full normal equations of the surface data up

to a certain d/o of spherical harmonics; for the higher degrees the surface data

contribute in form of block-diagonal normal equations to the solution.

2. Combine the satellite-data with full normal equations of the surface

data only.

The reason for following two different approaches is the amount of computer

memory (RAM and hard disk) needed for a full normal equation up to a certain

degree. For instance, one full normal equation up to d/o 180 (32,761 unknowns)

has the size of 4.3 GB (if stored as triangle matrix), while one full normal equation

up to d/o 360 (130,321 unknowns) has the size of 67.9 GB. The computation time

for setting up these normal equations behaves similarly exponential. Due to the fact

that each surface data set is handled as an individual normal equation system, the

necessary computer memory is 10–20 times larger (depending on the number of

different data sets). This is quite enormous for high d/o-solutions. Until the recent

past there was just no other way than to use block-diagonal techniques for these

solutions, since hard disk memory and RAM was too limited. With the use of block-

diagonal techniques, the size of the normal equation is reduced drastically, that is to

0.008 GB (for d/o 180) and 0.063 GB (for d/o 360), so that the higher d/o could be

easily calculated.

However, in order to make use of block-diagonal techniques, the used surface

data set must fulﬁll some requirements (Gruber, 2000):

• They must be complete and of the same data type.

• They must be gridded regularly and in equatorial symmetry.

• The weighting of the data must be constant for a certain latitude.

• The weighting of the data must be equator-symmetric.

If all these conditions are fulﬁlled, the resulting normal equations fall apart

into blocks due to the orthogonality of trigonometric functions. In particular, all

correlations are zero for (Gruber, 2000):

• all C- and S-coefﬁcients (e.g., C

• all coefﬁcients of different order (e.g., C

• all coefﬁcients of the same order if one coefﬁcient is of even degree while the

other is of uneven degree (e.g., C

3.2 Combination Models Derived from Full and Block-Diagonal

Normal Equations

The block-diagonal techniques offer the opportunity to calculate a high-resolution

model, whereby correlations between coefﬁcients are not neglected (thus it is out-

performing the numerical integration approach). On the downside there are the

mentioned requirements for the input data, which for instance do not allow merging

GOCE and Its Use for a High-Resolution Global Gravity Combination Model 237

gravity anomalies with geoid undulation as input data. Also, an individual weight-

ing is not possible. Therefore it must be the goal to use f ull normal equation for as

high degrees as computer resources allow.

Until the GFZ/GRGS-model EIGEN-5C (complete up to d/o order 360) block-

diagonal techniques had been used. In the following the strategy for a block-

diagonal combination is explained, using EIGEN-5C (Förste et al., 2008a) as pars

pro toto. Of course, t he presented strategy is optimized for the used GRACE data.

Its basic principle, however, can be adapted to GOCE very simple.

First, the satellite-only normal equation system (consisting of LAGEOS (com-

plete to d/o 30) and GRACE (complete to d/o 150)) was added to the full surface

data normal equation system (complete to d/o 280) in such a way that the surface-

and satellite-based coefﬁcients up to degree 70 were kept separate (see Fig. 4). That

means the s urface data did not contribute to the solution for n, m = 1...70; there-

fore the mentioned surface data problems in the long wavelengths (low degrees)

are avoided. The relative weights for the individual surface data normal equations

were chosen in accordance with the standard deviations of their underlying terres-

trial data. Furthermore, the surface data used for the complete normal equations had

been ﬁltered for frequencies higher than d/o 260 in order to avoid truncation errors,

which means that they do not contribute to the solution for n, m = 261...280.

In the overlapping spectral band, where both satellite and surface data are con-

tributing to the solution, there is a degree range where t he information of both

data types is completely used (from n, m = 70...90), followed by a degree range

where the satellite data are kept separate but bound together with the surface data

by introducing degree-dependent weighted pseudo-observations

C/S

nm_Surf.

− C/S

nm_Sat.

= 0 with weight p = 1/(0.8696 · 10

−10

· e

(n−90)/24.6

)

so that a smooth transition in the degree amplitude spectrum within the overlapping

part of the surface and satellite normal equations is assured, while the information

Fig. 4 Details of the EIGEN-5C combination

238 R. Shako et al.

of the satellite-only gravity ﬁeld model is kept up to its maximum d/o 150. The opti-

mum degree dependence of the weights in the above formula was found empirically

by testing a dedicated range of the therein given numbers after evaluating the for-

mal errors of the spherical harmonics of previous satellite-only models. The ﬁnally

resulting combined normal equation system was solved by Cholesky decomposition.

Last, t he block diagonal system was solved separately and the resulting solution

was used to extend the spherical harmonic coefﬁcients from degree 261 to degree

359, and the degree 360 coefﬁcients (obtained through numerical integration) were

added for completeness.

3.3 The GOCE-Model: Combination with Full Normal

Equations Only

As mentioned in the previous subsection, the use of block-diagonal techniques is

less preferable than the use of full normal equations. Due to the enormous develop-

ment in computer technology, which led to the availability of huge and affordable

hard disk and RAM memory, the handling of a dozen of 70 GB sized normal equa-

tions is not a too big problem anymore. The problem of the long processing time

has vanished due to the much higher speed of computer processors. Furthermore,

the use of clusters offer an alternative for massive parallel processing.

Furthermore, the software existing at GFZ has been optimized signiﬁcantly.

Table 1 demonstrates the huge speed-up. While the old software version, which

used internal compiler parallelization, was not able (due to technical reasons caused

by the code) to set up normal equations up to d/o 250, the new mpi-parallelized

software can in a reasonable time frame. When using 15 processors on a SUN V890

workstation, the computation of a complete global set of normal equation up to d/o

359 takes about 2 days.

In the following, the up-to-date s trategy for the combination of GOCE-SGG

observations with GRACE, LAGEOS, and surface data sets is discussed. This

strategy is targeting at a geopotential model complete up to d/o 360.

Table 1 Computation times for setting up one normal equation for 3,800 observations on a SUN

V890 workstation using 6 processors

Maximum degree

Normal equation

size (GB as

triangle matrix)

Old software

version (s) New version (s)

Increase of speed

(new/old)

100 0.4 2,036 91 2,237%

200 6.5 34,950 1,287 2,716%

250 15.8 Technically

impossible

3,248 ∞

300 32.8 Technically

impossible

6,348 ∞

359 67.9 Technically

impossible

9,862 ∞

GOCE and Its Use for a High-Resolution Global Gravity Combination Model 239

Fig. 5 Combination scheme for GOCE-SGG data with GRACE, LAGEOS, and surface data

Figure 5 shows the planned combination scheme. On the left hand side the pro-

cessing of surface data is displayed. As discussed in Sect. 3.1, the surface data are

preprocessed, homogenized and gridded in 30



×30



-resolution. Following this, they

will be ﬁltered for degrees higher than the resolution of the targeted geopotential

model (d/o 360). This is necessary for avoiding truncation errors and aliasing effects

of higher frequencies. Then, complete normal equations are computed up to d/o 359

(one normal equation for each individual surface data set). Since a 30



×30



-grid can-

not provide an invertible normal equation for the coefﬁcients of degree 360, these

coefﬁcients are calculated by numerical integration.

On the right hand side of Fig. 5 the processing of GOCE-SGG observations

is displayed. After the preprocessing and the elimination of time variable signals,

the GOCE normal equations are obtained (complete up to d/o 250), which are

then combined with both the surface normal equations and the GRACE/LAGEOS

normal equations (see below for further details). An almost complete gravity

ﬁeld model up to d/o 359 will result. The ﬁnal step is to extend this model

with the results of the numerical integration in order to make it complete up to

d/o 360.

Figure 6 shows the details of the combination step in the current stage of plan-

ning. First of all, it does not seem necessary to include CHAMP- and GOCE-SST

data. Internal studies have shown that their contribution to the quality of the solu-

tion does not justify their use. LAGEOS data, however, are very advantageous for

the accurate estimation of the C

-coefﬁcient. In analogy to the combination of the

EIGEN-5C model, LAGEOS will contribute to the solution from d/o 2 up to d/o 30

with complete normal equations, and GRACE will contribute with complete normal

equations up to d/o 150. The information of GOCE and surface data is eliminated

up to d/o 30 and d/o 120, respectively. This is done by keeping their low degree