Accurate knowledge of ocean surface currents at high spatial and temporal resolutions is crucial for a gamut of applications. The altimeter observing system, by providing repeated global measurements of the sea surface height, has been by far the most exploited system to estimate ocean surface currents over the past 20 years. However, it neither permits the observation of currents moving away from the geostrophic balance nor is it capable of resolving the shortest spatial and temporal scales of the currents. Therefore, to overcome these limitations, in this study the ways in which the high-spatial-resolution and high-temporal-resolution information from sea surface temperature (SST) images can improve the altimeter current estimates are investigated. The method involves inverting the SST evolution equation for the velocity by prescribing the source and sink terms and employing the altimeter currents as the large-scale background flow. The method feasibility is tested using modeled data from the Mercator Ocean system. This study shows that the methodology may improve the altimeter velocities at spatial scales not resolved by the altimeter system (i.e., below 150 km) but also at larger scales, where the geostrophic equilibrium might not be the unique or dominant process of the ocean circulation. In particular, the major improvements (more than 30% on the meridional component) are obtained in the equatorial band, where the geostrophic assumption is not valid. Finally, the main issues anticipated when this method is applied using real datasets are investigated and discussed.
Ocean surface currents are vitally significant for better understanding the complex ocean–atmosphere–biosphere interactions of the earth system, and are gaining greater interest for a wide range of applications (ship routing, search and rescue, pollution monitoring, etc.). Consequently, over the last decade, the necessity for accurate estimates of spatial and temporal current variations in the open ocean and the coastal waters at high resolution (<5 km) has emerged.
While ocean surface currents can be measured from in situ platforms, the basic requirements (global coverage and high spatiotemporal resolution) can be achieved only from space.
At present, no satellite observing system is available that can provide direct observations of the ocean surface currents. However, a large variety of active and passive remote sensing instruments have been put into orbit over the last few decades that provide continuous global information regarding the ocean, from which indirect estimates of the surface currents can be drawn. This includes altimeters, scatterometers, synthetic aperture radars (SAR), imaging radiometers operating at different wavelengths (microwave, infrared), and spectrometers.
Among the remote sensing satellites cited above, altimetry, by providing global, accurate, and repetitive measurements of the sea level, has been by far the most exploited system employed in the study of the variability of ocean surface currents over the last two decades. This is because as suggested by observations and scaling arguments, flows in the ocean interior (away from the boundary layers) and away from the equator are to the first order in the geostrophic balance, implying that the ocean surface velocity field can be readily obtained from the sea surface height (SSH) gradients. The effective temporal and spatial resolutions of the altimeter maps of surface currents obtained depend upon the satellite constellation. For example, while considering a 10-day temporal resolution, the spatial scales of the order of 150 km (100 km) are resolved with two (three) satellites (Pujol et al. 2012). Recent work by Ubelmann et al. (2015) uses the potential vorticity conservation equation to propagate the SSH information forward in time and hence improve the 10-day temporal resolution of the existing altimeter maps. Spatial resolution is another issue with which to contend. The expected spatial scales of the geostrophic structures in the ocean range from 200 km at the equator to 10–20 km at high latitudes (Chelton et al. 1998), so that the mapping capability of the present altimeter constellations fails to resolve the full geostrophic flow in the mid- to high latitudes. Finally, the geostrophic currents are just one component of the ocean total surface currents: other components include the Ekman current, which is set up by the wind field, tidal current, and several other ageostrophic currents. The altimeter observing system therefore presents two major limitations for monitoring the ocean surface currents: the geostrophic component of the currents alone can be derived and only for a limited range of the spatiotemporal scales.
To obtain more realistic ocean surface currents, the Ekman currents may be estimated utilizing knowledge of the wind field (Rio and Hernandez 2003; Rio et al. 2014; Lagerloef et al. 1999) and added to the altimeter-derived geostrophic currents (Bonjean and Lagerloef 2002; Sudre et al. 2013; Rio et al. 2014). Improved surface currents may also be obtained locally by blending the altimeter velocities and the surface velocities deduced from surface drifter trajectories (Berta et al. 2015). On the other hand, over the last two decades, information on tracer distribution at the ocean surface as sea surface temperature (SST) and ocean color (OC) has increased greatly because of the availability of several imaging radiometers and spectrometers. This information is connected with the ocean surface currents in two ways. First, as tracers, the evolution in time of the SST and OC patterns is partly due to advection by the surface currents. Second, since the early 1990s numerous studies have investigated the correspondence at various spatial scales between the SST structures and the altimeter sea surface heights, and thus the geostrophic component of the surface currents (Jones et al. 1998; Leuliette and Wahr 1999; Casey and Adamec 2002). Consequently, the estimation of ocean currents with the use of tracer information, alone or combined with the altimeter SSH, has become a very dynamic research domain and different methodologies are available in the literature.
One approach, based on surface quasigeostrophy (SQG) theory (Held et al. 1995), relies on the fact that, under favorable environmental conditions, upper-layer ocean dynamics is mostly driven by the distribution of surface density anomalies (Klein et al. 2008). As such, the upper-ocean dynamics may be drawn from high-resolution SST data (and the SSH measurements) with additional climatological data of the vertical stratification (Isern-Fontanet et al. 2008). Demonstrations using existing microwave SST images have already been performed with encouraging results (Isern-Fontanet et al. 2006, 2014). Information on the direction of the geostrophic currents may also be obtained by applying the singularity analysis on the SST images (Turiel et al. 2005, 2009). Another recent approach tested on model outputs (Gaultier et al. 2014) involves identification of the optimal correction to the altimeter-derived geostrophic velocities in order to minimize the distance between the maximum lines of the finite-size Lyapunov exponents (FSLE) calculated on the corrected velocity field and the normalized SST gradients calculated from high resolution SST images.
Finally, total near-surface velocity fields may be estimated from consecutive tracer images where the tracer is assumed to be transported by the currents. This constitutes one of the most challenging inverse problems in the study of ocean circulations, and different approaches have been investigated over the last few decades. The first approach, called the maximum cross correlation (MCC) method, estimates the velocity by identifying ocean patterns on tracer images. The patterns are tracked from one image to the next, assuming that the water parcels are passively advected by the flow field (Emery et al. 1986; Kelly and Strub 1992; Bowen et al. 2002). An alternative approach is to require the velocity field to obey a known equation, such as the heat conservation equation in two dimensions, and then invert this equation for the velocity vector (Kelly 1989). The main challenge encountered is that while the cross-isotherm velocity information can be retrieved from the tracer distribution at subsequent times, the along-isotherm component cannot be directly inferred (there is no contribution to tracer advection from the along-isotherm velocity or from either component in regions of negligible gradients). Assuming the field is varying slowly, this can be handled by formulating the problem as a statistical estimation of unknown parameters (Ostrovskii and Piterbarg 1995, 1997, 2000). This assumption, however, is not valid for mesoscale studies. Another approach is to estimate the cross-isotherm velocity component using the heat conservation equation employing various constraints or regularizations (horizontal divergence, vorticity, energy) to reconstruct the two-component velocity field (Kelly 1989; Vigan et al. 2000a,b; Zavialov et al. 1998). When the tracer is advected several pixels between two observation time steps, an overdetermined set of equations may also be obtained by dividing the image scene in many subarrays (Chen et al. 2008). To handle the high-noise sensitivity of the linear methods and their limitations to retrieve large displacement motion, the nonlinear solution for the heat equation has also been proposed by various authors, such as Chen (2011). Finally, based on an innovative approach first advanced by Piterbarg (2009, hereafter PIT09), Mercatini et al. (2010), followed by Piterbarg and Ivanov (2013), proposed a method where the uncertainty of the along-isotherm velocity is removed using background velocity information. The satellite tracer information is therefore used to obtain an optimized “blended” velocity. This is a rather interesting general method that includes a tracer equation with sources and sinks besides basic advection and diffusion. Mercatini et al. (2010) successfully checked the implementation of the method for the specific application of oil spill monitoring within the framework of a twin experiment. The background velocities used were issued from a numerical model. In conclusion, the strong potential of the method for future applications targeting sea surface temperature or ocean color datasets was emphasized.
In this study, the feasibility of applying this methodology to improve the low-resolution altimeter-derived geostrophic velocities by using the high-resolution information from sea surface temperature data is investigated. After presenting the method (section 2) and the data (section 3), the feasibility of the approach is investigated based on an observing system simulation experiment (OSSE) in which the sea surface height, velocity, and temperature outputs from an ocean numerical model are used to simulate a reference ocean. The results are presented in section 4. Considering the perfect knowledge of the SST data and the source and sink terms, the scores obtained in this OSSE study provide an upper limit for the results anticipated by applying the method on real datasets. However, unlike the work by Mercatini et al. (2010), the potential of the method is tested in a realistic situation, including the use of large-scale realistic models and the prescription of the source and sink terms. Besides, the dynamics and the geographical regions most affected are investigated. The application of this method on real altimeter and SST datasets will be the topic of a future paper. However, as preparatory work, two major expected issues related to the accurate estimation of the source and sink terms in the heat conservation equation, and to the spatiotemporal characteristics of the presently available sea surface temperature images measured from space are discussed, in section 5. The main conclusions are recorded in section 6.
Temporal variations of the ocean temperatures integrated over the mixed layer depth H are governed by the heat conservation equation [Eq. (1)]. They are mainly driven by atmospheric heat fluxes, diffusion (we note κx, κy, κz as the horizontal and vertical eddy diffusivity coefficients), and advection by ocean currents ,
where Tm is the mean temperature of the mixed layer, Tb is the temperature below the mixed layer, ρ is the density of seawater, and Cp is the specific heat of the seawater.
In Eq. (1) , where Qsw is the shortwave radiation flux (+150 W m−2 on average, worldwide), Qlw is the net infrared radiation flux (−50 W m−2 average), Ql is the latent heat flux due to evaporation and condensation (−90 W m−2 average), and Qs is the sensible heat flux due to conduction (−10 W m−2 average). While Qlw, Ql, and Qs are confined to the top surface layer, Qsw may penetrate to a depth, resulting in the warming up of the ocean surface layers. The amount of shortwave radiation QPEN(H) penetrating at depth H is a function of the water turbidity.
The temporal variations of the ocean temperature integrated over the mixed layer depth (dT) depend on the horizontal (Vt) and vertical (Wt) advective terms, the diffusion term (Dt), the heat flux term (Qt), and an entrainment velocity term (Et) depending on the entrainment velocity , due to the deepening of the mixed layer over time (we is zero in the case of shoaling).
The first part of this paper considers as well-recognized terms the vertical advection, entrainment velocity, diffusion, and heat fluxes (a discussion on the impact of forcing errors on our results is made in section 5), and F is defined as a generic forcing term [Eq. (2)],
Next, we consider the sea surface temperature measured from space as an estimate of the mean temperature in the mixed layer depth (SST = Tm). The heat conservation equation for the SST is thus expressed as
where u0 and υ0 are the zonal and meridional ocean surface velocities, respectively.
The method proposed by PIT09 involves inverting the SST conservation equation [Eq. (3)] for the horizontal surface velocity. The uncertainty on the along-isotherm velocity mentioned in the introduction is reduced using background velocity information (ubck,υbck). The optimized velocities (uopt,υopt) are then obtained [Eqs. (4) and (5)] as a function of the background velocities F and the spatial and temporal derivatives of the SST,
The approach is easily explained using geometrical considerations (see Fig. 1): A local reference frame (sst) linked to the SST gradient is defined by placing the xsst axis in the direction of the SST front and the υsst axis perpendicular to the SST front. In this rotated frame we get and , so that the across-front velocity υsst is readily obtained from Eq. (3),
To test the feasibility of applying this methodology to improve altimeter-derived velocities using high-resolution SST data, we employed an OSSE approach, where a reference ocean state is given by the outputs from a high-resolution () ocean numerical model. The modeled sea surface heights are used to produce simulated maps of the altimeter geostrophic velocities (calculation details available in section 3b) that are used as background velocities. The tracer field is the high-resolution modeled sea surface temperatures (SSTm1/12) and F is derived from the heat conservation equation
where are the high-resolution model surface velocities. The forcing field F therefore includes the heat fluxes, vertical advection, vertical entrainment, and diffusion.
a. Model data
In this study outputs from two different configurations (one free run and the other with data assimilation) of the Mercator Ocean system have been used. Both configurations employed the ORCA12 model (Le Galloudec et al. 2008; Hurlburt et al. 2009) with a homothetic tripolar ORCA type of grid (Madec and Imbard 1996; Timmermann et al. 2005). The horizontal resolution is 9.25 km at the equator, 7 km at the midlatitudes, and 1.8 km in the Ross and Weddell Seas. The ocean/ice code used is based on the Nucleus for European Modelling of the Ocean (NEMO) code (Madec 2008) developed and maintained in the framework of the European “NEMO consortium” (CNRS, Mercator Ocean, NERC, Met Office, INGV, and CMCC). The forcing fluxes come from ERA-Interim (Dee et al. 2011). The solar penetration is computed with the two bands’ extension algorithm (Jerlov 1968) and depends on the chlorophyll-a concentration monthly climatology drawn from the SeaWiFS data.
The OSSE results described in section 4 are based on one year (2013) of sea surface height, temperature, and velocity outputs from the resolution operational data assimilating system (Dombrowsky et al. 2009; Drévillon et al. 2008; Drillet et al. 2008). This system assimilates both in situ temperature and salinity profiles from the Coriolis data center; remotely sensed data, that is, sea level anomalies (SLA) calculated at Collecte Localisation Satellites (CLS) in the framework of the SSALTO/DUACS project and distributed by Archiving, Validation, and Interpretation of Satellite Oceanographic Data (AVISO); and SST data from NCEP [real-time global sea surface temperature (RTG_SST) data]. The assimilation system is a reduced-order Kalman filter using a 3D multivariate modal decomposition of the forecast error covariance (Tranchant et al. 2008). This simulation has been forced by the operational analyses of the ECMWF.
In section 5 the relative contribution of the different terms in the heat conservation equation based on model outputs is discussed. Data assimilating models are not adapted for such analyses because a nonphysical adjustment is done at each time step of the model to minimize the misfit between the model analysis and the observations. For this specific analysis, we have therefore used temperature and velocity outputs from a simulation without data assimilation for two time steps (1 and 2 August 2002). Finally, it must be noted that although it is defined on a resolution grid, the effective resolution of the two models used in this study may be considered to be around 5–10 times that value (Soufflet et al. 2016).
b. Simulated altimeter maps
To evaluate the potential impact of the method in order to improve the altimeter-derived geostrophic velocities, maps of the altimeter surface geostrophic currents were simulated from the modeled SSH using the identical procedure applied to generate the AVISO altimeter surface velocity maps calculated at CLS. Over the year in consideration (2013), four different altimeter satellites were available: Jason-2, CryoSat-2, Jason-1 (until 21 June 2013), and the Ka-band Altimeter (AltiKa; from 14 March 2013). First, the model SLA was calculated from the model SSH by subtracting the 2013 model mean SSH. The SLA grids thus obtained were then interpolated along the different altimeter satellite tracks. They were further filtered and subsampled in accordance with Pujol et al. (2016). The filtering lengths ranged from 120 km at the equator to 30 km at the pole and one point was kept over two (i.e., a 14-km distance between two along-track measurements). Then, these simulated along-track SLA observations were used to daily generate global ¼° maps of sea level anomalies using the optimal interpolation (OI) method. The OI procedure applied here has the same parameters as those used to construct the AVISO altimeter maps (Pujol et al. 2016)—that is, the same estimate of the altimetric noise budget, which varies depending on the mission; and the same correlation scales, which vary depending on latitude, both temporally (10–40 days) and spatially (100–500 km). Then, the 2013 model mean SSH was added back to the gridded SLA to obtain the SSH grids. Finally, again following the same procedure used at CLS, the geostrophic velocities were calculated using the nine-point stencil width method (Arbic et al. 2012). At the equator a beta-plane approximation is used, as described by Lagerloef et al. (1999).
a. One-year method skill statistics
In this section, the results of the OSSE study are reported. The PIT09 method is applied using as background the simulated altimeter geostrophic velocities, as described in section 3b. One year (2013) of optimal velocities are thus obtained and compared with the reference velocities (the model velocities). The efficiency of the method is tested by calculating a score δ [Eq. (8)] as proposed by PIT09:
Improvements (%) are calculated as 100(1 − δ) and are presented in Fig. 2 for the zonal (left panel) and meridional (right panel) components of the velocity.
Improvements are clearly obtained everywhere, with values ranging between 10% and 40% (the global mean improvement is 21.5% for the zonal component of the velocities and 22.5% for the meridional component). These values are somewhat less than the theoretical 30% improvement calculated by PIT09 under idealistic assumptions, but they concur well with the values obtained with realistic flows (12%–32%). As they are based on model results, improvements are also obtained in the seasonal sea ice zone at high latitude. It is noteworthy that in these areas, the data coverage will be limited in real datasets of both the SST and altimetric currents. A large spatial variability in improvement is observed in Fig. 2. Stronger values are obtained in areas where the greater error is made in deriving the total surface velocities from the altimeter data, giving the method more leeway for improvement. This is the case, for instance, in the equatorial area, where greater errors arise due to the failure of the geostrophic approximation. Also, stronger improvements are observed in the tropics when compared with the higher latitudes, which may be explained by the larger spatial correlation scales used (150–200 km) compared with higher latitudes (60 km) in the altimeter mapping procedure (Pujol et al. 2016). Consequently, in the tropics, the spatial scales resolved by the simulated altimeter maps are larger than at the higher latitudes, implying that the short-scale SST may provide more significant information. Another driver of the observed spatial variability is the SST front direction. For instance, the least improvement (around 10%) is obtained for the meridional component of the velocity in the strong upwelling systems west of the African and South American continents. The SST gradients in these areas being strongly zonal, , and υopt ≈ υbck in [Eq. (5)], so that only minimal information is brought in by the SST on the meridional velocity component. Inversely, it is in these areas that the strongest zonal improvements are obtained (30%–40%). In all strong zonal currents [Gulf Stream, Kuroshio, Antarctic Circumpolar Current (ACC), and equatorial current] the opposite is observed. Here and uopt ≈ ubck [Eq. (4)], implying that little information is generated by the SST on the zonal velocity component and that the improvement obtained is lesser on the zonal component (10%–20%) than on the meridional component (20%–30%) of the velocity. More generally, stronger improvements are usually observed on the meridional component than on the zonal component (Fig. 2). Greater uncertainty on the meridional component of the geostrophic velocity compared with the zonal velocity is anticipated because of the meridional orientation of the altimeter tracks (the across track; i.e., zonal geostrophic velocities are better resolved). Consequently, the mapping of the meridional component of the geostrophic velocity is less accurate than the mapping of the zonal component by 10%–20% in all latitudes (Le Traon and Dibarboure 1999).
The improvements calculated also feature strong seasonal variability. Hemispheric monthly averages (Fig. 3) reveal higher improvements in winter than in the summer for the zonal component of the velocity and the opposite for the meridional component. The regional variability of this seasonal cycle (Fig. 4 shows the mean improvements obtained in 2° boxes for the summer and winter months) clearly reveals that the strong regional seasonal cycle observed in the tropics and upwelling systems dominate the mean values. In other areas, a reversed seasonal cycle may be observed. Two main opposite mechanisms explain these regional behaviors.
The first driver is that in winter, when the ocean surface waters are well mixed, the surface temperature is much more correlated to the SSH (Jones et al. 1998) and to the corresponding background geostrophic velocities that respond to depth-integrated currents. In areas where the surface velocities are mostly geostrophic (in the center of the basins, for example), the background velocities are then almost similar to the reference velocities such that the SST information conveys relatively little additional information (i.e., the improvement of the method decreases). Accordingly, Fig. 4 shows that higher improvement is obtained during the summer in the zonal component in the strong western boundary currents (Gulf Stream, Kurshio) and in the latitudinal band between 30° and 40° (both hemispheres). Similarly, greater improvement is obtained in the summer in the meridional component of the velocity in the midlatitudes. The second opposite driver of these regional seasonal variations is the wind intensity seasonal variability. This effect is most often visible in the tropical band and in the ACC, where the stronger winds during winter produce stronger Ekman velocities, an ageostrophic component of the surface currents not measured by altimetry. Therefore, higher improvement is obtained in the winter for both components of the velocity, in these two areas. In the southeastern part of the Pacific subtropical gyre (latitudes between 60° and 40°S, longitudes between 80° and 160°W) the improvement is increased in winter, a phenomenon that is not due to stronger Ekman currents. This is more due to a drop in the correlation between the SSH and SST (from 0.8 in summer to 0 in winter, not shown), probably because of strong wind-driven barotropic currents.
b. Snapshot analysis in the Benguela upwelling system
In Fig. 5 an example of the efficacy of the method on 16 February 2013 in the Benguela upwelling system is seen. In this area the winds very often blow northward along the West African coast, generating strong westward-flowing Ekman currents (visible on the reference velocity image shown in Fig. 5c) that are not measured by the simulated altimeter data (Fig. 5a). The restitution of these currents is of prime importance in this area, as they induce the upwelling of the underlying colder (see the corresponding SST field on the same day in Fig. 5d) and nutrient-rich waters that sustain the very active Benguela ecosystem. The optimal velocities obtained are presented in Fig. 5b. It is clear that the method enables the partial recovery of the ageostrophic component of the surface current (mostly visible at the top part of the image, for latitudes equatorward of 30°S).
c. Snapshot analysis in the tropical Atlantic
In the eastern part of both the Atlantic and Pacific Oceans, water upwelling from the ocean depths generates a tongue of cold surface water that usually extends westward along the equator from the continents. These “cold tongues” are surrounded by warmer surface water in both hemispheres.
The velocity shear produced between the westward-flowing South Equatorial Current and the eastward-flowing Equatorial Undercurrent and Equatorial Countercurrent causes the cold tongues to deflect to the north and south, resulting in “tropical instability waves” characterized by a wavelength of about 1100 km and a propagation period of about 30 days (Von Schuckmann et al. 2008).
The surface structure of the tropical instability waves is most clearly evident in the satellite pictures of sea surface temperature and ocean color. A good example is seen in Fig. 6a, which shows the model sea surface temperature in the tropical Atlantic Ocean on 27 July 2013. The advection of the colder SST patterns toward the warmer, higher latitudes by the ocean surface currents associated with the tropical instability waves is clearly visible in this plot. Three main waves are visible, one centered at 16°W and 4°N, another at 25°W and 4°N, and the third at 35°W and 4°N. The speed of the reference field is shown in Fig. 6b. These waves are not resolved by the simulated altimeter data (Fig. 6c). The optimized velocity map (Fig. 6d) reveals the signature of the three waves. Although the intensity is lower than the reference velocities, the method successfully manages to reproduce a portion of these significant surface structures of the ocean.
An even more striking example is seen in Fig. 7 for the tropical Atlantic Ocean on 1 January 2013. The velocity maps clearly highlight strong differences between the background (Fig. 7a) and the reference velocities (Fig. 7c). In particular, centered between the longitudes of 28°W, 24°W and latitudes between 0°N and 5°N, the reference field features the signature of a tropical instability wave that is not resolved at all in the background field (where the flow is mainly zonal) but partially recovered in the optimized velocities. The background velocity field also appears to be much smoother than the reference velocity field, as is evident from the zonal power density spectra (Fig. 8): The simulated altimeter map resolution (background, in black) and the model resolution (“truth,” in orange) velocity spectra separate at scales ~150 km, concurring well with the expected resolution of the simulated AVISO maps in the area (correlation scales used for the mapping being of the order of 150–200 km). The velocity spectra of the optimized field follow the velocity spectra of the reference field until spatial scales of around 70 km and then continue with the same slope, while the drop in the power density of the reference field is much steeper. This highlights the presence of spurious short-scale variability in the optimized velocities. This short-scale noise is also visible in the examples shown in Figs. 5–7. This suggests that further filtering may be required in the optimized field. We calculated different filtered solutions using a low-pass filter with scales increasing in length from 10 to 50 km and further calculated the log-spectral distance between the reference field spectra (SR) and the optimized filtered spectra SO,
In the tropical Atlantic area as shown in Fig. 7, the minimal spectral distance was obtained for a 30-km filter length (red plain line spectra in Fig. 8). The two red dashed lines in Fig. 8 show the spectra obtained using a 20-km filter length (dashed line above the plain line, too noisy) and a 40-km filter length (dashed line below the plain line, too smooth). The same procedure was applied to the other parts of the ocean, focusing either on the high-variability areas (Kuroshio, Gulf Stream, Agulhas Current) or low-variability areas (subtropical gyres in the Atlantic, Pacific, and Indian Oceans). In all these cases, the minimal spectral distances were obtained using filtering lengths ranging between 20 and 30 km. We are therefore confident that this method can be effectively used to add a significant signal to the altimeter velocity maps in the range of 30 to 100–200 km (the effective resolution of the altimeter maps, depending on the geographical area and the altimeter constellation).
The efficacy of the method may be also characterized via calculation of the spectral coherence (Fig. 9). The spectral coherence between the background zonal velocities and the reference zonal velocities, close to 0.6 for very large scales, rapidly decreases to below 0.2 for spatial scales shorter than 500 km, whereas much higher coherence values (up to 0.6) are obtained for scales larger than 60–70 km, between the optimal zonal velocities and the reference zonal velocities. This indicates that the methodology can enhance the altimeter velocities at spatial scales not resolved by the altimeter system (i.e., below 150 km) but also for the larger scales, where the geostrophic equilibrium might not be the unique or dominant process of the ocean circulation. However, at scales shorter than 60–70 km, although the coherence is improved, the values remain low. Slightly higher coherences are obtained for the short scales using the 30-km low-pass filtered optimal velocities.
5. Discussion: Applicability of the method on real datasets
The OSSE results presented in the previous section have highlighted the potential of the methodology to improve the altimeter geostrophic velocities by adding short-scale and ageostrophic information to the large-scale geostrophic signal. This demonstration, however, is based on the a priori perfect information on both the high-resolution sea surface temperature and F of the heat conservation equation [Eq. (3)]. The next step of our work will be to apply the methodology using real observations. Although the practical application of the method on real observations is outside the scope of this paper, two main issues are anticipated that are discussed in this section.
a. F in the heat conservation equation
The large uncertainty on the forcing term is expected to be the more important issue to handle for this methodology to be applied to real datasets. To check the sensitivity of our results to forcing errors, the calculation of the optimized velocities was repeated by adding to the true forcing [Eq. (7)] a white noise of increasing RMS (from 10% to 50% of the forcing signal). In Fig. 10 it is clear that the improvements obtained (calculated for the global ocean) decrease as the noise intensity increases, with the meridional component of the velocity being more impacted than the zonal component. The improvement for the zonal component is 22% for the perfectly known forcing, which drops to 15% when 20% white noise is added, and to 7% with the addition of 30% white noise. This shows that error forcing should be less than 20% of the forcing intensity for the results to remain significant. Besides, if the noise is too large, then the method can actually degrade the background currents (negative improvements in Fig. 10).
As discussed in section 2, the forcing term required for our method is the sum of the different contributing terms [Eq. (2)], that is, atmospheric heat fluxes, vertical advection, entrainment velocity, and diffusion. To understand the relative contribution of these three terms and to decide whether any of them can be discarded for our application, we calculated them using the ORCA12 free-run simulation described in section 3.
We first checked that the heat budget is effectively closed in the model. The terms dT, Qt, Vt, Wt, and Dt in Eq. (1) were calculated together with a residual field R = dT − Qt − Vt − Wt − Dt using the hourly outputs from the Mercator Ocean free-run model configuration (the model time step is 6 min). The mixed layer depth may be calculated using different criteria. Here we discuss the results obtained using the model turbocline depth MLDturb, which corresponds to the depth at which the vertical diffusivity coefficient drops below 5 cm s−2 (R. Bourdalle-Badie 2015, unpublished manuscript). Other mixed layer depth (MLD) criteria have been checked, all of which resulted in higher residuals. The heat fluxes are those used to force the model, and QPEN(H) is calculated using the following formulation (Murtuggude et al. 2002):
where Chl is taken from the SeaWiFS chlorophyll-a concentration monthly climatology.
Because of the insufficient vertical spatial sampling of the model, we were unable to calculate accurate entrainment velocity values from the hourly MLD values. Therefore, this term has been neglected in the first approximation. It may, however, become significant (exceeding Qt) during fall at latitudes polarward of 40°N (Alexander et al. 2000). The residual field obtained is shown in Fig. 11 (bottom plot on the right) along with the other terms, dT, Qt, Vt, Wt, and Dt. It features both large-scale and small-scale structures that are small when compared with dT. Apart from the missing entrainment term, the discrepancies that arise are mainly due to nonlinearities in the temperature evolution equation and to the approximate evaluation of QPEN, which is linked to the mixed layer estimate, the Chl-a climatology, and the choice of the empirical formulation [Eq. (10)]. The results in Fig. 11 clearly identify Qt and Vt as the two main contributors of the temporal variations of the temperature in the mixed layer. The values of Wt and Dt are negligible compared with these two main components, and we might consider neglecting them when applying the methodology to real observations. On the other hand, the accurate knowledge of the QNET, QPEN, and Et terms will be essential, together with an appropriate choice of the depth H over which the heat budget is integrated. All of these items pose great challenges.
First, the different terms in QNET are very difficult to measure directly and are mostly calculated using bulk formulas. Various products are available based on numerical weather prediction reanalyses, ship-based analyses, satellite-based analyses, or blended analyses. However, at present, the global imbalance in the gridded products available ranges between 2 and 30 Wm−2, which is greater than the accuracy required for many applications (Yu et al. 2013).
On the other hand, an accurate estimation of QPEN is also very difficult. As it is turbidity dependent, it may be calculated using the information from ocean color images, together with an empirical formulation [e.g., the one given by Eq. (10)].
Finally, dT is the time difference in the temperature averaged over the mixed layer depth. Sea surface temperature measured from space corresponds to different depths depending on the radiometer operating wavelength. Precise definitions have been given by the Group for High Resolution Sea Surface Temperature (GHRSST). It is about 10–20 μm for infrared radiometers (the skin temperature) and about 1 mm for microwave radiometers (the subskin temperature). On the other hand, most gridded SST products resulting from the combination of both microwave and infrared data usually correspond to the foundation temperature (free from diurnal temperature variability, ~10 m). Therefore, care must be taken when using the observed SST products because the SST measured may differ from the mean mixed layer depth temperature.
The three issues listed above highlight the difficulties that will be encountered while attempting to apply the methodology using observational datasets. The cumulative uncertainties will significantly reduce the accuracy of the resulting optimal velocities.
However, a striking feature in Fig. 11 is the relative spatial scales of the Vt and QT terms, the latter exhibiting much larger scales than the former. Apart from the clear impact visible in Fig. 11 of the diurnal variation in the shortwave flux (due to the 1-h time step choice), this spatial-scale difference is mainly a result of the resolution of the forcing fields and atmospheric synoptic scales (Rossby radius in the atmosphere being much greater than the ocean one). This characteristic may enable us to circumvent the difficulties mentioned earlier. One solution could be to use statistical image analysis techniques applied to successive SST images, in order to separate the shorter spatial scales due to advection from the larger spatial scales due to atmospheric forcing. One could also benefit from the different time scales of the two terms (the atmospheric forcing temporal scales being much smaller than the temporal scales of the advection terms). This approach will avoid the use of additional and presently inaccurate observations.
b. Investigating the required spatiotemporal resolution of the observations
The second issue of applying the methodology to real SST observations is related to the spatiotemporal scales of the observed SST products available. As evident in the previous section, using the hourly fields of -resolution temperature, the heat budget equation integrated over the mixed layer depth is almost closed. A similar calculation performed with a daily time step leads to a much stronger short-scale residual signal (Fig. 12), mainly in the strong current areas. This additional signal is due to the improper calculation of the horizontal advection term. In one day a water parcel will travel many grid cells of a grid; therefore, calculating the SST differences due to advection will not match the SST time difference. For advective signals to be correctly resolved, the Courant–Friedrichs–Lewy (CFL) condition must yield
where U is the current speed, Δt is the temporal resolution, and Δx is the spatial resolution. In this section we discuss the required temporal resolution for a given spatial resolution, or inversely the required spatial resolution for a given temporal resolution. We approximate the velocity term U with the mean geostrophic velocities from the CNES-CLS13 mean dynamic topography solution (Rio et al. 2014) and obtain the following CFL conditions:
For Δx = 10 km (left plot of Fig. 13), Δt will be shorter than one day in strong currents (equatorial currents, Antarctic Circumpolar Currents, western boundary currents…). For Δx = 25 km (right plot of Fig. 13), the 1-day resolution is sufficient almost everywhere (white everywhere, that is, Δt > 24 h, except in the core of the strong currents).
If observations are available with a temporal resolution Δt = 1 day (Fig. 14), then the spatial resolution will be greater than 25 km in strong currents (up to 70–80 km in the core of the western boundary currents).
These results illustrate the fact that, when this methodology is applied to the real SST dataset, a careful examination of the effective spatial and temporal resolutions of the available products will be crucial and that different products may be selected depending on the oceanic regime where the method is to be applied: For low-variability areas, daily SST maps may be used (in that case, infrared images from radiometers embarked on polar-orbiting satellites are suitable), while for the highly variable areas daily SST maps are useful using only ¼° spatial resolution products (this corresponds to the spatiotemporal resolution of microwave SST–based products). Alternatively, in the highly variable areas, 3–6 hourly, spatial resolution SST may be used. This high spatiotemporal resolution may be provided by infrared radiometers embarked on geostationary satellites.
In the case of infrared SST images, however, another major issue is the potentially large amount of missing data due to cloud coverage. All of these issues will have to be carefully considered when applying the methodology to real SST datasets.
The objective of this work was to investigate the feasibility of improving the low-resolution altimeter geostrophic velocities using higher-resolution information from the sea surface temperature data by applying the method first proposed by PIT09. Using an OSSE approach based on the high-resolution model outputs from the Mercator Ocean system, we have simulated altimeter observations of the surface geostrophic currents and demonstrated that the method is capable of enhancing this background information everywhere with a skill of up to 30%–35%. It is quite interesting to note that major improvements are obtained in the equatorial area where the use of the altimeter data to calculate ocean surface currents is very limited due to the failure of the geostrophic approximation. Spectral analyses have revealed that improvements may be obtained on both large scales, by retrieving information regarding the missing ageostrophic current components (as for instance the Ekman current component), and smaller scales, totally unresolved by conventional altimetry. In the future context of high-resolution altimetry provided by the Surface Water and Ocean Topography (SWOT) satellite, such an approach will be fully relevant in helping to refine the temporal resolution of the velocity fields derived and to provide further information on the missing ageostrophic surface currents. However, the application of this methodology on the observed SST dataset is challenging, as discussed extensively in the last section of this paper. Several issues remain to be tackled, including the extraction of the short spatial scales from the SST images due to current advection and the careful selection of the SST images among the various products available in order to optimize between spatiotemporal resolution and spatial coverage. The application of this method on real datasets will be the subject of a future paper.
This study was performed within the framework of a Marie-Curie Fellowship cofunded by the European Union under the FP7-PEOPLE–Co-funding of Regional, National and International Programmes Grant Agreement 600407 and the RITMARE Flag Project. We thank the two anonymous reviewers for their very constructive comments.