Benefits of radar-derived surface current assimilation for South of Africa ocean circulation

Abstract

The oceanic circulation south of Africa is characterised by a complex dynamics with a strong variability due to the presence of the Agulhas current and numerous eddies. This area of interest is also the location of several natural gas fields under seafloor which are targeted for drilling and exploitation. The complex and powerful ocean currents induces significant issues for ship operations at the surface as well as under the surface for deep sea operations. Therefore, the knowledge of the state of the currents and the ability to forecast them in a realistic manners could greatly enforce the safety of various marine operation. Following this objective, an array of HF radar systems were deployed to allow a detailed knowledge of the Agulhas currents and its associated eddy activity. It is shown in this study that assimilation of HF radar allow to represent the surface circulation more realistically. Two kind of experiments have been performed, a one month analysis and nine consecutive forecast of two days each. The one month 4DVAR experiment have been compared to geostrophic currents issued from altimeters and highlight an important improvement of the geostrophic currents. Furthermore despite the restricted size of the area covered with HF radar, we show that the solution is improved almost in the whole domain, mainly upstream and downstream of the HF radar’s covered area. We also show that while benefits of the assimilation on the surface current intensity is significantly reduced during the second day of forecast, the correction in direction persists after 48 h.

Introduction

The oceanic circulation south of Africa is characterised by a complex dynamics with a strong variability due to the presence of the Agulhas current and numerous mesoscale eddies from the Mozambique Channel (Penven et al. 2006; Halo et al. 2014). More recently, high resolution modeling study by Tedesco et al. (2019) has highlighted the existence of numerous submesoscale eddies along the Agulhas cyclonic front.

Lutjeharms et al. (2003) observed the presence of cyclonic eddies embedded in the landward border of the southern Agulhas Current. These eddies have a diameter of about 50 km and are associated with a surface warm signature. Simulations suggest that those eddies remain trapped in the Agulhas Bank shelf bight and that eddies that travel downstream of the current represent leakages from the resident shear eddy. This occurs at a roughly 20 days occurrence frequency. The intensity of the meso-scale activity in this key region for the retro-flexion modulate the exchanges of heat and salt between oceans (Lutjeharms 1981; Reason et al. 2003; Van-Aken et al. 2013; Guerra et al. 2018) as well as towards the atmosphere (Messager and Stuart 2016).

This region exhibits furthermore a dynamical upwelling induced by the Agulhas Currents (Arnone et al. 2017) as observed by Goschen et al. (2015) during Natal Pulses. This upwelling, as been shown by Lutjeharms et al. (2000) to occurs on the landward side of the Agulhas Current and have an effect on the nutrient availability, stratification and primary productivity in the eastern Agulhas Bank. It as also been shown by Meyer and Niekerk (2016) that implementing an ocean current power plant in this region would outperforms onshore wind power plants and could increase the load carrying capacity of the country.

The area of interest of this paper, represented on Fig. 1 is also the location of several natural gas fields under seafloor which are targeted for drilling and exploitation. The complex and powerful ocean currents induces significant issues for ship operations at the surface as well as under the surface for deep sea operations. Strong ocean currents can also modify the height and direction of ocean waves, causing dangerous sea states (Quilfen et al. 2018). The risk of extreme waves is an important hazard for the shipping activity and off shore industry when crossing the main current systems. Therefore, knowledge of the currents state and the ability to forecast it in a realistic manners could greatly enforce the safety of various marine operations.

Following this objective an array of HF radar was deployed along the coast to allow a detailed knowledge of the Agulhas currents and its associated eddy activity. The purpose of the present document is to present and evaluate the impact of the 4DVAR assimilation of those radar data on ocean model simulation and forecast of the sea surface currents.

Data used for assimilation and validation are described in the following section. The model setup and the assimilation procedure are described in a third section while results are presented in section four and further discuss in the conclusion.

Fig. 1

Area of interest of this study. Main interest Focus on area 11b/12b and Brulpadda point. Credit: www.total.com

DATA

To monitor the variability of the Agulhas currents during offshore operations, three WERA HF radars, manufactured by Helzel Messtechnik GmbH, were installed by ACTIMAR and LWANDLE companies on the south coast of South Africa. The location of the radar system and the averaged area of measurement during April 2020 is represented on Fig. 2. The radial velocities are estimated by using the conventional method of Beam Forming with an extra filtering of the residual artefacts. Then the radial velocities are combined on a Cartesian grid at 6km resolution using the method describe by Barth et al. (2010) and made available every 30 min.

Comparisons with mobile and fixed ADCP measurements have been performed (cf. Fig. 3). For the fixed ADCP, differences intensity are observed for weak current ($\le$1.3 m/s) and a better matching is observed for stronger values. Current directions derived from radar are well correlated with ADCP measurements. The differences in intensity are explained by the low angular resolution of the BeamForming compared to the grid resolution (by a factor of about 4) at the position of the ADCP. For the mobile ADCP, the differences in intensity are lower compared to the mobile ADCP, while the differences in direction may be due to a poor calibration of the hull ADCP. Therefore, the currents provided with the Beam Forming method seems to be robust enough to be assimilated in ROMS simulations. Nevertheless, to overcome some inaccuracies, a hybrid Beam Forming/Direction Finding method developed by ACTIMAR called HYDDOA (cf. patent : FR 1562550) has been used for marine operations with better performances. Unfortunately, these data could not be used for this study.

Fig. 2

Black squares represent the emplacement of the HF radar installed to monitor the area delimited by the black contour (cf Fig. 1). The colored area represents the intensity of the averaged current measured by the radar during April 2020 and the white arrows are representative of the averaged direction of the surface currents during the same period

Fig. 3

Quantile-Quantile plot of speed and direction for comparison with fixed ADCP (a, c) and mobile ADCP (b, d)

In addition, Altimeters data were generated by a processing system including data from several altimeter missions: Sentinel-3A/B, Jason-3, HY-2A, Saral[-DP]/AltiKa, Cryosat-2, OSTM/Jason-2, Jason-1, Topex/Poseidon, Envisat, GFO, ERS-1/2 and delivered by E.U. Copernicus Marine Service Information. Being at a significant lower resolution than both model experiment those data were excluded from the assimilation process (although there are somehow assimilated in the Mercator ocean simulation used as boundary forcing) and may be considered as an independent source of observations for our validation process. Nonetheless an import bias in the representation of the Agulhas current by the altimeters data has been highlighted by Rouault et al. (2010).

Model setup & assimilation method (4DVAR)

The ocean circulation model used in this study was the Regional Oceanic Modeling System described in detail in Shchepetkin and McWilliams (2003, 2005). ROMS is a split-explicit, free-surface and terrain-following vertical coordinate oceanic model with 4DVAR capabilties (Moore et al. 2011c). Tracers momentum advection use a third order upstream biased advection scheme with no additional explicit horizontal dissipation/diffusion while on the vertical a GLS scheme is used to determine vertical mixing coefficients (Warner et al. 2005). The model grid, the atmospheric forcing, the initial and boundary conditions were all built using pyroms package freely available at http://www.myroms.org (doi:https://doi.org/10.5281/zenodo.3727272). The bottom topography is derived from Etopo1 (doi:https://doi.org/10.7289/V5C8276M). To ensure an usefull resolution in the upper ocean, 35 vertical levels with stretched s-coordinates improved double stretching function (Shchepetkin and McWilliams 2005, 2009) were used, using surface and bottom stretching parameters $\theta_s$ = 4, $\theta_b$ = 1 respectively. ROMS was initialized and forced at the lateral boundaries, by temperature, salinity and velocities profiles extracted from Mercator ocean global_analyses_forecast_phy_001_024 which provide weekly analyses and daily forecast. Atmospheric fluxes (heat and water) were extracted from ERA5 (fifth generation of ECMWF atmospheric reanalyses of the global climate—Copernicus Climate Change Service (C3S) (2017)) and introduced in the ocean model through a bulk formulae (Fairall et al. 2003). The model domain extends from $21^{\circ }\hbox { E}$ to $26^{\circ }\hbox { E}$ and from $-37^{\circ }\hbox { S}$ to $-33^{\circ }\hbox { S}$, on a 1.8 km regular grid. The ocean model has been run without data assimilation from January 2019 to April 2020 (called FREE experiment hereafter) and 4DVAR data assimilation of HF radar surface currents was performed during April 2020 only and compared to April issued from FREE.

Detailed description and evaluation of the 4DVAR data assimilation can be found in Di Lorenzo et al. (2007); Powell et al. (2008); Powell and Moore (2009); Broquet et al. (2009, 2011); Moore et al. (2011a, 2011b, 2011c); Song et al. (2016). In the present work, the dual formulation approach (Moore et al. 2011b; Gürol et al. 2014; Levin et al. 2019) has been used with 6 inner-loops and 2 outer-loops where the inner-loop correspond to the iterative linear minimization of the cost function and the outer-loop correspond to an updated non linear estimate of the circulation (for a complete description refer to Moore et al. (2011a)).

This setting has been determined after several experiments to reach an optimum between accuracy of the results and computational time. It is furthermore coherent with Levin et al. (2019) who used 7 inner-loops and 2 outer-loops in their Mid-Atlantic Bight configuration. The data were assimilated using 1-day assimilation windows. The 4D-Var analysis produced at the end of each day was used as initial condition for the next assimilation cycle. Computational cost of this choice is around 90 min by days (using 224 cpus)

Results

Hindcast

As detailed in “Model setup & assimilation method (4DVAR)” section, the 4DVAR simulation have been made for April 2020 and our analysis therefore target this periods. In the following the ROMS geostrophic velocities have been derived from surface elevation following Eq. 1 (where g represent the acceleration due to gravity, f the coriolis parameter and $\eta$ the surface elevation) for both 4DVAR and FREE simulations and are compared to altimeter derived geostrophic currents after being interpolated on the ROMS grid.

$$\begin{aligned} u=-\frac{g}{f} \frac{\partial \eta }{\partial y} v=\frac{g}{f} \frac{\partial \eta }{\partial x } \end{aligned}$$

(1)

Spatially averaged Root Mean Square Error (RMSE) between both simulations and altimeters derived geostrophic currents are shown on Fig. 4. Top panel of Fig. 4 represents the spatial average of the RMSE over the whole domain (hereafter global area), while bottom panel represents the RMSE averaged over the area corresponding to the HF radar observation zone (hereafter local area). Both panels show a significant decrease of the RMSE for the geostrophic currents with data assimilation, but while the local improvements are almost immediate it takes about 15 days to propagate them over the whole domain.

Figure 5 represents the temporal and spatial evolution of the RMSE’s differences between FREE, 4DVAR simulations during April 2020. April has been split in three periods; days 1 to 10 (panel a), days 11 to 20 (panel b) and days 21 to 30 (panel c). As expected from Fig. 4, panel (a) shows an improvement (blue color) mostly located in the HF radar area. Panel (b) shows an extension of the RMSE decrease to the south west and to the east and panel (c) a strengthening of this improvement when compared to altimetry data and highlight the positive impact of the data assimilation outside the area where HF radar data are assimilated.

While previous RMSE times series (Fig. 4) and maps (Fig. 5) illustrated a global improvement of the geostrophic currents over the whole April month, Fig. 6 focuses on a daily averaged (23 of April 2020). Panels (a) and (b) show the daily averaged geostrophic currents for FREE and 4DVAR respectively. Panels (c) and (d) show the averaged geostrophic currents derived from altimeter data and the surface currents measured by the HF radar and assimilated in the model during this day. This figure illustrates the improvement between FREE and 4DVAR experiments. Indeed, west of $24^{\circ }\hbox { E}$ the current is curving southward in FREE while 4DVAR is able to reproduce the northward bent seen by altimeters (panel (c)). Also, the area east of $24^{\circ }\hbox { E}$ and north of $36^{\circ }\hbox { S}$ is characterised by the presence of an anticyclonic eddy matching the eddy detected by altimeters. Finally 4DVAR also reproduces the intensity of the current in the southern branch of the eddy. It is furthermore interesting to note once again that beside the scarcity of the assimilated data (depicted on panel (d)) the 4DVAR assimilation is able to correct the circulation almost in the whole simulated domain.

Fig. 4

FREE and 4DVAR geostrophic velocities RMSE (against altimeter derived velocities) for the whole domain (a) and over the area of the HF radar observations (b)

Fig. 5

Maps of differences of FREE and 4DVAR geostrophic velocities RMSE (against altimeter derived velocities) during April 2020. a corresponds to the 1$^{st}$ to 10$^{th}$ of April, b corresponds to the 11$^{th}$ to the 20$^{th}$ of April and c corresponds to the 21$^{th}$ to the 30$^{th}$ of April. Blue color means RMSE reduction in 4DVAR solution. Dotted contour shows where HF radar data are assimilated

Fig. 6

Daily averaged geostrophic currents for the 23 April 2020 for: a the reference simulation, b the assimilated simulation, c from altimeters and d the assimilated HF radar data. white arrow represents current direction

Forecast

In the previous section, it has been shown that assimilating HF radar currents allows to improve the geostrophic circulation when compared to satellite-derived velocities. However, those altimeter data are at low resolution (25 km) with daily data only, while HF radar are available at 6 km resolution every 30 min. Since HF radar currents are assimilated in the model and therefore cannot be used for further validations, some forecast have been made starting from assimilated initial condition and FREE condition. This allow to validate the forecast against the HF radar data and explore the benefits of the assimilation on smaller scales and on the forecast capabilities of the current configuration.

To achieve this goal 48 h forecasts were performed during 9 consecutive days from the 21 to the 29 of April and were compared with the same forecasts initiated from FREE run. RMSE time series of both forecast against HF radar data are presented on Fig. 7. They show a strong improvement in intensity (top panel) and direction (bottom panel) when the forecast is initiated from 4DVAR simulation with a better performance of the first 24 hr of forecast (yellow curve) than the next ones (green curve).

When considering the first day of forecast, a reduction of more than 20% for both surface current intensity and direction is achieved for 87% and 99% of the time respectively. For the second day of forecast, those values change to 51% and 91% for intensity and direction respectively. When considering the averaged improvement of the current intensity RMSE, it is shown to move from 39 to 19% from the first to the second day of forecast. For directions, this RMSE improvement change from 42% to 28% from the first to the second day of forecast. Therefore, while the averaged forecast improvement is equivalent and around 40% during the first day of forecast, during the second day the averaged improvement is greater for the currents direction.

Figure 8 represents the maps of RMSE differences between forecast issued from 4DVAR and FREE for surface current intensity (Fig. 8a, b) and direction (Fig. 8c, d) along the 10 days of forecast. Right panels represent the first day of forecast and the left panels the second day of forecast. It confirms both forecast improvement of intensity and direction of the surface current. It also shows that current intensity and direction were significantly improved in the center of the area covered by the HF radar measurement. By comparing forecast day 1 and day 2, this figure also confirms that the faster degradation of the forecast is related to the current intensity rather than the direction.

Figure 9 shows the surface circulation of FREE, 4DVAR and HF radar currents averaged for two distinct first day of forecast as a superposition of arrows. It illustrates that the use of an assimilated initial condition allows to dramatically correct the path of the surface currents. Indeed on the forecast of the $21^{th}$ april 2020 (Fig. 9, left panel) while the FREE forecast is characterized by a southward deviation and an cyclonic circulation between $23^{\circ }\hbox { E}$ and $24^{\circ }\hbox { E}$, the 4DVAR forecast surface circulation is almost the opposite. Indeed it is characterized by a northward shift inducing a anti-cyclonic circulation between $22^{\circ }\hbox { E}$ and $24^{\circ }\hbox { E}$. While the whole eddy cannot be depicted by the HF radar measurement, the northward shift of the Agulhas corresponds to what it is observed by the HF radar. On the forecast of the $25^{th}$ of april 2020 (Fig. 9, right panel) while the eddy previously depicted seems to be dissipating, the northward bending of the 4DVAR forecast is once again coherent with the HF radar data and opposite to the southward bending of the FREE solution. This southward shift, seen in the FREE forecast is therefore an artefact of the model which can be corrected by using an assimilated initial condition.

Fig. 7

FREE and 4DVAR surface velocities (a) and direction (b) RMSE (against HF radar). 4DVAR-d1 and 4DVAR-d2 represent the first and the second day of forecast respectively

Fig. 8

RMSE maps between the two forecast initiated from FREE and 4DVAR and the HF radar data. Top panels (a, b) show surface currents intensities RMSE differences for the first (a) and the second (b) day of forecast. Surface currents directions RMSE differences for the first (c) and the second (d) day of forecast. Blue color means that the RMSE is reduced in the case of forecast issued from 4DVAR initial condition

Fig. 9

Daily averaged surface currents for the first day of forecast for the 21$^{th}$ and the 25$^{th}$ of April 2020. Blue arrows correspond to the FREE experiment. Red arrows correspond to the forecast and the green arrows correspond to the HF radar measured currents

Conclusion

In this study, the benefits of the 4DVAR assimilation of surface currents issued from HF radar in one of the most highly dynamic region of the world is presented. While the intense dynamics of the region make difficult for most of the numerical oceanic models to realistically reproduce the position and intensity of the Agulhas current and associated eddies, it has been shown that a 4DVAR assimilation of HF radar allow to represent the surface circulation more realistically. Two kind of experiments have been performed, a one month analyses (April 2020) and nine consecutive forecasts of 48 h each (21 to 29 of April 2020). The one month 4DVAR experiment have been compared to geostrophic currents issued from altimeters and highlight an important improvement of the geostrophic currents in 4DVAR when compared to FREE. Furthermore despite the size of the area covered with HF radar, it has been shown that the solution is improved almost in the whole domain, mainly upstream and downstream of the HF radar’s covered area.

To evaluate the forecast capability of such a configuration and to be able to use the HF radar surface currents as an independent set of data, nine consecutive forecast of 48 h each, starting from the analysed simulation (4DVAR) have been realised and compared to FREE. It has been shown that during the first day of forecast the averaged RMSE improvement is around 40% for both the surface current intensity and its direction. The second day of forecast has shown those improvements to be reduced by roughly 50% for the current and by 25% for the direction. This highlights a stronger persistence of the correction in direction than in intensity of the surface currents.

This improvement of the simulation, thanks to 4DVAR assimilation of HF radar surface currents, could have many impacts at all scales. Indeed long term reanalyses could provide better insight of the position of the Agulhas retro-flexion and the resulting ocean leakage between the Indian to the Atlantic ocean or the carbon uptake by the biological activity due to better representation of the up-welling areas and their variability. Furthermore, in this region of strong maritime activity, a realistic forecast of surface currents would increase marine safety and allow drilling campaign to be planed or suspended depending on the future surface current. Moreover although being out of the scope of this study, small scale currents having a strong impact on wave height variability (Ardhuin et al. 2017), the use of assimilated surface currents is also expected to improve wave forecast and therefore marine safety.