Computational Bayesian inference for active sonar localization under uncertainty in sound speed profile

Barros, Abner C.; Gendron, Paul J.

doi:10.1121/10.0035808

A computational Bayesian approach is presented to address active sonar localization under the challenges of small receive aperture, uncertain sound speed profile (SSP), and limited coherence time. The approach draws inference on wavevectors associated with the closely spaced angle/Doppler spread arrivals, characterizing the scattered acoustic field. The wavevector posterior density is mapped to the scattering body's location and speed under an uncertain SSP using eigenray interpolation and marginalization. SSP uncertainty is captured by a multivariate Gaussian and a low-dimensional subspace mode representation. A case study using SSPs from the Mediterranean Sea demonstrates the efficacy of this approach.

1. Introduction

Localization of underwater mobile objects using active sonar systems with small vertical aperture arrays presents significant challenges due to their limited angular resolution. Furthermore, the temporal integration duration is constrained by the limited coherence time associated with the motion of the object and the dynamics of the environment (Yang, 2006), resulting in reduced frequency resolution. Here, coherence time refers to the period over which a parametric model holds. These limitations hinder a sonar system's ability to resolve closely spaced wavevector components associated with the multipath arrivals of the measured pressure field. As a result, inferring the object location and state of motion from these arrivals becomes extremely challenging. This is further exacerbated in refractive undersea environments, as the solution must account for both the spatial dependence of angular refraction and the attendant path-dependent Doppler frequencies imparted at the object, which are also depth-dependent. Barros and Gendron (2019) developed a computational Bayesian approach that can resolve closely spaced multipath phase fronts in a horizontally stratified underwater environment to infer the location and speed of a submerged object in a bounded region of interest. The probability density function of the wavevector components of the arriving phase fronts is mapped to the range, depth, and speed of the scattering body using an efficient eigenray interpolation method. This approach was applied to scenarios in a channel characterized by a known sound speed profile (SSP). However, in many underwater acoustic applications, the SSP is not known precisely. Limited spatiotemporal measurements and incomplete knowledge of the driving forces of ocean circulation imply uncertain knowledge of the SSP. This uncertainty persists even when the data are taken with a conductivity–temperature–depth (CTD) profiler (Raiteri , 2018). For environments where measurements may not be available, prediction models such as the Hybrid Coordinate Ocean Model (HYCOM) can provide forecasts that can be used with empirical SSP models to generate profiles characterizing the channel. These profiles also carry a degree of uncertainty associated with the underlying prediction model.

This paper extends the Bayesian active localization approach to address uncertainty in the SSP by employing a computationally efficient solution based on the Laplace approximation and marginalization over uncertainty in sound speed. Uncertainty in the SSP is well captured by a multivariate Gaussian model with the mean SSP and a low-dimensional subspace mode representation. This model is constructed here using a HYCOM dataset. Similar expansions of the SSP using empirical orthogonal functions have been successfully applied in passive source localization (Collins and Kuperman, 1991; Huang , 2008; Lin and Michalopoulou, 2014). While computational Bayesian methods have been used to address uncertainties in environmental parameters, they have primarily been applied in the context of low-frequency passive localization (Dosso and Wilmut, 2008; Richardson and Nolte, 1991). This paper addresses this issue for object localization using high-frequency active sonar under the challenging constraints of a small-aperture array and a short processing window. The paper is organized as follows. Section 2 defines the active sonar signal model, which is used in Sec. 3 to develop the Bayesian approach that accounts for the uncertainty in the SSP. Section 4 then validates the method using a HYCOM dataset.

2. Active sonar signal modeling

The sonar system considered here consists of a small-aperture vertical linear array (VLA) in the water column. The VLA is positioned at such a distance from the region of interest that the far-field assumption of the scattered acoustic field holds between the scattering body and the receiver array. Furthermore, the directivity of the sonar projector and receiver enables the excitation and reception of the two most dominant acoustic paths, namely the direct and surface-reflected paths (See supplementary material). Therefore, the backscattered plane wave signal from the narrowband continuous active sonar (CAS) transmission, as observed by the mth element of the VLA receiver, is modeled as

r_{m} (t) = \sum_{p = 1}^{P = 2} a_{p} e^{j 2 π (f_{c} + Δ f_{p}) t} e^{- j k_{p} \cdot p_{m}} + n_{m} (t),

(1)

where P denotes the total number of acoustic paths, a_p is the attenuated complex amplitude, f_c is the known projector's center frequency,

Δ f_{p}

is the Doppler frequency,

p_{m}

is the array element position, and

n_{m} (t)

is temporally and spatially uncorrelated complex Gaussian white noise, i.e.,

n \sim C N (0, σ^{2})

⁠. The wavevector

k_{p}

and the Doppler frequency are given, respectively, by

k_{p} = \frac{2 π (f_{c} + Δ f_{p})}{c (z_{R})} {[cos (θ_{p}), sin (θ_{p})]}^{T}, Δ f_{p} = \frac{2 v_{T} f_{c}}{c (z_{T})} cos (ϕ_{p}),

(2)

where

c (z_{R})

is the sound speed at the receiver's depth, θ_p is the vertical angle of arrival, v_T is the object's speed,

c (z_{T})

is the sound speed at the object's depth, and

ϕ_{p}

is the angle of incidence at the object.

Using vector notation, the baseband array observation is expressed compactly as

r = a_{1} v_{1} + a_{2} v_{2} + n

⁠. Here, the vector

v_{p} (θ_{p}, Δ f_{p}) = a (θ_{p}) \otimes e (Δ f_{p}),

for p = 1, 2, is the Kronecker product of the plane wave array manifold

a (θ_{p}) = exp (- j k_{p} \cdot p_{m}),

where

m = - (M - 1) / 2, \dots, (M - 1) / 2,

and the temporal phase vector

e (Δ f_{p}) = exp (j 2 π Δ f_{p} n T_{s}),

where T_s is the sampling period and

n = - (N - 1) / 2, \dots, (N - 1) / 2.

The noise process n models ambient acoustic noise and can include modeling errors, reverberation, and unmodeled arrival paths. For simplicity, here n models only surface-generated ambient noise. At high frequencies (15–30 kHz), underwater ambient acoustic passband noise at low vertical grazing angles

(θ < - 50 °)

and relatively significant depths

(z_{R} > 50 m)

is well approximated as a Gaussian process (Heitmeyer, 2006). Thus, the baseband noise process is modeled as a symmetric complex, additive, zero-mean Gaussian with covariance matrix

Σ = σ^{2} I

⁠. Consequently, the likelihood for the baseband array observation is a multivariate complex Gaussian,

f (r | μ, Σ) = \frac{1}{π^{L} | Σ |} exp {- {(r - μ)}^{H} Σ^{- 1} (r - μ)},

(3)

where

L = M \times N

is the number of observations and

μ (k_{1}, a_{1}, k_{2}, a_{2}) = a_{1} v_{1} + a_{2} v_{2} .

3. Computational Bayesian approach for object localization

This section describes the computational Bayesian approach for wavevector inference, the transformation of the wavevector posterior to the object state posterior, and the extension of this approach to efficiently account for uncertainty in the SSP.

3.1 Wavevector posterior factorization

For the present estimation problem, the parameters of interest are

{k_{1}, a_{1}, k_{2}, a_{2}, σ^{2}},

with

k_{p} = {θ_{p}, Δ f_{p}},

forming a nine-dimensional parameter space, where

{a_{1}, a_{2}, σ^{2}}

are nuisance parameters. The objective is to draw inferences regarding these parameters using a Bayesian approach that constructs the following joint posterior probability density (PPD),

p (k_{1}, a_{1}, k_{2}, a_{2}, σ^{2} | r) \propto f (r | k_{1}, a_{1}, k_{2}, a_{2}, σ^{2}) π (k_{1}, a_{1}, k_{2}, a_{2}, σ^{2}),

(4)

where

f (\cdot)

is the likelihood function given by Eq. and

π (\cdot)

is the prior on the parameters. Since this PPD is over a high-dimensional parameter space, it does not admit a closed-form solution and must be solved numerically. Barros and Gendron (2019) developed a Markov chain Monte Carlo (MCMC) Gibbs sampling scheme that iteratively samples from each full conditional density of the parameters to construct the joint PPD in Eq. . Figure 1(a) illustrates the MCMC iteration and shows an example of the wavevector PPD for a scenario with the Munk SSP. This approach reduces the high-dimensional problem of constructing the PPD in Eq. into more manageable two-dimensional (2-D) and one-dimensional (1-D) problems. The selection of conjugate priors—based on a priori simulations and historical data—for the complex amplitudes and ambient acoustic noise power results in analytic conditional posteriors that arise from the normal-inverse-gamma family (Barros, 2020). However, the conditional densities of the wavevectors do not admit an analytic solution and must be sampled numerically. This is facilitated by recognizing and exploiting ordering features in their a priori domain (constructed from ray tracing) and using a 2-D quantile sampling approach.

Fig. 1.

View large Download slide

Illustration of a computational Bayesian approach for wavevector inference and object localization. (a) MCMC Gibbs sampling iteration. Here, sound speed uncertainty (denoted by λ_c) is illustrated with the deviation of wavevectors $k_{1}$ and $k_{2}$ ⁠. (b) Depiction of the transformation of the wavevector posterior (constructed under the mean SSP) to object state under an uncertain SSP.

3.2 Transformation of wavevector posterior to object state posterior

The wavevector PPD $p (k_{1}, k_{2} | r)$ is mapped to the joint PPD over range, depth, and speed using a numerical acoustic ray interpolation approach [see Fig. 1(b)]. This is a finite difference method that involves interpolation of the scatterer's range and depth on a coarsely sampled grid constructed a priori using the Gaussian beam ray tracing implemented by the bellhop program (Porter and Bucker, 1987). This approach is computationally efficient, allowing the entire joint density to be specified with a single ray tracing.

Let the transformation be denoted as

g : y \to x,

where

y = {[θ_{1}, θ_{2}, Δ f_{1}, Δ f_{2}]}^{T}

represents the wavevector parameters and

x = {[r_{T}, z_{T}, v_{T}]}^{T}

denotes the object state, with range r_T, depth z_T, and speed v_T. This formulation assumes that all samples of the wavevector PPD map to a bounded region of interest in the state space. The inversion

x^{(i)} = g (y^{(i)})

for the ith MCMC sample is formulated from a first-order approximation,

g (y^{(i)}) = g (y_{0}^{(i)} + Δ y^{(i)}) \approx g (y_{0}^{(i)}) + J^{+} Δ y^{(i)},

(5)

where

Δ y^{(i)} = y^{(i)} - y_{0}^{(i)}

and

y_{0}^{(i)} = argmin | | y^{(i)} - y_{m, n, k} | |^{2}

for all m, n, k corresponding to the range, depth, and speed indices of the pre-computed grid. Hence,

g : y_{0}^{(i)} \to x_{0}^{(i)}

is a point in the pre-computed grid that is closest to the sample

y^{(i)}

⁠. In Eq. ,

J^{+} = A {(A J^{T} J A)}^{- 1} A J^{T}

is the pseudoinverse of the Jacobian matrix

J = \partial y / \partial x |_{y_{0}},

which is computed numerically. The matrix

A = diag (Δ r, Δ z, Δ v)

is used for numerical stability, with diagonal elements corresponding to the grid spacing.

3.3 Posterior marginalization over sound speed uncertainty

Figure 1(a) provides a conceptual illustration of how uncertainty in sound speed manifests in the angle–Doppler space, represented by additional uncertainty λ_c in the wavevectors

k_{1}

and

k_{2}

⁠. Taking into account for this uncertainty in the wavevector sampling scheme is computationally expensive. Similarly, computing the marginal PPD

p (x | r) = \int_{C} p (x, c | r) d c

in the state space, which properly accounts for the uncertainty in the SSP c, is also demanding. The approach considered here provides a computationally fast solution using posterior factorization, Laplace approximation, and a low-dimensional subspace representation of SSP uncertainty. This approach is illustrated in Fig. 1(b). In this method, the wavevector joint PPD is constructed by running the MCMC scheme once with the mean SSP, denoted as

μ_{c}

⁠. We consider the mean-centered profile

δ_{c} = c - μ_{c}

and represent the lower-dimensional subspace uncertainty through its modes as

δ_{c} = Ψ ζ,

corresponding to the Karhunen–Loéve decomposition. Here, matrix

Ψ = [ψ_{1}, \dots, ψ_{J}]

consists of orthogonal mode functions, with J denoting the number of modes and

ζ = {[ζ_{1}, \dots, ζ_{J}]}^{T}

representing the mode amplitudes that scale each mode function accordingly. The prior SSP uncertainty can be accurately characterized by the mode set

{ψ_{j}}_{j = 1}^{J}

and the corresponding amplitudes

{ζ_{j}}_{j = 1}^{J},

provided the selected modes capture most of the energy or uncertainty in the profile. Sound speed uncertainty is modeled here as a multivariate Gaussian (MVG) distribution,

ζ \sim N (0, Λ),

where

Λ = diag (λ_{c, 1}^{2}, \dots, λ_{c, J}^{2})

and

λ_{c, j}^{2}

denotes the variance or eigenvalue associated with the jth mode. To approximate the solution for

p (x | r),

we consider the following posterior marginalization:

p (x | r) = \int_{ζ} q (x | r, c = μ_{c} + Ψ ζ) N_{ζ} (0, Λ) d ζ .

(6)

The mode set characterizing the SSP uncertainty can be determined a priori from available measured or predicted SSPs. Here, the mode set is computed from a dataset of profiles compiled using salinity, temperature, and depth data from the HYCOM database and Leroy's SSP model (Leroy , 2008). The curated dataset consists of 420 profiles [solid gray lines in Fig. 2(a)] for a region in the Mediterranean Sea near the Strait of Sicily (See supplementary material), covering an area of

20 \times 20 {km}^{2}

⁠. These profiles were generated using HYCOM data spanning five years (2019–2023), extracted for the first week of January each year and sampled daily at 0900 UTC at 12 evenly spaced geographical locations within the region. The sample covariance matrix (SCM) for the data is computed as

\hat{R} = \frac{1}{K} \sum_{k = 1}^{K} δ_{c, k} δ_{c, k}^{T},

where K is the number of profiles in the dataset. The eigendecomposition of the SCM is given by

\hat{R} = Ψ Λ Ψ^{- 1},

yielding the set of modes. Figure 2(b) presents the three most dominant eigenvectors (or mode functions) from this decomposition, while Fig. 2(c) shows the first 10 ordered eigenvalues out of 32, highlighting the prominence of the first three modes. Figure 2(a) displays 15 profiles (solid blue lines) sampled using the first five modes, J = 5, with the MVG model described above. These samples show good agreement with the HYCOM SSPs, underscoring the effectiveness of the MVG model and the low-dimensional mode set representation of the SSP uncertainty. It is important to note that the conditional variance in sound speed

λ_{c | HYCOM}^{2},

derived from the HYCOM predictor, is not readily available. Therefore, the variability observed in the SSP dataset shown in Fig. 2(a), representing the unconditional variance

λ_{c}^{2}

in the sound speed of a sample from the HYCOM dataset, serves as a surrogate for

λ_{c | HYCOM}^{2}

⁠. Specifically, the conditional variance

λ_{c | HYCOM}^{2} (t_{2}; x)

at time t₂ and geographical location x is replaced by

λ_{c}^{2} (t_{1}; x) [1 - ρ^{2} (t_{1}, t_{2}; x)],

where

ρ (t_{1}, t_{2}; x) = Cov [δ_{c} (t_{1}; x), δ_{c} (t_{2}; x)] / \sqrt{λ_{c}^{2} (t_{1}; x) λ_{c}^{2} (t_{2}; x)}

captures the spatiotemporal correlation between two profiles predicted for location x at times t₁ and t₂. This surrogate represents a worst-case scenario, as the unknown variances for salinity, temperature, and depth associated with the HYCOM predictor are expected to yield smaller conditional variance for the sound speed.

Fig. 2.

View large Download slide

Sound speed profiles (SSPs) from the Mediterranean Sea and their decomposition. (a) Gray solid lines denote SSPs computed using HYCOM data. The yellow dashed line indicates the mean SSP, and the blue solid lines are profiles sampled using the multivariate Gaussian model capturing the SSP uncertainty with the first five modes, J = 5. (b) The three most dominant mode functions from the eigendecomposition of the profiles. (c) The first 10 ordered eigenvalues (out of 32).

The posterior

q (x | r, c = μ_{c} + Ψ ζ)

in Eq. is assumed to belong to the exponential family of probability density functions and is approximated using only the quadratic terms, following the Laplace approximation:

N_{x | r, c} (m_{x} (c), Γ (c))

⁠. For simplicity, we assume

Γ (c) \approx Γ (μ_{c})

and approximate the posterior mean using the linear terms of a first-order Taylor expansion,

m_{x} (c) \approx m_{x} (μ_{c}) + S Ψ ζ,

where

S = \partial m_{x} (c) / \partial c |_{c = μ_{c}}

is the Jacobian matrix. Marginalization over the uncertainty in sound speed is achieved using the following posterior approximation,

p (x | r) \approx \int_{ζ} N_{\tilde{x} | r, μ_{c}} (m_{x} (μ_{c}), Γ (μ_{c})) N_{ζ} (0, Λ) d ζ = \int_{ξ_{x}} N_{\tilde{x} | r, μ_{c}} (m_{x} (μ_{c}), Γ (μ_{c})) N_{ξ_{x}} (0, P) | G | d ξ_{x},

(7)

where

\tilde{x} = x - ξ_{x},

with

ξ_{x} = S Ψ ζ .

Here, the transformation

ϕ : ζ \to ξ_{x}

maps the uncertainty in sound speed to the object state space, represented by the distribution

N_{ξ_{x}} (0, P),

where

P = S Ψ Λ Ψ^{T} S^{T}

is constructed using the low-dimensional mode set (i.e., with J modes). The Jacobian determinant is given by

| G | = | S Ψ Ψ^{T} S^{T} |^{- 1 / 2}

⁠. Equation represents a convolution of the two densities, i.e.,

p (x | r) \approx \int q (x - ξ_{x} | r, μ_{c}) h (ξ_{x}) d ξ_{x},

which can be efficiently computed. This approach replaces the posterior

q (x | r, c = μ_{c} + Ψ ζ)

in Eq. , which would otherwise require construction over many profiles, with the Laplace approximation

q (x | r, μ_{c}) \approx N_{\tilde{x} | r, c} (m_{x} (μ_{c}), Γ (μ_{c})),

which is constructed by running the MCMC scheme once with the mean SSP and inverting the resulting PPD to the object state space. A similar analysis using the law of total covariance provides a second-order lower-bound approximation, decomposing the posterior covariance as follows:

Cov [x | r] = E_{C} {Cov [x | r, c = μ_{c} + Ψ ζ]} + {Cov}_{C} [E {x | r, c = μ_{c} + Ψ ζ}] .

(8)

Since the posterior covariance

Cov [x | r, c = μ_{c} + Ψ ζ]

remains relatively constant across different SSPs and varies more significantly with respect to the signal-to-noise ratio (SNR) and array aperture (Barros, 2020), the first term on the right-hand side (RHS) of Eq. is well approximated as

E_{C} {Cov [x | r, c = μ_{c} + Ψ ζ]} \approx Cov [x | r, μ_{c}]

⁠. This is equivalent to the covariance

Γ (μ_{c})

of the Laplace approximation in Eq. . Similarly, the second term on the RHS approximates the covariance of the density of

ξ_{x}

in Eq. , i.e.,

{Cov}_{C} (E {x | r, c = μ_{c} + Ψ ζ}) \approx P .

One approach to approximate the density of

ξ_{x}

is by mapping the wavevector posterior mean,

E {k_{1}, k_{2} | r, μ_{c}},

to the object state space mean under a sample SSP,

{\hat{m}}_{x, k} = E {x | r, c_{k} = μ_{c} + Ψ ζ_{k}},

for

k = 1, \dots, K,

and then performing kernel density estimation using these samples.

4. Case study and discussion

To lend credence to the computational Bayesian approach presented in Sec. 3, this section presents a case study focusing on a short-range, low-SNR scenario in a channel characterized by the mean SSP shown in Fig. 2(a). Figure 3 depicts the scenario under consideration, showing ray tracing and the eigenrays for the ground truth location, indicated by the red marker at a range of $r_{T} = 4.17 km$ and a depth of $z_{T} = 167.3 m,$ with a velocity of $v_{T} = 5.3 {ms}^{- 1}$ ⁠. The angle and Doppler differences between the two paths are $Δ θ = 5 °$ and $Δ η = 0.74 Hz,$ respectively. Priors on the object state are uniformly distributed within the shaded surveillance region depicted in Fig. 3, with $r_{T} \sim U (2, 6) km, z_{T} \sim U (80, 250) m,$ and $v_{T} \sim U (3, 8) {ms}^{- 1}$ ⁠. The grid spacings used in the eigenray interpolation approach are $Δ r = 40 m, Δ z = 10 m,$ and $Δ v = 0.2 {ms}^{- 1}$ ⁠. The sonar projector's center frequency is $f_{c} = 15 kHz$ ⁠. At the receiver, the processor employs a window duration of $T_{d} = 81 ms$ (yielding a Rayleigh resolution of $Δ η_{R} = 12.3 Hz),$ and the ${SNR}_{i n}$ is $0 dB,$ which is the SNR measured at the output of a single hydrophone element. The VLA receiver consists of 19 elements centered at $z_{R} = 350 m$ (shown by the green marker in Fig. 3) with a standard element spacing of $λ_{f_{c}} / 2 = 0.05 m,$ calculated assuming a nominal sound speed of $c = 1500 {ms}^{- 1},$ resulting in an aperture of $L_{a} = 0.9 m$ (with a Rayleigh resolution of $Δ θ_{R} = 6 °) .$

Fig. 3.

View large Download slide

Scenario for a case study with the mean SSP from a region in the Mediterranean Sea (displayed for depths up to $450 m)$ ⁠. A ray tracing illustrates the object–receiver geometry, with eigenrays indicating the ground truth object location. The cyan shaded box outlines the prior surveillance region of interest.

Figure 4(a) shows the joint PPD, $q (x | r, μ_{c}),$ over the range and depth of the object (i.e., marginalized over speed) obtained from the wavevector inversion under the mean SSP. This result demonstrates the computational Bayesian method's ability to resolve closely spaced multipath arrivals even with a limited receive aperture and accurately localize the object in both range and depth, as indicated by the ground truth dashed lines. Results with speed inference are detailed in Barros (2020) and Barros and Gendron (2019) and are omitted here for brevity. The posterior over range, depth, and speed provides credible intervals that are essential for submersible localization and tracking applications. For instance, this PPD or its moments can serve as input to Kalman filters and particle filters for higher-level tracking.

Fig. 4.

View large Download slide

Probability densities over range and depth. (a) Joint PPD constructed under the mean SSP from the Mediterranean Sea dataset. The dashed lines denote ground truth. (b) Posterior obtained by convolving the densities in panels (a) and (c), i.e., approximated by Eq. (7). (c) Probability density $h (ξ_{x})$ depicting uncertainty in range $ξ (1)$ and depth $ξ (2)$ due to the SSP uncertainty.

To enable the computationally efficient solution given by Eq. (7), note that the PPD shown in Fig. 4(a) is now approximated via Laplace approximation as follows: $q (x | r, μ_{c}) \approx N_{\tilde{x} | r, c} (m_{x} (μ_{c}), Γ (μ_{c}))$ ⁠. Figure 4(c) shows the density of $ξ_{x},$ capturing uncertainty in range and depth due to the SSP uncertainty. As detailed in Sec. 3.3, this density is constructed by mapping the mean of the wavevector PPD, $E {k_{1}, k_{2} | r, C = μ_{c}},$ to the state space under K different profiles using the efficient eigenray interpolation approach from Sec. 3.2. The resulting density is then mean-centered. Furthermore, we consider a dataset comprising of 2000 profiles, of which 420 are the HYCOM profiles and 1580 were sampled using the MVG model capturing the SSP uncertainty, $ζ \sim N (0, Λ),$ with J = 5 modes. Figure 4(b) shows the density resulting from the convolution of the two densities depicted in Figs. 4(a) and 4(c). This density corresponds to the posterior approximation given in Eq. (7), incorporating the uncertainty in range and depth due to uncertainty in sound speed. Similarly, the covariance of the density shown in Fig. 4(b) approximates the covariance in Eq. (8). A detailed comparison of these densities reveals that the covariance in Fig. 4(b) increases by approximately a factor of 3 compared to the covariance in Fig. 4(a). The posterior shown in Fig. 4(b) is promising, not only because it demonstrates an accurate localization solution but also because the conditional variance in sound speed used here reflects a worst-case scenario.

Now consider the $100 (1 - α) %$ highest posterior density (HPD) interval, denoted by the lower and upper bounds $(θ_{L}, θ_{U})$ on the parameter θ (Chen and Shao, 1999). These bounds are computed by first marginalizing over the range–depth joint PPD. For the PPD shown in Fig. 4(a), the 90% HPD intervals computed for the marginalized range and depth densities are $(4102.5, 4394.3) m$ and $(157.95, 177.22) m,$ respectively. The same HPD intervals computed for the marginals of the density shown in Fig. 4(b) are $(4080, 4458.7) m$ for range and $(159.71, 175.29) m$ for depth. Although this case study focuses on a single scenario, it is noteworthy that SSP uncertainty has a greater influence on range than on depth, as observed in Figs. 4(b) and 4(c) and reflected by the HPD intervals. One possible explanation is that the object's depth lies below the mixed layer, which typically exhibits less variability. Overall, this case study demonstrates the effectiveness of the computational Bayesian approach in addressing uncertainty in environmental parameters such as sound speed. The approach produces an interpretable PPD over the object location and provides informative credible intervals, which are useful for other processing tasks such as tracking. To summarize, we have extended the posterior density analysis to provide a more accurate measure of uncertainty in localization associated with SSP uncertainty.

5. Conclusions

A computational Bayesian approach for active sonar localization is extended here to incorporate uncertainty in the SSP. Posterior inference under an uncertain SSP is demonstrated using a second-order variational Bayesian method that employs a multivariate Gaussian model of the SSP. A case study using simulated acoustic fields with profiles from the HYCOM database demonstrates the efficacy of the approach and lends further credence to the computational Bayesian approach.

Supplementary Material

See the supplementary material for figures illustrating the acoustic field and a map of the location where the HYCOM data were sampled in the Mediterranean Sea.

Acknowledgments

This research was supported by the Naval Undersea Warfare Center (NUWC) Division Newport's In-house Laboratory Independent Research (ILIR) program. P.J.G. is supported by the Office of Naval Research under Contract Nos. N000142412763, N000142212012, and N000142412231.

Author Declarations

Conflict of Interest

The authors have no conflicts to disclose.

Data availability

The data that support the findings of this study are openly available in the HYCOM database at https://hycom.org.

References

1.

Barros

,

A. C.

(

2020

). “

A computational Bayesian approach to localize a mobile scatterer from sub-Rayleigh resolvable acoustic wavefronts in refractive media

,” M.S. thesis,

University of Massachusetts Dartmouth

,

Dartmouth, MA

.

Google Scholar

2.

Barros

,

A. C.

, and

Gendron

,

P. J.

(

2019

). “

A computational Bayesian approach for localizing an acoustic scatterer in a stratified ocean environment

,”

J. Acoust. Soc. Am.

146

(

3

),

EL245

–

EL250

.

https://doi.org/10.1121/1.5126351

Google Scholar

Crossref

PubMed

3.

Chen

,

M.-H.

, and

Shao

,

Q.-M.

(

1999

). “

Monte Carlo estimation of Bayesian credible and HPD intervals

,”

J. Comput. Graph. Stat.

8

(

1

),

69

–

92

.

https://doi.org/10.1080/10618600.1999.10474802

Google Scholar

Crossref

4.

Collins

,

M. D.

, and

Kuperman

,

W. A.

(

1991

). “

Focalization: Environmental focusing and source localization

,”

J. Acoust. Soc. Am.

90

(

3

),

1410

–

1422

.

https://doi.org/10.1121/1.401933

Google Scholar

Crossref

PubMed

5.

Dosso

,

S. E.

, and

Wilmut

,

M. J.

(

2008

). “

Uncertainty estimation in simultaneous Bayesian tracking and environmental inversion

,”

J. Acoust. Soc. Am.

124

(

1

),

82

–

97

.

https://doi.org/10.1121/1.2918244

Google Scholar

Crossref

PubMed

6.

Heitmeyer

,

R. M.

(

2006

). “

A probabilistic model for the noise generated by breaking waves

,”

J. Acoust. Soc. Am.

119

(

6

),

3676

–

3693

.

https://doi.org/10.1121/1.2188916

Google Scholar

Crossref

7.

Huang

,

C.-F.

,

Gerstoft

,

P.

, and

Hodgkiss

,

W. S.

(

2008

). “

Effect of ocean sound speed uncertainty on matched-field geoacoustic inversion

,”

J. Acoust. Soc. Am.

123

(

6

),

EL162

–

EL168

.

https://doi.org/10.1121/1.2908406

Google Scholar

Crossref

PubMed

8.

Leroy

,

C. C.

,

Robinson

,

S. P.

, and

Goldsmith

,

M. J.

(

2008

). “

A new equation for the accurate calculation of sound speed in all oceans

,”

J. Acoust. Soc. Am.

124

(

5

),

2774

–

2782

.

https://doi.org/10.1121/1.2988296

Google Scholar

Crossref

PubMed

9.

Lin

,

T.

, and

Michalopoulou

,

Z.-H.

(

2014

). “

Sound speed estimation and source localization with linearization and particle filtering

,”

J. Acoust. Soc. Am.

135

(

3

),

1115

–

1126

.

https://doi.org/10.1121/1.4864787

Google Scholar

Crossref

PubMed

10.

Porter

,

M. B.

, and

Bucker

,

H. P.

(

1987

). “

Gaussian beam tracing for computing ocean acoustic fields

,”

J. Acoust. Soc. Am.

82

(

4

),

1349

–

1359

.

https://doi.org/10.1121/1.395269

Google Scholar

Crossref

11.

Raiteri

,

G.

,

Bordone

,

A.

,

Ciuffardi

,

T.

, and

Pennecchi

,

F.

(

2018

). “

Uncertainty evaluation of CTD measurements: A metrological approach to water-column coastal parameters in the Gulf of La Spezia area

,”

Measurement

126

,

156

–

163

.

https://doi.org/10.1016/j.measurement.2018.05.058

Google Scholar

Crossref

12.

Richardson

,

A. M.

, and

Nolte

,

L. W.

(

1991

). “

A posteriori probability source localization in an uncertain sound speed, deep ocean environment

,”

J. Acoust. Soc. Am.

89

(

5

),

2280

–

2284

.

https://doi.org/10.1121/1.400918

Google Scholar

Crossref

13.

Yang

,

T. C.

(

2006

). “

Measurements of temporal coherence of sound transmissions through shallow water

,”

J. Acoust. Soc. Am.

120

(

5

),

2595

–

2614

.

https://doi.org/10.1121/1.2345910

Google Scholar

Crossref

2021

U.S. Government

Computational Bayesian inference for active sonar localization under uncertainty in sound speed profile

1. Introduction

2. Active sonar signal modeling

3. Computational Bayesian approach for object localization

3.1 Wavevector posterior factorization

3.2 Transformation of wavevector posterior to object state posterior

3.3 Posterior marginalization over sound speed uncertainty

4. Case study and discussion

5. Conclusions

Supplementary Material

Acknowledgments

Author Declarations

Conflict of Interest

Data availability

References

Supplementary Material

Citing articles via

Submit your article

Sign up for alerts

Resources

Explore

pubs.aip.org

Connect with AIP Publishing

Computational Bayesian inference for active sonar localization under uncertainty in sound speed profile Free

1. Introduction

2. Active sonar signal modeling

3. Computational Bayesian approach for object localization

3.1 Wavevector posterior factorization

3.2 Transformation of wavevector posterior to object state posterior

3.3 Posterior marginalization over sound speed uncertainty

4. Case study and discussion

5. Conclusions

Supplementary Material

Acknowledgments

Author Declarations

Conflict of Interest

Data availability

References

Supplementary Material

Citing articles via

Submit your article

Sign up for alerts

Related Content

Resources

Explore

pubs.aip.org

Connect with AIP Publishing

This Feature Is Available To Subscribers Only

Computational Bayesian inference for active sonar localization under uncertainty in sound speed profile