# On the improvement of wave and storm surge hindcasts by downscaled atmospheric forcing: application to historical storms

Emilie Bresson, Philippe Arbogast, Lotfi Aouf, Denis Paradis, Anna

Kortcheva, Andrey Bogatchev, Vasko Galabov, Marieta Dimitrova, Guillaume

Morvan, Patrick Ohl, et al.

► To cite this version:

Emilie Bresson, Philippe Arbogast, Lotfi Aouf, Denis Paradis, Anna Kortcheva, et al.. On the improvement of wave and storm surge hindcasts by downscaled atmospheric forcing: application to historical storms. NATURAL HAZARDS AND EARTH SYSTEM SCIENCES, 2018, 18, pp.997-1012. 10.5194/nhess-18-997-2018. insu-04833693

HAL Id: insu-04833693

https://insu.hal.science/insu-04833693v1

Submitted on 12 Dec 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# On the improvement of wave and storm surge hindcasts by downscaled atmospheric forcing: application to historical storms

Émilie Bresson¹,ᵃ, Philippe Arbogast¹, Lotfi Aouf², Denis Paradis², Anna Kortcheva³, Andrey Bogatchev³, Vasko Galabov³, Marieta Dimitrova³, Guillaume Morvan², Patrick Ohl², Boryana Tsenova³, and Florence Rabier⁴

¹Centre National de Recherches Météorologiques – Groupe de Modélisation et d'Assimilation pour la Prévision, Toulouse, France

²Direction des Opérations pour la Prévision, Département Marine et Océanographie, Météo-France, Toulouse, France

³National Institute of Meteorology and Hydrology, Sofia, Bulgaria

⁴European Centre for Medium-Range Weather Forecasts, Reading, UK

ᵃnow at: Research Institute on Mines and Environment, Université du Québec en Abitibi-Témiscamingue, Rouyn-Noranda, Québec, Canada

**Correspondence:** Émilie Bresson (emilie.bresson@gmail.com)

Received: 5 April 2017 – Discussion started: 15 May 2017

Revised: 14 December 2017 – Accepted: 16 February 2018 – Published: 4 April 2018

**Abstract.** Winds, waves and storm surges can inflict severe damage in coastal areas. In order to improve preparedness for such events, a better understanding of storm-induced coastal flooding episodes is necessary. To this end, this paper highlights the use of atmospheric downscaling techniques in order to improve wave and storm surge hindcasts. The downscaling techniques used here are based on existing European Centre for Medium-Range Weather Forecasts reanalyses (ERA-20C, ERA-40 and ERA-Interim). The results show that the 10 km resolution data forcing provided by a downscaled atmospheric model gives a better wave and surge hindcast compared to using data directly from the reanalysis. Furthermore, the analysis of the most extreme mid-latitude cyclones indicates that a four-dimensional blending approach improves the whole process, as it assimilates more small-scale processes in the initial conditions. Our approach has been successfully applied to ERA-20C (the 20th century reanalysis).

## 1 Introduction

One of the most vulnerable areas affected by winter storms are coastal regions, as their soils are often easily eroded and their population density is high (Barredo, 2007; Clarke and Rendell, 2009; Ferreira et al., 2009; Ciavola et al., 2011; An-

dré et al., 2013). Such storm events are frequently responsible for severe damage, significant economic losses and many casualties. In Europe, sensitive regions include the Atlantic, Mediterranean and Black Sea coasts; in particular, storm surges as high as 2.5 m have been recorded along the Atlantic coasts and 1.5 m along the western Black Sea coasts (Marcos et al., 2009; Ryabinin et al., 1996). These extreme events are often associated with winter low pressure systems; those that affect western Europe are principally mid-latitude cyclones that originate in the Atlantic Ocean (Klawa and Ulbrich, 2003; Della-Marta et al., 2009; Usbeck et al., 2010), and the Bulgarian coasts are hit by cyclones generated in the Mediterranean region (Bocheva et al., 2007). The amplification of wind-generated waves and surges by equinox tides within deep low pressure systems can also produce a significant rise in sea level, resulting in coastal flooding.

For example, during Cyclone Xynthia, which hit the French Atlantic coast on 27 February 2010, a coastal flooding scenario occurred as a result of a tide coefficient of 102 that coincided with a highest astronomical tide between 0.96 and 1.15 m and wind gusts of 160 km h⁻¹ over coastal regions and about 120 km h⁻¹ over land (Rivière et al., 2012). As a result of these conditions, a damaging storm surge crested above 1.60 m at La Rochelle and Les Sables d'Olonne. This example demonstrates that better knowledge of the variability of these extreme coastal events is needed to improve high
surf and storm surge warning systems. In addition, evaluating
the frequency and severity of these events within the frame-
work of ongoing climate change is equally critical. Conse-
quently, a 20th century climatology of wave and storm surge
would provide a useful baseline for coastal protection and
risk management.

The lack of long-term wave records based on in situ mea-
surements and surge archives prevents the development of a
completely observational 20th century climatology for waves
and storm surges. Therefore, reconstructing wave and storm
surge by hindcast using numerical models represents an al-
ternative approach toward establishing a climatology. One
straightforward method for hindcasting involves using global
atmospheric reanalyses as the atmospheric forcing condi-
tions in wave and storm surge numerical models (Reistad
et al., 2011). Several weather forecast centres produce these
global atmospheric reanalyses, including the European Cen-
tre for Medium-Range Weather Forecasts (ECMWF).

The ECMWF Re-Analyses (ERA) include different prod-
ucts that have various date ranges, spatial resolutions and as-
similated datasets (Table 1, Uppala et al., 2005; Dee et al.,
2011; Poli et al., 2016). Although we can use the finer-scale
reanalysis as initial conditions for a given period, dynamical
downscaling of the global reanalyses is also necessary, since
they are too coarse to force the regional wave and storm surge
models. Furthermore, certain mesoscale processes related to
the formation of strong surface winds, such as sting jets
(Hewson and Neu, 2015), are absent even in ERA-Interim,
one of the higher-resolution reanalyses available from the
ECMWF. Therefore, in order to better resolve mesoscale fea-
tures associated with mid-latitude cyclone development and
their interaction with locally complex coastal topography,
dynamical downscaling can be applied to these reanalyses
using a high-resolution numerical model (e.g., Reistad et al.,
2011; Li et al., 2016).

In this study, we apply two different downscaling methods
on ERA datasets. The first one is a simple dynamical down-
scaling approach beyond the reanalysis truncation, whereas
the second is more complex. We evaluate to what extent the
mesoscale features resolved by the first downscaling tech-
nique impact our surge and wave reconstruction over the
French and Bulgarian coasts, followed by an examination of
the added value of the second downscaling method against
the first, simpler one. As observations are spatially and tem-
porally scattered in these regions, we focus on 30 extreme
events between 1924 and 2012 that targeted the French and
Bulgarian coasts. The selected cases offer a large panel of
observed extreme events with various affected areas (in par-
ticular, the French Atlantic and Mediterranean coasts and the
Bulgarian Black Sea coast), including cases with more or less
extended impacted zones, different cyclone trajectories and
amplitudes and varied highest astronomical tides (Table 2).
In the present paper, we first describe the methodology and
data used for the downscaling strategies (Sect. 2.1) and then
the wave and surge models' configurations (Sect. 2.2). In

Sect. 3, we first compare the results from the two downscal-
ing techniques on reconstructing an intense cyclone’s devel-
opment, then we evaluate wave hindcasts and storm surge
model skill, followed by an analysis of our early 20th cen-
tury cases. Finally, Sect. 4 summarizes our conclusions.

## 2 Methodology

### 2.1 Dynamical downscaling of reanalyses

The general method of dynamical downscaling uses a coarse-
resolution dataset, like global atmospheric reanalysis data, as
initial conditions for a numerical atmospheric model. Three
ECMWF reanalyses are selected for this study: ERA-20C,
ERA-40 and ERA-Interim (Table 1). They are all produced
by older versions of the Integrated Forecasting System (IFS),
the ECMWF’s operational forecasting coupled model sys-
tem. ERA-40 includes conventional observations (e.g., sur-
face stations, buoys, radiosondes), polar satellites and geo-
stationary satellites. ERA-Interim datasets benefit from im-
provements in assimilation methods and a large expansion of
available data, with observation quantity and quality increas-
ing over time. In order to mitigate this inhomogeneity in the
20th century reanalysis, only observations of surface pres-
sure and surface marine winds are assimilated in the ERA-
20C dataset. In order to provide the best possible atmospheric
conditions for wave and storm surge hindcast, the following
ERA datasets are downscaled for each event: ERA-20C for
cases before 1957, ERA-40 for the 1957–1978 period, and
ERA-Interim for storms occurring in 1979 and thereafter (Ta-
ble 2). The designator “ERA-x” is used in this manuscript to
describe a group of cases where more than one ERA reanal-
ysis product is applied.

Hereafter, this study focuses on the advantages of down-
scaling global atmospheric reanalysis for the development of
wave and storm surge hindcasts. Over both the French and
Bulgarian domains, numerical weather prediction (NWP)
models require high horizontal and temporal resolution, es-
pecially for the storm surge model hindcast. For French
events, the selected model, ARPEGE (Action de Recherche
Petite Echelle Grande Echelle), is the operational global
primitive-equation NWP system used at Météo-France and
is based on the ARPEGE-IFS software developed in collab-
oration with ECMWF (Table 3; Courtier et al., 1991). A
stretched grid allows for a finer horizontal resolution over
France (around 10 km). The version used here has 70 hy-
brid vertical levels from 17 m to 70 km height. The Bulgarian
events are hindcast from ALADIN (Aire Limitée, Adaptation
dynamique, Développement InterNational) model, which is a
limited-area model based on the ARPEGE system (Radnóti
et al., 1995). The model’s core characteristics are the same as
for ARPEGE.

Two dynamical downscaling methods are examined here,
hereafter referred to as D1 and D2, where D2 represents an
<table><caption>**Table 1.** Characteristics of ERA-20C, ERA-40 and ERA-Interim reanalyses. 4D-Var (3D-Var): four-dimensional (three-dimensional) variational analysis; VarBC: variational bias correction of surface pressure observations.</caption><thead><tr><th></th><th>ERA-20C</th><th>ERA-40</th><th>ERA-Interim</th></tr></thead><tbody><tr><td>Time period</td><td>1900–2010</td><td>1957–2002</td><td>1979–present</td></tr><tr><td>IFS version</td><td>Cy38r1</td><td>Cy23r4</td><td>Cy31r2</td></tr><tr><td>Data assimilation system</td><td>24 h 4D-Var; VarBC</td><td>6 h 3D-Var</td><td>12 h 4D-Var; VarBC</td></tr><tr><td>Spectral resolution</td><td>T159 (~ 125 km)</td><td>T159 (~ 125 km)</td><td>T255 (~ 80 km)</td></tr><tr><td>Number of vertical levels</td><td>91</td><td>60</td><td>60</td></tr><tr><td>Vertical scale (from the surface up to)</td><td>0.01 hPa (~ 80 km)</td><td>0.1 hPa (~ 64 km)</td><td>0.1 hPa (~ 64 km)</td></tr><tr><td>Pressure levels</td><td>37</td><td>23</td><td>37</td></tr><tr><td>Reference</td><td>Poli et al. (2016)</td><td>Uppala et al. (2005)</td><td>Dee et al. (2011)</td></tr></tbody></table>

<table><caption>**Table 2.** List of the 30 cases selected for this study. Coast: Atl–Med for Atlantic and Mediterranean. Tide gauges: number of available and useful tide gauges. Storm surge (metres): maximum storm surge recorded. Asterisk is for unknown information.</caption><thead><tr><th>Coast</th><th>Date</th><th>Tide gauges</th><th>Storm surge</th><th>Downscaling</th><th>ECMWF reanalyses</th></tr></thead><tbody><tr><td rowspan="10">Atlantic</td><td>8 October 1924</td><td>*</td><td>*</td><td>D1</td><td>ERA-20C</td></tr><tr><td>14 March 1937</td><td>*</td><td>*</td><td>D1</td><td>ERA-20C</td></tr><tr><td>31 January–1 February 1953</td><td>*</td><td>&gt; 3</td><td>D1</td><td>ERA-20C</td></tr><tr><td>13 February 1972</td><td>10</td><td>1.83</td><td>D1</td><td>ERA-40</td></tr><tr><td>30 November–2 December 1976</td><td>12</td><td>1.36</td><td>D1</td><td>ERA-40</td></tr><tr><td>11–13 January 1978</td><td>7</td><td>1.65</td><td>D1</td><td>ERA-40</td></tr><tr><td>15–16 October 1987</td><td>12</td><td>1.72</td><td>D1</td><td>ERA-Interim</td></tr><tr><td>26 February–1 March 1990</td><td>6</td><td>1.67</td><td>D1</td><td>ERA-Interim</td></tr><tr><td>2–4 January 1998</td><td>5</td><td>1.60</td><td>D1/D2</td><td>ERA-Interim</td></tr><tr><td>6 November 2000</td><td>8</td><td>1.00</td><td>D1/D2</td><td>ERA-Interim</td></tr><tr><td>17 December 2004</td><td>7</td><td>1.30</td><td>D1/D2</td><td>ERA-Interim</td></tr><tr><td>9 November 2007</td><td>2</td><td>2.20</td><td>D1/D2</td><td>ERA-Interim</td></tr><tr><td>10 March 2008 (Johanna)</td><td>7</td><td>1.30</td><td>D1/D2</td><td>ERA-Interim</td></tr><tr><td>23–24 January 2009 (Klaus)</td><td>10</td><td>1.29</td><td>D1</td><td>ERA-Interim</td></tr><tr><td>28 February 2010 (Xynthia)</td><td>8</td><td>&gt; 1.60</td><td>D1/D2</td><td>ERA-Interim</td></tr><tr><td rowspan="4">Mediterranean</td><td>6 November 1982</td><td>*</td><td>*</td><td>D1</td><td>ERA-Interim</td></tr><tr><td>6–7 February 2009</td><td>7</td><td>0.60</td><td>D1/D2</td><td>ERA-Interim</td></tr><tr><td>24–25 December 2009</td><td>6</td><td>0.50</td><td>D1/D2</td><td>ERA-Interim</td></tr><tr><td>19 February 2010</td><td>6</td><td>0.50</td><td>D1/D2</td><td>ERA-Interim</td></tr><tr><td>Atl–Med</td><td>27 December 1999 (Martin)</td><td>4</td><td>1.60</td><td>D1/D2</td><td>ERA-Interim</td></tr><tr><td rowspan="10">Bulgarian</td><td>5–21 October 1976</td><td>2</td><td>1.00</td><td>D1</td><td>ERA-40</td></tr><tr><td>16–21 January 1977</td><td>1</td><td>0.60</td><td>D1</td><td>ERA-40</td></tr><tr><td>13–23 February 1979</td><td>3</td><td>1.43</td><td>D1</td><td>ERA-Interim</td></tr><tr><td>7–10 January 1981</td><td>0</td><td>*</td><td>D1</td><td>ERA-Interim</td></tr><tr><td>24–31 December 1996</td><td>2</td><td>1.00</td><td>D1</td><td>ERA-Interim</td></tr><tr><td>15–19 December 1997</td><td>1</td><td>1.30</td><td>D1</td><td>ERA-Interim</td></tr><tr><td>20–27 January 1998</td><td>2</td><td>0.90</td><td>D1</td><td>ERA-Interim</td></tr><tr><td>1–3 July 2006</td><td>2</td><td>0.60</td><td>D1</td><td>ERA-Interim</td></tr><tr><td>8–11 March 2010</td><td>2</td><td>0.90–1.00</td><td>D1</td><td>ERA-Interim</td></tr><tr><td>7–9 February 2012</td><td>2</td><td>*</td><td>D1</td><td>ERA-Interim</td></tr></tbody></table>

improved version of D1. For D1, the necessary data from the global fields of ERA-x are interpolated to the plane model domain both on the horizontal and vertical scale for each NWP system, ARPEGE and ALADIN. The upper-air initialization step uses the spectral coefficients of ERA-x data.

Then we apply the Schmidt transformation, which is well defined in spectral space to project the fields into the ARPEGE stretched grid. The land-surface initialization is not straightforward, since there are many differences between the ERA reanalyses and the NWP models in terms of the applied land-
<table><caption>**Table 3.** Outline of the numerical models required for wave and storm surge hindcasts.</caption><thead><tr><th>Purpose</th><th>Model</th><th>Resolution</th><th>Coupling-initial conditions data</th><th>Domain</th></tr></thead><tbody><tr><td rowspan="3">Atmosphere</td><td>ARPEGE D1</td><td>T798 (~ 10 km)</td><td>ERA-x</td><td>global</td></tr><tr><td>ARPEGE D2</td><td>T798 (~ 10 km)</td><td>ERA-x + ARPEGE</td><td>global</td></tr><tr><td>ALADIN</td><td>10 km</td><td>ARPEGE D1</td><td>Bulgaria</td></tr><tr><td rowspan="2">Wave</td><td>MFWAM</td><td>0.1°</td><td>ARPEGE D1/D2</td><td>western Europe</td></tr><tr><td>SWAN</td><td>0.1°</td><td>ALADIN</td><td>Bulgaria</td></tr><tr><td rowspan="3">Surge</td><td>HYCOM</td><td>1 km</td><td>ARPEGE D1/D2 + bathymetry</td><td>ATL</td></tr><tr><td>HYCOM</td><td>1 km</td><td>ARPEGE D1/D2 + bathymetry</td><td>MED</td></tr><tr><td>MF model</td><td>0.0333°</td><td>ALADIN + bathymetry</td><td>Black Sea</td></tr></tbody></table>

surface parameterizations and physiographic databases. For instance, the Tiled ECMWF Scheme for Surface Exchanges over Land (TESSEL) scheme of ERA-x uses four soil layers with fixed thicknesses, each layer having its own water content. The land-surface scheme of ARPEGE, however, only uses two layers in our experiments; the top layer has a fixed size of 1 cm, and the second layer overlaps the first one and has a variable depth. Furthermore, for a given grid point, soil types are often very different in the two land-surface schemes. Therefore, using the raw land-surface datasets from ERA-x as initial conditions would be troublesome, since the water saturation fraction depends on the soil type. Thus, we interpolate the surface fields so as to preserve as much as possible the ERA-x surface heat and momentum fluxes (Boisserie et al., 2016). The procedure is based on the conservation of the soil wetness index (a relevant indicator for soil water availability) during the interpolation process, since soil water availability is supposed to regulate the partition of latent and sensible heat fluxes, which, in turn, influence energy and water exchanges between the atmosphere and the land surface. The resulting files are initial conditions (IC-1) for the NWP forecasts (Fig. 1, top). Then, hourly forecasts are produced twice a day, starting at 00:00 and at 12:00 UTC, and run for 18 h. Only hourly forecasts from +6 h to +18 h are used. The first 6 h are not taken into account to prevent model spin up, and after h+18, the next forecast time is considered (Fig. 1, top). Forecasts are produced from a week (d-7) before to 2 days (d+2) after the day (d) that the storm impacted the coastline. The D2 method is more complex than D1 (Fig. 1, bottom). The D2 method also uses hourly forecasts produced twice a day, at 00:00 and at 12:00 UTC, starting from h+06 to h+18, and the forecast starts 9 days (d-9) before and continue until 2 days after (d+2) the day (d) that the storm impacted the coastline. Instead of using independent initial conditions (IC-1) like in D1 for the 00:00 and 12:00 UTC forecasts, the initial conditions for D2 (IC-2) include information from the last 6 h forecast (Fig. 1, bottom). Consequently, the D2 method allows us to evaluate the importance of taking into account small wavelengths beyond the reanalysis truncation that are not considered in D1. Fur-

thermore, after a short period of time (3 h), non-linearities trigger small-scale processes which are consistent with the large scale. This small-scale information provided by the 6 h forecast is blended with the large-scale information given by the interpolated reanalysis (IC-1; Fig. 1, bottom). This procedure was cycled 4 times in 2 days before the first 00:00 UTC forecast used as forcing for the wave and storm surge models. Therefore, the determination of one single initial condition (IC-2) uses four reanalyses. The D2 technique is applied to 10 recent French coastal flooding events (Table 2). These 10 cases represent a diverse panel of events affecting different coastlines with adequate observational data (satellite altimeters and tide gauges) to evaluate the reconstruction of the wave and storm surge observations and to enable a comparison between D1 and D2.

## 2.2 Description of wave and storm surge models

In order to ensure consistency in our case studies, the selected wave and storm surge models share similar general characteristics, despite being adapted specifically to either the French or Bulgarian coasts.

### 2.2.1 Wave models

The French coast extreme wave events are hindcast with the Meteo-France WAve Model (MFWAM), a third-generation model of the operational wave forecasting system of Météo-France (Table 3). This model is based on the IFS-CY36R4 of the European wave model (ECWAM) with modified source terms for the dissipation by wave breaking and the air friction dedicated to swell damping as described in Ardhuin et al. (2010). The MFWAM model uses the wind input term as defined in Bidlot et al. (2005). The dissipation by wave breaking is directly related to the wave spectrum with a saturation rate of dissipation. The source term is a combination of an isotropic component and a direction-dependent component that controls the directional spread of the resulting wave spectra. It also includes a cumulative effect describing the smoothing of big breakers on small breakers. The term additionally uses a wave turbulence interaction com-
<figure><img src="image_3.png" /><figcaption>**Figure 1.** Schematic representation of D1 and D2 techniques. Energy spectra are within the small images. The red parts of forecast are the forecast data used as input forcing in the wave and storm surge models.</figcaption></figure>

ponent, which, as indicated in Ardhuin et al. (2010), is of
secondary importance. The MFWAM model uses a quadru-
plet non-linear interaction term based on the discrete inter-
actions approximation as defined in the ECWAM model. In
this study, a nested MFWAM model is implemented with a
grid size of 0.1° for western Europe, including the Mediter-
ranean Sea. The domain boundaries are 20–72° N, 32° W–
42° E (EURAT01 domain in Fig. 2). The wave spectrum is
discretized in 24 directions and 30 frequencies starting from

0.035 to 0.58 Hz. This regional model is forced by bound-
ary conditions provided by the global MFWAM model with
a grid size of 0.5°. The global MFWAM model is driven by
6 h ERA-x winds. The SWAN (Simulating Waves Nearshore)
model is used for the Bulgarian cases (Table 3). It is a third-
generation wave model that is especially designed to simu-
late waves in nearshore waters and is often applied to en-
closed and semi-enclosed seas, estuaries and lakes (Booij
et al., 1999). The model computes random short-crested,
<figure><img src="image_4.png" /><figcaption>**Figure 2.** Locations of EURAT01 (black), ATL (blue), MED (green) and BUL (red) domains used in the study, respectively, for European 0.1° resolution grid and Atlantic, Mediterranean and Bulgarian domains.</figcaption></figure>

wind-generated waves in coastal regions and inland waters. SWAN accounts for wave propagation and transitions from deep to shallow water at finite depths by solving the spectral wave action balance equation, which includes source terms for the wind input, non-linear interactions, whitecapping, bottom friction and depth-induced breaking. The model performance, the parameterizations of the wave generation and dissipation processes and other aspects of SWAN applied to the Black Sea basin have been addressed in previous studies (Akpinar et al., 2012; Arkhipkin et al., 2014; Rusu et al., 2014). The model domain that is used for the simulations of our historical Black Sea storms is based on a numerical grid covering the entire Black Sea area (40–47° N and 27–42° E; hereafter named BUL; Fig. 2) with a mesh size of 0.0333° in latitude and longitude. The spectral discretization is based on 36 directions and 30 frequencies logarithmically spaced from 0.05 to 1.00 Hz. The wind input parameterization follows Komen et al. (1984), and whitecapping is based on Hasselmann (1974), with the $\delta$ coefficient (which determines the dependency of whitecapping on wave number) set to 1 (following Rogers et al., 2003). This specific set of parameterizations is chosen to have the lowest bias, root mean square error (RMSE) and scatter index when compared to results

from the model and the along-track satellite altimetry data. The bathymetry data for the wave model are obtained by the digitalization of proprietary maps provided by the Bulgarian military's hydrographic service.

### 2.2.2 Storm surge models

The operational surge model of Météo-France (Daniel et al., 2001) is a barotropic two-dimensional version of the HYBRID Coordinate Ocean Model (HYCOM) implemented by SHOM (Service Hydrographique et Océanographique de la Marine) from the three-dimensional version (Table 3; Bleck, 2002; Baraille and Filatoff, 1995). The HYCOM code is managed by an international consortium, including COAPS (Center for Ocean-Atmospheric Prediction Studies, USA), NRL (Naval Research Laboratory, USA), SHOM (France), DMI (Danish Meteorological Institute, Denmark) and NERSC (Nansen Environmental and Remote Sensing Center, Norway). The model is run on two domains (as shown in Fig. 2): ATL corresponds to the northeast Atlantic area (Bay of Biscay, English Channel and North Sea) from 43 to 62° N and from 9° W to 10° E, and MED defines the Mediterranean Sea domain from 30 to 46° N and from 9° W to 37° E. In both domains, the model runs on a grid size of approximately 1 km on the French coast (curvilinear grid). The tides imposed at the marine boundaries are computed according to the 17 harmonic components from the COMAPI (COastal Modelling for Altimetry Product Improvement) project regional atlas implemented in the northeast Atlantic Ocean area (Cancet et al., 2010). The bottom friction coefficient is spatially variable and has been optimized to properly reproduce the propagation of tides. Tides are discarded in the storm surge computation, for which another computation of the tides, based on harmonic components obtained from measurements by SHOM, is added to the storm surge in order to more accurately represent the sea level at specific locations. The bottom friction coefficient is constant and taken as equal to 0.002. For both HYCOM configurations (ATL and MED), the drag coefficient used to compute the wind stress follows the Charnock (1955) scheme with a constant Charnock parameter of 0.025.

The simulations of storm surges for Black Sea cases are based on the storm surge model of Météo-France (Daniel et al., 2001), which was adapted for the Black Sea in Mungov and Daniel (2000) (Table 3). The model is depth integrated, and tides are not taken into account, as their amplitude is less than 9 cm in the Black Sea. The model grid for the Black Sea is a regular spherical grid with a spatial resolution of 0.0333° that covers the entire Black Sea. The bottom friction coefficient is $1.5 \times 10^{-3}$ over the shelf. In addition, the depth of the Black Sea mixed layer is considered as a liquid bottom given the very stable stratification of the Black Sea waters and the shallowness of the mixed layer depth, and as such, the bottom friction coefficient is defined as $1.5 \times 10^{-5}$ over the liquid bottom. Data about the seasonal variations of
<figure><img src="image_5.png" /><figcaption>**Figure 3.** Mean-sea level pressure (hPa) from observations (a) and ERA-Interim reanalysis at 06:00 UTC, 26 December 1999 (b), from 12 h forecast using the D1 (c) and D2 (d) downscaling methods at 18:00 UTC, 25 December 1999.</figcaption></figure>

the Black Sea mixed layer depth are taken from the study by Kara et al. (2009). Without this liquid bottom setup, the depth-integrated models for the Black Sea fail to simulate any surge, even if strong, constant winds are used as input. The bathymetry data for the storm surge model were ob-

tained by digitizing proprietary maps provided by the Bulgarian military hydrographic service.
<figure><img src="image_6.png" /><figcaption>**Figure 4.** Scatter plots of significant wave heights (SWHs) of model MFWAM and altimeters (ENVISAT and Jason-1) for the 2004, 2007, 2008 and 2010 French storms. (a) and (b) stand for runs with interpolated ERA-Interim and D1 wind forcing, respectively.</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>**Figure 5.** Comparison of the simulated significant wave heights (SWHs) with downscaled wind input and ERA-Interim wind input with the data from the ENVISAT track crossing the western Black Sea at 20:00 UTC on 7 February 2012. Purple and green stand for ERA-Interim and D1 forcing, respectively. Red line stands for ENVISAT observations.</figcaption></figure>

<caption>**Figure 6.** Comparison of the simulated significant wave heights using the two wind inputs (downscaled wind input D1 and ERA-Interim) with the data by ADCP located on the western Black Sea coast at 20 m depth during the storm of 7–8 February 2012. ADCP location coordinates: 43°04′49″ N, 28°01′40″ E. Purple and green stand for ERA-Interim and D1 forcing, respectively. The red line represents the ADCP observations.</caption>

# 3 Results

## 3.1 Impact of the two downscaling techniques on a deep cyclone development

The effects of the two downscaling techniques on the reconstruction of intense storms are presented for the case of the Lothar storm, an extreme cyclogenesis event (occurring a few hours before the Martin storm described further in Sects. 3.2 and 3.3) in December 1999. It is the most severe storm in terms of pressure gradient, surface winds and displacement velocity to hit France within the observational record (Wernli et al., 2002; Rivière et al., 2010). This storm did not produce extreme wave and storm surge, and thus it was not selected for hindcasts. Nevertheless, it is interesting to look at the behaviour of both downscaling strategies for this particular case due to its uniquely tight horizontal pressure gradient.

For this storm, the D1 method slightly improves the ERA-Interim reanalysis fields, but the D2 downscaling better reproduces the cyclone structure over northern France (Fig. 3). Statistical analysis using the mean, the bias, the root mean square error (RMSE) and the standard deviation (SD) error is performed with the 12 meteorological stations available in an area encompassing the low pressure system (48–50° N, 2–4° E). This analysis confirms that the use of D1 forcing is an improvement compared to using an ERA-Interim reanalysis with respect to surface observations. The use of D2 slightly improves the reconstruction of the observations (Table 4).
<table><caption>**Table 4.** Statistics for mean sea-level pressure from ERA-Interim reanalysis at 06:00 UTC, 26 December 1999, 12 h forecast using the D1 and D2 at 18:00 UTC, 25 December 1999, versus observations at 06:00 UTC, 26 December 1999. Mean (hPa), standard deviation (SD) error (hPa), bias (hPa), root mean square error (RMSE; hPa). Calculations are done for the nearest point. Small domain corresponds to 48–50° N, 2–4° E and includes 12 pairs of data and model values.</caption><thead><tr><th></th><th>Mean</th><th>SD</th><th>Bias</th><th>RMSE</th></tr></thead><tbody><tr><td>Obs</td><td>973</td><td>2</td><td>–</td><td>–</td></tr><tr><td>ERA-Interim</td><td>993</td><td>10</td><td>12</td><td>18</td></tr><tr><td>D1</td><td>980</td><td>1</td><td>6</td><td>6</td></tr><tr><td>D2</td><td>977</td><td>1</td><td>5</td><td>5</td></tr></tbody></table>

## 3.2 Wave hindcasts

For the wave reconstruction evaluation, simulated significant wave heights (SWHs) are compared against observations from satellite altimeter data and in situ observations. Several satellites operated over the French and Bulgarian coasts during the storms: TOPEX-Poseidon (1992–2005), ERS2 (1995–2011), ENVISAT (2002–2012) and Jason-1 (2002–2013). In addition, buoys and acoustic Doppler current profiler (ADCP) provide in situ SWH information. The limited scope of each of these observational datasets, together with the coarse resolution of altimeter measurements, preclude a comprehensive validation for all the selected cases. For an initial evaluation of our modelling approach, the results from the wave model driven by ERA-x and D1 data are compared to available altimeter data. The simulated wave heights are collocated with the altimeter tracks within a time window of 3 h. For the 2004, 2007, 2008 and 2010 French Atlantic coast storms and the 2012 Bulgarian storm, data are collected from two satellite altimeters, Jason-1 and ENVISAT. The scatter plots between model and altimeter wave heights indicate that the use of D1 winds provides a better fit to the data (Fig. 4). In particular, when compared to the results for the wave model driven by ERA-Interim initial conditions, the use of D1 data reduces the normalized root mean square error (NRMSE) from 17.1 to 13.1 %, largely owing to a significant reduction of bias from –35 to –4 cm (Fig. 4). The D1 downscaling also leads to a better fit for high SWHs, providing an important validation for extreme wave events. For the 1998, 1999 and 2000 storms, altimeter wave heights from TOPEX and ERS2 are also used for the evaluation of the modelled SWHs, and the same tendency is found, with an improvement of the reconstruction of SWHs using D1 winds over ERA-Interim winds (not shown).

As satellite altimeters provide data along a track, these observations can be useful for mapping the spatial distribution of the SWH. For further examination, we present the 2012 Bulgarian storm as an example of a more detailed evaluation of the reconstruction against observations. The wave model

<figure><img src="image_8.png" /><figcaption>**Figure 7.** Time series of significant wave heights (SWHs) for the storm on February 2010 near Nice (43°24′0″ N, 7°48′0″ E) in the Mediterranean Sea. Purple and green stand for ERA-Interim and D1 forcing, respectively. The red line shows the time series of the Nice buoy observations.</figcaption></figure>

outputs using ERA-Interim or D1 initial conditions are first compared to the 214 along-track data points measured by the Jason-1 and ENVISAT satellite altimeters on 7 and 8 February 2012. The wave reconstruction given by D1 forcing more closely matches the satellite observations, especially in terms of wave intensity over the southern part of the satellite track (Fig. 5). However, the maximum observed SWH value is not reached by the model for both the ERA-Interim winds and the D1 winds. Regarding the temporal evolution of the 2012 Bulgarian storm, we use in situ ADCP to check if the peak SWH occur at the same time in the observations and the reconstruction. In Fig. 6, we compare the SWH data from the ADCP located at Pasha Dere beach at 20 m depth provided by the Bulgarian Institute of Oceanology (Valchev et al., 2014) to our wave model outputs. The use of D1 generally overestimates the measured SWHs, while the use of ERA-Interim underestimates the wave heights. However, the use of D1 winds leads to better matching of the temporal structure of the wave. The overall improvement of the SWH reconstruction by using D1 is confirmed by the statistical analysis in Table 5. The temporal evolution of a storm can also be evaluated with in situ buoys. For example, for the 2010 Mediterranean storm, we compare the time series of SWHs from model and buoy data (43.4° N and 7.8° E) off the coast of Nice, France, at the peak of the storm (Fig. 7). The results show that the SWH induced by using D1 data more closely matches the buoy observations when compared to the ERA-Interim data forcing. Given our validation of the D2 approach discussed in Sect. 3.1, the D2-driven SWH hindcasts of the 2004, 2007, 2008 and 2010 French Atlantic storms are also compared to satellite altimeter data. The statistical analysis (bias and NRMSE) reveals that the use of D2 winds leads
<table><caption>**Table 5.** Comparison of SWAN wave model significant wave heights (SWHs; metres) and altimeter data from ENVISAT and Jason-1 satellites for the 2012 case over the Bulgarian coast.</caption><thead><tr><th rowspan="2">Time of satellite track</th><th rowspan="2">Pairs</th><th colspan="3">Mean</th><th colspan="3">Bias</th><th colspan="3">RMSE</th><th colspan="3">Scatter index</th></tr><tr><th>Obs</th><th>ERA-Interim</th><th>D1</th><th>ERA-Interim</th><th>D1</th><th>ERA-Interim</th><th>D1</th><th>ERA-Interim</th><th>D1</th><th>ERA-Interim</th><th>D1</th></tr></thead><tbody><tr><td rowspan="3">7 February 2012</td><td>08:00 UTC</td><td>44</td><td>3.9</td><td>3.5</td><td>4.1</td><td>-0.43</td><td>0.21</td><td>0.60</td><td>0.37</td><td>0.15</td><td>0.10</td><td></td></tr><tr><td>14:00 UTC</td><td>76</td><td>3.6</td><td>3.2</td><td>3.8</td><td>-0.41</td><td>0.15</td><td>0.66</td><td>0.57</td><td>0.18</td><td>0.16</td><td></td></tr><tr><td>20:00 UTC</td><td>51</td><td>6.4</td><td>5.3</td><td>6.3</td><td>-1.08</td><td>-0.09</td><td>1.14</td><td>0.37</td><td>0.18</td><td>0.06</td><td></td></tr><tr><td>8 February 2012</td><td>14:00 UTC</td><td>43</td><td>5.6</td><td>4.4</td><td>4.7</td><td>-1.22</td><td>-0.94</td><td>1.37</td><td>1.16</td><td>0.24</td><td>0.21</td><td></td></tr></tbody></table>

<figure><img src="image_9.png" /><figcaption>**Figure 8.** Variation of the bias (a) and the normalized root mean square error (NRMSE; b) of significant wave heights (SWHs) from the model MFWAM in comparison with the altimeters (ENVISAT and Jason-1) for the 2004, 2007, 2008 and 2010 French storms. Purple, green and blue stand for ERA-Interim, D1 and D2 forcing, respectively</figcaption></figure>

<table><caption>**Table 6.** Number of observations used for calculations of WNOE for each region and each forcing.</caption><thead><tr><th></th><th>ERA-x</th><th>D1</th><th>D2</th></tr></thead><tbody><tr><td>ATL</td><td>34</td><td>34</td><td>15</td></tr><tr><td>MED</td><td>13</td><td>13</td><td>13</td></tr><tr><td>BUL</td><td>9</td><td>9</td><td>0</td></tr></tbody></table>

<table><caption>**Table 7.** Portion of cases (in percent) with $||WNOE|| < 20\%$ for each coast (ATL: Atlantic; MED: Mediterranean Sea; BUL: Bulgarian; common cases: cases using D1 and D2 forcing).</caption><thead><tr><th></th><th>ERA-x</th><th>D1</th><th>D2</th></tr></thead><tbody><tr><td>ATL</td><td>21</td><td>63</td><td>80</td></tr><tr><td>MED</td><td>0</td><td>54</td><td>38</td></tr><tr><td>BUL</td><td>33</td><td>100</td><td>-</td></tr><tr><td>Common cases</td><td>18</td><td>64</td><td>61</td></tr></tbody></table>

to better results than the use of D1 winds (Fig. 8). Biases of SWHs are slightly improved using D2 winds over D1 winds; however, D2 winds slightly increase the NRMSE of SWHs for the 2004, 2007 and 2008 storms. The D2 method only slightly improves the NRMSE of SWH for Cyclone Xynthia (February 2010). While the application of the D2 method winds does not lead to an improved result over D1 in all

cases, D2 appears to show better skill for events with higher wind speeds, such as the ones observed during the Lothar storm.

## 3.3 Storm surge hindcasts

Storm surge hindcasts can be evaluated by tide gauge measurements. A network of 25 tide gauges along the French coasts is maintained to validate the surge model implemented at Météo-France. Furthermore, an additional 12 hydro-meteorological stations are located along the Bulgarian coasts for validation purposes. Depending on the storm extent and instrument condition, the number of available data points is different for each storm (Table 2). For a global hindcast evaluation of tide gauges, all available measurements with a peak in storm surge are selected. Weighted normalized observation error (WNOE) is calculated to highlight the overestimation and underestimation of the simulated maximum storm surges with respect to available measurements, and it is defined in Eq. (1).

$$WNOE = 100 \cdot \alpha(t_{\text{sim}}) \cdot \left( \frac{X_{\text{sim}} - X_{\text{mea}}}{X_{\text{mea}}} \right) \quad (1)$$

In this simple calculation, $t_{\text{sim(mea)}}$ is the time related to the simulation outputs (measurements, in hours), $X_{\text{sim(mea)}}$ is the simulated (measured) value of maximum storm surge (in centimetres) and $\alpha$ is the weighting coefficient. The value of $\alpha$ is equal to 0.9 if the simulated maximum of storm surge
<figure><img src="image_10.png" /></figure>

**Figure 9.** Storm surges (centimetres) at St Malo (a) and Dunkirk (b) from 14 December 2004 at 15:00 UTC to 19 December 2004 at 06:00 UTC. The measured surge (red line), the reconstructed surge by using the ERA-Interim forcing (purple line), the D1 forcing (green line) and the D2 forcing (blue line) are superimposed. The oscillatory dotted line in the lower part of the graph is used to indicate the time of high and low tides.

<figure><img src="image_11.png" /></figure>

**Figure 10.** Storm surges (centimetres) at Dunkirk from 7 November 2007 at 15:00 UTC to 11 November 2007 at 06:00 UTC. The measured surge (red), the reconstructed surge by using the ERA-Interim forcing (purple), the D1 forcing (green) and the D2 forcing (blue) are superimposed. The oscillatory dotted line in the lower part of the graph is used to indicate the time of high and low tides.

falls within a time window of ±3 h with respect to the observed peak time; if it is sooner or later, the weighting coefficient is set equal to 1.1 to reflect greater bias. For some cases, when no time information is available, no weighting is applied, and thus $\alpha = 1$. When $||WNOE|| < 20\%$, we consider errors to be low or moderate. Moreover, the values are evaluated regarding the number of samples (Table 6). First,

we evaluate the impact of using wind and mean sea-level pressure data from D1 instead of from ERA-x. The storm surge outputs using ERA-x forcing have a tendency to underestimate maximum storm surge compared to D1 forcing (Figs. 9, 10 and 11). Cases with low or moderate errors represent a larger proportion of storm surge events when D1 data are used. In particular, 63 % of storm surge events were associated with low and moderate error in the ATL basin, 54 % for BUL and 100 % for the MED domain. This represents a general improvement over the ERA-x data, which had low/moderate errors for 21 % of storm surge events for ATL, 0 % for BUL and 100 % for the MED domain (Table 7).

Second, the D2 method is applied to two examples of storm surge reconstruction (the Atlantic 2004 and 2007 storms in France) with corresponding statistical analysis. For the December 2004 storm, a deep low of 980 hPa crossed the northern French coasts from west to east, generating high waves and surge along the British Channel and the North Sea coasts due to strong northwesterly winds wrapping behind the system. The maximum observed surge exceeded 1 m at St Malo and Dunkirk during a period of below-average tide (Fig. 9). Over the course of this event, the application of ERA-Interim winds result in an underestimation of the surge by roughly 60 cm at St Malo and 20 cm at Dunkirk (Fig. 9). However, the use of D1 forcing successfully captures the peak of the surge in St Malo and Dunkirk. The use of D2 winds induces an overestimation of the surge of 20 cm at St Malo and roughly the same surge as D1 at Dunkirk. The second example of storm surge hindcast is provided by the November 2007 storm. This event affected the whole North Sea (including Dunkirk and Calais on the French coast) and parts of the eastern British channel. It was associated with a
<figure><img src="image_12.png" /><figcaption>**Figure 11.** The percentage of cases depending of their WNOE range when using ERA-x (purple), D1 (green) or D2 (blue) forcing. All the available observations with a maximum storm surge measurement are taken into account.</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>**Figure 12.** Surface pressure chart (hPa) at 06:00 UTC on 1 February 1953 from http://www.metoffice.gov.uk.</figcaption></figure>

strong northwesterly wind on the North Sea and lasted nearly 24 h. At the peak of the storm event, a surge of 2.30 m was recorded at Dunkirk (Fig. 10). While the ERA-Interim forcing significantly underestimates the surge by 80 cm (Fig. 10), a good fit is obtained by the model with both the D1 and D2 data forcing. For this particular storm, the D2 winds give slightly better surge results on 11 November 2007 at 00:00 UTC. These two storms are examples of the various responses of the storm surge hindcast with both types of downscaling: no significant trend could be highlighted. Overall,

the dispersion of WNOE values for the D2 results is larger than for D1 (Fig. 11), and Atlantic cases are better hindcasted with D2 forcing data (Table 7). The ability of D2 to simulate very deep cyclones could explain this point, since the mesoscale processes involved in strong winds are better described with the D2 approach.

### 3.4 Evaluation of early 20th century cases hindcast using ERA-20C

The 20th century extreme events that occurred before 1957 can be hindcast by using ERA-20C, the 20th century reanalysis ECMWF project (Poli et al., 2016). For these cases, even if there were no available wave observations, a storm surge evaluation is possible due to the availability of reliable sea-level observations.

To validate the concept of downscaling using ERA-20C reanalyses, we concentrated on the major storm that occurred in the North Sea in February 1953 (Fig. 12). It caused severe damage along the Dutch, Belgian and English coasts. Wind intensity around force 10 on the Beaufort scale (around $90 \text{ km h}^{-1}$) were measured in Scotland and northern England. The winds and the low atmospheric pressure combined with exceptional equinox tides were responsible for the surge, which was additionally exacerbated by the funnel shape and shallowness of the North Sea. The Netherlands were the worst affected, resulting in 1836 deaths and widespread property damage (Gerritsen, 2005). Most of the casualties occurred in the southern province of Zeeland; an additional 307 people were reported killed in England, 19 in Scotland and 28 in Belgium as a result of the storm. The most striking feature along the Dutch coast was a long swell with a peak period of 20 s, which induced wave flooding. In our reconstruction of the event, the MFWAM results using the D1 winds indicate SWH exceeding 16 m in the western part of the North Sea at 00:00 UTC on 1 February 1953 (Fig. 13). The storm surge hindcast produces a high surge which is unusual for this area; in particular, along the Dutch and Belgian coastlines storm surges exceeded 3 m either with ERA-20C or D1 data forcing (Fig. 14). The improvement of storm surge reconstruction induced by D1 forcing was particularly marked at IJmuiden, Ostend, Brouwershaven and
<figure><img src="image_14.png" /></figure>

**Figure 13.** Significant wave heights (metres; a) and peak wave period (seconds; b) from the wave model MFWAM with D1 winds outputs on the peak of the storm at 00:00 UTC on 1 February 1953. Mean wave direction is shown with black arrows in (a) when significant wave height are greater than 1.5 m.

<figure><img src="image_15.png" /></figure>

**Figure 14.** The highest simulated storm surges (centimetres) obtained for the period from 30 January to 2 February 1953, with the ERA-20C forcing (a) and with the D1 forcing (b) along the southern North Sea coast.

Dieppe, where the recorded peaks of storm surge are better represented than for ERA-20C (Fig. 15).

## 4 Conclusions

ECMWF reanalyses data are widely used for many climatological studies. However, due to the coarse spatial resolution and the limited temporal resolution of reanalysis model output, there is significant bias for high wind speeds associated

with extreme mid-latitude cyclones. To overcome this problem, dynamical downscaling techniques are implemented and applied to reproduce high-resolution historical atmospheric fields. ERA-20C, ERA-40 and ERA-Interim data are used to encompass the studied period of 1924–2012. Very short range forecasts using 10 km resolution and hydrostatic models initialized with ERA-x analyses provide the down-scaled data, which are used in turn to force wave and storm surge numerical models. This approach was already tested for the North Sea coast for a long period using only ERA-
<figure><img src="image_16.png" /><figcaption>**Figure 15.** The storm surges (centimetres) at (a) IJmuiden, the Netherlands, (b) Ostend, Belgium, (c) Brouwershaven, the Netherlands, and (d) Dieppe, France, from 18:00 UTC on 30 January to 18:00 UTC on 2 February 1953. Two surges are represented: those resulting from ERA-20C forcing (purple) and from the D1 outputs (green). The maximum observed storm surge is added (horizontal plain black line). The tide level is indicated by the dashed black line (at a reduced scale).</figcaption></figure>

40 data. In order to evaluate such downscaling technique on different initial conditions, 30 cases are selected over French and Bulgarian coastlines to offer a diverse selection of storm characteristics in terms of location, intensity, highest astrotomic tide and meteorological context. Some early 20th century cases generating extreme storm surge and waves are part of this selection due to the recent availability of ERA-20C. This study shows a significant and quasi-systematic improvement of wave and storm surge hindcast when using downscaled winds. The evaluation with independent wave observations (such as wave heights from altimeters) shows the strong reduction of bias and improved RMSE of significant wave height for extreme waves events. The downscaling techniques are also well suited for storm surge extreme

events, such as the 1953 storm, since the storm surge reconstruction using the presented approach fits with the recorded data from the Belgian and Dutch coasts. The D2 method, generally leads to an improvement in comparison with D1, especially for cases with small-scale, intense mid-latitude cyclones. Dynamical downscaling is a promising technique for providing an accurate reconstruction of waves and storm surges for the 20th century. After evaluation and calibration with observations, these model outputs can be useful for analyzing the interannual variability of coastal wind storms and for improving the thresholds used in the wave submersion warning system. Regional climate modelling in future studies is expected to address the response of extreme wave and surge variability to storm track modifications due to global
climate change. A further step towards this objective would be to use interactive models of wave and storm surge to enhance the hindcast. We expect that these approaches for reconstructing extreme events will prove valuable for coastal protection and risk management.

*Data availability.* Members of the ECMWF can access the MARS archive for the SYNOP weather station data used in Sect. 3.1. ERA-20C, ERA-40, ERA-Interim reanalysis data, ERS-2 and ENVISAT data, and TOPEX-POSEIDON and Jason-1 data can be obtained from the public server of, respectively, the ECMWF (http://apps.ecmwf.int/datasets/), the ESA (https://earth.esa.int/web/guest/data-access/browse-data-products) and the NASA (https://podaac.jpl.nasa.gov/datasetlist). The other data are available on request from the authors.

*Competing interests.* The authors declare that they have no conflict of interest.

*Acknowledgements.* The research was carried out as part of the IncREO (Increasing Resilience through Earth Observation) project with funding from the European Union Seventh Framework Programme under grant agreement no. 312461. The authors would like to thank the European Commission for its financial support through the 7th Framework Programme. We are also most grateful to Françoise Taillefer for her unconditional technical support and François Bouyssel for his constructive and valuable advice about the second downscaling method. Special thanks go to Philippe Dandin for his involvement in setting up this project. We also thank SHOM for providing the French storm surge measurements. Christophe-Thomas Simmons is warmly thanked for helping to improve the manuscript.

Edited by: Ricardo Trigo

Reviewed by: two anonymous referees

## References

Akpinar, A., van Vledder, G. P., Kömürcü, M. I., and Özger, M.: Evaluation of the numerical wave model (SWAN) for wave simulation in the Black Sea, Cont. Shelf Res., 50, 80–99, https://doi.org/10.1016/j.csr.2012.09.012, 2012.

André, C., Monfort, D., Bouzit, M., and Vinchon, C.: Contribution of insurance data to cost assessment of coastal flood damage to residential buildings: insights gained from Johanna (2008) and Xynthia (2010) storm events, Nat. Hazards Earth Syst. Sci., 13, 2003–2012, https://doi.org/10.5194/nhess-13-2003-2013, 2013.

Arduuin, F., Rogers, E., Babanin, A., Filipot, J.-F., Magne, R., Roland, A., Westhuysen, A. V. D., Queffeulou, P., Lefevre, J.-M., Aouf, L., and Collard, L.: Semi-empirical dissipation source functions for ocean waves: Part 1, definition, calibration and validation, J. Phys. Oceanogr., 40, 1917–1941, 2010.

Arkhipkin, V. S., Gippius, F. N., Koltermann, K. P., and Surkova, G. V.: Wind waves in the Black Sea: results of a hind-

cast study, Nat. Hazards Earth Syst. Sci., 14, 2883–2897, https://doi.org/10.5194/nhess-14-2883-2014, 2014.

Baraille, R. and Filatoff, N.: Modèle shallow-water multicouches isopycnal de Miami, Tech. Rep. 003/95, SHOM/CMO, 1995.

Barredo, J. I.: Major flood disasters in Europe: 1950–2005, Nat. Hazards, 42, 125–148, 2007.

Bidlot, J., Janssen, P., and Abdalla, S.: On the importance of spectral wave observations in the continued development of global wave models, in: Proc. 5th. Int. Symposium on Ocean Wave Measurements and Analysis (WAVES 2005), 2005.

Bleck, R.: An oceanic general circulation model framed in hybrid isopycnic – Cartesian coordinates, Ocean Modell., 4, 55–88, https://doi.org/10.1016/S1463-5003(01)00012-9, 2002.

Bocheva, L., Georgiev, C. G., and Simeonov, P.: A climatic study of severe storms over Bulgaria produced by Mediterranean cyclones in 1990–2001 period, Atmos. Res., 83, 284–293, 2007.

Boisserie, M., Decharme, B., Descamps, L., and Arbogast, P.: Land surface initialization strategy for a global reforecast dataset, Q. J. Roy. Meteor. Soc., 142, 880–888, https://doi.org/10.1002/qj.2688, 2016.

Booij, N., Ris, R. C., and Holthuijsen, L. H.: A third-generation wave model for coastal regions: 1. Model description and validation, J. Geophys. Res.-Oceans, 104, 7649–7666, 1999.

Cancet, M., Lux, M., Pénard, C., Lyard, F., Birol, F., Lemouroux, J., Bourgogne, S., and Bronner, E.: COMAPI: New regional tide atlases and high frequency dynamical atmospheric correction, in: Ocean Surface Topography Science Team meeting, 2010.

Charnock, H.: Wind Stress on a Water Surface, Q. J. Roy. Meteor. Soc., 81, 639–640, 1955.

Ciavola, P., Ferreira, O., Haerens, P., Van Koningsveld, M., Armaroli, C., and Lequeux, Q.: Storm impacts along European coastlines. Part 1: The joint effort of the MICORE and ConHaz Projects, Environ. Sci. Pol., 14, 912–923, 2011.

Clarke, M. L. and Rendell, H. M.: The impact of North Atlantic storminess on western European coasts: a review, Quat. Int., 195, 31–41, 2009.

Courtier, P., Freydier, C., Geleyn, J.-F., Rabier, F., and Rochas, M.: The ARPEGE project at Météo-France, ECMWF Annual Seminar, Eur. Cent. for Medium-Range Weather Forecasts, Reading, UK, 1991.

Daniel, P., Josse, P., and Ulvoas, V.: Atmospheric forcing impact study in Meteo-France storm surge model, in: Coastal engineering V, computer modeling of seas and coastal regions, Wessex, UK: WIT Press, 2001.

Dee, D. P., Uppala, S. M., Simmons, A. J., Berrisford, P., Poli, P., Kobayashi, S., Andrae, U., Balmaseda, M. A., Balsamo, G., Bauer, P., Bechtold, P., Beljaars, A. C. M., van de Berg, L., Bidlot, J., Bormann, N., Delsol, C., Dragani, R., Fuentes, M., Geer, A. J., Haimberger, L., Healy, S. B., Hersbach, H., Hólm, E. V., Isaksen, L., Källberg, P., Köhler, M., Matricardi, M., McNally, A. P., Monge-Sanz, B. M., Morcrette, J.-J., Park, B.-K., Peubey, C., de Rosnay, P., Tavolato, C., Thépaut, J.-N., and Vitart, F.: The ERA-Interim reanalysis: configuration and performance of the data assimilation system, Q. J. Roy. Meteor. Soc., 137, 553–597, https://doi.org/10.1002/qj.828, 2011.

Della-Marta, P. M., Mathis, H., Frei, C., Liniger, M. A., Kleinn, J., and Appenzeller, C.: The return period of wind storms over Europe, Int. J. Climatol., 29, 437–459, 2009.
Ferreira, Ó., Ciavola, P., Armaroli, C., Balouin, Y., Benavente, J., Río, L. D., Deserti, M., Esteves, L., Furmanczyk, K., Haerens, P., Matias, A., Perini, L., Taborda, R., Terefenko, P., and Trifo, E.: Coastal storm risk assessment in Europe: examples from 9 study sites, J. Coast. Res., 1632–1636, 2009.

Gerritsen, H.: What happened in 1953? The Big Flood in the Netherlands in retrospect, Philosophical Transactions of the Royal Society of London A: Mathematical, Phys. Eng. Sci., 363, 1271–1291, https://doi.org/10.1098/rsta.2005.1568, http://rsta.royalsocietypublishing.org/content/363/1831/1271, 2005.

Hasselmann, K.: On the spectral dissipation of ocean waves due to white-cap, Bound.-Lay. Meteorol., 6, 107–127, 1974.

Hewson, T. D. and Neu, U.: Cyclones, windstorms and the IMI-LAST project, Tellus A: Dynamic Meteorology and Oceanography, 67, 27128, https://doi.org/10.3402/tellusa.v67.27128, 2015.

Kara, A. B., Helber, R. W., Boyer, T. P., and Elsner, J. B.: Mixed layer depth in the Aegean, Marmara, Black and Azov Seas: Part I: General features, J. March Syst., 78, 169–180, https://doi.org/10.1016/j.jmarsys.2009.01.022, 2009.

Klawa, M. and Ulbrich, U.: A model for the estimation of storm losses and the identification of severe winter storms in Germany, Nat. Hazards Earth Syst. Sci., 3, 725–732, https://doi.org/10.5194/nhess-3-725-2003, 2003.

Komen, G. J., Hasselmann, S., and Hasselmann, K.: On the existence of a fully developed wind sea spectrum, J. Phys. Oceanogr., 14, 1271–1285, 1984.

Li, N., Cheung, K. F., Stopa, J. E., Hsiao, F., Chen, Y.-L., Vega, L., and Cross, P.: Thirty-four years of Hawaii wave hindcast from downscaling of climate forecast system reanalysis, Ocean Model., 100, 78–95, 2016.

Marcos, M., Tsimplis, M. N., and Shaw, A. G.: Sea level extremes in southern Europe, J. Geophys. Res.-Oceans, 114, C01007, https://doi.org/10.1029/2008JC004912, 2009.

Mungov, G. and Daniel, P.: Storm-surges in the Western Black Sea. Operational forecasting, J. Med. March Sci., 1, 45–51, 2000.

Poli, P., Hersbach, H., Dee, D. P., Berrisford, P., Simmons, A. J., Vitart, F., Laloyaux, P., Tan, D. G., Peubey, C., Thépaut, J., Trémolet, Y., Hólm, E. V., Bonavita, M., Isaksen, L., and Fisher, M.: ERA-20C: An Atmospheric Reanalysis of the Twentieth Century, J. Climate, 29, 4083–4097, https://doi.org/10.1175/JCLI-D-15-0556.1, 2016.

Radnóti, G., Ajjaji, R., Bubnová, R., Caian, M., Cordoneanu, E., Von Der Emde, K., Gril, J., Hoffman, J., Horányi, A., Isara, S., Ivanovici, V., Janousek, M., Joly, A., Le Moigne, P., and Malardel, S.: The spectral limited area model ARPEGE-ALADIN, PWPR report series, 111–117, 1995.

Reistad, M., Breivik, Ø., Haakenstad, H., Aarnes, O. J., Furevik, B. R., and Bidlot, J.-R.: A high-resolution hindcast of wind and waves for the North Sea, the Norwegian Sea, and the Barents Sea, J. Geophys. Res.-Oceans, 116, C05019, https://doi.org/10.1029/2010JC006402, 2011.

Rivière, G., Arbogast, P., Maynard, K., and Joly, A.: The essential ingredients leading to the explosive growth stage of the European wind storm Lothar of Christmas 1999, Q. J. Roy. Meteor. Soc., 136, 638–652, 2010.

Rivière, G., Arbogast, P., Lapeyre, G., and Maynard, K.: A potential vorticity perspective on the motion of a mid-latitude winter storm, Geophys. Res. Lett., 39, https://doi.org/10.1029/2012GL052440, 2012.

Rogers, W. E., Hwang, P. A., and Wang, D. W.: Investigation of wave growth and decay in the SWAN model: three regional-scale applications, J. Phys. Oceanogr., 33, 366–389, 2003.

Rusu, L., Butunoiu, D., and Rusu, E.: Analysis of the extreme storm events in the Black Sea considering the results of a ten-year wave hindcast, J. Environ. Prot. Ecol., 15, 445–454, 2014.

Ryabinin, V. E., Zilberstein, O. I., and Seifert, W.: Storm surges, World Meteorological Organization, 1996.

Uppala, S. M., Källberg, P. W., Simmons, A. J., Andrae, U., Bechtold, V. D. C., Fiorino, M., Gibson, J. K., Haseler, J., Hernandez, A., Kelly, G. A., Li, X., Onogi, K., Saarinen, S., Sokka, N., Alvan, R. P., Andersson, E., Arpe, K., Balmaseda, M. A., Beljaars, A. C. M., Berg, L. V. D., Bidlot, J., Bormann, N., Caires, S., Chevallier, F., Dethof, A., Dragosavac, M., Fisher, M., Fuentes, M., Hagemann, S., Hólm, E., Hoskins, B. J., Isaksen, L., Janssen, P. A. E. M., Jenne, R., McNally, A. P., Mahfouf, J.-F., Morcrette, J.-J., Rayner, N. A., Saunders, R. W., Simon, P., Sterl, A., Trenberth, K. E., Untch, A., Vasiljevic, D., Viterbo, P., and Woollen, J.: The ERA-40 re-analysis, Q. J. Roy. Meteor. Soc., 131, 2961–3012, https://doi.org/10.1256/qj.04.176, 2005.

Usbeck, T., Wohlgemuth, T., Pfister, C., Volz, R., Beniston, M., and Dobbertin, M.: Wind speed measurements and forest damage in Canton Zurich (Central Europe) from 1891 to winter 2007, Int. J. Climatol., 30, 347–358, 2010.

Valchev, N., Andreeva, N., Eftimova, P., and Trifonova, E.: Prototype of Early Warning System for coastal storm hazard (Bulgarian Black Sea Coast), CR Acad Bulg Sci, 67, 977, 2014.

Wernli, H., Dirren, S., Liniger, M. A., and Zillig, M.: Dynamical aspects of the life cycle of the winter storm “Lothar” (24–26 December 1999), Q. J. Roy. Meteor. Soc., 128, 405–429, 2002.