Synthetic aperture sonar (sas) systems are designed to observe stationary scatterers located near the sediment interface. Less commonly, a sas system may be used to observe scattering features located above the sonar in the water column. The Undersea Remote Sensing (USRS) project, sponsored by the Office of Naval Research, was a collaborative Directed Research Initiative (DRI) focused on studying dynamic estuarine water column features. During the USRS DRI, researchers from multiple institutions gathered to observe tidal features at various estuaries along the coast of the United States using both in situ and remote sensing techniques, including sas. The first studied estuary was the mouth of the Connecticut River (CTR). Data captured by a sas system deployed during a tidal event were post-processed to create three-dimensional observations of the structure of the leading edge of the CTR's ebb plume front. From these observations, lobed structures similar in scale to previously reported instabilities are revealed, with the present observations providing additional insight regarding the structure of the bubble distribution behind the front. Velocity estimates of plume features were also determined from sas data and shown to compare favorably with concurrent marine radar estimates.

From 2017 to 2021 the Undersea Remote Sensing Directed Research Initiative (USRS DRI), sponsored by the Office of Naval Research (ONR), funded the deployment of an autonomous underwater vehicle (AUV) mounted synthetic aperture sonar (sas) system in three different estuaries [Connecticut River (CTR) James River, and Mobile Bay]. The purpose of the deployments was to ascertain how the tidal features and fresh and salt-water interactions in these environments affect sonar imaging—and to determine if sas systems could sense hydrodynamic features in a way that may provide oceanographic insight. Thus, sas derived inferences associated with the spatial scales of plume features and propagation speeds are compared to estimates determined from other in situ and remote sensing techniques. In this context, favorable comparisons to widely adopted measurement techniques reveal the potential of sas systems as an additional tool in physical oceanographic studies.

The usage of a synthetic aperture system for making observations of mid-to-upper water column hydrodynamic features diverges significantly from most previous sas applications. Synthetic aperture sonar imaging exploits the motion of a transmitter and hydrophone array along a trajectory to create, via coherent post-processing, an artificially long “synthetic” aperture.1,2 This synthetic aperture can be used to create imagery with resolution that is, allowing for certain caveats, both range and frequency independent. This functionality, however, is closely linked to stationarity and medium homogeneity assumptions that are violated by the phenomena linked to estuarine hydrodynamics that this study wishes to observe. Furthermore, the signal-to-noise ratio (SNR) limiting effects of multipath in shallow water environments have had a major impact on the design of sas arrays,3–5 and most systems, including the one used in the present study, are specifically designed to attenuate the type of acoustic backscattering the present study wishes to observe, namely, backscattering originating from features located in the mid-to-upper water column regions.

If these challenges can be overcome, however, the sas system deployed during the USRS DRI has much to offer: it has a planar array that can be used to estimate azimuthal and vertical angle-of-arrival (AOA) with sufficient accuracy to create meaningful three-dimensional (3D) “point clouds” (xyz scatter plots) from each ping. The broad field of view allows multiple observations to be made of the same locations, potentially allowing the evolution of hydrodynamic features to be observed. In subsequent sections of this paper, the method used for generating point clouds from the deployed sas system will be described, along with a method by which the velocities of point clouds corresponding to hydrodynamic features were determined. The results of applying the developed methods to observe the features and velocities of two different tidal plume fronts encountered during a survey operation at the mouth of the CTR will be discussed and compared to concurrent measurements made via a land-based radar.

Many synthetic aperture systems use a pair of horizontal real-aperture arrays with a vertical offset to make vertical (elevation) angle-of-arrival (AOA) estimates for incoming signals. Examples include the HiSAS 10326,7 and Kraken MinSAS.8 The vertical baseline between arrays creates a time-difference of arrival for backscattering echoes that can be used to create an interferogram containing bathymetric information. One of the most common difficulties faced by interferometric processing is the resolving of angle ambiguities caused by phase wraps. Much literature has been devoted to the topic of interferogram phase unwrapping, with most papers being published in the related field of synthetic aperture radar,9–11 though literature specifically directed toward unwrapping interferograms created by synthetic aperture sonar systems also exists.12–16 Sonar systems often have a wider fractional bandwidth than radar systems, and several resources discussing various methods for exploiting bandwidth to resolve phase wrap ambiguities are available.17–20 

Unlike the more common two-stave interferometric systems, the sas system used in the present study has a planar array with more than two elements in the vertical dimension. Nevertheless, the spacing of the elements in the vertical dimension is still multiple wavelengths, and phase ambiguities must be dealt with for AOA estimation. Furthermore, the fractional bandwidth of the system (i.e., the signal bandwidth divided by the center frequency) is just at the edge of being sufficiently broad to leverage bandwidth for ambiguity resolution using the techniques described in Saebo18 or Synnes et al.:20 ambiguity resolution can be quite good when the SNR is very high but tends to wain when the SNR drops. SNR can be increased at the cost of range-resolution using the cross correlation strategy described in Saebo17 or Cook.19 However, this strategy is ordinarily reserved for stave pairs, whereas the current system employs a vertical array and data is available from a larger number of staves. To leverage the SNR-enhancing properties of correlation but retain compatibility with the vertical array of the sas system, a cross correlation beamforming framework21 was adopted for AOA estimation.

Consider a hypothetical sonar system analogous to the deployed sonar, with a center frequency f0 of 150 kHz, a bandwidth fbw of 45 kHz, and a planar array with five elements in the vertical direction. The vertical width wrcv of the elements is 3λ0, where λ0 is the wavelength (0.01 m) at the center frequency f0 given a sound-speed of c = 1500 m/s. The vertical spacing of the elements is 1.5wrcv. Furthermore, assume the system has a rectangular transmitter wtrn = 4.5 λ0 tall, and that the transmitter and receiver array are both oriented with a depression angle θ dep = –5 degrees relative to horizontal. To illustrate the correlation beamforming process, a data collection from this system was simulated for a sonar six meters above a distribution of random scatterers located along an interface, as shown in Fig. 1.

FIG. 1.

(Color online) The geometry of the sonar simulation showing the sediment (black line at −10 m), the sonar ( –4 m), and an illustration of the aggregate beam pattern formed by the transmitter and one of the array elements (red curve). The simulated sonar has a −5 degree depression angle which steers the mainlobe of the beam toward the sediment, and not much energy is directed above the sonar. In the simulation, the sediment is extended to 200 m in the horizontal direction.

FIG. 1.

(Color online) The geometry of the sonar simulation showing the sediment (black line at −10 m), the sonar ( –4 m), and an illustration of the aggregate beam pattern formed by the transmitter and one of the array elements (red curve). The simulated sonar has a −5 degree depression angle which steers the mainlobe of the beam toward the sediment, and not much energy is directed above the sonar. In the simulation, the sediment is extended to 200 m in the horizontal direction.

Close modal
The signal model for this simulation is similar to the point scatterer simulation tool described in Brown et al.22 i.e.,
P n ω = S ω m = 1 M b t θ g , m , ω b r θ g , m , ω sin θ g , m e i k R f , m + R b , m , n + ϕ m R f , m + R b , m , n + η ω ,
(1)
b t θ g , m , ω = sinc sin θ dep + θ g , m ω w trn / 2 c ,
(2)
b r θ g , m , ω = sinc sin θ dep + θ g , m ω w rcv / 2 c ,
(3)
where ω is angular frequency in radians per second, S ω is the frequency-domain signal amplitude, which in the present case is a Tukey-window with a Tukey parameter of 0.25, b t , r are the transmitter and receiver directivity weights calculated via Eqs. (2) and (3), θ g , m is the grazing angle relative to horizontal of the acoustic ray from the transmitter to the mth scatterer, R f , m is the distance between the transmitter and the mth scatterer, and Rb,m,n is the distance between the mth scatterer and the nth receive element in the vertical receive array. Last, ϕ m is a random phase rotation applied to each scatterer, k is the acoustic wavenumber ω/c, where c is modelled as a homogeneous 1500 m/s in the water column, and η ω is complex uniform Gaussian noise added to evaluate the performance of correlation beamforming as the SNR degrades. Note that the sine term in Eq. (1) is included to simulate a reflectivity that obeys Lambert's law,23 which is normally expressed as cosine of the incident angle, but the angular convention of Eq. (1) is with respect to the horizontal sediment plane (i.e., grazing) rather than vertical.
Applying Eqs. (1)–(3) results in a set of N frequency-domain signals representing the signal acquired at each of the vertical staves in the array. A basic multi-stave, narrow-band, far-field AOA estimate for each range-cell of the signals can be determined from a shear-average24 using Eqs. (4)–(6),
P 1 N r = F 1 P 1 N ω ,
(4)
ϕ intf r = arg 1 N 1 P n r P n + 1 * r ,
(5)
θ AOA r = sin 1 λ 0 ϕ intf r 3 π w rcv .
(6)
In Eq. (4), F 1 represents the inverse Fourier transform operator operating along the ω dimension, and Eq. (4) defines the transformation from a frequency domain signal to a signal that is a function of range, where the substitution r = (tc/2) has been made to the independent variable. In Eq. (5), the * symbol represents the complex conjugate, and ϕ intf r is the average interferometric phase between adjacent pairs across the array. Equation (6) converts the interferometric phase to AOA using the vertical array spacing value. The AOA estimated for the simulated data and the SNR of the signal vs range, determined by the added noise term η ω, are shown in Fig. 2.
FIG. 2.

(Top): The AOA estimated from the simulated data via Eqs. (3)–(5) and the SNR as a function of range (bottom). Angles are defined relative to the vector normal to the array, which is tilted −5 degrees relative to horizontal. The SNR curve in the lower plot was computed by determining the range-dependent root-mean-squares of the separate signal and noise terms in Eq. (1) using a 3 meter sliding window, and multiplying the log of the ratio by 20.

FIG. 2.

(Top): The AOA estimated from the simulated data via Eqs. (3)–(5) and the SNR as a function of range (bottom). Angles are defined relative to the vector normal to the array, which is tilted −5 degrees relative to horizontal. The SNR curve in the lower plot was computed by determining the range-dependent root-mean-squares of the separate signal and noise terms in Eq. (1) using a 3 meter sliding window, and multiplying the log of the ratio by 20.

Close modal
Combined with spreading loss the noise term η ω results in a range of SNR values, highlighting the rapid breakdown in performance as SNR deteriorates. The visible discontinuity in the AOA estimates near the 30 meter range mark is the result of phase wrapping in the value of ϕ intf r. This is not unexpected given the array geometry: at 4.5λ spacing, the AOA corresponding to the first phase wrap can be calculated by plugging π into Eq. (6), i.e.,
θ wrap = ± sin 1 λ 0 3 w rcv .
(7)
For the simulated system, θ wrap = ∼6.4 degrees, whereas the first null of the aggregate transmitter and beamformer beampattern lies at ∼13 degrees, placing the first wrap well within the mainlobe. As illustrated in Fig. 2, many near-range scatterers require unwrapping; this is due to the steep grazing angles at which backscattering occurs for near range objects below the sonar in this simulation. More generally, sonar systems will be operating in a variety of environments and it is not necessarily possible to guarantee that data from a particular range span will not require unwrapping.
Correlation beamforming, the technique used in the current study to improve the SNR, is performed by sectioning the data from each array into range-bins containing multiple samples and performing cross-correlations between each pairwise combination of staves for each range bin, including the autocorrelation terms, i.e.,
C i , j r k = P i r k P j r k , for all i , j 1 N ,
(8)
where ⋆ represents the cross correlation operation and r k represents the kth bin of range samples. Including autocorrelation, the total number of correlations in C i , j r k is N2; however, the total number of unique pairwise combinations is only (N + 1)N/2 = 15 for the simulated case of N = 5. Beamforming is performed by determining a set of candidate arrival angles θ cdt and, for the range r ¯ k corresponding to the center of the range bin r k, calculating and applying the delays associated with each angle, and summing across correlations. For the transmitter location (xtrn, ytrn) and N array locations (xrcv,1…N, yrcv,1…N),
x cdt , y cdt = x trn + r ¯ k cos θ cdt , y trn + r ¯ k sin θ cdt ,
(9)
τ i , j θ cdt = ( x cdt x rcv , i 2 + y cdt y rcv , i 2 x cdt x rcv , j 2 + y cdt y rcv , j 2 ) / c ,
(10)
C ¯ t k , θ cdt = F 1 i , j F 1 C i , j t k e j ω τ i , j θ cdt .
(11)
Note that in Eq. (11), the independent variable substitution t k= 2 r k / c has been made. C ¯ t k , θ cdt is a matrix of values that spans the length of the signal range bins in one dimension and the range of candidates in the other dimension. The AOA's associated with the dominant scatterers are the values of θ cdt that maximize C ¯ t k , θ cdt. In practice, three-point parabolic interpolation leveraging the peak and adjacent samples is used to estimate the correlation peak location to sub-sample precision. An example of this peak localization process can be found in Sec. III B and Eq. (5) from Michael et al.25 Fig. 3 shows an example of C ¯ t k , θ cdt, calculated via Eqs. (8)–(11) with range-bin lengths of 16 samples (24 cm), and an overlap factor of 50% between bins.
FIG. 3.

The correlation beamforming result for the simulated sonar data. The brightest curve is the mainlobe of the beamformer; adjacent curves lying above and below the mainlobe are ambiguity lobes. These are both less focused and attenuated relative to the mainlobe. The signal has been normalized by the maximum value at each range.

FIG. 3.

The correlation beamforming result for the simulated sonar data. The brightest curve is the mainlobe of the beamformer; adjacent curves lying above and below the mainlobe are ambiguity lobes. These are both less focused and attenuated relative to the mainlobe. The signal has been normalized by the maximum value at each range.

Close modal

In Fig. 3(a), a clear signal trace is visible both at close ranges and steep grazing angles, and out to the maximum range despite the low SNR. Determining the peak AOA value and comparing the resulting bathymetry estimates with those of the narrow-band technique, it is clear that the noise is reduced and the phase-wrap ambiguities have been resolved.

Note that the division of data into bins results in a reduction of range resolution. With a bin length of 16 samples (24 cm) and overlap of 50%, there is an eightfold reduction in the number of data-points in the correlation beamformer results shown in Fig. 4. The benefits of the beamformer, e.g., phase unwrapping and improvements in accuracy, outweigh the consequences of resolution reduction in the current study, in which hydrodynamic features having scales on the order of meters to tens of meters are being observed.

FIG. 4.

(Color online) The correlation-beamforming bathymetric height estimates (black dots) vs the narrow-band results (red dots). Note that, due to phase-wrap issues, many of the height values determined using the narrowband estimator are significantly wrong, and some even map above the air/water interface. The correlation beamformer, however, resolves the phase wrap issues and provides significant improvement both in accuracy and variance.

FIG. 4.

(Color online) The correlation-beamforming bathymetric height estimates (black dots) vs the narrow-band results (red dots). Note that, due to phase-wrap issues, many of the height values determined using the narrowband estimator are significantly wrong, and some even map above the air/water interface. The correlation beamformer, however, resolves the phase wrap issues and provides significant improvement both in accuracy and variance.

Close modal

Several additional notes need to be made regarding the usage of this technique with the deployed sas system: prior to performing the previously described AOA estimation procedure, the real-aperture array data from each horizontal row of elements were first beamformed into polar coordinates (range and azimuth). Following azimuthal beamforming, correlation beamforming is performed and peak values are extracted to create a point cloud in spherical coordinates: range, azimuth, and elevation. Finally, to geo-rectify the point cloud, the spherical coordinate representation is transformed into Cartesian coordinates and vehicle-to-earth coordinate transformations supplied by the gyroscope of the AUV’s inertial navigation system (INS) are applied to the point cloud data, alongside translations based on INS estimates for the AUV location. The rapid oscillations frequently experienced by the AUV, particularly roll, were found to cause serious point cloud misalignments if the earth-coordinate transformations were applied uniformly to each point within a ping. This was found to be due to the failure to account for the rapidly changing orientation of the array during the reception of ping echoes. This problem was resolved by sub-sampling the INS orientation state estimates to achieve unique corrections for the point cloud data corresponding to each range bin of the point cloud data. Figure 5 is an example of a point cloud generated using the deployed system, showing sediment backscattering (magenta, near −10 m), a fish cloud (yellow), and bubbles near the surface (green).

FIG. 5.

(Color online) An example of a single-ping point cloud created using correlation beamforming. The points are colorized based on height above sediment; at longer ranges, the dominant features are located near the water surface. Point luminosity is determined by backscattering strength, (bright luminosity = high backscattering, low luminosity = low backscattering), with dynamic range compression applied using the square-root function.

FIG. 5.

(Color online) An example of a single-ping point cloud created using correlation beamforming. The points are colorized based on height above sediment; at longer ranges, the dominant features are located near the water surface. Point luminosity is determined by backscattering strength, (bright luminosity = high backscattering, low luminosity = low backscattering), with dynamic range compression applied using the square-root function.

Close modal

Note that in Fig. 5, backscattering intensity information has been retained in the point cloud using luminosity. In this figure and subsequent point clouds shown in this paper, luminosity is determined by the magnitude of the correlation beamforming result after undergoing dynamic range compression via the square-root operator. This method of visualization was found to improve the visibility of features with high backscattering intensity in point cloud visualizations.

This section summarizes inferred plume properties based on sas-derived measurements and provides quantitative comparisons to similar measurements derived from other in situ and remote sensing techniques to demonstrate the effectiveness of sas systems in this context. To do this we provide relevant information about the CTR plume and supporting sampling methods including optical drone, X-band radar, echosounder, and current velocity and water property profiles. Combined, these are used to make direct comparisons of spatial-scales corresponding to plume features and estimates of propagation speeds both for the expanding plume front and hydrodynamic instabilities propagating along the front.

In 2017, the research teams associated with the USRS DRI performed a set of surveys and measurements of the ebb plume fronts that occur at the mouth of the CTR. The tidally pulsed freshwater outflow from the CTR initially discharges into a crossflow corresponding to the ebbing waters of Long Island Sound, meeting at a strongly convergent plume front a few meters thick.26–28 The location and water properties in the vicinity of the front vary temporally and with river discharge conditions, but the ebb plume front consistently provides conditions that result in elevated acoustic scattering from physical oceanographic features. First, the convergence of the river discharge with coastal water causes waves to steepen and break near the front, as described by Baschek et al.29 Measurements from the CTR and elsewhere show that this bubble entrainment occurs in the immediate vicinity of the front.29–31 Directly corresponding to the horizontal velocity convergence at the CTR ebb plume front are strong downwelling currents that have been measured to reach 0.20 to 0.25 m/s.28,32 This downward vertical flow can draw bubbles deeper into the water column, whereby advection driven by ambient currents can create diffuse bubble clouds extending tens of meters or more behind the front. Frontal features and processes that can entrain bubbles are associated with strong salinity gradients that can also be notable sources of backscattering at near-normal incidence.33 Bubbles, however, dominate the scattering in the near-field of the frontal features both in the CTR and in other estuarine environments.29,31

During the experiment, radar coverage of the CTR mouth was constantly maintained. The radar was stationed atop a telescoping tower located 1 km inside the river mouth. The CTR plume front was visible in radar images as it formed at the mouth and as it expanded westward. Figure 6 shows a map of the estuary as visible in Google Earth, and a corresponding example of a radar image. A tidal front is visible within the radar image. The 360-degree radar footprint updated at 0.8 Hz (48 RPM), and 80-s averages (i.e., 64-snapshot averages) were used to create de-speckled image frames. By measuring the temporal changes in the manually extracted location of the tidal front, the frontal velocities could be ascertained from the radar data.34 Uncertainty in these velocities varies with manual ridge-finding accuracy and with the smearing induced by time averaging. A conservative estimate of front position uncertainty is 13.6 m, half the 3 m antenna beam width at the measurement location. Smearing uncertainty must be calculated a posteriori, since it varies with front speed. Combined, the uncertainty can be estimated as 2 σ x front / Δ t, where σ x front = 13.6 + 80 V plume ( RADAR ) / 2 and Δ t = 900 s between front sample times.

FIG. 6.

(Color online) The mouth of CTR opening into the Long Island Sound. The figure on the left is a Satellite Image of the region surveyed by the DRI, gathered from Google Earth, (Google, map data: SIO, NOAA, U.S. Navy, NGA, GEBCO). The axes definitions are North = y = up, East = x= right. The figure on the right shows an 80-s average of image data captured by the Oregon State radar. The manifestation of a tidal front has been annotated in the radar image.

FIG. 6.

(Color online) The mouth of CTR opening into the Long Island Sound. The figure on the left is a Satellite Image of the region surveyed by the DRI, gathered from Google Earth, (Google, map data: SIO, NOAA, U.S. Navy, NGA, GEBCO). The axes definitions are North = y = up, East = x= right. The figure on the right shows an 80-s average of image data captured by the Oregon State radar. The manifestation of a tidal front has been annotated in the radar image.

Close modal

Representative examples of the physical oceanographic conditions at the CTR ebb plume front, in addition to upward-looking echosounder measurements from a REMUS 100 vehicle35 are shown in Fig. 7. This example, which is characteristic of the CTR front, shows that a large plume of bubbles is injected into the water column near the plume front and that these bubbles remain in the upper water column for 10 to 100+ meters behind the plume front.

FIG. 7.

(Color online) (A) An echogram from an upward-looking 200 kHz echosounder during a crossing of the plume front. The line around 5 m is associated with the transducer's blanking distance and the saturated data near the surface is from the processing sidelobes from the surface echo. The scattering is attributed to bubbles as increases in scattering are seen closer to the plume front. (B) Acoustic Doppler current profiler measurements of cross-front velocity overlaid with salinity contours. (C) The sound velocity profiles, calculated from a towed array of conductivity, temperature, and depth sensors, during a crossing of the plume front.

FIG. 7.

(Color online) (A) An echogram from an upward-looking 200 kHz echosounder during a crossing of the plume front. The line around 5 m is associated with the transducer's blanking distance and the saturated data near the surface is from the processing sidelobes from the surface echo. The scattering is attributed to bubbles as increases in scattering are seen closer to the plume front. (B) Acoustic Doppler current profiler measurements of cross-front velocity overlaid with salinity contours. (C) The sound velocity profiles, calculated from a towed array of conductivity, temperature, and depth sensors, during a crossing of the plume front.

Close modal

Relevant to AOA estimation, the sound velocity profile shown in Fig. 7 is not as extreme as might be expected given the strong salinity gradient. This is due to the thermal profile offsetting some of the effects of salinity: the freshwater layer at the time of measurement was warmer than the salt water in the Long Island Sound. The sound velocity that is present is upward refracting within the plume, which will bias the height estimation of scatterers within the plume to positions below their true location. This is because the current method used for AOA estimation assumes a homogenous medium and linear propagation. In an upward refracting medium the vertical angle of the ray between the sonar and scatterer will be at a minimum at the face of the sonar: the linear ray extending from the face of the sonar at the observed AOA lies beneath the actual ray, which is upward refracted during propagation inside the plume. The dense bubble cloud present at the interface between the plume and saltwater was generally observed to occlude scatterers deep within the plume, however, and the effect of this gradient on the geometry of the plume does not appear to strongly affect the reconstruction of observed plume features, which are mostly scattered from bubbles entrained near the plume boundaries.

Figure 8 shows a 3D point cloud generated from a single ping in which the trajectory of the sas system was oblique to the front and on course to intersect it. As in Fig. 5, the color scheme used to visualize the point cloud data encodes height information via hue and magnitude information using luminance. The color-bar on the right side of Fig. 8 shows the hue and luminosity values that may be used to describe height and magnitude. A ground-plane projected photograph of the front captured by an aerial drone during the same experiment but from a different tidal cycle is included to show the surface manifestation of the lobe-and-cleft features observed by the sonar.

FIG. 8.

(Color online) (A) A top-view of a point cloud generated at the mouth of the CTR showing a tidal plume front above sand ripple bed-forms. Heights are relative to the air/water interface. Near surface scatterers are colored green. Scatterers with positive height values, though few in number, are present and the reason for their presence is discussed in the text. The front has a lobed structure, and several prominent lobes have been identified in the image. Fifteen-meter lines in the northern and eastern directions are shown for scale and direction reference. (B) Shows a photograph of the front (from a different tidal cycle) captured by an aerial drone; the image has been projected onto the water surface, is 30 meters on edge, and is scaled for consistency with the sonar image. (C) Shows a ∼1-degree slice near the central angle of the sonar field-of-view, revealing the profile of a section of Lobe 1.

FIG. 8.

(Color online) (A) A top-view of a point cloud generated at the mouth of the CTR showing a tidal plume front above sand ripple bed-forms. Heights are relative to the air/water interface. Near surface scatterers are colored green. Scatterers with positive height values, though few in number, are present and the reason for their presence is discussed in the text. The front has a lobed structure, and several prominent lobes have been identified in the image. Fifteen-meter lines in the northern and eastern directions are shown for scale and direction reference. (B) Shows a photograph of the front (from a different tidal cycle) captured by an aerial drone; the image has been projected onto the water surface, is 30 meters on edge, and is scaled for consistency with the sonar image. (C) Shows a ∼1-degree slice near the central angle of the sonar field-of-view, revealing the profile of a section of Lobe 1.

Close modal

Several characteristics of the under-water features of the tidal plume front are visible in the point cloud of Fig. 8. For example, the plume exhibits a prominent lobed structure: two lobes, both ∼11 m long, are shown in the figure, while other similarly sized or smaller lobes are also present. Optical imagery captured from a drone that was deployed during the study [e.g., Fig. 8(B)] also shows the presence of lobes in wave breaking and bubble patterns on the water surface at the front. The presence of these structures was primarily attributed to instabilities arising from cross-front shear in the along-front velocity.36 The along-front length scales revealed in the sas derived point cloud are on the low end of those reported by Simpson et al.36 (roughly 20 m), but not inconsistent with the spatial and temporal variability observed using the drone.

The lobes appear to have bubble scattering signatures trailing behind them and these trailing signatures are not oriented perpendicular to the front, but rather trail at an angle in a south-easterly direction. This is consistent with the velocity vectors of the Long Island Sound current, which drives the south-easterly advection of the bubbles in the mixing layer at the base of the plume. The depth of the bubble plume as a function of range along the direction of the central beam is visible in Fig. 8(C). The distribution qualitatively appears to be linear near the leading edge indicating that a spatial descent rate may be calculated using a fit with a linear slope model. A plume front slope estimation was made by projecting bubbles along a portion of the front onto a vertical plane oriented perpendicular to the front, and culling points greater than 1 meter above the estimated water surface (note that points above z = 0 exist due to system biases, which will be discussed later), and less than 4 meters below the interface. The first 5 meters of the front were used for slope estimation. Slope estimation was performed using total-least squares via a singular value decomposition. The projected points used for the fit and the resulting best-fit line are shown in Fig. 9. The estimated slope is –0.62 meters in altitude for every meter behind the front, approximately −32 degrees, which is a value consistent with plume front structure observed using echo sounder measurements such as the one shown in Fig. 7(A), as well as others captured during the experiment.

FIG. 9.

Bubble “descent rate” behind the plume: the dots correspond to the green plume-front points shown in Fig. 8(B). The solid line is a best fit linear model through the points (R2 = 0.79) with a downward slope of −32 degrees.

FIG. 9.

Bubble “descent rate” behind the plume: the dots correspond to the green plume-front points shown in Fig. 8(B). The solid line is a best fit linear model through the points (R2 = 0.79) with a downward slope of −32 degrees.

Close modal

Figures 8 and 9 were generated using data from a single ping; data from multiple pings were also combined to achieve a larger scale picture of frontal features. Because the front is travelling, successful multi-ping mosaicking requires estimation of the relative shifts between point clouds prior to combination. Fortuitously, estimating the relative translations between pings can also be used to estimate the speed of travel of the features along the tidal front. The alignment process begins with combining point clouds of pings containing tidal plume front features in non-overlapping sets; in the present case, sets of five pings were combined to form dense, independent point cloud observations of the front. Following grouped point cloud creation, two-dimensional (2D) histograms of the x and y distributions of the points were created to make images of the front. To reject stationary sediment contributions and accentuate the frontal features only, the points having height values within the vertical span –4 m < z < 2 m relative to the air/water interface were used to create the histograms. The positive values in the span, corresponding to heights above the air/water interface, have been included due to the uncertainty in AOA measurements. Ambient acoustic noise, thermal noise in the sonar hardware, ray bending caused by sound-velocity profile inhomogeneities, vehicle depth, or residual roll-correction errors can erroneously position points above the interface. This effect is visible, for example, in Fig. 9, in which the leading edge of the front appears to be about half a meter above the interface. Culling points above the interface would result in the loss of the leading edge of the front, as indicated by the presence of coherent, high intensity scatterers, so the rejection threshold has been extended above the front by two meters. Another source of points that can appear near the air/water interface is multipath,3,4 particularly the bottom-surface path described in Belettini and Wang3 which, due to the reflection off of the surface and the Lloyd Mirror effect,37 will cause bottom scatterers to appear above the air/water interface. Similarly, the reciprocal surface-bottom path described in Belettini and Wang appears below the sediment interface. As shown in Belettini and Wang, the bottom-surface path defines the lower-bound (in the height dimension) of all multipath contributions for which the air/water interface is the last interface interaction. Likewise, the surface-bottom multipath, which maps below the sediment, is the upper bound of all multipath arrival angles for multipath rays which have the sediment as the last interface reflection. The principle can therefore be inferred that, in an isovelocity medium, multipath signals never map between bounding water column interfaces. This principle is useful for ruling out observed phenomena within the water column as multipath scattering.

The histogram images created from these points were intensity weighted, i.e., for each histogram bin (pixel) Hi,j,
H i , j = 1 N n 1 N P n 2 ,
(12)
where N is the number of points contributing to the spatial histogram bin and Pn is the magnitude of the nth contributor to the bin. Subsequently, the histogram images are filtered to improve the cross correlation peak prominence between adjacent scans: gaps in the histogram are filled using the median of the surrounding values, the image is thresholded to remove lower-level features, a spatial high-pass filter is applied to accentuate edges, and the dynamic range is compressed using the hyperbolic tangent as a soft-clipping function. Figure 10 shows an example intensity weighted histogram of a five-ping set before and after this pre-correlation image conditioning has been applied.
FIG. 10.

Plot (A) shows the intensity-weighted spatial histogram of the point-set for all points with height values –4 < z < 2. The lobe-and-cleft structure of the front is clearly visible in this set of pings. Plot (B) shows the data in plot (A) following image pre-filtering for pair-wise image cross-correlations.

FIG. 10.

Plot (A) shows the intensity-weighted spatial histogram of the point-set for all points with height values –4 < z < 2. The lobe-and-cleft structure of the front is clearly visible in this set of pings. Plot (B) shows the data in plot (A) following image pre-filtering for pair-wise image cross-correlations.

Close modal

In Fig. 10(A), the dominant features in the histogram image are the prominent lobe and cleft structures at the leading edge of the front and the current-driven bubble plume that persists downstream of the front. Notably, as discussed later in greater detail, regions of relatively weak intensity trail behind the cleft features. Several band-like regions corresponding to dropouts are visible in the histogram. These dropout regions correspond to locations in which bright crests of sand-ripples were the dominant scatterers rather than the plume front. As described previously, the values in these dropout bands were “inpainted” using the median of the surrounding pixels to prevent the large discontinuities they create from biasing correlations between image pairs.

Following generation of the pre-filtered and compressed front representations, a complete set of pair-wise correlations are performed between all sets that have field-of-view overlap of greater than 50%. As previously done for AOA estimation, parabolic interpolation is applied to the resulting correlation peaks to make displacement estimates with sub-pixel precision. These displacement estimates are then used to determine the relative northern and eastern translations between histogram images. Displacement measurements found using correlation peaks can be unreliable when correlation peaks have poor SNR; therefore, a quality value based on local correlation peak prominence is used to weight the displacement estimates created from each correlation pair. The quality metric used in the present study is the ratio of the correlation peak height to the median of the absolute value of the cross correlation function for a 15 meter radius surrounding the correlation peak. Similar to the sub-aperture alignment method described in Marston and Kennedy,38 the measured image shifts and quality estimates are combined in an over-determined system of equations and a solution for the locations of each histogram is found via iterative-reweighted least squares (IRLS), i.e.,
x j = A T W j A + γ 1 1 W j A T x ̃ ,
(13)
W j + 1 = γ 2 W 1 γ 2 + x ̃ A x j 2 ,
(14)
where superscript (j) represents iteration number, A is the matrix relating the aperture offset x to the vector of relative aperture shift estimates x ̃, W is a diagonal matrix of weights, initialized in W 1 using the previously described quality values based on the correlation peak local prominence, and γ 1 , 2 are regularization parameters: γ 1 is a very small value ( γ 1 = 10 6 in the present study) which forces the solution x j to correspond to the zero-mean solution, and γ 2 functions as a measurement outlier penalizer, where smaller values of γ 2 result in stronger outlier rejection and a sparser weighting matrix. In the present implementation, γ 2 = 0.1. Equations (13) and (14) are iterated to solve for the x- and y- offsets that align the histograms and fit the displacement observations in a weighted least squares sense. The present implementation limited the total number of IRLS iterations to five. From the resulting x- and y- displacements, the average plume-front feature speed and direction are found via
V plume ( SAS ) = x N x 1 N 1 2 + y N y 1 N 1 2 τ N τ 1 N 1 ,
(15)
θ plume ( SAS ) = π 2 tan 1 y N y 1 x N x 1 ,
(16)
and their associated uncertainties are characterized by the standard deviations of their zero-mean inter-set differences,
σ v _ plume ( SAS ) = ± 1 N 2 n = 1 N 1 V n V plume 2 ,
(17)
σ θ _ plume ( SAS ) = ± 1 N 2 n = 1 N 1 arg e i θ n θ plume 2 ,
(18)
where
V n = x n + 1 x n 2 + y n + 1 y n 2 τ n + 1 τ n ,
(19)
θ n = π 2 tan 1 y n + 1 y n x n + 1 x n .
(20)
In Eqs. (15) and (19), τ n is the average time-stamp of the pings contributing to the nth point cloud. Note that in this paper, the convention used for the plume direction of motion θ plume is converted to the standard map angle convention: North = 0 degrees, East = 90 degrees. Correspondingly, the convention for the x-axis is north and the y-axis is east.

Finally, a new fully aligned point cloud is generated from the entire set of pings visualizing the front by interpolating the offset locations x τ 1 N , y τ 1 N to the time-stamps corresponding to the individual pings. In the present implementation, linear interpolation is used during this step. The point clouds from all contributing pings, each of which is now compensated for the motion of the plume features, are then combined. A combined set of visualizations in a stationary, geo-rectified frame should undergo defocusing due to “motion blur” induced by the non-stationarity of the front. The motion-compensated version, an image created in the travelling coordinate system of the local plume observations, should have improved focus and enable better resolution of plume features. As shown in Fig. 11, this is in fact what is observed.

FIG. 11.

(Color online) Plot (A) shows a 32-ping weighted histogram image of the point cloud points spanning –4 < z< 2 meters relative to the surface. The image was created in the fixed, geo-rectified coordinate frame. Plot (B) shows the same data but with a geo-rectification strategy that compensates for the moving features of the front. Superimposed on plot (B) are the velocity vectors corresponding to the average motion of the plume front features as measured by the sas [Vplume(sas)= 0.20 m/s] and the frontal velocity in the front-normal direction as measured by the radar [Vplume(RADAR) = 0.14 m/s].

FIG. 11.

(Color online) Plot (A) shows a 32-ping weighted histogram image of the point cloud points spanning –4 < z< 2 meters relative to the surface. The image was created in the fixed, geo-rectified coordinate frame. Plot (B) shows the same data but with a geo-rectification strategy that compensates for the moving features of the front. Superimposed on plot (B) are the velocity vectors corresponding to the average motion of the plume front features as measured by the sas [Vplume(sas)= 0.20 m/s] and the frontal velocity in the front-normal direction as measured by the radar [Vplume(RADAR) = 0.14 m/s].

Close modal
It is evident from Fig. 11(B) that many fine-scale lobe and cleft features that were originally obscure or defocused in the stationary frame become visible when combined in the moving frame of the plume. Figure 11(B) shows two velocity values: Vplume(sas), the average velocity of the features determined by frame shifts in the sonar imagery, and Vplume(RADAR), the velocity of the front in the front-normal direction detected by the radar. As estimated by the sonar, Vplume(sas) = 0.20 m/s [ σ v _ plume ( SAS )= ±0.04] at θplume(SAS)= 200 degrees [ σ θ _ plume ( SAS )= ±9]. The radar front-normal velocity is Vplume(RADAR) = 0.14 m/s [ σ v _ plume ( RADAR ) = ± 0.04] at θplume(RADAR) = 247 degrees. Note that the directions appear to be very different; however, the sonar is tracking the velocity of the individual cusp features which, in addition to a velocity component due to the expansion of the front, include an along-front velocity component as well.36 Projecting the velocity of the plume features in the direction of the front-normal determined by the radar results in
V normal = V plume cos θ radar θ plume .
(21)
When plugging in the observed values, V normal = 0.14 m/s, which indicates good agreement between the sonar and radar observations. Figure 12 places the sonar measurement in the larger context of the radar observations for the rest of the tidal front structure, annotated with front-normal velocities and the sonar position.
FIG. 12.

(Color online) Oregon State Radar (OSU) observations of the front at the time of sonar data collection. Front normal velocities are shown in the figure along with the position of the sonar (blue circle). The image represents an average of 64 radar scans collected over 80 s. Date of data capture is June 26, 2017.

FIG. 12.

(Color online) Oregon State Radar (OSU) observations of the front at the time of sonar data collection. Front normal velocities are shown in the figure along with the position of the sonar (blue circle). The image represents an average of 64 radar scans collected over 80 s. Date of data capture is June 26, 2017.

Close modal
Though the radar gives an accurate estimate of the front-normal velocity, it does not have a high enough resolution to determine the along-front velocity of the cusp features. Using sonar data, the along-front component of the velocity can be determined via computing the velocity perpendicular to the normal,
V along _ front = V plume cos θ radar θ plume π / 2 .
(22)
Plugging in the observed values results in a speed of: V along _ front = 0.15 m/s. Simpson et al.36 report that the propagation speeds of these features along the front, as estimated via optical footage captured by a drone, had a range of velocities from nearly stationary to almost 1 m/s, with a mean of 0.35 m/s and wide peak in the distribution from 0.3 to 0.5 m/s.

A 3D visualization of the set of points along the plume was created by generating a 3D intensity-weighted histogram and rendering using alpha-blending, a technique that allows for occlusion by assigning each voxel a transparency value.38 Voxel transparency was linearly mapped to the intensity-weighted histogram value of the voxel. Figures 13 and 14 show the spatial scales of the motion-compensated plume features. The scalloping of the injected bubble-plume travelling behind and underneath the front is especially visible in Fig. 14, which shows the bubble-plume from underneath, after having removed the sediment contributions to prevent occlusion. A notable feature is that the depths of the bubbles extending from the cleft structures appear to be shallower, as evidenced by the green (i.e., higher altitude) veins trailing behind the clefts. The mechanisms for this observation are unclear but it suggests that areas extending from the clefts may be areas of localized upwelling, or at least weaker than downwelling right at the front. The lobe and cleft structures observed at CTR vary from another type of ebb plume instability referred to as “lobe-cleft instability.” Previous observations made by Horner-Devine and Chickadel using infrared images from the Merrimack River ebb plume showed unique thermal signatures extending from similar-looking cleft like features.39 These features are related to localized upwelling drawing cooler waters from below. The instabilities driving the lobes and clefts at the CTR are hypothesized to be associated with the cross-front shear, which varies from those in Horner-Devine and Chickadel40 where there is no significant shear across the front. Simpson et al.36 hypothesize a different mechanism for the generation of the features in CTR. Further insight into the source of these lobe and cleft structures is beyond the scope of this manuscript.

FIG. 13.

(Color online) Renderings from different perspectives of the intensity-weighted 3D histogram of the motion-compensated point cloud containing 32 pings worth of observations, or ∼13 s, of the tidal plume front. The reference lines are 10 meters in the North, East, and Z dimensions, with ticks spaced every 2.5 meters for scale reference. The 3D histogram bin size was 25 cm on edge.

FIG. 13.

(Color online) Renderings from different perspectives of the intensity-weighted 3D histogram of the motion-compensated point cloud containing 32 pings worth of observations, or ∼13 s, of the tidal plume front. The reference lines are 10 meters in the North, East, and Z dimensions, with ticks spaced every 2.5 meters for scale reference. The 3D histogram bin size was 25 cm on edge.

Close modal
FIG. 14.

(Color online) A view looking upwards from below the plume. The scalloping of the bubble plume beneath the lobes is more visible from this viewpoint. The reference lines are identical in size and position to the reference lines shown in Fig. 13.

FIG. 14.

(Color online) A view looking upwards from below the plume. The scalloping of the bubble plume beneath the lobes is more visible from this viewpoint. The reference lines are identical in size and position to the reference lines shown in Fig. 13.

Close modal

The same plume motion estimation and alignment procedure was applied to an observation of the sonar front on the 28th of June 2017, two days after the previous example. In this scan, the AUV was travelling along to the front at a standoff distance of ∼50 meters. Figure 15 contains volume renderings of data from this scan. The velocity magnitude of the features in Fig. 15 is much higher than in the previous encounter (0.4 ± 0.08 m/s); however, the observed front normal velocity was much smaller: ∼0.06 ± 0.01 m/s, and once again was similar to the radar-estimated front-normal velocity of 0.08 ± (0.04) m/s.

FIG. 15.

(Color online) Twenty-ping volume renderings of a second encounter with the ebb plume front at the CTR, rendered in the moving frame of the front. (A) Shows a view from above; the same 10-meter reference lines are visible in this plot as in Figs. 13 and 14. Note that motion blur causes the sediment to defocus because the image has been created in the moving frame of the plume front features. (B) A view from the underside of the plume front feature. The lobe features are visible; however, unlike in Figs. 13 and 14, the scalloping is not as visible, and it appears that the height distribution of the bubbles under the plume is more homogeneous. (C) Shows a maximum-intensity projection through the northings-dimension of the volume. Unlike in the previous encounter with the plume, the bubbles appear to ascend after an initial descent; the region in which the bubbles are ascending is annotated via the white arrow in Fig. 15(C).

FIG. 15.

(Color online) Twenty-ping volume renderings of a second encounter with the ebb plume front at the CTR, rendered in the moving frame of the front. (A) Shows a view from above; the same 10-meter reference lines are visible in this plot as in Figs. 13 and 14. Note that motion blur causes the sediment to defocus because the image has been created in the moving frame of the plume front features. (B) A view from the underside of the plume front feature. The lobe features are visible; however, unlike in Figs. 13 and 14, the scalloping is not as visible, and it appears that the height distribution of the bubbles under the plume is more homogeneous. (C) Shows a maximum-intensity projection through the northings-dimension of the volume. Unlike in the previous encounter with the plume, the bubbles appear to ascend after an initial descent; the region in which the bubbles are ascending is annotated via the white arrow in Fig. 15(C).

Close modal

The vertical scalloping effect readily visible behind the clefts of The front in Figs. 13 and 14 appears mostly absent in Fig. 15. Not only are the depths of the bubbles much more homogenous, but the slope of the plume as it descends appears more gradual and is followed by a region in which bubbles ascend back to the surface. This is consistent with bubble kinematics, in that bubbles are brought to depth by downwelling at the front, but the downwelling signal only extends a few meters behind the front. Beyond this range, bubbles are advected by local currents and will gradually rise to the surface driven by their own buoyancy or will dissolve. As shown in Fig. 7, similar features are observed using echosounders. Another feature of importance in Fig. 15(C) is the strong amount of scattering that appears to be originating from significant distances above the surface. Though some of the near-surface contributions to this region may be due to the previously mentioned uncertainties in arrival angle and vehicle roll, these are not significant enough to account for the vertical spread shown in the figure. In fact, these vertical features are more likely due to the higher order multipath contributions described in Belettini and Pinto,4 specifically, the multipath rays in which the last interaction is with the surface or, in the present case, with near-surface bubbles.

Though the ray-path bending effects of hydrodynamic features such as internal waves have been observed and studied using sas systems,41–43  sas systems have not previously been used as a tool to make direct observations of hydrodynamic features in the mid-to-upper water column. Despite being an unorthodox usage of a sas array, the field data collected at the mouth of the CTR during the USRS DRI revealed that sas systems are capable of capturing near-surface estuarine-specific hydrodynamic features such as the tidal-plume fronts shown in this paper. Furthermore, these features were captured without any extra efforts oriented at suppressing sediment features; rather, backscattering from bubbles in the tidal plume front was frequently the dominant backscattering feature for the portions of imagery containing their backscattering contributions.

The challenges associated with coherently processing non-stationary scatterers in a standard synthetic aperture multi-ping coherent post-processing framework were circumvented by operating in an alternative post-processing framework that leveraged point cloud data created using physical array processing applied to real-aperture data collected from individual pings. Dense multi-ping point clouds and 3D histograms were synthesized from combinations of these single-ping point clouds. To leverage the planar array of the deployed sas system, solve angular ambiguities caused by phase wraps, and reveal the 3D structure of the observed environment features, a correlation-beamforming algorithm was exploited to determine the elevation angles corresponding to the dominant scatterer observed at any given range and azimuth. Using this approach, point clouds could be created from every single ping transmitted by the sonar. Matching point cloud features between pings enabled the velocities of features to be determined and improved the focus of multi-ping two- and 3D histogram renderings of the plume front.

A variety of tidal plume front features were observed with the deployed sas system and the post-processing described in this paper: a lobe-and-cleft structure was present in both encounters of the front, and the lobe features had length scales ranging from several meters to just over ten meters. In both observations, the along-front speed of the plume features matched or exceeded the front-normal velocity. Furthermore, the front-normal velocities in both cases nearly matched independent estimates derived from land-based radar observations.

The first Encounter with the front showed an initial descent of –0.62 m meters in height per meter traveled behind the front, and this tapered off when the bubbles reached approximately 3 meters in depth. The measured plume depth is consistent with other sensors and the estimated descent rate is also reasonably consistent with most measurements obtained using other acoustic sensors.32,35 The lobe-and-cleft structure of the plume front was accompanied by a scalloping in the vertical distribution of the bubble-plume which was visible when viewed from below (Fig. 14). In contrast, the second encounter did not show a steep initial descent, and bubbles seemed to be distributed more homogeneously and more closely to the surface. Furthermore, approximately 15 meters after the leading edge of the front, the bubbles appear to rise in elevation. This is likely driven by buoyant bubbles slowly rising in the water column in the absence of downwelling currents away from the front. Ray bending is an unlikely explanation in this scenario because the observed pattern would require downward refraction, whereas when crossing the front, an upward refracting sound speed profile was observed within the buoyant plume waters.

An argument based on the Lloyd's Mirror effect was put forward that the angles of arrival corresponding to multipath arriving from the surface never map beneath the surface in an iso-velocity environment. Combined with the similar reality that, as shown in Belletini et al.,3,4 multipath arriving from the sediment never maps above the sediment, the general rule for multipath arrivals in an isovelocity medium is that they never map between the sediment and water interfaces. This principle can be used to reject multipath as a potential source for many of the bright near-surface features visible in the sonar point clouds and 3D histograms shown in this paper, but it also explains the distribution of energy farther above the surface that was particularly noticeable in the second observation of the front and is shown in Fig. 15(C).

The potential benefits to oceanographic studies based on sas are derived from the ability to resolve much larger areas, in 3D, than alternative sonar and in situ sampling techniques. For example, echosounder or acoustic Doppler velocity current profiles, which have been used in estuarine hydrodynamics studies for decades (Famer,44 Trump,30,32 Kilcher,45 Lavery33) can achieve much higher resolutions than physical sampling techniques and can resolve scattering from features bubbles, turbulent microstructure, stratification, and suspended sediment. Because their sampling is limited to the relatively small ensonified beam, temporal and spatial variability must either be assumed or inferred from differences between transects occurring at different times. Additional spatial coverage can be obtained through the use of multibeam or sidescan sonars (Thorpe46) but overall coverage is still small when compared to the sas approach described here. The ability to capture 3D features of an environment over a large region in a single snapshot enables feature tracking that is much more difficult with systems that sample more sparsely. Herein, the ability to effectively identify and track features of interest associated with a tidal plume from using a sas system is demonstrated and results generally agree with observations derived from commonly employed sampling techniques. This suggests that new methods for observing features using sas could potentially lead to new insights into the hydrodynamic processes that are observed by the system.

This work was supported as part of the Undersea Remote Sensing Directed Research Initiative, and funded via Office of Naval Research grant numbers: N00014-16-1-2854, N00014-16-1-2948, N00024-17-D-6421, and N00014-19-1-2593. The authors would additionally like to thank Rocky Geyer, David Ralston, Andone Lavery, and Joe Jurisa for providing the supporting physical oceanographic and acoustic scattering measurements presented here.

1.
R. E.
Hansen
,
Introduction to Synthetic Aperture Sonar
(
Intech
,
London, UK
,
2011
), pp.
3
28
.
2.
M. P.
Hayes
and
P. T.
Gough
, “
Synthetic aperture sonar: A review of current status
,”
IEEE J. Oceanic Eng.
34
,
207
224
(
2009
).
3.
A.
Bellettini
,
M.
Pinto
, and
L.
Wang
, “
Effect of multipath on synthetic aperture sonar
,” in
Proceedings of the World Congress on Ultrasonics
,
Paris, France
(
September 7–10
,
2003
). pp.
531
534
.
4.
A.
Bellettini
and
M.
Pinto
, “
Design and experimental results of a 300-kHz synthetic aperture sonar optimized for shallow-water operations
,”
IEEE J. Oceanic Eng.
34
,
285
293
(
2009
).
5.
A. E. A.
Blomberg
,
A.
Austeng
,
R. E.
Hansen
, and
S. A. V.
Synnes
, “
Improving sonar performance in shallow water using adaptive beamforming
,”
IEEE J. Oceanic Eng.
38
,
297
307
(
2013
).
6.
T. O.
Saebo
,
H. J.
Callow
,
R. E.
Hansen
,
B.
Langli
, and
E. O.
Hammerstad
, “
Bathymetric capabilities of the HISAS interferometric synthetic aperture sonar
,” in
Proceedings of OCEANS 2007
(
IEEE
,
New York
,
2007
).
7.
R. E.
Hansen
, “
Mapping the ocean floor in extreme resolution using interfermetric synthetic aperture sonar
,”
Proc. Mtgs. Acoust.
38
,
055003
(
2019
).
8.
D.
Shea
,
D.
Dawe
,
J.
Dillon
, and
S.
Chapman
, “
Onboard real-time SAS processing—Sea trials and results
,” in
Proceedings of 2014 Oceans–St. John's
(
IEEE
,
New York
,
2014
).
9.
D. C.
Ghiglia
and
M. D.
Pritt
,
Two-Dimensional Phase Unwrapping: Theory, Algorithms, and Software
(
Wiley
,
New York
,
1998
).
10.
R.
Gens
, “
Two-dimensional phase unwrapping for radar interferometry: Developments and new challenges
,”
Int. J. Remote Sens.
24
,
703
710
(
2003
).
11.
H.
Yu
,
Y.
Lan
,
Z.
Yuan
,
J.
Xu
, and
H.
Lee
, “
Phase unwrapping in InSAR: A review
,”
IEEE Geosci. Remote Sens. Mag.
7
,
40
58
(
2019
).
12.
S. M.
Banks
,
T.
Sutton
, and
H.
Griffiths
, “
Interferometric synthetic aperture sonar
,” in
Proceedings of the London Communication Symposium
,
London, United Kingdom
(
September 18–21
,
2000
). pp.
172
176
.
13.
S. R.
Silva
,
S.
Cunha
,
A.
Matos
, and
N.
Cruz
, “
Shallow water height mapping with interferometric synthetic aperture sonar
,” in
Proceedings of OCEANS 2008
,
Kobe, Japan
(
April 8–11
,
2008
).
14.
H.
Zhong
,
J.
Tang
,
S.
Zhang
, and
M.
Chen
, “
An improved quality-guided phase-unwrapping algorithm based on priority queue
,”
IEEE Geosci. Remote Sens. Lett.
8
,
364
368
(
2011
).
15.
H.
Zhong
,
J.
Tang
,
S.
Zhang
, and
X.
Zhang
, “
A quality-guided and local minimum discontinuity based phase unwrapping algorithm for InSAR/InSAS interferograms
,”
IEEE Geosci. Remote Sens. Lett.
11
,
215
219
(
2014
).
16.
B.
Thomas
,
A.
Hunter
, and
S.
Dugelay
, “
Phase wrap error correction by random sample consensus with application to synthetic aperture sonar micronavigation
,”
IEEE J. Oceanic Eng.
46
,
221
235
(
2021
).
17.
T. O.
Saebø
,
S. A. V.
Synnes
, and
R. E.
Hansen
, “
Wideband interferometry in synthetic aperture sonar
,”
IEEE Trans. Geosci. Remote Sens.
51
,
4450
4459
(
2013
).
18.
T. O.
Saebø
, “
Seafloor depth estimation by means of interferometric synthetic aperture sonar
,” Ph.D. thesis,
University of Tromsø
,
Tromsø, Norway
(
2010
).
19.
D. A.
Cook
, “
Synthetic aperture sonar motion estimation and compensation
,” M.S. thesis,
Georgia Institute of Technology
,
Atlanta, GA
(
2007
).
20.
S. A. V.
Synnes
,
T. O.
Saebo
, and
R. E.
Hansen
, “
Interferometry using phase slope estimation
,” in
Proceedings of EUSAR 2016: 11th European Conference on Synthetic Aperture Radar
,
Hamburg, Germany
(
June 6–9
,
2016
), pp.
761
764
.
21.
E.
Ruigrok
,
S.
Gibbons
, and
K.
Wapenaar
, “
Cross-correlation beamforming
,”
J. Seismol.
21
,
495
508
(
2017
).
22.
D. C.
Brown
,
S. F.
Johnson
, and
D. R.
Olson
, “
A point-based scattering model for the incoherent component of the scattered field
,”
J. Acoust. Soc. Am.
141
,
EL210
EL215
(
2017
).
23.
E.
Coiras
,
Y.
Petillot
, and
D. M.
Lane
, “
Multiresolution 3-D reconstruction from side-scan sonar images
,”
IEEE Trans. Image Process.
16
,
382
390
(
2007
).
24.
M.
Hayes
and
A.
Hunter
, “
Synthetic aperture bathymetry estimation with a multielement array
,” in
Proceedings of the 27th Conference on Image and Vision Computing
,
Wellington, New Zealand
(
November 25–27,
2012
), pp.
364
369
.
25.
T.
G-Michael
,
B.
Marchand
,
J. D.
Tucker
,
T. M.
Marston
,
D. D.
Sternlicht
, and
M. R.
Azimi-Sadjadi
, “
Image-based automated change detection for synthetic aperture sonar by multistage coregistration and canonical correlation analysis
,”
IEEE J. Oceanic Eng.
41
,
592
612
(
2015
).
26.
R. W.
Garvine
, “
Physical features of the Connecticut River outflow during high discharge
,”
J. Geophys. Res.
79
,
831
846
, https://doi.org/10.1029/JC079i006p00831 (
1974
).
27.
R. W.
Garvine
and
J. D.
Monk
, “
Frontal structure of a river plume
,”
J. Geophys. Res.
79
,
2251
2259
, https://doi.org/10.1029/JC079i015p02251 (
1974
).
28.
J.
O'Donnell
,
G. O.
Marmorino
, and
C. L.
Trump
, “
Convergence and downwelling at a river plume front
,”
J. Phys. Oceanogr.
28
,
1481
1495
(
1998
).
29.
B.
Baschek
,
D. M.
Farmer
, and
C.
Garrett
, “
Tidal fronts and their role in air-sea gas exchange
,”
J. Mar. Res.
64
,
483
515
(
2006
).
30.
G. O.
Marmorino
and
C. L.
Trump
, “
High-resolution measurements made across a tidal intrusion front
,”
J. Geophys. Res.
101
(
C11
),
25661
25674
, https://doi.org/10.1029/96JC02384 (
1996
).
31.
C.
Bassett
and
A. C.
Lavery
, “
Observations of high-frequency acoustic attenuation due to bubble entrainment at estuarine fronts
,”
Proc. Mtgs. Acoust.
45
,
005001
(
2021
).
32.
C. L.
Trump
and
G. O.
Marmorino
, “
Mapping small-scale along-front structure using ADCP acoustic backscatter range-bin data
,”
Estuaries
26
,
878
884
(
2003
).
33.
A.
Lavery
,
W. R.
Geyer
, and
M.
Scully
, “
Broadband acoustic quantification of stratified turbulence
,”
J. Acoust. Soc. Am.
134
,
40
54
(
2013
).
34.
D. A.
Honegger
,
M. C.
Haller
,
W. R.
Geyer
, and
G.
Farquharson
, “
Oblique internal hydraulic jumps at a stratified estuary mouth
,”
J. Phys. Oceanogr.
47
,
85
100
(
2017
).
35.
C.
Bassett
,
A. C.
Lavery
,
R.
Petitt
, and
S.
Loranger
, “
Autonomous platforms for measuring broadband backscatter
,”
Proc. Mtgs. Acoust.
45
,
005002
. (
2021
).
36.
A. J.
Simpson
,
F.
Shi
,
J. T.
Jurisa
,
D. A.
Honegger
,
T.-J.
Hsu
, and
M. C.
Haller
, “
Observations and modeling of a buoyant plume exiting into a tidal cross-flow and exhibiting along-front instabilities
,”
JGR Oceans
127
,
e2021JC017799
(
2022
).
37.
G. V.
Frisk
,
Ocean and Seabed Acoustics: A Theory of Wave Propagation
(
Pearson Education
,
Upper Saddle River, NJ
,
1994
).
38.
T.
Porter
and
T.
Duff
, “
Compositing digital images
,” in
Proceedings of the 11th Annual Conference on Computer Graphics and Interactive Techniques
,
Minneaplois, MN
(
January 23–27
,
1984
), pp.
253
259
.
39.
T. M.
Marston
and
J. L.
Kennedy
, “
Spatially variant autofocus for circular synthetic aperture sonar
,”
J. Acoust. Soc. Am.
149
,
4078
4093
(
2021
).
40.
A. R.
Horner-Devine
and
C. C.
Chickadel
, “
Lobe-cleft instability in the buoyant gravity current generated by estuarine outflow
,”
Geophys. Res. Lett.
44
,
5001
5007
, https://doi.org/10.1002/2017GL072997 (
2017
).
41.
R. E.
Hansen
,
A. P.
Lyons
,
T. O.
Saebø
,
H. J.
Callow
, and
D. A.
Cook
, “
The effect of internal wave-related features on synthetic aperture sonar
,”
IEEE J. Oceanic Eng.
40
,
621
631
(
2015
).
42.
R. E.
Hansen
,
A. P.
Lyons
,
D. A.
Cook
, and
T. O.
Saebø
, “
Detection of Internal Waves Using Multi-Aspect Processing in Synthetic Aperture Sonar
,” in
Proceedings of EUSAR 2016: 11th European Conference on Synthetic Aperture Radar
,
Hamburg, Germany
(
June 6–9
,
2016
), pp.
757
760
.
43.
N.
La Manna
,
A. P.
Lyons
, and
R. E.
Hansen
, “
The effect of internal waves on synthetic aperture sonar resolution
,”
Proc. Mtgs. Acoust.
44
,
070013
(
2021
).
44.
D. M.
Farmer
and
J. D.
Smith
, “
Tidal interaction of stratified flow with a sill in Knight Inlet
,”
Deep Sea Res. Part A. Oceanogr. Res. Papers
27
,
239
254
(
1980
).
45.
L.
Kilcher
and
J.
Moum
, “
Structure and dynamics of the Columbia River tidal plume front
,”
J. Geophys. Res.
115
,
C05S90
, https://doi.org/10.1029/2009JC006066 (
2010
).
46.
S. A.
Thorpe
,
M. S.
Cure
,
A.
Graham
, and
A. J.
Hall
, “
Sonar observations of langmuir circulation and estimation of dispersion of floating particles
,”
J. Atmos. Oceanic Technol.
11
,
1273
1294
(
1994
).