WP 3.3.5. Influence on Ice-Mass Balance

Gabriel Strykowski, Geodynamics Division, KMS

1. Introduction

The description of this work package, WP 3.3.5., in the original technical proposal is very short and very general. For practical reasons it was necessary to narrow the scope of the investigation. This was done in corporation with prof. Carl Christian Tscherning. In this work package we are concerned with the question of whether the improved accuracy of the estimated Spherical Harmonic Coefficients, SHC, by future satellite missions (e.g. GOCE), can ensure a better aerial estimation of depths to the bedrock below the ice, in particular below the Greenland Inland Ice. Gravity information obtained from the satellite missions is viewed as additional, but indirect, source of information to e.g. depths obtained from the airborne radar echosoundings (which were conducted in Greenland in late 1970's and early 1980's). The relation to the ice-mass balance is indirect; a better knowledge of the bedrock topography below the ice results in a better glaciological model for the ice-mass balance. To fulfill the work package, it was agreed to make an error propagation study. The investigated area covers the whole of Greenland (55oN - 85oN and 0oW - 90oW).

2. The error propagation study

The idea is simple. In a local area, the error in the gravity field masks the information about parts of the gravity signal below the noise level. If the noise is spatially uncorrelated and centred (i.e. the mean value for each point in space is zero), it is still possible (e.g. by Wiener filtering or even by a simple moving window averaging technique), from noisy data, to recover gravity features of certain spatial extension. This statement is valid even if the gravity signal associated with such features is below the noise level for each point in space. Thus, the "visibility" of certain features of the signal masked by the noise depends not only on the signal to noise ratio for a particular point in space, but also on the relation between the spatial extension of the signal feature and the correlation properties of the noise.

In physical geodesy it is customary to view the local gravity signal as a realization of a 2nd order stochastic process with certain, idealized properties (ergodicity, zero mean, azimuthal isotropy). From measured and centred data an empirical covariance function can be estimated and, if necessary, fitted by a model function C(s) where s denotes a spherical distance. The model function yields a variance C0 and a correlation length s1 (i.e. spherical distance for which C(s1 ) = C0). A similar analysis can be conducted for the noise e, yielding a noise variance C(e) and a noise correlation length s1(e).

How to interpret the masking property of the correlated additive noise? A detailed answer can be quite complicated. The theoretical background from Signal Analysis is the theory of matched filters and the theory of Wiener filters with correlated noise (see e.g. Papoulis, 1984), except that in the present case the whole power spectrum of the correlated additive noise is not known. (In principle, a noise characterization by its variance and its correlation length is weak with respect to constraining the power spectrum. In other words, different covariance models, each with a distinctive power spectrum, can fit such weak noise characteristics. The design of an optimal filter in the above methods depends on the particular power spectrum.)

Consequently, instead of the cumbersome theoretical derivations, one can circumvent the question by asking; at what level will the signal features be visible above the noise level of the additive correlated noise? It was decided to answer this question empirically. The procedure is as follows:

(1) Computing realizations of noise. The computer program HARMEXG.F (see WP 2.1) is used to generate 20 realizations of noise caused by the errors in SHC, and for the area of investigation (see Sec. 1). In details, SHC of EGM96 have been perturbed 20 times according to the expected SHC error statistics of the GOCE mission. This yields 20 realizations of the field contaminated by the noise in different height levels H. (In practise, the direct evaluation of the noisy gravity field in different height levels was chosen as an alternative to the harmonic downward continuation of the noisy gravity data from the flight level, H = 300 km). Subsequently, the reference gravity field (EGM96) was computed at each height level and subtracted from the 20 realizations of the error contaminated gravity field. For the area of investigation, and for each height level H, this yields 20 realizations of the noise according to the expected SHC error statistics of the GOCE mission. Fig. 1 shows the reference gravity signal (EGM96) for H=3000m. Fig. 2 shows the same gravity signal contaminated by the SHC noise according to the apriori error statistics of the GOCE mission. Fig. 3 shows the corresponding noise.

(2) Empirical covariance properties. The noise realizations described in (1) are given on a regular grid in geographical coordinates (spacing: dlat x dlon= 0.25o x 0.50o). These data can be used to estimate the empirical covariance functions as it is commonly done in physical geodesy, e.g. in GRAVSOFT software package [Tscherning et al., 1992]. The gravity error covariance properties are expressed as variance C0(eg) and the correlation length s1(eg). Another problem is to transform these expressions to the formal error covariance properties for the undulations of the topography. The technique is directly related to the notion of the linear approximation for topographic effects, see e.g. [Forsberg, 1984, Sec. 7.3]. It is assumed, that the transition from the ice to the bedrock in Greenland is associated with a fixed mass density contrast drho=(2670 kg/m3 - 920 kg/m3] = 1750 kg/m3. In the spherical case, the linear approximation yields:

dgH = (4 pi G drho) dhH (1)

where G is the gravitational constant, dgH denotes a harmonic downward continued gravity signal from the satellite height to the height level H (above or below the zero-level), dhH denotes the undulation of the rock-ice interface just below the level H.

Eq.(1) can be viewed as a linear functional relating dgH and dhH. Thus, the law of covariance propagation applies (Moritz, 1980). (It is important to emphasize, that we do not postulate that the sources of the gravity field are located just below the height level H. Eq. (1) will only be used for the purpose of covariance propagation, i.e. in order to determine the "gravi- equivalent" covariance properties in different depths for the undulation error.) Using the functional expressed by eq. (1) the variance and the correlation lengths of the gravity error can be directly translated to the equivalent error of the interface undulation in some depth H:

C0(eg ) = (4 pi G drho)2 C0(eh ) (2a)

and the correlation length

s1(eg) = s1(eh ) = s1 (2b)

The following investigations of the noise covariance properties (as suggested by prof. C.C.Tscherning) were followed in this investigation:

Latitude bands of width =3o were chosen.

(3) Visibility of the undulations of the ice-rock interface. As explained above, it is not easy to relate the formal undulation error parameters C0(eh) and s1(eh), obtained from eqs. (2a)-(2b), to the equivalent masking property of the correlated noise. Instead, a purely empirical method was used. A small subarea of the investigated area has been chosen. A rectangular initial prism of height h0 = [C0(eh)] and a horizontal width s0=dx0=dy0=s1 (both in EW- and NS-directions) has been chosen. The centre of mass of the prism coincides (in horizontal) with the location of the gravity grid point in the centre of the chosen subarea. In the vertical, the bottom side of the prism is located 100 m below the height level H. The mass density contrast to the surroundings is, drho = 1750 kg/m3, thus representing the undulation of the rock topography below the ice in the height level H. Subsequently, the thickness (i.e. the height) of the prism is changed in steps dh= 1 m, so that h = h0 + Ndh, where N is a whole number (positive/negative), until the gravity response from the prism is clearly seen above the gravity noise (and for all 20 realizations of the noise). The obtained value of h is denoted hlim.

Starting with the prism with parameters (hlim, s0, s0) the horizontal extension of the prism are changed gradually and simultanously in both directions. The steps are ds = 0.1 km, until the gravitational attraction from the prism is no more visible above the noise level. Thus, for each depth level H the specific parameter slim is obtained. The results of the investigation show, that the question of whether certain features are visible depends mainly on h. For h<hlim the prism response is not visible above the noise level, irrespectible of the horizontal extension. However, for h > hlim, the features are clearly visible, even for features with relatively small horizontal extensions as compared to the correlation length, see Table 4. The precise meaning of the parameters hlim and slim is, thus, that for h>hlim and s>slim the undulations of the ice-rock interface with spatial resolution (h,s,s) will be visible in the gravity signal.

 Fig. 4 shows the gravity response in height H = 0 m of the rectangular prism (h0,dx0, dy0) = (11.4 m, 75 km, 75 km) with the bottom side located in a height level H=-100 m. Fig. 5 shows the noise in the same subarea. Fig. 6 shows the gravity response of the prism shown on Fig. 4 with the background noise. Clearly, the prism gravity response is hidden in the noise. Fig. 7 shows the gravity response of a prism (hlim,dx0, dy0) and the background noise and Fig. 8 shows the gravity response of a prism (hlim,slim, slim) and the background noise. As explained above, Fig. 7 and Fig. 8 demonstrate the concept of the visibility of the undulations of the ice-rock interface in the presence of the additive correlated noise.

 3. The results

In Sec. 2. the background for the results presented in the tables below is explained in some details. Furthermore, the computations were repeated for the error statistics of EGM96 to illustrate the improvement. The results shown in Table 1 indicate an improvement by a factor 4.4 in the formal standard deviation of the covariance function. However, Table 4 shows that the improvement in "visibility" (see Sec. 2) yields a factor 6.4 for hlim and a factor 7-10 for slim.

Table 1. The whole area (55oN - 85oN and 0oW - 90oW). Covariance properties of noise for different heights.

H (metres)

Gravity noise, standard deviation (mgal) correlation length, s1 (km) Rock-ice interface undulation noise, standard deviation (metres)

300000 0.0085 0.1832 135 250
3000 1.55 6.93 85 60
200 1.68 7.43 75 60 11.42 50.64
100 1.67 7.43 75 60 11.41 50.64
0 1.67 7.43 75 60 11.40 50.64
-100 1.68 7.45 75 60 11.43 50.75
-200 1.68 7.47 75 60 11.46 50.88
-300 1.69 7.48 75 60 11.49 50.99
-400 1.69 7.50 75 60 11.52 51.11
-500 1.69 7.52 75 60 11.55 51.23

The undulations of the ice-rock interface in Greenland are located in height levels between some -500m and 200m.The height level H = 300000 m is that of the GOCE satellite, and H = 3000 m is the height just above the highest point of the Inland Ice. Computation of the formal equivalent height undulation error statistics is not relevant for these two heights.

Table 1. shows the formal covariance parameters (here the standard deviations) of the error empirical covariance function. For the height levels between -500m and 200m there is a slight but systematic change of the variance (i.e. the standard deviation) and the corresponding change in variance (standard deviation) of the ice-rock interface undulation error. The apparently constant correlation length should be understood in the following way: to determine values of the empirical covariance function the data were grouped in intervals of width =20km. Subsequently, the value of the correlation length is determined from the noisy empirical covariance estimates.

Table 1 show, that, for the area as a whole, and except for the weak change of variance (or equivalently the standard deviation), the average covariance properties do not change significantly in height levels between -500 m and 200 m. However, as shown in Table 2, there is a clear and significant change of covariance properties with latitude. There is a clear increase of error variance with latitude, and a small decrease of noise variance with height.

Table 2. Change of covariance properties with latitude (H=-500 m and H=200 m)

H (metres)

Gravity noise, standard deviation (mgal) correlation length s1 (km) Rock-ice interface undulation noise, standard deviation (metres)
55oN - 58oN -500 1.18 7.18 85 70 8.04 48.49
55oN - 58oN 200 1.16 7.13 85 70 7.92 48.68
58oN - 61oN -500 1.18 7.31 85 65 8.04 49.85
58oN - 61oN 200 1.17 7.26 70 70 7.97 49.45
61oN - 64oN -500 1.31 7.32 70 60 8.93 49.91
61oN - 64oN 200 1.30 7.27 70 65 8.87 49.60
64oN - 67oN -500 1.37 7.46 80 60 9.34 50.86
64oN - 67oN 200 1.36 7.41 80 65 9.24 50.34
67oN - 70oN -500 1.44 7.43 70 65 9.81 50.61
67oN - 70oN 200 1.42 7.38 75 50 9.71 50.46
70oN - 73oN -500 1.56 7.37 80 60 10.63 50.24
70oN - 73oN 200 1.54 7.32 80 60 10.49 49.86
73oN - 76oN -500 1.71 7.62 85 65 11.65 51.91
73oN - 76oN 200 1.69 7.57 80 60 11.51 51.56
76oN - 79oN -500 1.92 7.75 80 60 13.08 52.78
76oN - 79oN 200 1.90 7.69 80 60 12.93 52.33
79oN - 82oN -500 2.09 7.85 85 65 14.24 53.46
79oN - 82oN 200 2.07 7.79 85 60 14.14 53.21
82oN - 85oN -500 2.54 7.82 85 65 17.31 53.28
82oN - 85oN 200 2.51 7.76 85 60 17.09 52.84

Table 3 below shows the covariance properties for the whole area and in one height level (H = 0 m), but for different realizations of noise. The idea is to demonstrate how robust the estimated error covariance properties are.

Table 3. Error covariance properties for different realizations of noise, H = 0 m.

Perturbation number

Gravity noise, standard deviation (mgal) correlation length, s1 (km) Rock-ice interface undulation noise, standard deviation (metres)
1 1.54 7.29 85 65 10.49 49.66
2 1.75 7.20 85 65 11.91 49.00
3 1.48 7.35 85 65 10.11 50.21
4 1.65 7.52 85 60 11.21 51.09
5 1.59 7.41 90 65 10.82 50.43
6 1.60 7.22 85 65 10.93 49.32
7 1.69 7.77 85 65 11.51 52.92
8 1.53 6.95 85 60 10.41 47.29
9 1.79 7.17 90 65 12.19 48.83
10 1.71 7.35 80 60 11.66 50.12
11 1.72 7.39 80 65 11.70 50.27
12 1.50 7.54 85 70 10.20 51.27
13 1.86 7.60 90 65 12.70 51.89
14 1.81 7.95 90 70 12.32 54.11
15 1.79 7.29 85 60 12.21 49.36
16 1.60 7.75 85 70 10.88 52.70
17 1.57 7.42 85 65 10.71 50.62
18 1.89 7.64 90 70 12.88 52.07
19 1.60 7.42 90 65 10.88 50.46
20 1.72 7.32 90 65 11.69 49.75

 Finally, Table 4. shows the visibility parameters hlim and slim (see Sec. 2) in different heights. It should be mentioned, that a decision about which parameter to choose is somehow subjective. The difference in the effect of changing parameters of the prism is not sharp and depends on the chosen location (a grid point). For this reason, the results displayed in Table 4. should be used with great caution. For comparison, the excersise was repeated for EGM96 error statistics. These results were obtained using a rougher resolution for the parameters (+/-0.5km for slim and +/-5m for hlim).

Table. 4. The parameters slim and hlim related to the concept of the visibility of the undulations of the ice-rock interface, see Sec. 2.
H (metres) slim (km) hlim (metres)
-500 0.2 1.5 43 270
-400 0.2 1.5 43 270
-300 0.2 1.5 42 265
-200 0.2 1.5 41 260
-100 0.3 2.0 40 250
0 0.4 3.0 37 240
100 0.6 5.0 34 220
200 0.8 8.0 32 205

 4. Figures

Figure 1. Greenland, reference gravity anomalies (EGM96, H=3000m).

Figure 2. Greenland, reference gravity anomalies contaminated by a realization of the expected SHC-noise of GOCE mission (EGM96, H=3000).

Figure 3. Greenland, one realization of the expected SHC-noise of GOCE mission (EGM96, H=3000).

Figure 4. Greenland, gravitational attraction of a prism of size (dx, dy, dh) = (75 km, 75km, 11.4 m), drho = 1750 kg/m3

Figure 5. Greenland, one realization of the expected SHC-noise of GOCE mission, see Figure 3, H=0 m.

Figure 6. Greenland, gravitational attraction of a prism shown on Figure 4 and the additive correlated noise shown on Figure 5, H = 0 m.

Notice that although the prism is (almost) visible on the figure the criterion of visibility is that it is also visible in other realizations of noise.

Figure 7. Greenland, as in Figure 6 except that the size of the prism is changed to (dx, dy, hlim) = (75 km, 75 km, 37 m), H= 0 m. Notice that the

gravitational attraction of the prism is clearly seen above the correlated additive noise.

Figure 8. Greenland, as in Figure 7 except that the size of the prism is changed to (slim, slim, hlim) = (400 m, 400 m, 37 m), H= 0 m.


Forsberg, R.: A Study of Terrain Reductions, Density Anomalies and Geophysical Inversion Methods in Gravity Field Modelling. The Ohio State University, Dept. of Geodetic Science and Surveying, Rep. No. 355, 1984.

Moritz, H.: Advanced Physical Geodesy. Herbert Wichmann Verlag, Karlsruhe, 1980.

Papoulis, A.: Signal Analysis. McGraw-Hill Electrical and Electronic Engineering Series, 1984.

Tscherning, C.C., R. Forsberg, and P. Knudsen: The GRAVSOFT Package for Geoid Determination. In: Holota P. and M. Vermeer (Eds.) Proc. IAG, First Continental Workshop on the Geoid in Europe. ISBN 80-901319-2-1, pp. 327-334, Prague, Czechoslovakia, May 11-14,1992. 


Last update 1999-05-05 by cct.