NEW WAVE HINDCAST FOR THE RIO DE LA PLATA ESTUARY

This paper present a wave hindcast for the Río de la Plata Estuary. It was based on the version 5.16 of the numerical wave model WAVEWATCH III ®, forcing with CFSR winds and considering non-stationary water levels and currents obtained from the hydrodynamic model TELEMAC 2D. A multi-mission altimetry database processed by IFREMER was considered as the ground truth. It was taken as reference to calibrate and validate the wave model, and also to do a local validation of CFSR winds. Tuning the Γ coefficient of the JONSWAP parametrization of dissipation by bottom friction (Sbot) and the BETAMAX coefficient of the energy input by wind (Sin), it was possible to correct the large negative BIAS in significant wave height (Hs) obtained with the default configuration. The hindcast cover the period 1985–2016 providing wave parameters with a spatial resolution of 1 km and a time step of 1 hour and full spectra at a 10 km offshore distance along the coast.


INTRODUCTION
The Río de la Plata (RDP) is a large estuary formed by the confluence of Paraná and Uruguay rivers that discharge into the Atlantic Ocean.As shown in Figure 1, it is 290 km long and has a NW-SE orientation.An imaginary line that connects Montevideo and Punta Piedras is the upstream limit of the outer zone of the estuary.This zone is wide (O(200 km)), and water depth varies between 10 and 20 m.The wave conditions in this zone are a combination of wind-seas and swells that come from the south and dissipate as they travel into the estuary.On the other hand, the intermediate and inner zone is narrower (O(50 km)) and shallower (O(5 m)) and wave-conditions are dominated by sea-winds.
As a large and shallow water body, wave-bottom interaction processes are relevant and the composition of the bottom become an important fact.Fine sediment carried by the two large tributaries mostly covers the bottom conforming different patches that differ from one another according to the proportion of clay, silt and sand (Wells and Daborn, 1998).Hence, it is a heterogeneous bottom where fine sediment predominates.The availability of reliable wave data is useful for multiple activities.Currently, the state of the art of wave modeling, the computing power and the availability of reliable wind data from atmospheric reanalysis (e.g.CFSR), make it possible to improve wave information at any place by performing a wave hindcast.However, its implementation for the RDP has many difficulties.Most of them are common for any estuary where water level variations, currents and bottom composition strongly influence wave transformation processes.This is evident in Alonso et al. (2015), where the results of a high-resolution wave hindcast for Uruguayan waters show good agreement with measurements on the Atlantic region but not so good agreement for the RDP.Particularly, a systematic underestimation of significant wave height (Hs) is observed at the intermediate and inner region of the estuary.
For that reason, a new wave hindcast for the RDP was developed and it is presented in this paper as follows.On Materials and Methods section, the wave model configuration is presented and the input data as well as the altimetry database used as ground-truth are described.This section finish with the presentation of the methodologies used to validate input winds and to calibrate the wave model.After that, the obtained results are shown and discussed, to end with the main conclusions and the identification of necessaries future steps to keep improving wave data for this part of the globe.

Wave model configuration
The model was implemented on the wave modeling framework WAVEWATCH III ® version 5.16 (WWDG, 2016).The parametrizations named ST4 was considered.It incorporates the latest results of the physics of energy input by wind and whitecapping dissipation (Sin+Sds), that are exposed on Arduhin et al. (2009 and2010), Leckler et al. (2013) and Rascle and Arduhin (2013).The choice of ST4 was based on the results of Stopa et al. (2016), where for different global hindcast, the one that use ST4 has the best performance on the western side of the south Atlantic Ocean.
Regarding to dissipation by bottom friction (Sbot), for simplicity and as an initial step, the linear JONSWAP parametrization (Hasselmann et al. 1973) was used.
Five grids were used on a multi-grid two-way nesting mode.Starting with a coarse global grid (Grid 1) and reaching the Uruguayan coast and the intermediate and inner RDP (Grid 5) with a resolution of 40"•(~1km).The main characteristics of the grids can be found on Table 1 and Figure 2. The spectral discretization is the same for whole grids.Spectra were discretized in 36 directions uniformly distributed and 25 frequencies distributed as a logarithmic grid covering the range 0.0418 Hz -0.4114 Hz.

Model inputs
Bathymetry.ETOPO1 (Amante and Eatkins, 2009) was used on grids 1 and 2, while local and more detailed nautical charts were used on grids 4 and 5.In the case of grid 3, the nautical charts do not cover all the domain, so there were complemented with ETOPO1.
Shoreline.It is important for the definition of the sub-grid obstacles, such as islands.The GSHHG data set (Wessel and Smith, 1996) was used on grids 1, 2 and 3.The high resolution data set was used on grid 1 while the full resolution one was used on grids 2 and 3.For grids 4 and 5, the shoreline was obtained from the nautical charts.
Wind.The same wind data were used for all the grids.It was obtained from the atmospheric reanalysis CFSR (Saha et al. 2010) and its extension CFSv2 (Saha et. al. 2014).The data were interpolated to a regular grid with 0.3125º x 0.3122º spatial resolution and 1h time step.
Non-stationary water levels and currents.There were obtained from the hydrodynamic model TELEMAC 3D (Hervouet, 2007) implemented for the RDP domain taking into account fluvial discharges of Paraná and Uruguay rivers, tides at the oceanic boundary (astronomical and meteorological from a regional model), and wind and sea level pressure from CFSR reanalysis.The model was calibrated and validated using in situ measured data.The results show good agreement with the measured data, satisfactorily reproducing the main features of the RDP dynamics (see Santoro et al. 2016).

Altimetry database.
It is a multi-mission database processed in IFREMER (Queffeulou and Croizé-Fillon, 2013).It covers the period 1991-2013.The Figure 3 shown the spatial distribution of the data for the area between latitudes 38 S -32 S, and longitudes 58W -52W.It was considered the reference data, at first, to validate the CFSR wind, and then to calibrate the wave model.

CFSR winds validation
Each altimetry data was paired with one from CFSR.To obtain a CFSR value for a specific location and time, bi-lineal interpolation on space and lineal interpolation on time was done.For each of the sub-regions of Figure 3, standard error metrics were calculated.In addition, 20 quantiles equispaced on a Gumbel scale were compared superimposed on dispersion diagrams.

Wave model calibration
The year 2002 was chosen as the calibration period, because it has more altimetry data than an average year and also more extreme events than an average year.To define an extreme event, the 95 th percentile was considered.
Firstly, a simulation of 4 month (i.e.January 2002-April 2002) with the default parameter were performed.Based on the results of the default simulation, the parameter BETAMAX of the ST4 energy input by wind parametrization and the parameter Γ of the JONSWAP parametrization of dissipation by bottom friction were selected as the calibration parameters.
The nine simulations described on Table 2 were considered.The objective was make the BIAS of Hs tent to zero and decrease the scatter index (SI) for the four sub-regions defined on Figure 4. Based on the results of these nine simulations, the value of one of the two calibration parameters was fixed and then, more simulations were done tuning the parameter that was left undefined.

CFSR winds validation
The difference between the wind velocity at 10 m height (U10) obtained from CFSR reanalysis and the obtained from altimetry radars is shown on Figure 5 on terms of BIAS, root mean square error (RMSE) and correlation coefficient.Dispersion diagrams superimposed with q-q plots corresponding to the regions of Figure 4 are presented on Figure 6.
From a general point of view, CFSR U10 shows a good performance on the whole study area with a small positive BIAS and RMSE around 2 m/s.As expected, the worst performance occur on the RDP, particularly on the inner zone.Special attention was paid to the underestimation observed on the higher quantiles of the intermediate and inner RDP (top-left plot of Figure 6), arriving to the conclusion that it is due to only a few events and different causes (i.e.time lag of the peaks, not feasible altimetry data, meteorological phenomena of smaller scale than the reanalysis, and events not captured by the reanalysis) were identified.Hence, no kind of correction was attempted and CFSR U10 was used directly.

Wave model calibration
The results obtained with the default configuration of the model are sown on Figure 7.It is observed a small negative BIAS on the Atlantic and a greater one on the RDP.Hence, the strategy was increase Hs from BETAMAX increase combined with ||Γ|| decrease.The obtained results in terms of BIAS and SI for each region defined on Figure 4 are presented on tables 3 to 6.The results are much more sensitive to Γ than to BETAMAX.Particularly in the RDP, the results evidence the relevance of wave-bottom interaction processes.BETAMAX was set on 1.55 because a better fit on high quantiles in the Atlantic were observed.Result consistent with Pereira et al. (2017).After that, Γ was continue tuned reaching the best results for Γ=-0.012m 2 s -3 .These results can be appreciated on Figure 8.A good agreement between modeled and altimetry Hs is observed.The main persistent problem is an overestimation of Hs for the highest waves in the intermediate and inner RDP.It is concluded that with JONSWAP parametrization is not possible to simulate properly moderate and high waves with the same Γ.Increasing || Γ|| it is possible to correct the BIAS of higher waves, but it will induce a negative BIAS on most frequent waves.Hence it is considered that BETAMAX = 1.55 and Γ=-0.012m 2 s -3 is the best configuration using JONSWAP parametrization and CFSR winds.
Finally, the calibrated model was used to wave hindcasting the period 1985-2016, providing a database of wave parameters with a spatial resolution of 1 km and a time step of 1 hour and full spectra at a 10 km offshore distance along the coast, with a distance between them also of 10 km approximately.

CONCLUSIONS AND FUTURE WORK
A new wave hindcast for Uruguay was developed.It has a better performance than the previous one, particularly for the Río de la Plata.The improvement was obtained from better winds and bathymetry, from a higher model resolution, from better physics of source terms, but mainly from including non-stationary water levels and from considering the dissipation by bottom friction to calibrate the model.
The results clearly show the key role of bottom friction on wave dissipation, particularly on the intermediate and inner zone of the estuary.Due to its shallowness, water level variation due to tides and surges induce significant changes on water depths, so consider non-stationary water levels is a first mandatory step to proper model wave-bottom interaction in this zone.Calibrate the parameter Γ was the other step that was taken, and it is as far as it can be reached using JONSWAP parametrization.
An overestimation of Hs is observed for the highest waves in the intermediate and inner RDP.It is concluded that with JONSWAP parametrization is not possible to simulate properly moderate and high waves with the same Γ.So, the others parametrizations available on WWIII are going to be tested, i.e., SHOWEX (Arduhin et al. 2003), Liu & Dalrymple (1978) and Ng (2000), and probably a new configuration will be necessary.It seems that for moderate waves, the mud of the bottom is consolidated being a smoother surface than a sandy bottom, but for higher waves the mud liquefies, triggering a wave damping process.So to represent this behavior, it is necessary to include a source term corresponding to wave damping by fluid mud, that have to be activated from a certain threshold.The dependence on wave condition and bottom composition of the threshold will be topic of future research.
On the other hand, although the use of CFSR winds was validated, a component of the errors of the wave model are inherited from the wind fields.So improving them is the other parallel way to keep improving wave hindcasting for this region.

Figure 1 .
Figure 1.Location and bathymetry of the study area.

Figure 3 .
Figure 3. Spatial distribution of the altimetry database.

Figure 5 .
Figure 5.Comparison between CFSR data and altimetry measures of U10.Spatial distribution of BIAS (left).RMSE (center) and Correlation coefficient (right).

Figure 7 .
Figure 7. Dispersion diagram and q-q plot for the intermediate and inner RDP (left), and the Atlantic (right).Default simulation.