Model-data comparison of sound propagation in a glacierized fjord with a simulated brash ice surface^a)

The cryosphere, the Earth's collection of frozen water as snow, sea and lake ice, ice sheets, glaciers, and frozen ground, is undergoing a period of significant and rapid change (Vaughan et al., 2013). Glaciers are masses of perennial ice that flow downhill due to their own weight and deformability. Many glaciers reside at or near their melting points and, as such, the retreat or advance of their terminus serves as an indicator of a changing climate as a result of their inherent sensitivity to fluctuations in temperature and precipitation (Zemp et al., 2008). Marine terminating glaciers are among the most dynamic glaciological systems as they are governed by the combination and interaction of simultaneous glacial and oceanographic processes (Howat et al., 2007; Joughin et al., 2012; Motyka et al., 2003; Nick et al., 2009). However, many of the physical processes taking place at the ice-ocean boundary are poorly understood and warrant further study.

Glacial fjords are among the loudest natural oceanic noise environments as a combination of calving (Glowacki, 2020; Glowacki et al., 2015; Pettit, 2012), ice cracking, and bubbles escaping from entrained air from melting glacier ice (Lee et al., 2013; Urick, 1971) have resulted in annual average sound pressure levels of 120 dB re 1 μPa (Pettit et al., 2015). These ambient sounds provide insight into the intensity, regularity, daily/weekly/seasonal variability, and directionality of these physical processes (Deane et al., 2014). However, the ambient noise field is strongly dependent on the propagation environment, which is partially dependent on the prevalence and distribution of ice at the ocean surface, necessitating acoustic propagation measurements and analyses (Glowacki et al., 2016). Additionally, the ice at the ocean surface acts as a unique boundary with complex acoustic reflection and scattering characteristics, as well as consisting of geometries capable of obstructing signal transmission paths. The combination and interaction of inflowing warm, deep salt water with freshwater from submarine melt along the glacier face driven by subglacial discharge (Jackson et al., 2020; Motyka et al., 2003) and additional melt from icebergs distributed throughout the fjord results in a complex, range-, depth-, and seasonally dependent sound speed profile (SSP; Glowacki et al., 2016).

To better recognize, quantify, and locate acoustic source mechanisms at the ice–ocean interface within glacial fjords, acoustic propagation within these environments must be more appropriately characterized (Zeh et al., 2019). While experimental and modeling work has been conducted to analyze propagation under sea ice (Alexander et al., 2016; Ballard, 2019; Diachok, 1976), the accumulation and formation of brash ice, itself materially different than sea ice and within the confines of a glacial fjord, necessitate propagation work in this environment. We present an analysis of an acoustic propagation experiment within a tidewater glacial fjord, focusing on inferring an effective acoustic model of the brash ice surface observed during the experiment. Through comparison of the propagation data with a ray-based acoustic propagation model of the environment, an effective brash ice surface model, including geometry and acoustic reflection characteristics, was determined.

An acoustic propagation experiment was performed in September 2017 in Hornsund Fjord in southwestern Svalbard, a Norwegian archipelago in the Arctic Ocean, as shown in Fig. 1(a). To the authors' best knowledge, this is the first active acoustic experiment conducted in a glacial fjord. The experiment used two stationary receiver HTI-96-min hydrophones (High Tech, Inc., Long Beach, MS), hung at nominal depths of 8 and 17 m to a buoy moored to the ocean floor, roughly 500 m from the terminus of Hansbreen glacier. An ITC-1007 (International Transducer Corporation, Santa Barbara, CA) source hung from a small boat at a nominal depth of 6.7 m was allowed to drift away from the two receivers, yielding propagation paths ranging from 20 to 200 m. The water depth along the entire propagation path was approximately constant at 67 m. A schematic diagram of the experimental geometry is shown in Fig. 2. The source signal was an m-sequence, 255 bits long, single cycle, phase modulated with an 11 kHz carrier and transmitted with a source level of 150 dB re 1 μPa at 1 m. The sequences were sent continuously for 10 min and then paused for several seconds to record the ambient environment. The receiver signals were high-pass filtered at 8 kHz and digitized at 96 kHz with 16 bits of resolution. Recordings of the receivers were each cross correlated with a replica of the transmitted signal to increase the gain and temporal resolution of the received signal by factors of 255 and 255, respectively. These gain values, calculated with the transmitted voltage signal as measurements of the in-water acoustic source signal, bearing the effect of the transducer response, were not performed. From the correlator output shown in Fig. 3, two distinct arrivals were observed and are considered to be the direct and surface reflected arrivals from the source.

The source and receiver signals were not time-synchronized; hence, the time difference of arrival (TDOA) between the direct and reflected signals was computed. The arrivals were originally noted visually through tracking peaks in scaled images of the cross correlation, similar to those shown in Fig. 3. The visually noted maxima were then used as a guide to window around ±10 samples, or roughly 0.1 ms, of the tracked peaks to determine the maximum arrival amplitudes and specific TDOAs of the two maxima within these windows. The noise, defined as the average sound pressure level within each 100 ms of the recorded signal outside of the 0.1 ms window of the direct arrival, was also saved. The direct and surface reflected arrivals that were less than 10 dB greater than the noise were discarded.

The range from the glacier and mooring float to the boat was documented regularly throughout the experiment using a Bushnell Tour Z6 JOLT laser rangefinder (Bushnell Corporation, Overland Park, KS), which has an accuracy of ±0.45 m between 4.5 and 114 m and ±0.91 m between 115 and 1189 m. To supplement and verify the rangefinder measurements, a handheld Global Positioning System (GPS) unit (Garmin GPSMap^® 64, Schaffhausen, Switzerland) was periodically used during the experiment, and a strong agreement was found between the two. The distance measurements from the rangefinder were correlated with the time of the recording to document the signal as a function of the distance from the source. The transmission loss (TL) of the direct and surface reflected arrivals was calculated as the source level minus the amplitude of the correlation process outlined above. The TDOA for the direct and surface reflected paths and the amplitude of these arrivals as a function of the distance from the source for the recorder at a 17-m depth (nominal) is shown in Fig. 4. The TDOA and arrival intensity data were then filtered using the matlab^® medfilt1 function (The MathWorks, Natick, MA) so as not to disproportionately weight outliers by filtering the data to the median observed across each 50 range data points or roughly 0.35 m. A curve fit with evenly spaced range values was then interpolated from the median-filtered TDOA for ease of comparison with the subsequent model.

A Valeport miniCTD profiler (Totnes, Devon, UK) was used to perform conductivity, temperature, and depth (CTD) measurements at the beginning and end of the propagation track, yielding SSPs at the receiver location (r = 0 m) and end of the drifting track (r = 200 m), that is, 500 and 700 m from the glacier terminus, as shown in Fig. 5. Small-scale perturbations in the sound speed, particularly those within the first 20-m depth of the original profiles, are considered to be local effects and unlikely to be present at all profile points along the experimental range. However, the dominant aspects of the sound speed environment, as depicted in the smoothed profiles in Fig. 5, are reasoned to be consistent along the propagation path and are similar to those measured two years prior in the same environment (Glowacki et al., 2016). The range-dependence of the sound speed environment was integrated into the eventual model through use of the quadrilateral SSP interpolation available in BELLHOP, resulting in a smoothly varying sound speed environment between the smoothed SSPs at 0 and 200 m.

Among the most straightforward propagation models used to predict the amplitude of and TDOA between ray paths is the two-ray reflection model (Rappaport, 2002). This model is analogous to this experiment when considering the direct and surface reflected ray paths, as outlined in Fig. 6. Although typically generalized for a stationary source and receiver, the geometric relationships from this model hold when allowing the source to drift along the experimental track at a speed far below the speed of sound in the water, as performed for this experiment.

Referring to Fig. 6, z_R is the depth of the receiver, z_S is the depth of the source, r is the horizontal distance between the source and receiver, d is the distance traveled by the direct ray, and s and $s ′$ are the distance traveled by the surface reflected ray for the ice-free and ice-covered surface cases, respectively. By using the method of images, several important relationships between these quantities may be determined. Assuming a constant sound speed c and flat surface at the air–water interface, the TDOA $Δ t$⁠, as a function of the horizontal distance from the source to the receiver r, may be expressed as

Also of importance is the range at which the surface reflected ray encounters the top surface, r₁ or r₂, which may be calculated as

Assuming spherical spreading along the direct and surface reflected path lengths, d and s are the TL along these ray paths, $TL d$ and $TL s$⁠, respectively, and may be approximated as

Through the addition of a top surface of variable height, h(r), within the model as shown in Fig. 6, when influenced by the top surface, the TDOA of the direct and surface reflected rays, $Δ t surf$⁠, yields a similar result to that in Eq. (1),

With the original $Δ t$ having no surface and the new $Δ t surf$⁠, the difference between these two cases, $Δ T$⁠, is

where $Z = z S + z R$⁠.

Or solving Eq. (5) for h(r) with the requirement that the surface height be less than z_S and z_R yields

While Eq. (1) is exact, Eq. (4) is an approximation. This approximation is carried through to Eqs. (5) and (6). Equation (4) provides a solution for $Δ t surf$ at each range value, r, using the variable height of the surface, h(r). However, this solution is approximate as previous and subsequent values for h(r) are ignored when calculating $Δ t surf$⁠. That is, the possibility of shading due to the variable height of the surface is disregarded for this calculation. This limitation is not ideal but is required to employ the straightforward calculation using the two-ray reflection model presented here.

Although an isovelocity SSP was not observed during the test, the simplification of the propagation environment to an isovelocity SSP for this step was made to better qualify the variability of the surface heights and make the calculations straightforward. However, to more appropriately account for the difference between the calculated isovelocity case and SSP observed within the tests, a correction was made to the experimentally encountered $Δ t surf$ data. The corrections were calculated through tracking the differences in the TDOAs as a function of the range in two, geometrically identical BELLHOP modeled environments with the same experimental source and receiver depths in an isovelocity SSP $Δ t ISO$ and the experimental range-dependent environment $Δ t SSP$⁠. The range-dependence of the sound speed environment was integrated in the BELLHOP modeled TDOA, $Δ t SSP$⁠, through a quadrilateral interpolation between the two smoothed SSPs in Fig. 5. The direct and surface reflected ray paths from the source to receiver were tracked along the range of the experiment in both SSP environments, and the ratios between the TDOAs between the isovelocity and range-dependent cases were determined as $C ( r ) = Δ t SSP / Δ t ISO$⁠. The correction was approximately unity near the beginning of the experimental track at 20 m for the deeper hydrophone, gradually increasing to 1.05 for $Δ t surf$ at 90 m, the maximum experimental range used to infer the surface in this work. This correction was used to effectively remove the effects of the experimentally encountered SSP by dividing the $Δ t surf$ of Eq. (4) by C(r).

The source and receiver depths are susceptible to surface wave motion caused by wind, calving at the glacier face (Massel and Przyborska, 2013), and capsizing icebergs (Levermann, 2011; MacAyeal et al., 2011). Additionally, the source and receiver were deployed via hanging cables from the receiver mooring and the drifting source boat, and as such, the source-to-receiver distance was subject to the movement of these cables from motion of the ocean and boat. As travel time difference measurements are sensitive to small deviations in the geometry of the system, a better understanding of this variability throughout the experiment is necessary. Further, as there were no depth sensor measurements taken alongside the acoustic source and receiver, a deviation from the reported and inferred depths is considered.

Toward this end, the TDOA data were used to refine the nominal source and receiver depths listed in Sec. II A. Again, the TDOA data were used, as opposed to absolute travel time of the signal from the source to the receiver, as the source and receiver instruments were not synchronized. Solving Eq. (1) for z_S for the deep receiver $z R 1$ and shallow receiver $z R 2$ with the requirement that the source depth be positive results in

Additionally, although the source and receivers are not synchronized, the pair of receivers was synchronized, recording simultaneously on two channels of the same device. As such, noting the TDOA between the direct arrivals, $Δ t dir$⁠, of the shallow and deep receivers results in

Solving for z_S yields

With Eqs. (7) and (8) determined, the source, z_S, and receiver depths, $z R 1$ and $z R 2$⁠, may be solved as a function of the range. This calculation assumes a constant, flat surface at the air–water interface and an isovelocity SSP.

The rays were tracked within the BELLHOP modeled environments through monitoring the arrivals at the calculated receiver depths along the experimental track. As a program, BELLHOP assumes a stationary source at a range of 0 m with receivers at prescribed ranges in the environment. However, invoking the principle of reciprocity, the source and receiver locations may be swapped to model a moving source, as was done in this work. 9999 rays launched at angles from −85° to +85° were tracked at 1001 equally spaced receivers at ranges from 0 to 250 m. The bottom reflected rays were discarded through windowing the depth of interest to disregard the bottom interactions. The travel time difference and arrival amplitudes of the first two arrivals were determined at ranges along the experimental track. The ray launch angles were documented, and the paths of the two rays along the transmission path were tracked to confirm that the rays were the direct and surface reflected arrivals (Zeh et al., 2019).

In BELLHOP, Gaussian beam shapes were used. Gaussian beams more often result in accurate TL calculations by eliminating shadows and caustics often observed using hat-shaped beams (Porter and Bucker, 1987). Gaussian beams often result in the perception of many simultaneous arrivals of varying amplitude at each receiver location due to their zone of influence. However, internally, BELLHOP incorporates logic to group together and sum simultaneous ray arrivals into a single ray arrival. The BELLHOP-reported ray arrival, whose amplitude is often the sum of several rays, is reported in this work.

Initially, an ice-free BELLHOP model was considered, using refined source and receiver depths in an isovelocity sound speed environment. The TDOA between the direct and surface reflected ray arrivals, $Δ t ( r )$⁠, was then calculated for this case. Comparing this initial, ice-free BELLHOP model to that of the measured data, corrected for the influence of the SSP as discussed in Sec. III A, the difference of the two cases' TDOAs as a function of the range, $Δ T ( r )$⁠, was calculated. $Δ T ( r )$ was then used to approximate the inferred surface height as a function of the range throughout the experiment through use of Eq. (6). The range at which the surface reflected ray path encounters the top surface was also recorded. The inferred surface height as a function of the range, h(r), was then rounded, effectively “squaring off” the surface, ensuring individual rays will only interact with flat surfaces. Although the interactions of the rays with the brash ice within the environment are not likely of this type, as an initial step in understanding the relative depth of penetration of these obstructions, this simplification was deemed appropriate.

The inferred surface height, h(r), was then implemented as the top surface in the range- and depth-dependent sound speed environment from the experiment. BELLHOP was then run again using this environment, and the TDOAs between the direct and surface reflected arrivals were compared with the measurements to evaluate the effects of the inferred surface on the TDOA, as discussed in Sec. IV B and shown in Fig. 9.

With the difference in the TDOAs between the data and modeled cases used to model the effective height and shape of the surface, additional efforts were undertaken to better understand the reflection characteristics of this surface through analyzing the TL/arrival amplitude from the model, which incorporates the inferred surface and the data. Initially, the inferred top surface shape from Sec. III C 1 was modeled as a pressure-release (p = 0) boundary, and the amplitudes of the modeled direct and surface reflected ray arrivals at the receiver were determined. Noting the grazing angle at the surface of the modeled surface reflected ray and the amplitude of the surface reflected ray from the data, the angle-dependent reflection coefficient, $R ( θ )$⁠, of the surface was calculated. This inferred reflection coefficient is comprehensive, including intrinsic loss due to the material properties of the surface as well as scattering losses dependent on the frequency and angle of the incident wave on the rough ice surface.

Comparing the amplitude of the surface reflected arrivals from the data to the modeled pressure-release top surface, a ratio between the two, equivalent to the reflection coefficient, R, of the top boundary, was determined as

where $p S , data$ and $p S , model$ are the pressures of the surface reflected arrivals in the experimental and BELLHOP modeled cases, respectively. For the instances that the pressure of the surface reflected ray from the data exceeded that of the model, the reflection coefficient was assumed to be one.

The angle made between the surface reflected ray and top surface is not recorded internally within BELLHOP. However, the launch angles of the rays from the source are recorded, and given the experimental geometry and lack of curvature observed in the s₁ segment of the surface reflected path, these launch angles are assumed to be approximately equal to the incident angle that the surface reflected ray makes with the top surface. At a 90 m range, for example, the ray launch angle at the source location is 15.3° and the actual grazing angle that the ray makes with the surface is 15.9°. The reflection coefficient, R, is recorded in terms of the grazing angle, θ, of the surface reflected ray paths at each range approximated under the assumption of straight rays.

The inferred reflection coefficient, $R ( θ )$⁠, was then integrated within the range- and depth-dependent sound speed environment from the experiment. BELLHOP was run again in this environment, and the amplitudes of the data and modeled surface reflected rays were noted. These modeled results were compared with the measured data to evaluate the modeled reflectivity of the surface in predicting the TL data.

Given the potential for source and receiver depth uncertainty as discussed in Sec. II C, as well as the lack of absolute depth measurements throughout the test, Eqs. (7) and 8 were used to refine the source and receiver depths. A constant sound speed was assumed with no top surface. Only the first 50 m of the experimental track were used to infer the source and receiver depths due to the minimal effect that slight deviations in the SSP from fine-scale salinity or temperature gradients will have on the acoustic propagation at these ranges and, as such, an isovelocity SSP was assumed for this calculation. A possible error of two samples (⁠ $± 2 / f s$⁠) from the addition of a possible error in the direct arrival time data, $Δ t dir$⁠, of $1 / f s$ and the surface reflected time data, $Δ t 1 , 2$⁠, also $1 / f s$⁠, was considered. As a majority of the variation within the inferred source and receiver depths is assumed to be due to the top boundary surface and not the variations in the source or receiver depth, the deepest values of the inferred source and receiver depths within the first 50 m of the experimental track are used in the modeling of the top surface below. The source depth was calculated to be 7.00 m (±0.01 m) deep, the shallow receiver was at 8.93 m (±0.11 m), and the deep receiver was at 17.75 m (±0.12 m).

These calculated values are 0.3–0.9 m deeper than the nominal measurements made at the start of the test from on board the boat. As discussed in Sec. III B, this discrepancy is likely partially attributable to natural environmental factors throughout the experiment, including surface waves that may or may not simultaneously effect the depths of one or both of the source and receivers, as well as possible swinging of the source or receivers from movement of the boat or mooring. Additional human error may also have occurred because of possibly inappropriately noting the lengths of the source and receiver cables, as these measurements were made only at the beginning of the experiment. The lack of certainty regarding these depths due to the lack of onboard depth sensors is not ideal. However, given the combination of natural variability present in this environment, as well as possible human error in the cable measurement, a discrepancy of up to 0.9 m is reasonable.

With the appropriate source and receiver depths calculated and equivalent uncertainty noted, the environment was modeled as outlined in Sec. III C 1 and used to determine an inferred top surface shape. Although the shallow and deep receivers are used to inform the source and receiver depths, only the deeper receiver was considered when inferring the top surface as a result the clarity of the surface reflected ray arrival over a greater range of the experiment. As the surface reflected arrival is difficult to regularly distinguish beyond ranges greater than 90 m, as shown in Fig. 4, only this portion of the propagation range is considered in modeling the top surface. With the surface shape determined and modeled in BELLHOP, the inferred reflection coefficient, $R ( θ )$⁠, was then calculated and incorporated into the model as outlined in Sec. III C 2. The inferred surface and subsequently determined surface loss as a function of the range are shown in Fig. 7. The surface loss as a function of the range was then used to inform the surface loss as a function of the grazing angle, shown in Fig. 8. As noted in Sec. III C 2 and Eq. (9), the reflection coefficient was assumed to be one, or a surface loss of 0 dB, for instances where the amplitude of the surface reflected ray exceeded that of the direct arrival.

With the surface height and reflection characteristics modeled, BELLHOP was then run again with the inferred surface and reflection characteristics. The resulting TDOA curves and TL curves are shown in Fig. 9.

The inferred surface, as shown in Fig. 7, varies between 0 and 30 cm with protuberances ranging from 0.5 to 5 m wide from peak-to-peak. Using pictures taken from the experimental site, the ice within the fjord ranges in height above the water surface from nearly 0 cm, where the surface is or is nearly ice-free, to portions with larger icebergs reaching 20 cm above the water surface. The majority of the ice along the propagation path appears to lie flush with the water surface, approximately 1–5 cm above the air–water interface. Due to buoyancy, roughly 88% of floating icebergs are submerged in the water column with the remaining ice necessarily at or above the water surface, and given a range of 0–20 cm ice heights above the water surface from the experimental pictures, it should be expected that 0–1.5 m of this ice would penetrate into the water column with the majority of ice penetrating from 10 to 35 cm (International Ice Patrol, U. S. C. G., 2015). Additionally, given ice fragments of this approximate height, the corresponding length scales are expected to be <5 m wide given the International Ice Patrol classifications (International Ice Patrol, U. S. C. G., 2009). The approximate height and width dimensions of the experimentally encountered ice fragments, estimated from pictures taken during the experiment, correspond well with those calculated using the method above, indicating that the experimentally derived surface is reasonably close to that physically encountered in the environment, as shown in Fig. 1(b).

With the inferred surface incorporated into the modeled environment using BELLHOP, as discussed in Sec. III C 1 and shown in Fig. 7, the TDOA between the direct and surface reflected ray paths observed in the data were compared to that from the model, as shown in Fig. 9(a). The tuned surface accurately reproduces the main features of the TDOA curve within the 20–90 m of the experimental track considered, corresponding to the 15–65 m range of the top surface. At greater ranges, the surface reflected arrival is often dominated by noise making the inferred surface difficult to predict and is not shown. Again, this inferred surface assumes ray interactions only with flat surfaces, eliminating unknown scattering effects at the water–ice interface that may influence the surface reflected path time of arrivals. The modeled direct and surface reflected arrivals with the reflective surface are shown in Fig. 9(b). Whereas the modeled surface reflection TL strongly mirrors that of the data, the modeled direct TL is typically 1–3 dB greater than that witnessed within the data. The increased intensity of the direct arrival witnessed in the data compared to that of the model could be due to the simultaneous arrival of multiple rays that perceptively act as a single, more intense individual arrival. Additionally, as the documented source level (150 dB re 1 μPa at 1 m) of the source was measured prior to the field experiment, the increased intensity of the direct arrival in the data could possibly be attributed to a slightly higher source level from the source during deployment.

As outlined in Sec. III C 2, differences in the surface reflected ray intensity between a perfectly reflective version of the inferred surface and that observed in the data were used to determine an angular-dependent reflection coefficient. The angular-dependent surface loss, shown in Fig. 8, may be compared to those derived from models that integrate the physical properties of glacier ice, including compressional wave speed, c_L = 3447 m/s, and compressional wave attenuation, $α L = 0.161$ dB/m, as measured at 11 kHz near 0 °C at Langenferner glacier (Meyer et al., 2019). Because of the lack of shear wave speed and shear attenuation data, these values were calculated based on relationships to the experimentally derived compressional wave properties, resulting in a shear wave velocity, $c S = c L / 2 = 1724$ m/s, and shear wave attenuation, $α S = 6 α L = 0.968$ dB/m (Clee et al., 1969; McCammon and McDaniel, 1985). Considering an elastic layer bounded by a fluid and vacuum, that is, glacier ice between water and air in this case, the surface losses from several different thicknesses of ice versus the grazing angle were determined for an 11 kHz signal, as shown in Fig. 10. The ice thicknesses considered encompass those witnessed during the experiment. Grazing angles <10° and >50° are not shown for the inferred surface in Sec. IV B as they are not encountered within the BELLHOP modeling due to the experimental geometry and corresponding paths between the source and receiver.

The inferred surface reflection from this experiment (Fig. 8) does not compare favorably with the modeled surface reflection in structure but does exhibit surface losses on the same order of magnitude as the model. There are several possible origins of this discrepancy. In particular, the elastic properties of ice are strongly dependent on the temperature, reasoned to be the result of relaxation and slipping at grain boundaries (Clee et al., 1969). As the temperature increases from below freezing and approaches melting at 0 °C, the compressional and shear wave speeds in the ice decrease, lending to an increased loss at lower grazing angles. Although the values used for the modeled surface loss in Fig. 10 are taken from measurements performed on a temperate glacier at roughly 0 °C, the ice encountered during this test was actively melting due to exposure to nonfreezing water in Hornsund Fjord. Further, although the relationship is assumed to be linearly related in this work, the compressional wave speed and attenuation in ice are more sensitive to increasing temperatures than their shear counterparts (Vaughan et al., 2016). Additionally, the variability in these physical properties across glaciers varying in geographic location, internal temperature, and structure is not well known (Meyer et al., 2019). The inferred surface implemented in the model also consists only of squared-off, flat surfaces. As such, ray interactions, particularly those at shallow grazing angles, do not take into account the finer structure and roughness of the ice. As discussed in Sec. III C 2, the inferred surface losses determined in this work are inclusive of scattering and intrinsic material losses, whereas the modeled surface loss in Fig. 10 only considers intrinsic losses. Scattering losses, dependent on the topography of the under-ice surface, can dominate the intrinsic material losses (Ballard et al., 2020), likely accounting for some of the discrepancy between the modeled and inferred reflection losses calculated in this work.

While treated sequentially in this work, the sound speed, source and receiver depths, and top surface jointly influence the direct and surface reflected sound propagation pathways. Given the experimental limitations discussed above, including the lack of absolute depth measurements of the source and receiver, undersampling of the water column, and poor understanding of the shape and reflection characteristics of brash ice, the solutions to these undeniably integrated factors often required independent investigation to achieve meaningful results. Although the solution presented here is not unique, the approach detailed presents a straightforward documentation of the dominant effects considered and provides a solution informed by the collected data.

Through calculating the amplitudes and TDOA of the direct and surface reflected rays from the experiment conducted in Hornsund Fjord, an effective top surface shape with reflection characteristics was inferred through tracking the ray paths within simulated environments in BELLHOP and use of the straightforward two-ray reflection model. When modeled in BELLHOP, inclusion of the inferred surface and reflection characteristics, shown in Fig. 7, results in ray TL and a TDOA that closely approach those observed in the experiment, as shown in Fig. 9. Inclusion of this inferred brash ice surface yields a significant improvement in the data-model agreement as compared to the model with a flat pressure-release surface (Zeh et al., 2019).

However, the brash ice present within this and other glacial fjords is a function of a variety of glacial and oceanographic conditions. Whereas the simulation results compare favorably with those observed within the TDOA data for this particular data set, additional fieldwork is necessary to determine the efficacy of this method to characterize the brash ice surface in similar environments over a variety of time spans. To most effectively use the two-ray reflection model as discussed in this work, future experiments must appropriately constrain the experimental geometry, including the source and receiver depths and range along the propagation path, limiting the propagation variability to that of the environment and collection of brash ice at the ocean surface. Continuous monitoring of the sound speed environment along the propagation path must be performed to effectively capture the spatiotemporal variability of the water column. While performing these field experiments, coincident data that yields the shape, density, and regularity of the brash ice surface from aerial drones and/or underwater cameras could be used to verify the coustically derived results. Better documenting the shape of brash ice from these observations will provide insight toward the influence that local slopes of the ice have on ray interactions, a factor disregarded in this work. Additionally, laboratory experiments analyzing the reflection characteristics of glacier ice are essential as unique material features, particularly entrained pressurized air in the ice, may significantly influence reflection and scattering tendencies when compared to other forms of naturally formed ice.

To best exploit the rich acoustic environment in tidewater glacial fjords, the acoustic propagation within these environments must be better understood. Through determining the role that brash ice plays in sound propagation, researchers may improve the usefulness of passive acoustic monitoring and more accurately predict underwater acoustic communication by anticipating the effects of this unique top surface. This work presents an initial step toward this end through analysis of modeled ray interactions with an inferred brash ice surface. Future work further accounting for the influence of ice in both forward and backward acoustic propagation models and recordings may allow for more effective localization, characterization, and quantification of notable glaciological processes, including melting and calving at the glacier face.

The authors foresee that acoustic monitoring stations may be deployed in tidewater glacial fjords. These stations will continuously record the ambient noise environment while moored acoustic sources will periodically recharacterize the propagation environment. These installations will allow researchers to track acoustic changes from the hourly and daily to seasonal and yearly timescales, and paired alongside visual data, such as time-lapse photography and satellite images, will aid the glaciological and climate research communities in monitoring trends in these dynamic environments.

M.Z. would like to acknowledge the National Defense Science and Engineering Graduate Fellowship (NDSEG) provided by the Department of Defense (DoD) and The University of Texas Graduate Continuing Fellowship for funding this work. The authors also wish to acknowledge the Office of Naval Research (ONR) and the Polish National Science Centre for funding the data collection efforts performed by O.G. and G.D. of the Scripps Institution of Oceanography. This work was supported by a subsidy from the Polish Ministry of Education and Science for the Institute of Geophysics, Polish Academy of Sciences. We gratefully acknowledge the support of the staff at the Polish Polar Station Hornsund. G.B.D. gratefully acknowledges financial support through ONR Grant No. N00014-21-1-2316.