Along-Track Geopotential Difference and Deflection of the Vertical from GRACE Range Rate: Use of GEOGRACE

We present a theory and numerical algorithm to directly determine the time-varying along-track geopotential difference and deflection of the vertical at the Gravity Recovery and Climate Experiment (GRACE) satellite altitude. The determination was implemented using the GEOGRACE computer program using the K-band range rate (KBRR) of GRACE from the Level-1B (L1B) product. The method treated KBRR, GPS-derived orbit of GRACE and an initial geopotential difference as measurements used in the least-squares estimation of the geopotential difference and its formal error constrained by the energy conservation principle. The computational procedure consisted of three steps: data reading and interpolation, data calibration and estimations of the geopotential difference and its error. The formal error allowed removal of KBRR outliers that contaminated the gravity solutions. We used the most recent models to account for the gravity changes from multiple sources. A case study was carried out over India to estimate surface mass anomalies from GEOGRACE-derived geopotential differences. The 10-day mass changes were consistent with those from the MASCON solutions of NASA (correlation coefficient up to 0.88). Using the geopotential difference at satellite altitude avoids the errors caused by downward continuation, enabling the detec-tion of small-scale mass changes.


InTRODUCTIOn
The Gravity Recovery and Climate Experiment mission (GRACE; Tapley et al. 2004) has been used for mapping the Earth's time-varying gravity field for more than a decade. Most of the regular products released to the scientific community are based on series expansions of global functions such as spherical harmonics at one-month interval, which have limited temporal and spatial resolutions. Because of the loss of higher-resolution gravity change signals in both space and time in a global solution, several scientific centers have produced improved versions of the GRACE gravity solutions that are different from the regular products in the spatial and temporal resolutions. Depending on the processing method, the spatial resolution of GRACE gravity from a global solution is about 250 km (half of the Nyquist wavelength). Thus, enhancing the GRACE spatial resolution to resolve small-scaled features is one of the key efforts in the new gravity products. Some sample enhanced global gravity solutions are provided in the form of geopotential coefficients to certain degrees and order with varying nominal temporal resolutions, e.g., the 10-day solution from the Groupe de Recherches de Geodesie Spatiale (GRGS) (Bruinsma et al. 2010) and the 1-day solution from the IGG (Institut für Geodäsie und Geoinformation) (Kurtenbach et al. 2012). Because gravity variations in a particular region may not be extracted effectively when the global solution is used, a regional solution has also been attempted in the last half decade. Currently, the only available regional gravity product is the 10-day/4° solution from the Goddard Space Flight Center (GSFC) (Rowlands et al. 2005). A regional solution may capture local mass change signatures not seen in a global solution as the short wavelength components are normally limited (e.g., maximum harmonic degree/order 60) in the latter.
Most of the GRACE gravity solutions are provided as a ready-to-use product in the form of, e.g., spherical harmonic coefficient, equivalent water thickness, for which several geophysical corrections (except hydrological correction) have been applied. The comparison between the solutions is not very consistent because different data centers applied different methodologies and background models. Such methods also prevent the users from applying their own geophysical corrections (background models) to regionally extract the gravity signatures. Therefore, it is believed that an alternative gravity product is needed for users interested in improving GRACE gravity solutions at both the global and regional scales. For example, Han et al. (2011) showed that GRACE along-track K-band range rate (KBRR) measurements could be used to detect major earthquakes that could not be seen in the monthly, global-function-based gravity products. Ramillien et al. (2011), Tangdamrongsub et al. (2012) used GRACE measurements to determine geopotential differences, which were in turn used to improve regional gravity signatures.
The objective of this paper is to illustrate the GEOGRACE (GEOpotential from GRACE) computational procedure of a computer program giving the basic theories for the program. GEOGRACE allows users to obtain a wide range of gravity solutions depending on the options (e.g., resolution, signal target). This opens an opportunity for improving the gravity solution in both theoretical and practical aspects. GEOGRACE was developed based on Han et al. (2006) and Tangdamrongsub et al. (2012) with a state-of-art modification. Extending from Tangdamrongsub et al. (2012), a numerical method and a clear procedure to estimate satellite-to-satellite geopotential differences and its formal error at the GRACE satellite altitude using the public GRACE L1B (Level-1B) product (Case et al. 2010) is described. The along-track geopotential difference is determined at the finest observation interval, i.e., 5 s. The background models are stated clearly for flexible model modifications. The surface mass anomalies derived from geopotential differences over India were computed to demonstrate the potential of our numerical method.

Along-Track Geopotential Difference
The primary GRACE measurements are the K-band range (KBR) and range rate (KBRR) observed between satellites A and B. We used KBRR to determine the geopotential differences. The theories used in GEOGRACE are based largely on the works of Jekeli (1999) and Han et al. (2006). A new least-squares parameter estimation technique and a new accelerometer calibration technique are presented in this paper and in GEOGRACE, along with other features.
The convention used in the paper is the convention used for matrix operations. For example, this paper uses x T y for the inner product of vectors x and y instead of x·y. The observation equation for a KBRR, AB t o , can be expressed as (Han et al. 2006) where f is a function of the relative position and velocity, e is the residual of AB t o following the normal distribution with a zero mean and an error variance of KBRR 2 v . The subscript AB indicates the relative (differentiated) quantity between the A and B satellites, ( , , ) r r r r T In order to compute geopotential differences from GRACE we used KBRR observations with constraints, based on the energy conservation concept. Under the ideal condition, potential and kinetic energies are conserved, but in reality surface perturbing forces (e.g., atmospheric drag and solar radiation) are commonly break the conservation state. According to this idea, Jekeli (1999) derived the following approximation equation for the satellite-to-satellite geopotential difference: where F is a function containing the relative position, velocity and geopotential difference between satellites A and B, E kin is the kinetic energy, E rot is the energy due to the Earth's rotation, E diss is the dissipation energy and C AB(0) is an integration constant. The term E AB rot can be improved using a geopotential model, but this is a subject for future study and not pursued here. The expressions E AB kin , E AB rot , E AB diss can be found in Appendix A. In Eqs. (2) and (A2), all input vectors were obtained from the GRACE L1B product. As such (see also section 3), the relative positions and velocities were from GNV1B and accelerometer data were from ACC1B (GNV1B and ACC1B are the terms given to the orbit and accelerometer products at ftp://podaac-ftp.jpl.nasa.gov/all-Data/grace). However, the non-gravitational acceleration [vector a in Eqs. (A1) and (A2)] measured by the onboard accelerometers contains instrument-induced errors in the satellite reference frame (SRF, Wu et al. 2006), which must be removed (see section 2.2).
Equations (1) and (2) are the basic equations for parameter estimation and are linearized in the least-squares parameter estimation with constraints given below. The parameters to be estimated in this study are included in the vector ( , , ) In practice, we use where x 0 is a vector containing a priori values of the parameters to be estimated and Δx is its correction vector. However, it is not possible to estimate the 7 parameters simultaneously in x using only Eqs.
(1) and (2). Therefore, we used additional observation equations from the GPS-derived relative position and velocity vectors (total 6) and the geopotential difference of GRACE A and B (total 1). The geopotential difference observation was computed using Eq. (A2) in Appendix A, based on the GPS-derived position and velocity vectors and the assumption that C AB(0) is zero. Note that C AB(0) was determined during the accelerometer data calibration (see below). With a KBRR observation at a given epoch the linearized observation equation is where 0 e and A are the observation residual and design matrix (a row vector) and Δy 0 is the observation. The matrices in A elements are shown in Appendix B. In order to solve all 7 parameters (6 for GPS and one for the geopotential), we introduced 7 additional pseudo observations as The term ps e is a residual vector of the GPS-geopotential observations, I ps is a 7 × 7 identity matrix, the vector Y ps is computed as Y X X ps ps ps 0 = -, where X ps and X ps 0 are vectors containing the 6 GPS-derived observations and one geopotential difference observation, and their initial values, respectively. Combining Eqs. (4) and (5), we have where P is the weight matrix (diagonal) containing 8 diagonal elements (P i ), including three diagonal elements for the three elements in r o and one for V AB . In GEOGRACE the numerical values for the elements in P are: The magnitudes of v for P i were based on the results of Jekeli (1999) and Wu et al. (2006). See also Tangdamrongsub et al. (2012) for specific values. Note that a zero value or a large value of VAB v implies that the parameter V AB is fixed or completely unknown. In addition to the 8 observation equations we added one constraint equation to the solution, based on the energy integral in Eq. (2): , where c e and B are the residual and design matrix (a row vector) of the constraint (see below) equation, respectively. F 0 is the a priori value of the constraint. Appendix B shows the elements in B. By using the least-squares method with a small variance for the constraint, the x D U parameters can be estimated as where s 2 is a small error variance associated with the constraint equation (Koch 1987). In GEOGRACE we use s = 10 -20 . The a posteriori variance of the unit weight is estimated using where e is the residual vector associated with x D U , n is the total number of observation equations, u is the number of parameters to be estimated and q is the number of constraint equations. For each epoch we have n = 8, u = 7, and q = 1. The error covariance matrix of the estimated parameters is computed as As an example, Fig. 1 shows the standard geopotential differences errors estimated from GEOGRACE at a 5-s interval using 30 days of GRACE data in April 2007. In Fig. 1 the background model corrections were not applied to the raw GRACE KBR measurements. The magnitudes of the standard errors were 6 -7 orders smaller than the estimated geopotential differences. Figure 1 shows that some of the tracks (formed by GRACE A and B) contain large standard errors that require further processing such as outlier removal and filtering. Such large standard errors were due to the large 0 2 v W in Eq. (10). Errors shown in Fig. 1 provide clues on the quality of the estimated geopotential differences. Specifically, the standard errors from Eq. (11) could be used to identify large errors originating from the KBR measurements. This is particularly useful for KBR quality check for a regional solution aimed at looking for gravity change signals due to, e.g., seismic events and small-scaled hydrological events.

Calibrating Accelerometer Data
The accelerometers on-board GRACE satellites are important for removing the non-gravitational accelerations acting on the GRACE satellites. To demonstrate the importance of accelerometer data, Fig. 2 shows the residual geopotential (the estimated geopotential minus a reference gravitational potential generated based on GGM03S to degree and order 180; see Table 1). The residual geopotentials of GRACE A were at the level of 10 6 m 2 s -2 , which was too large. With the raw accelerometer measurements applied, the geopotentials of GRACE A were about the same order of magnitude as before and contain a negative trend. This drives the need to calibrate the accelerometer data separately for GRACE satellites A and B. In this paper the accelerometer data calibration was carried out in SRF by estimating the following three components (bias, scale, and drift) where E diss (cal) and E diss(raw) are the calibrated and raw dissipative energy, respectively. The term E diss(raw) can be computed using the dissipative energy part in Eq. (A1) when a is obtained from the raw accelerometer data. The last term in Eq. (12) represents the error part where the row vector K = (K 1 , K 2 , K 3 ) is the coefficient vector and the column vector a e = (a bias , a scale , a drift ) T is the accelerometer error vector. The elements of K are expressed as where e 3 is a unit row vector with 3 elements, a SRF raw is the raw measurement directly from the accelerometer, and r SRF o is the satellite velocity in SRF.
In the first step we applied the correction values recommended by Bettadpur (2009) to the raw accelerometer data. We then inserted Eqs. (12) to (A1) and subtracted a reference gravitational potential (U, based on GGM03S to degree and order 180) from both sides of Eq. (A1). The reference geopotential should be the one that best represents the Earth static gravity field (such as GGM03S model). We further assumed that the left side of the equation (V -U) was close to zero. The parameters a bias , a scale , and a drift and C 0 can then be determined as follows where the bold-faced symbols indicate the vectors containing the values (or coefficient) in the underlying period. In our case 3 accelerometer parameters were solved for every day, which led to a vector length of approximately 1 day / 5 s = 17280. The result of using the calibrated accelerometer data is shown in Fig. 2 (bottom). After the accelerometer data was calibrated the epoch wise time-varying geopotential differ- Because the computed VAB D included the effects of multiple mass change sources, the desired source can be extracted by removing the contributions of all other sources. For example, if the targeted signal is a mass anomaly resulting from the variation in terrestrial hydrology, the contribution from, e.g., N-body, ocean tides, solid Earth tides, atmosphere and ocean barotropic response, solid Earth pole tides and ocean pole tides must be removed from the final product. These contributions to geopotential differences ( V AB c D ) could be generated using the models described in Table 1. After these contributions are removed the time-varying geopotential differences caused by the hydrology signal are However, VAB hydro D as the estimation from the above procedure may still contain systematic errors from the KBR measurements. Such systematic errors are modeled in this study as follows where p is the empirical parameter to be estimated, ~ is the frequency associated with one cycle per revolution (CPR) of the satellite around the Earth (about 5400 s) and N is the maximum number of sinusoidal term pairs. We used N = 2, resulting in 6 terms in Eq. (15). Following the recommendations of Kim (2000) and Han et al. (2006), for very satellite revolution p 1 and p 2 were estimated once, while p 3p 6 were estimated twice. The targeted geopotential difference was obtained as (in this case, hydrology-induced mass change).

Along-Track Deflection of the Vertical
Because of GRACE's polar orbit the along-track geopotential difference from GRACE (at the satellite altitude) can be translated to the north-south component of the deflection of the vertical (DOV) direction as where g is the geoidal height, c is the normal gravity, U AB is the geopotential difference of satellites A and B due to a normal ellipsoid, all evaluated at satellite altitude (Torge 1989, pp. 37 -39), and Δs is the inter-satellite distance (about 220 km). In general, any surface phenomenon  originating from time-varying mass change can be associated with GRACE-derived DOVs using the potential theory (Heiskanen and Moritz 1967), making it possible to calibrate the parameters associated with the phenomenon. The GRACEderived DOV will be most sensitive to surface mass changes that have a large gradient in the north-south direction. Equation (16) shows the evident weakness of the GRACE orbit configuration: the twin-satellites allow only for the north-south DOV component to be resolved because of its polar orbits. Therefore, the GRACE mission yields no information about the west-east DOV component. Using GRACE-derived DOVs for gravity solutions is like using high-inclination altimeter data for marine gravity determination using DOVs from sea surface heights (Hwang and Parsons 1996). Both cases (DOV from GRACE and from a high-inclination altimeter) will lead to stripe-like artifacts in the resulting gravity fields. This weakness in GRACE can be mitigated by adding data from GRACE-like missions that have lower inclination angles. The surface mass change will be determined in this study only from along-track geopotential differences (section 3.2), rather than DOVs. However, if two GRACE-like missions of varying inclination angles are available at about the same altitude, their along-track geopotential differences can be used to determine gravity anomalies at the satellite altitude, which can be conveniently linked to phenomena producing time-varying gravity changes. Figure 3 shows the numerical procedure for direct determination of geopotential differences using GEOGRACE, which was coded in FORTRAN. The computational steps in GEOGRAE (Fig. 3) are explained as follows:
(2) The data were decimated and interpolated at a given interval (e.g., 5 s) using a cubic-spline interpolation. This step also filled in missing data and was useful in later steps, especially when a time series analysis was needed. (3) GNV data provided in the Terrestrial Reference Frame (TRF) was converted into Celestial Reference Frame coordinates (CRF, Seeber 2003). The long-term Earth orientation data based on IAU1980 framework were used to form the rotational matrices. (4) ACC data were corrected for the systematic errors (e.g., bias, trend) using the values provided by Bettadpur (2009). SCA data provided in the form of quaternion elements were used to form the rotational matrices for transforming the satellite coordinates between CRF and SRF. The calibration method described in section 2 was then applied. KBR data were also corrected using the correction values provided in the KBR product. (5) The gravitational potentials originating from several sources were computed using the background models described in Table 1. Some remarks are given as follows:  6) The reference gravitational potential was generated based on the static gravity model (GGM03S) to degree and order 180. (7) All available data were used to form the design matrix A and Y , and the geopotential differences were computed based on Eq. (7). (8) The final time-varying geopotential differences were obtained by removing the contributions from the dissipative energy, non-desired gravitational potential differences and the empirical errors. Figure 4 shows the geopotential differences in 3 different stages using the GRACE data between 5 and 8 January 2007. The geopotential differences estimated by Eq. (9) clearly show the 1-CPR signature throughout the time series. With the reference geopotential differences (step 6) removed, the remaining geopotential differences contain a 2-CPR component, suggesting that the 2-CPR terms should also be included in the empirical error model. The time series also showed that the oscillation increases with the length of satellite arc. This was most likely caused by the unmodeled errors in the accelerometer measurements, whose cumulative effect increased with time. Finally, the time varying geopotential differences caused by hydrological variations were obtained when the calibrated dissipative energy, non-desired (nonhydrological related) gravitational potential differences and the empirical errors were removed. The geopotential difference variation range was between ±0.02 m 2 s -2 . For simplicity we used the term time-varying geopotential differences to indicate the time-varying geopotential differences caused by hydrology variations below. Figure 5 shows the time-varying geopotential differences in four selected months in 2007, computed at a 5-s interval, with about 51000 observations per month. Although the observations were very dense, Fig. 5 suggested that the data were not equally distributed, with gaps at around 40°S latitude in every month. This uneven data distribution (in space) caused an aliasing effect, which could be reduced only when the orbital configuration was re-designed (Wiese et al. 2012). The time-varying geopotential differences showed clear signatures due to hydrology changes. For ex-ample, evident hydrological signatures were seen over the South African and South American regions, especially in April and October. A clear seasonal variation is seen in these 2 months as the 2 neighboring basins had different hydrological signatures in the spring and fall in the southern hemisphere. Furthermore, the variations in geopotential differences in different months over Tibet in Fig. 5 were likely caused by snow-induced mass changes. The evident geopotential signatures over Greenland and the West Antarctic may be partially caused by relatively large errors in the background models at the land-ocean boundaries. Because the post glacial rebound (PGR) effect was very small at the monthly or shorter time scales, it was not accounted for in Fig. 5. PGR correction implementation is considered in the further GEOGRACE version.
We also investigated the possibility of using the timevarying geopotential differences to improve the GRACE temporal resolution. We experimented with solutions at time scales 3 -4 times shorter than one month. After intensive offline experiments we found that it is possible to resolve the surface mass anomalies every 5 days or 1 week when there are not too many missing data records. Although the data gaps were filled at the early step (see Fig. 3) certain gravity signatures could not be properly resolved because of interpolation errors. Therefore, we decided to compute 10day solutions and the interpolated data were used only in the empirical parameter estimation stage (step 3). The MAS-CON solution (Rowlands et al. 2005) provided an independent result to validate the geopotential differences estimated in this paper. Figure 6 shows the GEOGRACE output, 10-day geopotential differences over India in the four months in 2007. The size of the selected land region in Fig. 6 was about 3 -4 times smaller than the size of the hydrologically active region like South Africa and South America. The GRACE observations in the region shown in Fig. 6 were normally contaminated by high level of noise, leading to difficulty in identifying small-scaled mass change features. Therefore, the experiment over India was a good example to verify the quality of our mass anomaly solution. The rightmost column in Fig. 6 shows the time-varying geopotential differences in one complete month. The spatial data coverage for the 10-day solutions shown in Fig. 6 are different because of bad data removal. In Fig. 6 significant variations in geopotential differences occurred at around 20 -30°N in all 4 months, while the variation in October was extended to 15°N. Note that different 10-day solutions exhibited different variations. This suggested that the one-month gravity solution was not able to reveal terrestrial water variations at time scales smaller than one month (e.g., weekly).

Regional surface Mass Variation from Geopotential Difference
We demonstrated the mass variation computation based on the GEOGRACE result. The direct Newtonian representation could be used as a connecting function be-tween the time-varying geopotential differences and the surface mass anomalies. It is well-known that without applying the weight matrices or filtering to the potential-mass system equations the resulting mass anomalies would be greatly contaminated by the north-south (N-S) stripe errors, which were demonstrated in Fig. 7. The N-S stripe Fig. 6. Sample 10-day time-varying geopotential differences: day 1 -10 in column 1, day 11 -20 in column 2 and day 21 -30 (or 31) in column 3. The one-month solutions are shown in column 4. (Color online only) errors were caused by misalignments between the GRACE measurements and background models in both the spatial and temporal domains, so that the surface mass anomalies were obscured by the N-S stripe errors and could not be resolved without regularization or filtering. Similar to the method described in Tangdamrongsub et al. (2012), the space-time covariance function representing the surface mass anomalies distribution over India was constructed based on the 3-hours Global Land Data Assimilation System (GLDAS, Rodell et al. 2004). In this paper the derived covariance function followed the Gaussian form ( , ) (2 150 10 ) where } is the spherical (angular) distance between the geographical grids and Δt is the time difference between the solutions. The regularization parameters used to stabilize the solutions were estimated based on the L-curve criterion method (Hansen and O'Leary 1993). It was difficult to estimate the standard time-varying geopotential differences errors caused by the hydrology signal due to the difficulty in modeling the background model errors. It is emphasized here that the standard geopotential differences errors given in Fig. 1 were derived from KBRR only and should not be used directly until the error terms from the background models are applied. We decided to use the L-curve criterion as discussed in Tangdamrongsub et al. (2012) to adjust an initial error covariance that was always an identical matrix, relative to the signal covariance matrix. This adjustment ensures the optimal surface anomaly determination from KBRR. Figure 8 shows the 10-day surface mass anomalies for the 4 months in Fig. 6. Due to the downward continuation process in the surface mass anomalies computation from the geopotential differences, we found large edge effects near the boundary. Therefore, the resulting surface mass anomalies close to the boundary were masked and not discussed here. In order to compare our results with the MASCON solutions the surface mass anomalies outside the Indian continent were excluded in the comparison. Figure 8 shows that large surface mass anomalies occurred at around 20 -30°N in January, April, and July and at around 20°N in October 2007. These surface mass anomalies concurred with the large geopotential difference variations seen in Fig. 6. It turned out that the surface mass anomaly magnitude depended on the regularization parameter choice. However, since the surface mass anomaly was beyond the scope of this paper, we would not perform spectral analysis and further comparisons for this quantity. Rather, we determined the correlation coefficients and Root Mean Square (RMS) differences between the surface mass anomalies from this paper and from the MASCON solutions, as shown in Table 2. Despite the different 10-day batch windows used in this paper and in MASCON, the two solutions had correlation coefficients reaching 0.88 with the RMS difference of approximately 5 cm. The comparison with MASCON verified the reliability of geopotential differences derived by GEOGRACE.

DIsCUssIOns
A direct use of geopotential differences at the GRACE altitude can avoid the error due to downward continuation commonly seen in surface mass anomalies estimated from GRACE. With the undesired background temporal signals properly modeled, the along-track geopotential difference can be used to detect episodic events such as earthquakes. For example, for a sufficiently large earthquake, the alongtrack geopotential differences (also geoidal difference) at satellite altitudes around the epicenter before and after the earthquake will be different and one can use such a difference (at the GRACE altitude) to constrain the fault parameters associated with this earthquake. Future GRACE-like missions will greatly enhance the sensor performance to detect small-scaled signatures caused by mass changes. Our method greatly increases the sensitivity of gravity signal detection, especially for signals that occur along satellite tracks.
Although the standard KBRR errors shown in Fig. 1 are given only for the demonstration, the final standard errors can be computed when the background model contributions,  V AB C l are accounted for. By introducing V AB C , the partial derivatives described in Eqs. (B3) -(B6) will be altered and the geopotential differences in Eq. (9) will need to be reestimated. Using a full error covariance matrix will lead to the optimal signal-to-noise ratio in the GRACE solution, and it is a subject for future research.
Several improvements in geopotential estimation by GEOGRACE can be made over the current result, e.g., updating the GRACE satellite orbits and the background models (either globally or regionally), and improving the potential rotational term (Guo et al. 2015). The benefits of such improvements in resolving temporal mass change are seen in publications such as Han et al. (2011). In addition, the systematic errors in the KBR measurements can be reduced by pre-processing (Ramillien et al. 2011), and the estimation of empirical errors in the last step (section 3) will be no longer required.
The line-of-sight KBR is sensitive only to the Earth's time-varying gravity field horizontal component and is insensitive to the radial gravity variation. This weakness in KBR measurements in sensing gravity changes can be removed by using measurements in multiple dimensions and from multiple sensors. For example, the GPS-derived orbital perturbations of COSMIC (Lin et al. 2012) are shown be able to estimate low-degree gravity changes that agree well with the satellite ranging-only and GRACE-only results. Combining GRACE KBR measurements with GPS measurements from low-earth orbiting satellites such as CHAMP and COSMIC may also enhance the mass variation temporal and spatial resolutions from GRACE.

COnClUsIOns
We showed a clear numerical procedure and a computer program, GEOGRACE, to compute along-track geopotential differences and their formal errors at the GRACE altitude. In the numerical examples given in this paper we used the state-of-art background models and a best static geopotential model. The results show strong hydrological signals over some hydrologically active regions like South Africa and South America and signals related to ice mass variations over Central Asia. Using 5-s geopotential differences enables a user to obtain temporal gravity changes at short spatial and temporal scales. The high correlation coefficients between the MASCON solutions and our solutions for mass anomaly suggest that our approach for the alongtrack geopotential difference determination is adequate and can resolve surface mass anomalies at time scales less than 1 month.

AppEnDIx A: EnERGy InTEGRAl EqUATIOns fOR sInGlE AnD TwIn sATEllITEs
The energy integral equation of one satellite can be expressed as