DETERMINISTICALLY-MODIFIED INTEGRAL ESTIMATORS OF GRAVITATIONAL TENSOR

The Earth’s global gravity field modelling is an important subject in Physical Geodesy. For this purpose different satellite gravimetry missions have been designed and launched. Satellite gravity gradiometry (SGG) is a technique to measure the second-order derivatives of the gravity field. The gravity field and steady state ocean circulation explorer (GOCE) is the first satellite mission which uses this technique and is dedicated to recover Earth’s gravity models (EGMs) up to medium wavelengths. The existing terrestrial gravimetric data and EGM scan be used for validation of the GOCE data prior to their use. In this research, the tensor of gravitation in the local north-oriented frame is generated using deterministicallymodified integral estimators involving terrestrial data and EGMs. The paper presents that the SGG data is assessable with an accuracy of 1-2 mE in Fennoscandia using a modified integral estimatorby the Molodensky method. A degree of modification of 100 and an integration cap size of 2.5 for integrating 5 5    terrestrial data are proper parameters for the estimator.


INTRODUCTION
Satellite gravity gradiometry (SGG) is a space method for recovering the Earth's gravity field from the second-order derivatives of the field.Delivering Earth's gravity models (EGMs) with higher resolution than those recovered from former satellite gravimetry techniques is expected from these data.The gravity field and steady-state ocean circulation explorer (GOCE) mission (BALMINO et al. 1998, 2001, ESA 1999, ALBERTELLA et al. 2002, RUMMEL et al. 2002, DRINKWATER et al. 2003) was successfully launched on 17 th March 2009 as the first satellite mission which used the SGG technique.So far GOCE has delivered EGMs todegree and order 250 corresponding to a resolution of 60 km by 60 km.The second-order derivatives of the gravity field is measured by a gradiometer mounted on the GOCE spacecraft and validation of them is very important prior to their use for any purpose.Having erroneous observations leads to the deviation of the recovered gravity field from the true one and wrong interpretation for the geophysical phenomena.
Numerous studies have been done in the Earth's gravity field modelling from SGG for years, but the idea of validating the SGG data was presented in a few different ways.In the following we divide the validation methods into two different categories: a) Validation of GOCE products Today, the concentrations of the geodetic researchers are on the evaluation of the GOCE products determined by three different methods of time-wise (TIM), spacewise (SPW) and direct (DIR); see Pail et al. (2011).Hirt et al. (2011) compared some of GOCE EGMs with terrestrial gravimetric data over Switzerland, Austria and astrogeodetic deflections over Europe.They observed some improvements between degrees 160-165 and 180-185. Gruber et al. (2011) compared some of the GOCE EGMs for reproducing the orbit of the Gravity recovery and Climate Experiment (GRACE) (Tapley et al. 2005) and concluded that they do not outperform the GRACE orbit; therefore, combination of GRACE and GOCE data is useful.They found significant improvements between degrees 50 and 200 for geoid computation goal.Janak and Pitonak (2011) evaluated the GOCE products in central Europe and Slovakia.They mentioned that TIM2 and SPW2 to degree 210 are much better than the previous releases to the same degree and GOCO02S (GOIGINGER et al. 2011) has a significant improvement comparing to GOCO01S (PAIL et al. 2010).Sprlak et al. (2012) did the same study in Norway and mentioned that the direct solutions are highly affected by a priori information and time-wise solution is more reliable.Abdalla et al. (2012) evaluated the GOCE EGMs in Sudan and concluded that the SPW1, SPW2, TIM1, TIM2 and GOCO01S are consistent with the local data.Abdalla and Tenzer (2012) validated EGMs in New Zealand.Guimaraes et al. (2012) tested the EGMs in Brazil and found out that TIM3 is much better than the previous ones, as expected.Eshagh and Ebadi (2013) also investigated different EGMs and evaluated them over Fennoscandia.So far evaluation of the GOCE EGMs was done based on comparison of the EGM products with external sources of data, like gravity anomaly, disturbing gravity, geoid and/or astrogeodetic deflections.However, it should be considered that the errors of EGMs have been also presented.Wanger and McAdoo (2012) noticed that the errors of GOCE EGMs are not realistic and tried to calibrate them based on EGM08 (Pavlis et al., 2008(Pavlis et al., , 2012)).Eshagh (2013) studied the reliability and calibration of GOCE EGMs as well.Eshagh and Ebadi (2014) presented a method for error calibration of some GOCE EGMs based on condition adjustment models.
b) On board validation of GOCE data Bouman et al. (2003) has set up a calibration model based on the instrument (gradiometer) characteristics to validate the SGG data.Bouman and Koop (2003) presented an along-track interpolation method to detect the outliers.Their idea is to compare the along-track interpolated gradients with measured gradients.If the interpolation error is small enough, the differences should be predicted reasonably by an error model.Also, Bouman et al. (2004) concluded that the method of validation using high-low satellite-to-satellite tracking data fail unless a high-resolution EGM is available.Kern and Haagmans (2004) and Kern et al. (2005) presented an algorithm for detecting the outliers in the SGG data in the time domain.

c) Validation of GOCE data by external sources of data
The simplest methodis the direct comparison of the observed SGG data with the generatedones using an existing EGM; (see ESHAGH and ABDOLLAHZADEH 2010, 2011). Also, Haagmans et al. (2002) and Kern and Haagmans (2004) used the extended Stokes and Hotine formulae for using the terrestrial gravimetric data for this purpose.Mueller et al. (2004) used the terrestrial gravity anomalies to generate them, and after that Wolf (2007) investigated the deterministic approaches to modify the integrals and validation.In fact, the spectral weighting scheme (SJÖBERG 1980and 1981and WENZEl 1981) was used by Wolf (2007).Least-squares collocation (LSC) can be used for the same purpose and Tscherning et al. (2006) considered this method and concluded that the SGG data are predictable with an error of 2-3 mE in the case of an optimal size of the collection area and optimal resolution of the data.Zielinski and Petrovskaya (2003) proposed a balloon-borne gradiometer to fly at 20-40 km altitude simultaneously with satellite mission and proposed downward continuation of satellite data and comparing them with balloon-borne data.Pail (2003) proposed a combined adjustment method supporting high quality gravity field information within the well-surveyed test area for the continuation of the local gravity field upward and validating the SGG data.Bouman et al. (2004) stated that there were some limitations in generating the SGG data using terrestrial gravimetry data and the EGMs.When such a model is used, higher degrees and orders of EGMs should be taken into account and the recent ones seem to be able to remove the greater part of the systematic errors.In their regional approach, they concluded that the bias of the SGG data can be recovered very well using LSC.Toth et al. (2005) investigated the generation of the SGG data using the Torsion balance data by LSC.Jarecki et al. (2006) did a similar study but by LSC and integral formulae without any modification.Sprlak and Novak (2014a) used the deflections of vertical for generating the gravity gradients at satellite level, also Sprlak and Novak (2014b) found integral relations between the GOCE and GRACE types of data which can be used for validation purpose of their data.The stochastic modification was used by Eshagh (2010a) to modify the second-order radial derivative of the extended Stokes formula.He proposed two methods of a) modification prior to derivative and b) derivative prior to modification.The former method is similar to the work done by Mueller et al. (2004) and Wolf (2007) but not in deterministic way of modification.Eshagh (2010b) also modified the second-order radial derivative of the Abel-Poisson formula in a least-squares sense to generate the second-order radial gradient at satellite level using an EGM and geoid model.The least-squares modification of the vertical-horizontal and horizontal-horizontal derivatives of the extended Stokes formula was done by Eshagh andRomeshkani (2011) andRomeshkani (2011) based on the theoretical study done for the possibility of this idea by Eshagh (2009a). d) The present work This paper is very similar to the work presented by Eshagh and Romeshkani (2011) with the difference of investigating the deterministic methods of modifying the integral estimatorsfor generation of the SGG data instead of the stochastic ones.This study is important as the deterministic methods are much simpler than stochastic onesas they are solely dependent on the integration domain and not the quality of the data.However, this method is not optimal in statistical point of view.Wolf (2007) has done some studies about validation of the SGG data using integral formulae modified deterministically, but she considered limited number of methods for this goal.Here, we consider methods of Molodensky (1962), Vanicek-Kleusberg (1987), Meissl (1971), Heck andGrunningar (1987), Featherstone et al, (1998), Wong and Gore (1969) methods for modifying the integral estimators for generating the SGG data from the gravity anomalies at sea level.

SATELLITE GRAVITY GRADIOMETRY OBSERVABLES
The gravitational tensor contains second-order derivatives of the gravitational field.Geocentric frame, local north-oriented frame (LNOF), orbital frame, gradiometer frame are some well-known ones for studying this tensor.
Here and after, we use the LNOF in our mathematical presentations, which is defined as the frame whose z-axis is pointing upwards in the geocentric radial direction, the x-axis towards the north and the frame is right-handed, implying that they-axis is directed to the west.Having considered harmonicity of the gravitational potential, the gravitational tensor contains 5 independent elements.GOCE measured the second-order derivatives of gravitational potential, but since our goal is to use the gravity anomalies at sea level, we consider the derivatives of the disturbing potential (T).However, by removing the contribution of the normal gravity field from GOCE data, we can derive the second-order derivatives of T. In this study,we use T,and therefore, zz T is named vertical-vertical (VV) gradient as it is the second-order derivative of T in the direction of the z-axis.T are called horizontal-horizontal (HH) gradients as there is no derivative in the direction of the z-axis.
A gradient estimator is an integral formula connecting the gravity anomaly at sea level to gradients at satellite level.Due to the limited coverage of the terrestrial data, such integral formulaes hould be modified in such a way that the contribution of the far-zone data is minimised.Here, we use some deterministic approaches to modify the estimators and test theirs uccessfulness for the validation of VV, VH and HH gradients.Below we present these integral estimatorsand we call them VV, VH and HH estimators, respectively (see e.g.ESHAGH and ROMESHKANI 2011): where is the geocentric angle between the integration and the computation points, 0  the integration domain, g  the gravity anomaly at sea level and at the integration points, d the surface integration elements  is the azimuth between the integration and computation points., i = 0, 1 and 2 (4) is the kernel of the integral terms and i is a coefficient to specify the type of the integral estimator, i.e. when 0 i  the estimator is related to Tzz, when 1 i  to Txz and Tyz and when 2 i  to Txx, Tyy and Tyx.Also,
In Table 1, we summarise some of these methods which are developed for the SGG data and the estimators presented in Eqs. ( 1)-(3).In fact, this table presents the mathematical formulae of the modified kernels based on each deterministic approach and the truncation coefficient of the corresponding integral formulae.Modified Kernel , for 0 , 0 for where and , i =0, 1 or 2, (7) are the Paul coefficients of order 0, 1 and 2 and Paul (1973) and finally we have:

NUMERICAL INVESTIGATIONS
In order to perform a numerical investigation into the estimation of the SGG data at the 250 km level, the EIGEN-51C geopotential model (Bruinsma et al. 2010) is used as a reference EGM for generating the gravity anomaly and gravity gradients for our simulation study.The aim is to use the gravity anomalies and the deterministically-modified integral estimators to produce the SGG data and compare their corresponding ones to those generated by the EGM.Through this study, the Tscherning-Rapp's (1974) model is used for generating the signal spectra of the gravity anomaly.Eshagh (2009c) showed that the results of the modification of the Stokes formula in the case of using this model and EGM08 (Pavlis et al. 2008) more or less the same and there is no need to use EGM08 for generating the signal spectra.The truncation coefficients of the estimators are simply generated by their spectral forms.Eshagh (2010a) showed that because of the inclusion of the parameters  the series will be convergent due to the elevation of the satellite, but the series should be used to high degrees.He showed that a maximum degree of 5000 should be enough for this purpose.He also discussed that the recursive formulae presented by Shepperd (1982) are not suitable when the data are of satellite type.Shepperd (1982) mentioned that his formulae might be used for the altitudes below 20 km.
The cap size of integration and the degree of modification are two important factors which should be investigated for applying any integral formula and should be specified through different numerical studies.Therefore, we divide our numerical investigations into two parts, in the first part, the behaviours of the deterministicallymodified kernels are presented, in the second one, the SGG data are generated using the gravity anomalies at sea level by using the modified integral estimators.

Behaviours of Isotropic Parts of the Modified Kernels
Now, the isotropic parts of the kernels of the integral estimators ( 1), ( 2) and ( 3) are presented.The significance of the far-zone gravity anomalies depends on these parts of the kernels.Plotting these isotropic functions shows that whether the modification has been done successfully or not.Also, it can somehow give an idea about the significance of the data being integrated and the cap size of integration.Eshagh (2014) mentioned that the contribution of the far-zone data was also dependent on the type of the data being integrated. plays a less significant role, but the opposite is true for the degree modification.
The modified kernels by the VK and Mol modifications are almost coincidedandshow that the terrestrial data have the same contribution in their corresponding estimators.Featherstone (2003) presented that the fluctuation of the spheroidal Stokes kernel is because of the oscillations of the low-degree Legendre polynomials, and this fluctuation grows with the degree of spheroidal modification.Large oscillations of the kernels cause difficulties in their numerical implementations of the integral part of the estimators.Since the value of the non-modified kernels at the rim of integration cap that is close to zero, therefore, the modified kernel by the Meissl method almost coincidesthe original kernel.
The values of the non-modified kernels of the HH and VH estimators are not small at the rim of 0  therefore the modified kernels by the Meissl method deviate from the non-modified ones.In fact, the Meissl modification needs to determine 0  by minimising the mean squared errorof the estimator in a least-square sense.Because of the limited coverage of the terrestrial data this is not always useful.To achieve the minimum mean squared error, we take its derivative respect to   0 , Sr  and equate the result to zero (SJÖBERG and HUNEGNAW 2000).For minimising the contribution of the far-zone terrestrial data in the modified estimators, the kernel values should be close to zero outside the integration domain.As can be seen, Figure 1 shows this issue for the modified kernels, by comparing them to the original ones, especially for the modified kernels of the VH and HH estimators.This can somehow show that the modification processes have been successful.The system of equations from which the modification parameters are obtained is ill-conditioned for the VH and HH estimators (ESHAGH, 2010a).This is due to inclusion of i k u in the mathematical expressions of the elements of the coefficient matrix of the systems.i k u increases unboundedly and for large degree of modifications, the system will be even more unstable.However, this instability is harmless and modification of the estimators will be successful if a simple regularisation method is used for solving the system Eshagh (2010a).

Generation of Satellite Gravity Gradiometry Data Using Modified Integral
Estimators Over Fennoscandia A grid of the terrestrial gravity anomalies is produced using the EIGEN-51C model to degree and order 360 with a resolution of 55   at sea level in Fennoscandia, limited between the latitudes 50 ° N and 75 ° N and the longitudes 0° E and 35 ° E. The SGG data in the LNOF are generated using the same model and the nonsingular formulae of the gravity gradients presented by Eshagh (2009a,b).Here, our idea is to regenerate them from the simulated gravity anomalies by using the modified integral estimators.It should be stated that, the use of EIGEN-51C to degree and order 360 is reasonable and any other EGM can be used in this respect.However, we cannot expect to recover the all frequencies to degree and order 360 due to the satellite elevation.The GOCE EGMs have been computed to degree and order 250.
We have generated the gravity anomalies with a resolution of 55   to reduce the discretisation error in the integration process but we know that recovering higher frequencies than 250 is not meaningful.Table 2 shows the statistics of the mentioned gravity anomalies and the SGG data at 250 km level and Figure 3 shows their maps over Fennoscandia.T (E)  The differences between the SGG data derived from the EIGEN-51C and those from the gravity anomalies are regarded as the external errors of the estimators.In the following, different 0  and L, are considered and tested to find the best conditions for the successful generation of the SGG data.

Test of Cap Size of Integration and Degree of Modification
Table 3 (see appendix) presents the statistics of the differences between the generated SGG data using the integral estimators modified to degree L = 150 and different 0  .Generally, it shows that the modified estimators can reproduce the SGG data with an error of about 1-2 mE.This level of error reduces insignificantly by increasing 0  meaning that 0  = 2 o is sufficient for integrating the gravity anomalies.
The large value of L = 150 causes that the estimators become insensitive to the choice of 0 This means that the series terms of their estimators are stronger than the integral ones and have more contribution to the estimates.The estimators, modified by the Meissl, HG and WG methods have this property.The errors of the estimates, derived by the modified estimators by MOL, increase withincreasing 0  .The reason is the dependence of these components to the second terms of their corresponding estimators so that by increasing 0  insignificant changes are seen in the results.
Here, we test different L to see if the choice of a smaller L is successful according to the selected 0  .Having small L decreases the dependence of the estimators to the EGM and increases the contribution of the gravity anomalies in the integration domain.
Table 4 (see appendix) presents the statistics of the differences between the generated SGG data when L = 75.
It presents that the produced SGG data are more sensitive to 0  than the case where L = 150, as expected.By selecting L = 75, the estimator will consider higher weight for the gravity anomalies than the EGM.Also, when 0 Meissl methods have this property.In the case of using the MOL method the error increases with increasing 0  due to the dependency of the estimators to the EGM.0  = 3 o seems to be suitable for this purpose.Consequently, we keep this fixed and reduce L to find the best agreement between 0  = 3 o and L. Table 5 shows the results of this investigation. 3o and L = 125 can be suitable for some of the SGG data.In this table, the methods of WG and HG and Meissl exhibit best results when the contribution of the second terms increase with increasing L to 175.Now, we can test 0  versus L> 75.The results are presented in Table 6 (see appendix).0  is decreased to 2 o and 2.5 o and L = 100 and L = 125 are considered.
From the table, one can reach to accuracies of about 1-2 mE for xz T and yz T ; and 2 mE for zz T .Therefore, their modified integral estimators with 0  = 2.5 o and L = 100 work properly.Also, when L = 125 the SGG data an error of 2-3 mE can be produced with 0  = 2 o and 0  = 2.5 o for the VV and VH gradients.This means that L = 125 is suitable and the estimator is sensitive enough to 0  and subsequently to the gravity anomalies.

CONCLUDING REMARKS
The cap size of integration ( 0  ) and the degree of modification (L), which are two important parameters to be selected before using the estimators, should beL=125 and 0  = 2.5 o , the VV gradient can begenerated with an error less than 2 mE over Fennoscandia.The choices of L=150 and 0  = 2.5 o lead to the same level of accuracy for the VH gradients and L=125 and 0  =3 o for the HH gradients.Insignificant variations are seen in the quality of the generated gradients due to changes in the resolution of the terrestrial gravity anomalies.Our numerical studies show that, in most cases, the Molodenski modification delivers better results than the rest of the deterministic modification methods.  o and 2.5 o with L = 100 and 125 from 30 30   gravity anomalies.Unit: 1 mE.Mean (M).
the parameters which are estimated according to the modification type.n g  is the Laplace harmonics of g  at the computation point located at sea level. and  are the longitude and the co- latitude of the computation point.Furthermore,


where R is the radius of the reference sphere, r the geocentric distance at is the associated Legendre function of degree n.

Figure 1 -
Figure 1 -Isotropic parts of kernels of the (a) VV,(b) VH and (c) HH integral estimator before and after modification

Table 2 -
Statistics of gravity anomalies and SGG data in Fennoscandia  .The generation of xx yy TT  is successful with an error of less than 2 mE by some estimators and again the anomalies seem to have no impact on the quality of the with 5 mE error, but the rest of the SGG data with a lower accuracy level.One cannot find a method, which delivers the best results for all SGG data.As the table shows, the Molodensky modified estimators deliver almost the best estimates for zz yz T xy T is too small therefore, the best results derived by those estimators whose series terms are more precise and accurate.The estimators modified by the HG, WG and

Table 5 (
see appendix) shows the statistics of the errors of the produced SGG data using the estimators modified with 0 = 3 o at different L .As thetable presents, by increasing L from 75 to 100, the errors of zz T , xz T and yz T reaches to 3 mE and xy T to 7 mE.By considering L= 125 even smaller errors are seen in the SGG data, which means that considering larger L is not necessary.However, for generating some of the SGG data L = 150 is required.For reproducing the HH components L = 125 and 0  = 3 o should be enough to reach to 2-3 mE error.As the table shows, the increase of L to 175 leads to errors less than 1 mE.In short, for producing the SGG data,

Table 4 -
Statistics of errors of generated gradients from 55   gravity anomalies at different 0  and L = 75.Unit: 1 mE.

Table 5 -
Statistics of errors of generated gradients for 0 3  o  at different degrees of modification.Unit: 1 mE.

Table 6 -
Statistics of errors of generated gradients using the modified estimator with 0 2