Effects of New Level-1B Data on GRACE Temporal Gravity Field Models and Precise Orbit Determination Solutions

Abstract

The quality of Gravity Recovery and Climate Experiment (GRACE) observation is the prerequisite for obtaining the high-precision GRACE temporal gravity field model. To study the influence of new-generation GRACE Level-1B Release 03 (RL03) data and the new atmosphere and ocean de-aliasing (AOD1B) products on recovering temporal gravity field models and precise orbit determination (POD) solutions, we combined the global positioning system and K-band ranging-rate (KBRR) observations of GRACE satellites to estimate the effect of different data types on these solutions. The POD and monthly gravity field solutions are obtained from 2005 to 2010 by SHORDE software developed by the Shanghai Astronomical Observatory. The post-fit residuals of the KBRR data were decreased by approximately 10%, the precision of three-direction positions of the GRACE POD was improved by approximately 5%, and the signal-to-noise ratio of the monthly gravity field model was enhanced. The improvements in the new release of monthly gravity field model and POD solutions can be attributed to the enhanced Level-1B KBRR data and the AOD1B model. These improvements were primarily due to the enhanced of KBRR data; the effect of the AOD1B model was not significant. The results also showed that KBRR data slightly improve the satellite orbit precision, and obviously enhance the precision of the gravity field model.

1. Introduction

The Gravity Recovery and Climate Experiment (GRACE) satellite was launched in 2002 [1] and has greatly contributed to restoring the accuracy and spatial resolution of the Earth’s gravity field. The GRACE satellites represent a research hotspot in the field of satellite gravimetry. The significant contributions of GRACE satellites to the Earth’s gravity field are due to the satellite payload, including a global positioning system (GPS) receiver, K-band ranging (KBR) system, accelerometer, and star camera assembly. Among these, the KBR system is the primary instrument for recovering the Earth’s gravity field models and obtaining the phase difference of K/Ka microwave signals. The GPS receiver is used for precise orbit determination (POD), and the accelerometer measures the non-conservative forces on the satellite. Further, the star camera provides data on the orientation of the GRACE satellites. These payloads can provide high–low (e.g., on-board GPS) and low–low satellite-to-satellite tracking (SST) observations (e.g., KBRR), and they have been employed in the recovery of the Earth’s gravity field models and the associated variation information [1]. GRACE temporal gravity fields have been widely used in geophysics, geodesy, oceanography, and glaciology, enhancing our understanding of critical global issues, including global and regional water storage capacity variations [2,3,4,5].

Due to the significance of GRACE gravity field models in obtaining time-variable gravity field signals (e.g., hydrological signals), several research institutions have dedicated tremendous efforts to obtain better time-variable parts of the gravity field model. In addition to the official GRACE data processing centers (i.e., the Center for Space Research (CSR), GeoForschungsZentrum (GFZ), and Jet Propulsion Laboratory (JPL)), 15 other research institutes publish monthly gravity field solutions via the International Centre for Global Earth Models (ICGEM) (http://icgem.gfz-potsdam.de/ (accessed on 5 July 2021)). The methods for utilizing GRACE Level-1B data to recover gravity field models predominantly follow a dynamic approach, short-arc approach, acceleration approach, celestial mechanics approach, and energy-balance approach. These methods are further distinguished based on the type of observation. One method based on precise orbit and KBRR observations uses the dynamic approach [6,7,8,9,10], short-arc approach [11,12,13], and energy balance approach [14,15,16] to develop the gravity field model. Another method joins GPS and KBRR observations using the dynamic approach [17,18,19] or celestial mechanics approach [20,21] to obtain the gravity field model and POD solutions. The second method directly uses the original observations to simultaneously determine the precise satellite orbit and gravity field model.

Precise satellite orbits are a prerequisite for obtaining a high-quality gravity field model, and Earth’s gravity field models helps in improving the accuracy of satellite orbits. Many methods (e.g., the kinematic approach, dynamic approach and reduced-dynamic approach) can be employed to obtain centimeter-level precision GRACE orbits. In addition, the dynamic and reduced-dynamic approach [22,23] can provide continuous, smooth and stable low earth orbits (LEOs). Kang et al. [24] have used GRACE acceleration data to obtain dynamic orbits, with an accuracy of approximately 2 cm in each component for the GRACE orbits and determined that accurate orbits can be obtained using accelerometer data during high solar activity. Jäggi et al. [25] have used undifferenced or doubly differenced GPS data to obtain GRACE orbits with accuracies between 2 cm and 2.5 cm, which were validated by satellite laser ranging observations. Additionally, Li et al. [26] have compared the GRACE precision orbits obtained using the zero-difference kinematic, dynamic and reduced-dynamic approaches. Although these studies used different methods to obtain precise orbits of GRACE satellites, none of them introduced KBR inter-satellite observation data for GRACE satellite orbit determination. KBR inter-satellite observations have only been used as external data to verify the orbit accuracy of the satellite [24,25]. The impact of different release KBRR data on GRACE orbit accuracy needs to be studied. The issue would be discussed in this study.

With the improvements in the GRACE data processing algorithms and a deepened of understanding of observation error sources [27], the accuracies of the satellite orbit and the Earth’s gravity field model have been greatly enhanced. The improved data processing algorithms have typically focused on refining background force models, GRACE Level-1B data, gravity field recovery methodologies and noise models [13]. With improved understanding of the error sources of GRACE observations, JPL released the GRACE Level-1B Release 03 (RL03) dataset in 2018. Compared with the GRACE Release 02 (RL02) data, the GRACE Level-1B RL03 dataset publishes the reprocessed spacecraft attitudes and KBR data. The primary difference exists in the spacecraft attitude data, in which the time-resolution of the satellite attitude data was improved from 5 s to 1 s. As attitude data are used to prepare KBR Level-1B data, the KBR Level-1B data were also updated [28]. KBR data have been reprocessed to enhance the accuracy of the solution [29]. However, accelerometer and GPS observations remained unchanged. In addition, the AOD1B RL06 (Release 06) model was released in 2017 [30], which contains more data types, a higher time-resolution (3 h rather than 6 h), and a higher degree and order (d/o; maximum d/o of 180 rather than 100) compared with previous releases 05 (RL 05). The AOD1B product is independent and not among the GRACE level-1B data; however, it is one of the background model data necessary to obtain the GRACE Level-2 time-variable gravity field model using Level-1B data.

With the improvements in the Level-1B data accuracy, the next generation of GRACE Level-2 RL06 monthly gravity field solutions have been developed [31,32,33]. Processing centers have obtained the new-release GRACE monthly gravity field models (e.g., CSR RL06, GFZ RL06 and JPL06) based on the GRACE Level-1B RL03 and AOD1B RL06 data, and a few agencies have been updating the temporal gravity field model [13,34].

The improvements in the quality of Level-1B data benefit to improve the precision of GRACE gravity field model and POD solutions. The new release of Level-1B data provides a key opportunity to study the influence of different products on the estimation of the GRACE gravity field model and POD solutions. Quantifying the effects of these improvements on temporal gravity field models and satellite orbit can continuously improve the accuracy of these solutions. Therefore, we combine on-board GPS and KBRR observations using the dynamic approach to recover the new release GRACE monthly gravity field model and POD solutions. We also estimate the impact of each type of observation on the accuracy of the monthly gravity field model and POD solutions. This study can serve as a basis for selecting and utilizing GRACE data for various applications.

2. Methods and Background Models

2.1. Functional Model for Dynamic Approach

GRACE satellites are equipped with a GPS receiver and KBR system, providing simultaneous on-board GPS and KBRR measurements. On-board GPS measurements are used to determine the precise orbit of LEO satellites and gravity field model. However, there are limitations in recovering middle and high degrees of gravity field models using only GPS observations [35]. Combining GPS with KBRR measurements can improve the accuracy of mid- and high-degree gravity field models [1]. We produced GRACE POD and monthly gravity field solutions from two types of measurements by the dynamic approach [31,32,36,37,38,39,40]. This method can simultaneously estimate satellite orbit and time-variable gravity field models. Additionally, this method is more rigorous and can achieve a greater degree of assimilation between the observational data and model. Shanghai astronomical observatory Orbit Determination software developed by the Shanghai Astronomical Observatory was used for LEO POD and gravity field model determination [26,40].

2.1.1. GPS Observation Equation

The dynamic approach used in this study to recover gravity field models and POD solutions involves applying Newton’s equations of motion to the satellites [36]. On-board GPS measurements include dual-frequency pseudo-range and phase measurements. To eliminate ionospheric delays, we formed an ionosphere-free carrier-phase linear combination to perform the basic measurements to achieve the POD solutions and recover the Earth’s gravity field model. The zero-difference phase observation equation used in this study can be found in the literature [41]. Most of the errors in the zero-difference phase observation can be corrected. The phase ambiguity and receiver clock error were considered the estimated parameters. The observation equations for GRACE-A and GRACE-B are the same. We linearized the zero-difference ionospheric-free phase observation equation for the GRACE-A/B satellite, which can be expressed as Equation (1):

$$\delta {\rho }^G=\frac{\partial {\rho }^G}{\partial \overrightarrow{r}(t)}\frac{\partial \overrightarrow{r}(t)}{\partial \overrightarrow{r}(t_0)}\delta \overrightarrow{r}(t_0)+\frac{\partial {\rho }^G}{\partial \overrightarrow{r}(t)}\frac{\partial \overrightarrow{r}(t)}{\partial K}\delta K+\frac{\partial {\rho }^G}{\partial \overrightarrow{r}(t)}\frac{\partial \overrightarrow{r}(t)}{\partial U}\delta U+\frac{\partial {\rho }^G}{\partial t^r}\delta t^r+\frac{\partial {\rho }^G}{\partial b_3}\delta b_3$$

(1)

where $\delta {\rho }^G$ represents the correction of the zero-difference ionosphere-free GPS carrier phase, $\overrightarrow{r}(t)$ represents the GRACE-A/B position vector of the satellite at time $t$, and $\delta \overrightarrow{r}(t_0)$ represents the correction of the GRACE-A/B position vector at the initial time. $\delta K$ refers to the correction of the gravity field coefficients, which includes $\left\{C_{nm},S_{nm}\right\}$. $\delta U$ is the correction of the bias and scaling factor in the accelerometer, $\delta t^r$ refers to LEO satellite clock correction, and $\delta b_3$ is the ambiguity correction of the zero-difference ionosphere-free GPS carrier phase measurements.

We can obtain the normal equation of GPS phase observations using the least-squares method. If $K$ represents the gravity field model parameters and $Q$ represents the other parameters, the normal equation can be simplified as Equation (2):

$${\left[\begin{array}{l}{N}_{QQ}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}{N}_{QK}\\ N_{KQ}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}{N}_{KK}\end{array}\right]}_{PHA}\left(\begin{array}{l}\Delta Q\\ \Delta K\end{array}\right)={\left(\begin{array}{l}{R}_Q\\ R_K\end{array}\right)}_{PHA}$$

(2)

where $N=B^{T}PB$ represents the coefficient matrix of normal equation, and $R=B^{T}PL$ represents the residual matrix of normal equation. $B$ represents the coefficient matrix of observation equation, $P$ represents the weight of observations, and $L$ represents the residual of observations.

2.1.2. KBRR Observation Equation

In this study, we used the KBR inter-satellite range rate as the basic measurement for recovering the gravity field model. The observation equation for KBRR observations can be formed using the dynamic approach, as illustrated by Equation (3):

$$\begin{array}{ll}\delta \dot{\rho }& =\frac{\partial \dot{\rho }}{\partial {\overrightarrow{r}}_1(t)}\frac{\partial {\overrightarrow{r}}_1(t)}{\partial {\overrightarrow{r}}_1(t_0)}\delta {\overrightarrow{r}}_1(t_0)+\frac{\partial \dot{\rho }}{\partial {\dot{\overrightarrow{r}}}_1(t)}\frac{\partial {\dot{\overrightarrow{r}}}_1(t)}{\partial {\dot{\overrightarrow{r}}}_1(t_0)}\delta {\dot{\overrightarrow{r}}}_1(t_0)+(\frac{\partial \dot{\rho }}{\partial {\overrightarrow{r}}_1(t)}\frac{\partial {\overrightarrow{r}}_1(t)}{\partial U_1^{}}+\frac{\partial \dot{\rho }}{\partial {\dot{\overrightarrow{r}}}_1(t)}\frac{\partial {\dot{\overrightarrow{r}}}_1(t)}{\partial U_1^{}})\delta U_1\\ & +\frac{\partial \dot{\rho }}{\partial {\overrightarrow{r}}_2(t)}\frac{\partial {\overrightarrow{r}}_2(t)}{\partial {\overrightarrow{r}}_2(t_0)}\delta {\overrightarrow{r}}_2(t_0)+\frac{\partial \dot{\rho }}{\partial {\dot{\overrightarrow{r}}}_2(t)}\frac{\partial {\dot{\overrightarrow{r}}}_2(t)}{\partial {\dot{\overrightarrow{r}}}_2(t_0)}\delta {\dot{\overrightarrow{r}}}_2(t_0) +(\frac{\partial \dot{\rho }}{\partial {\overrightarrow{r}}_2(t)}\frac{\partial {\overrightarrow{r}}_2(t)}{\partial U_2^{}}+\frac{\partial \dot{\rho }}{\partial {\dot{\overrightarrow{r}}}_2(t)}\frac{\partial {\dot{\overrightarrow{r}}}_2(t)}{\partial U_2^{}})\delta U_2\\ & +(\frac{\partial \dot{\rho }}{\partial {\overrightarrow{r}}_1(t)}\frac{\partial {\overrightarrow{r}}_1(t)}{\partial K}+\frac{\partial \dot{\rho }}{\partial {\dot{\overrightarrow{r}}}_1(t)}\frac{\partial {\dot{\overrightarrow{r}}}_1(t)}{\partial K}+\frac{\partial \dot{\rho }}{\partial {\overrightarrow{r}}_2(t)}\frac{\partial {\overrightarrow{r}}_2(t)}{\partial K}+ \frac{\partial \dot{\rho }}{\partial {\dot{\overrightarrow{r}}}_2(t)}\frac{\partial {\dot{\overrightarrow{r}}}_2(t)}{\partial K})\delta K+\frac{\partial \dot{\rho }}{\partial Z}\delta Z\end{array}$$

(3)

where ${\overrightarrow{r}}_1,{\overrightarrow{r}}_2,{\dot{\overrightarrow{r}}}_1,{\dot{\overrightarrow{r}}}_2$ are the position and velocity vectors of GRACE-A and GRACE-B, respectively; $\delta \dot{\rho }$ refers to the correction of the KBRR observations; $\delta K$ refers to the correction of the gravity field coefficients, which include $\left\{C_{nm},S_{nm}\right\}$; and $\delta U$ refers to the corrections of the bias and scaling factor in the accelerometer. $\delta Z$ stands for the corrections of the KBRR nuisance parameter. We can obtain the normal equation of the KBRR observations based on the least-squares method. If $K$ represents the gravity field model parameters and $Q$ represents the other parameters, then the normal equation can be simplified as Equation (4).

$${\left[\begin{array}{l}{N}_{QQ}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}{N}_{QK}\\ N_{KQ}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}{N}_{KK}\end{array}\right]}_{KBRR}\left(\begin{array}{l}\Delta Q\\ \Delta K\end{array}\right)={\left(\begin{array}{l}{R}_Q\\ R_K\end{array}\right)}_{KBRR}$$

(4)

where $N,R$ represent the same as in Equation (2).

2.1.3. Combining the Normal Equations

Combining the normal equations of the GPS phase and KBRR observations with a priori accuracies of ${\sigma }_{PHA}^{}$ = 0.8 cm for the phase observations and ${\sigma }_{KBRR}^{}$ = 0.2 μm·s⁻¹ for the KBRR observations, we can obtain the combined normal equation, as illustrated by Equation (5):

$$\left(\frac{1}{{\sigma }_{PHA}^2}{\left[\begin{array}{l}{N}_{QQ}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}{N}_{QK}\\ N_{KQ}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}{N}_{KK}\end{array}\right]}_{PHA}+\frac{1}{{\sigma }_{KBRR}^2}{\left[\begin{array}{l}{N}_{QQ}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}{N}_{QK}\\ N_{KQ}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}{N}_{KK}\end{array}\right]}_{KBRR}\right)\left(\begin{array}{l}\Delta Q\\ \Delta K\end{array}\right)=\left(\frac{1}{{\sigma }_{PHA}^2}{\left(\begin{array}{l}{R}_Q\\ R_K\end{array}\right)}_{PHA}+\frac{1}{{\sigma }_{KBRR}^2}{\left(\begin{array}{l}{R}_Q\\ R_K\end{array}\right)}_{KBRR}\right)$$

(5)

where $N,R$ represent the same as in Equation (2). $\Delta K$ represents the correction of gravity field model parameters, and $\Delta Q$ represents correction of the other parameters (including bias and scaling factor in the accelerometer, LEO satellite clock, zero-difference ionosphere-free GPS carrier phase ambiguity, and KBRR nuisance parameter).

2.2. Background Models and Processing Strategy

The background model configurations and processing strategy used to determine the temporal gravity field and LEO POD solutions are critical. It includes perturbations caused by the N-body perturbation, solid Earth tide, solid Earth pole tide, ocean tide, ocean pole tide, relativistic perturbation, changes in atmospheric and oceanic mass distributions, and non-conservative forces. During data processing, the perturbations caused by the N-body perturbation, solid tides, solid pole tides, and relativity perturbation, were calculated according to IERS2010 [42]. The ocean tide and ocean pole tide perturbations were computed using the d/o 80 FES2014b model [43] and the d/o 30 Desai model [44], respectively. The non-tidal atmospheric and oceanic influences were calculated using the AOD1B model [30,45]. EIGEN-6C [46] is used as an a priori gravity field, and the static part is up to d/o 180, and the trend, annual and semi-annual part is up to d/o 50. The a priori gravity field model used same as the GFZ RL05a model [32]. Furthermore, perturbations due to non-conservative forces were calculated based on the accelerometer and satellite attitude observations. The accelerometer and empirical KBRR parameters were estimated using our proposed method. The bias and scales of the accelerometer were assessed at 1 h and 1-day intervals, respectively, whereas the empirical KBRR parameters were estimated once every 1.5 h using a polynomial fit. The fit parameter included the range-rate bias, drift, and cos/sin coefficient of the range bias.

To compare different GRACE monthly gravity field models and POD solutions, identical background models and processing strategy were employed to obtain both the solutions. The models and parameters are summarized in

Table 1. Overview of background models and processing strategy.

Forces and Parameters	Description	Notes
Forces
A priori gravity model	EIGEN-6c model	d/o 180 × 180
N-body perturbation	Sun and Moon et al.	JPL DE/LE 421
Solid Earth tide	IERS2010
Solid Earth pole tide	IERS2010
Ocean tide	FES2014b
Ocean pole tide	Desai model
Relativity perturbation	IERS2010
Non-conservative forces	Level-1B ACC and ATT
Atmosphere and oceanic variability	AOD1B
Solved parameters
Department parameters
Initial state vector	3D positions and velocity	One group per arc
Accelerometer		One group bias per hour One group factor per day
KBRR empirical parameters		One group per 90 min
LEO receiver clock		Every epoch
Ambiguities		Dependent on cycle
Global parameters
Geopotential coefficients parameters	d/o 2–60

.

3. Results

First, we combined on-board GPS and KBRR observations using the dynamic approach [32,40] to produce the GRACE monthly gravity field model and POD solutions from 2005 to 2010. The datasets from 2005 to 2010 were the most complete, with few missing data; thus, these data were used in this study. The monthly gravity field and POD solutions obtained from the Level-1B RL02 and RL03 datasets were labeled Solution 1 (SOL01) and Solution 2 (SOL02), respectively. Subsequently, we analyzed the effects of each type of Level-1B data on monthly gravity field models and POD solutions. Compared with the RL02 dataset, the Level-1B RL03 dataset updated the AOD1B model data and KBRR observations. This study analyzed the effects of the KBRR observations and AOD1B model data on the precision of the monthly gravity field and POD solutions. The results were evaluated via the following aspects: the PSD (Power Spectral Density) of KBRR data, the post-fit residuals of the GPS phase and KBRR data, precision of GRACE POD solutions, degree variances of the temporal gravity solutions, spatial distribution of the global terrestrial water storage variations, and time series of the changes in the mean mass of the Amazon Basin and the Pearl River Basin.

3.1. PSD of KBRR Observations

The difference in the frequency bands of KBRR data is important for gravity field retrieval. To analyze the differences in KBRR RL02 and RL03 data, the PSD (Power Spectral Density) of KBRR RL02 and RL03 data on 1 January 2006 was estimated. Figure 1 shows the PSD of KBRR data before estimating the monthly gravity model and precise orbit. The difference between KBRR RL02 and RL03 data is obvious, and the PSD of KBRR RL03 is lower in high frequency.

3.2. Post-Fit Residuals of KBRR Observations

KBRR observations are vital for the recovery of temporal gravity field models, and their post-fit residuals can reflect internal estimation accuracy of temporal gravity field solution. In this study, the post-fit residuals of the KBRR observations were calculated from the Level-1B RL02 and RL03 datasets. The daily RMS of the KBRR post-fit residuals were analyzed. Furthermore, the post-fit residuals of the 2006 KBRR data shown in Figure 2 revealed significantly lower and stable values for SOL02 than for SOL01, with average values of 0.218 μm·s⁻¹ and 0.199 μm·s⁻¹, respectively. This decrease in the post-fit residuals of KBRR data was attributed to the combined changes in the AOD1B and KBRR data. However, the magnitude of influence of each dataset remained undetermined. Therefore, we calculated post-fit residuals of the KBRR using the AOD1B RL05 and RL06 data (represented as SOL02-05 and SOL02-06, respectively) based on the GRACE Level-1B RL03 data, as shown in Figure 3. The post-fit residuals of the 2006 KBRR observations for the SOL02-05 and SOL02-06 exhibited average values of 0.200 μm·s⁻¹ and 0.199 μm·s⁻¹, respectively. The AOD1B data exhibited minimal influence on these solutions. The following results illustrated that the decrease in the KBRR post-fit residuals can be attributed to the improved AOD1B and KBRR data, albeit with the KBRR data exhibiting a dominant influence.

3.3. Post-Fit Residuals of GPS Phase Observations

The internal estimation accuracy of POD and monthly gravity field solutions obtained from on-board GPS and KBRR observations can be determined from the fitting residuals the observations. We calculated the GPS phase post-fit residuals of SOL01 and SOL02, respectively. We could not display all the results in the figures. As the results from 2005 to 2010 exhibited the same trends as the result in 2006, so the results in 2006 are shown. The daily root mean square (RMS) of the phase post-fit residuals of SOL01 and SOL02 (including GRACE-A and GRACE-B satellites) were shown in Figure 4. Daily RMS of phase post-fit residuals were lower than 7 mm, with mean values of 5.66 and 5.64 mm for SOL01 and SOL02, respectively. The RMS of the phase post-fit residuals from SOL01 and SOL02 did not exhibit a significant difference. These results indicate that the post-fit residuals of GPS phase observations are not significantly improved by KBRR and AOD1B data.

3.4. POD Precision

The precision of the GRACE satellite orbits is vital in the recovery of Earth’s gravity field models. The POD solutions from the Level-1B RL02 and RL03 datasets were labeled SOL01 and SOL02, respectively. To analyze the effects of each dataset on the precision of the GRACE POD solutions, the GRACE reduced-dynamic orbits published by JPL (GRACE Level-1B navigation RL 02 solution) [27] were used as the reference orbits. Daily RMS values of the 3D position differences (3D RMS) between the reference JPL orbits and GRACE-A/B orbits in the Conventional Terrestrial System are shown in Figure 5. The SOL02 orbits were more accurate than the SOL01 orbits. Relative to the reference JPL orbits, the 3D RMS differences were 3.67 and 3.46 cm for SOL01 and SOL02 of the GRACE-A orbits, respectively, and 3.54 and 3.39 cm for SOL01 and SOL02 of the GRACE-B orbits, respectively. These results demonstrate that the positions of the 3D orbit are approximately 5% more accurate when the Level-1B RL03 dataset is used; the improved orbit precision was attributed to the improvements in AOD1B and KBRR data.

To estimate the contributions of the AOD1B model and KBRR data in improving the orbit precision, we calculated the precise orbit using the AOD1B RL05 and RL06 data (represented as SOL02-05 and SOL02-06, respectively) based on GRACE Level-1B RL03 data, as shown in Figure 6. The SOL02-06 orbits were similar to SOL02-05 orbits. Relative to the reference JPL orbits, the 3D RMS were 3.47 cm and 3.46 cm for SOL02-05 and SOL02-06 of the GRACE-A orbits, and 3.38 cm and 3.39 cm for SOL02-05 and SOL02-06 of the GRACE-B orbits, respectively. The results demonstrated negligible effects of AOD1B RL05 and RL06 datasets on the orbit precision. The improvement in the POD precision was primarily due to enhanced KBRR data.

3.5. Comparison of Monthly Gravity Field Model Degree Variances

The accuracy of Earth’s gravity field model is typically expressed by the degree variance and geoid error between different models. The degree variances in this study represent variance of differences between monthly models and an average gravity model. The average of monthly gravity field solutions over the time span (2005–2010) as their reference solution. The mean gravity model for different solutions (SOL01/02, GFZ05a/06, JPL05/06, CSR05/06) is individually computed for each type of solution. The degree variance is the total intensity of the time-variable signal and noise of monthly gravity field model. Therefore, the consistency of temporal gravity field solutions can be inferred by comparing the degree variances of the gravity field models. The d/o 60 monthly gravity field solutions calculated from the Level-1B RL02 and RL03 datasets were labeled SOL01 and SOL02, respectively. Additionally, the official RL05/RL06 models from the CSR, JPL, and GFZ were denoted as CSR05/06, JPL05/06, and GFZ05a/06, respectively. The precision of monthly gravity field solutions corresponding to SOL01 and SOL02 was determined by comparing the official RL05 and RL06 models published by CSR, JPL, and GFZ. The degree variances of eight monthly gravity field solutions calculated for August 2006 and September 2006 included SOL01/SOL02, CSR05/CSR06, GFZ05a/GFZ06, and JPL05/JPL06. GRACE satellites are intrinsically inadequate for accurately measuring the C20 coefficient; thus, the calculated C20 coefficient associated with these models exhibits the variations. The degree variances of the monthly gravity field solutions from d/o 3 to d/o 60 are compared in Figure 7. The degree variances for SOL01 and the official RL05 models for August 2006 and September 2006 are shown in Figure 7a,d, respectively. The temporal signal degree variances corresponding to the SOL01 model showed consistency with those corresponding to the official RL05 models. Furthermore, the degree variances for the SOL02 and official RL06 models for August and September 2006 are shown in Figure 7b,e, respectively. The degree variances for the SOL02 model were consistent with the official RL06 models, thereby validating the reliability of the method used in this study. In contrast, the degree variances corresponding to SOL01 and SOL02 in August 2006 and September 2006 (Figure 7c,f, respectively) revealed a significant decrease after degree 25. The GRACE time-variable gravity field model signals are predominantly concentrated before the degree 25, and as the degree increases, the noise gradually increases. Figure 7c,f show that the degree variance of a new-release, time-variable gravity field model is significantly reduced after the degree 25. Therefore, based on the definition of the degree variance, it can be inferred that the sum of the time-variable signal and noise is significantly reduced. Although noise after the degree 25 is dominant, we cannot judge whether the time-variable signal is affected. This issue is discussed further in Section 3.6 and Section 3.7.

The degree variances corresponding to SOL02 in Figure 7c,f were significantly lower than those corresponding to SOL01. This decrease in the degree variance of the temporal signal can be attributed to the combined changes in the AOD1B model and KBRR data. However, the magnitude of influence from each dataset remained undetermined. Therefore, we calculated the degree variances of the temporal gravity field models using the AOD1B RL05 and RL06 models (i.e., SOL02-05 and SOL02-06, respectively) based on Level-1B RL03 data, with corresponding data for August 2006 and September 2006, as shown in Figure 8a,b. The degree variances for temporal gravity field models differed slightly. Therefore, considering the current precision level of the observations the AOD1B RL06 model exhibits a negligible impact on the overall precision of the recovered temporal gravity field solutions.

3.6. Spatial Distribution of Temporal Signals and Noise in the Monthly Gravity Field Solutions

To analyze the distribution of signals and noise in the temporal gravity field models, the temporal gravity signals corresponding to SOL01 and SOL02 were represented as equivalent water heights (EWHs). First, the C20 coefficient calculated from satellite laser ranging was used to replace the monthly gravity field solutions. Then, the P4M6 decorrelation filter method [3] was used to remove the correlations between the Stokes coefficients of odd and even degrees of the same order. This method is to subtract a fourth-order polynomial from the odd and even degree components for the residual coefficients with the same order from six to the maximum degree. Finally, the 300 km Gaussian smoothing method was applied to the temporal gravity field models to calculate the monthly global water storage signal. The comparisons of the spatial distributions of temporal signals for the SOL01 and SOL02 solutions are shown in Figure 9, with Figure 9a,c displaying SOL01 and SOL02 solutions for August 2006, and Figure 9b,d displaying SOL01 and SOL02 solutions for September 2006. From the result in Figure 9, the north–south stripes for the SOL02 solutions were evidently significantly weaker than those for the SOL01 solutions; however, the spatial distribution of temporal signals was consistent. This conclusion is also shown in the official RL06 models [31,32,33].

3.7. Distribution of Temporal Signals and Noise in Specific Regions

To quantify the differences between SOL01/02 monthly gravity field solutions and RL06 products released by the CSR, GFZ, and JPL, we calculated the time series of the changes in the mean masses of the Amazon Basin and the Pearl River Basin from 2005 to 2010. To compare local water storage changes at different spatial scales, the mean mass changes were based on monthly gravity field models for these regions (Figure 10). The monthly gravity field solutions were decorrelated and smoothed, and the EWH values corresponding to the mass changes for these regions were obtained from the weighted averages of their latitudes.

Figure 10a,b show the average water storage capacity of the Amazon Basin and the Pearl River Basin during 2005–2010 based on monthly gravity field solutions. The time series of changes in the mean masses of these regions was consistent. Regarding the Amazon Basin, the correlations between SOL01/02 model and the official RL06 model typically exceeded 0.99. For the Pearl River Basin, the correlation between the SOL02 model and the official RL06 model was greater than 0.98, whereas that between the SOL01 model and the official RL06 model was greater than 0.96. EWH values obtained from SOL01 and SOL02 models exhibited sufficient consistency with those based on the official RL06 model. Furthermore, temporal signals the SOL01 and SOL02 models were essentially identical in amplitude; however, the degree variance for the SOL02 model was significantly lower than that for the SOL01 model (Figure 7c,f). The lower power of the degree variances for degrees above 25 is interpreted to be caused by a lower noise level in the SOL02 solution as compared to SOL01.

4. Conclusions

Recovering high precision models of the Earth’s gravity field is the core objective of satellite gravity measurements. In this study, we obtained monthly gravity field models and POD solutions using the different release Level-1B data using the dynamic approach. The monthly gravity field solutions obtained in this study exhibited accuracy comparable to those of the official GRACE RL06 solutions. We also simultaneously obtained high precise GRACE orbits. To quantify the effects of the Level-1B dataset on the monthly gravity field models and POD solutions, the SOL01 and SOL02 from this research were compared. According to the post-fit residuals of observations (GPS and KBRR), POD precision, and monthly gravity field model precision, the following conclusions were reached.

(1) By comparing the monthly gravity field models and POD solutions simultaneously computed from the different Level-1B datasets, it was observed that the post-fit residuals of GPS-phase observations exhibited no apparent decrease. In contrast, utilizing the new Level-1B dataset, the post-fit residuals of the KBRR data decreased by 0.019 μm·s⁻¹ (approximately 10%) on average and improved the GRACE POD 3D position precision by 0.19 cm (approximately 5%) on average. The improved precision of these solutions was attributed to improved KBRR data, while the use of different releases of AOD1B products does not significantly change the results.

(2) By comparing monthly gravity field models from SOL01 and SOL02, the level of the degree variances of the SOL01 solution is clearly above that of SOL02 above degree 25. Furthermore, the differences between the degree variances of these models increased as the degree increased. Comparing local water storage changes at various spatial scales revealed no significant difference between the amplitude of the temporal signals from SOL01 and SOL02. Due to the significantly lower degree variance of SOL02 compared to that of SOL01, the SNR of the former was significantly higher than that of the latter.

(3) The improved precision of SOL02 (i.e., monthly gravity field model and POD solutions) was associated with the improvements in KBRR data and AOD1B model, which were primarily due to the enhanced KBRR data; the effect of the AOD1B model was less than that of the KBRR data. The enhanced KBRR data precision significantly impacted the accuracy of the temporal gravity field models and consequently enhanced the POD precision of the LEOs. However, the various AOD1B model releases do not significantly impact these solutions.