2.1.1. Comparison of Offshore Wave Data

The assessment of the whole dataset available for the Ponza wave buoy and for ECMWF time (point E3) is graphically represented with polar diagrams assembled in Figure 2. The wave dataset obtained from the Ponza wave buoy by means of the "geographic transposition of wave gauge data" to point O is also reported. The geographic transposition has been applied according to the method originally formulated by Contini and De Girolamo [96].

**Figure 2.** Wave climate referenced to different wave height classes (in legend): (**a**) at Ponza buoy; (**b**) at point O, as obtained by transposition of Ponza wave buoy data; (**c**) at ECMWF point E3.

The method is based on the following hypotheses:


*Water* **2020**, *12*, 1936

Under the above conditions, the spectral significant wave height *Hm0* and the peak period *Tm* can be estimated using the Sverdruv Munk Bretschneider (SMB) method [98]:

$$\frac{\text{g} \times H\_{m0}}{\text{l}I\_A^2} = 1.6 \times 10^{-3} \times \left| \frac{\text{g} \times F}{\text{l}I\_A} \right|^{\frac{1}{2}} \tag{1}$$

$$\frac{g \times T\_{\rm m}}{l I\_A^2} = 2.857 \times 10^{-1} \times \left| \frac{g \times F}{l I\_A} \right|^{\frac{1}{3}} \tag{2}$$

where *g* is the gravity acceleration, *UA* represents the wind-stress factor and *F* is the effective fetch.

Equations (1) and (2) can be written at real and virtual sites and, under the assumption that wind conditions are the same for both stations, the following equations can be derived:

$$\frac{H\_{m0,V}}{H\_{m0,R}} = \left| \frac{F\_V}{F\_R} \right|^{\frac{1}{2}} = K\_H \tag{3}$$

$$\frac{T\_{m,V}}{T\_{m,R}} = \left| \frac{F\_V}{F\_R} \right|^{\frac{1}{2}} = K\_T \tag{4}$$

where the subscripts *R* and *V* denote the variables referring to the real and virtual station, respectively.

The "transposition coefficients" *KH* and *KT* allow us to calculate from the real wave buoy records the transposed wave gauge dataset at the virtual station. For the sake of completeness, the virtual station O was selected between the ECMWF grid points E3 and E4 before the hindcast data analysis. In this way, once we selected the best ECMWF reference point from one of the two points, a sufficient comparison with the transposed dataset was ensured. More than 65% of the annual wave energy comes from the sector 220◦–280◦, in accordance with the long fetch facing the Gulf and with the mesoscale climate conditions, with swells approaching from distant storms coming from the NW sector of the Mediterranean Sea.

The comparison between the ECMWF dataset and buoy records (both real and transposed) shows some differences. In particular, the rate of waves coming from the east is significantly reduced. Then, the lowest values of energy in the ECMWF points are noticeable, particularly in the highest power class. It can be seen that the average wave power moves from 3.85 kW/m computed at the Ponza wave buoy (4.73 kW/m for point O) to 2.19 kW/m at point E3 (Table 2).


**Table 2.** Main wave climate parameters (based on the whole datasets) at Ponza wave buoy, point O and ECMWF grid point E3.

Moreover, Table 3 shows the average differences between the measurements carried out by the buoy (at Ponza and after transposition at point O, respectively) and the hindcast data. Such values of Hs,PONZA BUOY / Hs,ECMWF and Hs,POINT O / Hs,ECMWF are organized by wave class (in terms of *Hs* ranges).


**Table 3.** Average differences between the buoy records (at Ponza and transposed at point O, respectively) and ECMWF hindcast data for point E3.

It is possible to note that for calm conditions (*Hs* < 0.5 m), the datasets are very similar. The highest discrepancy was found for 0.5 m < *Hs* < 0.75 m, where values of *Hs* recorded by the Ponza wave buoy (point O, respectively) were on average 1.75 times (1.77, respectively) higher than those reported for the E3 hindcast data. The overall mean discrepancy between the Ponza buoy and the ECMWF data was 1.37, while between point O and ECMWF, it was 1.42. A comparison of the wave height time series obtained from the different datasets is highlighted in Figure 3.

**Figure 3.** Comparison of the wave height obtained from buoy records, transposition of buoy records at point O and ECMWF data for the reference point E3.

The bulk of these differences can be attributed to the dissimilar measurement conditions. The smaller sampling frequency for the hindcast data involves peak attenuation, acting as a band-pass filter and smoothing the signal. The underestimate of the ECMWF data was previously highlighted within the WW-Medatlas projects [28]. Moreover, through intercomparison with NCEP (National Centers for Environmental Prediction) Climate Forecast System Reanalysis [26] and with wave buoy data, the lowest values of energy in the ECMWF points were detected, especially in the highest power class (e.g., [6,10,27,36]). Hence, the use of the ERA-Interim dataset could be considered adequate for slightly conservative wave power potential and studying long-term variations in wave height [10] but, at the same time, should be examined carefully during detailed resource assessments or for arriving at the design wave condition or to build a detailed nearshore wave model. The main parameters of the wave climate at each grid point are reported in Table 4.


**Table 4.** Main wave climate parameters (based on 39-year average) at ECMWF grid points.

A tentative contour map (based on interpolation of power rate at 12 grid points) has been provided in Figure 4, where wave power isolines are depicted, ranging from 2.5 to 5 kW/m.

**Figure 4.** The 18-year averaged energy flux for the 12 ECMWF grid points and contour lines of the estimated mean wave power flux per unit crest on the central and southern Tyrrhenian Sea.

#### *2.2. Study Site and Nearshore Wave Instrumentation*

The nearshore study site is represented by the Bagnoli-Coroglio Bay, located within the Gulf of Naples, a natural semi-enclosed embayment within the Gulf of Pozzuoli (also known as the Gulf of Puteoli). Its mean depth is ca 60 m, with a maximum depth of 110 m and a surface area of 33 kmq.

Due to the proximity to the city of Naples, the whole area historically represents one of the best studied coastal areas of the Mediterranean Sea [99].

Thanks to the presence of the Stazione Zoologica "Anton Dohrn" (SZN) since 1872, marine investigations have been carried out for more than a century and half [100]. Recently, two Monitoring and Environmental Data Units (MEDA) of the SZN have been installed in the Gulf of Naples and in the Gulf of Pozzuoli (Figure 5). These MEDA units are mainly used for the chemical, biological and environmental monitoring of the marine ecosystem and both are equipped with an ADCP [61]. This shallow marine area is also famous, as it is the most highly active volcanic district in the coastal zone of SW Italy [101]. The geographical coordinates and water depth for MEDA A (Gulf of Pozzuoli) and MEDA B (Gulf of Naples) are indicated in Table 5. Close to MEDA B, a DWSD wave buoy (DWSD-B hereafter), provided by the University of Campania, was installed during a field campaign which took place from May to June 2016 [93]. Instead, in the period of February–March 2017, a DWSD buoy (DWSD-A) was placed close to MEDA A.

**Figure 5.** (**a**) Map of the Gulf of Naples and location of the study site. (**b**) Zoom on the study area. The positions of point W (offshore Pozzuoli's Gulf), Monitoring and Environmental Data Units, MEDA-A and MEDA-B, are also depicted. The brown contour line defines the remediation site boundaries.



### 2.2.1. Acoustic Doppler Current Profiler

The ADCPs used for the test campaign are one of the most widely used instruments in oceanographic research for measuring the wave velocity profile. Such instruments are also able to provide wave information. ADCP-A and ADCP-B are part of the aforementioned MEDA A and MEDA B, respectively. The ADCP (Figure 6) is a bottom-mounted upward-looking instrument which takes the measurements of the waves basically using three independent techniques.

**Figure 6.** Stazione Zoologica "Anton Dohrn" (SZN) instruments: (**a**) MEDA A; (**b**) MEDA B; (**c**) acoustic Doppler current profilers (ADCPs).

**(a) (b) (c)** 

The first method is wave measuring using the basic principle of Doppler shifting to evaluate the orbital velocities of waves, ensonifying the entire water column along four inclined beams. The orbital velocity measured by the ADCP along each distant beam provides information above the directional and non-directional wave spectrum. In addition to wave orbital velocity measurement, the ADCP also measures the non-directional spectra through echo ranging (surface track) and bottom pressure with a pressure transducer, providing redundant measurements of wave height and water depth.

In Tables 6 and 7, the ADCP specification and the parameters used for the spectral analysis are described. Details of the ADCP wave measurements are described in [102–106].


**Table 6.** ADCP configuration.


#### 2.2.2. Directional Wave Spectra Drifter-Derived Wave Buoy

The DWSD buoy uses the GPS sensor package in order to measure *w(t)*, *u(t)* and *v(t),* which represent respectively the vertical, horizontal E-W and horizontal S-N buoy velocity components, from changes in the frequency of the GPS signal [86,93,107,108]. The measurements are made for a 17 min long sample of *u(t)*, *v(t)* and *w(t)* every hour, divided into overlapping 4-min segments with 1 Hz of sampling frequency. The power spectral density, co-spectra and quadrature-spectra parameters are derived from a Fourier transform of the correlation functions related to each pair of the three aforementioned signals, giving the first five independent Fourier coefficients (a0, a1, a2, b1, b2) and thus the wave spectra for each hourly sea state. For each measured sea state, the three velocity components, the computed first five Fourier coefficients and the main wave data parameters are transmitted in real time through the Iridium satellite system. All these wave data, including data on battery voltage and the pressure, temperature and humidity of the hull, are accessible in real time from a dedicated website.

The DWSD buoy (Figure 7) has a simple spherical geometry with a diameter of 0.39 m and weight of 12 daN, reducing in this way the installation and maintenance costs, being very easy to handle and to install.

**Figure 7.** The directional wave spectra drifter (DWSD) buoy: (**a**) ashore, prior to launching; (**b**) right after deployment at MEDA-B. The orange float required for the mooring system is also shown.

#### *2.3. Wave Propagation and Model Calibration*

In order to consider the intricate variations in wave energy density occurring from offshore of the Gulf of Naples to the Bagnoli-Coroglio Bay, the nearshore energetic patterns have been studied by means of the numerical suite MIKE 21 SW spectral wave model, developed by DHI Water and Environment [109]. The model takes into account the effects of refraction and shoaling due to varying depths and local wind generation and energy dissipation due to bottom friction and wave-breaking, and it has been validated by comparison with data from buoys and satellites by several authors [8,110–113]. Moreover, several scientific papers (e.g., [114–117]) discussed the overall satisfactory agreement between MIKE 21 SW and SWAN, TOMAWAC and STWAVE. In particular, Ilia and O'Donnell [118] found that both MIKE 21 SW and SWAN were largely consistent in their observations during storms, even if MIKE 21 SW predicted some of the storm peaks slightly better than SWAN. Therefore, the results from this study can be of interest for applications with other spectral wave models.

The basic equations in the model are derived from the conservation equation for the spectral wave action density *Z*, based on the approach proposed by the authors of [119]. In fact, in the presence of currents, wave action is conserved whilst the wave energy is not [120]. The source/sink term that represents all physical processes which generate, dissipate or redistribute energy, *Stot*, can be written as:

$$S\_{tot} = S\_{in} + S\_{surf} + S\_{dw} + S\_{bvt} + S\_{nl} \tag{5}$$

where *Sin* represents the energy transfer from wind to waves, *Ssurf* is the dissipation of wave energy due to depth-induced breaking, *Sdw* is the dissipation of wave energy due to whitecapping, *Sbot* is the dissipation due to bottom friction and *Snl* is the energy transfer due to non-linear triad (three-wave) interactions. The following approaches/models are used in the model:


Operatively, the models compute the evolution of *Z* by solving the action balance equation [127], which in the Cartesian co-ordinates can be written as:

$$\frac{\partial Z}{\partial t} + \nabla\_{\mathbf{x}, \mathbf{y}} \left[ (\mathbf{C}\_{\mathcal{S}} + \mathcal{U}) \mathbf{Z} \right] + \frac{\partial}{\partial \sigma} (\mathbf{C}\_{\sigma} \mathbf{Z}) + \frac{\partial}{\partial \theta} (\mathbf{C}\_{\theta} \mathbf{Z}) = \frac{\mathbf{S}\_{\text{tot}}}{\sigma} \tag{6}$$

where *Z* = *V*/σ, *V* being the variance density and σ the relative angular frequency, θ is the mean wave direction measured clockwise from true north, *Cg* is the group velocity, *U* is the current velocity vector and *C*<sup>σ</sup> and *C*<sup>θ</sup> are the propagation velocities in spectral space (σ*,*θ). The left-hand side of the above equation represents the local rate of change of the wave energy density, propagation in geographical space and shifting of frequency and refraction due to the spatial variation of the depth and current.

For wave propagation over slowly varying depths *h*, σ can be written by means of the linear dispersion relation

$$
\sigma = \sqrt{\text{g } k \tanh(kl)}\tag{7}
$$

in which *k* is the wave number.

The magnitude of the group velocity *Cg* is given by

$$C\_{\mathcal{S}} = \frac{\partial \sigma}{\partial k} = \frac{1}{2} \left[ 1 + \frac{2kh}{\sinh(2kh)} \right] \sqrt{\frac{\mathcal{S}}{h} \tan h(kh)}\tag{8}$$

The implicit assumption of these equations is that properties of the medium (water depth and current) as well as the wave field itself vary on time and space scales that are much larger than the variation scales of a single wave.

The model takes into account diffraction by using the approximation proposed by Holthuijsen et al., [128], based on the revised version of the mild slope equation model of Berkhoff [129] that was proposed by Porter [130].

It is worth noting that the source functions *Sin*, *Snl* and *Sds* in MIKE 21 SW are similar to the source functions implemented in the WAM Cycle 4 model [131,132]. The latter provides the basis for the ECMWF wave hindcast dataset [16,133,134]. One of the main restrictions of the model is that when propagation leads to waves moving nearly parallel and close to the coast, there is an unrealistic loss of energy caused by the large second-order diffusion error [135]. In this case, moreover, the main assumption that the source integration time step has to be shorter than or equal to the propagation time step is at fault. Hence, an intrinsic sensibility to direction can be detected, representing a warning if significant diffraction-reflection conditions can be found.

The model solves the governing equation by means of finite element-type methods to discretize geographical and spectral space. A parameterization of the conservation equation in the frequency domain is performed by introducing the zeroth and the first moment of the action spectrum as dependent variables.

The computational domain was discretized using an unstructured grid with meshes based on linear triangular elements (Figure 8) and performed using the cell-centered finite volume method. The seabed was performed by interpolating at the grid nodes the information provided by the General Bathymetric Chart of the Oceans (GEBCO) database [136]. The grid resolution was assumed to be variable linearly between the maximum depth to 150 m for depths in the range of 500 m to 100 m. Constant values of 150 m and 1000 m of the grid resolution have been assumed for water depth shallower than 100 m and deeper than 500 m, respectively.

**Figure 8.** Zoom on the Bagnoli-Coroglio Bay with focus on the computational mesh implemented in MIKE 21 SW:

The wave model was run as forced wave-by-wave with data from the ECMWF internal WAve Model (WAM) with the ERA-Interim dataset related to source point E3. The basic data necessary to fulfill the offshore requirements are the significant wave height (*Hm0*), mean wave period (*Tm*) and mean wave direction (θ), provided by 6-h hindcast wave data. Wave power series was calculated from the resulting dataset provided by the transformation model. For natural sea states, where waves are random in height and period (and direction), the spectral parameters have to be used. The wave energy flux can be defined as:

$$P = \rho \text{g} \int\_0^{2\pi} \int\_0^\infty \mathbb{C}\_{\mathcal{S}}(f, h) \mathbb{S}(f, \theta) df d\theta \tag{9}$$

where is the sea water density, *S*(*f*,θ) denotes the 2D wave spectrum as a function of the spectral wave frequency *f* and mean wave direction θ and *Cg*(*f*,*h*) denotes the wave group velocity, expressed by Equation (8).
