Sonar observation of diffuse hydrothermal flows

Abstract

Previous studies have used multibeam sonar data to map diffuse hydrothermal flows by exploiting the decorrelation of successive seafloor echoes caused by temperature fluctuations. The present work extends this approach to quantify temperature fluctuations integrated along the sonar line of sight. The method is illustrated using data from the Cabled Observatory Vent Imaging Sonar deployed at the Grotto Vent complex on the Endeavour Segment of the Juan de Fuca Ridge. Inversion results are presented in the form of the structure function for integrated temperature fluctuations. A three‐parameter fit to the structure function provides time series encapsulating statistics of these fluctuations.

Introduction

The diffuse low‐temperature (tens of °C) flow of a hydrothermal vent field discharges from a variety of settings in the surroundings of high temperature (>250°C) sites including the sides of sulfide mounds, isolated cracks, lava tubes, and generalized areas of permeable seafloor [Bemis et al., , Figure ; Barreyre et al., ]. Although often unimpressive except for the biological cover, diffuse flow is an important part of the overall hydrothermal system, as the cumulative heat output through areas of a vent field covered by diffuse flow discharge may equal or exceed that of focused flow from associated smokers [Barreyre et al., ; Rona and Trivett, ; Schultz et al., ]. Direct estimates of diffuse heat transport have generally combined spot field measurements of temperature and discharge rates with mapping of the areal extent of diffuse discharge [e.g., Schultz et al., ; Barreyre et al., ; Rona et al., ], amongst a multitude of relevant studies. For example, in 1988, Schultz et al. [] used an electromagnetic‐based flowmeter and several thermistors on the sulfide mound called Peanut in the southern half of Main Endeavour Field (MEF) to simultaneously measure temperature and vertical velocity in diffuse discharge for about 2 months. Average diffuse discharge temperatures (7–13°C) and flow rates (0.07–0.15 m/s) observed over a 0.07 m² area beneath the sensor imply a diffuse heat flux of 2.91±0.23 MW/m² which was extrapolated to suggest total diffuse heat transfer of approximately 60 MW for Peanut (20 m² summit) and 9000 MW for MEF (3227 m² venting area). More recently, acoustic estimates of areal extent were combined with limited temperature and flow rate measurements to estimate heat transport by diffuse discharge at North Tower of the sulfide mound Grotto in MEF as 33–380 MW or a heat flux density of 0.33–3.8 MW/m² over the 100 m² summit of the North Tower [Rona et al., ]. These studies (as well as many others) document considerable spatial and temporal variability in diffuse discharge, observing both tidal modulation and localized episodic excursions in temperatures, [Barreyre et al., , ; Rona et al., ; Schultz et al., ] and distinct styles of discharge between mounds, cracks, tubes, and mats [e.g., Barreyre et al., ; Meyer et al., ; Sohn, ]. Thus, simple extrapolation over a mapped area is rendered inaccurate by (1) local temporal and spatial variations in temperature and discharge rates and (2) systematic variations related to distinct settings (i.e., cracks and mounds).

Fig. 1. The incidence (red) of correlation below an arbitrary threshold of 0.95 is shown for those locations visible (red and green) to COVIS on a bathymetric map of Grotto mound. This figure is an average over data for April 2015. Areas of high decorrelation (red) are generally located in areas geologically likely to have diffuse discharge. COVIS's location is shown on the map for reference (yellow diamond). The interval between bathymetric contour lines is 1 m, and COVIS is at a low point so that any contour crossed moving away from COVIS indicates a decrease in water depth. Inversions were performed at a and c on a slope several meters above COVIS. Inversions were also performed and averaged over the region between a and b. Bathymetry from Clague et al. [].

Multibeam sonar observations of hydrothermal flows began with 3‐D imaging of focused flows (“black smokers”) [Rona et al., ] and 2‐D images of diffuse flows [Rona et al., ]. Methods have been subsequently developed to obtain quantitative information on focused flows, including interaction with ambient current [Rona et al., ], volume flux [Xu et al., ], and heat flux [Xu et al., ]. The present effort seeks to extend this work, following Ivakin et al. [], providing an inversion method for quantitative observation of diffuse flows. The inversion method is applied to data obtained by the Cabled Observatory Vent Imaging Sonar (COVIS) at the Endeavour node of the NEPTUNE observatory operated by Ocean Networks Canada (ONC). The deployment of COVIS at this site at a depth of 2200 m over the period 29 September 2010 to 10 September 2015, with a 1 year interruption, provided a time series of 4 year length. In section 2, the sonar data are described and the inversion algorithm used to estimate geophysically interesting parameters is developed. Section 3 gives inversion results obtained by applying the algorithm to the data, and section 4 offers conclusions and suggestions for future work.

Methods

Sonar System

COVIS is a high‐frequency (200 and 400 kHz) sonar system, based on the Reson 7125 multibeam sonar. The source and receiver transducers are mounted on a tripod 4 m high on a novel triaxial rotator. COVIS provides near‐real‐time images and monitoring of hydrothermal flow [Bemis et al., ]. COVIS was positioned and oriented to scan in azimuth and elevation a volume encompassing the Grotto Vent cluster and associated buoyant plumes for maximum ranges of about 75 m. The raw acoustic data were transmitted via the cable to a shore station and to the University of Victoria where raw data and sonar metadata are archived as compressed files by the ONC's Data Management and Archival System (DMAS). DMAS provides processing programs and processed data to the community.

COVIS uses a wide‐vertical angle transmitter (22° 3 dB full beam width at 200 kHz). The transmitter directivity is wide in azimuth (130° 3 dB full beam width), and this sector is subdivided by multiple, narrow receiver beams. The receiver array, with digital beamforming, provides 128 beams having 3 dB azimuthal width of 1.0° to cover a sector of width 128°. In acquiring diffuse flow data, COVIS functions in much the same way as a conventional sonar, with echo time series obtained for each of the 128 receive beams. The echoes are due to scattering by the seafloor, and a two‐dimensional map of the echo intensity would yield an image of the scattering level of the seafloor, with azimuth determined by beam number and range determined by time of flight and the speed of sound. As will be noted shortly, the rugged topography at Grotto causes some portions of the seafloor to be shadowed acoustically. Although the sonar echoes are due to scattering by the seafloor, the algorithm to be detailed quantifies the changes in seafloor echoes due to sound speed changes. This “acoustic scintillation” is slight, but it is possible, using the known relation between water temperature and sound speed, to remotely sense temperature change.

The COVIS team has adopted the term “sweep” to denote a complete set of pings required to carry out a particular type of measurement. Unlike plume imaging sweeps, which require stepping of the transducer pointing direction in elevation, the diffuse flow sweeps use a single elevation setting, with the sonar pointed horizontally. In this mode, azimuthal coverage is provided by the 128 digitally formed beams covering a 128° sector, and range is determined as in most sonars by echo time of flight out to ranges of about 50 m. As in earlier acoustic mapping efforts [Rona et al., ], COVIS exploits the expected constancy of ping returns from a solid object (like the seafloor) to detect the phase changes caused by travel time variations due to the variable water temperature between COVIS and the seafloor [Rona and Jones, ]. As there is time jitter between the instant of transmission and beginning of digitization, the transmitted signal is digitized along with the echo data and used to align the phase of successive echo time series. COVIS uses digital beamforming in the horizontal with a Hamming window to provide azimuthal resolution of 1°. At a typical range of 20 m, this corresponds to cross‐range resolution of 0.35 m. The processing algorithm to be defined sets the range resolution at 0.4 m. Improvement of these resolution values would require narrower beam width and wider bandwidth, neither of which are available in the sonar employed in this work. In the diffuse flow mode of operation, in a single sweep, COVIS transmits a burst of 10 pings at a rate of 4 Hz, with sweeps repeated every 3 h. The transmitted pulse has an approximately rectangular envelope of width 0.3 ms and a nominal source level of 200 dB re 1 μPa@1m. In previous work the correlation lag variable was set to 1 s (e.g., ping No. 5 compared to ping No. 1), and Figure  shows a map (for 10 May 2015 at ~03:50 UTC) of the areas interpreted by this method as having diffuse flow, as well as the areas visible to COVIS. It should be understood that this is not a conventional sonar image in which image intensity would be determined by echo intensity. The black regions are those shadowed by topography, while the red and green regions are those visible to the sonar. The red regions show a “ping‐to‐ping correlation” (defined in the next section) of 0.95 or less and are interpreted as regions of diffuse flow. In normal operation, COVIS provides lags up to 2.25 s by comparing echoes from a single 10 ping sweep. To provide greater lags, special sweeps were performed 19–20 May 2015. In these sweeps, bursts of 1000 pings were transmitted at a rate of 1 Hz. Normal and special sweeps were performed every hour for a total of 24 normal‐special sweep pairs. The normal sweeps provide lags of 0.25 s to 2.25 s in 0.25 s steps, while the special sweeps provide lags of 1 s to 999 s in 1 s steps.

Inversion Algorithm

The basis of the inversion algorithm is that changes in sound speed cause changes in the sonar echo from the seafloor. These changes are quantified in terms of ping‐to‐ping correlation, and decorrelation (the departure of correlation from unity) is related to temperature change between pings, as there is a simple linear relationship between temperature change ΔT and sound speed change Δc. For this work, dc/dT = 3.5 ms⁻¹ K⁻¹ is used, appropriate for nominal temperature 10°C, salinity 35 practical salinity unit (psu), and depth 2200 m. Departure from the linearity in the sound speed‐temperature relation can be quantified by the second derivative, d²c/dT² =  − 0.04 ms⁻¹ K⁻². The resulting error in the assumed slope would be 0.4 ms⁻¹ K⁻¹ at 20°C. Salinity fluctuations have a negligible effect, as a salinity change of about 15 psu would be required to match the sound speed change caused by a (plausible) 5° temperature change, whereas the salinity measured by the Remote Access Sampler at a diffuse flow site within 5 m of the one examined here has a standard deviation of 0.58 psu over a 10 month period beginning July 2013 (D. A. Butterfield, personal communication, 2016). Current fluctuations have no effect on the sonar signal, as time shifts due to this mechanism are cancelled over the round‐trip acoustic path.

The inversion proceeds in three steps. The first step involves signal processing: a given sonar echo is compared with a later echo using a correlation algorithm (ping‐to‐ping correlation). Second, this correlation is related to an integral measure of temperature change (“path‐averaged temperature change,” to be defined) via an acoustic model. Finally, path‐averaged temperature change is expressed statistically in terms of a structure function.

The ping‐to‐ping correlation function is formed using two (complex) echo time series s₁(t) and s₂(t)from the same sonar pointing direction, transmitted at “true” times t₁ and t₂ = t₁ + τ and with “acoustic travel time” time arguments t measured from the center point of each transmission. True time can be regarded as the usual clock time (e.g., UTC), while acoustic travel time resets to zero with each sonar transmission. The Reson sonar provides travel time series in “baseband” form with magnitude and phase equal to the time‐varying envelope and phase of the real signal (the real signal is the real‐valued function oscillating at the carrier frequency (200 kHz) as it appears at the output terminals of the receiving transducer and the first stages of gain of the receiver electronics). The ping‐to‐ping correlation estimator is $$\tilde{\rho }=\frac{K_{12}}{\sqrt{K_{11}{K}_{22}}},$$ where $$K_{12}=\underset{-\infty }{\overset{\infty }{\int }}{s}_1\left(t\right)s_2^{*}\left(t\right)w\left(t-t_0\right)\mathrm{d}t,$$ and with K₁₁ and K₂₂ being given by similar integrals with the squared magnitude of each echo time series. If the two echo time series are identical, $\tilde{\rho }$ has a value of unity. If the two are substantially different, $\tilde{\rho }$ will have a magnitude much smaller than unity and may be complex. If the two time series are identical except for a phase difference, $\tilde{\rho }$ will have unit magnitude and phase equal to the phase difference. The integral in is over acoustic travel time with weighting by a window function w(t) whose position time, t₀, determines path length (“slant range” in sonar parlance), r = ct₀/2, where c is the speed of sound (about 1495 ms⁻¹ the case of interest). The window width is set by a trade‐off between range resolution and the desire to reduce statistical fluctuations. The integral is performed as the sum over samples in the time domain with a sampling frequency of 34 kHz.

Key to the inversion process is an acoustic model giving a relation between ping‐to‐ping correlation and sound speed change. As noted earlier, there is a linear relationship between sound speed change and temperature change. Accordingly, the discussion in the next several paragraphs will center on sound speed change rather than temperature change. The acoustic travel time for a sonar signal to travel to a point on the seafloor and back to the sonar is $$t_1=2\underset{0}{\overset{r}{\int }}\frac{\mathrm{d}l}{c_1}$$ for ping 1, with the integral being along the nearly straight path of length r shown in Figure . At a later time, the position‐dependent sound speed will have changed from c₁ to c₂ = c₁ + Δc, and the acoustic travel time for ping 2 will be altered to $$t_2=t_1\left(1-\alpha \right),$$ where the fractional change in acoustic travel time is $$\alpha \approx \frac{1}{\mathit{cr}}\underset{0}{\overset{r}{\int }}\mathrm{\Delta }c\mathrm{d}l.$$

In deriving this expression ray bending is neglected, based on ray tracing simulations showing that refractive effects are negligible in the circumstances of interest. The parameter c is a nominal, position‐independent sound speed, and $$\mathrm{\Delta }c=c_2-c_1$$ is dependent on position along the path and is assumed to have absolute value much smaller than c. Using the value dc/dT = 3.5 ms⁻¹ K⁻¹ given above, typical temperature changes of order 5 K cause sound speed changes of order 20 ms⁻¹. The model neglects dependence upon coordinates parallel to the seafloor within the sonar resolution cell; thus, Δc must be interpreted as an average over these coordinates.

Fig. 2. Geometry for sonar measurement of a sloping diffuse flow region. The line‐of‐sight path length r and layer thickness b along this line are indicated as are the unit vectors pointing toward the sonar and normal to the seafloor.

Next, consider how the time scaling caused by sound speed change affects the second echo signal s₂(t) relative to the first, s₁(t). Treating the sonar and environment as a time‐invariant system, one can write $$s_1\left(t\right)=\underset{-\infty }{\overset{\infty }{\int }}{s}_0\left(t^{\prime }\right)e^{i{\omega }_{0}t\prime }g\left(t-t^{\prime }\right)\mathrm{d}{t}^{\prime },$$ where g(t) is the “impulse response” of the combined sonar‐environment system at the time of the first transmission, s₀(t) is the envelope of the transmitted signal, and ω₀ is the center frequency of the sonar in rad/s. The impulse response incorporates the complexity of scattering by the seafloor as well as the effects of acoustic transmission and sonar properties. It will be possible to ignore these complications by making the assumption that, considering a formal ensemble of statistically equivalent seafloors, $$<g\left(t_1\right)g^{*}\left(t_2\right)>=G_0\left(t_1\right)\delta \left(t_1-t_2\right).$$

This assumption posits that the scattering interface has no regular structure and is frequently used in applications of the point‐scattering model [Faure, ]. As a counterexample, a seafloor whose roughness is due to regular ripples would not satisfy . Due to the time scaling shown in , the second echo is $$s_2\left(t\right)=\underset{-\infty }{\overset{\infty }{\int }}{s}_0\left(t^{\prime }\right)e^{i{\omega }_{0}t\prime }g\left(t-\mathit{\alpha t}-t^{\prime }\right)\mathrm{d}{t}^{\prime }.$$

This picture differs from that of Jackson and Dworski [] who assumed as an approximation that sound speed change simply stretches or shrinks the echo itself, when, in fact, only the impulse response is scaled.

It is not feasible to obtain a closed‐form model for the expected value of the estimator $\tilde{\rho }$, as it is an irrational functional of the two time series. Instead, one can readily obtain a formal result for the following expression $$\rho =\frac{\mid <K_{12}>\mid }{\sqrt{<K_{11}><K_{22}>}}.$$

Inversion will be performed by using as the estimate of correlation and comparing with , which is treated as the model. To obtain the required formal result for , first combine and – to obtain $$<K_{12}>=\underset{-\infty }{\overset{\infty }{\int }}\underset{-\infty }{\overset{\infty }{\int }}{s}_0\left(t_1\right)s_0^{*}\left(t_1-\mathit{\alpha t}\right)e^{i{\mathit{\alpha \omega }}_{0}t}{G}_0\left(t-t_1\right)w\left(t-t_0\right)\mathrm{d}{t}_1\mathrm{d}t.$$

Assuming that shifts due to time scaling are much smaller than the pulse length (αt much smaller than the pulse length) and that G₀(t) is slowly varying compared to both the window function w(t) and the transmitted envelope, can be approximated as $$<K_{12}>=e^{i{\mathit{\alpha \omega }}_0{t}_0}\underset{-\infty }{\overset{\infty }{\int }}\mid s_0\left(t_1\right)\mid {}^2\mathrm{d}{t}_1{G}_0\left(t_0\right)\underset{-\infty }{\overset{\infty }{\int }}{e}^{i{\mathit{\alpha \omega }}_{0}t}w\left(t\right)\mathrm{d}t.$$

Finally, can be written as follows: $$\rho =\frac{\mid \underset{-\infty }{\overset{\infty }{\int }}{e}^{i{\mathit{\alpha \omega }}_{0}t}w\left(t\right)\mathrm{d}t\mid }{\underset{-\infty }{\overset{\infty }{\int }}w\left(t\right)\mathrm{d}t}.$$

The inversion equates expression and the measured correlation to determine α. This yields path‐averaged sound speed change via , which is then divided by dc/dT to obtain path‐averaged temperature change. The model for correlation is proportional to the Fourier transform of the window function. In previous work (including that exhibited in Figure ), a rectangular window function was used, but this has the potential for ambiguity, as more than one value of the unknown α may correspond to a given value of correlation if the measured correlation is less than 0.217 (Figure ). Any one of several bell‐shaped window functions having low Fourier sidelobes would offer improvement. In the present work a Hamming window is used, and a fast Fourier transform yields the correlation function shown in Figure . It can be seen that ambiguity is not a problem for correlations larger than about 0.0075. This limit is not approached by the data employed in the present work. For example, in the processing to be described covering the year 2013, 0.3% of the ping‐to‐ping correlation values fell below 0.217, and only 0.02% fell below 0.0075.

Fig. 3. Model correlation function resulting from Hamming and rectangular analysis window functions of width 1 ms with sonar center frequency 200 kHz.

The path‐averaged sound speed change is $$\mathrm{\Delta }{c}_{\mathrm{pa}}=\mathit{\alpha c}=\frac{1}{r}\underset{0}{\overset{r}{\int }}\mathrm{\Delta }c\mathrm{d}l.$$

With reference to Figure , the largest Δc_pa that can be accommodated is about 15 m/s. Returning to the main theme of temperature change, it follows that the inversion method allows an estimate of path‐averaged temperature change, $$\mathrm{\Delta }{T}_{\mathrm{pa}}=\mathrm{\Delta }{c}_{\mathrm{pa}}/\left(\mathrm{d}c/\mathrm{d}T\right).$$

Given that path‐averaged temperature change, rather than path‐averaged temperature itself is estimated, a natural statistical characterization is provided by the structure function defined as $$S\left(\tau ,t_s\right)=<{\left[T_{\mathrm{pa}}\left(t_s+\tau \right)-T_{\mathrm{pa}}\left(t_s\right)\right]}^2>,$$ where t_s is a true‐time variable. The variable τ is a time difference, so identification with either acoustic travel time or true time is unnecessary. The brackets <> denote a formal average over a hypothetical infinite ensemble but will be implemented as a sample average. Note that the structure function is simply the second moment of temperature change occurring after an elapsed time τ, which will be referred to as “lag.” Although the structure function for path‐averaged temperature is estimated rather than the structure function for temperature itself, the following discussion often will refer to “temperature” rather than “path‐averaged temperature” for brevity.

The structure function is related to the more commonly used variance and correlation function as follows: $$S\left(\tau ,t_s\right)=2{\sigma }_T^2\left(t_s\right)\left[1-{\rho }_T\left(\tau ,t_s\right)\right],$$ where ${\sigma }_T^2\left(t_s\right)$ is the temperature variance and ρ_T(τ, t_s) is the correlation function, having a maximum magnitude of unity. Because of the linear relationship between temperature change and sound speed change, both temperature and sound speed will have the same correlation function, and the temperature variance will be the sound speed variance divided by a known factor (the square of dc/dT). It should be noted that the correlation function ρ_T(τ, t_s) is not the ping‐to‐ping correlation . To summarize, the ping‐to‐ping correlation function is used to determine path‐averaged temperature change, from which the structure function is obtained.

The structure function is dependent upon sonar geometry and so is not an intrinsic property of the observed flow. The following discussion outlines a method for obtaining an intrinsic measure. While is an integral along the sonar line of sight, if the thermal layer is, on average, stratified and does not extend to heights above the sonar, the integral can be taken along any straight path through the layer. In particular, the path may be vertical. Thus, the inversion output can be taken to be proportional to $$\mathrm{\Theta }=\underset{0}{\overset{\infty }{\int }}\mathrm{\Delta }T\mathrm{d}z,$$ where ΔT = Δc/(dc/dT). The advantage of converting path‐averaged temperature change to vertically integrated temperature change is that the latter is independent of sonar‐seafloor geometry, while the former is not. In order to understand the relationship between ΔT_pa and Θ, it is helpful to consider a thermal layer having thickness b measured along the line of sight as in Figure . For simplicity, ΔT will be assumed constant over the layer. Then ΔT_pa = bΔT/r. The thickness of the layer in the direction normal to the seafloor is equal to b multiplied by the sine of the angle that the line of sight makes with the seafloor. Thus, the perpendicular thickness is b n • e, where n is the unit vector normal to the seafloor and e is the unit vector pointing opposite the sonar line of sight. Finally, the thickness in the vertical direction is the perpendicular thickness divided by the cosine of the angle that the seafloor normal makes with the vertical, which is the vertical component of n, denoted n_z. These considerations give Θ = b ΔT n · e/n_z. Comparing the expressions for ΔT_pa and Θ, $$\mathrm{\Theta }=d_{\text{scale}}\mathrm{\Delta }{T}_{\mathrm{pa}},$$ where the “scale factor” has the dimensions of length and is $$d_{\text{scale}}=r\mathbf{n}·\mathbf{e}/n_z.$$

While the argument leading to assumed a uniform thermal layer, it is applicable to layers that are stratified in the direction normal to the seafloor. To convert the structure function to the structure function for vertically integrated temperature (an intrinsic property of the flow), it must be multiplied by $d_{\text{scale}}^2$.

Inversion Results

The inversion technique has been applied to locations within the field of view of COVIS near the Grotto Vent complex (Figure ). Data from the normal and special sweeps of 19–20 May 2015 were inverted using a Hamming window of width 1 ms, and the second moment of path‐averaged temperature change was obtained as an average over 24 sweeps spanning 24 h. This second moment is the desired structure function estimate and is displayed in Figure  for two locations, denoted a and c in Figure . Note that the structure function displayed is for path‐averaged temperature, not for vertically integrated temperature. The error bars in the figure represent the RMS error, assuming that all 24 samples are independent. Using the fact that the asymptote of the structure function is twice the variance, the standard deviation of path‐averaged temperature is about 1 K.

Fig. 4. Structure function for path‐averaged temperature fluctuations obtained using COVIS. The inversion was performed for 24 sweeps at locations a and c shown in Figure . The error bars represent ±RMS error in the estimates.

Thermistor measurements conducted during the September 2015 cruise during which COVIS was recovered indicate layer thicknesses less than 1 m, similar to results provided by Hautala et al. []. The thermistor data were obtained with a fast‐response SBE 39plus and show a standard deviation near the seafloor of about 2.5 K and a thickness of about 20 cm. These measurements were too sparse in time and space and too uncertain in location (about 2 m navigational error) to serve as ground truth.

An attempt to convert the measured structure function to that for vertically integrated temperature was unsuccessful. Considering possible navigational error, the scale factor was computed for a set of locations within 1.5 m of “c” in Figure  and found to range from 0.17 m to 29 m, a variation of over 2 orders of magnitude that is due to extreme variations in seafloor depth and slope. It would be better to perform inversions at sites with less severe bathymetry and to obtain temperature measurements over much wider intervals of time and space than in the present data set.

Inversion for the structure function over a wide range of lags using data from special runs does not exploit the lengthy deployment of COVIS, and it is a challenge to extract key parameters suitable for forming time series. In normal operation COVIS provides only a sparse set of lags, up to 2.25 s using single sweeps and multiples of 3 h by comparing different sweeps. As one method of extracting time series from the voluminous COVIS data set, consider the following three‐parameter fit to the structure function for path‐averaged temperature fluctuations: $$S\left(\tau \right)=2{\sigma }_T^2\left[1-e^{-{\left(\tau /{\tau }_0\right)}^{\mu }}\right],$$ where ${\sigma }_T^2$ is the variance of path‐averaged temperature fluctuations, τ₀ is a time scale, and μ is an exponent giving the slope at small lags in log‐log plots such as Figure . Comparing with , it can be seen that the correlation function for path‐averaged temperature fluctuations is $${\rho }_T\left(\tau ,t_s\right)=e^{-{\left(\tau /{\tau }_0\right)}^{\mu }}.$$

Note that the three fitting parameters are functions of true time t_s. Using data from “normal” sweeps, the slope is estimated by linear regression of the logarithm of the structure function with respect to the logarithm of the lag evaluated at the first several lags (up to 1.25 s). The variance is estimated by comparing pings separated by the 3 h sweep interval, and the time scale parameter is estimated as a function of the regression parameters and variance. Figure  shows time series for these parameters for a 1 year segment of COVIS data with the structure function averaged over an active region composed of an arc at fixed range extending between points a and b in Figure  using the Hamming analysis window defined previously and including eight contiguous sonar beams of width 1°. In order to restrict the analysis to a single ambient current regime, data were rejected if the coefficient of determination (R²) for the linear regression was less than 0.9, or ambient current was sufficient to cause at least 10° bending of the North Tower plume, as determined using COVIS imaging [Xu et al., ]. The exponent and standard deviation of path‐averaged temperature change appear to be rather stable as functions of time, while the time scale parameter is rather variable.

Fig. 5. Time series for structure function parameters for 2013.

It is hoped that the exponent will be useful in constraining the spectrum of diffuse flow temperature fluctuations via the Taylor hypothesis as in Duda and Trivett []. For fully developed turbulence in the inertial region, the structure function increases as the 5/3 power of lag [Ishimaru, ]. If the inversion were a point measurement rather than an average over the processing volume, one might expect the exponent in Figure  to be 5/3 ≃ 1.67. The inversion gives exponents with a mean of 1.1 and a standard deviation of 0.3. While the processing volume is finite, it is rather small, being limited in the azimuthal direction by sonar directivity to about 30 cm and in the downrange direction by the analysis window to about 40 cm. The vertical extent is equal to the downrange extent multiplied by n · e_r/n_z. This factor is not known to any useful accuracy, but it is likely smaller than unity, so the vertical dimension of the analysis volume is most likely less than 40 cm. If the Taylor hypothesis is to be invoked, one must ask whether heterogeneity is being advected by ambient current or by thermally induced vertical flow. The data selection criteria remove data with plume bending greater than 10°. Given North Tower plume rise velocities of order 20 cm/s at heights at which bending is apparent [Xu et al., ], this should exclude data obtained when ambient current is greater than about 3.5 cm/s. The measurements made at MEF by Schultz et al. [] showing diffuse flow velocity ranges from 7 to 15 cm/s. It seems likely, then, that thermally induced vertical flow is responsible for the observed fluctuations at small lag times. If this is the case, the time scale parameter may provide a measure of vertical flow rate. In conjunction with the standard deviation, this might provide the basis for another inversion step: estimation of heat flux.

Conclusions

This work is intended to provide quantification of diffuse flow activity using multibeam sonar data and is also intended as a first step toward estimation of heat flux from diffuse flow sites, analogous to previous work on focused flows [Xu et al., ]. One important element is to obtain better geometric data, and this requires less severe bathymetry than seen at the COVIS Grotto site. The largest hurdle is in constructing a geophysical model for diffuse flow suitable for use in an acoustic model. This may be difficult, in part, owing to the variety of diffuse flow types.

Ping‐to‐ping correlation of multibeam sonar data has been used to quantify diffuse flow at the Grotto Vent complex at the Endeavour Segment of the Juan de Fuca Ridge. An inversion algorithm converts ping‐to‐ping correlation to values for change in sound speed integrated along the sonar line of sight. Using the approximate linear proportionality between sound speed change and temperature change, one obtains the spatially integrated temperature change. The second moment of the temperature change is its structure function, and the structure function obtained from multibeam data may provide information on the turbulent structure of diffuse flow. To facilitate development of time series, a three‐parameter fit has been made to the structure function. An effort should be directed toward coupling this inversion technique with a geophysical model for diffuse flow, with the ultimate goal of remotely measuring heat flux. Future field work should invest considerable resources toward obtaining ground truth and would benefit from deployment of COVIS at a site with less severe bathymetry than Grotto.

Acknowledgments

The authors acknowledge the contributions of the late Peter Rona, who initiated the application of acoustics to the diffuse flow problem and led the collaboration for two decades. The APL‐UW team of Russ Light (leader), Vernon Miller, Michael Kenney, and Chris Jones developed COVIS and associated software with invaluable help from Eric Shug, Paul Jubinski, and Mike Mutschler of Reson. We thank Ocean Networks Canada (ONC) for their superlative assistance in deploying COVIS and for dive time and staff assistance in subsequent years. The support provided by the ROPOS and Jason ROV teams was essential to the success of this effort. Data acquired by COVIS are available on the ONC DMAS.