Wind-wave characteristics and extremes along the Emilia-Romagna coast

. This study examines the wind-wave characteristics along the Emilia-Romagna coasts (northern Adriatic Sea, Italy) with a 10-year wave simulation for the period 2010– 2019 performed with the high-resolution unstructured-grid WAVEWATCH III (WW3) coastal wave model. The wave parameters (signiﬁcant wave height, mean and peak wave period, and wave direction) were validated with the in situ measurements at a coastal station, Cesenatico. In the coastal belt, the annual mean wave heights varied from 0.2–0.4 m, and the seasonal mean was highest for the winter period ( > 0 . 4 m). The Emilia-Romagna coastal belt was characterized by wave and spectra seasonal signals with two dominant frequencies of the order of 10 and 5–6 s for autumn and winter and 7–9 and 4 s for spring and summer. The wavelet power spectra of signiﬁcant wave height for 10 years show considerable variability,


Introduction
The wind-induced stress on the sea surface gives rise to wind waves that affect human activities on the coasts (Armaroli et al., 2019). The prevailing wind waves of a region determine the defence performance of coastal and offshore structures, and therefore precise information on wind waves is crucial for coastal operations and defence systems. During extreme events, the wind waves modify the total water-level elevation, leading to a higher risk of overtopping which can damage infrastructures. The Intergovernmental Panel on Climate Change (IPCC, 2007) has also highlighted the need for a long-term evaluation of wind-wave climate trends for coastal resilience (Hemer et al., 2012). The IPCC (2021) indicates the necessity of a regional evaluation of climate change, with various target factors that can aid in risk management and policy-making. The report points out that over the 21st century, nearshore regions will encounter sea level rise, thereby adding to more persistent coastal flooding (across low-lying regions) and associated coastal erosion.
Numerous studies have been reported for the Adriatic Sea, using numerical models to demonstrate the wind-wave climate characteristics. In the Adriatic there are many windwave forecast systems, including the Henetus forecast system described in Bertotti et al. (2011). Other state-of-the-art models include the Nettuno system as reported in Bertotti et al. (2013), which combines the atmospheric model COSMO (Steppeler et al., 2003) and the wave model WAM (Komen et al., 1994), and SWAN-MEDITARE, which combines COSMO with SWAN (Valentini et al., 2007;Russo et al., 2013). Donatini et al. (2015) have also implemented highresolution model chains for wind-wave forecasting in the Mediterranean and Adriatic seas, which use a combination of the atmospheric model WRF and wave model MIKE 21 (DHI, 2017). In a study over the Gulf of Taranto in southern Italy, a multi-nesting approach was adopted to evaluate coastal wave dynamics and hydrodynamics (Gaeta et al., 2016). In the Adriatic Sea, Sikiric et al. (2018) implemented the unstructured grid WAVEWATCH III (WW3) (WW3DG, 2016) with 2 km wind forcings from ALADIN forecasts (Farda et al., 2007). The study showed a good match with satellite measurements (SARAL) as compared to CryoSat-2 and Jason-2. The results were in agreement with the studies by Sepulveda et al. (2015), which showed that SARAL estimates of wave heights were far better than CryoSat-2 and Jason-2.  also reported on the application of SARAL data, producing good results.
In a study of the northern Adriatic, Lionello et al. (2012) used the WAM model to predict extreme wind waves and associated storm surge effects. In the Adriatic a modelling combination of WAM + SHYFEM (Komen et al., 1994;Umgiesser et al., 2014) forced with analysis and forecast ECMWF winds was used to forecast the 29 October 2018 storm (Cavaleri et al., 2019) conditions in northern Italy. The application of corrected forecast winds (ECMWF) within these models provided consistent results in line with measurements. High waves in the northern Adriatic Sea were reported in a recent study by Cavaleri et al. (2021).
Studies by Katalinić et al. (2015) reported that in the Adriatic basin, the wind speed and wave height increase from the northern to the southern areas with a maximum mean (annual) H s of 0.68 m. These results are underestimated as compared with the findings of Queffeulou and Bentamy (2007) resulting from a 14-year (1992-2005) satellite mission that revealed a mean H s of 0.85 m. Queffeulou and Bentamy (2007) also showed that in the Adriatic Sea, 80 % of the H s values were lower than 1.10 m. An intercompar-ison of WAM and WW3 models in the Adriatic and North Sea, based on testing various input physics, was reported by Benetazzo et al. (2021). The analysis aided in investigating the processes that lead to the generation of higher waves in the context of storms.
In the light of several hazardous and extreme events in the Emilia-Romagna (ER) coastal area, several studies have investigated the following: (i) coastal risk and vulnerability to flooding and erosion (Armaroli et al., 2009(Armaroli et al., , 2019Sekovski et al., 2015;Armaroli and Duo, 2018;Sanuy et al., 2018;Ferrarin et al., 2020); (ii) sea level rise, land subsidence, and littoral hydrodynamics (Perini et al., 2017;Gaeta et al., 2018); (iii) identification of storm thresholds (Armaroli et al., 2012); and (iv) forecasting of coastal flooding (Biolchi et al., 2021(Biolchi et al., , 2022. To the best of our knowledge, no studies have been carried out to date on the wind and wave characteristics and extremes in the ER coastal belt with high-resolution wind-wave models. Our study focusses on the prevailing wind-wave climatology in the coastal belt of the ER (northern Adriatic Sea) for a period of 10 years (2010-2019), the characterization of the wind-wave regimes, and the study of extreme wave conditions along the coastal belt to quantitatively determine the extreme wave hazard. We use a specific probability distribution function fitting procedure with the wind-wave model data and thereby extract hazard indices for different coastal points. We believe that our 10-year model simulation with appropriate validation at a coastal location will be useful for hazard estimations along the ER coastal area. For the first time we discuss the probability distribution of waves that are essential to quantify the extremes and their hazard.
The paper is organized as follows. Section 2 outlines the study area. Section 3 describes the wind-wave model used in the study, the model forcing, and the validation buoy data used. Section 4 describes the wind and wave climate in the ER coastal belt together with the wave spectra characteristics and wavelet analysis. Section 5 presents the analysis of the probability density distribution and the hazard index for extreme events. Finally, Sect. 6 summarizes the key findings from the study with a brief conclusion.

The study area -the Emilia-Romagna coast
The study area is the coastal waters of Emilia-Romagna, situated in northern Italy along the Adriatic Sea, with a coastline including natural zones and dunes to long stretches sheltered by groynes and breakwaters (Armaroli et al., 2012). The coastline is 130 km long (Harley et al., 2016) with the Po Delta as the northern boundary and the town of Riccione at the southernmost point. Figure 1 shows the study area in the ER coastal belt.
There are two main wind patterns in this region -the bora and sirocco winds (Pandzic and Likso, 2005;Umgiesser et al., 2021). Severe wind storms occur from the east-northeast, i.e. the prevailing direction of the bora winds. The sirocco winds are associated with low-pressure systems over the Italian Peninsula and the Ionian Sea. Owing to the restricted fetch, i.e. limited extension of the wave generation area, the bora winds generate young and steep waves that break frequently , while the sirocco winds generate longer-fetch waves across the Adriatic Sea (Cavaleri, 2000). Thus the swell seas are controlled by the sirocco winds and the seas are dominated by the bora winds (Bonaldo et al., 2017).
The ER coastal area is subdivided into three major zones ( Fig. 1c) which correspond to different coastal trophic conditions (Fiori et al., 2016). The station locations are the land town locations perpendicular to which the environmental agency monitoring transects are carried out monthly and weekly to monitor the marine ecosystem environmental status. Thus, knowing the prevailing winds and waves at these locations could be of importance for the management of this important coastal area.
The prevailing hydrodynamics show that the region is microtidal with spring tides (80-90 cm) and neap tides (30-40 cm), with strong diurnal and semi-diurnal components (Armaroli and Duo, 2018). A low-energy wave climate (IDROSER, 1996;Ciavola et al., 2017) has been reported along the coastal belt of ER, i.e. 60 % H s < 1 m. Armaroli et al. (2012) reported that waves originating from the east correspond to a proportion of 91 % H s < 1.25 m, owing to the controlled fetch.

Numerical wave model set-up
In this study, the third-generation unstructured-grid spectralwave model, WW3 (version 5.16; WW3DG, 2016) was used to evaluate nearshore waves. WW3 is a universally accepted wave model (Tolman et al., 2002) with continuous updates of ocean wave physics. The model is formulated by solving the action-density balance equation: (1) The left-hand side of Eq. (1) denotes the changes in wave action density (i.e. local rate), generation in physical space, shifting of action density (frequency/direction) owing to spatio-temporal changes in depth, and current. λ denotes longitude, ϕ latitude, θ direction of wave propagation, and k wave number, and σ and t represent the intrinsic angular frequency and time respectively. The source term used in this paper, S in Eq. (1), is the wind input and dissipation source package ST4 (Ardhuin et al., 2010) or ST6 (Zieger et al., 2015;Rogers et al., 2012;Babanin, 2011), the bottom friction JONSWAP (Joint North Sea Wave Project) parameterization (Hasselmann et al., 1973), or SHOWEX (Shoaling Waves Experiment) formulation (Ardhuin et al., 2003) for sandy bottoms. In the section on sensitivity experiments we used a combination of these source terms. The WW3 model grid ( Fig. 1) is divided into 15 392 elements, linked with 8148 nodes, with a resolution of about 300 m at the coast and 2.5 km at the open boundary (Fig. 1a). The merged European Marine Observation and Data Network (EMODnet) data (250 m resolution) and multibeam high-resolution measurements from Arpae (Regional Agency for Prevention, Environment and Energy of Emilia Romagna) serve as the bathymetry of the ER domain (Fig. 1b). The model spectrum is sampled in 24 directions and 30 frequencies (0.0500-0.7932 Hz), with an increment factor of 1.1. The model time steps are set as the (i) maximum global time step, 200 s; (ii) maximum CFL (Courant-Friedrichs-Lewy) time step X-Y , 50 s; (iii) maximum CFL time step k-θ , 50 s; and (iv) minimum source term time step, 10 s. The source term for linear input and wind input uses the parameterization formulated by Cavaleri and Malanotte-Rizzoli (1981) and Donelan et al. (2006). The Generalized Multiple DIA (GMD) was used to simulate the non-linear interactions (Tolman, 2010(Tolman, , 2013(Tolman, , 2014; the dissipation physics were based on Rogers et al. (2012); and the SHOWEX formulations by Ardhuin et al. (2003) were used to simulate the bottom friction. The SHOWEX parameterization is ripple-induced bottom friction, which considers the formation of sand ripples on the bottom. Breaking (depth-induced) is activated using Battjes and Janssen (1978) physics.
The WW3 Table 1 showed correlations of the order of 0.7, with bias of −0.2 m s −1 indicative of underestima-   (Yamaguchi, 1984) based on the CMEMS wave parameters (significant wave height, peak period, and mean direction).

Wave observational data set and validation method
In order to validate the model hindcasts, we used the wave buoy Nausicaa in Cesenatico (44.21 • N, 12.47 • E; station 6) as shown in Fig. 1c. This station is situated away from the coast of Cesenatico municipality and is supported with a Datawell Directional Waverider (MkIII 70 wave) buoy called Nausicaa (https://www.arpae.it/it/temi-ambientali/ mare/dati-e-indicatori/dati-boa-ondametrica, last access: 20 March 2022) which has been maintained by Arpae since 23 May 2007. The location of the buoy is 8 km offshore of Cesenatico, at a depth of approximately 10 m, in a region inaccessible to fishing, navigation, and moorings. Wave data such as height (H s ), period, and the direction of waves every 30 min constituted the basic validation data set for the modelling period from January 2010 to December 2019.
Wave model parameters such as wave height, period, and direction were extracted and analysed for eight control points as shown in Fig. 1c. The details of the control points are described in Table 2. The model-simulated 1D wave spectra are extracted and analysed based on seasons. The skill of the model to reproduce the observations at the Nausicaa buoy location was assessed by standard statistics, namely the correlation coefficient (R), bias, and root mean square error (RMSE): where model estimates are denoted by P , O represents observational data, n indicates number of data points, and the overbar denotes mean values.

Sensitivity experiments for wave model parameterizations
Three sets of sensitivity experiments using WW3 were executed using a combinations of wave physics: for the representative months of February (winter) and September (autumn) 2018. The combination of ST6 with JONSWAP is not considered because this bottom friction is not suitable for sandy beaches as EXP1 will show. The H s comparison with the Nausicaa buoy is shown in Table 3, which highlights that the best physics is given by EXP3. The mean buoy H s is 0.92 and 0.38 m for February and September 2018 respectively. The comparison of the mean wave period, T m (not shown), for the three experiments showed a higher performance using the combination of ST6 + SHOWEX. The sensitivity study produced sufficient confidence in using the ST6+SHOWEX physics for the ER region. This wave physics was thus adopted for the 10year simulation.

Validation of wave hindcasts
The model outputs, such as significant wave height (H s ), mean wave period (T m ), peak wave period (T p ), and mean wave direction (θ m ), were compared with the buoy observations for the 10-year period 2010-2019. Figure 2 shows the 10-year comparison of H s , which qualitatively demonstrates that the overall model H s also followed the buoy values in peak events at the Cesenatico station (station 6 in Fig. 1c). This is a consistency check of model against observations as required for "goodness" indicators in numerical weather predictions (Murphy, 1993). The model also captures the seasonal variations at the coastal location. In general, the lower H s values are slightly overestimated, while higher H s values are underestimated. Table 4 shows the validation statistics for each year. The mean of model/buoy estimates for H s , T m , and T p are 0.40/0.45 m, 3.02/3.23 s, and 4.17/4.36 s respectively. On average the model underestimates the measurements (as seen from the negative bias for most of the years). A high correlation is shown, ranging from 0.81 to 0.93 for 2010-2019, with the highest correlation for 2017. The T m comparison revealed a lower correlation of the order of 0.72 to 0.81, compared to H s . The negative bias (−0.371 to −0.018 s) indicated an underestimation of T m , with a corresponding RMSE of the order of 0.79 to 0.91 s. Similarly, the T p also showed a lower correlation (0.53 to 0.70) in comparison to H s and T m . T p also showed underestimations as revealed from the bias of the order of −0.382 to 0.151 s, with an RMSE varying from 1.48 to 1.78 s. Figure 3 represents the observations-model scatter plot of H s for the period 2010-2019 (Fig. 3a) and the seasonal scatter as shown in Fig. 3b-e for station 6. The dashed red line denotes the best data fit for the comparison. The comparison of H s for 2010-2019 (Fig. 3a) shows that there is relatively good agreement between model H s and measurements with a high correlation of 0.90. There is a slight underestimation (bias = −0.05 m), with RMSE = 0.21 m. The seasonal scatters for winter, spring, and autumn (Figs. 3b, c, and e) showed high correlations, with a slight underestimation in relation to buoy observations. The summer seasons (Fig. 3d) showed a comparatively lower correlation with an underestimation of H s . In general, the model H s underestimates the buoy data, specifically the higher H s , and similar underestimations have been reported in many past studies such as Ardhuin et al. (2007), Korres et al. (2011), andClementi et al. (2017).
The comparison of T m for 2010-2019 is shown in Fig. 3f and for the seasons in Fig. 3g-j, revealing a larger scatter in comparison to H s . During 2010-2019 (Fig. 3f), the simulated T m is lower than the buoy measurements and shows a lower performance (R = 0.75) in comparison to H s . The winter, spring, and autumn seasons (Fig. 3g, h,  Below we present the wind climatology in the ER region based on the ECMWF analysis winds over a period of 10 years. The seasons are presented as winter (December-January-February), spring (March-April-May), summer (June-July-August), and autumn (September-October-November).

Climatology of wind speed and direction
Analysis of wind speed and direction over the ER coast for the period 2010-2019 is presented in Fig. 4. The annual mean characteristics showed a very precise pattern, with the winds reaching the coast from the east-northeast. The annual mean wind speeds were of the order of 0.5-2 m s −1 , with a large standard deviation (SD) of 1.6-3 m s −1 . The lowest wind speeds were observed during spring and summer (1.5 and 1.8 m s −1 ), followed by autumn (2.4 m s −1 ), with the highest wind speeds (2.9 m s −1 ) during winter. Overall, for the winter and spring the approaching wind is easterly, related to the bora wind climatological direction. In the summer, the mean wind direction is from the southeast, owing to sirocco events. The spatial distribution of the seasonal and annual SD of wind speed from 2010-2019 is shown in Fig. 4f-j. The annual SD varies from 1.6 to 3 m s −1 in the entire domain (Fig. 4j), and the annual maximum is further offshore from the ER coastal belt. During the winter (Fig. 4f), the SD varies from 2.2 to 3 m s −1 and, in spring, from 1.2-3 m s −1 (Fig. 4g). While in summer and autumn, the SDs are 2.2-2.6 and 1.6-3 m s −1 respectively. To better study the wind characteristics along the ER coast, the wind rose diagrams are shown for the eight control points in Fig. 5a to h. Points 1 to 5, belonging to Zone A and Zone B, have the highest wind speeds, approaching at an angle of 45 to 135 • . The wind speed ranging from 3 to 4 m s −1 is more frequent at these control points at an approaching angle ranging from 45 to 90 • . The average coastal angles of Zone A and Zone B are nearly 45 • . Points 6 to 8 fall along the concave side of the coastal area, i.e. in Zone C. Along these control points, the maximum wind speed approaches from W to NNW. The wind speeds up to 3.5 m s −1 show a marked increase in frequency. The frequent wind speeds approach from NW and ENE for station 6, NNE for point 7, and NNW for point 8. Moving from point 1 to 8, there is a gentle shift in the maximum wind speed approaching from NNE to ENE.

Wave height and direction climatology
Figure 6 (top panels) depicts the annual mean H s for the ER coast and the seasonal H s means for winter, spring, summer, and autumn. The SD for each event is illustrated in the bottom panels from 2010 to 2019. The waves converge at the southern and the northern part of the study domain due to the shape of the coastline. There is divergence in wave energy in the middle region of the coastal domain (i.e. Zone B as reported in Fig. 1c). The annual H s mean (Fig. 6e) in the domain varied from 0.08-0.6 m. The annual average H s is higher (0.5-0.7 m) off the ER coast and at the boundary in the open ocean, and in the central ER domain H s is of the order of 0.5-0.6 m. However, in the ER coastal belt, the annual mean H s is < 0.4 m owing to the bathymetric features.
The seasonal climatology of H s in the winter season ( Fig. 6a) indicates higher waves offshore of the order of 0.1-0.9 m, where the ER coastal belt has H s < 0.5 m. In spring (Fig. 6b) and summer (Fig. 6c) the H s values are comparatively lower and vary in the range of 0.1-0.5 and 0.1-0.33 m respectively. The autumn H s mean in the ER coastal belt is < 0.4 m. The spatial H s field structure and direction approximately resemble the bathymetric contour lines (Fig. 1b). The annual SD (Fig. 6j) varies from 0.09-0.71 m in the ER domain. The summer season (Fig. 6h) shows the lowest SD (0.1-0.38 m) compared to all other seasons.
The detailed features of the model in the coastal zone are shown by means of wave rose diagrams (Fig. 7) for the eight points in Fig. 1c. The waves at control point 1 fall in the Lido di Pomposa region where the coast is sheltered and exposed to winds and marine currents. The bathymetric contour enables the waves to converge in control point 1, where the maximum wave heights approach from E to SE. From points 2 to 7 along Porto Garibaldi to Rimini, the approaching angles of wave heights are from NE to SE. The maximum waves approach from ENE to E for points 2 to 4 and NE to E for points 5 to 8. The maximum wave activity is observed at point 3. Point 1 is a relatively calm area compared to the other control points, perhaps because of the shadow zone. The concave shape of the coast, well represented by the highresolution unstructured-grid model, and bathymetric patterns   are key to understanding the prevailing wave characteristics in the ER coastal belt. The wave energy converges at the end points and diverges at the middle points.
Based on available buoy data for the Cesenatico station, the observed wave roses are compared with the model estimates for selected years as shown in Fig. 8. Overall, the modelled wave roses (Fig. 8b, d, and f) show reasonable correspondence with the observed data (Fig. 8a, c, and e), even with some difference in magnitudes. An underestimation of model wave heights in the lower ranges is noted. Comparing the directional distribution of waves, the directions are comparable and in the same sectors but there exist higher differences in their magnitudes. A similar wave climate by the Nausicaa buoy located offshore of Cesenatico is reported in studies by Armaroli et al. (2012) and Romagnoli et al. (2021), which shows that this is the representative wave climate of the Emilia-Romagna coast. This qualitative comparison shows that at the Cesenatico station, overall characteristics of waves are fairly reproduced by the model. Figure 9 shows the offshore wave climate, presented as wave rose diagrams at the control points along the boundaries of the study domains (control points 9 to 14). In Fig. 9a at point 9, the waves approach from NE to SE with maximum H s approaching from ENE to ESE. At point 10, the predominant waves are at an angle of 30 to 150 • where the maximum H s approaches from NE and SE directions (see Fig. 9b). For points 11 to 14, the predominant wave directions are from 30 to 150 • , where the maximum H s approaches from NE and SE directions. Deep-water control points 10 to 14 receive waves from all directions.

Wave spectra characteristics
In the ER region, there are hardly any studies on the spectral characteristics of the waves. Cavaleri et al. (2019) analysed model spectra for the event of 29 October 2018 in the northern Adriatic Sea and compared them with measurements on the Venice coastline. The simulated wave spectra on the 25th of the months corresponding to winter (February), spring (May), summer (August), and autumn (November) at 12:00 LT are represented in Fig. 10a to d for station 6 for 2010-2019. Figure 10a shows the simulated instantaneous spectra in February (25th of the month, 12:00 LT) with the highest peak energy of 2.0234 m 2 Hz −1 for 2018 and the lowest of 0.0008 m 2 Hz −1 for 2012 and 2014. February, which is a representative month of the winter season, shows a combination of single-peaked and double-peaked spectra with swell dominance at the coastal location. In all the seasons, the swell dominates the spectral energy with a peak at around 0.11 Hz (9 s). The shorter wave peaks range from 0.21 to 0.54 Hz (1.8 to 4.7 s). The spectra vary considerably over the years, and in general, the wave spectra at the Cesenatico coastal location showed signatures of single-and double-peaked spectra for the period 2010-2019 (Table 5). The wave spectra were prominently double-peaked during all seasons (45 %-53 %), along with single-peaked spectra but with a lower percent-  age of occurrences (27 %-33 %). Double peakedness was highly prominent in the summer season (53 %), while winter, spring, and autumn showed dominance of single-peaked spectra (31 %-33 %). As evident from Table 5, the percentages of the number of peaks (single/double) in the Cesenatico location clearly depict the co-existence of sea-swell characteristics in the study domain. The monthly mean wave spectra for winter, spring, summer, and autumn corresponding to the typical months of February, May, August, and November for 2010-2019 are represented in Fig. 11a to d. During February (Fig. 11a), the averaged spectra showed prominent single peaks for most of the years with peak energies of the order of 0.1615-0.722 m 2 Hz −1 . The highest peak energies were in 2012 (0.701 m 2 Hz −1 ) and 2014 (0.722 m 2 Hz −1 ), and during the 10-year period the peak frequencies ranged from 0.0974 to 0.1726 Hz. Figure 11b shows the averaged spectral characteristics for May (spring). As seen from Fig. 11b, 2019 had the highest peak energies of 0.253 m 2 Hz −1 , and the spectra also highlight a few secondary peaks in some of the years with the peak frequency ranging from 0.1072 to 0.2297 Hz. During the summer season (August), the spectra show single/double peaks with peak energies varying from 0.0102 to 0.0686 m 2 Hz −1 (Fig. 11c). The maximum peak energy was for 2016 (0.0686 m 2 Hz −1 ) with comparatively low energies for the rest of the years and with peak frequencies varying from 0.1427 to 0.278 Hz. Similarly, during autumn (Fig. 11d), the averaged spectra were mostly single-peaked with peak energies of the order of 0.1362 to 0.740 m 2 Hz −1 . The highest peaks with energies of 0.740 m 2 Hz −1 were in 2019, with the lowest energy in 2015. The dominant frequencies corresponding to the peak energies were 0.0974-0.2089 Hz.
Overall, the highest and lowest spectral peaks are in winter and summer, with energies of 0.722 and 0.0686 m 2 Hz −1 , as shown in Fig. 11a and c. The mean wave spectra for 2010 to 2019 exhibit a peak in variance for 2014, 2019, 2016, and 2019 for winter, spring, summer, and autumn respectively. The spectra show more or less similar characteristics for spring and autumn. There is also a reversal of spectrum curves for winter and spring as swells clearly dominate the coastal location. The spreading of spectra, which is dependent upon the wind characteristics and the prevailing fetch, is variable during all seasons.

Wavelet analysis
Wavelet analysis is an important tool to analyse spectral components and occurrence time (Torrence and Compo, 1998). The wavelet considers spectral components' time localization, as well as time-frequency rendering of signal into realization, such that the frequencies in the wavelet analysis are associated with the time domain. Thus, wavelet analysis (based on the Morlet mother wavelet) provides an understanding of spectral characteristics and their variability in time.
In this study the wavelet transform for H s (Fig. 12) was applied to the coastal location of Cesenatico for 2010-2019, using the mean model estimates. Figure 12a represents the wavelet power, with the x axis representing the time and y axis denoting the component periods. Figure 12b represents the global wavelet power spectrum, i.e. time-averaged power spectrum, which uses the same y axis. Cesenatico was selected as it was the station where the wave parameters were validated with the model estimates. The idea of presenting the wavelet transform is to accurately represent the variance in the spectrum. In the power wavelet (Fig. 12a), the real signals can be observed enclosed in the black contours with a 95 % confidence level, while the region below the dashed magenta line indicates the cone of influence, in which the time series analysis edge effects are significant. In the global spectrum, the peaks indicate the combined signal throughout the analysis. The dashed red line in the global spectrum corresponds to a confidence level of 95 %.
In Fig. 12, the largest signal occurs in the 256-512 d band, which contains the seasonal frequency, and sporadic signals can be identified by comparatively short times (2-3 months). Figure 12a indicates that over the 10-year period, intermittent oscillations are in band 16-128 in the years 2011, 2012, 2014, 2016, 2018, and 2019. Figure 12c

Extreme-wave analysis
In this study, the statistical characteristics of H s were analysed using the methodology of fitting a probability distribution function (PDF) to the wave time series at the control points of the ER coastal belt. Many studies have indicated that the probability distribution used to model longterm distributions of wave heights is well represented by the two-parameter Weibull distribution (Muraleedharan et al., 1993(Muraleedharan et al., , 1998(Muraleedharan et al., , 1999. The PDF of a random variable x with the Weibull distribution (Weibull, 1951) is defined for positive values, x > 0, as where κ and λ (> 0) are the shape parameter (dimensionless) and scale parameter (m) respectively. It is clear that when κ = 1, the PDF reduces to an exponential distribution. Fitting this PDF to the data enables the hazard index to be calculated, which is the probability that the waves will exceed a threshold, let us say x C in H s . The hazard index is then defined as To compute the best-fit shape and scale parameters of the Weibull distribution for each of the eight control points, the maximum likelihood method (MLM) was used. This is the most widely used technique among parameter estimations and finds a value of the parameter that maximizes the likelihood function. The values of the Weibull parameters for each control point are presented in Table 6, which shows that the mean, standard deviation, and skewness computed from the model data are very similar to those estimated from the Weibull fit parameters. This thus highlights that the Weibull distribution well represents the behaviour of the H s model data. The mean value and the corresponding variance of H s at the Cesenatico station are larger than for the other control points, as the station is far from the coast with the highest water depth. The analysis results show that the fitted Weibull distributions have positive kurtosis, which indicates that the distribution has fat tails.
To evaluate the goodness of fit of the Weibull distribution, the classical chi-squared (χ 2 ) test was used. This test determines how well the theoretical distribution fits the given model data distribution. If the chi-squared value is lower than a critical χ 2 value, we retain the null hypothesis and conclude that there is no significant difference between the observed  Table 6. The best-fit Weibull scale and shape parameters for H s (columns 3-4) at the eight control points. Column 5 shows the estimated χ 2 values. Mean, variance, skewness, and kurtosis of H s (columns 6-9) computed from the model data (left sub-columns) and from the Weibull fit parameters (right sub-columns) and the wave height hazard index (with threshold value X c , H s = 1.08 m) calculated with (Eq. 6) for the eight control points (indicated in the last column) along the Emilia-Romagna coastal strip. and the expected distributions. The estimated χ 2 values for each control points are given in Table 6. The decision rule for the χ 2 test depends on the level of significance (set to 0.05) and the degrees of freedom, defined as df = N − np, where N is the number of bins (set to 30) and np is the number of distribution parameters (i.e. 2), so the critical value of χ 2 is 41.34 (taken from the χ 2 distribution table). Table 6 high-lights that the two-parameter Weibull distribution fits the H s data well. In all the eight locations, since the shape value κ is close to 1, the fitted Weibull distributions behave like the exponential distribution. Figure 13 compares the Weibull distribution fit (red line) and the histogram of the model data for three relevant locations (Porto Garibaldi, Lido Adriano, and Ce- for location) using mean model estimates. (a) Wavelet power spectrum for H s . The colour bar stands for the formation of H s variation. Power spectra intensity is represented by colours varying from navy blue colour (i.e. weak) to dark yellow (i.e. strong). The contours represent the total variance percentage, and the black contours indicate amplitude significance (greater than 95 % level). The dashed magenta line is the cone of influence (region of spectrum with the significant edge effects), where zero padding has reduced the variance. Figure 11a shows that power is concentrated in the 256-512 d band, which is a strong signal. (b) Global wavelet power spectrum, where the blue curve indicates the fast Fourier transform of the complete data. The dashed red line is the significance (95 %) for the global wavelet spectrum, assuming the same significance level and background spectrum as in wavelet power spectra. (c) Scaled averaged time series over a 16-512 d band showing variance of H s . The dashed red line is the 95 % confidence level for H s . senatico). The two-parameter Weibull distribution appears to fit the data well in the coastal study area.
After evaluating the Weibull distribution fit and the statistical moments, we estimated the hazard index as shown in Table 6 (column 10) for a threshold value X c (i.e. H s = 1.08 m, 3 times the mean standard deviation). The hazards were shown to increase 7-fold from the northern control points (Lido di Pomposa) to Cesenatico and then to decrease again. In the future it will be interesting to compare the hazards for different coastal areas around the Adriatic Sea.

Summary and conclusions
To accurately simulate the wind-wave climate in the Emilia-Romagna coastal belt, a high-resolution numerical modelling study using unstructured-grid WW3 was executed for a 10year period. The WW3 model was driven by the ECMWF analysis winds, and the model was validated with available wave buoy data at a coastal location. The sensitivity tests have shown the good accuracy of ST6+SHOWEX physics for wave hindcasts in the study area. The results of a comparison of model estimates with measurements were promising. An H s correlation of 0.86 to 0.93 was found for the 10year simulations with observed data. The underestimations in H s were indicative of a negative bias (−0.076 to −0.016 m) with an RMSE of 0.19 to 0.25 m. The comparison of T m and T p with buoy measurements revealed a correlation of 0.70 to 0.80 and 0.53 to 0.70 respectively, for the 10 years. Our database was then used to study and characterize the present climate of waves for the region, and a hazard index for extreme events was defined and computed. The following conclusions were drawn: -The spatial mean wind speed for winter, spring, summer, and autumn varied in the range 1.1-2.9, 0.5-1.5, 0.5-1.8, and 0.5-2.4 m s −1 respectively. The lowest wind speeds were observed during spring and summer (1.5 and 1.8 m s −1 ) considering the study domain, followed by autumn (2.4 m s −1 ), and the highest wind speeds (2.9 m s −1 ) were observed during winter seasons.
-The annual H s mean in the ER domain varied from 0.08-0.6 m, and in the ER coastal belt the annual mean H s was < 0.4 m owing to bathymetric features. In the ER coastal belt, the seasonal climatology of H s in the winter showed a mean H s < 0.5 m, while in spring and summer the H s values were comparatively low (H s < 0.38 and H s < 0.21 m). The autumn H s mean was < 0.4 m. It should be noted that there was more wave activity in winter and autumn than in spring and summer. -The analysis of instantaneous spectra showed that during all seasons the spectra exhibited bi-modal characteristics (double-peaked), with a dominant occurrence percentage of 53 % in summer, while during spring, winter, and autumn the spectra also showed single-peaked spectra (31 %-33 %). The average spectra analysis showed peak frequencies of the order of 0.097 to 0.172 Hz, 0.107 to 0.229 Hz, 0.142 to 0.278 Hz, and 0.097 to 0.208 Hz for winter, spring, summer, and autumn seasons respectively.
-With the aid of a wavelet analysis tool, the power features (time-frequency) of H s data showed substantial variability in H s for a 10-year period, with the occur-rence of monthly and seasonal periods. The 256-512 d band showed a higher power concentration which represents the seasonal frequency.
-The coastal control point time series was well fitted by a Weibull PDF. The Weibull at all control points was congruent with an exponential distribution. Using the Weibull PDF fit, we calculated a hazard index which indicated that for waves higher than 3 standard deviations from the mean, i.e. the highest hazard is reached at the Cesenatico station (hazard index of 0.07 (Fig. 13c) instead of 0.06 (Fig. 13b) and 0.04 (Fig. 13a)).
The limitation in the present study is the non-availability of wave spectra measurements at the coastal locations for val-idation. Future study will aim to consider data assimilation (de Rosnay et al., 2022) and to have higher-resolution winds as forcings for the wave model. In the context of heavy-tailed data sets, the Weibull distribution may not represent the best description of the peak and tail. These limits can be overcome by adding more parameters such as the four-parameter exponentiated Weibull distribution (Mahmoudi et al., 2018) such that the extra shape parameter can provide more versatility regarding the distribution in the shape of the tails. Another limitation is that a 10-year period would generally not be enough to bring out the complete climatological wave response to winds. ECMWF winds were too low resolution before 2010, and no downscaled limited-area meteorological forcing is available. Hence this 10-year period could be the first reference database for the prevailing wind-wave characteristics in the coastal belt for researchers and coastal engineers/designers. Future works should definitely deal with longer simulation periods and also higher-resolution winds.
Our analysis highlights the importance of long-term wave databases, which can aid in the design requirements of coastal engineering applications. It also demonstrates the useful application of PDFs to the estimate of hazards along coastal belts. The study also highlights the need for extensive wave spectra comparisons (Lobeto et al., 2021) with measurements for selected locations on the coastal belt, which will update the coastal wave database. The early detection of hazards such as coastal erosion and associated shoreline changes are still challenging (Le Cozannet et al., 2020) due to the non-availability of long-term observations. As reported by Vousdoukas et al. (2018), by the end of the century, the community encountering marine flooding is estimated to rise from 1.52 to 3.65 million, and considering global vulnerability (Luijendijk et al., 2018), low-lying nearshore regions (one-fourth) are retreating, and eroded land  remains twice what is acquired. Better knowledge of the prevailing wave characteristics on the ER coastal belt will aid in predicting coastal impacts.
Code and data availability. The data and codes used in this study can be accessed at the following Zenodo archive: https://doi.org/10.5281/zenodo.6360348 (Pranavam Ayyappan Pillai et al., 2022).
Author contributions. UPAP contributed in terms of conceptualization, data curation, investigation and analysis, methodology, validation, visualization, results interpretation, and writing the original draft; NP contributed in terms of project supervision, conceptualization, methodology, results interpretation, and writing -review and editing; IF contributed in terms of data curation, methodology, and writing -review and editing; SC contributed in terms of data curation, methodology, and writing -review and editing; FT contributed in terms of data curation, methodology, results interpretation, and writing -review and editing; SU contributed in terms of data curation and writing -review and editing; AV contributed in terms of data curation and writing -review and editing.