Assessing the impact of corrections included in the GRACE Level-3 data on gravimetric polar motion excitation estimates

Between 2002 and 2017, the Gravity Recovery and Climate Experiment (GRACE) mission provided datasets of temporal variations of the Earth’s gravity field, among others in the form of maps of terrestrial water storage (TWS) changes (Level-3 datasets). This paper examines the impact of several corrections included in the GRACE Level-3 data on the estimated series of hydrological plus cryospheric angular momentum (HAM/CAM). We tested the role of removing the glacial isostatic adjustment (GIA) signal, adding degree-1 coefficients of geopotential (DEG1 correction), replacing the degree-2 zonal coefficient with a more accurate estimate from satellite laser ranging (C20 correction), and restoring ocean bottom pressure geopotential coefficients (GAD). The contribution of improved separation of land and ocean signals by using the Coastal Resolution Improvement (CRI) filter was also assessed. We examined the change in agreement between HAM/CAM and the hydrological plus cryospheric signal in geodetically observed excitation (geodetic residuals, GAO) when the corrections are applied. The results show that including GIA, DEG1, C20, and GAD corrections in the GRACE data increases HAM/CAM trends and reduces overall HAM/CAM variability. The exploitation of corrections slightly heightens consistency between HAM/CAM and GAO for χ1 and χ2 in the non-seasonal spectral band and for χ1 in the seasonal spectral band. The results from this study demonstrate how the different corrections combine to make the overall improvement in agreement between HAM/CAM and GAO and which corrections are most valuable.


Introduction
Between 2002 and 2017, the Gravity Recovery and Climate Experiment (GRACE) mission collected over 15 years of time-variable gravity data, which provide valuable information on large-scale mass redistribution and its variation within the Earth system (Tapley et al. 2004(Tapley et al. , 2019. The successor of the mission, GRACE Follow-On (GFO), was launched in 2018 and has since delivered more than 2 years of observations . The principle of GRACE/GFO measurements is based on tracking the distance between two satellites located in the same orbit that are positioned approximately 220 km apart (Tapley et al. 2004;Kornfeld et al. 2019). The appropriate processing of these measurements, which is routinely performed by many institutes around the world, allows to obtain quasi-monthly models of temporal changes in the Earth's gravity field.
These unique data from GRACE/GFO have been used in a number of geophysical tasks related to studying processes in the solid Earth and atmosphere, ocean, hydrosphere, and cryosphere. Such studies included the analysis of mass changes in total water storage components such as groundwater, soil moisture, and surface water (e.g. Chen 2012; Rodell et al. 2018;Scanlon et al. 2016), estimating global sea-level variation (Jeon et al. 2018), and quantifying ice mass loss and gain in polar regions (e.g. Loomis et al. 2021;Velicogna et al. 2020). GRACE/GFO temporal gravity field models are delivered by various computing centres around the world, mainly in the form of spherical harmonic (SH) coefficients of geopotential (Level-2 data) and maps of distribution of terrestrial water storage (TWS) anomalies (Level-3 data). TWS variations in GRACE/GFO Level-3 datasets can be obtained using either the SH approach or the mascon method.
The mascon approach of transforming GRACE/GFO observations into mass anomalies has gained its popularity in recent years (e.g. Lenczuk et al. 2020;Jing et al. 2019;Ran et al. 2018;Scanlon et al. 2016). The word "mascon" is a contraction of the term "mass concentration block". It was first used in the 1960s in the studies of changes in the gravity field on the Moon to describe regions of excess gravitational attraction (Müller and Sjögen 1968). Currently, mascons are used in GRACE/GFO data processing to determine mass changes, which in this case are concentrated in blocks (mascons) located on the surface of the Earth. There are several classes of mascon solutions, but in most of them, in contrast to the SH-based models, there is a direct relationship between intersatellite range-rate measurements and mass variations (Watkins et al. 2015). The advantage of mascon solutions over SH-based datasets is that there is no need to filter TWS maps that could remove part of the actual geophysical signal. Mascon solutions are also characterized by better separation between land and ocean areas. To date, three data centres have processed and delivered GRACE/GFO mascon solutions, namely the Jet Propulsion Laboratory (JPL) (Watkins et al. 2015;Wiese et al. 2016), Center for Space Research (CSR) , and the Goddard Space Flight Center (GSFC) (Loomis et al. 2019a, Luthcke et al. 2013. For geodetic purposes, GRACE/GFO datasets have been widely exploited to analyse and interpret variations in polar motion (PM) induced by changes in mass redistribution of Earth's surficial fluids, especially in the hydrosphere and cryosphere (Brzeziński et al. 2009, Cheng et al. 2011, Göttl et al. 2018, Nastula and Sliwińska 2020, Seoane et al. 2009,Śliwińska et al. 2020a, b, 2021a, b, Wińska et al. 2016. Such contribution can be described with hydrological plus cryospheric angular momentum (HAM/CAM). It has been shown that the use of gravimetric data from GRACE enables the study of PM excitation in interannual (Meyrath and van Dam 2016;Seoane et al. 2011), intraseasonal (Seoane et al. 2011), seasonal (Brzeziński et al. 2009;Chen et al. 2012;Seoane et al. 2011), and non-seasonal  spectral bands as well as the examination of PM excitation trends (Zotov et al. 2021). Previous research has also shown that the use of GRACE solutions in determination of HAM/CAM allows for a higher compliance with the sum of hydrological and cryospheric signal in the geodetically observed excitation of PM than the use of hydrological models (Jin et al. 2010(Jin et al. , 2012Nastula et al. 2019;Wińska et al. 2017). The advantage of GRACE/GFO datasets over geophysical models of the hydrosphere is that they provide information on changes in total water storage, including components such as groundwater and ice, which are usually not covered by the models. PM excitation determined from gravimetric data is often referred to as gravimetric PM excitation.
In order to obtain the highest possible accuracy of TWS and to remove unwanted signals, a number of corrections are routinely incorporated in mascon data by processing centres. The main modifications are as follows: (1) the removal of the impact of postglacial rebound (PGR) through the application of an appropriate glacial isostatic adjustment (GIA) model Peltier et al. 2015Peltier et al. , 2018, (2) adding degree-1 (DEG1) SH coefficients that are related to the geocentre position in relation to the Earth-fixed reference frame (not measured by GRACE/GFO satellites) (Swenson et al. 2008;Sun et al. 2016;Wu et al. 2012Wu et al. , 2017, (3) replacing the degree-2 zonal SH coefficient (C 20 ) determined by GRACE/GFO with a more precise estimation obtained from Satellite Laser Ranging (SLR) (Cheng et al. 2013, Loomis et al. 2019b, and (4) restoring ocean bottom pressure geopotential coefficients (GAD) (Bettadpur 2018). The latter correction is essential in the studies of ocean bottom pressure and should be zero over continents. However, this is not the case, especially in land areas near the ocean boundary, which may be due to model errors. The above-mentioned corrections are introduced by all GRACE/GFO main data centres to the Level-3 data (both from SH coefficients and from mascons). However, there are also corrections that are only included by some institutes. For example, the JPL team applied the Coastline Resolution Improvement (CRI) filter in their latest mascon solution to improve the separation of land and ocean areas.
The objective of this study was to examine the influence of various corrections included in GRACE mascon solutions on estimated HAM/CAM. As a result of the short length of the GFO data and the one-year gap in observations between the end of GRACE and start of GFO, we focused on solutions from GRACE only. We tested the impact of applying corrections both on the general characteristics of the HAM/CAM series (time series variability, trends) and on the consistency between HAM/CAM and hydrological plus cryospheric signal in geodetically observed PM excitation called here geodetic residuals or GAO. We considered the corrections included in all GRACE mascon solutions (GIA, DEG1, C 20 , and GAD), but also analysed the contribution of the CRI filter used in the solution provided by JPL. The influence of the selected GIA model on the estimated HAM/CAM trends was also examined. We studied the impact of applying the corrections on the HAM/CAM series by examining trends, seasonal changes, and non-seasonal variations in HAM/CAM. The results of this research can be helpful in understanding how the different corrections included in GRACE mascon solutions combine to make the overall improvement in consistency between HAM/CAM and GAO. This study also showed the contribution of each correction individually and indicated which corrections are most important.

Determination of the sum of hydrological and cryospheric signal in observed PM excitation-geodetic residuals (GAO)
A common method of HAM/CAM validation obtained from either GRACE data or hydrological models is to compare them with the sum of hydrological and cryospheric signal in geodetically observed PM excitation. This signal is determined from the geodetic angular momentum (GAM) series, computed from precise measurements of pole coordinates, after removing atmospheric angular momentum (AAM) obtained from models of atmospheric pressure and winds, and oceanic angular momentum (OAM) derived from models of ocean bottom pressure and currents. The resulting series are denoted as GAM-AAM-OAM, GAO, or geodetic residuals and mainly reflect the contribution of the continental hydrosphere and cryosphere to PM excitation . Several geophysical models have been used to determine AAM and OAM needed for GAO computation Seoane et al. 2011;Wińska andŚliwińska 2019). In Nastula et al. (2019), different GAO estimates were compared, and it was concluded that the choice of AAM and OAM model should not noticeably affect the correlation between GAO and HAM/CAM. In the current study, we used the following datasets for GAO computation: • Time series of χ 1 and χ 2 components of GAM computed from the Earth Orientation Parameters (EOP) 14 C04 solution (Bizouard et al. 2019) and provided by the International Earth Rotation and Reference Systems Service (IERS) Earth Orientation Center (Brzezinski 1992;Eubanks 1993); • Time series of χ 1 and χ 2 components of AAM based on the ECMWF (European Centre for Medium-Range Weather Forecasts) model (Dobslaw et al. 2010) and provided by GeoForschungsZentrum (GFZ); • Time series of χ 1 and χ 2 components of OAM based on MPIOM (Max Planck Institute Ocean Model) (Jungclaus et al. 2013) and provided by GFZ.
Mass and motion components of AAM and OAM developed by GFZ have been processed consistently with HAM (hydrological angular momentum) and SLAM (sea-level angular momentum) as a part of GFZ's Earth System Model. The determination of AAM, OAM, HAM and SLAM was performed with global mass balance being kept. This means that the total mass in atmosphere, oceans, and all land storages should be constant at any given time. To approximate barystatic sea-level variations and to account for the global mass balance among AAM, OAM and HAM, all atmospheric and terrestrial masses over the whole globe were integrated at any time-step (Dobslaw and Dill 2019). It was shown that the preservation of the global mass balance is of particular importance in the analysis of length-of-day variations (Dobslaw and Dill 2019; Yan and Chao 2012).
We selected AAM and OAM series based on ECMWF and MPIOM models, because the same models were applied in GRACE/GFO Atmosphere and Ocean De-Aliasing Level-1B Release-6 (AOD1B RL06) dataset (Dobslaw et al. 2017), which was used by the computing centres to remove non-tidal atmospheric and oceanic contributions from GRACE/GFO gravity fields. Exploiting the same models in GRACE dealiasing data and in GAO allows for the highest possible consistency between HAM/CAM and GAO.

CSR RL06M solution
We used the CSR RL06M v1 solution for the study of the impact of GIA, DEG1, C 20 , and GAD corrections on HAM/CAM, which was obtained from the CSR website (http://www2.csr.utexas.edu/grace/RL06_mascons.html). Details of this solution can be found in the works of Save et al. (2016) and . The advantage of the mascon solution provided by CSR over data delivered by other data centres is that it offers not only TWS maps with all the relevant corrections applied, but also each correction to TWS in the same form (TWS maps) and units (centimetres) separately. Hence, it is possible for users to remove the corrections processed by CSR and include their own corrections or to analyse each correction separately.
In this study, we used the corrected and non-corrected TWS as well as all corrections (given in the same format and units as the TWS data) to obtain time series of corrected HAM/CAM, non-corrected HAM/CAM and corrections to HAM/CAM resulting from GIA, DEG1, C 20 , and GAD. The relation between corrected and non-corrected TWS is described with the following equation: It should be kept in mind that only corrected data should be used for the proper interpretation of TWS changes and the resulting HAM/CAM variations. Nevertheless, the approach used in our research allows the assessment of the impact of individual corrections on HAM/CAM.

JPL RL06M solution
Compared with the previously released JPL mascon solution (JPL RL05M), the most recent JPL RL06M dataset improved the separation of land and ocean mascons. This was achieved through the introduction of the CRI filter that helps in separating the land and ocean portions of mass in each mascon lying on the land-ocean boundary at the post-processing stage (Watkins et al. 2015;Wiese et al. 2016). This operation helped to reduce leakage effects across coastlines. The JPL offers TWS data with and without the CRI filter applied. We use both variants to quantify the impact of the CRI on HAM/CAM. However, it should be remembered that for the correct interpretation of TWS changes and the resulting HAM/CAM variations, the data with the CRI filter applied should be used. The JPL RL06M solution used in this study was accessed from Physical Oceanography Distributed Active Archive Center (PO.DAAC) (https://podaac-tools. jpl.nasa.gov/drive/files/GeodeticsGravity/tellus/L3/mascon/ RL06/JPL/v02) managed by JPL.

GSFC v02.4 and GSFC RL06 solutions
It is known that PGR, which is the rise of land masses after the removal of ice cover accumulated during the last glacial period, creates measurable effects on vertical and horizontal crustal motion, global sea levels, gravity field, crustal stress, and Earth rotation (Peltier et al. 2015). We showed in our previous research (Śliwińska et al. 2020b, 2021a, b) that signals from PGR have a particular impact on HAM/CAM trends. Therefore, removing the impact of PGR by applying a GIA model is crucial to correctly assess HAM/CAM trends. Three main GIA models have been used by data centres in the most recent GRACE mascon solutions to account for the impact of PGR, namely the ICE-6G_D model (Peltier et al. 2018) used in the JPL RL06M, CSR RL06M, and GSFC RL06 solutions; the model developed by Geruo A (A et al. 2013) applied in the GSFC v02.4 solution; and the ICE-5G model used in GSFC v02.4 (Peltier et al. 2015).
The newest mascon solution from GSFC, GSFC RL06, uses the same GIA model (ICE-6G_D) as the mascon datasets from CSR and JPL (Loomis et al. 2019a). However, the previous release, GSFC v02.4, offers mascon data with two different GIA models included, either the model provided by Geruo A or ICE-5G model (Luthcke et al. 2013). We used all three available solutions from GSFC (i.e. GSFC RL06 with the ICE-6G_D model applied, GSFC v02.4 with the ICE-5G model applied, and GSFC v02.4 with the model provided by Geruo A applied) to study the impact of the chosen GIA model on HAM/CAM trend. It should be kept in mind, however, that while the two GSFC v02.4 solutions differ from each other only in terms of GIA model used, the GSFC RL06 dataset, in addition to different GIA model, also uses other background models, such as the mean static gravity field model, atmosphere and ocean dealiasing data, and data for C 20 replacement (Loomis et al. 2019a). Therefore, HAM/CAM trends from the newest GSFC mascon solution may differ from the others not only because of the different GIA model. All the GSFC data used in this study were obtained from GSFC website (https://earth.gsfc. nasa.gov/geo/data/grace-mascons).

Determination of HAM/CAM from GRACE mascon data
The excitation of PM induced by the continental hydrosphere and cryosphere is usually described with time series of two equatorial components of HAM/CAM, χ 1 (oriented along the 0°meridian) and χ 2 (oriented along the 90°E meridian). Information about TWS distribution over land is needed for their computation. The relationship between (χ 1 , χ 2 ) and TWS anomalies can be described with the following formulas (Eubanks 1993;Gross 2015): where R e is Earth's mean radius (R e = 6.371·10 6 m); C and A are Earth's principal moments of inertia (C − A = 2.6398·10 35 kg m 2 , Gross 2015); ϕ, λ, and t are latitude, longitude, and time, respectively, and dS is the surface area. The resulting HAM/CAM series are given in the units of milliarcseconds (mas). The method described above is also applicable to TWS derived from other than GRACE datasets such as hydrological or climate models. GRACE mission provides information of mass changes for both land and ocean areas, but we are particularly interested in signals from continental hydrosphere and cryosphere. In order to consider only signals from continents, in the case of mascon solutions from CSR and JPL, we applied masks to the ocean provided by those data centres. For the solutions from GSFC, the HAM/CAM series were computed excluding the mascons defined as ocean. Nevertheless, all the analyses performed in this paper could be repeated for the data without placing a mask over the oceans. Preliminary study of the impact of masking oceans on HAM/CAM showed that this operation reduces general variability of series and decreases trends in χ 1 (see S1 in supplementary material). It was also noticed that GIA and GAD corrections has visibly higher impact on HAM/CAM when mask is not used. In the remainder of this paper, however, all analyses performed were limited to land areas only. We believe that such an approach will allow to isolate, to a large extent, only signals coming from the terrestrial hydrosphere.
We applied Eqs. 2 and 3 to TWS provided by the CSR RL06M, JPL RL06M, GSFC v02.4, and GSFC RL06 solutions as well as to corrections to TWS provided by the CSR RL06M solution. All these datasets were provided in the form of maps of TWS anomalies given in centimetres. As a result of this computation, we obtained time series of HAM/CAM and corrections to HAM/CAM given in milliarcseconds (mas). Table 1 summarizes all the GRACE data used in this study and the computed HAM/CAM series.

Time series processing
The HAM/CAM and GAO series have a different temporal resolution, therefore to make them consistent with each other, we filtered them using Gaussian filter with full width at half maximum (FWHM) parameter equal to 60 days and interpolated them into the same moments for the period between January 2003 and July 2016. In the following sections, we analysed trends, seasonal variations (sum of annual, semiannual, and terannual changes), and non-seasonal variations. Trends and seasonal oscillations were computed together by fitting a model of seasonal variations to the series. This model consisted of a sum a first-degree polynomial (trend) and sinusoids and co-sinusoids with periods of 12 (annual oscillation), 6 (semiannual oscillation), and 4 months (terannual oscillation). The model of seasonal variations was fitted to the series with least squares method. Before analysing seasonal variations, we removed trends from the series and analysed them separately. Non-seasonal oscillations were determined from filtered and interpolated HAM/CAM and GAO series by removing trends and seasonal oscillations.
The processed HAM/CAM series were evaluated by comparing them with GAO using the following parameters: standard deviation (STD) of the time series, error of STD (STDerr), correlation coefficient (Corr) between HAM/CAM and GAO, percent of variance of GAO explained by HAM/CAM (relative explained variance, Var exp ), and root-mean-square deviation (RMSD). These metrics are described in detail inŚliwińska et al. (2021b).

Results
The results are divided into three parts: (1) analysis of corrections to HAM/CAM determined from corrections introduced in GRACE mascon data and their impact on consistency between HAM/CAM and GAO-performed with the use of the CSR RL06M solution (Sect. 3.1), (2) analysis of the impact of including the CRI filter on HAM/CAM series and their consistency with GAO-performed with the use

Impact of including corrections in the GRACE CSR RL06M solution on HAM/CAM
The first part of this section aims to analyse the overall size of the individual corrections and their contribution to HAM/CAM variability. This analysis is performed for overall time series, i.e. without filtering and separating into trends, seasonal, and non-seasonal oscillations. Figure 1 shows a comparison of the time series of χ 1 and χ 2 for HAM/CAM obtained with and without various corrections applied as well as for each correction to HAM/CAM. The plot is supplemented by Table 2, which presents STD of series, trend rates, and Corr between each series and HAM/CAM (without corrections). The purpose of presenting Corr is to explore how the individual corrections differ from each other and from HAM/CAM in terms of phases. As HAM/CAM and each correction describe different phenomena, high correlations between all series should not be expected. Therefore, Corr should not be treated as a reliable measure of the quality of various corrections to HAM/CAM, but only as a parameter describing the phase compatibility between individual series. Figure 1 and Table 2 show that of all corrections, the highest variability and the biggest amplitudes are observed in the time series obtained from degree-1 coefficients (DEG1). The other corrections seem to have a smaller effect on amplitudes of HAM/CAM. Correction by DEG1 is also characterized by higher trends than those for the other corrections. Negative Corr between DEG1 correction and HAM/CAM indicate differences in phases between the series, which should not be a surprise as these series reflect different phenomena. It is known that degree-1 coefficients of geopotential (C 10 , C 11 , S 11 ), which are added to GRACE data with the use of the DEG1 correction, are proportional to geocentre motion. The geocentre motion is defined as the motion of the centre of mass (CM) of the Earth system (including solid Earth,  oceans, and atmosphere) with respect to the centre of figure (CF) of the solid Earth surface (Sun et al. 2016;Swenson et al. 2008;Wu et al. 2012Wu et al. , 2017. The degree-1 coefficients can be observed directly by tracking satellites, which are orbiting the CM, while the positions of ground stations determine the CF. The GRACE satellites cannot measure the geocentre motion, so degree-1 coefficients have to be determined by the use of other methods. A common way to derive their values is to use ocean and atmospheric models in combination with GRACE coefficients for degrees 2 and higher (Sun et al. 2016;Swenson et al. 2008). Table 2 shows that the series obtained from C 20 corrections (both from SLR and from GRACE) are noticeably correlated with HAM/CAM, and the correction based on SLR is characterized by higher Corr with HAM/CAM. Both SLR and GRACE estimates of C 20 have rather minor effect on HAM/CAM trends. The trends in C 20 from SLR are only slightly smaller than the trends in C 20 from GRACE (Table  2). In general, variations in C 20 are proportional to changes in the Earth's oblateness and usually are used to compute excitation of length-of-day caused by surface mass variations (Chen and Wilson 2008). The C 20 is the largest component of the Earth's time-varying gravitational field and is essential in recovering regional mass variability (Loomis et al. 2019b). Previous studies that compared the C 20 changes obtained from SLR and GRACE (Chen and Wilson 2008;Cheng et al. 2011Cheng et al. , 2013Cheng and Ries 2017) showed that the SLR provides the most reliable estimates of this coefficient, which proves the validity of replacing the C 20 determined from GRACE with an estimate obtained from SLR by GRACE data centres in Level-3 GRACE solutions.
At this point, it should be mentioned that C 20 correction used in this study was based on C 20 series described in GRACE Technical Note 11 (TN 11, https://podaactools.jpl.nasa.gov/drive/files/allData/grace/docs/TN-11_ C20_SLR.txt) and provided by Cheng and Ries (2017). These series were exploited in CSR RL06M v01 version of mascon solution used in this paper. However, currently GFO project recommends to use GRACE Technical Note 14 (TN 14, https://podaac-tools.jpl.nasa.gov/drive/files/allData/ gracefo/docs/TN-14_C30_C20_GSFC_SLR.txt), produced by Loomis et al. (2020). In TN 11, C 20 coefficient was determined from SLR observations to five geodetic satellites, namely Laser Geodynamics Satellite (LAGEOS)-1 and 2, Starlette, Stella and Ajisai. In TN 14, in addition to the satellites mentioned, Larets (starting in February 2004) and Laser Relativity Satellite (LARES) (starting in March 2012) were also used. Nevertheless, the choice of either TN 11 or TN14 to replace the value of the C 20 coefficient obtained from GRACE with the value determined from SLR measurements should not affect HAM/CAM determination (see S2 in supporting information).
The series obtained from GAD correction are not correlated with HAM/CAM and have little impact on trends, but have non-negligible contribution to HAM/CAM variability as indicated by the STD of the series (Table 2). In general, GAD data provide monthly averages of non-tidal ocean bottom pressure variations, which are added back to the GRACE solution to restore the full ocean bottom pressure variations because in GRACE data products, the gravity coefficients are estimated relative to this model (Chambers and Bonin 2012). Because of smoothing, GRACE solutions deal with leakage of mass from the continental hydrosphere and ice sheets into the ocean, which is especially evident across land-ocean boundaries. In this study, we put a mask into the ocean areas to consider only the signals from land. The fact that GAD has a relatively high impact on HAM/CAM despite the use of a mask is a bit surprising. We believe that there is no single reason why the impact of GAD is greater than expected. One of the causes might be leakage at the coasts. The GAD data were originally distributed as SH coefficients up to degree and order 180. Thus, one could not expect sharp boundaries between land and ocean on a half-or quarter-degree latitude-longitude grids. The size of GAD over continents may also be influenced by the mask used in this study to eliminate signals from the oceans. The mask used has been shown in previous research to be important in the estimates of ocean mass (Uebbing et al. 2019). In the current study, we used quarter-degree mask provided by CSR for their mascon solution. For this mask, certain areas of the ocean, such as bays (the largest of which is Hudson Bay), waters around islands and peninsulas, are treated as land, which may result in an overestimation of the GAD correction in land areas (see S3 in supporting information). The fact that some part of the GAD signal is present on continents may indicate the imperfections of the GAD correction itself, which also affects the land in some way.
As can be seen from Fig. 1 and Table 2, GIA correction, which represents a long-term phenomenon, has little impact on amplitude variation of HAM/CAM (see values of STD in Table 2). It is also clear that the GIA removal impacts only HAM/CAM trends, and its role is stronger for the χ 2 component. This is because the rise of land masses after the removal of the huge weight of ice accumulated during the last Ice Age is especially prominent in North America and to a lesser extent in Europe and Antarctica (Peltier 2004). The definition of χ 1 and χ 2 components makes χ 2 more sensitive to mass changes in these areas than χ 1 , which clearly makes the χ 2 trend more sensitive to removing the GIA signal than the trend in χ 1 . Similar conclusions were drawn in our previous work (Śliwińska et al. 2021b) in which we analysed the impact GIA correction on trends in GAO. Nevertheless, the contribution of GIA signal to the HAM/CAM trend is not as prominent as the impact of DEG1, which relates to geocentre position (Table 2). Śliwińska et al. (2021b) suggested that apart from these two factors, the ice mass changes can also affect HAM/CAM trends. The relatively small contribution of the GIA correction to HAM/CAM trends, especially in the case of the χ 1 , might be surprising. To explain this result, we analysed regional maps of the GIA correction to the χ 1 and χ 2 components of HAM/CAM (see S4 in supporting information). It was proved that the strongest signals from GIA are from Laurentia and to a lesser extent from Scandinavia. However, for the area of oceans, which are masked in our study, the GIA contribution is still non-negligible. It is undeniable that GIA signals over the oceans are weaker than those over continents, but the ocean regions where the GIA signal occurs are more extensive than the land areas where the GIA is observable. Therefore, summing up these effects might show that globally, the influence of GIA from the oceans has a non-negligible impact on the HAM/CAM trends. Another case is that the trend rates in GIA correction to HAM/CAM are not uniform in terms of sign for all land areas. A comparison of trends in GIA correction to HAM/CAM for lands only (-0.05 mas/year for χ 1 and + 0.30 mas/year for χ 2 ) with corresponding values computed for data without masking the oceans (+ 1.02 mas/year for χ 1 and -4.60 mas/year for χ 2 ) proves unquestionable impact of applying a mask on the trends in HAM/CAM (see S4 in supporting information). In general, including the sum of all corrections in GRACE data increases HAM/CAM trends for both χ 1 and χ 2 , reduces HAM/CAM variability for χ 1 and has no effect on HAM/CAM variability for χ 2 , as shown by the STD values. The high correlation between HAM/CAM without corrections and HAM/CAM with corrections included (0.96 and 0.94 for χ 1 and χ 2 , respectively) shows that the application of corrections to the GRACE data has no meaningful impact on the HAM/CAM phase. Of all the corrections, the DEG1 series has the most noticeable impact on the amplitudes and trends in the HAM/CAM; however, it is negatively correlated with HAM/CAM, which requires further research.
After analysing the contribution of each correction to HAM/CAM, we focused on the total impact of applying all these corrections on the consistency between HAM/CAM and GAO. This analysis was conducted in terms of trends, seasonal, and non-seasonal oscillations. We examine the HAM/CAM series computed without including corrections and that obtained after introducing all corrections and compare these data with GAO. Figure 2 shows the trends alongside the overall filtered time series. It should be kept in mind that the impact of GIA is not removed in the geodetic observations of pole coordinates; therefore, these signals are also preserved in GAO. For this reason, we eliminated this contribution from GAO using the ICE-6G_D model provided by Peltier et al. (2015), which was also applied in the GRACE CSR RL06M solution. For this purpose, we used C 21 and S 21 coefficients of GIA distributed by Richard Peltier. The trend rates shown in Fig. 2 suggest that including corrections in GRACE data decreases trend differences between HAM/CAM and GAO by 0.74 mas/year for χ 1 and 0.51 mas/year for χ 2 (see Table 3). Nevertheless, despite applying corrections, GRACE-based HAM/CAM underestimates the trend rates observed for GAO and this inconsistency is particularly evident for χ 2 . This trend discordance between the series may arise because residual ocean signals may not be completely removed from GAO as the ocean model used to derive OAM is not ideal. In the case of GRACE-based HAM/CAM, ocean signals were fully eliminated by masking ocean areas.
The plots of seasonal oscillations in GAO and GRACEbased HAM/CAM (corrected and non-corrected) shown in Fig. 3 demonstrate that including corrections in GRACE solution decreases amplitudes of seasonal oscillations of the obtained series without a noticeable change in phase. In order to analyse the impact of the GRACE corrections on the compliance of HAM/CAM with GAO in seasonal spectral band in more detail, we calculated STD, STDerr, Corr, Var exp , and RMSD (Table 4). These values suggest that for seasonal variations, the inclusion of corrections results in decreased STD of HAM/CAM, as previously suggested by Fig. 3. Notably, for χ 2 , the consistency between seasonal HAM/CAM and seasonal GAO is lower for corrected data than for non-corrected data. In the case of χ 1 , all parameters (except Corr) are improved after applying corrections to the GRACE solution. The values of change in consistency between HAM/CAM and GAO that improved after applying the corrections are italicized The motion of the pole is a combination of two circular movements, namely prograde, which is a circular motion in a counter-clockwise direction, and retrograde, which is a motion in a clockwise direction. The sum of prograde and retrograde terms gives an elliptic motion. Therefore, in order to study seasonal oscillations in more detail, we decompose them into retrograde and prograde terms separately for annual and semiannual variations and present them in phasor diagrams (Fig. 4). The plots in Fig. 4 show the amplitudes and phases of seasonal oscillation in GAO and HAM/CAM in the form of vectors whose length corresponds to the magnitude of the amplitude and their direction reflects the phase of the oscillation. Figure 4 shows that for annual variation, including corrections in GRACE data decreases the amplitude of oscillation and heightens phase differences between HAM/CAM and GAO for both prograde and retrograde terms. For semiannual oscillation, including corrections in GRACE data increases HAM/CAM amplitudes for both prograde and retrograde terms and reduces the phase differences between HAM/CAM and GAO for retrograde oscillation without remarkably affecting the phase of the prograde term (Fig. 4).
The finding that the use of corrected GRACE data may reduce consistency between HAM/CAM and GAO for annual oscillation requires a more thorough analysis. Therefore, we computed amplitudes and phases of seasonal oscillations for each correction to HAM/CAM (except the GIA correction, which affects only trends in HAM/CAM) to determine which of them contributes most to the deterioration in compatibility between HAM/CAM and GAO (Table  5). Table 5 shows that of all corrections considered, DEG1 has the highest amplitudes of annual oscillation, but its phases clearly differ from the phases of GAO, which has certainly reduced the compatibility between HAM/CAM and GAO after applying corrections to the GRACE data. The correction with the second largest amplitude of annual oscillation is GAD, which also exhibits phase inconsistency with GAO. The corrections relating to C 20 have the smallest contribution Seasonal variations are considered. The parameters that improved after applying the corrections are italicized. For seasonal variations, the critical value of Corr for 68 independent points within a 95% confidence level is 0.20, while the standard error of the difference between two correlation coefficients is 0.18   to the HAM/CAM amplitude of annual variation. Similar to annual variation, DEG1 has the highest impact on the size of the amplitudes of semiannual oscillation in HAM/CAM. The vector of DEG1 is closer to the GAO vector than noncorrected HAM/CAM for semiannual retrograde oscillation; therefore, DEG1 generally improves the consistency between HAM/CAM and GAO for this spectral band. We now come to an analysis of non-seasonal oscillations in HAM/CAM and GAO. The comparison of corrected and non-corrected HAM/CAM series shows that corrections included in GRACE data do not clearly change the nonseasonal variations in HAM/CAM (Fig. 5). Overall, the use of corrections increases the variability of HAM/CAM, as shown by STD values (Table 6). The increase in consistency between HAM/CAM and GAO is detected for both χ 1 and χ 2 components when corrected GRACE data are used, as indicated by the values of STDerr, Corr, Var exp , and RMSD. Nevertheless, the differences between corrected and non-corrected HAM/CAM are not prominent. Despite improving results for non-seasonal oscillation after applying GIA, DEG1, C 20 , and GAD to the GRACE solution, the consistency between HAM/CAM and GAO remains unsatisfactory, especially if we take into account the variability of the time series. Even after applying corrections, GRACEbased HAM/CAM visibly underestimates the amplitudes of GAO, which may indicate that GRACE satellites are unable to capture some effects that have an impact on GAO. It may be relevant that both the geodetic observations used to determine GAM and the geophysical models exploited to calculate AAM and OAM are available at daily resolution, while GRACE data are only monthly. It is possible that GAO filtration did not fully remove effects that contribute to GAO variability at sub-monthly time periods. It is also apparent that the mismatch in the magnitude of amplitudes between HAM/CAM and GAO is stronger for χ 1 , which is more sensitive to mass changes within the oceans than χ 2 , possibly indicating a non-negligible role of residual ocean signals in GAO variability.

Impact of including the CRI filter in the JPL RL06M solution on HAM/CAM
In this section, we test the impact of applying the Coastline Resolution Improvement (CRI) filter used by JPL in their latest mascon solution (JPL RL06M) to improve the separation of land and ocean areas. First of all, we consider overall variability of HAM/CAM series and their consistency with GAO. As we focused on the contribution of the CRI filter, the corrections relating to the application of GIA, DEG1, C 20 , and GAD were included in HAM/CAM (both with and without the CRI applied). This analysis was performed for the overall time series, i.e. without filtering and separating into trends, seasonal, and non-seasonal oscillations. Table 7 presents STD of series, their trends as well as Corr between them. High Corr values indicate that the application of the CRI filter has little effect on the HAM/CAM phases (Table  7). However, it increases the variability of series, especially for χ 2 , as demonstrated by the STD values. After applying the CRI filter, the trend increased for χ 1 and decreased for χ 2 , and the change is more evident for χ 1 . The overall filtered time series of GAO, HAM/CAM without the CRI, and HAM/CAM with the CRI included are shown in Fig. 6 together with values of trend rates. Trend differences between HAM/CAM and GAO are presented in Table 8. Figure 6 and Table 8 show that including the CRI filter in GRACE data reduces trend differences between HAM/CAM and GAO by 0.93 mas/year for χ 1 but increases this disagreement for χ 2 (by 0.28 mas/year). This shows that the correct separation of land and ocean signals is more important for trends in the χ 1 component of HAM/CAM, which is as expected because χ 1 is more sensitive to mass changes over oceans than χ 2 , which is more responsive to mass variations over continents. The comparison of the   Non-seasonal variations are considered. The parameters that improved after applying the corrections are italicized For non-seasonal variations, the critical value of Corr for 36 independent points within a 95% confidence level is 0.28, while the standard error of the difference between the two Corr values is 0.25   We now focus on the impact of the CRI filter inclusion on seasonal variations in the HAM/CAM series and their consistency with GAO (Fig. 7). Figure 7 shows that applying the CRI filter in the GRACE JPL RL06M solution increases amplitudes of seasonal oscillations for χ 2 without any noticeable change for χ 1 , which is also indicated by the STD of the series (Table 9). The phases of the HAM/CAM series after the application of the CRI filter are almost unchanged. A detailed evaluation of HAM/CAM series with the use of GAO as a reference indicates that in general, applying the CRI filter only slightly improves the consistency between HAM/CAM and GAO for χ 1 . For χ 2 , the inclusion of the CRI filter results in improved STDerr, Var exp , and RMSD values, without any change in Corr. This shows that the CRI filtration is especially useful in improving seasonal amplitudes of the χ 2 component.
A study of phasor diagrams of annual oscillation demonstrates that including the CRI filter in the GRACE data increases amplitudes of oscillation for both the prograde and retrograde terms, and decreases the phase difference between HAM/CAM and GAO for the prograde term (Fig. 8a, b). The use of the CRI filter in the GRACE data has little effect on the semiannual oscillation-only the amplitude of the prograde term is slightly increased after applying the CRI filter (Fig. 8c, d).
An analysis of non-seasonal oscillations in HAM/CAM shows that including CRI filter in GRACE data has a minor effect on time series variability and on consistency with GAO (Fig. 9). The agreement between GAO and HAM/CAM is only slightly worsened when the CRI filter is applied to the GRACE data, as indicated by the Corr, Var exp , and RMSD values. The improvement of the separation between land and ocean signals slightly increases the STD of series and decreases the STDerr (Table 10). Nevertheless, the differences between the two variants of HAM/CAM are negligible in the non-seasonal spectral band.

Impact of the different GIA models used in the GSFC v02.4 and GSFC RL06 solutions on trends in HAM/CAM
A detailed analysis of the impact of GIA correction on HAM/CAM series is conducted in Sect. 3.1, and it was generally concluded that GIA only affects HAM/CAM trends. Because mascon solutions processed by GSFC are available with three different GIA models applied, we now study the impact of GIA model used on the trends in HAM/CAM. Figure 10 shows the overall filtered time series of HAM/CAM determined from GSFC mascon solutions with three GIA determinations, i.e. ICE-5G, ICE-6G_D, and model provided by Geruo A, and the values of the trend rates. We also present trend differences between GAO and each HAM/CAM estimate (Table 11). Figure 10 shows that all HAM/CAM series differ from each other mainly in terms of trend rates. Moreover, the choice of GIA model affects the trend in χ 2 almost twice as much as it affects the trend in χ 1 as the maximum difference in trend magnitude between particular HAM/CAM series reaches 1.14 mas/year for χ 2 and 0.55 mas/year for χ 1 . Despite applying the GIA model in each HAM/CAM series, there are still differences with respect to the GAO trend. It is likely that other processes affect trend rates in PM excitation, such as residual ocean signals or ice mass changes (Śliwińska et al. 2021b). For both χ 1 and χ 2 components, the ICE-6G_D model provides the smallest trend differences between HAM/CAM and GAO (3.14 mas/year for χ 1 and 2.51 mas/year for χ 2 ), but this should not be surprising as the same model was also applied in GAO. The non-negligible trend differences between various HAM/CAM estimates shows that the choice of the GIA model has a noticeable impact on the trends in HAM/CAM. The trends in PM excitation require further study, especially regarding the causes of inconsistency between HAM/CAM and GAO. Phase (ϕ) in (a, b) is defined by the annual term as sin(2π(t − t 0 ) + φ), ϕ in (c, d) is defined by the semiannual term as sin(4π(t − t0) + ϕ), t 0 is the reference epoch (1 January 2004)

Summary and conclusions
In this study, we examined for the first time the impact of various corrections routinely introduced in GRACE Level-3 data on the computed hydrological and cryospheric signal in PM excitation expressed as the HAM/CAM series. The analyses were performed for trends, seasonal, and non-seasonal oscillations. Our results showed that the inclusion of GIA, DEG1, C 20 , and GAD corrections in GRACE data mainly affects trends in HAM/CAM, while seasonal and non-seasonal oscillations are only marginally affected by this operation. Of the four corrections, DEG1 (for both χ 1 and χ 2 ) and GIA (for χ 2 only) have the highest impact on HAM/CAM trends; the trends in DEG1 correction are as high as + 0.53 mas/year for χ 1 and + 0.63 mas/year for χ 2 , while for GIA correction, these values are as high as -0.05 mas/year and + 0.30 mas/year for χ 1 and χ 2 , respectively. The DEG1 correction also has the highest impact on overall variability of HAM/CAM. In general, including GIA, DEG1, C 20 , and GAD corrections in GRACE data increases HAM/CAM trends for both χ 1 and χ 2 , reduces overall HAM/CAM variability for χ 1 , and has no effect on HAM/CAM variability for χ 2 . The application of those corrections to the GRACE data does not have a meaningful impact on the HAM/CAM phases. The use of GIA, DEG1, C 20 , and GAD corrections does not have a noticeable effect on the consistency between HAM/CAM and GAO in terms of seasonal and non-seasonal  Non-seasonal variations are considered. The parameters that improved after applying the CRI filter are italicized. For non-seasonal variations, the critical value of Corr for 36 independent points within a 95% confidence level is 0.28, while the standard error of the difference between the two Corr values is 0.25 oscillations; it slightly improves results for χ 1 and χ 2 in non-seasonal spectral band and for χ 1 in seasonal spectral band, and slightly worsens results for seasonal oscillation in χ 2 . The fact that adding corrections deteriorates the amplitude and phase agreement between HAM/CAM and GAO for annual oscillation and improves this consistency for semiannual oscillation requires further explanation. However, the use of corrections visibly increases trend consistency between HAM/CAM and GAO. The analyses performed in this study indicate that the CRI filter, which is applied by JPL for better separation of land and ocean signals, is especially useful for improving the trend agreement between HAM/CAM and GAO in χ 1 (improvement by 0.93 mas/year) as well as for increasing consistency between HAM/CAM and GAO in terms of amplitudes and phases of annual oscillation. However, the contribution of the CRI filtration to non-seasonal variability of HAM/CAM is almost negligible. It is likely that the inclusion of the CRI filtration is more important in the case of regional TWS analyses, especially in coastal areas, where the separation of land and ocean signals is essential than in the case of global analyses such as the HAM/CAM examination.
The comparison of HAM/CAM with different GIA models applied indicates that the choice of GIA model does not affect the trend in χ 1 , but does have an impact on the HAM/CAM trend in χ 2 . Of the three GIA models considered, ICE-6G_D provided the highest trend consistency of HAM/CAM with GAO. Nevertheless, non-negligible discrepancies remain between HAM/CAM and GAO in terms of trends, and it is likely that other processes, such as residual ocean signals or ice mass changes, affect trends in PM excitation.
The sources of inconsistencies between HAM/CAM and GAO may be not only related to inaccuracies in GRACE data processing, but also errors related to geodetic observations of the pole coordinates as well as errors in geophysical  HAM/CAM (ICE-6G_D) ± 3.14 ± 2.51 models used to determine AAM and OAM. Future research should focus on finding the causes of the discrepancies between HAM/CAM and GAO, especially in terms of trends and amplitudes of non-seasonal changes. One of the steps towards finding higher compatibility between HAM/CAM and GAO is to test alternative sources of observational data for GAM computation such as the combined EOP solution provided by JPL. Another option is to use other atmosphere and ocean models for GAO determination; however, Nastula et al. (2019) showed that the choice of the AAM and OAM series did not have a significant impact on the trend in GAO or correlation between GAO and HAM/CAM. It is clear that the observations of the new GFO mission will enhance the quality of the HAM/CAM determination; nevertheless, at the moment the length of these data is about 4 and a half years, which is insufficient to draw objective conclusions.

Conflict of interest The authors declare no conflict of interests.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecomm ons.org/licenses/by/4.0/.