Improving Data Quality for the Australian High Frequency Ocean Radar Network through Real-Time and Delayed-Mode Quality-Control Procedures

Quality-control procedures and their impact on data quality are described for the High-Frequency Ocean Radar (HFR) network in Australia, in particular for the commercial phased-array (WERA) HFR type. Threshold-based quality-control procedures were used to obtain radial velocity and signal-to-noise ratio (SNR), however, values were set through quantitative analyses with independent measurements available within the HFR coverage, when available, or from long-term data statistics. An artifact removal procedure was also applied to the spatial distribution of SNR for the first-order Bragg peaks, under the assumption the SNR is a valid proxy for radial velocity quality and that SNR decays with range from the receiver. The proposed iterative procedure was specially designed to remove anomalous observations associated with strong SNR peaks caused by the 50 Hz sources. The procedure iteratively fits a polynomial along the radial beam (1-D case) or a surface (2-D case) to the SNR associated with the radial velocity. Observations that exceed a detection threshold were then identified and flagged. After removing suspect data, new iterations were run with updated detection thresholds until no additional spikes were found or a maximum number of iterations was reached.


Introduction
The Australian Ocean Radar facility at the University of Western Australia is part of the Integrated Marine Observing System (IMOS), a national collaborative research infrastructure tasked with collection and dissemination of ocean data.The radar facility measures met-ocean parameters from an array of High-Frequency Radar (HFR) systems deployed around the coastline of Australia.Commercial direction-finding (SeaSonde) and phased-array (WERA) HFR systems, provided respectively by Codar Ocean Sensors (COS) and Helzel MessTechnik, are in operation in the Rottnest Shelf and Turquoise Coast (Western Australia, WA; Figure 1), Coffs Harbour and Newcastle (New South Wales, NSW, Australia), and the South Australian Gulf regions (SAG).Each HFR node is configured primarily to sample ocean currents with a maximum range of over 200 km; however, at selected locations where the phased-arrays systems are installed, waves and winds can also be obtained.Radar data, freely available from the IMOS portal (https://portal.aodn.org.au/), are used for scientific research, operational modeling, coastal monitoring, fisheries, and other applications [1][2][3][4][5][6][7].
This paper focuses on the quality-control (QC) procedures that have been implemented for the phased array systems in Australia and are applied operationally on both the near real-time (NRT) and delayed-mode (DM) products.The main focus is on radial currents and surface current maps, Remote Sens. 2018, 10, 1476 2 of 15 but some of the proposed methodologies have the potential to be applied to wind and wave data available from the phased-array systems.
Routinely, QC tests for HFR data are applied at different processing levels, including Doppler spectra, radial velocity maps, and surface current vectors [8].Most of the tests on Doppler spectra (e.g., noise floor and computation, first-order Bragg peak detection, and the detection and removal of radio-frequency interference and ship echo [9,10]) are embedded within the manufacturer software.Other tests (for example, on radial velocity speed or signal-to-noise ratio (SNR)) rely on the definition of threshold values which should be set on the basis of historical knowledge or statistics derived from more recently acquired data [8].Examples of operational procedures applied in RT mode for HFR in general and WERA systems in particular can be found in [11,12], while advanced methodologies capable of handling complex cases of external radio frequency interferences are described in [13].
The paper is organized as follows.Section 2 describes the network set up and some of the present limitations.Section 3 provides details on the choice of the data quality proxy and the thresholds, and documents the accuracy levels of the HFR systems deployed across the country.Additional details on the NRT and DM QC procedures are also provided in this section.Section 3 focuses mostly on the HFR systems deployed in WA, due to the large number of mooring deployments in the region that can be used for validation of the HFR data.However, results and methodologies are applied operationally to the HFR systems across Australia and can be easily extended to the global HFR network.Discussions and conclusions are presented in Sections 4 and 5, respectively.

Materials and Methods
The Australian Ocean Radar network includes the compact cross-loop systems direction-finding HFR system provided by Codar Ocean Sensors (US) and the WERA phased-array HFR systems manufactured by Helzel MessTechnik (Germany).Operational settings (frequency, bandwidth, integration times, and output rates) are detailed in Table 1.Since 2014, all systems have operated within the ITU frequency bands in compliance with the 2012 World Radiocommunication Conference, ITU Resolution 612, which designated specific bands between 3 and 50 MHz to support oceanographic radar operations.The WERA systems operate at 9.335 MHz with 33.4 kHz bandwidth in the Rottnest shelf (WA) and the South Australia Gulfs region (SAG), and at 13.5 MHz with 100 kHz bandwidth at Coffs Harbor (COF, NSW).Within these frequency bands, a general increase in the average noise levels was documented [11] and a more pronounced diurnal variability was observed with respect to the frequency bands previously utilized by the Australian HF radar network.In general, operating range decreased after the frequency change due to stronger radio frequency interference.At each WERA node, the two systems operate in alternate mode with short 5 min acquisition cycles to avoid mutual interference and for monitoring purposes.The receive array is composed of equally spaced 16 element linear arrays, and the standard 4 element configuration is used for the transmit array.When operated in real-time mode, radial current maps are transmitted over a 3 G connection to a central server for QC and upload to the data portal after conversion to netCDF format.Every hour, radial maps from the two stations are averaged and converted to NRT-QC vector maps on a regular grid with 4 km (ROT and SAG) and 1.5 km (COF) horizontal resolutions.Vector maps in netCDF format are also uploaded every hour to the data portal and made available in near real time.Raw data are stored on site and retrieved on a quarterly basis (every three months) for offline processing.Proprietary software provided by the manufacturer is used at site and in offline mode to process the raw data; real-time and offline conversion, quality control, and vector mapping are performed through processing tools developed in Python and MATLAB languages.
Quality-controlled subsurface mooring data for the Rottnest Shelf region (Figure 1) are used for validation and fine-tuning purposes, primarily for the large number of sustained moorings and the wide range of oceanographic processes in the area [5][6][7].The data set used here spans a time period between 01 February 2015, 00:00 UTC, and 30 June 2015, 23:59 UTC.Subsurface current data, collected by the IMOS Australian National Mooring Network (ANMN) Facility, are available for download through the data portal [14][15][16][17][18][19][20][21][22], and include Nortek 190 KHz Continental current profilers at WATR20 and WACA20 moorings, a 600 KHz RDI Workhorse ADCP at WATR04 mooring, a RDI Workhorse Sentinel ADCP at NRSROT, a RDI Long Ranger ADCP at WATR50, and a Nortek Aquadopp Pro at WATR10 mooring sites.Details on the mooring design and the QC processing can be found in [23].
equally spaced 16 element linear arrays, and the standard 4 element configuration is used for the transmit array.When operated in real-time mode, radial current maps are transmitted over a 3 G connection to a central server for QC and upload to the data portal after conversion to netCDF format.Every hour, radial maps from the two stations are averaged and converted to NRT-QC vector maps on a regular grid with 4 km (ROT and SAG) and 1.5 km (COF) horizontal resolutions.Vector maps in netCDF format are also uploaded every hour to the data portal and made available in near real time.Raw data are stored on site and retrieved on a quarterly basis (every three months) for offline processing.Proprietary software provided by the manufacturer is used at site and in offline mode to process the raw data; real-time and offline conversion, quality control, and vector mapping are performed through processing tools developed in Python and MATLAB languages.
Quality-controlled subsurface mooring data for the Rottnest Shelf region (Figure 1) are used for validation and fine-tuning purposes, primarily for the large number of sustained moorings and the wide range of oceanographic processes in the area [5][6][7].The data set used here spans a time period between 01 February 2015, 00:00 UTC, and 30 June 2015, 23:59 UTC.Subsurface current data, collected by the IMOS Australian National Mooring Network (ANMN) Facility, are available for download through the data portal [14][15][16][17][18][19][20][21][22], and include Nortek 190 KHz Continental current profilers at WATR20 and WACA20 moorings, a 600 KHz RDI Workhorse ADCP at WATR04 mooring, a RDI Workhorse Sentinel ADCP at NRSROT, a RDI Long Ranger ADCP at WATR50, and a Nortek Aquadopp Pro at WATR10 mooring sites.Details on the mooring design and the QC processing can be found in [23].The QC metrics used here involve the calculation of correlation coefficient (R) and the root-mean-square difference (RMSD) between subsurface currents and HFR radial velocity time series at the current meter bin closest to surface and free from surface interference, and the radar radial velocity at the grid cell closest to the mooring locations.The following equations are used for R, RMSD, and the standard deviation σ: Subsurface current data are projected onto the direction between the mooring location and the radar receiver before any comparisons.A synthetic SNR field that represent the typically observed effects of the 50 Hz contamination is used to test the capabilities of the iterative identification algorithm in a controlled environment and to test the effects of different background noise levels on the detection rate.

Velocity Threshold Optimization
Data QC on radial velocity is primarily based on the application of radial velocity thresholds.The conventional approach involves setting threshold velocity values a priori, using either long-term statistics of HFR data at each radar grid cell or quantitative analyses and comparisons with near-surface current data.Data can be derived from moorings, drifters, or even glider observations, provided that they are available within the radar footprint and that differences and limitations of each individual sampling strategy are properly accounted for.In the specific deployments across Australia, both approaches have proven to be reliable and complementary: under the assumption of normally distributed radial velocity data, 99% confidence levels (3 times the deviation from the mean, σ) can be set from the velocity distributions and anomalous radial velocity data can be identified and either flagged or removed.Examples of radial current distributions from the radar stations and the current meter data for the FRE HFR system are reported in Figure 2. Regardless of the relative distances between near-surface radar data and mooring observations, which are 13-56 m below the surface, there is a surprisingly good agreement in terms of mean radial velocity and 99% confidence levels across all moorings.Before applying any velocity threshold, R values are in the range R = 0.15 to R = 0.86, with RMSD in the range 10.3 cm/s to 40.8 cm/s (Table 2).After applying a 1.5 m/s (absolute value) on radial velocity, comparison metrics improved in general for both FRE and GUI HFR stations, with R increasing up to 25% or more and RMSD decreasing by up to 50% (Table 2) Comparison metrics are in good agreement with previously reported findings in different ocean regions [24][25][26].Significant differences and poor statistical comparison results can be observed for the radar-mooring pairs FRE-WATR04 (distance from surface 13 m), FRE-WATR10 (distance from surface 43 m), GUI-NRSROT (distance from surface 15 m), and GUI-WACA20 (distance from surface 32 m).Differences arise primarily from radar radial velocity often exceeding 1.5 m/s, which are inconsistent with both HFR data at nearby locations, or in comparison to mooring velocity data.At these locations, comparison metrics improve significantly when 99% thresholds are applied to the HFR observations with increased (decreased) R (RMSD) values (Table 2).Comparison metrics are in good agreement with previously reported findings in different ocean regions [24][25][26].Significant differences and poor statistical comparison results can be observed for the radar-mooring pairs FRE-WATR04 (distance from surface 13 m), FRE-WATR10 (distance from surface 43 m), GUI-NRSROT (distance from surface 15 m), and GUI-WACA20 (distance from surface 32 m).Differences arise primarily from radar radial velocity often exceeding 1.5 m/s, which are inconsistent with both HFR data at nearby locations, or in comparison to mooring velocity data.At these locations, comparison metrics improve significantly when 99% thresholds are applied to the HFR observations with increased (decreased) R (RMSD) values (Table 2).

SNR Threshold Optimization
The SNR of the Doppler lines in a Doppler spectra is in general considered to be a good proxy for HFR data quality and, in combination with velocity threshold, are more than adequate to identify

SNR Threshold Optimization
The SNR of the Doppler lines in a Doppler spectra is in general considered to be a good proxy for HFR data quality and, in combination with velocity threshold, are more than adequate to identify and flag the majority of suspect anomalous data.For direction-finding systems such as the SeaSonde systems, it was shown in particular that low SNR values are useful for identifying anomalous data and that data quality generally improves as SNR increases [27].Analyses of the five-month data set used here suggest that similar assumptions hold for phased-array radars such as the WERA systems.Distribution of radial velocities for different classes of SNR (not shown here) suggest that suspect radial velocities exceeding 1.5 m/s (absolute value) magnitude tend to be clustered between 5 and 10 dB, while radial velocities between [−1.5, −1] m/s and [1, 1.5] or +/− 1 m/s span similar ranges for SNR, including values below the 10 dB threshold.As for radial velocity, SNR thresholds are usually defined a priori under the rule of thumb that data quality improves with increasing SNR values.Thresholds are also generally assumed to be time-or space-invariant; i.e., it is assumed they are constant with time or space, and regional, spatial, or temporal variations due to such factors as noise levels or the transmit patterns are not taken into account.Restrictive thresholds can have the undesirable effect of removing valid observations; on the other hand, the opposite would happen with poorly constrained thresholds.
An optimization of the SNR threshold is possible and desirable using independent datasets: different threshold values are used, poorly SNR-constrained radial velocities are removed, and R and RMSD metrics are computed with the aim of finding that particular value for which a significant change in (R, RMSD) is obtained with an acceptable data loss.Results for the Rottnest HFR deployment are summarized in Figure 3. Figure 3 in particular shows that the (R, RMSD) for SNR thresholds in the 0-40 dB range, while Table 3 shows results only for the default threshold used in the proprietary software (6 dB) and the 10 dB threshold for which no more statistically significant changes are observed.
and flag the majority of suspect anomalous data.For direction-finding systems such as the SeaSonde systems, it was shown in particular that low SNR values are useful for identifying anomalous data and that data quality generally improves as SNR increases [27].Analyses of the five-month data set used here suggest that similar assumptions hold for phased-array radars such as the WERA systems.Distribution of radial velocities for different classes of SNR (not shown here) suggest that suspect radial velocities exceeding 1.5 m/s (absolute value) magnitude tend to be clustered between 5 and 10 dB, while radial velocities between [−1.5, −1] m/s and [1, 1.5] or +/− 1 m/s span similar ranges for SNR, including values below the 10 dB threshold.As for radial velocity, SNR thresholds are usually defined a priori under the rule of thumb that data quality improves with increasing SNR values.Thresholds are also generally assumed to be time-or space-invariant; i.e., it is assumed they are constant with time or space, and regional, spatial, or temporal variations due to such factors as noise levels or the transmit patterns are not taken into account.Restrictive thresholds can have the undesirable effect of removing valid observations; on the other hand, the opposite would happen with poorly constrained thresholds.
An optimization of the SNR threshold is possible and desirable using independent datasets: different threshold values are used, poorly SNR-constrained radial velocities are removed, and R and RMSD metrics are computed with the aim of finding that particular value for which a significant change in (R, RMSD) is obtained with an acceptable data loss.Results for the Rottnest HFR deployment are summarized in Figure 3. Figure 3 in particular shows that the (R, RMSD) for SNR thresholds in the 0-40 dB range, while Table 3 shows results only for the default threshold used in the proprietary software (6 dB) and the 10 dB threshold for which no more statistically significant changes are observed.In general, (R, RMSD) values are consistent with other radar validation analyses performed in different environments [24][25][26].Little or no effect is introduced for SNR values below 4 dB, in agreement with the minimum threshold applied during the inversion of the Doppler spectra to radial currents.For both HFR stations used here (FRE and GUI), the most statistically significant changes occur between 5 and 10 dB.In this range, RMSD values drop from 40 to 10-12 cm/s, R values increase from R < 0.2 to R > 0.7, while data loss is within 10-15%.The analysis also shows that some locations are less sensitive than others to variations in the SNR thresholds, and no statistically significant changes can be detected regardless of the applied SNR thresholds.

Advanced Artifact Removal
When correctly chosen or fine-tuned, thresholds on radial velocity or SNR (Sections 3.1 and 3.2) are in general capable of handling the majority of anomalous velocities in a radial map.There may be cases, however, when artifacts appear in the radial velocity field, such as unrealistic velocities at fixed and known range cells, i.e., distances, from the HFR receiver.Most of the time, similar artifacts are introduced by either external radio frequency interference (RFI) or by modulations of the 50-60 Hz 220 V power line.As a consequence, strong signals appear in the range-Doppler spectra at multiples of the 50 Hz frequency, for instance at approximately 20, 40, 60, 80, 100, and 120 km offshore (Figure 4a); at higher harmonics, this signal also tends to spread over frequency and potentially interfere with the detection of the Doppler peaks, thus introducing spurious radial currents.
There are several possible ways of dealing with this contamination.A first approach, presently implemented and operational across the Australian Ocean Radar HFR network, makes use of an iterative method that fits a 1-D or 2-D reference signal to the radial SNR map, and identifies and flags anomalous data that clearly stand out from the background SNR values.Other approaches, currently under investigation, act on the I/Q time series of the range-gated signals at each antenna, and subsequently identify and filter any anomalous signal in the frequency domain in a manner similar to [28][29][30][31].The conventional beam-forming and the Doppler spectra inversion steps are then applied on the filtered data.This section focuses on the first approach, which is operationally implemented across the network for both NRT and DM operations.Other methods still need to be In general, (R, RMSD) values are consistent with other radar validation analyses performed in different environments [24][25][26].Little or no effect is introduced for SNR values below 4 dB, in agreement with the minimum threshold applied during the inversion of the Doppler spectra to radial currents.For both HFR stations used here (FRE and GUI), the most statistically significant changes occur between 5 and 10 dB.In this range, RMSD values drop from 40 to 10-12 cm/s, R values increase from R < 0.2 to R > 0.7, while data loss is within 10-15%.The analysis also shows that some locations are less sensitive than others to variations in the SNR thresholds, and no statistically significant changes can be detected regardless of the applied SNR thresholds.

Advanced Artifact Removal
When correctly chosen or fine-tuned, thresholds on radial velocity or SNR (Sections 3.1 and 3.2) are in general capable of handling the majority of anomalous velocities in a radial map.There may be cases, however, when artifacts appear in the radial velocity field, such as unrealistic velocities at fixed and known range cells, i.e., distances, from the HFR receiver.Most of the time, similar artifacts are introduced by either external radio frequency interference (RFI) or by modulations of the 50-60 Hz 220 V power line.As a consequence, strong signals appear in the range-Doppler spectra at multiples of the 50 Hz frequency, for instance at approximately 20, 40, 60, 80, 100, and 120 km offshore (Figure 4a); at higher harmonics, this signal also tends to spread over frequency and potentially interfere with the detection of the Doppler peaks, thus introducing spurious radial currents.
There are several possible ways of dealing with this contamination.A first approach, presently implemented and operational across the Australian Ocean Radar HFR network, makes use of an iterative method that fits a 1-D or 2-D reference signal to the radial SNR map, and identifies and flags anomalous data that clearly stand out from the background SNR values.Other approaches, currently under investigation, act on the I/Q time series of the range-gated signals at each antenna, and subsequently identify and filter any anomalous signal in the frequency domain in a manner similar to [28][29][30][31].The conventional beam-forming and the Doppler spectra inversion steps are then applied on the filtered data.This section focuses on the first approach, which is operationally implemented across the network for both NRT and DM operations.Other methods still need to be optimized in order to account for the spreading over Doppler for the higher-order harmonics (Figure 4a).
Remote Sens. 2018, 10, x FOR PEER REVIEW 8 of 15 implemented across the network for both NRT and DM operations.Other methods still need to be optimized in order to account for the spreading over Doppler for the higher-order harmonics (Figure 4a).An example of heavily contaminated range-Doppler spectra and corresponding radial map is given in Figures 4a,b, respectively.Data shown here refer to the HFR station at Red Rock (RRK), NSW, where the 50 Hz contamination is particularly significant, although it is commonly observed at all the installations with the exception of Cape Wiles (CWI), SA.The range-Doppler spectrum-derived following [32,33] (Figure 4a) shows the expected Bragg peaks from the backscattering ocean waves, but also strong signals at 20, 40, 60, 80, 100, 120, and 140 km offshore, with SNR values comparable in magnitude to the dominant Bragg peaks.While at closer ranges this feature is well separated from the true Doppler peaks, the spreading over frequency that is observed with range biases the Doppler peak detection at 60, 80, 100, 120, and 140 km offshore.Effects are evident in both the radial and SNR maps, where well-defined rings are clearly detected and inconsistent radial current patterns are found.Since the SNR values associated with radial currents at these range cells typically exceed 30-40 dB, any spike identification based on the conventional SNR threshold would fail detection, remove potentially valid observations, and propagate anomalous observations through to the vector mapping step.Similarly, a conventional threshold on radial velocity would possibly remove the obviously erroneous data points; however, the method would fail to identify the problematic range cells.For this specific example, poorly SNR constrained radial currents are still within realistic ranges ([−1.5, 1.5] m/s).Radial velocities exceeding a maximum value of 1.5 m/s are flagged; some of them are associated with poor SNR values, whilst others exceed the 10 dB threshold.Peaks at 80 and 120 km cannot be identified using standard constraints on SNR or radial velocity.In contrast, most of the observations associated with the peak at the 140 km range may be flagged and removed with SNR and speed combined.
The artifact removal procedure developed at the Ocean Radar facility fits either a 1-dimensional polynomial to the SNR distribution along each radial direction, or a 2-D surface on the SNR map and, on the basis of a "distance" between the model and the data, identifies and flags suspect data.The procedure is iterative and continuously updates the distance between the fit model and the observations until no suspect data points are found, or a maximum number of iterations is reached.Detailed steps are provided below: An example of heavily contaminated range-Doppler spectra and corresponding radial map is given in Figure 4a,b, respectively.Data shown here refer to the HFR station at Red Rock (RRK), NSW, where the 50 Hz contamination is particularly significant, although it is commonly observed at all the installations with the exception of Cape Wiles (CWI), SA.The range-Doppler spectrum-derived following [32,33] (Figure 4a) shows the expected Bragg peaks from the backscattering ocean waves, but also strong signals at 20, 40, 60, 80, 100, 120, and 140 km offshore, with SNR values comparable in magnitude to the dominant Bragg peaks.While at closer ranges this feature is well separated from the true Doppler peaks, the spreading over frequency that is observed with range biases the Doppler peak detection at 60, 80, 100, 120, and 140 km offshore.Effects are evident in both the radial and SNR maps, where well-defined rings are clearly detected and inconsistent radial current patterns are found.Since the SNR values associated with radial currents at these range cells typically exceed 30-40 dB, any spike identification based on the conventional SNR threshold would fail detection, remove potentially valid observations, and propagate anomalous observations through to the vector mapping step.Similarly, a conventional threshold on radial velocity would possibly remove the obviously erroneous data points; however, the method would fail to identify the problematic range cells.For this specific example, poorly SNR constrained radial currents are still within realistic ranges ([−1.5, 1.5] m/s).Radial velocities exceeding a maximum value of 1.5 m/s are flagged; some of them are associated with poor SNR values, whilst others exceed the 10 dB threshold.Peaks at 80 and 120 km cannot be identified using standard constraints on SNR or radial velocity.In contrast, most of the observations associated with the peak at the 140 km range may be flagged and removed with SNR and speed combined.
The artifact removal procedure developed at the Ocean Radar facility fits either a 1-dimensional polynomial to the SNR distribution along each radial direction, or a 2-D surface on the SNR map and, on the basis of a "distance" between the model and the data, identifies and flags suspect data.The procedure is iterative and continuously updates the distance between the fit model and the observations until no suspect data points are found, or a maximum number of iterations is reached.Detailed steps are provided below: 1.
fit a polynomial model to the SNR distribution along each radar beam (1-D case) or fit a 2-D surface to the SNR map; 2.
estimate a "distance" between the data and the polynomial fit and set confidence levels as n * σ(data-fit), where n = 2, 3; 3.
remove suspect SNR data outside of the confidence levels; 4.
move to the next radial beam (1-D case); 5.
repeat Steps 1-4 until (i) no more anomalous data are detected or (ii) a maximum number of iterations is reached.
For the 1-D case, the actual implementation of the algorithm allows for the choice of a second-order, third-order, or exponential function, while for the 2-D case two options that fit either a third-or fifth-order surface to the SNR maps are available.An example of the results of the 1-D iterative identification and removal along one specific radial beam (92 • NCW) is provided in Figure 5 for the cubic fit (Figure 5a), the quadratic fit (Figure 5b), and exponential decay (Figure 5c).No other thresholds are applied to this example.
1. fit a polynomial model to the SNR distribution along each radar beam (1-D case) or fit a 2-D surface to the SNR map; 2. estimate a "distance" between the data and the polynomial fit and set confidence levels as n * σ(data-fit), where n = 2, 3; 3. remove suspect SNR data outside of the confidence levels; 4. move to the next radial beam (1-D case); 5. repeat Steps 1-4 until (i) no more anomalous data are detected or (ii) a maximum number of iterations is reached.
For the 1-D case, the actual implementation of the algorithm allows for the choice of a second-order, third-order, or exponential function, while for the 2-D case two options that fit either a third-or fifth-order surface to the SNR maps are available.An example of the results of the 1-D iterative identification and removal along one specific radial beam (92° NCW) is provided in Figure 5 for the cubic fit (Figure 5a), the quadratic fit (Figure 5b), and exponential decay (Figure 5c).No other thresholds are applied to this example.In order to test the effectiveness of the 1-D and 2-D formulation and optimize the choice of the fit function, a series of sensitivity tests have been performed on a synthetic SNR map with realistic patterns, for instance matching the one that generated the radial current shown in Figure 4b.The SNR model used here is a simple linear decay over range for all bearings, which starts from 60 dB at the grid cell closest to the receiver site and decreases to 4 dB at 180 km offshore.A Gaussian-shaped noise function is added to range cells 90-95 and 120-128 with a peak amplitude of 50 dB, and to range cells 150-159 with a peak amplitude of 15 dB, so as to reproduce a typical pattern in the data.An additional noise term is included in the sensitivity tests, with variance levels of 0 (no noise), 1, 2, 5, and 10 dB.No assumptions are made in regard to directional patterns.Sensitivity analysis aimed to determine the speed, the number of iterations required to flag and remove the known artifacts, and the number of "false-flag" i.e., incorrect detections, as a function of the noise level.
Results for the 1-D case for both the noise-free and randomly distributed noise cases can be summarized as follows.In the absence of noise, the quadratic, cubic, and exponential fit equally succeed in detecting and flagging more than 70% of the anomalies, while minimizing the false-flag detection (0% in this specific realization).On the other hand, the exponential fit requires a significantly larger number of iterations to achieve results comparable to the quadratic and cubic function; however, there is no major improvement in the detection rate as the number of iterations increases.Detection rates decrease when noise is added to the SNR maps, but the lower success rate is related to the missed identification of the weaker SNR peaks at the outer ranges rather than an increasing rate of false-detection.
When applied to field data, results suggest that a third-order polynomial (Figure 6c) is more effective in identifying and removing anomalous observations from the radial maps than the quadratic (Figure 6b) or exponential fit (Figure 6d).For this specific example, the artifact removal is applied to radial currents from one station (RRK), and additional thresholds have been applied to radial SNR and radial velocity before the polynomial fit.
Results for the 2-D fit improve with respect to the 1-D case when noise is added to the synthetic SNR test case, but no improvement is obtained when a noise term is excluded.In general, the detection rate decreases with the increasing noise level similarly to the 1-D case, but the false detection rate is significantly lower in the 2-D case for the same noise level.Processing speed also increases slightly in respect to the 1-D case when noise is added to the SNR map.When applied to a real radial velocity map, no major differences can be found between the third-or fifth-order surfaces, yet processing speed greatly improves-up to a factor of 10 compared to the 1-D case.In order to test the effectiveness of the 1-D and 2-D formulation and optimize the choice of the fit function, a series of sensitivity tests have been performed on a synthetic SNR map with realistic patterns, for instance matching the one that generated the radial current shown in Figure 4b.The SNR model used here is a simple linear decay over range for all bearings, which starts from 60 dB at the grid cell closest to the receiver site and decreases to 4 dB at 180 km offshore.A Gaussian-shaped noise function is added to range cells 90-95 and 120-128 with a peak amplitude of 50 dB, and to range cells 150-159 with a peak amplitude of 15 dB, so as to reproduce a typical pattern in the data.An additional noise term is included in the sensitivity tests, with variance levels of 0 (no noise), 1, 2, 5, and 10 dB.No assumptions are made in regard to directional patterns.Sensitivity analysis aimed to determine the speed, the number of iterations required to flag and remove the known artifacts, and the number of "false-flag" i.e., incorrect detections, as a function of the noise level.
Results for the 1-D case for both the noise-free and randomly distributed noise cases can be summarized as follows.In the absence of noise, the quadratic, cubic, and exponential fit equally succeed in detecting and flagging more than 70% of the anomalies, while minimizing the false-flag detection (0% in this specific realization).On the other hand, the exponential fit requires a significantly larger number of iterations to achieve results comparable to the quadratic and cubic function; however, there is no major improvement in the detection rate as the number of iterations increases.Detection rates decrease when noise is added to the SNR maps, but the lower success rate is related to the missed identification of the weaker SNR peaks at the outer ranges rather than an increasing rate of false-detection.
When applied to field data, results suggest that a third-order polynomial (Figure 6c) is more effective in identifying and removing anomalous observations from the radial maps than the quadratic (Figure 6b) or exponential fit (Figure 6d).For this specific example, the artifact removal is applied to radial currents from one station (RRK), and additional thresholds have been applied to radial SNR and radial velocity before the polynomial fit.
Results for the 2-D fit improve with respect to the 1-D case when noise is added to the synthetic SNR test case, but no improvement is obtained when a noise term is excluded.In general, the detection rate decreases with the increasing noise level similarly to the 1-D case, but the false detection rate is significantly lower in the 2-D case for the same noise level.Processing speed also increases slightly in respect to the 1-D case when noise is added to the SNR map.When applied to a real radial velocity map, no major differences can be found between the third-or fifth-order surfaces, yet processing speed greatly improves-up to a factor of 10 compared to the 1-D case.

Discussion
Routinely, QC procedures for HFR radial velocity maps rely on the definition of threshold values that are applied to SNR or radial velocity speed.Examples of operational procedures applied in RT mode for HFR in general and WERA systems in particular can be found in [28,29], while advanced methodologies capable of handling complex cases of external radio frequency interferences are described in [13].In most cases, these values are defined "a priori" without any assumption on the real range of variability or the physical processes in the region.This work suggests that a fine-tuning of the thresholds is beneficial to improving the overall quality of the HFR data set.Velocity thresholds can be set through comparisons with available data from independent sampling devices or, alternatively, using long-term statistics of radar observations.A proper definition of the radial velocity threshold significantly improves the comparison metrics at several locations across the HFR domain for both FRE and GUI radar stations (Table 2).

Discussion
Routinely, QC procedures for HFR radial velocity maps rely on the definition of threshold values that are applied to SNR or radial velocity speed.Examples of operational procedures applied in RT mode for HFR in general and WERA systems in particular can be found in [28,29], while advanced methodologies capable of handling complex cases of external radio frequency interferences are described in [13].In most cases, these values are defined "a priori" without any assumption on the real range of variability or the physical processes in the region.This work suggests that a fine-tuning of the thresholds is beneficial to improving the overall quality of the HFR data set.Velocity thresholds can be set through comparisons with available data from independent sampling devices or, alternatively, using long-term statistics of radar observations.A proper definition of the radial velocity threshold significantly improves the comparison metrics at several locations across the HFR domain for both FRE and GUI radar stations (Table 2).
In agreement with previous studies, SNR can be considered a valid proxy for HFR data quality.However, similarly to radial velocity, the SNR threshold needs to be tuned for each specific deployment, in order to avoid removing potentially valid observations.Tests performed across several mooring locations also suggest that a single threshold value may not be adequate for the entire domain and that this threshold may vary over time as a result of varying environmental factors, external noise, interference sources, or other problems.The analysis performed in Section 3.2 suggests that small variations in the SNR threshold can greatly improve data quality without significantly compromising the data availability.For example, setting a minimum SNR of 7-8 dB instead of the default threshold results in significant improvements (Table 3).Operationally, this SNR threshold was increased to 10 dB for the Australian HFR systems, and the combined effects of SNR and radial velocity thresholds are documented in Table 4. Thresholds on radial velocity or SNR may both be inadequate when more complex cases arise, such as the contamination of the 50 Hz signal, that can significantly bias direction and speed of the resulting current vector when propagated to the step of vector mapping.For this purpose, a fast and effective iterative approach has been developed that iteratively identifies and removes anomalous radial observations from a radial map based on the spatial distribution of the SNR maps and the basic assumption that its spatial distribution can be modeled as a cubic polynomial (1-D case) or a third-degree surface in a Cartesian (X, Y) plane.The first method can be considered an improvement to the method proposed in [10] in its application to ship detection and tracking.In combination with recent advances in the waveform design and implementation [13], the proposed method has the potential to improve radial data quality for WERA HFR systems.The actual implementation of the procedure uses a Python code for the NRT applications where computational speed needs to be optimized to run every 5 min across the network.An equivalent version in MATLAB language, which also includes the 2-D fit, can be run both in NRT and in DM for debugging, testing, or development purposes.

Conclusions
NRT and DM QC procedures developed for ocean HFR systems can greatly improve the quality of HFR station data and resulting vector maps, thus facilitating data uptake and usage.This study describes new methods that optimize data quality using a combination of radial velocity and SNR thresholds.Thresholds can be set through quantitative analyses and comparisons with independent data sets, when available, or a self-validation approach that uses distributions and statistics from HFR data.Methodologies can be applied both in NRT and DM operations, and, although designed for radial velocity data collected from commercial-type WERA HF systems, they can easily be extended to any other type or HFR for oceanographic use.

Figure 1 .
Figure 1.Coverage of the HFR systems in the Rottnest Shelf deployment, WA, with the dominant bathymetric features.Red (black) markers are used for the WERA (SeaSonde) FRE and GUI stations (GHED and LANC); red (black) areas show the typical coverage of each HFR station.Mooring locations are marked as black diamonds inside the WERA HFR footprint.

Figure 1 .
Figure 1.Coverage of the HFR systems in the Rottnest Shelf deployment, WA, with the dominant bathymetric features.Red (black) markers are used for the WERA (SeaSonde) FRE and GUI stations (GHED and LANC); red (black) areas show the typical coverage of each HFR station.Mooring locations are marked as black diamonds inside the WERA HFR footprint.

Figure 2 .
Figure 2. Examples of radial velocity distributions for FRE HFR station, and corresponding subsurface radial velocity components, at two deep-water mooring locations within the radar domain: (a) 56 m from the surface; (b) 32 m from the surface.Units for radial velocity are m/s.

Figure 2 .
Figure 2. Examples of radial velocity distributions for FRE HFR station, and corresponding subsurface radial velocity components, at two deep-water mooring locations within the radar domain: (a) 56 m from the surface; (b) 32 m from the surface.Units for radial velocity are m/s.

Figure 3 .
Figure 3. Correlation (R; black solid lines), RMSD (dashed grey lines), and percent data available (solid grey line) as a function of different SNR threshold values for (a) GUI and (b) FRE HFR stations, for all the mooring deployments in Figure 1.Units for RMSD are cm/s.

Figure 4 .
Figure 4. Range-Doppler spectra for Receive Antenna #09 at the Red Rock (RRK) HFR station (a) showing both the contamination from the 50 Hz across range (distance) from the receive array in the first order Bragg peaks and the spectral spreading over Doppler; (b) the radial velocity map extracted from this data set, documenting the spatial inconsistencies in radial velocity maps at fixed annular ranges from the antennas.Peak SNR values in these annular regions typically exceed 40 dB.

Figure 4 .
Figure 4. Range-Doppler spectra for Receive Antenna #09 at the Red Rock (RRK) HFR station (a) showing both the contamination from the 50 Hz across range (distance) from the receive array in the first order Bragg peaks and the spectral spreading over Doppler; (b) the radial velocity map extracted from this data set, documenting the spatial inconsistencies in radial velocity maps at fixed annular ranges from the antennas.Peak SNR values in these annular regions typically exceed 40 dB.

Figure 5 .
Figure 5. Example of the effects of the proposed 1-D QC methodology on radial velocity map from RRK HFR station along one radial beam: (a) cubic 1-D fit; (b) quadratic 1-D fit; (c) exponential decay 1-D fit.Units are km for offshore range, dB for SNR, and m/s for radial velocity (color scale).

Figure 6 .
Figure 6.Effects of the proposed QC methodology on radial velocity map from the RRK HFR station: (a) no QC applied; (b) quadratic 1-D fit; (c) cubic 1-D; (d) exponential decay 1-D fit.(b-d) also include threshold on radial velocity and SNR.Units for radial velocity are m/s.

Figure 6 .
Figure 6.Effects of the proposed QC methodology on radial velocity map from the RRK HFR station: (a) no QC applied; (b) quadratic 1-D fit; (c) cubic 1-D; (d) exponential decay 1-D fit.(b-d) also include threshold on radial velocity and SNR.Units for radial velocity are m/s.

Table 1 .
Operational settings for the High-Frequency Radar (HFR) systems across Australia; locations of the HFR systems for the Western Australia (WA) node are given in Figure1.

Table 2 .
Effects of velocity thresholds on radar data: a , no threshold applied; b , 1.5 m/s (absolute value) threshold applied.

Table 3 .
Effects of SNR thresholds on radar data: a no threshold applied; b threshold applied.

Table 4 .
Effects of SNR and velocity thresholds on radar data: a no threshold applied; b threshold applied.