Global marine gravity gradient tensor inverted from altimetry-derived deﬂections of the vertical: CUGB2023GRAD

. Geodetic applications of altimetry have largely been inversions of gravity anomaly. Previous studies of Earth’s gravity gradient tensor mostly presented only the vertical gravity gradient (VGG). However, there are six unique signals that constitute the gravity gradient tensor. Gravity gradients are signals suitable for detecting short-wavelength topographic and tectonic features. They are derived from double differentiation of the disturbing potential and hence are susceptible to noise ampliﬁcation which was exacerbated by low across-track resolution of altimetry data in the past. However, current generation of altimetry observations have improved spatial resolutions, with some better than 5 km. Therefore, this study takes advantage of current high-resolution altimetry datasets to present CUGB2023GRAD, a global (latitudinal limits of ± 80°) 1 arcmin model of Earth’s gravity gradient tensor over the oceans using deﬂections of the vertical as inputs in the wavenumber domain. The results are ﬁrst assessed via Laplace’s equation, whereby the resultant residual gradient is virtually zero everywhere. Further analysis at local regions in the Arctic and south Indian oceans showed that T xy , T xz and T yz are the most dominant gravity gradients for bathymetric studies. This proves that bathymetric signatures in the non-diagonal tensor components are worth exploiting. Bathymetric coherence analysis of T zz over the Tonga Trench showed strong correlation with multibeam shipboard depths. This study proves that current generation of altimetry geodetic missions can effectively resolve Earth’s gravity gradient tensor. The CUGB2023GRAD model data can be freely accessed at https://doi.org/10.5281/zenodo.10511125 (Annan et al., 2024).


Introduction
It is now 50 years since Skylab, the first satellite altimetry mission, was launched in 1973.This was shortly followed by the GEOS-3 (Geodynamic Experimental Ocean Satellite) and Seasat missions, which spanned 1975Seasat missions, which spanned -1979Seasat missions, which spanned and 1978, respectively. , respectively.Satellites following these first-generation missions have shown improvements on the knowledge acquired and reflect improvements in the technologies developed during the lifespans of their respective predecessors (Escudier et al., 2018).
Developments in satellite altimetry over the years -such as the improved range accuracy from the Ka-band of SAR-AL/AltiKa -have resulted in more accurate sea surface heights (SSHs) (Verron et al., 2021(Verron et al., , 2018)).This and better spatial resolution and other advancements from the Ku-band missions (i.e., the Jason series, HY-2 series, Cryosat-2, Sentinel series and the recently launched SWOT mission) have enabled diverse applications of satellite altimetry in geodesy, geophysics, glaciology, oceanography and hydrology.
For an altimetry satellite's observations to be considered for gravity field recovery, the observations ought to have been acquired during the geodetic mission (GM) phase of the satellite (i.e., in a long repeat orbit).Most satellites begin life in the exact repeat mission (ERM) phase, where they repetitively observe the same tracks of ocean surface in a short period, resulting in better temporal resolution at the expense of spatial resolution.The GM phase is usually considered the end-of-life of the satellite although some altimetry missions' GM phase is the start or middle phase of the satellite; this yields higher across-track spatial resolution at the expense of temporal resolution.The higher across-track spatial resolution enables the mapping of short-wavelength features in the gravity field (Andersen, 2013;Andersen et al., 2021).It helps to map out finer details of mean sea surface (MSS), which is used to improve sea level anomalies for the ERM, whereas the GM phase also benefits from long-term MSS modeled through the ERM phase.The MSS is used to reduce SSH measurements from the GM phase to obtain the geoid -the surface of equilibrium potential.A description of these two satellite phases is presented well in Andersen et al. (2021).
Although the geoid is the base gravity field signal recovered through satellite altimetry, it is sensitive to longwavelength features.Conversely, short-wavelength features, which are of more interest to researchers, are better revealed through derivatives (i.e., deflections of the vertical, gravity anomaly and gravity gradient tensor) of the disturbing potential.Deflections of the vertical and gravity anomaly can be computed as first-order derivatives of the disturbing potential in the horizontal and vertical directions, respectively.Gravity gradients are the second-order derivatives and are better at revealing bathymetric and tectonic signatures.Gravity anomalies and gravity gradients can be recovered from geoid heights directly (i.e., through the inverse Stokes formula and double differentiation) or from deflections of the vertical (i.e., through the inverse Vening Meinesz formula and Laplace's equation).Previous studies have indicated that the use of deflections of the vertical is more accurate as it minimizes long-wavelength errors (Olgiati et al., 1995;Andersen, 2013).
Even though there are numerous studies about Earth's marine gravity field, most of them are focused on gravity anomaly and, to some degree, on deflections of the vertical.Studies in which gravity gradients have been studied usually discuss only the vertical component (often denoted as T zz ) of the gradient tensor although the tensor comprises six unique components.It suffices to conclude that more research has been conducted on marine gravity anomaly and its vertical derivative than full tensor of gravity gradients.Evidently, only Scripps Institution of Oceanography (SIO) releases publicly available gravity gradient models; even these are models of T zz only.One of the reasons for the limited literature on marine gravity gradient tensors is that methods for inverting them from altimetry data are comparatively few, unlike those for inverting gravity anomaly.Another signifi-cant justification for this has been the low spatial resolution of altimetry observations in the past.This is because higher differentiation of the disturbing potential results in amplification of high frequencies, which unfortunately includes noise in the signal (Sideris, 2016;Bouman et al., 2011;Wan et al., 2023).However, current datasets from the GMs of Jason-1, Jason-2, HY-2A, SARAL/AltiKa and Cryosat-2 are more accurate and densified enough to instigate a revisit to altimetryderived full tensor of gravity gradients.Generally, observations with 8 km across-track spatial resolution are deemed acceptable for gravity field recovery.With the exception of SARAL/AltiKa, which has variable across-track spatial resolution (i.e., 1-15 km) due to its drifting phase (Verron et al., 2021), the spatial resolutions of these other satellites are all better than 8 km (Andersen et al., 2021;Annan and Wan, 2021).
Therefore, this study takes advantage of the abovementioned highly densified datasets to develop CUGB2023GRAD (Annan et al., 2024), a global marine gravity gradients product consisting of all six components of the tensor.We compute the gravity gradients in the wavenumber domain through the remove-compute-restore method by using the north-south and east-west components of deflections of the vertical as input signals.
The remove-compute-restore approach demands the removal of long wavelengths in the form of an initial gravity signal (i.e., deflections of the vertical).It is later restored in the form of the desired signal (i.e., gravity gradient tensor) after computations.To this end, the global geopo-tential model EGM2008 (Pavlis et al., 2012) was used to construct the required reference gravity signals.EGM2008 was obtained as spherical harmonic coefficients from the International Centre for Global Earth Models (ICGEM: http: //icgem.gfz-potsdam.de/,last access: 1 February 2023).We used EGM2008 because it is the most widely used global geopotential model for studies involving marine gravity field (Sandwell et al., 2019;Zhang et al., 2017;Zhu et al., 2020Zhu et al., , 2019;;Andersen and Knudsen, 2019).The reference signals were simulated at a maximum degree of 2190 using the GrafLab program developed by Bucha and Janák (2013).

Derivation of gravity gradient tensor
Marine gravity gradient tensor is derived from the secondorder differentiation of the disturbing potential in the horizontal and vertical directions.It is a tensor with nine components, of which three are redundant; therefore, there are six unique tensor components.The gravity gradient tensor in the local north-oriented reference frame is defined as (Petrovskaya and Vershkov, 2006;Bouman et al., 2011) T is the disturbing potential, which according to the Bruns formula, is related to the geoid, N , via the normal gravity in geoid, γ ; i.e., T = γ N. r is the mean radius of Earth, λ is longitude and ϕ is colatitude.(x, y, z) is the coordinate in the local reference frame.x points in the latitudinal direction towards north, y points in the longitudinal direction towards west and z points in the radial direction outside of Earth.
In order to implement the remove-compute-restore approach as illustrated in Fig. 1, we compute the residual components of deflections of the vertical, ξ and η, by subtracting reference ξ and η derived from EGM2008, from the altimetry-derived ξ and η.
The relationship between components of deflections of the vertical and T can be expressed as (2) Therefore, using residual components of deflections of the vertical, it is obvious to infer from Eq. ( 2) that The first derivative of T in the vertical direction produces the radial disturbing gravity gradient, T r .Its residual form, T r , is computed in the wavenumber domain using the residual components of deflections of the vertical as inputs.
where γ is the normal gravity.k = k 2 x + k 2 y such that k x and k y are defined as 1 λ x and 1 λ y , respectively; λ x and λ y are the wavelengths in the horizontal direction.F and F −1 are the Fourier transform and inverse Fourier transform, respectively.
Having computed the residual signals T λ , T ϕ and T r from Eqs. ( 3) and ( 4), the derivative property of the Fourier transform is then applied on them to obtain (Wan et al., 2023) (5) By substituting Eqs. ( 3)-( 5) into Eq.( 1), residual components of the gravity gradient tensor can now be computed as Finally, EGM2008-simulated components of the gravity gradient tensor are then added to the residual tensor components to obtain the gravity gradient tensor.

Altimetry-derived gravity gradient tensor
The input deflections of the vertical, and the inverted gravity gradient tensor are presented in Figs. 2 and 3, respectively.Gravity gradients are known to be sensitive to topographic variations and, as such, they are good at revealing shortwavelength bathymetric and tectonic features.Even though some tectonic features can be seen in the deflections of the vertical (Fig. 2), they are better depicted in the various components of the gravity gradient tensor.For instance, the outline of the Mid-Atlantic Ridge is well revealed in Fig. 2, whereas its spreading is perfectly exposed in addition to its outline in Fig. 3. Furthermore, the boundaries of the African and South American tectonic plates can be more clearly seen in Fig. 3    Nazca tectonic plate (which borders the Juan Fernández and Easter microplates in the eastern Pacific Ocean) can be seen in the east deflection component; however, it is better exposed in T yy , T zz and T yz .These observations are attestations to one key characteristic of the gravity potential field: higher differentiations reveal high frequencies.
In order to further substantiate the short-wavelength nature of gravity gradients, Fig. 4 presents the multibeam bathymetry over the Tonga Trench near Fiji in juxtaposition with the inverted gravity field signals.The multibeam depth data were obtained through the AutoGrid web tool of the National Centers for Environmental Information (NCEI: https://www.ncei.noaa.gov/maps/autogrid/,last access: 5 December 2023).From Fig. 4, one can observe bathymetric signatures in the various gravity field signals, including the two components of deflections of the vertical.It is obvious that the bathymetric signatures resolved by the deflections of the vertical have longer wavelengths than those resolved by the gravity gradients.Additionally, this clearly proves that deflections of the vertical also contain valuable bathymetric information that is worth exploiting in the absence of the widely used gravity anomaly and vertical gravity gradient (VGG) (Annan and Wan, 2022).
To check the accuracy of the gravity gradient tensor, we test the Laplacian equation on the gravity gradient tensor, i.e., T xx +T yy +T zz = 0. Apart from its ability to tell how accurate the inverted gradient tensor is, the result from the Laplacian equation is also an indication of the effectiveness of the inversion method used to derive the signals.The residual gradient signal shown in Fig. 5 is the result of the Laplacian operation.It can be seen that the residuals are practically zero everywhere.The average residual gradient is −0.0012E, with a standard deviation of 0.0472 E. The high accuracy reported in Fig. 5 is an alternative interpretation of the accuracy of the altimetry observations.Also, it consequently serves as an indicator of the accuracy of the deflections of the vertical.This is because each component of the gravity gradient tensor is computed from the same north-south and east-west components of the deflections of the vertical.
Additionally, the coherency between the inverted T zz and multibeam bathymetry of the Tonga Trench was computed.The result is juxtaposed with corresponding coherency values computed from the DTU21GRA-derived T zz and SIO's VGG product (i.e., curv_32.1.nc)as shown in Fig. 6.The curves in Fig. 6 are nearly identical, with low coherency values seen at low wavelengths.The small coherency values at the low wavelengths are caused by upward continuation of gravity field from the seafloor to the sea surface (Smith and Sandwell, 1994).Analysis of Fig. 6 shows that with a minimum coherency of 0.5, the inverted T zz and the VGGs from DTU and SIO would poorly detect bathymetric features with wavelengths less than 25 km.Bathymetric features with wavelengths greater than 45 km would be detected with higher accuracy.This is because the coherency values of these wavelengths are greater than or approximately equal to 0.60 in each of the three vertical gravity gradients.It can be seen that the inverted T zz is slightly closer to the signals from DTU than those from SIO.To answer this question, we adapted the deep learning method of bathymetry inversion developed by the authors in Annan and Wan (2022) to assess the bathymetric significance of each tensor component.We assume that the most accurate bathymetric model would be obtained from the combined use of all six gravity gradients.Therefore, bathymetry inversion https://doi.org/10.5194/essd-16-1167-2024Earth Syst.Sci.Data, 16, 1167-1176, 2024 .Table 1 gives a summary of the bathymetric influence of the gravity gradients.In the Arctic Ocean region, the three most influential gravity gradients ranked in increasing order are T xz , T xy and T yz , whereas in the south Indian Ocean region, the order was T xy , T xz and T yz .It is worth mentioning that T xz and T yz possess infor-  mation in both vertical and horizontal directions; thus, it is possible that in addition to the vertical information, the horizontal information in these two gradients also contribute to refining bathymetric prediction.Indeed, there are previous works that have inverted bathymetry from T zz (Hu et al., 2021(Hu et al., , 2014(Hu et al., , 2015;;Tozer et al., 2019); however, as shown in this section, T zz is not the most dominant gravity gradient for bathymetric prediction.The results from this analysis are consistent with findings in our previous study (Wan et al., 2023), in which a shallow neural network was adapted for bathymetric predictions.
Apart from bathymetry inversion, gravity gradients can also be used for identifying seamounts (Kim and Wessel, 2015) and for studying their evolution (Wessel et al., 2022).Another interesting application of gravity gradients was conducted by Harper et al. (2021), in which seafloor spreading was studied by analyzing the distribution and tectonic importance of "see-saw" ridges.It is worth mentioning that the findings of Kim and Wessel (2015), Wessel et al. (2022) and Harper et al. (2021) were all derived from T zz only, which is in the vertical direction only.Therefore, by having access to the full gravity gradient tensor from the present study, in addition to the vertical perspective from T zz , it would be more interesting to analyze their findings from different directions if the other five tensor components are used.
In summary, the gravity gradients presented in this paper prove that the high spatial density and SSH accuracy of currently available GM datasets are capable of resolving the various components of Earth's gravity gradient tensor over the oceans.The results from this study further substantiate a statement of Sandwell et al. (2013), who recently asserted that gravity field signals inverted from current generation of altimetry datasets are becoming more superior in quality than most of the publicly available shipborne gravimetry datasets.Therefore, if the geoscience community would invest similar efforts in the techniques of inverting full tensor of gravity gradients like has been invested in gravity anomaly and VGG, the accuracy of future models of gravity gradient tensors would improve.We say this in light of the high range accuracy of the Ka-band mission (i.e., SARAL/AltiKa), as https://doi.org/10.5194/essd-16-1167-2024 Earth Syst.Sci.Data, 16, 1167-1176, 2024 well as the high across-track sampling from Cryosat-2 and the recently launched SWOT mission, which incorporates interferometric technology.

Data availability
The global marine gravity gradient tensor model CUGB2023GRAD is available at the ZENODO repository: https://doi.org/10.5281/zenodo.10511125(Annan et al., 2024).The dataset consists of GMT-readable geospatial grids in NetCDF file format (i.e., vector of latitudes, vector of longitudes and matrix of gravity gradients).

Conclusion
Components of deflections of the vertical, inverted from altimetry-derived SSHs, have been used as input signals to invert marine gravity gradient tensor over the globe.The resultant gravity gradient tensor was assessed via the Lapla-cian equation, with the corresponding residual gradient having magnitudes close to zero across the globe.Assessment of the inverted T zz through bathymetric coherence analysis showed that it correlates well with multibeam depths.Further analysis at local regions in the Arctic and south Indian oceans proved that the frequently used vertical gravity gradient is not the most dominant tensor component for bathymetric prediction.Instead, the most influential tensor components are the three non-diagonal gravity gradients, thereby proving that the bathymetric and tectonic information in the other five gravity gradients are worth exploiting.The average across-track sampling of current generation of altimetry observations is better than the 8 km minimum required for gravity field inversion.Therefore, with the anticipated higher accuracy and better spatial resolution of the recently launched SWOT mission and upcoming Ka-band altimetry missions, coupled with an increase in research interest and investment, the accuracy of future gravity gradient tensor models should improve.

Figure 1 .
Figure 1.Illustration of the remove-compute-restore approach used.

Figure 2 .
Figure 2. Altimetry-derived deflections of the vertical: (a) north and (b) east components.

Figure 5 .
Figure 5. Result of the Laplacian operation: (a) map view and (b) histogram of residual gradient signal.

Figure 6 .
Figure 6.Coherency between vertical gravity gradient and multibeam bathymetry of the Tonga Trench.

Figure 7 .
Figure 7. (a) Vertical gravity gradient and (b) predicted bathymetry of Arctic Ocean region; (c) vertical gravity gradient and (d) predicted bathymetry of south Indian Ocean region.

Table 1 .
Analyzing the bathymetric influence of the gravity gradients (unit: m).