© The Author(s), Published by EDP Sciences, 2022

1 Introduction

A common method to estimate abundance of fish in the oceans is based on measurement of acoustic backscatter from the water column [1, 2]. While 38 kHz is a main measurement frequency for quantification, simultaneous measurements in the range 18 kHz to 400 kHz are typically used to support species identification alongside trawl sampling. Multifrequency echosounder systems can operate at discrete frequencies and lately also across multiple continuous frequency intervals to improve range resolution and frequency response measurements. Some of the higher frequencies are also recommended for abundance estimation, e.g., 200 kHz for Atlantic mackerel [3, 4].

Quantitative measurements of acoustic backscattering are made by calibrated echosounders. Elastic spheres of materials such as copper (Cu) or tungsten carbide (WC) are suspended in the sound beam, and the echosounders’ measurements of their target strengths are adjusted against theoretical values (e.g. [5, 6]). Calibration data thus obtained for single targets on the sound beam axis are extended to the measurement of volume backscatter by applying a factory-calibrated beam pattern for each echosounder transducer [6, 7].

Linear sound propagation is assumed in fisheries acoustics and many other applications, with spherical spreading, sound absorption, and scattering as the only transmission loss mechanisms [8]. This approximation applies to small-amplitude signals. Above an application-dependent amplitude threshold, nonlinear effects must be considered, and the sound is said to be of finite amplitude [9]. As higher frequencies are employed for species identification and increasingly efficient echosounder systems have been taken into use, the assumption of linear sound propagation is not always valid for transmit power levels accessible in fisheries acoustics. Signal distortion due to finite-amplitude sound propagation has been observed to cause additional transmission loss when the transmit power is increased for measurement at longer distances (several hundred meters). Investigations have demonstrated, first experimentally [10], and later by numerical simulations [10–14] and further experiments [13, 14], that finite-amplitude effects can be significant for echo sounders used for fish abundance estimation and species identification. For available signal transmission powers (up to a few kW) and frequencies, the predominant finite-amplitude effects are convection of the medium and dependence of sound speed upon pressure [9]. Amplitude-dependent waveform distortion is observed as increased transmission loss at the measurement frequency [10] (herein referred to as “nonlinear loss”) and generation of harmonics that can interfere with simultaneous measurements at higher frequencies [10]. As multifrequency systems have evolved from operating at discrete frequencies to continuous frequency intervals, it has become more difficult to avoid such interference [15].

Nonlinear loss can lead to error in the evaluation of the backscattering cross-section σ_bs and the volume backscattering coefficient s_v used in fish abundance estimation and species identification. The waveform distortion caused by finite-amplitude sound propagation accumulates over propagation distance through regions where the sound pressure amplitude is sufficiently high, and depends upon transmit power, signal frequency, transducer directivity, medium parameters, range, and angle off the sound beam axis. The resulting nonlinear loss manifests itself as excess axial transmission loss and flattening of the beam pattern [11–13].

Advice on maximum power settings was given for some echosounder transducers and operating frequencies in 2002 [16] and 2008 [17], to reduce finite-amplitude effects to a presumably insignificant level. However, the recommendations differ and overly simplified extrapolations are used, based on nominal power per transducer surface area [17]. The transmit signal power is adjustable in many of today’s echosounder systems. A recent system also makes it possible to limit the maximum user-selectable power setting for each connected transducer. The available recommendations do not indicate how much nonlinear loss may remain when the recommended power settings are used. Further work is needed to manage finite-amplitude effects at a quantitative level, including more accurate calculations specific to the relevant echosounder transducers. A method should be established formally, to quantify finite-amplitude effects so that related measurement errors can be avoided or compensated for. The present work does not replace the existing recommendations directly, but intends to improve the empirical and methodological basis for revisiting such recommendations and their form.

A set of field data, acquired in 2004, shows systematic differences in measured s_v in experimental acoustic survey measurements when alternating between “low” and “high” transmit powers. The data presented here are improved in volume and quality compared to what has been reported earlier [13]. Numerical simulations indicate that the observed systematic differences can be ascribed to finite-amplitude sound propagation in seawater. A scheme for model-based correction for finite-amplitude sound propagation effects was proposed in [13], building on the conventionally used power budget equations (1) and (2). Formal derivations supporting the validity of these equations in acoustics were however not found in scientific literature at the time of the original study. This lack of foundation was an obstacle to publishing the experimental and numerical results and thus led to a series of works [18–21] to establish the required formal derivations. The present paper presents, firstly, a first comparison of experimental volume backscattering field data with numerical simulation results within this newly founded mathematical framework [21]. It is demonstrated that the framework established in refs. [18–21] is functional for describing the measurement of volume backscattering under conditions of significant finite-amplitude sound propagation. Secondly, a method to account for nonlinear loss in echosounder calibration is proposed within the same mathematical framework, to arrive at a complete model for quantifying and correcting for finite-amplitude sound propagation in fisheries research echosounders. A “finite-amplitude correction factor” is introduced to account for nonlinear loss upon both echosounder calibration and measurements of volume backscattering.

The paper is organized as follows: Section 2 describes the theoretical models used. This includes (i) the small-amplitude power budget equations conventionally used in fisheries acoustics, (ii) recent extensions of these power budget equations to finite-amplitude conditions, (iii) a further amendment to account for finite-amplitude effects during echosounder calibration and field operation, and (iv) theory for numerical simulation of finite-amplitude sound propagation in a fluid. The oceanic survey measurements and numerical simulations are described in Sections 3 and 4, respectively. Experimental and numerical results are presented and discussed in Sections 5 and 6, with conclusions given in Section 7. Boldface letters are used to indicate complex-numbered quantities.

2 Theory

2.1 Small-amplitude power budget equations (linear sound propagation)

A spherical coordinate system (r, θ, φ) is used with origin centred at the face of the echosounder transducer and the coordinate axis (r, 0, 0) coincident with the sound beam axis. θ ∈ [0, π] and φ ∈ [0, 2π] are the polar and azimuthal angles, respectively [18].

Under assumption of small-amplitude signal propagation, the electroacoustic power budget equations used conventionally in fish abundance estimation and species identification are [18, 22]

$${\sigma }_{\mathrm{bs}}\left(r,\theta ,\phi \right)=\frac{{\mathrm{\Pi }}_{\mathrm{R}}^{\mathrm{st}}}{{\mathrm{\Pi }}_{\mathrm{T}}^{\mathrm{st}}}\frac{16{\pi }^2 r^4{e}^{4\alpha r}}{{\lambda }^2 F_{\mathrm{\Pi }} G^2\left(\theta ,\phi \right)},$$(1)

$$s_v\left(r\right)=\frac{{\mathrm{\Pi }}_{\mathrm{R}}^{\mathrm{v}}}{{\mathrm{\Pi }}_{\mathrm{T}}^{\mathrm{v}}}\frac{32{\pi }^2 r^2{e}^{4\alpha r}}{G_0^2\psi {\lambda }^2 c_0{\tau }_{\mathrm{p}} F_{\mathrm{\Pi }}},$$(2)

$$s_a\equiv \underset{r_1}{\overset{r_2}{\int }}{s}_v\left(r\right)\mathrm{dr}.$$(3)

Here, σ_bs (units of m²) is the backscattering cross-section of a single target in the far field of the transducer; s_v (units of m⁻¹) is the volume backscattering coefficient for a thin far-field spherical shell of thickness $\frac{1}{2}{c}_0{\tau }_{\mathrm{p}}$, centred at the transducer front; and s_a (dimensionless) is the area backscattering coefficient for a far-field interval of observation ranges [r₁, r₂].

In equation (1), r is the distance from the transducer to the single target. In equation (2), r is the mid-range to the thin spherical shell. τ_p is the duration of the transmitted voltage signal (τ_p = 1 ms may often be used, corresponding to a spherical shell thickness $\frac{1}{2}{c}_0{\tau }_{\mathrm{p}}\approx 0.75 \mathrm{m}$). α is the acoustic absorption coefficient in seawater, expressed in Np/m. λ = c₀/f is the acoustic wavelength at the measurement frequency f, and c₀ is the small-signal sound velocity of the fluid propagation medium (seawater). The electrical termination factor F_Π accounts for the electrical impedance load at the transducer’s terminals upon signal reception [13].

G(θ, φ) and G₀ (both non-dimensional) are the transducer gain and the axial transducer gain, given as [18]

$$G\left(\theta ,\phi \right)=G_0{\left|{\mathit{B}}_{\mathrm{i}}\left(\theta ,\phi \right)\right|}^2, \hspace{1em}{G}_0=\eta D_0,$$(4)

respectively, where η (non-dimensional) is the electroacoustic conversion efficiency of the transducer. B_i(θ, φ) = P_i(r, θ, φ)/P_i(r, 0, 0) is the beam pattern of the incident sound pressure field, P_i being the incident sound pressure amplitude assuming small-signal (linear) sound propagation.

$$D_0=\frac{4\pi }{{\int }_{4\pi }^{}{\left|{\mathit{B}}_{\mathrm{i}}\left(\theta ,\phi \right)\right|}^2 d\mathrm{\Omega }}$$(5)

is the axial directivity factor [18, 23]. G(θ, φ) and G₀ represent the transducer’s one-way electroacoustic conversion efficiency per unit solid angle, in the (θ, φ) and axial directions, respectively, for lossless sound propagation conditions in the fluid [21].

For reciprocal transducers, having the same beam pattern for transmission and reception, ψ (in steradians, sr) is the transducer’s equivalent two-way beam solid angle [18],

$$\psi ={\int }_{4\pi }^{}{\left|{\mathit{B}}_{\mathrm{i}}\left(\theta ,\phi \right)\right|}^4 d\mathrm{\Omega }=\frac{1}{G_0^2}{\int }_{4\pi }^{}{G}^2\left(\theta ,\phi \right) d\mathrm{\Omega },$$(6)

where dΩ ≡ sinθdθdφ is a solid angle element. ψ represents an effective beam width of the transducer’s intensity field, expressed in terms of a solid angle, accounting for the combined effect of transmission and reception [2, 13, 18, 21].

${\mathrm{\Pi }}_{\mathrm{T}}^{\mathrm{st}}$ and ${\mathrm{\Pi }}_{\mathrm{T}}^{\mathrm{v}}$ are the average electrical transmit signal powers during transmission of single-frequency tone bursts for single-target and volume backscattering measurements, respectively. ${\mathrm{\Pi }}_{\mathrm{R}}^{\mathrm{st}}$ and ${\mathrm{\Pi }}_{\mathrm{R}}^{\mathrm{v}}$ are the corresponding average electrical powers received after two-way sound propagation, backscattering, and two-way electroacoustic transduction [21]. “Average” here refers to time averaging over one cycle of the monochromatic signal waveform with frequency f. In conventional operation ${\mathrm{\Pi }}_{\mathrm{T}}^{\mathrm{st}}={\mathrm{\Pi }}_{\mathrm{T}}^{\mathrm{v}}$ is often used, and F_Π is commonly omitted from equations (1) and (2) [13, 18, 21, 22].

The assumptions underlying equations (1)–(3) are summarized in, e.g., [18–21], with references therein to relevant literature on the subject.

2.2 Finite-amplitude power budget equations (nonlinear sound propagation)

Since the amplitude of the backscattered wave is much weaker than that of the incident wave, finite-amplitude effects are considered most influential in the forward-propagating (incident) sound beam while the backscattered wave is considered to propagate without nonlinear loss [21]. Amendments to equations (1) and (2) to account for finite-amplitude sound propagation effects in measurement of σ_bs and s_v were proposed in [13, 24] and more rigorously derived in [21]. Their form as given in [21] is

$${\sigma }_{\mathrm{bs}}(r,\theta ,\phi )=\frac{{\mathrm{\Pi }}_{\mathrm{R}}^{\mathrm{st}}}{{\mathrm{\Pi }}_{\mathrm{T}}^{\mathrm{st}}}\frac{16{\pi }^2 r^4{e}^{4\alpha r}}{{\lambda }^2 F_{\mathrm{\Pi }} G^2\left(\theta ,\phi \right)}\cdot \frac{1}{{\left|{\mathit{C}}_{\mathrm{i}}^{\mathrm{n}}\left(r\right)\right|}^2\cdot {\left|{\mathit{B}}_{\mathrm{rel}}^{\mathrm{n}}\left(r,\theta ,\phi \right)\right|}^2},$$(7)

$$s_v\left(r\right)=\frac{{\mathrm{\Pi }}_{\mathrm{R}}^{\mathrm{v}}}{{\mathrm{\Pi }}_{\mathrm{T}}^{\mathrm{v}}}\frac{32{\pi }^2 r^2{e}^{4\alpha r}}{G_0^2\psi {\lambda }^2 c_0{\tau }_{\mathrm{p}} F_{\mathrm{\Pi }}}\cdot \frac{1}{{\left|{\mathit{C}}_{\mathrm{i}}^{\mathrm{n}}\left(r\right)\right|}^2\cdot {\psi }_{\mathrm{rel}}^{\mathrm{n}}\left(r\right)}.$$(8)

Three factors in equations (7) and (8) to account for finite-amplitude sound propagation effects. The “axial finite amplitude factor”, defined as [21]

$${\mathit{C}}_{\mathrm{i}}^{\mathrm{n}}\left(r\right)\equiv \frac{{\mathit{P}}_{\mathrm{i}}^{\mathrm{n}}(r,\mathrm{0,0})}{{\mathit{P}}_{\mathrm{i}}\left(r,\mathrm{0,0}\right)},$$(9)

accounts for nonlinear loss on the sound beam axis in the incident complex sound pressure amplitude due to finite-amplitude sound propagation, ${\mathit{P}}_{\mathrm{i}}^{\mathrm{n}}\left(r,\theta ,\phi \right)$. The absolute value ${|\mathit{C}}_{\mathrm{i}}^{\mathrm{n}}\left(r\right)|$ is in the interval 〈0, 1〉.

The “beam pattern finite amplitude factor”, defined as

$${\mathit{B}}_{\mathrm{rel}}^{\mathrm{n}}\left(r,\theta ,\phi \right)\equiv \frac{{\mathit{B}}_{\mathrm{i}}^{\mathrm{n}}\left(r,\theta ,\phi \right)}{{\mathit{B}}_{\mathrm{i}}\left(\theta ,\phi \right)},$$(10)

accounts for finite-amplitude effects on the beam pattern (beam flattening) in the forward-propagating (incident) wave. ${\mathit{B}}_{\mathrm{i}}^{\mathrm{n}}\left(r,\theta ,\phi \right)={\mathit{P}}_{\mathrm{i}}^{\mathrm{n}}\left(r,\theta ,\phi \right)/{\mathit{P}}_{\mathrm{i}}^{\mathrm{n}}\left(r,\mathrm{0,0}\right)$ is the beam pattern of the forward-propagating pressure wave under finite-amplitude sound propagation conditions. $\left|{\mathit{B}}_{\mathrm{rel}}^{\mathrm{n}}\left(r,\theta ,\phi \right)\right|\approx 1$ for small-amplitude (i.e., linear) sound propagation, and $\left|{\mathit{B}}_{\mathrm{rel}}^{\mathrm{n}}\left(r,\theta ,\phi \right)\right|\ge 1$ under finite-amplitude conditions. Since ${\mathit{B}}_{\mathrm{i}}\left(\theta ,\phi \right)$ and ${\mathit{B}}_{\mathrm{i}}^{\mathrm{n}}\left(r,\theta ,\phi \right)$ are both unity on the sound beam axis, ${|\mathit{B}}_{\mathrm{rel}}^{\mathrm{n}}\left(r,\mathrm{0,0}\right)|=1$.

The “beam solid angle finite-amplitude factor” ${\psi }_{\mathrm{rel}}^{\mathrm{n}}\left(r\right)$, defined as [21]

$${\psi }_{\mathrm{rel}}^{\mathrm{n}}\left(r\right)\equiv \frac{{\psi }^{\mathrm{n}}\left(r\right)}{\psi },$$

$${\psi }^{\mathrm{n}}\left(r\right)={\int }_{4\pi }^{}{\left|{\mathit{B}}_{\mathrm{i}}^{\mathrm{n}}\left(r,\theta ,\phi \right)\right|}^2{\left|{\mathit{B}}_{\mathrm{i}}\left(\theta ,\phi \right)\right|}^2 d\mathrm{\Omega }=\frac{1}{G_0^2}{\int }_{4\pi }^{}{G}^2\left(\theta ,\phi \right) {\left|{\mathit{B}}_{\mathrm{rel}}^{\mathrm{n}}\left(r,\theta ,\phi \right)\right|}^2 d\mathrm{\Omega },$$(11)

accounts for a change in the shape of the radiated sound beam, while the beampattern for reception of the backscattered signal remains unchanged. ψⁿ(r) is the equivalent two-way beam solid angle under conditions of finite-amplitude propagation of the forward-propagating (incident) sound.

The finite-amplitude factors ${\mathit{C}}_{\mathrm{i}}^{\mathrm{n}}\left(r\right)$, ${\mathit{B}}_{\mathrm{rel}}^{\mathrm{n}}\left(r,\theta ,\phi \right)$, and ${\psi }_{\mathrm{rel}}^{\mathrm{n}}\left(r\right)$ can be obtained from hydrophone measurements for relevant echosounder transducers [21].

2.3 Influence of finite-amplitude effects in echosounder calibration

Equations (7) and (8) can be amended to account for possible nonlinear loss when the echosounder is calibrated against a reference target with known backscattering cross-section σ_bs,ref. The reference target is positioned on the sound beam axis at a range r = r_cal in the far field of the transducer (see, e.g., [6, 22]). Calibration is commonly viewed as determining the axial transducer gain G₀ [6] so that equation (1) is satisfied with σ_bs = σ_bs,ref and (r, θ, φ) = (r_cal, 0, 0). The value for G₀ obtained by the calibration measurement is here denoted $G_0^{\mathrm{meas}}$, given by

$${\left(G_0^{\mathrm{meas}}\right)}^2=\frac{{\mathrm{\Pi }}_{\mathrm{R}}^{\mathrm{st}} 16{\pi }^2}{{\mathrm{\Pi }}_{\mathrm{T}}^{\mathrm{st}} {\lambda }^2{F}_{\mathrm{\Pi }} {\sigma }_{\mathrm{bs},\mathrm{ref}}}\cdot r_{\mathrm{cal}}^4 e^{4\alpha r_{\mathrm{cal}}}.$$(12)

Subsequent measurements of single-target backscattering are carried out using $G^{\mathrm{meas}}\left(\theta ,\phi \right)=G_0^{\mathrm{meas}}\cdot {\left|{\mathit{B}}_{\mathrm{i}}\left(\theta ,\phi \right)\right|}^2$ for G(θ, φ) in equation (1) or (7), and volume backscattering measurements are made using $G_0^{\mathrm{meas}}$ for G₀ in equation (2) or (8).

If finite-amplitude sound propagation effects are significant upon calibration, the propagation loss becomes higher than accounted for by equation (12), leading to error in $G_0^{\mathrm{meas}}$. The correct value for G₀ in presence of nonlinear loss is found by applying equation (7) in place of equation (1) to the calibration measurement:

$$G_0^2=\frac{{\mathrm{\Pi }}_{\mathrm{R}}^{\mathrm{st}} 16{\pi }^2}{{\mathrm{\Pi }}_{\mathrm{T}}^{\mathrm{st}} {\lambda }^2{F}_{\mathrm{\Pi }} {\sigma }_{\mathrm{bs},\mathrm{ref}}}\cdot \frac{r_{\mathrm{cal}}^4 e^{4\alpha r_{\mathrm{cal}}}}{{\left|{\mathit{C}}_{\mathrm{i},\mathrm{cal}}^{\mathrm{n}}\left(r_{\mathrm{cal}}\right)\right|}^2}=\frac{{\left(G_0^{\mathrm{meas}}\right)}^2}{{\left|{\mathit{C}}_{\mathrm{i},\mathrm{cal}}^{\mathrm{n}}\left(r_{\mathrm{cal}}\right)\right|}^2}.$$(13)

The subscript “cal” indicates that the axial finite amplitude factor ${\mathit{C}}_{\mathrm{i},\mathrm{cal}}^{\mathrm{n}}\left(r_{\mathrm{cal}}\right)$ is evaluated for the source level and propagation medium upon calibration.

Measurements based on a calibration for which linear sound propagation has been inaccurately assumed, thus inadvertently using a calibration result $G_0^{\mathrm{meas}}$ that is different from the actual transducer gain G₀, can be compensated for the aggregated effects of nonlinear loss by substituting equation (13) for G₀ in equations (7) and (8), yielding

$${\sigma }_{\mathrm{bs}}=\frac{{\mathrm{\Pi }}_{\mathrm{R}}^{\mathrm{st}}}{{\mathrm{\Pi }}_{\mathrm{T}}^{\mathrm{st}}}\frac{16{\pi }^2 r^4{e}^{4\alpha r}}{{\lambda }^2{F}_{\mathrm{\Pi }}{\left(G^{\mathrm{meas}}\left(\theta ,\phi \right)\right)}^2}\cdot \frac{{\left|{\mathit{C}}_{\mathrm{i},\mathrm{cal}}^{\mathrm{n}}\left(r_{\mathrm{cal}}\right)\right|}^2}{{\left|{\mathit{C}}_{\mathrm{i}}^{\mathrm{n}}\left(r\right)\right|}^2 {\left|{\mathit{B}}_{\mathrm{rel}}^{\mathrm{n}}\left(r,\theta ,\phi \right)\right|}^2},$$(14)

$$s_v=\frac{{\mathrm{\Pi }}_{\mathrm{R}}^{\mathrm{v}}}{{\mathrm{\Pi }}_{\mathrm{T}}^{\mathrm{v}}}\frac{32{\pi }^2 r^2{e}^{4\alpha r}}{{\left(G_0^{\mathrm{meas}}\right)}^2 \psi {\lambda }^2 c_0{\tau }_{\mathrm{p}} F_{\mathrm{\Pi }}}\cdot \frac{{\left|{\mathit{C}}_{\mathrm{i},\mathrm{cal}}^{\mathrm{n}}\left(r_{\mathrm{cal}}\right)\right|}^2}{{\left|{\mathit{C}}_{\mathrm{i}}^{\mathrm{n}}\left(r\right)\right|}^2 {\psi }_{\mathrm{rel}}^{\mathrm{n}}\left(r\right)}.$$(15)

The transducer gain G^meas (θ, φ) in equation (14), based on the calibration measurement, and assuming linear sound propagation in both calibration and field operation, is

$$G^{\mathrm{meas}}\left(\theta ,\phi \right)=G_0^{\mathrm{meas}}\cdot {\left|{\mathit{B}}_{\mathrm{i}}\left(\theta ,\phi \right)\right|}^2=G_0 \left|{\mathit{C}}_{\mathrm{i},\mathrm{cal}}^{\mathrm{n}}\left(r_{\mathrm{cal}}\right)\right| {\left|{\mathit{B}}_{\mathrm{i}}\left(\theta ,\phi \right)\right|}^2.$$(16)

Equations (14) and (15) account for possible finite-amplitude effects upon both calibration and subsequent field measurements of σ_bs and s_v, when the calibration has been performed assuming linear sound propagation.

2.4 Finite-amplitude sound propagation modelling

The Khokhlov–Zabolotskaya–Kuznetsov (“KZK”) equation [9, 25, 26] is used here to model nonlinear sound propagation and thus obtain numerical values for ${\mathit{C}}_{\mathrm{i}}^{\mathrm{n}}\left(r\right)$, ${\mathit{C}}_{\mathrm{i},\mathrm{cal}}^{\mathrm{n}}\left(r_{\mathrm{cal}}\right)$, and ${\psi }_{\mathrm{rel}}^{\mathrm{n}}\left(r\right)$ for use in equations (14) and (15). The KZK equation is parabolic and suitable for forward-marching solution by finite differences to the second order in nonlinearity, with validity inside the paraxial region of a sound beam.

A Cartesian coordinate system (x, y, z) is defined with origin and z-axis coinciding with the origin and axis of the previously defined spherical coordinate system. The KZK equation can be written, for the sound pressure p = P − P₀ [9],

$$\left(\frac{{\partial }^2}{\partial z \partial t^{\mathrm{\prime }}}-\frac{c_0}{2}\left(\frac{{\partial }^2}{\partial x^2}+\frac{{\partial }^2}{\partial y^2}\right)-\frac{\delta }{2c_0^3}\frac{{\partial }^3}{\partial t^{\mathrm{\prime }3}}\right)p=\frac{\beta }{2{\rho }_0{c}_0^3}\frac{{\partial }^2{p}^2}{\partial t^{\mathrm{\prime }2}}.$$(17)

P₀ and P are the hydrostatic and total pressures, respectively, and t′ = t − z/c₀ is the retarded time assuming constant velocity c₀. β is the coefficient of nonlinearity [9], and ρ₀ is the ambient density of the fluid propagation medium. Classical absorption of a thermo-viscous fluid is accounted for by the diffusivity $\delta =c_0^3{\alpha }_2/2{\pi }^2$, where the absorption coefficient α = α₂ f² and α₂ is frequency-independent [23]. α₂ can be adjusted to specify correct absorption at the operating frequency. The chemical relaxation mechanisms in seawater [27] however make the frequency-squared power law inaccurate for the harmonic frequency components generated through finite-amplitude signal distortion.

Equation (17) is solved in the frequency domain using the “Bergen Code” finite-difference algorithm [29]. The frequency domain approach allows insertion of the Francois–Garrison seawater absorption coefficient [28] for each harmonic component. This yields somewhat improved agreement with the experimental data. 50 harmonics are included in the calculations.

The KZK equation and the Bergen Code have been tested extensively against experimental measurements over the recent 3–4 decades. For an overview of such work in published literature, see, e.g., [9, 13]. Favourable agreement between simulations and experiments has been reported in general, over a wide range of simulation parameters.

3 Measurements

In-situ acoustic measurements on fish schools were performed from the research vessel G.O. Sars in the North Sea [13]. A Simrad EK60 echo sounder system was operated for simultaneous measurements at 38 kHz, 120 kHz, and 200 kHz acoustic frequencies. The three transducers are designed with the same far-field beam width 2θ_−3dB = 7° (at main thickness resonance) and mounted close together on a protrudable instrument keel. Effects of air bubbles were minimized by extending the instrument keel 3 m below the hull to 8.5 m depth. The overlap between the 38 kHz sound beam and each of the 120 kHz and 200 kHz sound beams has been estimated to 85% at 30 m depth, assuming that the beams were parallel [17]. Signal transmission was simultaneous between the frequencies so that the ship’s motion did not increase the offset between the ensonified volumes. The echosounder system was calibrated before the survey cruise with reference targets at range r_cal = 22 m (about 30 m depth), using electrical power settings of 250 W for the 120 kHz echosounder and 120 W for the 200 kHz echosounder. These are in the following referred to as the “low” power settings, while 1000 W is referred to as the “high” power setting for both operating frequencies. The measurements at 38 kHz were made with constant power setting of 2000 W. Nonlinear loss at this reference frequency is thus constant throughout the measurement series, while factors such as fish density or school size cause similar variability in the measurements at all frequencies.

The vessel crossed repeatedly over an area with several schools while the transmit power for the 120 kHz and 200 kHz transducers were alternated between the “low” and “high” settings for each crossing (Fig. 1). The schools extended predominantly from 20 m to 110 m depth, with more than 90% of the backscattered energy from depths between 40 m and 100 m. The depth distribution was similar throughout the measurement series. There were essentially three schools (purple blobs in Fig. 1) in the area, passed in repeated straight cruise tracks (blue curves in Fig. 1). The schools were remarkably stationary during the 11 h from 13:30 to 00:30 UTC on 6–7 November 2004. The cruise tracks covered 42 nautical miles (nmi) and 70,500 pings over 11 h. The ship followed the schools for another 2.5 h as the schools started to move slowly. Data from this final period were not included in the analysis. A pelagic trawl haul was carried out at depth 70 m–90 m and indicated 95% Atlantic mackerel (Scomber scombrus L.).

Figure 1 Cruise track (blue curves) for the acoustic survey and observed fish schools (purple blobs). The axes indicate longitude and latitude in decimal degrees. The values in the legend are the nautical area scattering coefficient sA [nmi−1] [33]. Note that the cruise tracks do not overlap completely. Each school could be crossed straight over in one leg and on a side in the next, even if the schools did not move.

Ambient noise was removed from the acoustic data using the method described in [30]. Acoustic data were analysed by means of Bergen Echo Integrator (BEI [31]) during the survey, and in a more thorough analysis immediately after the survey. Schools were encircled, species identified, and the result was stored to a database. The successor of BEI, Large Scale Survey System (LSSS [32]), was used to repeat the analysis and produce the results presented here. A total of 42 passes over schools were extracted from the acoustic survey data. The relative frequency response identified the species as Atlantic mackerel.

The measured s_v(r) were integrated with a threshold of S_v = 10 log₁₀(s_v/1 m⁻¹) = −82 dB re 1 m⁻¹ [33] at all frequencies. The removal of ambient noise reduced the average backscatter s_v by 0.15%–0.25% depending on frequency. Each of the encircled schools contained 350–700 pings over 0.2 nmi–0.4 nmi path lengths.

The common variability between measurements at 38 kHz, 120 kHz, and 200 kHz cancels partially in the relative frequency response [34],

$$r_a\left(f\right)\equiv \frac{s_a\left(f\right)}{s_a\left(38 \mathrm{kHz}\right)} .$$(18)

r_a(f) is thus a measurement result with less point-to-point variability, but with the same dependence upon transmit power at f = 120 kHz and 200 kHz, as the respective s_a(f). The experimental results for s_a(f) are calculated according to equations (2) and (3).

4 Numerical calculations

The transmitting echosounder transducer is modelled as a circular planar piston, vibrating with uniform axial velocity in a rigid baffle of infinite extent. Within the parabolic approximation, the velocity amplitude U_s of the piston can be expressed via a source condition pressure amplitude ${\mathit{P}}_{\mathrm{s}}={\rho }_0{c}_0\cdot {\mathit{U}}_{\mathrm{s}}$ [9, 29]. The source pressure amplitude and effective piston radius a were calculated from far-field hydrophone measurements with low amplitudes, using echosounder transducers of the same types as for the oceanic survey [13]. The source pressure amplitudes were scaled by measured voltage amplitudes across the electrical terminals of the transducers upon signal transmission [13].

Table 1 gives the parameters used in numerical simulations to compare with the in-situ survey measurements. The propagation medium parameters of seawater [28, 35, 36, 37] are calculated based on CTD (conductivity, temperature, and pressure) sampling performed regularly during the survey cruise, averaged over depths from 9 m to 120 m. The measured water temperature was between 9.2 °C and 10.5 °C, and the sound speed was between 1490 m/s and 1492 m/s throughout the water column [35]. The calculated coefficient of nonlinearity [36] was between 3.58 and 3.59.

Numerical calculations were made to produce pairs of corresponding simulated sound pressure fields ${\mathit{P}}_{\mathrm{i}}^{\mathrm{n}}(r,\theta ,\phi )$ of finite-amplitude sound propagation (corresponding to the “low” and “high” experimental power settings) and ${\mathit{P}}_{\mathrm{lin}}(r,\theta ,\phi )$ where finite-amplitude effects were avoided by setting the source condition pressure amplitude to $\left|{\mathit{P}}_{\mathrm{s}}\right|=\left|{\mathit{P}}_{\mathrm{s},\mathrm{lin}}\right|=1 \mathrm{Pa}$. From equation (9), simulated values for the axial finite-amplitude factor are calculated as

$${\mathit{C}}_{\mathrm{i}}^{\mathrm{n}}\left(r\right)=\frac{{\mathit{P}}_{\mathrm{s},\mathrm{lin}}}{{\mathit{P}}_{\mathrm{s}}}\cdot \frac{{\mathit{P}}_{\mathrm{i}}^{\mathrm{n}}\left(r,\mathrm{0,0}\right)}{{\mathit{P}}_{\mathrm{lin}}\left(r,\mathrm{0,0}\right)},$$(19)

$${\mathit{C}}_{\mathrm{i},\mathrm{cal}}^{\mathrm{n}}\left(r_{\mathrm{cal}}\right)=\frac{{\mathit{P}}_{\mathrm{s},\mathrm{lin},\mathrm{cal}}}{{\mathit{P}}_{\mathrm{s},\mathrm{cal}}}\cdot \frac{{\mathit{P}}_{\mathrm{i},\mathrm{cal}}^{\mathrm{n}}\left(r_{\mathrm{cal}},\mathrm{0,0}\right)}{{\mathit{P}}_{\mathrm{lin},\mathrm{cal}}\left(r_{\mathrm{cal}},\mathrm{0,0}\right)},$$(20)

where the subscript “cal” indicates the sound propagation conditions for the case of echosounder calibration.

The beam pattern finite-amplitude factor is calculated from equation (10) and the definitions of ${\mathit{B}}_{\mathrm{i}}$ and ${\mathit{B}}_{\mathrm{i}}^{\mathrm{n}}$, yielding

$${\mathit{B}}_{\mathrm{rel}}^{\mathrm{n}}\left(r,\theta ,\phi \right)=\frac{{\mathit{P}}_{\mathrm{lin}}(r,\mathrm{0,0})}{{\mathit{P}}_{\mathrm{i}}^{\mathrm{n}}(r,\mathrm{0,0})}\cdot \frac{{\mathit{P}}_{\mathrm{i}}^{\mathrm{n}}(r,\theta ,\phi )}{{\mathit{P}}_{\mathrm{lin}}(r,\theta ,\phi )}.$$(21)

The beam solid angle finite-amplitude factor ${\psi }_{\mathrm{rel}}^{\mathrm{n}}\left(r\right)$ is calculated by insertion of equation (21) into equation (11). Since the KZK equation is valid under the paraxial approximation and significant finite-amplitude effects are expected only in the main lobe, the integrals in equation (11) are evaluated for θ up to an integration limit θ_lim = 20°. Parameter variation shows that the reported results are unchanged with θ_lim down to 6°, while the first minimum of the beam pattern is at approximately θ = 8°. For rotational symmetry with respect to φ one has dΩ = 2π sinθ dθ, and the expression for the beam solid angle finite-amplitude factor becomes

$${\psi }_{\mathrm{rel}}^{\mathrm{n}}\left(r\right)\approx \frac{\underset{0}{\overset{{\theta }_{\lim }}{\int }}{G}^2\left(\theta ,\phi \right)\cdot {\left|{\mathit{B}}_{\mathrm{rel}}^{\mathrm{n}}\left(r,\theta ,\phi \right)\right|}^2\sin \theta \mathrm{d\theta }}{\underset{0}{\overset{{\theta }_{\lim }}{\int }}{G}^2\left(\theta ,\phi \right) \sin \theta \mathrm{d\theta }}.$$(22)

Since the nonlinear loss is range dependent, finite-amplitude effects in measurement of r_a(f) should be accounted for by using equation (15) to calculate s_v(r), before integrating over a range interval [r₁, r₂] to obtain s_a according to equation (3). The experimental measurements of r_a(f) were however processed using standard methods as described in Section 3, and therefore no correction for nonlinear loss was made prior to integration.

For clarity in the following, the actual values of s_v, s_a, and r_a are identified with superscript “correct”, while the corresponding (erroneous) measurement results, assuming linear sound propagation, are identified with superscript “meas”. $s_v^{\mathrm{meas}}$ is thus the measurement result obtained using equation (2) (assuming linear sound propagation in the survey measurement), with $G_0^{\mathrm{meas}}$ used in place of $G_0$ (assuming linear sound propagation in echosounder calibration), given as

$$s_v^{\mathrm{meas}}\left(r\right)=\frac{{\mathrm{\Pi }}_{\mathrm{R}}^{\mathrm{v}}}{{\mathrm{\Pi }}_{\mathrm{T}}^{\mathrm{v}}}\frac{32{\pi }^2 r^2{e}^{4\alpha r}}{{\left(G_0^{\mathrm{meas}}\right)}^2\psi {\lambda }^2 c_0{\tau }_{\mathrm{p}} F_{\mathrm{\Pi }}}.$$(23)

From equations (15) and (23) one can calculate a “finite-amplitude correction factor” for the total effect of nonlinear loss,

$$\frac{s_v^{\mathrm{correct}}\left(r\right)}{s_v^{\mathrm{meas}}\left(r\right)}=\frac{{\left|{\mathit{C}}_{\mathrm{i},\mathrm{cal}}^{\mathrm{n}}\left(r_{\mathrm{cal}}\right)\right|}^2}{{\left|{\mathit{C}}_{\mathrm{i}}^{\mathrm{n}}\left(r\right)\right|}^2\cdot {\psi }_{\mathrm{rel}}^{\mathrm{n}}\left(r\right)}.$$(24)

For comparison with the field data on the form of measured $r_a\left(f\right)$ under conditions of alternating power settings, the quantity of interest is $r_a^{\mathrm{meas}, \mathrm{high}}/r_a^{\mathrm{meas}, \mathrm{low}}$. By noting that ${\mathit{C}}_{\mathrm{i}}^{\mathrm{n}}\left(r\right)$ and ${\psi }_{\mathrm{rel}}^{\mathrm{n}}\left(r\right)$ are both approximately constant at long ranges, for given frequency and source level [21], one can approximate ${\left|{\mathit{C}}_{\mathrm{i}}^{\mathrm{n}}\left(r\right)\right|}^2\cdot {\psi }_{\mathrm{rel}}^{\mathrm{n}}\left(r\right)$ by its average over the full observation range interval, or for a thinner layer,$\left[r_1,r_2\right].$ From equations (3) and (24), one finds

$$\frac{s_a^{\mathrm{correct}}}{s_a^{\mathrm{meas}}}=\frac{\underset{r_1}{\overset{r_2}{\int }}{s}_v^{\mathrm{correct}}\left(r\right) \mathrm{dr}}{\underset{r_1}{\overset{r_2}{\int }}{s}_v^{\mathrm{meas}}\left(r\right) \mathrm{dr}}=\frac{{\left|{\mathit{C}}_{\mathrm{i},\mathrm{cal}}^{\mathrm{n}}\left(r_{\mathrm{cal}}\right)\right|}^2\underset{r_1}{\overset{r_2}{\int }}{s}_v^{\mathrm{correct}}\left(r\right) \mathrm{dr}}{\underset{r_1}{\overset{r_2}{\int }}{\left|{\mathit{C}}_i^{\mathrm{n}}\left(r\right)\right|}^2{\psi }_{\mathrm{rel}}^{\mathrm{n}}\left(r\right) s_v^{\mathrm{correct}}\left(r\right) \mathrm{dr}}\approx \frac{{\left|{\mathit{C}}_{\mathrm{i},\mathrm{cal}}^{\mathrm{n}}\left(r_{\mathrm{cal}}\right)\right|}^2 \underset{r_1}{\overset{r_2}{\int }}{s}_v^{\mathrm{correct}}\left(r\right) \mathrm{dr}}{\langle {\left|{\mathit{C}}_{\mathrm{i}}^{\mathrm{n}}\left(r\right)\right|}^2\cdot {\psi }_{\mathrm{rel}}^{\mathrm{n}}\left(r\right)\rangle \underset{r_1}{\overset{r_2}{\int }} s_v^{\mathrm{correct}}\left(r\right) \mathrm{dr}}=\frac{{\left|{\mathit{C}}_{\mathrm{i},\mathrm{cal}}^{\mathrm{n}}\left(r_{\mathrm{cal}}\right)\right|}^2}{\langle {\left|{\mathit{C}}_{\mathrm{i}}^{\mathrm{n}}\left(r\right)\right|}^2\cdot {\psi }_{\mathrm{rel}}^{\mathrm{n}}\left(r\right)\rangle },$$(25)

for $r_1$ and $r_2$ in the far field. Calibration was performed at one range $r_{\mathrm{cal}}$, using the “low” power setting, so that a single value of ${\left|{\mathit{C}}_{\mathrm{i},\mathrm{cal}}^{\mathrm{n}}\left(r_{\mathrm{cal}}\right)\right|}^2$applies to each measurement frequency. The average product $\langle {\left|{\mathit{C}}_{\mathrm{i}}^{\mathrm{n}}\left(r\right)\right|}^2\cdot {\psi }_{\mathrm{rel}}^{\mathrm{n}}\left(r\right)\rangle$ for a range layer, on the other hand, depends on the power setting during field measurements, so that one can write

$$s_a^{\mathrm{meas},\mathrm{high}}\approx \frac{s_a^{\mathrm{correct}}{\cdot \langle {\left|{\mathit{C}}_{\mathrm{i}}^{\mathrm{n}}\left(r\right)\right|}^2\cdot {\psi }_{\mathrm{rel}}^{\mathrm{n}}\left(r\right)\rangle }_{\mathrm{high}}}{{\left|{\mathit{C}}_{\mathrm{i},\mathrm{cal}}^{\mathrm{n}}\left(r_{\mathrm{cal}}\right)\right|}^2},$$(26)

$$s_a^{\mathrm{meas},\mathrm{low}}\approx \frac{s_a^{\mathrm{correct}}{\cdot \langle {\left|{\mathit{C}}_{\mathrm{i}}^{\mathrm{n}}\left(r\right)\right|}^2\cdot {\psi }_{\mathrm{rel}}^{\mathrm{n}}\left(r\right)\rangle }_{\mathrm{low}}}{{\left|{\mathit{C}}_{\mathrm{i},\mathrm{cal}}^{\mathrm{n}}\left(r_{\mathrm{cal}}\right)\right|}^2}.$$(27)

Since s_a (38 kHz) does not depend on the power settings for the 120 kHz and 200 kHz frequencies, equations (18), (26), and (27) can be combined to

$$\frac{r_a^{\mathrm{meas},\mathrm{high}}}{r_a^{\mathrm{meas},\mathrm{low}}}\approx \frac{s_a^{\mathrm{meas},\mathrm{high}}}{s_a^{\mathrm{meas},\mathrm{low}}}\approx \frac{{\langle {\left|{\mathit{C}}_{\mathrm{i}}^{\mathrm{n}}\left(r\right)\right|}^2\cdot {\psi }_{\mathrm{rel}}^{\mathrm{n}}\left(r\right)\rangle }_{\mathrm{high}}}{{\langle {\left|{\mathit{C}}_{\mathrm{i}}^{\mathrm{n}}\left(r\right)\right|}^2\cdot {\psi }_{\mathrm{rel}}^{\mathrm{n}}\left(r\right)\rangle }_{\mathrm{low}}}.$$(28)

5 Results

Figure 2 shows measurement results $r_a^{\mathrm{meas}}\left(f\right)$ using the “low” and “high” power settings for the two measurement frequencies 120 kHz and 200 kHz (relative to 38 kHz), used in the survey measurements on Atlantic mackerel schools, over the observation depth interval [30 m, 110 m]. Results are shown on a logarithmic (dB) scale. The calculated average measurement results for each power setting, 〈r_a(f)〉, are shown using dashed and full lines for the “low” and “high” power settings, respectively. The ratios between 〈r_a(f)〉 measured using the two power settings, $\langle r_a^{\mathrm{meas},\mathrm{high}}\rangle /\langle r_a^{\mathrm{meas},\mathrm{low}}\rangle$, are 0.75 (−1.2 dB) for the 120 kHz measurements, and 0.48 (−3.2 dB) for the measurements made at 200 kHz.

Figure 2 Measured relative frequency response, ra(f), obtained from oceanic in-situ acoustical survey measurements of sv on Atlantic mackerel schools. Data from within an observation depth interval of 30 m to 110 m are included. Each point (x or o) is based on measurements from a whole school. (a) 120 kHz measurement frequency, using 250 W and 1000 W electrical power settings in survey operation (alternating), and 250 W during calibration. (b) 200 kHz measurement frequency, using 120 W and 1000 W electrical power settings in survey operation (alternating), and 120 W during calibration.

Figure 3a shows numerical calculations of ${\left({\left|{\mathit{C}}_{\mathrm{i}}^{\mathrm{n}}\left(r\right)\right|}^2 {\psi }_{\mathrm{rel}}^{\mathrm{n}}\left(r\right)\right)}^{-1}$ for the “low” and “high” power settings. This represents the effect of nonlinear loss on volume backscattering measurement (Eq. (8)). Results are shown on a logarithmic (dB) scale vs. logarithmic range. For the 120 kHz measurements with the “low” power setting, the nonlinear loss effect grows from about −0.1 dB towards −0.3 dB (factors of 0.98–0.94 in natural units) over the ranges 10 m–300 m. With the “high” power setting the nonlinear loss effect goes from −0.4 dB to −1.0 dB (factors of 0.90–0.78 in natural units) over the same ranges.

Figure 3 (a) Numerically calculated values for ${\left({\left|{\mathit{C}}_{\mathrm{i}}^{\mathrm{n}}\left(r\right)\right|}^2\cdot {\psi }_{\mathrm{rel}}^{\mathrm{n}}\left(r\right)\right)}^{-1}$, representing the integrated nonlinear loss over the sound beam and thus the nonlinear loss in the sv measurement. The parameters in Table 1 were used for the calculations. (b) Numerical and experimental values for the ratio $r_a^{\mathrm{meas},\mathrm{high}}/r_a^{\mathrm{meas},\mathrm{low}}$ between volume backscattering measurements using the “high” and “low” power settings. The experimental data are divided in depth bins with 20 m width, indicating mean values and standard deviations for the 11-hour measurement series.

For the 200 kHz operating frequency with the “low” power setting, the nonlinear loss effect goes from −0.2 dB at 10 m range to −0.5 dB at 300 m (factors of 0.95 and 0.89, respectively). The effect increases to −2.1 dB and −3.5 dB (factors of 0.61 and 0.44 in natural units) at the same ranges when the “high” power setting is used.

Figure 3b compares the numerically calculated $r_a^{\mathrm{meas},\mathrm{high}}/r_a^{\mathrm{meas},\mathrm{low}}$ as a function of depth, to the measurement results separated into 20 m depth bins. The measurement results are shown as mean values and standard deviations at the centre depth of each bin. The numerical and experimental data agree within the experimental uncertainty.

Table 2 compares the measurements and numerical results for the full observation range interval. The first two rows give the average values $\langle {\left|C_{\mathrm{i}}^{\mathrm{n}}\left(r\right)\right|}^2 {\psi }_{\mathrm{rel}}^{\mathrm{n}}\left(r\right)\rangle$ calculated numerically using equations (19) and (22), for a depth interval [40 m, 100 m], for the “high” and “low” power settings, respectively. The third row in Table 2 gives the simulated $r_a^{\mathrm{meas},\mathrm{high}}/r_a^{\mathrm{meas},\mathrm{low}}$ from equation (28) and the first two rows of the table. The corresponding in-situ measured values $\langle r_a^{\mathrm{meas},\mathrm{high}}\rangle /\langle r_a^{\mathrm{meas},\mathrm{low}}\rangle$, cf. Figure 2, with relative standard deviations, are given in the next two rows. The numbers of datasets exported from LSSS and processed for mean value and standard deviation are given in the last line. The simulated $r_a^{\mathrm{meas},\mathrm{high}}/r_a^{\mathrm{meas},\mathrm{low}}$ at “high” and “low” power settings agree with the measured $\langle r_a^{\mathrm{meas},\mathrm{high}}\rangle /\langle r_a^{\mathrm{meas},\mathrm{low}}\rangle$ to within 0.7 dB.

Numerically calculated values for the correction factor $s_v^{\mathrm{correct}}∕s_v^{\mathrm{meas}}$ (Eq. (24)) are shown as an example in Figure 4. Results for the 1000 W power setting and the parameters in Table 1 are used for both calibration and subsequent measurements of s_v(r). Values for five ranges r_cal to the calibration reference target have been included to visualize the effect of nonlinear loss during calibration on the total measurement error. Due to the flattening of the main lobe, s_v(r) is always overestimated at the calibration range when nonlinear loss is significant but not accounted for. Further finite-amplitude sound propagation beyond the calibration range reduces the total measurement error until it may pass unity and the error from calibration is outweighed.

Figure 4 Numerically calculated example values for the finite-amplitude correction factor $s_v^{\mathrm{correct}}∕s_v^{\mathrm{meas}}$, with sound propagation parameters from Table 1 and 1000 W power setting for both calibration at range rcal and sv measurement at range r. The calibration range is indicated with a mark on each curve. Logarithmic scale, with corresponding ratios shown on the right-hand ordinate axis. (a) 120 kHz, (b) 200 kHz.

The same numerical results are used in Figure 5 to illustrate the combined result of axial nonlinear loss, flattening of the main lobe, and echosounder calibration. Dashed curves indicate the two-way intensity beam pattern |B_i(θ, φ)|⁴ for the 200 kHz echosounder without accounting for nonlinear loss, cf. equation (4). The full-drawn curves show the same two-way beam pattern divided by the ratio

$${\frac{{\sigma }_{\mathrm{bs}}^{\mathrm{correct}}}{{\sigma }_{\mathrm{bs}}^{\mathrm{meas}}}|}_{\left(r,\theta ,\phi \right)}=\frac{{\left|{\mathit{C}}_{\mathrm{i},\mathrm{cal}}^{\mathrm{n}}\left(r_{\mathrm{cal}}\right)\right|}^2}{{\left|{\mathit{C}}_{\mathrm{i}}^{\mathrm{n}}\left(r\right)\right|}^2 {\left|{\mathit{B}}_{\mathrm{rel}}^{\mathrm{n}}\left(r,\theta ,\phi \right)\right|}^2},$$(29)

where ${\sigma }_{\mathrm{bs}}^{\mathrm{meas}}$ is calculated using equation (1) and ${\sigma }_{\mathrm{bs}}^{\mathrm{correct}}$ is from equation (14) using the 1000 W power setting. Echosounder calibration is here assumed to take place with r_cal = 14 m. The leftmost plot assumes linear sound propagation upon calibration, i.e., ${\mathit{C}}_{\mathrm{i},\mathrm{cal}}^{\mathrm{n}}\left(r_{\mathrm{cal}}\right)=1$. The middle plot illustrates finite-amplitude effects during calibration, where the echosounder is calibrated to measure σ_bs correctly on the sound beam axis. The flattening of the main lobe leads to overestimated measurements elsewhere at the same range and thus for s_v. The rightmost plot illustrates the effect of nonlinear loss at 300 m range after the same calibration at r_cal = 14 m, visualizing how further nonlinear loss between ranges r_cal and r > r_cal can counteract the overestimate of s_v experienced at the calibration range.

Figure 5 Illustration of nonlinear loss in the sound beam, using the numerical results for 200 kHz, 1000 W. Left: Two-way intensity beam pattern at 14 m range, with and without correction for nonlinear loss in the forward propagating wave. Centre: Two-way intensity beam pattern at 14 m range, corrected for nonlinear loss upon both measurement and calibration at the same range. Right: Two-way intensity beam pattern at r = 300 m range, corrected for nonlinear loss upon measurement and upon calibration at rcal = 14 m range.

6 Discussion

Systematic and significant differences are observed between the in-situ r_a measurements obtained using “low” and “high” available power settings, for both of the studied operating frequencies. Fair agreement, within 0.5 dB, was found with simulations of finite-amplitude propagation in the forward-propagating sound beam. This is taken to indicate that the observed amplitude dependency in backscattering measurements is caused predominantly by finite-amplitude sound propagation effects in seawater. This result confirms earlier studies in which numerical data were compared to hydrophone measurements in freshwater [10, 13, 14] and seawater [13], and echosounder measurements on reference targets (calibration spheres) in seawater [13]. The systematic effects are of a magnitude that can cause significant errors in quantitative echosounder measurements at both 120 kHz and 200 kHz measurement frequencies. The 1000 W electrical power setting used here as “high” has been standard in the past, e.g., in the widely used echosounder Simrad EK500 with 120 kHz and 200 kHz transducers with 7° opening angles (2θ_−3dB). Figure 4 shows examples of finite-amplitude correction factors that could be used in such cases.

Improved quantitative agreement may be achieved by accounting more accurately for the seawater and operational parameters. Such parameters include the characteristics of the sound source, measurement frequency, source level, range, hydrostatic pressure (depth), and seawater parameters such as the parameter of nonlinearity and the absorption coefficient. The model uncertainty is also important. Factors of relevance in this context are discussed in the following.

Further field data for the same and other relevant operating frequencies would strengthen the experimental basis for quantifying measurement errors across echosounder systems and environmental conditions. Target strength measurements using calibration spheres, comparing results from multiple power settings, have shown fair agreement with simulations on and off the sound beam axis [13]. Absolute sound pressure measurements have been made using hydrophones in freshwater and seawater for distances up to 9 m–12 m [10, 13, 14]. Additional hydrophone measurements of finite-amplitude sound propagation in seawater at longer distances would be useful as a basis for further work if sufficient positioning and acoustic measurement accuracy can be achieved.

With respect to source parameters, scientific echosounders use advanced transducers designed with low sidelobe levels. Their beam patterns are different from that of a uniform circular piston with the same −3 dB beam width, in particular at large angles off the sound beam axis. The paraxial approximation also affects the validity of the numerical calculations at large angles. The importance of this discrepancy is limited since significant nonlinear sound propagation occurs predominantly in the central main lobe in the present application. Some work has been done on the use of alternative source conditions for the simulations [13], and further study can still be useful. Considering the total body of experimental and numerical uncertainties, a possible outcome of such work may however be that the generality of a simple piston source condition outweighs the achievable benefit of improved agreement with specific transducers.

For use of the numerical solution algorithm [29] in the context of the present work, calculations of the factors ${\mathit{C}}_{\mathrm{i}}^{\mathrm{n}}\left(r\right)$, ${\mathit{B}}_{\mathrm{rel}}^{\mathrm{n}}\left(r,\theta ,\phi \right)$, and ${\psi }_{\mathrm{rel}}^{\mathrm{n}}\left(r\right)$ are based mainly on the pressure amplitude ratio ${\mathit{P}}_{\mathrm{i}}^{\mathrm{n}}(r,\theta ,\phi )/{\mathit{P}}_{\mathrm{i}}^{}(r,\theta ,\phi )$. Their uncertainties are expected to be smaller than for simulations of absolute sound pressure amplitudes, since some systematic model uncertainties cancel each other, partially or completely, in the ratios between the numerical results.

The sensitivity of the numerical results to algorithm configuration parameters, source condition, and propagation medium parameters has been studied in [13]. Further sensitivity study on seawater parameters, including consideration of their spatial and temporal variability, could be a useful contribution. Air bubbles in the upper ocean layer could also be investigated as a possible cause of variability in relation to finite-amplitude effects. Such studies would help to determine how accurately one can estimate the effects of finite-amplitude sound propagation on fishery acoustic survey data, and thus to what extent model-based compensation would reduce their measurement uncertainty with respect to backscattering cross sections and volume and area backscattering coefficients. Frequency-dependent absorption coefficients for seawater were inserted for each harmonic component in place of the frequency-squared power law inherent in the KZK equation being solved. Possible effects of this ad-hoc adaptation could be investigated, perhaps by comparing to a model derived with relaxation terms.

7 Conclusions

Oceanic survey measurements on Atlantic mackerel have been carried out to investigate possible influence of acoustic signal amplitude on parameters used in fish abundance estimation and species identification. For the relative frequency response, r_a, 14% and 47% deviations have been measured at 120 kHz and 200 kHz, respectively, by alternating between “low” and “high” electrical transmit power settings. The “low” power settings were 250 W and 120 W for the 120 kHz and 200 kHz measurement frequencies, respectively. The “high” power settings were 1000 W for both frequencies. Echosounder calibration was done with the “low” power setting for both frequencies.

The simulations strongly indicate that the observed systematic deviations can be ascribed to finite-amplitude sound propagation effects in seawater. Simulations and measurements were compared for the ratio between r_a(f) obtained using high and low power settings. The agreement was found to be mainly within the experimental uncertainty. The formal framework established through [18–21] is demonstrated to be functional for describing the measurement of volume backscattering under conditions of significant finite-amplitude sound propagation. A method to account for nonlinear loss in echosounder calibration is proposed within the same framework. The result is a complete method for quantifying and correcting for finite-amplitude sound propagation in fisheries research echosounders or establishing power limits with known and acceptable nonlinear loss for an application at hand. Potentials for further improved agreement are discussed.

If not accounted for, the demonstrated effects can lead to significant bias in biomass abundance estimation and acoustic species identification. Finite-amplitude effects are not accounted for in the conventional methods used today within fishery acoustics, but minimized approximately by recommended limits on maximum transmitted power [16, 17], partially based on overly simplified extrapolations [17]. The effects are however well understood and can, as shown here, be accounted for using available methods.

The results provide an indication of typical errors that can be encountered in fisheries acoustic measurements when high transmit powers are used and finite-amplitude sound propagation effects are not accounted for. The approach demonstrated here may serve as a basis to refine previous recommendations and add quantitative information on nonlinear loss. Nonlinear loss and generation of harmonics can be significant even at relatively low source levels and frequencies. A well-chosen transmitted power thus depends on transducer directivity, efficiency, frequency, and on medium parameters, and should be rooted in a quantitative requirement on measurement uncertainty. Further parametric study could provide a fuller map of measurement uncertainty contributions and potential corrections for finite-amplitude effects in both historical and current acoustic survey data. If finite-amplitude effects are accounted for in echosounder signal processing, higher powers may be used when needed to improve signal-to-noise ratio or extend measurement range in echosounder and sonar measurements. The acoustic quantities accounting for finite-amplitude sound propagation effects, ${\mathit{C}}_{\mathrm{i}}^{\mathrm{n}}\left(r\right)$, ${\mathit{B}}_{\mathrm{rel}}^{\mathrm{n}}\left(r,\theta ,\phi \right)$, and ${\psi }_{\mathrm{rel}}^{\mathrm{n}}\left(r\right)$, can be calculated for relevant field conditions, and measured for relevant echosounder transducers in a laboratory tank setup.

Conflict of interest

The authors declared no conflict of interest.

Acknowledgments

The work was supported by The Research Council of Norway (MARE 152790) and Norwegian Seafood Research Fund, Norway (FHF 222042). Discussions with M. Vestrheim, University of Bergen, and A.C. Baker, Christian Michelsen Research AS, are acknowledged. The associate editor and two reviewers are thanked for their valuable feedback.

Data availability statement

The research data associated with this article is available in Zenodo (https://doi.org/10.5281/zenodo.6059530 [38]). The data include measurement data exported to ASCII by LSSS and simulation results in a MATLAB .mat file. Accompanying MATLAB scripts produce Figures 2–5 and the results in Table 2 and include guidance to the data formats.

References

All Tables

All Figures

Figure 1 Cruise track (blue curves) for the acoustic survey and observed fish schools (purple blobs). The axes indicate longitude and latitude in decimal degrees. The values in the legend are the nautical area scattering coefficient sA [nmi−1] [33]. Note that the cruise tracks do not overlap completely. Each school could be crossed straight over in one leg and on a side in the next, even if the schools did not move.

In the text

Figure 2 Measured relative frequency response, ra(f), obtained from oceanic in-situ acoustical survey measurements of sv on Atlantic mackerel schools. Data from within an observation depth interval of 30 m to 110 m are included. Each point (x or o) is based on measurements from a whole school. (a) 120 kHz measurement frequency, using 250 W and 1000 W electrical power settings in survey operation (alternating), and 250 W during calibration. (b) 200 kHz measurement frequency, using 120 W and 1000 W electrical power settings in survey operation (alternating), and 120 W during calibration.

In the text

Figure 3 (a) Numerically calculated values for ${\left({\left|{\mathit{C}}_{\mathrm{i}}^{\mathrm{n}}\left(r\right)\right|}^2\cdot {\psi }_{\mathrm{rel}}^{\mathrm{n}}\left(r\right)\right)}^{-1}$, representing the integrated nonlinear loss over the sound beam and thus the nonlinear loss in the sv measurement. The parameters in Table 1 were used for the calculations. (b) Numerical and experimental values for the ratio $r_a^{\mathrm{meas},\mathrm{high}}/r_a^{\mathrm{meas},\mathrm{low}}$ between volume backscattering measurements using the “high” and “low” power settings. The experimental data are divided in depth bins with 20 m width, indicating mean values and standard deviations for the 11-hour measurement series.

In the text

Figure 4 Numerically calculated example values for the finite-amplitude correction factor $s_v^{\mathrm{correct}}∕s_v^{\mathrm{meas}}$, with sound propagation parameters from Table 1 and 1000 W power setting for both calibration at range rcal and sv measurement at range r. The calibration range is indicated with a mark on each curve. Logarithmic scale, with corresponding ratios shown on the right-hand ordinate axis. (a) 120 kHz, (b) 200 kHz.

In the text

Figure 5 Illustration of nonlinear loss in the sound beam, using the numerical results for 200 kHz, 1000 W. Left: Two-way intensity beam pattern at 14 m range, with and without correction for nonlinear loss in the forward propagating wave. Centre: Two-way intensity beam pattern at 14 m range, corrected for nonlinear loss upon both measurement and calibration at the same range. Right: Two-way intensity beam pattern at r = 300 m range, corrected for nonlinear loss upon measurement and upon calibration at rcal = 14 m range.

In the text