Calibration, Validation, and Analysis of an Empirical Algorithm for the Retrieval of Wave Spectra from HF Radar Sea Echo

The accuracy of the wave products retrieved by a 12-MHz high-frequency (HF) phased-array radar is evaluated for a 5-month period. The two stations composing the systemwere deployed in 2011 to overlook the Wave Hub, a test site for marine renewable energy devices located on the southwestern coast of the United Kingdom. The system was conceived and configured to reduce the inaccuracies introduced by short time averaging and minimal overlap between stations, both associated with the most traditional HF radar deployments, whose primary activity is current measurement. Wave spectra were retrieved by an empirical algorithm distributed with Wellen Radars (WERA), which were calibrated using in situ measurements collected within the radar footprint. Evaluated through comparison against measurements acquired by three in situ devices, the results revealed estimates of significant wave height with nearly zero bias, linear correlations higher than 90%, and RMS errors that range from 29 to 44 cm. The relative error of wave energy period comparisons was within 10% for periods between 8 and 13 s, while both underand overestimations were observed above and below that range, respectively. The validation demonstrated that when locally calibrated, the algorithm performs better than in its original form in all metrics considered. Observed discrepancies are mainly attributable to single-site estimations, antenna sidelobes, and the effect of the secondharmonic peaks of the Doppler spectrum.


Introduction
Marine energy is expected to provide 1.3 GW of electricity to the U.K. grid by 2020 (DECC 2009), a projection that has been the stimulus for extensive research focused on developing efficient, economically feasible wave energy converters (WECs). As a result, different devices have been produced and are now being tested in offshore demonstration sites in order to assess their efficiency and survivability (Saulnier et al. 2012). With this step forward of the offshore marine renewable energy (MRE) industry, reliable wave data sources are now not only necessary for assessing the resource (Boudière et al. 2013) but also for providing continuous input throughout the operational stage of the wave energy projects (Venugopal et al. 2011b).
High-frequency (HF, 3-30 MHz) radars are shorebased remote sensing systems that can help to satisfy the ever-increasing need for wave and surface current data within the MRE sector, complementing the information obtained from satellites and point measurement in situ devices. This technology has the capability of measuring surface currents, waves, and wind direction from locations up to 200 km offshore, depending on the radar's operating frequency and geometry, and the variable being measured. Results are provided in near-real time with a spatial resolution that can range between hundreds of meters and tens of kilometers, depending on the allocated bandwidth and the antenna design. Such an ability to provide high-resolution data over a large area in near-real time could be extremely useful for MRE, which would benefit from spatial measurements when assessing the resources and monitoring the wave energy sites.
Remote sensing of the ocean by HF radar is based on the transmission of vertically polarized electromagnetic waves, which travel along and are scattered off the rough sea surface. To first order, the scatter occurs when the radio wave interacts with waves of exactly half its wavelength traveling toward and away from the radar. These waves, known as Bragg waves, are in resonance with the transmitted wave, and therefore produce coherent reflections that result in two strong peaks in the recorded backscatter spectrum (Paduan and Graber 1997), at a relative Doppler frequency determined by the phase velocity of the scattering ocean waves and the underlying current . Because the expected Doppler shift due to the phase speed of the Bragg waves is readily derived from the wave dispersion relation, any displacement from this theoretical value is attributed to the surface current, and this constitutes the basis for current measurement with HF radar. However, the recorded signal not only contains the returns from the fundamental first-order Bragg waves. Nonlinear wave interactions and double-scattering processes result in higher-order waves with a wavelength equal to onehalf the transmitted wavelength that will produce Bragg resonant scatter and contribute to the Doppler signature at different frequencies, forming the so-called secondorder continuum. Because waves of any wavelength may be involved in the second-order scatter, the continuum contains information about the entire wave directional spectrum.
The strong signal produced by the resonant first-order Bragg waves, along with the simple relation between current velocity and the Doppler shift of the first-order peaks, provides a robust and reliable HF radar measurement of the surface current, which is routinely used for oceanographic studies and operational applications (i.e., Fujii et al. 2013;Liu et al. 2014;Lorente et al. 2015). The second-order spectrum, on the other hand, is characterized by a weaker signal that can be sometimes buried in noise; thus, the inversion is in this case a further constrained problem, limited to about 50% of the range over which surface current measurements are possible. Additionally, this part of the spectrum arises from nonlinear wave interactions that make the interpretation of the radar signal an intricate process not as straightforward as the estimation of surface current radial velocities. Consequently, the inversion of the second-order part of the backscatter spectrum into meaningful wave data is an area of active research, and the validity of the radar's wave products is still under scrutiny. Current measurement is therefore the primary activity of the majority of HF radars installed worldwide, which have accordingly been optimized for this purpose (i.e., Haus et al. 2010;Voulgaris et al. 2008;Wyatt et al. 2009Wyatt et al. , 2013. Using data acquired by these systems for wave processing can have implications in the extent and quality of the estimates. For example, a substantial alongshore separation between radars, which is sometimes required to maximize coverage for current measurement, can result in a minimal overlap across the short ranges where wave measurements are possible (i.e., Voulgaris et al. 2008). This often entails the use of only one station to obtain wave estimates, which are consequently restricted to nondirectional wave products due to the intrinsic direction ambiguity associated with a technique based on radial measurements. Moreover, the integration time for measurements is normally restricted to a few minutes, in order to allow resolution of the rapidly changing conditions sometimes related to surface currents. This implies data samples that are not long enough to reduce sampling variability impacts on wave calculations (Wyatt et al. 2009), and it has been found to result in noisy estimates (Voulgaris et al. 2008).
To overcome the above-mentioned constraints, two phased-array Wellen Radars (WERA; Gurgel et al. 1999) specifically configured for wave measurement were installed in the southwest coast of the United Kingdom to overlook an offshore wave energy test site located in the area. This is the first application of this technology to monitor wave conditions for MRE installations, and it is intended to provide high-quality data that should allow for an appropriate evaluation of HF radar wave measurements. The arrangement of the two radars that compose the system is such that it maximizes their beam overlap across the ranges where wave measurements are possible. In addition, data are acquired over 17-min sampling intervals, an integration time similar to that of the in situ devices, and more appropriate for wave measurement (Wyatt et al. 2009). Results from this system are therefore expected to provide a good indication of the HF radar capabilities for wave measurement and to help prove the suitability of this technology for MRE.
The backscattered signals received at the radar described above are inverted into frequency ocean wave spectra using an empirical method (Gurgel et al. 2006; WERA algorithm) based on a simple linear relation between the ratio of second-to first-order Doppler energy and the ocean wave spectrum weighted by an angular spreading function. The method, in spite of being provided with all WERA radars, has received little attention in the literature; hence, very little is known about its performance. Furthermore, although it provides an estimation of the frequency wave spectrum, only validations of the wave height results have been reported in the literature (Gurgel et al. 2006;Savidge et al. 2011;Toro et al. 2014). Savidge et al. (2011) provide a comparison of radar-derived wave heights and in situ measurements in the form of scatterplots, but no quantitative measures for the goodness of fit are reported. Nevertheless, the figures indicate high scatter and positive bias on the radar results, somewhat poor results that were nonetheless attributed to a suboptimal configuration of the radar. Toro et al. (2014) report on the results of their own adaptation to the original WERA algorithm. Their modifications comprise the inclusion of a weighting function derived by Barrick (1977), a wind speed-dependent relation between Doppler and ocean wave frequencies, and four different parameterizations-also wind speed dependent-to scale the Doppler spectra. A linear correlation (R) of 0.73 and 0.39-m RMSE are reported for a 16-day comparison of buoy-and radar-derived wave heights, limited to cases when the wind was blowing toward the radar with a speed (U 10 ) higher than 8 m s 21 . Toro et al. (2014) also provide what is the only quantitative information related to the shape of the frequency spectrum retrieved by the WERA algorithm. This is represented by an estimation of the wind-sea spectral peak-a parameter that was found to agree very well with in situ measurements, showing a linear correlation R of 0.9 for a limited dataset corresponding to the same wind conditions referred to above.
Here, we intend to complement these rather scarce validations by providing a thorough evaluation of the results produced by the WERA software, using both the original calibration found in Gurgel et al. (2006) and a local calibration performed as part of this work. Using a longer dataset than most published validations of empirical inversion algorithms (i.e., Maresca and Georges 1980;Wyatt 1988;Heron and Heron 1998;Essen et al. 1999;Ramos et al. 2009;Toro et al. 2014), and adding surface gravity wave spectra and energy period (T e ) to the more typical significant wave height (H s ) validations, we hope to enhance the understanding of the algorithm's accuracy, strengths, and limitations. Assessing outputs extracted from three different locations within the radar's coverage is expected to help determine the spatial reliability of these measurements and contribute to the evaluation of the adequacy of HF radar for wave monitoring at wave farm installations.
This document is organized as follows. Section 2 describes the datasets used, as well as the area where they have been collected. In section 3, the empirical algorithm and the process followed to calibrate it are described. Validation of the results over a 5-month period is presented in section 4 and a discussion is given in section 5. Finally, a summary of the major conclusions drawn from the results presented hereby is provided in section 6.

a. HF radar
The dataset used in this study (D. C. Conley 2013, unpublished data) was obtained with two WERA operated by Plymouth University. Both are located on the north coast of Cornwall, United Kingdom, from where they overlook a test site for offshore renewable energy technologies. The site, known as Wave Hub, covers 8 km 2 and hosts an underwater electrical hub, located 16 km offshore at an average depth of 50 m and connected to the U.K. grid through an undersea cable. Wave Hub is intended to help demonstrate the commercial viability of wave energy by providing an offshore site for full-scale demonstration of devices. The radars were installed as part of the many projects that have been launched over the past years to support the activities at Wave Hub and assess the environmental impacts of wave energy installations.
The two radar stations that compose the system, located at Pendeen and Perranporth (Fig. 1), approximately 40 km apart, have been operational since 2011. Each site consists of a 16-element phased-array receiver and a square four-element transmitter located roughly The open circles represent grid points where dual data availability during 2012 was higher than 70%. Closed circles represent points with data availability higher than 70%, but in this case both dual-and single-radar data are considered. The in situ mooring devices are also depicted: wave buoy (star), ADCP-W (inverted triangle), and ADCP-E (triangle). The gray rectangle delimits the Wave Hub test site for MRE devices. parallel to the coast. The receiving array at Pendeen is orientated at 1138 from true north (clockwise); hence its boresight (beam direction perpendicular to the receivers) is directed 238, also from north. At Perranporth, the receiving array is 358 from north and its centroid beam is directed to 3058N, looking at the direction from which the prevailing westerly swell approaches the area. Independent measurements are synchronously acquired at the two stations for 17 min 45 s, every hour, at approximately 1-km range resolution and 78 angular resolution. Both radars operate in a ''listen before talk'' mode (Gurgel and Schlick 2007), determining the cleanest operating frequency for transmission within a 250-kHz bandwidth, and around a center frequency of 12 MHz. At this frequency, the transmitted waves are backscattered off ocean waves 12.5 m long, from distances up to 101 km. However, only surface current measurements are achievable over the full range, while wave products are limited to half this distance due to the lower signal-to-noise ratio (SNR) of the second-order returns compared to the first-order echo. Considering only grid points with data availability higher than 70%, this translates to approximately 50 km 3 50 km of spatial coverage for wave measurement (Fig. 1). This area is further reduced to approximately 19 km 3 19 km if we take into account only grid points where data from the two radars are available simultaneously more than 70% of the time. The information backscattered from the illuminated ocean area is then transformed into an orthogonal coordinate system and set into a uniform grid at 1-km spacing.

b. Wave buoy
The wave buoy used for this study is a Seawatch Mini II directional buoy, manufactured by Fugro Oceanor, deployed at 53-m depth, and 20 km offshore from Pendeen and 30 km from Perranporth. The buoy uses a combination of gyroscopes, accelerometers, and a compass to resolve its position in the north, east, and vertical directions. Measurements are acquired for 17 min 4 s, every 30 min, at 2 Hz. According to the manufacturer, these measurements are accurate to 0.1 m of the measured elevation, 2% of the wave period, and 0.38 of the measured direction. The recorded displacement time series were subsequently processed using Wave Analysis for Fatigue and Oceanography (WAFO), a toolbox of MATLAB routines for statistical analysis, and simulation of random waves and loads (Brodtkorb et al. 2000). Directional spectra were calculated using a 256-point cross-spectral analysis with 50% overlapping, and a Hanning window. The extended maximum entropy method (EMEM; Hashimoto 1997) was applied to estimate the directional distribution function over 90 directions. This process resulted in directional spectra at 0.0078-Hz frequency resolution and 48 directional resolution. These spectra were posteriorly resampled to 0.01 Hz, and wave parameters were calculated by integrating them over the range of 0.05-0.25 Hz, which is the frequency range of the radar outputs.

c. ADCP
Two upward-looking Teledyne RD WorkHorse broadband acoustic Doppler current profilers (ADCPs) were also deployed within the radar's measurement domain for validation purposes. The two ADCPs operate at 600 kHz and have a Janus configuration with four beams, each of them located at 208 from the vertical. Measurements are acquired for 20 min, every 2 hours, at 2 Hz. Both instruments were located at an average of 30-m depth over the two deployments that provide the data used here. The ADCP-E (east) was located at 24 km from Pendeen and 19 km from Perranporth, while the ADCP-W (west) was moored at 16 km from Pendeen and 29 km from Perranporth. The data collected were subsequently processed with WavesMon, the manufacturer's software. The latter calculates the directional wave spectrum over 64 frequencies at 0.016 Hz, and 90 directions, using the iterative maximum likelihood method (IMLM; Krogstad 1988). The raw data used to perform the calculations are derived from the orbital velocities measured by the sensor, which at 600 kHz have an accuracy of 60.3%. Wave parameters were calculated from the spectra over the frequency range of 0.05-0.25 Hz.

Calibration of the empirical algorithm
a. Description of the algorithm Following Hasselmann (1971), the algorithm presented in Gurgel et al. (2006) is based on the assumption that the continuous second-order sidebands surrounding the first-order peaks of the Doppler spectrum are proportional to the frequency ocean wave spectrum. This simple relationship, however, was modified by introducing an angular spreading function as follows: where a k are the coefficients resulting from regressing radar sidebands against buoy spectra, with k representing the index of the spectral frequencies (0.05-0.25 Hz). Variables S m and S p refer to negative (m of minus) and positive Doppler shift relative to the first-order Bragg peak, respectively. The ocean wave spectrum H k is measured in situ; F is an angular spreading function; f r is the radar look direction, which is defined as the angle between the radar boresight and a particular location within the radar coverage (the wave buoy in this case); and f k is the in situ-measured mean wave direction, both measured clockwise from north. Gurgel et al. (2006) analyzed Eq. (1) using different angular spreading functions, including a cos 2 and a sech 2 . The latter is reported to produce correlations of 0.75-0.9 for the positive sideband and of about 0.3-0.4 for the negative side when compared to buoy spectra. No correlation results are given for the cos 2 function, but the authors state that both resulted in a similar set of regression coefficients. The correlation of the unweighted sidebands and buoy spectra was also examined and was found to result in a 10% reduction with respect to the sech 2 correlations. The regression coefficients were obtained by applying a cos 2 function to approximate the angular distribution of the short scattering ocean waves. The solution of Eq. (1) for the cos 2 (0:5f) is H k 5 a k (S mk 1 S pk ) . (2) In the WERA software, the correction coefficients found by Gurgel et al. (2006) are applied to Eq.
(2) in order to calculate the ocean wave spectrum.

b. Calibration of the algorithm
In their paper, Gurgel et al. (2006) state that the algorithm and correction coefficients derived using radars operating at 27.65 MHz can be used for other radar frequencies if multiplied by the squared ratio of the latter frequency and the new radar operating frequency. This correction yields a new set of coefficients that provides reasonable results for the 12-MHz Plymouth University radar, mainly in terms of significant wave height. However, rather less satisfactory estimates of wave period motivated a recalibration using data collected locally, in an attempt to further improve the results.
The radar second-order sidebands appearing in Eq. (2), and necessary to calibrate the algorithm, were obtained using the WERA software. In the first stage of the processing performed by the WERA software, the received returns are range and azimuth sorted to grid points across the radar coverage area (see Gurgel et al. 1999 for details), and a Doppler spectrum is generated via FFT at each grid point. In the next stage, the firstorder peaks are identified for each spectrum, and the first-and second-order regions, used to estimate current and wave information, delineated. The stronger of the two Bragg peaks identified is then used to normalize the second-order sidebands surrounding it. Thereafter, Doppler frequencies of the normalized secondorder sidebands are linearly transformed into ocean frequencies by subtracting the Bragg frequency, and the associated spectral amplitudes are interpolated onto a vector ranging from 0.05 to 0.25 Hz, at a resolution of 0.01 Hz. These results are stored and provided as outputs. Using these variables as a starting point, a quality control procedure equivalent to that implemented in the WERA software was performed. As part of the process, data points with a signal-to-noise ratio lower than 25 dB in the first-order region, and lower than 10 dB in the second-order sidebands, were flagged as erroneous (see Fig. 2a). Furthermore, spectra with an abnormally high signal in the second-order bands were also flagged, as it is likely to find interference with the first-order region in this situation. Flags for these tests were raised for 21% and 11% of the training datasets recorded at Perranporth and Pendeen, respectively. Once all the spectra on the measuring grid had been assigned a quality flag, the result at each node was calculated as a nine-point spatial average, provided all neighbors satisfy the criteria referred to above. The end products were the S m and S p sidebands appearing in Eq. (1), which are averaged in space (over 9 km 2 approximately) and have a quality flag assigned to them.
As mentioned in the previous section, Gurgel et al. (2006) found a better correlation between Doppler sidebands and in situ-measured spectra when the latter were weighted by an angular spreading function, and therefore applied a cos 2 function prior to finding the empirical parameters that allow estimation of the ocean wave spectra from the Doppler sidebands. Here, we carried out a similar analysis, and studied the correlation between the Doppler sidebands and the buoy-measured spectra with and without the angular distribution weighting. As opposed to Gurgel et al. (2006), our results did not show a significant improvement in terms of correlation after applying an angular spreading function; hence, the unweighted spectra were used to obtain the correction coefficients. Linear correlations (R) above 0.7 were found between radar returns and ocean wave spectra at frequencies between 0.06 and 0.1 Hz. These numbers decreased toward the high frequencies to reach values around 0.4.
The statistically significant correlations suggest that there is indeed proportionality between the energy on the radar sidebands and that on the measured wave spectra, mainly at the low frequencies. Therefore, applying a correction such as that suggested by Gurgel et al. (2006) could be a good approach for obtaining the ocean wave spectrum from HF radar measurements, especially if a persistent swell propagates over the measurement area. Data collected over three monthsnamely, March, August, and September 2012-were selected as the inputs for the training dataset. A varied group of sea states, characterized by wave heights ranging from less than 0.5 m to more than 5 m, peak periods between 4 and 20 s, and wind speeds up to 18 m s 21 blowing from various directions (Fig. 3), suggests that the chosen period can provide a representative dataset for developing a relevant set of correction coefficients that can then be used to calculate the onedimensional ocean spectrum from the radar signal.
The fitting procedure was formulated as a least squares minimization problem between wave spectra measured in situ and the sum of the energy in the sidebands surrounding the most energetic side of the Doppler spectra [Eq. (2)]. The a k coefficients of Eq. (2) were found by minimizing the sum of the square of the errors between radar and in situ spectra. Two independent regression analyses were performed, using

250
data collected by each of the two radar stations comprising the system.

a. Result of the calibration
The regression analysis detailed in the previous section was conducted using data collected at three different points within the radar's field of view, where in situ measurements were available (Fig. 1). The results, depicted in Fig. 4, show that similar correction factors need to be applied at different ranges and azimuthal directions within the radar's coverage, suggesting that the Doppler spectra do not differ significantly from each other at those locations. At the low frequencies, the main difference among the three curves derived for Perranporth lies in the higher values obtained at the buoy node. This might be the result of the attenuation of the radar signal with distance, as this location is the farthest from the radar station. At Pendeen, the curves obtained for the ADCP-W and wave buoy are similar in shape, but the former shows lower values. The mooring is the closest to the radar station, and the angle between the prevailing swell and the radar beam is lower than at the buoy's position. The combination of these two factors results in a stronger signal backscattered from this location, hence the lower correction coefficients. Around the frequencies 0.14-0.15 Hz, a minimum appears as a common feature to all curves. Assuming a linear relation, such frequencies equate to a Doppler frequency of 62 1/2 f B (the frequency fluctuation would be due to the variation in transmitting frequency mentioned in section 2a), where the second harmonic peaks (SHP; Fig. 2b) appear in the Doppler spectrum (Ivonin et al. 2006). Resulting from the resonant scatter produced by the second harmonic of the Bragg waves, these peaks do not represent the actual shape of the ocean wave spectrum that generated the scatter (Barrick 1977), and the minima appear to compensate for this frequency localized problem. Finally, the higher variability of wind seas and the noisiness of the radar data at high frequencies, both represented by the broad confidence intervals shown in Fig. 4, are likely to be responsible for any differences between the correction coefficients above 0.17 Hz. Overall, the most noticeable discrepancy can be seen in the coefficients of the two radar stations, especially over the lower frequencies, where the values obtained using data collected by the Pendeen radar are approximately 60% higher. The prevailing westerly swell propagating across the area intersects Pendeen's beam almost perpendicularly at most azimuthal directions; a situation known to result in a lower signal as compared to that emerging from waves propagating toward or away from the radar (Wyatt 2002), the latter being the common situation at Perranporth.
Despite the differences found, which are primarily due to the known dependency of the radar spectral signature on the angle between wave direction and radio beam direction (Lipa and Barrick 1986), the similarity between the three sets of coefficients obtained for each station suggests that a constant curve can be a valid approximation to inverting the radar backscattered power over the measurement grid. The coefficients obtained using data acquired with the wave buoy were therefore chosen as the reference set to use for future calculations. This decision was made based on data quality, as the wave buoy was located at low angles from the boresight line of both radars, where the digital beamforming performs best.

b. Assessment of the calibration
The coefficients obtained following the process described above were used to calculate the one-dimensional ocean spectrum. The calculations were performed using the average of the normalized second-order power of the two stations and their respective mean coefficients when possible, while the station-specific coefficients were applied when dual data were not available. Wave parameters were calculated by integrating the latter spectra over the full range of available frequencies (0.05-0.25 Hz). Significant wave height (H s ) was calculated as 4 times the square root of the zeroth moment of the FIG. 4. Coefficients obtained regressing radar Doppler spectra to in situ-measured spectra recorded at Pendeen (thin lines) and Perranporth (thick lines). Three curves, corresponding to the mooring devices deployed in the area, are shown for each station: ADCP-E (dashed black line), ADCP-W (blue dotted line), and wave buoy (gray solid line). The coefficients found in Gurgel et al. (2006), used on the WERA software, are also shown (yellow dashed-dotted line). The error bars are the 95% confidence intervals of the coefficient estimation. spectrum (Venugopal et al. 2011a), and the energy period (T e ), a relevant parameter for wave energy applications, was obtained by dividing the first negative moment by the zeroth moment of the spectrum (Venugopal et al. 2011a). The results retrieved by the calibrated algorithm were then compared to in situ measurements in order to assess their accuracy. The agreement between methods is quantified in terms of the linear correlation coefficient (R), the bias of the estimations, the root-mean-square error (RMSE), the scatter index (SI), and the mean percentage error (MPE). The bias was obtained from the difference between the means of both methods. The SI was defined as the RMSE normalized with the averaged in situ value and the MPE was calculated as the percentage difference between the radar and in situ values divided by the magnitude of the latter, and normalized by the number of records N. The latter was computed for different wave height and period classes, in order to verify the performance of the calibration under different sea states. In calculating the above-mentioned statistical descriptors, the in situ measurements were assumed to be sea truth, and any differences between methods are accounted for as error in the HF radar estimate. Such an assumption, however, is not strictly true, as differences can arise from the distinct measurement principles, sampling variability, or systematic offsets (Krogstad et al. 1999). For our particular case, a potential source of variability would be that emerging from the different sampling volumes of each device, which range from approximately 1 km 2 in the case of the HF radar to some 100 m 2 in the case of the ADCPs, and down to the few meters of its excursion radius, in the case of the wave buoy.

1) RADAR-IN SITU COMPARISONS
Comparisons of the radar estimates against data collected by the three-point measuring devices located within its footprint are presented in Fig. 5 and summarized in Tables 1 and 2, where the results obtained with the original coefficients found in Gurgel et al. (2006) are also included. The calibrated algorithm retrieved estimations of H s with nearly zero bias, linear correlations higher than 90%, an average scatter index of 0.2, and RMSEs that range from 29 cm at the wave buoy coordinates to 44 cm at the ADCP-W. The higher errors found for the ADCP-W, which was deployed at a high angle from Perranporth's boresight (458), are probably associated with antenna sidelobes. The error is some 14% higher in the results retrieved by the original WERA algorithm, which nonetheless performed similarly to the calibrated scheme for this parameter. The improvement of the calibration relative to the unmodified algorithm is most noticeable in the T e estimations, which revealed up to a 30% error reduction.  Figure 6a shows the mean percentage errors for wave height bins of 0.4 m. While the error on the radar's significant wave height estimations is within 10% for most wave amplitudes, deviations occur in the range of 2.8-4-m waves for comparisons against the ADCP-E measurements, which present errors closer to 20% with respect to the in situ measurements. The overestimations in this range are associated with low-frequency waves (T e . 9 s) and wind speeds higher than 10 m s 21 , blowing from the west (Fig. 5c). The Doppler spectra recorded at both stations under these conditions (not shown) have their energy well distributed over the frequencies. However, they are peaked at the maximum of the correction coefficients and show some extra energy around 0.14 Hz. The combination of these two factors seems therefore to be responsible for the overestimations. On the other hand, the underestimation of wave heights higher than 4 m shown in Fig. 6a is also associated with high wind speeds (.14 m s 21 ) (not shown), but in this case, the energy-containing waves are traveling from the northnorthwest, aligned with the wind (see Fig. 5c). Under these conditions, the Doppler spectra recorded at Pendeen (not shown) are dominated by the SHP and do not represent the ocean gravity spectra that produced the scatter. At Perranporth, the radar-to-wind angle reaches the values for which the hydrodynamic coupling coefficient is low (Lipa and Barrick 1980), resulting in a relatively weak backscatter. As a consequence, H s is underestimated. An increase in the error can also be observed in Fig. 6a below 1 m and until the 0.4-m minimum wave height detectable at the radar's transmitted frequency of 12 MHz. The figure also reveals the benefit of using only dual data to calculate wave height-a situation resulting in a 20% error reduction in the radar's estimations with respect to the buoy's measurements at the upper tail of the wave height distribution (6-6.4 m). An exception to the latter can be observed on the low wave heights obtained at the ADCP-W, where single-site estimations show 30% less error than dual estimations. The reason for this is that, as it can be seen on Fig. 5e, such single-site estimations are usually derived from Pendeen's measurements which, as opposed to Perranporth's, are unaffected by antenna sidelobes at this location. The same analysis was conducted for the T e results, which were discretized into 1-s bins. The results show accurate estimations of wave periods between 8 and 13 s, while both under-and overestimations can be observed above and below that range, respectively (Fig. 6b). In addition, and as opposed to wave height estimations, this parameter does not show any significant improvement by using only dual-radar estimations.
Compared to the original WERA algorithm, the highest impact of the calibration was in the results calculated using data collected exclusively at Perranporth, which represent 28% and 6% of the ADCP-E and buoy validation datasets, respectively. For those cases, the results derived from the calibrated algorithm show a 26% and 31% average RMS error reduction in H s and T e , respectively. When using dual data, the estimations obtained with the new coefficients show a 15% reduction in H s RMSE and a 13% reduction in T e . On the other hand, the results obtained from the signal received at Pendeen only did not show a significant improvement with the new calibration.

2) SPATIAL COMPARISONS
Radar wave height maps normally show a high spatial variability over the instrument measurement region. Whether this is real or a product of radar artifacts is a question that remains partially unanswered, and although the results presented in the previous section have shown that noisy results near the edge of the radar's field of view may be in fact a consequence of antenna sidelobes, a complete spatial assessment of the results is essential. For that, in situ measurements acquired at two distant locations were compared to determine the variability captured by such devices. The same exercise was then conducted with radar estimates at equivalent positions, to allow for comparison against the in situ measurements. Unfortunately, whenever the ADCP-W is involved in the comparison, the results cannot be considered trustworthy due to the noise introduced by the antenna sidelobe issue mentioned earlier, so here we will focus the discussion on the buoy-ADCP-E pair, deployed 12 km apart. The in situ devices exhibited highly correlated measurements, which show almost no bias between the wave buoy and ADCP-W in either H s (Table 3) or T e (Table  4). Although still highly correlated, the measured conditions at the ADCP-E are smaller, both in terms of H s and T e . Such modification of the wave field is due to its interaction with the varying seabed and possibly the tidal flow, and might be in fact responsible for some of the error on the radar estimates at the ADCP-E location. The highest waves (.3 m) measured by the buoy were seen to suffer an amplitude reduction when measured by the ADCP-E, which is represented by a 0.2-m bias (Table 3). Although this is also the case in the radar-to-radar comparison, the difference is smaller (0.07 m bias). In general, the radar-to-radar comparisons of the buoy-ADCP-E pair revealed 6% more scatter (SI) and 12% less bias on the H s results than the in situ comparisons (Table 3). The T e estimations derived from the radar are 6% less correlated between the two locations (see Table 4) than the in situ pair, and unlike the latter they do not show smaller T e measured at the ADCP-E.

3) SPECTRAL COMPARISONS
This section addresses the accuracy of the spectral estimates obtained with the radar. Average radar and in situ spectra were computed for the entire validation period, and for two subsets corresponding to low-and high-frequency wave components. To obtain the latter, the spectra were crudely divided into two frequency bands, corresponding to swell and wind seas, establishing the frequency separation at 0.12 Hz. The results, depicted in Fig. 7, show that the calibrated algorithm results in an accurate estimation of the average spectral shape and amplitude of the entire validation dataset (Fig. 7a), as well as the swell cases (Fig. 7b), only with a slight overestimation of 0.6 m 2 Hz 21 at the peak of both spectra. A different picture is observed for the frequency band higher than 0.12 Hz (Fig. 7c), where the peak frequency is underestimated by 0.02 Hz, explaining the overestimation of the lowest periods seen in the previous section. The accuracy of the estimates is further examined in the second row of the same figure, which depicts individual realizations of three different sea states. Figure 7d is an example of a single-site estimate derived from Perranporth's signal and shows the improvement achieved using station-specific coefficients for the inversion in that situation (yellow spectrum in Fig. 7d). As in the average results, the lowfrequency component (Fig. 7d) is well represented by the radar, which reproduces both the shape and amplitude of the spectrum quite accurately. This translates into an H s estimate with only a 10-cm error and an agreement on the T e obtained with both techniques (11 s). Wind information obtained from a wave model (Boudière et al. 2013) indicates that a 8.2 m s 21 wind was blowing from the northeast direction, and the peak wave direction measured by the wave buoy is eastnortheast. As already seen in the average spectrum, the wind sea example (Fig. 7e) reveals a spectrum with its peak shifted to a lower frequency with respect to the buoy. However, the energy contained in both spectra is the same, and the significant wave height is accurately estimated by the radar. Finally, Fig. 7f shows a mixed sea represented by a double-peaked spectrum. Both wind and peak waves are from the west, with the former blowing with a speed of 11 m s 21 . The spectrum can be seen as a combination of the two previous cases, with the swell component accurately represented both in frequency and amplitude, and the high-frequency peak shifted to a lower frequency, which coincides with the SHP.

Discussion
The results presented in the previous section revealed deviations in the results, which were mainly related to single-site estimations, sidelobe contamination, or the occurrence of wind seas. In this section, these and other factors potentially affecting the performance of the algorithm are discussed. First of all, and bearing in mind that the Doppler signature depends on factors such as the radar-to-wind direction (Lipa and Barrick 1986), the distance from the radar station (Ramos et al. 2009), and the bathymetry (Lipa et al. 2008)-all of which vary across the radar footprint-a potential source of error arises from calibrating the algorithm exclusively using measurements collected at its center. To determine whether this introduces a significant error to the results, these were recalculated using the regression coefficients derived from the ADCP-E measurements. The results obtained were then compared against the validation dataset measured by that same device (ADCP-E), which revealed a 3% reduction in the H s error and 1% in T e , with respect to the estimates obtained with the correction coefficients derived from the buoy measurements. This suggests that in spite of the differences in azimuth (308 variation in Pendeen's look direction), range (10-km difference from Perranporth's station), and depth (;20 m) between these two locations, using a constant set of regression coefficients does not significantly increase the error in the estimates.
Another factor found to affect the quality of the results was the use of a single radar to perform the calculations. Among the single-site estimations, those derived from the Perranporth radar are often closer to the in situ measurements. Furthermore, the results of this station show differing quality depending on the direction of both wind and waves. When the two are aligned, traveling toward the radar and roughly parallel to its beam, the results provide a better fit to the in situ measurements, and the single-site estimations of T e are strongly correlated to the ADCP-E and buoy measurements (R 5 0.9). Violante-Carvalho et al. (2004) found a stronger nonlinear interaction between short scattering waves and swell in this situation. This may account for the better results, since such nonlinear interactions have been found to be of importance in defining the swell peaks in the Doppler spectrum (Shen et al. 2014); thus, conditions enhancing their occurrence are likely to increase the quality of the results, as observed.
The less satisfactory results obtained at Pendeen are related to the perpendicularity between this station's beam and the prevalent westerly swell waves throughout most of the radar coverage. As opposed to what was described in the previous paragraph, when the angle between the waves and radar beam approaches 908, the nonlinear interactions between short scattering waves and swell are weak (Lipa and Barrick 1980). As a consequence, the energy of the low-frequency part of the spectrum is underestimated. An attempt to improve the results was made by finding different scaling curves [alpha coefficients in Eq. (2)] depending on wind direction, a relevant factor in defining the Doppler spectrum, already reported to have an important influence on the results (Haus et al. 2010). Using a specific coefficient set for cases when the wind was blowing from the western sector (2258-3158N), which is prone to produce the highest underestimations, reduced the RMS error of the radar energy period estimates with respect to the wave buoy by 60%, while wave height was 10% closer to the in situ measurements. However, the shape of the spectra obtained after applying the new coefficients closely matches that of the new scaling curve, which is imposing a certain shape on the low-frequency part of the spectrum, where there was barely any energy. The effect of this correction therefore warrants further research and discussion, both of which are beyond the scope of this paper.
In addition to the above, there are spatial variations in the quality of the results, which was seen to degrade from the buoy to the ADCP-W, both in the original estimations retrieved by the WERA algorithm and the hereby proposed calibration. Because the ADCP-W was deployed at a high angle from Perranporth's boresight, the radar signal backscattered from that point is likely to be affected by antenna sidelobe contamination. Haus et al. (2010) have addressed similar geometries, and they found that at angles higher than 458 sidelobe returns may be interpreted as second-order peaks and produce spurious wave estimates. Restricting the outputs to grid nodes located at angles lower than 458 from the boresight of both stations seems therefore essential to avoid noisy estimates.
It has been shown that the WERA algorithm has a limit on the accuracy of period estimations at both ends of the frequency spectrum, either with its original settings or after calibration. The underestimation of the longest swells seems to be in part related to the use of a prescribed set of coefficients to scale the second-order frequency bands. This may impose a certain shape onto the spectra, which are preferentially peaked at the maximum of the correction factors, around 0.09-0.1 Hz. Additionally, the peaks generated by these very long waves tend to merge with the first-order region of the spectrum, ultimately causing the underestimation.
On the other hand, the overestimation of the lowest periods seems to be related to a general inefficacy of the WERA algorithm in describing wind seas. The algorithm assumes a linear dependency between Doppler and ocean wave frequencies, a simple relation that, even though it does not agree with theory (Forget et al. 1981), can be a good approximation for long ocean waves. Characterized by small wavenumbers, the relation describing their interaction with a higher-frequency component to give a Bragg scatterer with wave vector 22k 0 can be approximated as (Shen et al. 2013) where K s and K w are the wave vectors of the two interacting waves (swell and wind, respectively), and K B is the Bragg wave vector. The swell frequency can therefore be approximated as where v d is the Doppler frequency, and m 1 and m 2 represent the sign of the interaction. The above-mentioned approximation is not valid at higher frequencies; thus, the Doppler and ocean wave frequencies will not be equivalent. This may therefore be partially responsible for the inaccuracies found in the high-frequency band of the spectra calculated by applying the abovementioned linear relation between Doppler and ocean frequencies.
When wind seas do occur, there seems to be a tendency to produce spectra peaked at 0.14-0.15 Hz (i.e., Fig. 7e). As explained in section 4a, these peaks are related to the resonant scatter produced by the second harmonic of the 12.5-m Bragg ocean waves traveling toward and away from the radar station. This is illustrated in Fig. 8, where it can be seen that winds blowing from sectors centered at the north and south generate a peak at 0.14 Hz on the positive and negative Doppler sidebands measured at Pendeen. This situation is replicated at Perranporth but for winds blowing from the west and the east. The occurrence of these peaks is well described in the HF radar theory (Barrick 1977), and their effect is usually corrected for by a weighting function also developed as part of Barrick's (1977) work. In the approach used here, the SHPs are partially compensated for by a minimum value in the correction coefficients (Fig. 4). However, these do not retain the Doppler frequency dependence, while the weighting function does. Additionally, multiplication of the Doppler spectra by the correction coefficients is performed after the secondorder energy has been interpolated onto the vector of ocean frequencies obtained by subtracting the Bragg frequency, and that might be spreading the SHP energy onto adjacent frequencies. Toro et al. (2014) combined the inversion approach of Gurgel et al. (2006) with Barrick's weighting function. This modification, along with a wind speed-dependent transformation between Doppler and ocean frequencies, might have been partly responsible for the better representation of the high-frequency part of the spectrum obtained by Toro et al. (2014), which is exemplified by a correlation of 0.9 between the radar and in situ wind-sea spectral peak. Nevertheless, because of the limited range of wind conditions addressed by Toro et al. (2014), further research is needed to determine the suitability of their modifications for a broader spectrum of wind and wave conditions. Such an analysis, however, is beyond the scope of this publication. It is worth mentioning that in spite of the spectral shape discrepancy between the spectra measured in situ and with radar, accurate wave height estimations are routinely retrieved under the conditions shown in Fig. 8, when the spectral energy is mostly concentrated around 0.14-0.15 Hz. This seems to be in accordance with Zhou and Wen (2014), who found strong correlations between the ratio of the SHP to the first-order Bragg peaks, and the significant wave height. They fitted a linear regression model between in situ-measured wave height and the normalized SHP energy, and used it to retrieve estimates that presented a 0.34-and 0.26-m standard deviation for single-and dual-radar estimations, respectively. Here, using Doppler spectra that had their maximum and most of their energy at either 0.14 or 0.15 Hz, we found RMS errors of 0.25 and 0.55 m for the buoy and ADCP-E datasets, respectively.

Conclusions
An empirical algorithm used to retrieve wave spectra from HF radar returns has been presented and calibrated for the area where the radar is deployed. By conducting the longest quantitative validation of one of such empirical algorithms, and adding wave spectra and period to the more typical wave height comparisons, this work has provided a deeper insight into the algorithm's skill and the source of any inaccuracies. The validation of the results demonstrated that, locally calibrated, the algorithm performs better than in its original form in all metrics considered here. Additionally, the results show that the error associated with using a uniform set of calibration coefficients for the entire radar coverage is less than 9%, at the locations examined.
The improvement obtained with the new coefficients is most noticeable on the T e estimations, and more specifically on the results obtained using data only from Perranporth's station, which see their error notably reduced by 25% in the case of wave height, and by 30% for the wave period. Attempts to improve Pendeen's singlesite estimations have also been undertaken by using correction factors fitted with data collected by that radar only, as well as introducing a specific set for westerly blowing winds, associated with the highest underestimations of H s . However, the weak signal recorded and an almost invariable shape of the Doppler spectra limit the retrieval of reliable wave data from this station. The aforementioned invariable shape of the Doppler spectrum seems to be associated with a major effect of the second harmonic peaks, which under certain conditions can mask the rest of the Doppler spectrum structure. The inclusion of a weighting function, such as that proposed by Barrick (1977), might therefore be necessary to correct for this frequency localized problem.
In an attempt to determine the spatial reliability of the radar estimations, in situ measurements acquired at two locations 12 km apart were compared against each other, in order to determine the variability captured by such devices. The same exercise was then conducted with radar estimates at equivalent positions. The in situ devices revealed highly correlated measurements between the two points, and a reduction of the highest wave heights and periods from the deeper to the shallower location. Although this is also the case with the radar-to-radar H s comparison, the bias is 12% lower. The T e radar estimations are only 6% less correlated than the in situ pair, but the period reduction between points is not evident. Despite the differences mentioned, the point measurement devices and the radar seem to capture the same variability between two distant points, and it seems reasonable to assume that most differences can be attributed to the differing measuring technique (Krogstad et al. 1999).
The results presented hereby point to the setting of the Plymouth University radar as adequate for wave FEBRUARY 2016 L O P E Z E T A L . measurement in this area, especially for dual estimations. Single-site retrieval of wave spectra, on the other hand, has proven to be more challenging, particularly at the station looking roughly perpendicularly to the direction of travel of prevailing swell waves. However, the inclusion of its results was beneficial for dual functioning, as it diluted directional effects on the results. The main limitation of the algorithm has been found to be related to the estimation of both very long period waves (.13 s) and wind seas (,6 s), which are under-and overestimated, respectively. Significant wave height, on the other hand, has been proven to be a very robust measurement, especially when the returns from the two stations are combined, which resulted in correlations above 0.9 at the three locations examined, including one grid point affected by antenna sidelobe contamination.