## Abstract

This paper describes the instrumentation and techniques for long-term targeted observation of the centimeter-scale velocity structure within the oceanic surface boundary layer, made possible by the recent developments in capabilities of autonomous platforms and self-contained pulse-coherent acoustic Doppler current profilers (ADCPs). Particular attention is paid to the algorithms of ambiguity resolution (“unwrapping”) of pulse-coherent Doppler velocity measurements. The techniques are demonstrated using the new Nortek Signature1000 ADCP mounted on a Lagrangian float, a combination shown to be capable of observing ocean turbulence in a number of recent studies. Statistical uncertainty of the measured velocities in relation to the ADCP setup is also evaluated. Described techniques and analyses should be broadly applicable to other autonomous and towed applications of pulse-coherent ADCPs.

## 1. Introduction

Ocean current velocities vary on a wide range of horizontal, vertical, and temporal scales, particularly in the turbulent surface and bottom boundary layers, presenting a substantial observational challenge (D’Asaro 2014). Turbulent velocities are known to be patchy in space and intermittent in time; at the same time, direct velocity observations are comparatively sparse, mostly limited to moored and shipboard measurements and inference from the trajectories of drifters and floats. Accelerating development of autonomous observational platforms (Perry and Rudnick 2003) combined with the advances in the acoustic Doppler velocity measurement technology would enhance our ability to conduct persistent targeted observations of ocean currents worldwide.

There are two basic methods of oceanographic acoustic Doppler velocity measurements: incoherent (single pulse) and pulse-to-pulse coherent. Both methods essentially rely on detecting the phase shift of the reflected acoustic signals arising from the motion of the scatterers embedded in the water current. The first method, dating back to Pinkel (1979), estimates the Doppler shift using a single return of a relatively long pulse, allowing long-range profiling. Later (“broadband”) modifications of this method employed pseudorandom modulation of the sound pulse (Brumley et al. 1991), allowing substantial reduction of the relatively high noise levels characteristic of this method. The second method uses a pair of relatively short pulses that are separated in time but maintain coherence (Lhermitte and Serafin 1984); the first pulse is received and is allowed to die down before the second one is transmitted. This pulse-to-pulse coherent method achieves improved accuracy and range resolution, but it has a relatively shorter profiling range and suffers from phase ambiguity [see section 2c(1)]. We will refer to the single-pulse and pulse-to-pulse coherent Doppler velocity-profiling methods as long-range (LR) and high-resolution (HR), respectively.

The new Nortek Signature acoustic Doppler current profiler (ADCP) used in this study is capable of interleaving LR and HR profiling modes, allowing a uniquely rich view of the ocean current structure. When such an ADCP is mounted on an autonomous platform [a Lagrangian float (LF) in our case], unattended open-ocean turbulence observations over prolonged periods can be achieved.

Our work follows a number of prior implementations of predominantly broadband ADCPs on autonomous and semiautonomous platforms over the past decades, including self-propelled autonomous underwater vehicles (Fong and Jones 2006), underwater gliders (Todd et al. 2011; Ellis et al. 2015), drifting surface buoys (Thomson et al. 2015), moored profilers (Rusello et al. 2011), etc. The main obstacles in the broader implementation of HR ADCPs were the relatively low ranges (compared to LR broadband ADCPs) and low tolerance to platform motion and relative current speeds (RD Instruments 2002). A combination of LF platform capabilities with the new instrumentation and algorithms allowed us to relax these limitations substantially.

Here, we will focus specifically on processing and evaluating HR velocity observations, as well as the unique advantages offered by HR/LR interleaved sampling. Processing of autonomous LR observations and reconstruction of the overall velocity profile are conceptually similar to these tasks for the lowered ADCP (LADCP) (Fischer and Visbeck 1993; Visbeck 2002), an electromagnetic velocity-profiling float (Sanford et al. 2005), or an ocean glider (Todd et al. 2011; Ellis et al. 2015). Appropriate techniques have been discussed in detail in the respective publications.

Even though this article describes a particular instrument mounted on a particular platform, the methods of acquisition and postprocessing of HR Doppler data should be broadly applicable, especially to other autonomous vehicles equipped with Doppler profilers such as underwater gliders or surface drifters, as well as moored and towed systems.

## 2. Materials and procedures

### a. Signature1000

Nortek Signature is a new line of broadband acoustic Doppler current profilers, including 55-, 250-, 500-, and 1000-kHz versions. The instruments are based on the patented acoustic Doppler dual current profiler (AD2CP) platform (U.S. Patent 7,911,880) that introduced the ability to interleave acoustic sampling with different configurations of the transmit pulse and signal processing. Effectively, this feature allows quasi-simultaneous measurements within the same water volume with two “virtual ADCPs” (hence “dual” in AD2CP).

The 500- and 1000-kHz Nortek Signature ADCPs are five-beam self-contained instruments, with four beams in a standard symmetric Janus configuration pointing 25° off the ADCP axis and a fifth axial (vertical) beam. The instruments are equipped with a pressure sensor, a thermistor, a three-axis accelerometer, and a magnetometer. A smaller externally powered four-beam version of the Signature1000 is also available, designed to be integrated into underwater ocean gliders (Siegel and Rusello 2013).

A distinguishing feature of the highest-frequency Signature1000 instrument used in this study is the HR pulse-coherent sampling capability. In keeping with the AD2CP concept, HR and LR sampling can be interleaved [this feature is exploited in section 2c(5)]. HR sampling can be configured with up to 256 range cells, with the minimum cell size of 2 cm. The maximum useable range is roughly 4−8 m, limited by the amount of scatterers in the water and flow velocities (see section 2c).

### b. Lagrangian float platform

The LF (Fig. 1) is a proven versatile platform for upper-ocean observations developed and built at the Applied Physics Laboratory of the University of Washington (D’Asaro 2003). These 1.5-m-long instruments were originally designed to measure turbulence in the ocean mixed layer by accurately following the three-dimensional motion of water parcels through a combination of actively maintained neutral buoyancy and high drag provided by controllable flexible drogues. Accurate and highly adaptable automatic buoyancy control and a relatively heavy payload distinguish LFs from other types of profiling floats. The floats can provide uniform sampling of the mixed layer turbulence in fully Lagrangian water-following mode (D’Asaro et al. 2014). They can also operate in either isopycnal or isopycnal/Lagrangian mode in the pycnocline (D’Asaro 2008), allowing observations of turbulence associated with internal waves. Profiling across a given depth range is also a possibility. A typical mission combines these modes of operation to adapt to the changing stratification and forcing. Because of these capabilities, LFs have been increasingly used as platforms for targeted open-ocean observations with a wide range of experimental payloads, such as gas-tension sensors (D’Asaro and McNeil 2008), ambient-sound recorders (Zhao et al. 2014), a finescale profiling thermistor–fluorometer combo (Shcherbina et al. 2015), and multisensor biogeochemical payloads (Alkire et al. 2012).

LFs are particularly attractive platforms for high-resolution turbulence measurements because of their underwater stability, with the typical RMS tilt of <0.05° below the surface layer influenced by the waves. Because of the water-following properties of LFs in fully Lagrangian and isopycnal modes, nearby relative velocities are always close to zero, minimizing the range of velocities that need to be resolved by the onboard Doppler profilers. Relative vertical velocities may reach 20−30 cm s^{−1} during fast profiling, but this profiling speed can be readily obtained by differentiation of concurrent high-resolution pressure measurements.

### c. HR Processing

#### 1) Doppler phase shift estimation

Modern Doppler sonars, including the Nortek Signature1000, typically measure the radial velocity of scatterers in the water by detecting temporal change in the phase of the reflected acoustic signal. In pulse-coherent (HR) mode, Doppler velocity *υ* is obtained from the phase *ϕ* of the complex correlation of the echoes of the two pulses separated by time lag *τ*, which is typically much longer than the pulse duration—the so-called pulse-pair covariance method (Miller and Rochwarger 1972; Lhermitte and Serafin 1984; Zedel et al. 1996),

where *c* is the sound speed, *F*_{0} is the sonar carrier frequency, and *C*(*τ*) is the complex pulse-pair correlation, which is estimated as

where *Z*_{m} = *I*_{m }*+ iQ*_{m} are the complex representations of the *m*th pulse echo created from the signal components that are in phase (*I*_{m}) and in quadrature (*Q*_{m}) with the transmitted waveform, and the asterisk represents a complex conjugate.^{1} The implicit time variable *t* is determined by the range gate corresponding to a particular range cell; it will be omitted for brevity. The amplitude of the complex pulse-pair correlation |*C*(τ)| is a valuable metric for the uncertainty of the measurement (see section 3c).

In broadband incoherent (LR) mode, the lagged autocorrelation of a single modulated pulse is used in a similar way; in this case, *τ* represents the lag between the modulation sequences within the pulse or, equivalently, the pulse decorrelation time.

Finally, a hybrid of the two modes could be used, where each of the two HR pulses is itself broadband, which is achieved through either the coded pulse modulation (Pinkel and Sherman 1990; Brumley et al. 1991) or, as in the case of the Nortek Signature1000, the chirp modulation (Hay et al. 2015). Such broadband pulses have shorter decorrelation times, therefore increasing the number of quasi-independent subsamples that could be obtained from a single pulse pair. For a chirp with the duration *T* and the bandwidth Δ*F*, the number of subsamples is given by the bandwidth–time product *M =* Δ*F T.* The pulse duration is typically matched to the range cell size *h*, giving *T* = 2*h/c* and *M =* 2*h*Δ*F/c.* For multiple subsamples, the complex pulse correlation is estimated as

where angle brackets represent ensemble averaging over *M* subsamples.

An inherent limitation of Doppler methods arises from “phase wrapping”—a ±2π*n* ambiguity in the determination of the phase shift in (1) that results in ±2*nV*_{R} uncertainty of measured radial velocity, where *V*_{R} = *c*/(4*F*_{0}*τ*) is the maximum unambiguous measured velocity amplitude (or “velocity range”) and *n* is an arbitrary integer (Figs. 2, 3). Unless accounted for, phase wrapping is particularly detrimental for broadband pulses that can become completely decorrelated within the span of a few 2π phase wraps (Fig. 2b), leading to an unacceptable increase in the measurement noise level (Hay et al. 2015).

Phase wrapping is rarely a problem with LR measurements, because they are typically configured with sufficiently short modulation sequences, leading to comfortably large velocity range values (1−5 m s^{−1} for the chirp modulation in the Signature1000). On the other hand, HR measurements must have a pulse separation longer than the round-trip sound travel time to the furthest measurement cell. As a result, the maximum profiling range *L*_{max} of a pulse-to-pulse coherent Doppler sonar and the radial velocity range are linked with a fundamental trade-off relationship *V*_{R }*L*_{max }*= c*^{2}/(8*F*_{0}) ≈ 0.28 m^{2} s^{−1} for 1-MHz sonar (Lhermitte and Serafin 1984). Achieving a profiling range of 5–10 m therefore necessitates *V*_{R} on the order of 3–5 cm s^{−1}. The accuracy of the HR velocity measurements also tends to improve with increasing pulse separation (see section 3). Therefore, it is not optimal to set up an HR profiler with such short pulse separation that phase wrapping never happens. Instead, we can rely on a number of “phase unwrapping” methods that can be applied both during the signal processing and in subsequent analysis.

#### 2) Multicorrelation pulse-coherent processing

The Doppler phase-wrapping problem can be partially mitigated by using advanced signal processing techniques, combining the benefits of broadband and pulse-coherent processing (Lohrmann and Nylund 2008). The Nortek Signature1000 employs the multicorrelation pulse-coherent (MCPC) extended velocity range (EVR) method, taking advantage of the increased bandwidth of the chirp HR pulses (Hay et al. 2015). Because of the chirp frequency modulation, the amplitude of the pulse correlation |*C*(*τ* + *ζ*)| has a relatively sharp peak at an additional lag ζ = 2*υτc*^{−1} relative to the nominal pulse lag *τ* (Fig. 2a). While it is not practical to estimate the Doppler velocity from the location of the peak alone, this location can be used to resolve phase ambiguity.

In the MCPC EVR method, several “side correlations,”

are obtained in addition to the “central” correlation *C*_{0} = *C*(*τ*), as illustrated in Fig. 2 (*T*_{0} = *F*_{0}^{−1} is the carrier period). The lag of the correlation maximum is then estimated using a piecewise quadratic interpolation, giving a rough velocity estimate of , which is unaffected by phase wrapping. Finally, a more accurate value of the Doppler velocity is determined from the phase of the side correlation with maximum amplitude *C*_{m} as in (1),

except that the correct number of phase wraps *n* is determined by the minimization of . As a result of this procedure, the range of unambiguously measurable velocities is extended to *V*_{EVR} = 2*n*_{max}*V*_{R}, where *n*_{max} is the number of side correlations considered, which is constrained by the real-time computational load; the current implementation of the MCPC algorithm in the Signature1000 firmware defaults to *n*_{max} = 6. An added benefit of MCPC processing is that the Doppler phase is calculated where the correlation value is at its maximum, thus ensuring optimal performance of the pulse-pair method (relationship between the signal correlation and measurement uncertainty is discussed in section 3). The MCPC method appears to be conceptually similar to that employed by the “water profiling mode 11” of the RDI Workhorse ADCPs (RD Instruments 2002) but different from the EVR mode implemented on Nortek Aquadopp HR ADCPs (Rusello 2009).

#### 3) Correcting phase resolution errors

The EVR algorithm does not eliminate the ambiguity of the measured velocities entirely (Fig. 3b). If the correlation peak is incorrectly located (e.g., because of the noise), then the velocity reading may receive an erroneous offset equivalent to a random integer multiple of the velocity range, ±2*nV*_{R}. These errors are easily identified in individual velocity profiles obtained using the EVR method as isolated spikes with the amplitude exceeding *V*_{R} and reduced correlation (typically). For each of these spikes, we seek a correction equivalent to an integer number of phase wraps, such that the corrected profile is as smooth as possible. Choosing the sum of the squared finite differences of the velocity profile as a measure of roughness, we formulate an integer linear least squares (ILS) minimization problem and solve it to obtain the necessary spike correction offsets (see the appendix for details).

The ILS phase ambiguity correction algorithm is adequate for dealing with isolated phase resolution errors caused by noise (Fig. 3b) but not systematic wrapping caused by strong velocity gradients. In these cases, more advanced phase unwrapping techniques need to be applied in postprocessing, as described in the next section.

#### 4) Phase unwrapping in postprocessing

If the EVR method is not available (e.g., in case of a legacy pulse-coherent sonar), or if the Doppler shift exceeds the maximum number of phase wraps considered by the EVR method, then additional phase unwrapping may be necessary during postprocessing. Detection and correction of phase wraps in various signals is a common problem in optics, magnetic resonance imaging (MRI), synthetic aperture radar (SAR) interferometry, and other fields (Ghiglia and Pritt 1998). Basic unwrapping can be readily achieved with the Itoh (1982) algorithm: calculate phase differences Δ*ϕ* between the adjacent points, locate the “jumps” (i.e., |Δ*ϕ*| > π), correct the jumps by ±2π so that |Δ*ϕ*| < π, and reintegrate to obtain the unwrapped phase. In practice, there are two major obstacles to accurate phase unwrapping:

The change of phase in the “true” signal between the adjacent points cannot exceed π by absolute value, because these changes would be detected as jumps and unwrapped incorrectly.

The presence of noise in the phase data with amplitude exceeding π could result in spurious 2π discontinuities in the unwrapped phase; erroneous values are then propagated through the rest of the record (Fig. 3).

The first problem can be typically addressed by selecting a fine-enough resolution of phase sampling in space and time. A Nortek Signature AD2CP with velocity range *V*_{R} = 4 cm s^{−1} and vertical cell resolution *h* = 2 cm sampling at *f*_{s} = 1 Hz would allow correct unwrapping of velocity profiles with shears of up to *V*_{R}/*h* = 2 s^{−1} and accelerations *V*_{R }*f*_{s} = 40 mm s^{−2}, which should be sufficient for most ocean boundary layer observations (except perhaps surface wave orbital motions in the upper meters of the water column).

The noise problem is more serious, as it often prevents obtaining a unique solution to the phase unwrapping problem. It cannot be addressed by the standard noise mitigation methods such as spatial or temporal averaging, because they would also remove “legitimate” jumps in the wrapped phase and make unwrapping impossible (consider Fig. 3). The problem becomes tractable if a two-dimensional wrapped phase record is considered (e.g., sampled uniformly in time and space), as the phase unwrapping can often circumvent isolated erroneous phase jumps by considering various integration paths in 2D. Two-dimensional phase unwrapping algorithms can be readily applied to postprocessing of pulse-coherent ADCP observations, but, to our knowledge, this has never been done before.

A large number of 2D phase unwrapping algorithms exists, offering varying degree of noise robustness and features specific to particular applications (Ghiglia and Pritt 1998). The algorithms fall into two broad classes: path following (or branch cut) and path independent (or global) algorithms. The former heuristically seeks a suitable unwrapping path (e.g., Arevallilo Herráez et al. 2002). The latter class of algorithms solves the global optimization problem of finding the smoothest possible unwrapped phase field (e.g., Ghiglia and Romero 1994). Global algorithms tend to be computationally intensive but more robust to noise. Most modern algorithms of either class also allow segmentation of the large datasets into a number of tiles or regions.

We evaluated a number of 2D phase unwrapping algorithms and selected two that combine good performance, readily available implementation, and features useful for HR velocity unwrapping: The first algorithm is the two-dimensional sorting by reliability following a noncontinuous path (2D-SRNCP), a fast path-following algorithm (Arevallilo Herráez et al. 2002), with C implementation available on request from Liverpool John Moores University.^{2} The second is the statistical-cost, network-flow algorithm for phase unwrapping (SNAPHU), a global optimization algorithm originally developed for SAR interferometry (Chen and Zebker 2000), with C implementation available from Stanford University.^{3} Unlike path-following algorithms, SNAPHU can make use of the signal correlation to construct an appropriately weighted optimization cost function and to improve the algorithm’s robustness to noise. SNAPHU’s statistical cost functions can be further tuned with some knowledge of the noise statistics of the observed signal (Chen and Zebker 2001), although we have not fully explored this particular feature with respect to HR Doppler velocity unwrapping.

Both algorithms clearly outperform 1D Itoh unwrapping when applied to Signature1000 AD2CP HR velocity data (Fig. 4). In the low-shear, low-noise part of the record (0−100 s in Fig. 4), algorithms produce identical results. As the noise level increases (100−200 s), 2D-SRNCP becomes prone to creating patches with erroneous phase offsets (marked in Fig. 4); in these instances, SNAPHU is more robust.

#### 5) HR/LR cross validation

The HR velocity unwrapping discussed in the previous sections is inherently nonunique, allowing an arbitrary ±2*nV*_{R} shift of the whole unwrapped record to remain. If the record is split into the reliable data segments separated by periods of noisy or invalid data, then each segment may receive a ±2*nV*_{R} offset not recoverable by 2D phase unwrapping. We commonly encounter the latter situation when the float-mounted AD2CP periodically approaches the surface. In these cases, cross validation of HR and LR velocity observations allows detecting and correcting these offsets. Such cross validation is made particularly easy with the Signature AD2CP’s capability to interleave HR and LR measurements within each second. Velocity offset is established as the mean difference between HR and LR velocities for each of the AD2CP’s beams over the overlapping range, rounded to the nearest ±2*nV*_{R} value. Relatively low resolution and a high noise level of single-ping LR velocity measurements are not obstacles for such cross validation because only the averages over prolonged data segments (hundreds of seconds) need to be compared, and agreement needs to be accurate to only within half the velocity range of HR measurements (on the order of 1 cm s^{−1}).

#### 6) Fish detection

Potential contamination of acoustic Doppler velocity measurements by swimming fish have been long recognized as a problem (Freitag et al. 1992; Wilson and Firing 1992).^{4} ADCPs depend on the presence of scatters (including fish, plankton, and density discontinuities) to infer velocity. However, if these scatterers move on their own, then Doppler velocity measurements can no longer be attributed solely to water motion. Broadband (LR) ADCPs typically do not resolve single organisms, but they would incur a velocity bias associated with schools of swimming fish. These biases have been used to quantify movement of fish and zooplankton assemblages from ADCP records (see Simmonds and MacLennan 2007 for a review).

High-resolution, low-noise AD2CP measurements resolve the motion of what appears to be individual fish (Fig. 5). These motions are too widespread and pronounced to ignore, at least in some environments. Individual fast-moving targets present a challenge for Doppler phase unwrapping: Most of the time 2D unwrapping algorithms can work around the fish-induced velocity anomalies, but velocities of the fish themselves may remain aliased (e.g., see the observed Doppler velocity of target D in Fig. 5 change sign despite the steady downward progression of the echo). If the fish velocities are of interest, then they can be unwrapped separately using the ping-to-ping echo movement as a reference.

Fish contamination of an ADCP record can be diagnosed in postprocessing by analyzing the logarithmic echo intensity *E* recorded by the instrument alongside the velocity profile. The operational algorithm used in RD Instruments’ ADCPs detects fish presence by one of the beam’s echo intensity exceeding that of the rest of the beams (Freitag et al. 1993). Similarly, we use the echo intensity anomaly relative to the median echo intensity in the same range cell across the rest of the beams. Even though the fifth (vertical) AD2CP beam is operated asynchronously with the four slanted beams, the 0.125-s delay is short enough to ignore for the purposes of fish detection. We also neglect the fact that the slanted beams’ echo intensity falls slightly faster (<10%) with the range cell number compared to that of the vertical beam as a result of the longer travel path and higher absorption. Unlike the fixed predefined fish detection threshold used in the RD Instruments algorithm, we estimate the locally appropriate threshold value for each profile as

where is the standard deviation of for a given profile and is the “universal threshold” for range cells. Instead of a sample standard deviation, we use a more robust median absolute deviation (MAD) estimate where square brackets indicate the median value over a profile (Wahl 2003).^{5} This approach amounts to flagging the echo intensity values that are improbably high given the local statistics. The resulting threshold value is typically 5–10 dB, increasing to 20–30 dB near the surface.

This echo-intensity-based fish rejection algorithm detects only unusually strong scatterers; some of them can be perfectly neutral in the water and cause no velocity bias (e.g., target E in Fig. 5). Conversely, strong velocity biasing can occur without the echo intensity rising above the fish detection threshold (e.g., target C). Therefore, additional velocity despiking steps may be necessary (e.g., Goring and Nikora 2002).

## 3. Assessment

The main application of the float-mounted HR ADCP is obtaining a close-up view of turbulence in the ocean interior away from the surface and the bottom boundaries, provided it is energetic enough. This view complements both bottom-mounted and surface-based ADCP observations that cannot extend more than a few meters into the respective boundary layers because of the limited range of HR instruments.

With the vertical resolution on the order of a few centimeters, autonomous HR ADCP observations fill the gap between the long-range broadband ADCP velocity measurements typically operating at resolutions on the order of 1–10 m and the millimeter-scale turbulence microstructure profiling. These observations also complement the characterization of large-scale (10−100 m) energy-containing turbulent eddies that can be derived from the motions of the Lagrangian float itself (Lien and D’Asaro 2006).

A Lagrangian float equipped with a Nortek Signature1000 AD2CP was deployed in Tacoma Narrows and Colvos Passage, a pair of linked tidal channels in Puget Sound, Washington. The site is characterized by highly variable levels of turbulence, with energy dissipation rates spanning the range between 10^{−8} and 10^{−4} W kg^{−1} depending on location and tidal phase (Lien and D’Asaro 2006). During the 17–20 October 2016 deployment period, average vertical stratification, computed from the float CTD profiles, corresponded to buoyancy frequency *N* = 4.5 × 10^{−3} s^{−1}. Over 3 days, the float made three full transits of the channels; after each ~24-h transit, it was moved back to the starting point and redeployed around the time of the maximum tidal currents. The float was operating in fully Lagrangian drift mode most of the time, following the turbulent motions within the tidally driven boundary layer that typically spanned the entire 120-m depth of the channel (Fig. 6f).

### a. Turbulence observations

Examples of turbulent vertical velocity structure observed with the fifth (vertical) beam of the float-mounted Signature1000 AD2CP in HR mode during the Colvos Passage deployment are shown in Figs. 6a–e. HR sampling was conducted at 2-Hz rate, interleaved with 1-Hz LR sampling (not shown). The range resolution was set to 195 2-cm cells, for a total range of 4 m (see Table 1 for details on the other setup parameters). This setup corresponded to *V*_{R} = 7 cm s^{−1}. Relative vertical velocities reaching 40 cm s^{−1} were successfully unwrapped using a combination of in-instrument EVR processing and SNAPHU algorithm in postprocessing [see section 2c(4)]. Ranges of up to 8 m (corresponding to *V*_{R} = 3 cm s^{−1}) were achievable during the periods of weaker flow, but they resulted in poor correlation and problems with phase unwrapping during the episodes of stronger turbulence.

Each HR sample of the Nortek Signature1000 AD2CP allows computation of turbulent velocity spectra over roughly two decades of vertical wavenumbers directly from the measured along-beam velocity variability (“direct” method).^{6} Averaging the spectra over an ensemble of successive samples provides the necessary statistical stability of the spectral estimate without compromising wavenumber resolution. These spectra are not affected by the float motion in any way and are equivalent to an Eulerian view of the turbulence. Vertical wavenumber spectra, calculated from 120-sample 1-min bursts of HR observations (Fig. 7a) generally followed the *k*^{−5/3} slope of the canonical Kolmogorov spectra of turbulence in the inertial subrange, with the levels consistent with the turbulent kinetic energy (TKE) dissipation rate *ε* = 10^{−8}–10^{−4} W kg^{−1}. Spectral levels at higher wavenumbers (*k*_{z }5–10 cpm) show pronounced flattening as a result of unresolved turbulent velocity variance within the sampling volume (Veron and Melville 1999; see also section 3c).

Frequency spectra of vertical velocity observations from a float in fully Lagrangian drift mode (Fig. 7b) represent the Eulerian–Lagrangian (E-L) view of an observer moving with the large turbulent eddies. As shown by Fung et al. (1992), E-L frequency spectra in fully developed homogeneous isotropic turbulence are expected to have an inertial subrange with the same *f *^{−2} spectral slope as a Lagrangian frequency spectrum but at somewhat lower level for a given ε; however, appropriate Kolmogorov constants for either case are not well established (Lien and D’Asaro 2002). Observed frequency spectra are generally consistent with the theoretical E-L spectra corresponding to ε = 10^{−8}…10^{−4} W kg^{−1}_{}. At lower frequencies, the E-L spectra derived from ADCP observations merge with the Lagrangian spectra derived from the float pressure record (Fig. 7b), as expected.

It is obviously possible to fit the “universal” spectra to those observed (either in wavenumber or frequency) and estimate TKE dissipation rates over the period of the float deployment. This, however, is beyond the scope of the present study. It remains to be determined whether the deviations of the observed spectra from the *k*^{−5/3} and *f*^{−2} slopes indicate meaningful departures from the canonical isotropic turbulence model as a result of the presence of stratification, anisotropy, or nonequilibrium forcing. Here, we merely demonstrate the potential of the float-mounted HR ADCP observations of ocean turbulence that may lead to insights beyond the canonical turbulence models. The internal self-consistency of spatial and temporal spectra is demonstrated in the next section.

### b. Turbulence profiling

A float is typically scheduled to pause its Lagrangian drift on a regular basis to execute a vertical profile and surface for communication. During these maneuvers, the float-mounted HR ADCP becomes a turbulence profiler, conceptually similar to those operating on McLane moored profilers (Rusello et al. 2011) but without the complicating effects of the mooring line motion and with a different beam configuration.

Each of the range cells of the HR ADCP vertical beam sweeps through the turbulent velocity field during the profile. The time series of HR vertical velocity measurements *w*(*t*) can therefore be converted to spatial transects *w*(*z*) and interpreted with the help of the Taylor’s frozen turbulence hypothesis, as is customary for microstructure turbulence profilers (Lueck et al. 2002); streamwise wavenumber spectra can then be estimated (“indirect” method). The streamwise coordinate *z* is the length of the path of the sampling volume through the “frozen” turbulence field, which is obtained, to the first order, by integration of measured relative vertical velocity. Even though the pathlength *z* can also be calculated from the pressure time series, it would not capture the relative advection of the small-scale turbulence through the sampling volume by the large turbulent eddies or, in presence of stratification, by internal waves. To calculate vertical wavenumber spectra, measured velocities *w*(*z*) are interpolated onto a regular grid, with the spacing δ*z* equivalent to the mean first difference of *z*(*t*). With this procedure, each of the ADCP range cells produces a quasi-independent estimate of the vertical wavenumber spectrum of vertical velocity; these estimates are subsequently averaged to improve their stability.

The comparison of the direct and indirect methods of estimation of vertical wavenumber spectra of vertical turbulent velocities in profiling mode shows good agreement in the overlapping wavenumber range (Fig. 8). The indirect method extends the spectra to lower wavenumbers, although the applicability of the frozen turbulence hypothesis degrades as the scales increase (Lumley 1965). The typical vertical velocity of the float during the profile was 20–30 cm s^{−1}. At the 2-Hz sampling rate, this produced an equivalent average vertical resolution δ*z* = 10–15 cm, which was somewhat coarser than the typical HR range resolution (2–4 cm). Reducing the profiling speed and increasing the sampling rate (up to the maximum 16 Hz with a single beam) would improve the special resolution of the indirect method and expand its wavenumber range. It would, however, run into the same noise floor limitation as the direct method (see the previous section).

Even though indirect estimation of wavenumber spectra of turbulent velocities via frozen turbulence assumption could be formally performed in a fully Lagrangian drift mode (as done by Lien and D’Asaro 2006), a few problems may require special consideration: Unlike an acoustic Doppler velocimeter employed by Lien and D’Asaro (2006), ADCP cannot measure full 3D velocity within a given volume, relying instead on an approximation based on radial velocities along the three–five slanted beams diverging in space. When the measurement platform moves with the large turbulent eddies (as does an LF in drift mode), relative velocities advecting turbulence past the float are themselves turbulent, varying substantially on the scales comparable to the ADCP beam separation. This violates the core ADCP assumption of horizontal homogeneity of the flow that is necessary to estimate 3D velocity vectors (Gordon 1996). Furthermore, in absence of dominant motion, relative advection trajectories form convoluted and potentially looping paths, complicating the application of the frozen turbulence hypothesis. Finally, in order to circumvent the difference between the spectra of streamwise and spanwise turbulent velocity components, full (3D) turbulent kinetic energy spectra need to be considered (Lien and D’Asaro 2006) which is, again, not possible without the point measurements of 3D velocity. These problems are avoided during profiling because it induces a steady relative flow of turbulence past the ADCP, which is also conveniently aligned with the vertical ADCP beam and therefore can be measured accurately.

### c. Uncertainty statistics of HR ADCP velocity measurements

Uncertainty of HR ADCP velocity measurements is determined by the instrument’s ability to resolve the phase change between the pulses, which in turn depends on the sharpness of the signal autocorrelation function or, equivalently, the Doppler spectrum width (Miller and Rochwarger 1972; Lhermitte and Serafin 1984). The amplitude of pulse correlation *R*(*τ*) = |*C*(*τ*)| is therefore routinely calculated along with its phase and is saved as a part of ADCP diagnostics. With *C*(*τ*) given by (3), *R*(*τ*) represents a maximum-likelihood estimator for the true signal autocorrelation function amplitude at lag τ (Miller and Rochwarger 1972).

Broadening of the Doppler spectrum (or pulse decorrelation) can be caused by a number of extraneous factors, such as changing scatterer composition and sonar beam divergence, as well as unresolved laminar shear and turbulent velocity fluctuations within the sampling volume (Cabrera et al. 1987; Lhermitte and Lemmin 1994; Zedel 2014). The HR velocity noise level is therefore highly variable and difficult to predict (in contrast to LR measurements, where the spectral width and therefore the measured velocity uncertainty are primarily defined by the pulse duration).

Unresolved turbulence is a major contributor to the measured velocity uncertainty (except perhaps in low-turbulence, low-scattering environments). With the assumption of Kolmogorov similarity hypothesis, the unresolved turbulent velocity variance within the range cell can be expressed as

where = 0.707 is the transverse Kolmogorov constant (Yeung and Zhou 1997), *η* is the Kolmogorov length microscale, and *k*_{d} = *d*^{−1} is the wavenumber corresponding to the characteristic sampling volume size *d* ≫ *η* [wavenumbers are expressed in cycles per minute (cpm) here].^{7} The width of the sampling volume increases with range as a result of the finite beamwidth. Thus, the effective range resolution of the HR velocity measurements may be reduced compared to the nominal range cell size. With the setup used in the Colvos Passage deployment (Table 1), the beamwidth exceeded the cell size for ranges greater than 0.4 m, reaching 20 cm for the last cell. This loss of resolution is evident in the flattening of the TKE spectra at wavenumbers *k*_{z }5−10 cpm roughly corresponding to the inverse average beamwidth (Fig. 7a). With the assumption of an idealized Gaussian Doppler spectrum, the corresponding autocorrelation amplitude is given by

(Zedel et al. 1996). Combined with (6), this gives an estimate for the Doppler signal correlation in presence of turbulence,

Decorrelation of the signal with the increasing cell size is illustrated in Fig. 9 showing mean profiles of correlation for 5-min segments of the Colvos Passage record (same as used for calculating frequency spectra shown in Fig. 7). As expected from (8), stronger decorrelation was observed at higher turbulence levels. Deviations of the observed correlation amplitude profiles from those predicted by (8) may be explained by the spatial variability of turbulence or other factors causing the Doppler spectrum broadening. Low-correlation spikes in the near-field (range < 0.5 m) may be indicative of signal contamination by acoustic scattering from the float superstructure or by a hydrodynamic wake of the float. It is not clear whether there is a corresponding effect of either factor on the measured velocity; visual inspection of the velocity fields (Figs. 6a–e) is inconclusive, requiring further investigation.

To clarify the relationship between the measured velocity variance and the signal correlation for different Nortek Signature1000 settings, a set of tests with varying cell size (*h* = 2−20 cm), pulse lag (τ*c*/2 = 2−8 m), and signal bandwidth (Δ*F*/*F*_{0} = 6.25%−50%) were run in a quiescent natural channel (Fig. 10). The test results are presented in terms of the Doppler phase variance, which is related to measured velocity variance via

As expected, the tests show a general decrease of measurement uncertainty with increasing correlation and the number of subsamples. The observed phase variance appears to flatten at a minimum value of about 4 × 10^{−4} rad, possibly indicating the level of instrumental (numerical) noise of signal processing.

A number of theories that relate observed pulse-to-pulse correlation *R* to Doppler phase variance var(*ϕ*) have been developed in the Doppler radar and sonar literature for various signal configurations (Miller and Rochwarger 1972; Zrnić 1977; Lee et al. 1994). Some insights could be obtained from Monte Carlo simulations as well (Dillon et al. 2012).^{8} A major complication is that both *R* and *ϕ* are nonlinear estimators that produce non-Gaussian and, in the case of *R*, highly skewed distributions even when the measured velocity fluctuations are Gaussian (Lee et al. 1994). Because of this, using (7), which relates the true values of variance and autocorrelation amplitude, would not be appropriate.

Various approximations tend to converge for high values of correlations and the number of subsamples, as the distributions of *R* and *ϕ* become asymptotically Gaussian (Dillon et al. 2012). Miller and Rochwarger (1972) used the central moment expansion to develop a simple relationship for this limiting case,

that is further simplified here with the assumption of a high signal-to-noise ratio. Equation (10) describes the general tendencies observed in testing (Fig. 10), except for the cases with a low number of subsamples and low correlations, where it overestimates the phase variance by a factor of 3−4. Such disagreement is to be expected, since the assumptions of the central moment expansion are formally violated in these cases (Dillon et al. 2012). Moreover, for a low number of subsamples (*M*10), the statistical relationship between *R* and var(*ϕ*) retains a dependence on the underlying distribution of the true signal correlation , which is generally unknown. Therefore, *R* cannot be considered a universal predictor for var(*ϕ*) for a low number of subsamples even in principle. Nonetheless, the underlying statistics of signal correlation is expected to vary slowly in any given deployment, so the tendencies presented in Fig. 10 can be used as guidelines for anticipating the effect of the HR ADCP settings on the uncertainty of the HR velocity measurements. Combining (9)–(10), we obtain

Thus, for a given distribution of signal correlation, the noise variance can be expected to be inversely proportional to bandwidth Δ*F*, cell size *h*, and the square of the pulse lag τ. It is worth reiterating that while (11) allows a posteriori estimation of measurement uncertainty, it does not fully characterize the dependence of this uncertainty on the instrument setup, since the signal correlation itself may depend on the setup parameters, particularly cell size and bandwidth. Therefore, theoretical (Cabrera et al. 1987) and numerical (Zedel et al. 1996) predictions for sonar performance still need to be developed further.

## 4. Conclusions and recommendations

A combination of Nortek Signature’s MCPC signal processing and 2D phase unwrapping in postprocessing allows HR Doppler velocity measurements to be pushed to a new level of range and accuracy. It should be kept in mind that the temporal frequency of the sampling needs to be sufficiently high to resolve temporal changes of the velocity field, which is a prerequisite for the advanced 2D phase unwrapping techniques [see section 2c(4)]. This also puts constraints on the uncontrolled platform accelerations that translate to accelerations of measured relative velocities. Platform acceleration is typically not a problem for Lagrangian floats, which generally follow the lateral motion of the surrounding water and therefore have negligible relative horizontal accelerations; relative vertical accelerations during profiling are somewhat larger (up to 1 mm s^{−2} depending on the drogue design and operation), but they are still well resolved with 1-Hz HR sampling with typical ambiguity velocity of 4 cm s^{−1}. In one of our test deployments, however, the float was loosely tethered to a small surface buoy. The tether allowed the wave motions to be transferred to the float, causing periodic accelerations with up to 0.1 m s^{−2} amplitude. These accelerations were not well resolved with 2-Hz HR sampling with 3.5 cm s^{−1} velocity range, causing the phase unwrapping to fail. A newer version of Signature AD2CP will include an optional high-accuracy attitude and heading reference system (AHRS), which would replace the standard accelerometer–magnetometer combo and could potentially allow correction of platform motion.

## Acknowledgments

We thank Nortek AS and Atle Lohrmann for the opportunity to participate in beta testing of the Signature1000 HR mode. Fruitful discussions with Jim Thomson, Len Zedel, and Robert Pinkel were also extremely helpful. We are indebted to the MLF engineering team—Mike Kenney, Mike Ohmart, and Nick Michel-Hart—for making it all work. We also thank two anonymous reviewers for their thorough reading of our manuscript and thoughtful comments. The Colvos Passage deployment dataset is available on request from the authors. This work was supported, in part, by NASA SPURS Award NNX15AG63G, ONR LCDRI Awards N00014-14-1-0732 and N00014-14-1-0666, and a grant from the Gulf of Mexico Research Initiative (SA-15-15).

### APPENDIX

#### Matrix Formulation of Phase Correction Problem

We assume that a velocity profile over a range of *n* cells **v ***=* [*υ*_{i}], *i* = 1, ..., *n* obtained using EVR is mostly correct, with the exception of a subset of cells *s*_{j} = {*s*_{1}, …, *s*_{m}}, *m* ≪ *n,* identified as erroneous phase resolution spikes (or otherwise suspect). For each of these suspect cells, we seek a subtractive correction δ*υ*_{j}. Since we anticipate that the spikes are due to the erroneous phase ambiguity resolution, we constrain the corrections to be equivalent to an integer number of phase wraps , that is, δ*υ*_{j }*=* 2*k*_{j}*V*_{R}. The corrected velocity profile can be expressed as

where is a “spike embedding” operator—an *n* × *m* matrix constructed from an *n* × *n* identity matrix by selecting the columns corresponding to the suspect cells {*s*_{j}}.

We would like the corrected profile to be smooth and choose the sum of the squared finite differences of the velocity profile as a measure of roughness to be minimized. Here,

is the (*n*−1) × *n* first-difference matrix operator.

Thus, we obtain an ILS minimization problem,

where **v′= ****v** is the first difference of the velocity vector, = 2*V*_{R}(), and the matrix product () can be constructed by selecting the columns of corresponding to the suspect cells {*s*_{j}}.

A number of efficient algorithms of solving the ILS problem [(A1)] have been developed in various application fields (Xu et al. 2012). In our case of isolated ambiguity resolution errors, a suboptimal, yet very straightforward “rounding” solution (Teunissen 1999) appears to be sufficient: We replace the ILS problem [(A1)] with an ordinary least squares problem,

that is trivial to solve [e.g., using the left matrix division operator (“\”) in MATLAB] and then obtain an approximate ILS solution, where the overbar denotes rounding to the nearest integer. Even though this solution is more prone to errors than the optimal algorithms (Teunissen 1999), it works reasonably well when the phase resolution errors are sparse.

## REFERENCES

_{2}, NO

_{3}

^{−}, and POC through the evolution of a spring diatom bloom in the North Atlantic

*Proceedings of OCEANS’87: The Ocean—“An International Workplace*,” Vol. 1, IEEE, 93–97, https://doi.org/10.1109/OCEANS.1987.1160903.

*2015 IEEE/OES Eleventh Current, Waves and Turbulence Measurement (CWTM)*, St. Petersburg, FL, IEEE, 10 pp., https://doi.org/10.1109/CWTM.2015.7098120.

*Proceedings of OCEANS’92: Mastering the Oceans through Technology*, Vol. 2, IEEE,

*Proceedings of OCEANS’93: Engineering in Harmony with Ocean*, Vol. 2, IEEE, II394–II397, https://doi.org/10.1109/OCEANS.1993.326127.

*Two-Dimensional Phase Unwrapping: Theory, Algorithms and Software*. Wiley-Interscience, 512 pp.

*2015 IEEE/OES Eleventh Current, Waves and Turbulence Measurement (CWTM 2015)*, IEEE, 210–216, https://doi.org/10.1109/CWTM.2015.7098130.

*Proceedings of the IEEE/OES/CMTC Ninth Working Conference on Current Measurement Technology*, IEEE, 19–24, https://doi.org/10.1109/CCM.2008.4480837.

*2011 IEEE/OES/CWTM Tenth Working Conference on Current, Waves and Turbulence Measurements (CWTM)*, J. Rizoli White and A. J. Williams III, Eds., IEEE, 259–265, https://doi.org/10.1109/CWTM.2011.5759562.

*Proceedings of the IEEE/OES Eighth Working Conference on Current Measurement Technology*, J. Rizoli White and S. Anderson, Eds., IEEE,

*Fisheries Acoustics: Theory and Practice*. 2nd ed. Wiley, 256 pp.

*2015 IEEE*

*/OES Eleventh Current, Waves and Turbulence Measurement (CWTM)*, IEEE, 79–83, https://doi.org/10.1109/CWTM.2015.7098107.

*Proc. 2014 Oceans—St. John’s*, St. John’s, NL, Canada, 5 pp., https://doi.org/10.1109/OCEANS.2014.7003130.

## Footnotes

© 2018 American Meteorological Society. For information regarding reuse of this content and general copyright information, consult the AMS Copyright Policy (www.ametsoc.org/PUBSReuseLicenses).

^{1}

Note that a number of different normalizations of correlation amplitude have been used historically (Lhermitte and Serafin 1984; Zedel et al. 1996).

^{4}

We use the term *fish* to encompass all acoustic targets that move (swim) relative to the water. This would also include zooplankton, oil or gas bubbles, buoyant marine debris, research and fishing gear, etc.

^{5}

Universal threshold is an easily computable (although generous) upper bound on the expected maximum value among *n* samples drawn from a standard normal distribution.

^{6}

Second-order spatial structure functions (SF) could be equivalently used to characterize the observed turbulence velocities (Rusello 2009; Thomson et al. 2015).

^{7}

This expression is likely implied by Cabrera et al. (1987), although they refer to an unpublished manuscript for actual derivation.

^{8}

But note that Dillon et al. (2012) discuss is in terms of the unknown “true” correlation coefficient rather than the sample correlation *R*.