Improving the Altimeter-Derived Surface Currents Using Sea Surface Temperature (SST) Data: A Sensitivity Study to SST Products

Measurements of ocean surface topography collected by satellite altimeters provide geostrophic estimates of the sea surface currents at relatively low resolution. The effective spatial and temporal resolution of these velocity estimates can be improved by optimally combining altimeter data with sequences of high resolution interpolated (Level 4) Sea Surface Temperature (SST) data, improving upon present-day values of approximately 100 km and 15 days at mid-latitudes. However, the combined altimeter/SST currents accuracy depends on the area and input SST data considered. Here, we present a comparative study based on three satellite-derived daily SST products: the Remote Sensing Systems (REMSS, 1/10 ∘ resolution), the UK Met Office OSTIA (1/20 ∘ resolution), and the Multiscale Ultra-High resolution SST (1/100 ∘ resolution). The accuracy of the marine currents computed with our synergistic approach is assessed by comparisons with in-situ estimated currents derived from a global network of drifting buoys. Using REMSS SST, the meridional currents improve up to more than 20% compared to simple altimeter estimates. The maximum global improvements for the zonal currents are obtained using OSTIA SST, and reach 6%. Using the OSTIA SST also results in slight improvements (≃1.3%) in the zonal flow estimated in the Southern Ocean (45 ∘ S to 70 ∘ S). The homogeneity of the input SST effective spatial resolution is identified as a crucial requirement for an accurate surface current reconstruction. In our analyses, this condition was best satisfied by the lower resolution SST products considered.


Introduction
Oceanic currents are a key factor in modulating both the short-term and climatological dynamics of the ocean-atmosphere system. On a large scale, their monitoring is needed to evaluate the meridional transport of heat and salt and better predict ocean thermohaline circulation variability and changes [1]. At the oceanic meso-to sub-mesoscales, the characterization of the marine currents is also crucial. Mesoscale eddies are non-stationary and energetic recirculation features with horizontal scales of 10 to 100 km and can persist on timescales of weeks to months. They can migrate for several miles carrying heat, salt and nutrients, and their perturbations can drive intense vertical exchanges. Submesoscale features like eddies, fronts and filaments are characterized by present study, we tested both the CMEMS OSTIA (1/20 • resolution), and the Multiscale Ultra-High resolution SST (1/100 • resolution), both relying on the interpolation of several satellite and in-situ SST observations (see Section 2 for more details). We focused on the period ranging from 2014 to 2016, which yielded more accurate surface current reconstructions in RS18.

Data
The following datasets were used in our study, all covering the period ranging from 2014 to 2016: 1.
The altimeter-derived sea surface currents computed at CLS (the list of undefined abbreviations is available at the end of the manuscript) in the framework of the DUACS project and distributed by the CMEMS Sea Level Thematic Assembly Center: two different products were used, referred to as "2SAT" and "4SAT". Both products are gridded data provided on a regular 1/4 • grid. The 2SAT product is calculated merging observations from two altimeters: Jason-2 and AltiKa, with Jason 3 only from March 2016. The 4SAT product is obtained using a four altimeter constellation: Jason-2(3), Cryosat, Altika and HY-2A. The 4SAT dataset can be seen as the best altimeter-derived surface current estimate in the 2014-2016 period. On the other hand, the 2SAT version is less accurate, being based on observations from only two altimeters like for the altimeter-derived currents of the early altimetry era (the early 1990s) [36]. A 2SAT altimeter constellation is the minimum required for resolving the larger mesoscales circulation structures, providing spatial-temporal resolutions around 150-200 km and 10-15 days. However, merging information from four (or more) altimeters enables to improve the retrieval of mesoscale features missing in the 2SAT estimates, achieving effective spatial-temporal resolutions around 100 km and, at best, 7 days ( [17,19,37], https://www.aviso.altimetry.fr); 2.
The SST daily observations are the REMSS processing centre: we used the high resolution product based on the combination of microwave (from TMI, AMSR-E, AMSR2, WindSat and GMI) and infrared (Terra MODIS, Aqua MODIS, VIIRS) data. These SST observations are corrected using a diurnal model and represent a foundation SST ( 10 m depth) [38,39]. These data are calculated using an Optimal Interpolation scheme with 100 km and 4-day correlation scales and are provided on a 1/10 • regular grid ( [40], http://www.remss.com); 3.
The SST daily observations from the Operational Sea Surface Temperature and Sea Ice Analysis (OSTIA) system: OSTIA uses satellite data including AMSR-E, AVHRR (GAC+LAC), IASI, SEVIRI, TMI, GOES, SSMIS, SSM/I sensors together with in-situ observations to determine the sea surface temperature [41,42]. The analysis is performed using the optimal interpolation (OI) scheme described by [43]. The analysis is produced daily and is provided on a 1/20 • regular grid; 4.
The Multiscale Ultra-high Resolution global SST analyses (MUR): such SST dataset relies on high and low resolution satellite observations in the microwave and infrared bands (e.g., from AMSR-E, AMSR-2, WindSat, AVHRR, MODIS). The satellite data are also merged with in-situ SST estimates via a Multiresolution Variational Analysis Method [44]. The analysis is produced daily and is provided on a 1/100 • regular grid; 5.

Methods: The Optimal Reconstruction
We reconstruct the oceanic surface currents relying on an optimized blending of coarse-resolution geostrophic currents from altimetry measurements and higher resolution multi-sensor SST products. The method rationale and theoretical background is described in [34] and its application to satellite derived data is thoroughly described in [33,46,47] both on global and regional scales, including both qualitative and quantitative validations. For clarity, the algorithm to derive the optimal currents (OPC hereinafter) is reiterated in this section. The reader is also referred to [33,34,46,47] for further details. The sea surface currents are inferred from the SST dynamical evolution equation: In Equation (1) (u,v) are, respectively, the zonal and meridional components of the ocean surface flow, (x,y) are the zonal and the meridional directions and F is the forcing term, representing the source and sink terms for the SST dynamical evolution equation. Piterbarg 2009 derived the expressions of an optimized flow field accounting for the merged contribution of a large-scale, background flow with the dynamical information contained in high-resolution tracers. If the background flow is given by the geostrophic currents and the high-resolution tracer by SST, the optimized zonal and meridional flows can be expressed by Equation (2): where the subscript (GEO/OPC) stands for (geostrophic/optimal current), A = ∂ x SST, B = ∂ y SST, E = ∂ t SST − F. In Equation (2), the subscript "t" stands for temporal derivative and the subscripts "x,y" respectively indicate the spatial derivative in the zonal and meridional directions. Equation (2) expresses the Piterbarg 2009 method rationale: the geostrophic currents are corrected by means of a factor (u CORR , v CORR ) depending on the spatio-temporal derivatives of the high-resolution tracer observations and on the forcing term (F) regulating the tracer dynamical evolution. However, Equation (2) is only valid when F is assumed to be known perfectly, i.e., with a strictly null uncertainty.
In RS18 and in the present study, F is approximated by a low-pass filtered SST temporal derivative (cut-off wavelength at 500 km). Physically, this indicates that the major contributor to the SST evolution source/sink terms is identified with a large scale input, i.e., the atmospheric forcing. Such an approximation was firstly derived by Rio et al. 2016 [46] in a numerical experiment and successfully tested with space-based observations by RS18 and Ciani et al., 2019 [47]. When the F term is approximated, an uncertainty on its estimates has to be accounted for. RS18 proved this is an essential requirement to avoid the appearance of spurious surface current values in the reconstructed OPC, especially when the SST spatial gradient |∇SST| ranges from 0 to 1.5 × 10 −5 • C × m −1 (see also Section 2.3 for the formal expression of |∇SST|).
If F is not known perfectly, the Piterbarg method equations account for the uncertainties on both the background geostrophic currents and the forcing term F. The correction terms in Equation (2) are modified as indicated below: where: Moreover, the set of Equation (5) illustrates the functions f , g (expressed as functions of a generic variable γ) and the quantities p, q, α and β appearing in Equation (4): In Equation (5), σ u , σ v respectively indicate the uncertainties on the zonal and meridional background velocities while h is the uncertainty on the forcing term.
To compute σ u and σ v , the Altimeter-derived geostrophic currents (2SAT and 4SAT) have been interpolated along the trajectories of the drogued SVP drifters (in the 1993 to 2017 period). Then, the Root Mean Square Error (RMSE) between the geostrophic and the in-situ measured currents has been computed on 4 • × 4 • boxes (as in RS18).
Using a similar approach, the SST measured by the SVP drifting buoys was used to derive h. Recalling that F is approximated as the low-pass filtered satellite SST temporal derivatives (F sat ), the SVP drifters allow to compute an in-situ estimate of F (F i ) via the following expression (for every drifter): where dT = 1 day and the subscript i stands for in-situ measured. h is then computed interpolating F sat along the drifter trajectories and computing the RMSE in 4 • × 4 • boxes. This operation was repeated for each of the SST datasets involved in our study. The only difference with respect to RS18 is that a common 2002-2017 period was chosen for the OSTIA and MUR SST (constrained by MUR SST data availability). As an example, maps of h and σ u/v are provided as Supplementary Materials ( Figure S1). Finally, RS18 also pointed out the need of an empirical calibration of h to optimize the OPC reconstruction. Indeed, based on comparisons with in-situ measured currents (see also Equation (6)), they obtained the best global performances if h was increased by a factor of three. In our study, we rely on the Piterbarg 2019 equations calibrated as in RS18, while specific calibrations are derived for the OPC based on OSTIA and MUR SST. Our calibration factor (c h ) varies from 1.8 to 3 moving from the equator to ±80 • latitude degrees. This was achieved running several OPC reconstructions with constant c h varying from 1.5 to 4 and then evaluating the local improvements by means of Equation (6). Our choice maximizes the improvements of the Altimeter currents in both low-mid and high latitude areas.

Optimal Reconstruction Validation Metrics
The performances of the OSTIA, REMSS and MUR SST for the OPC reconstruction will be intercompared using the following metrics: the RMSE and Percentage of Improvement (PI) [33,47]. The PI is defined by Equation (6), In Equation (6), the subscripts U,V stand for zonal and meridional PI, respectively. The PI indicates the improvements of the Altimeter-derived currents after the optimal combination with the satellite SST data. The RMSE OPC/GEO computation relies on the knowledge of a reference surface current, provided by the SVP drifters estimates. In the 2014-2016 period, the drifters network provided 2,388,610 validation measurements over the global ocean. Their distribution is provided as Supplementary Materials ( Figure S2). Following the recommendations of RS18 and Ciani et al., 2019, the PI regional variability will be analyzed and discussed. In addition, the PI will be studied as a function of the local SST spatial gradient magnitude (|∇SST| = (∂ x SST) 2 + (∂ y SST) 2 ).

Results
In this section (as well as Sections 3.2 and 3.3) we present the PI regional variability in 20 • × 20 • boxes [33]. Such choice guarantees significant statistics in the different boxes, providing a number of validation points ranging from 100 up to 15,000 in all the areas of the global ocean [33]. The results are shown for local SST spatial gradients larger than approximately 1 × 10 −5 • C· m −1 , i.e., when the improvements brought by the optimal reconstruction become evident (see also Figure 4).

Reconstruction Based on REMSS Data
In a good agreement with RS18, the OPC computed using the REMSS data exhibit a PI exceeding 20% (up to 24% locally), with better performances for the meridional component of the motion and in the equatorial band. The results indicate that our method is more efficient when applied to low-quality altimetry products, as for the 2SAT case (see also Pascual et al., 2006 [37]). However, still in agreement with RS18, Figure 1 shows that even the 4SAT geostrophic currents do benefit from the optimal combination with SST. Indeed, local maximum improvements range from 10% at mid-latitudes to 20% in the equatorial band for both the zonal and the meridional currents. Therefore, combining the altimeter-derived currents with SST data is advantageous for both altimetry products based on information from two altimeters (like in the early altimetry era) and more recent versions combining observations from several altimeters. The number of available altimeters can indeed vary from two to four (sometimes five) in the 1993 to 2016 period [36]. On average, the PI of the zonal OPC is lower than in the meridional case. This is easily explained considering the north-south orientation of the nadir-looking altimeter ground tracks, allowing to derive more accurate zonal surface currents via the geostrophic approximation (see, e.g., [16]). However, as previously emphasized by RS18, the optimal reconstruction of zonal ocean currents exhibits poor performance in high latitudes, as indicated by negative PIs that can reach −10%.

Reconstruction Based on OSTIA Data
Here we present the OPC derivation based on the use of OSTIA SST. The results are given in Figure 2. The overall performances of the OSTIA OPC are consistent with the reconstructions obtained with the REMSS SST, as shown in Section 3.1. The regional variability of the PI is analogous, indicating that the optimal reconstruction brings larger benefits for the 2SAT case, particularly for the meridional component of the currents and in the equatorial band. However, the OSTIA SST data are smoothed to fulfill operational requirements in numerical weather prediction (https://data.noaa.gov/dataset/dataset/global-sst-sea-ice-analysis-l4-ostia-uk-met-officeglobal-0-05-daily-2013-present1). According to Equation (2), one may expect average reduced performances than in the REMSS case, mostly due to the smoothing of the spatial SST gradients. Studying the PI as a function of the local |∇SST| (Figure 4) partially confirms this expectation. When the local |∇SST| becomes larger ( 5 to 7×10 −5 • C·m −1 ), the meridional PI of the 2SAT OPC based on OSTIA SST barely exceeds 15% while, for the REMSS 2SAT case, this value can reach 23%.
When considering the zonal component, however, the averaged OPC improvements exceed 6% for the OSTIA case, compared to 4% found with the REMSS SST. An explanation for this behavior will be provided in Section 4.

Reconstruction Based on MUR Data
The OPC derived from the MUR SST are here presented for the 2SAT case. The 4SAT OPC, except for the averaged lower improvements, exhibit similar regional variability and dependence on the local SST spatial gradient magnitude (|∇SST|). The results of 4SAT OPC are provided as Supplementary Materials (Figures S3 and S4).
The MUR OPC showed reduced performances than in the REMSS and OSTIA cases. This is evident from Figures 3 and 4. The zonal degradation area found at high latitudes is even more evident than in the REMSS-based optimal reconstruction. The zonal OPC are degraded compared to the Altimeter currents below 45 • S as well as in the equatorial Pacific and Indian Oceans. On a global average, the zonal OPC degradation is also confirmed by a negative PI, whose minimum value can go down to −4%, as shown by Figure 4.
The regional and global averaged meridional PIs are lower compared to the OPC computed with the REMSS and OSTIA SSTs. Indeed, the maximum averaged meridional PI barely reaches 9%. This is verified in regions where the local |∇SST| is larger than 6.5 × 10 −5 • C· m −1 (Figure 4). However, there are regions of the ocean where the local PI can reach 10% to 20% for both components of the motion. This is found in the retroflection area of the Agulhas Current, the North Brazilian current and off Eastern Australia.

Discussion and Conclusions
The optimal combination of altimeter-derived currents and SST is a promising technique to improve past and present-day gridded altimetry products. The information contained in high-resolution SST data enables to overcome some of the limitations of the altimeter system for deriving the global ocean circulation. Such limitations are related to the altimeter along-track sampling (with a preferential north-south orientation) and the assumption of geostrophic balance, which fails in low latitudes areas. However, the corrections computed away from the equatorial band by the optimal combination method still lead to mostly geostrophic estimates. As pointed out by RS18, deriving more accurate ageostrophic motions would require the knowledge of the high-resolution source and sink term of the SST dynamical evolution equation, here approximated as a large scale SST temporal derivative. The optimal combination was tested with three global L4 SST datasets: the REMSS, OSTIA and MUR SSTs. These datasets led to different performances in deriving the global ocean OPC. The regional variability of the optimal reconstruction has been described in Section 3, while the global percentage of improvement is here discussed as a function of the local |∇SST| for the 2SAT OPC ( Figure 4). The 4SAT case is analogous and is provided as Supplementary Materials ( Figure S4).
On average, the PI depends on the magnitude of the local spatial SST gradient. For the meridional PI, all the reconstructions indicate a fairly linear increase up to |∇SST| 4×10 −5 • C· m −1 . For even larger SST gradients the PI can further increase, as for the REMSS and MUR cases or stabilize as for the OSTIA OPC.
For the zonal case, we get similar behavior in terms of the PI linearity, except a decreasing tendency in areas of highest SST gradients. Interestingly, we notice the enhanced performances of the OSTIA OPC in the 0.2 to 4 × 10 −5 • C· m −1 |∇SST| range, mostly due to improved local performances in the Southern Ocean (Figure 2). In the end, the MUR SST OPC constitute an exception, indicating averaged negative PI for the 2 to 5.6 ×10 −5 • C· m −1 |∇SST| range. This additional analysis, mostly driven by the knowledge of Equations (2) and (3) indicates how the optimal reconstruction method strongly depends on the local SST gradients. Indeed, in areas of low |∇SST| we can get little to zero improvements. Under these conditions, the optimal reconstruction method cannot correct the background surface currents, as no information related to the upper ocean circulation is available in the spatial patterns of the L4 SST field.
The anomalous behavior of the MUR SST was unexpected. Indeed, the very high nominal resolution of this product should in principle allow a more accurate description of the surface dynamical processes as also emphasized by [44,48] and references therein. In order to further investigate this point, we performed a spectral analysis of the three SST datasets used in our study ( Figure 5). Although a thorough quality assessment of the different SST products is out the scope of our paper, this analysis provides a possible explanation for the anomalous behavior of the MUR OPC. The mean spatial spectral content of the REMSS, OSTIA and MUR SST was analyzed in a test case during wintertime, focusing on scales up to 10 km. This was carried out computing the Power Spectral Density (PSD) via wavenumber spectral analysis in a land-free area of the North Atlantic, on 1 January 2016. Such spatial-temporal selection is due to the results of [49] indicating that, during winter, the submesoscale turbulence is more efficient than in summer. Hence, the expression of the smaller scale structures in SST should be favored under the chosen conditions. The spectral analyses are performed in two distinct boxes of the North Atlantic, both extending longitudinally from 60 • W to 30 • W: a first box with latitudinal extent from 30 • N to 36 • N and a second from 36 • N to 40 • N (the boxes separation is identified by the white line in Figure 5c). In both boxes, the baroclinic Rossby radius of deformation is around 30 km [50]. Investigating scales up to 10 km also enables to capture submesoscale features in the area, whenever they are present.
In the southern box (30 • N to 36 • N, see also Figure 5e), the spectra of the REMSS, OSTIA and MUR SST are initially superimposed and so they evolve until scales of 200 km. Afterward, the OSTIA spectrum begins to separate, indicating a decreased amount of spatial variance at finer scales, eventually flattening around 50 km. On the other hand, the REMSS and MUR SST exhibit lager PSD than the OSTIA SST, with a dominance of MUR until scales of 15 km, where REMSS is instead quite flat. If we consider the different nature of the three SST datasets described in Section 2, this is the expected spectral behavior of the input SSTs for our application.
If the same spectral analysis is performed in the Northern box (Figure 5c,d), i.e., where the MUR SST looks smoother, the spectral properties of the MUR and REMSS SST are inverted. Indeed, the REMSS spectrum is characterized by larger PSDs at scales lower than 100 km. Therefore, the REMSS SSTs are associated with the finer effective spatial resolution in the Northern Area (36 • N to 40 • N). This result is indicating that the global scale OPC based on MUR data can be negatively affected by unequal effective SST resolutions in different areas of the ocean. Visual inspection of the MUR SST fields evidenced that this non uniformity can happen throughout the year in different areas of the ocean. This spectral inconsistency can have drawbacks on the computation of the SST spatio-temporal gradients, in which non physical dynamical structures can be recognized. An example is provided in Figure 6, showing |∇SST| maps on three different dates in the Gulf Stream, Kuroshio Current and Agulhas Current areas. On the other hand, the REMSS and OSTIA SSTs did not evidence non-homogeneities in the effective resolution of the L4 fields, being coherent with the example provided in Figure 5 during the period 2014-2016.
However, selecting areas where the MUR SST are spatially uniform, this problem could be overcome and the fine resolution of the SST field can be fully exploited. In Figure 7, we perform a qualitative comparison between the MUR and REMSS OPC in the Agulhas Current. In the chosen area, the following three conditions are satisfied: (i) the SST spectra are comparable with the condition of Figure 5e (not shown); (ii) according to Figures 1 and 3, both the components of MUR and REMSS OPC are improved based on a three-year statistic; (iii) the trajectory of a drifter in correspondence of a SST gradient was available (see also [33]).
It is evident that the MUR SST describe a finer SST gradient than in the REMSS case. Moreover, the gridded current vectors provided by the MUR OPC exhibit a reduced angle with respect to the drifter trajectory in the north-eastern section of the domain (where the drifter time is closer to the date of the Eulerian fields, 25 March 2015).
Concluding, our study can be summarized by the following main items: • The derivation of the sea surface currents from space can benefit from the synergy between optimally interpolated space-based Earth observations from multiple platforms, today available at an operational level and with nominal daily and global coverage. Based on recent works [33,47] we attempted to optimize the surface currents retrieval by combining the altimeter-derived geostrophic currents with satellite SST from Infrared (IR) and Microwave (MW) sensors. Our study was based on three different L4 SST estimates: a dataset provided by Remote Sensing Systems, fully based on the optimal interpolation of satellite observations (IR and MW), and two additional datasets based on the optimal interpolation of satellite and in-situ data, i.e., the OSTIA and MUR SSTs.

•
The REMSS and OSTIA OPC exhibited the best performances, with maximum overall improvements equal to or larger than 15% with respect to the geostrophic estimates (for the meridional flow). This was achieved transferring the high resolution dynamical content of the satellite-derived SST into the coarser resolution geostrophic current estimates [33,34]. However, the OSTIA OPC are characterized by larger improvements than the REMSS OPC in the 0.2 to 4 ×10 −5 • C· m −1 |∇SST| range. This result is due to enhanced performances of the OSTIA OPC in the 45 • S to 70 • S latitudinal band. This is illustrated by Figures 1, 2 and 4 and summarized by Table 1: the OSTIA OPC partially solve the degradation of the zonal flow in the Southern Ocean, also yielding slightly larger averaged PIs for the meridional flow.
Most likely, this is due to the larger number of sensors used to build the L4 SSTs as well as the use of in-situ observations in the optimal interpolation procedure. This could optimize the SST estimates at high latitudes, where cloud coverage is often very high and where average intense surface winds degrade the satellite infrared and microwave SST retrievals [51,52]. An additional analysis is provided as Supplementary Materials ( Figure S5).
However, despite the improved performances of the OSTIA OPC compared to the REMSS case, occasional degradations in the polar regions can also occur with the OSTIA SST, as shown in Figure 2. This emphasizes the importance of providing very high quality, high resolution and synoptic SST measurements in those areas, highlighting the strong potential of future ESA satellite missions like the Copernicus Imaging Microwave Radiometer (CIMR) [53] (https://cimr.eu). CIMR will provide global-scale, all-weather, 15 km effective spatial resolution SST observations, also guaranteeing sub-daily coverage at latitudes higher than ±60 • . All the SST-based applications will benefit from the future CIMR remote sensing capabilities. Applying the RS18 method in high-latitude areas could also benefit from additional oceanic tracers as sea surface salinity. In polar areas, the contribution of salinity in determining the ocean dynamics can be relevant [54].

•
The OPC based on the use of MUR very high resolution data, though showing degraded performances at global scale, suggest that fine-scale SST gradients can be retrieved successfully. This is achieved when the MUR SST effective high resolution is homogenous in the study area, mostly indicating a potential for local applications. Therefore, in order to guarantee an accurate global OPC reconstructions via the method of Piterbarg 2009 and RS18, the L4 SSTs have to contain the signature of fine scale SST gradients and guarantee spatial homogeneity of the optimally interpolated field. In this way, the zonal and meridional SST gradients are not affected by the inhomogeneities like the ones reported in Figures 5 and 6. This work also identifies the RS18 method as a tool for indirect validation of the L4 SST products, allowing to test their dynamical content.  Figure S1: Uncertainties on the zonal geostrophic currents (σ u ) and on the Forcing term (h) used in Equation 3.2. Figure S2: Lagrangian buoy trajectories over the 2014-2016 period. Figure S3: Percentage of improvement of the Optimal Currents based on the MUR SST and the 4SAT geostrophic currents. Figure S4: Zonal and meridional percentages of improvement of the Optimal Currents based on the 4SAT geostrophic estimates (seen as a function of the local SST spatial gradients). Figure S5: Performances of the REMSS and OSTIA SST evaluated via in-situ SST match-up analysis.