1. Introduction
In underwater environments, electromagnetic signal is easily dissipated; thus, sound is mainly used for underwater target localization. As a fundamental function for sonar systems, underwater target localization has gained significant attention [
1,
2,
3,
4,
5,
6]. Target localization by measuring the bearing of target’s sound is an important issue in passive target tracking.
This paper handles localizing a non-cooperative target, such that many sonobuoys are located to measure the bearing of the target’s sound. In other words, we consider the case where multiple sonobuoys are located to measure the bearing of an underwater target’s sound. We consider a Directional Frequency Analysis and Recording (DIFAR) sonobuoy, which consists of an omni-directional hydrophone and two dipole sensors oriented orthogonally [
7]. The directional information provided by the dipoles makes it possible to estimate the bearing to an underwater acoustic source [
7]. A sonobuoy has very low bearing accuracy, such as 10 degrees [
8,
9]. The details of a sonobuoy system is addressed in [
8,
10].
The detection of the target’s sound is associated to the estimation of the target’s bearing using multiple hydrophones. How to estimate the bearing of the target sound using multiple hydrophones is discussed in [
7,
10,
11,
12,
13]. However, this bearing estimation (target detection phase) is not within the scope of our paper.
Each sonobuoy can measure the bearing (azimuth) of the target sound. However, a sonobuoy cannot measure the elevation angle of the target sound. Thus, the depth of a target cannot be estimated using sonobuoy measurements. This implies that 3D position of the target cannot be estimated using sonobuoy measurements only.
We assume that by processing the target sound, each sonobuoy can measure the bearing of the target sound in real time. The bearing measurements of each sonobuoy are transmitted to a base station, where the target’s 2D location, except for the target’s depth, is calculated in real time. By processing the bearing measurements of every sonobuoy, the base station estimates the underwater target’s 2D position, except for the target’s depth, in real time. This paper handles how to estimate the target’s 2D position using the bearing measurements of multiple sonobuoys.
In practice, we can use multiple heterogeneous bearing sensors. Thus, the variance of a sensor noise may be different from that of another sensor. In addition, the maximum sensing range of a sensor may be different from that of another sensor. For instance, the performance of a cylindrical array sonar may be different from that of another sonar, depending on the geometry of array segments [
14].
In the considered bearings-only tracking problem, a sonar sensor measures the bearing of the target, but not the elevation of the target. Bearings-only measurements have been widely used for target localization [
15,
16,
17,
18]. Bearings-only measurements can be further utilized for tracking a moving target [
19,
20,
21,
22,
23]. Multiple turn models, such as the Interacting Multiple Model (IMM), were utilized to track a maneuvering target based on bearings-only measurements [
19,
24].
Bearing measurements can be utilized for tracking a moving target [
17,
25], such as a target with a constant velocity [
26]. The Range-Parameterized Extended Kalman Filter (RPEKF) [
16,
27] is a Gaussian-sum filter with multiple EKFs, each initialized at an estimated target range. In this way, we can reduce the initial range estimation error. The target position estimate is then calculated by merging all EKF outputs. The RPEKF assumes that the maximum sensing range of the observer is known in advance. Inspired by the RPEKF, we also assume that the maximum sensing range of a sonobuoy is known in advance. The RPEKF considers tracking a target using measurements of a single bearing sensor. In our paper, we handle the case where multiple bearing sensors are used for locating a target.
Considering bearings-only tracking using a single bearing sensor, Ref. [
28] addressed the observability analysis. The authors of [
28] showed that sensor maneuver is required to get the observability on the target state. In our paper, we utilize multiple bearing sensors to locate a target. Thus, sensor maneuver is not required to observe the target state.
The authors of [
29] presented a formation control of multiple bearing sensors for tracking a target based on bearings-only measurements. The authors of [
29] assumed that each bearing sensor can measure the noisy target bearing in real time. Ref. [
29] showed that uniform angular distribution centered at the target is the optimal formation configuration. However, in practice, the target position is not known a priori. Thus, the bearing sensors cannot form uniform angular distribution centered at the target. Optimal configuration of bearing sensors is not within the scope of our paper.
Consider the case where multiple bearing sensors are used to localize a target. Least Square (LS) estimation with a closed-form solution, was designed in [
30,
31,
32] for target localization using bearings-only measurements. The bearings-only LS suffers from severe bias problems. In [
33], an improved LS with bias compensation strategy was developed for bearings-only passive target tracking. Moreover, Ref. [
33] assumed that the bearing noise is sufficiently small, such that the first order Taylor series on the bearing noise can be applied. Our paper considers the case where bearing accuracy of sonobuoy sensor is very low; thus, the first order Taylor series on the bearing noise cannot be applied.
The LS was extended to the Weighted Least Square (WLS) estimation in [
34]. In the WLS, measurement noise variance was used to enhance the estimation result, assuming that the bearing measurement noises are sufficiently small. The authors of [
31,
35] extended the LS to the Total Least Square (TLS) estimation algorithm. The TLS was developed based on orthogonal vectors with the advantage of simplicity and reduced bias in the presence of bearing noise and observer position errors. Refs. [
31,
35] argued that the TLS has better performance than the LS. However, LS [
30,
31,
32], WLS [
34], and TLS [
31,
35] did not consider the maximum sensing range of a sonar sensor. Using MATLAB simulations, we verify the outperformance (considering both location accuracy and computational load) of the proposed localization method by comparing it with LS, WLS, and TLS localization methods.
Several papers [
36,
37] addressed bearings-only localization using multiple passive sensors, which requires an initial guess for target position. Reference [
36] introduced the Gauss Newton (GN) algorithm applied to bearings-only target estimation. The GN method is an iterative method to derive the target estimation. The convergence of the GN estimator is sensitive to the initial guess and the step size [
38]. Assuming that bearing measurements and inter sensor bearing measurements are available, Ref. [
37] addressed a subspace algorithm which requires an initial guess for target position. In our paper, we do not compare the proposed localization approach with [
36,
37], since [
36,
37] require an initial guess for target location. Note that the localization approach proposed in our paper can be used as an initial guess required in [
36,
37].
We consider a scenario where only a few sonobuoy sensors exist and each sensor has very low bearing accuracy, such as 10 degrees [
8]. In order to cope with poor accuracy of sonobuoy sensors, this paper introduces a novel target localization approach based on multiple Virtual Measurement Sets (VMS). Here, each VMS is derived considering the bearing measurement noise of each sonar sensor. Among all feasible target estimates derived from each VMS, we select an estimate which has the minimum residual (the difference between the actual and predicted measurements), such that the estimate satisfies the sonar sensing range constraints. The selected estimation is then utilized as the final target estimation.
Inspired by the RPEKF, we assume that the maximum sensing range (e.g., 15 km which appeared in [
9]) of a sonobuoy is known in advance. As far as we know, our paper is novel in localization of a target using heterogeneous sonobuoy sensors with low accuracy, considering the sensing range constraints of a sonar sensor. Our paper is distinct from LS [
30,
31,
32], WLS [
34], and TLS [
31,
35], since our location approach considers the maximum sensing range of a sonar sensor. Using MATLAB simulations, we demonstrate the superiority (considering both time efficiency and location accuracy) of the proposed localization, by comparing it with LS, WLS, and TLS localization methods.
This paper is organized as follows.
Section 2 presents the bearings-only target localization. Using MATLAB simulations in
Section 3, we verify the performance of the proposed estimation method based on VMS.
Section 4 provides the conclusions.
2. Bearings-Only Target Localization
Suppose that there are
K heterogeneous sonobuoy sensors which measure the bearing of the target. We assume that by processing the target’s sound, each sonobuoy can measure the bearing of the target’s sound in real time. How to calculate the target’s bearing using a sonobuoy is discussed in [
7,
10,
11,
12,
13].
The bearing measurements of each sonobuoy are transmitted to a base station, where the target’s 2D location is calculated in real time. By processing the bearing measurements of every sonobuoy, the base station estimates the underwater target’s 2D position in real time. This paper handles how to estimate the target’s 2D position using the bearing measurements of multiple sonobuoys.
Let
denote the true target position. Let
denote the estimated target position. Let
denote the
k-th sonobuoy sensor position (
). We assume that each sonobuoy sensor is localized using Global Positioning Systems (GPS). Let
present the bearing angle of the target with respect to the
k-th sonar sensor.
Here,
represents the true bearing associated to true relative position. In (
1),
is used to represent the phase (angle) of a complex number
.
in (
1) has Gaussian distribution with zero mean and variance
. Here,
(
) is feasible, since we handle heterogeneous sensors. We assume that
is known a-priori using experiments with the
k-th sensor.
The empirical rule states that 99.7 percents of data observed following a normal distribution exists within 3 standard deviations of the mean. Thus,
in (
1) is within the following interval:
with 99.7 percents probability. In (
1),
can be considered as Unknown-But-Upperbounded (UBU) by
.
Many papers on bearing-only localization [
33,
39,
40] assumed that the bearing noise is sufficiently small, such that the first order Taylor series on the bearing noise can be applied. This implies that [
33,
39,
40] used
where
in (
1). Our article considers the case where bearing accuracy is very low; thus, the first order Taylor series on the bearing noise cannot be applied. Therefore,
cannot be ignored, as in (
3).
Let
denote the maximum sensing range of the
k-th sensor. Here,
(
) is feasible, since we consider heterogeneous sensors. Inspired by the RPEKF [
16,
27], we assume that
is known in advance. Distinct from the RPEKF, we consider the case where there are more than one sensor.
In MATLAB simulations, the maximum sensing range
is set as 15 km for all
. This sensing range information appeared in [
9].
2.1. Least Square (LS), Weighted Least Square (WLS), and Total Least Square (TLS) Methods
Ignoring the noise term in (
1), we have
(
4) leads to
where
We calculate the Least Square (LS) solution [
31,
32] as
The LS was extended to the Weighted Least Square (WLS) estimation in [
34]. In the WLS, measurement noise variance was used to enhance the estimation result, assuming that the bearing measurement noises are sufficiently small. We calculate the WLS solution [
34] as
where the weight matrix is
Here,
is a diagonal matrix whose diagonal terms are given by
in this order. In addition,
is a diagonal matrix whose diagonal terms are given by
in this order. Here,
is
Since is not available in practice, we use for estimation of .
The TLS in [
31,
35] is an extension of (
8). The TLS is given as
Here,
is the smallest singular value of matrix
. Refs. [
31,
35] showed that the TLS has better performance than the LS. From a numerical analyst’s point of view, (
12) indicates that the TLS solution is more ill-conditioned than the LS solution, since it has a higher condition number [
31]. The implication is that errors in the data are more likely to affect the TLS solution than the LS solution. This is particularly true for the worst case perturbations [
31]. However, from a statistician’s point of view, (
12) implies that the TLS method asymptotically removes the bias by subtracting the error covariance matrix (estimated by
) from the data covariance matrix
[
31].
2.2. Bearings-Only Target Localization by Simulating Multiple VMSs
VMS-Based Localization Method
In order to improve the location accuracy compared to the LS-based solutions (
Section 2.1), we develop a novel localization approach of using Virtual Measurement Sets (VMS). Suppose that we use
Q VMSs in total. Each VMS has
K virtual measurements in total. The
q-th VMS (
) has
K virtual measurements, say
. Here, the
k-th virtual measurement in the
q-th VMS is generated by adding a random value in the interval
(
2) to the bearing measurement of the
k-th sensor (
). Recall we assume that
is known a-priori using experiments with the
k-th sensor (
).
In the
q-th VMS,
(
k-th virtual measurement) is used to derive the true bearing
, which is
in (
1). By processing the virtual measurements associated to the
q-th VMS (
), we generate a feasible target estimate, say
.
Among all feasible target estimates (), we select an estimate which has the minimum residual (the difference between the actual and predicted measurements), such that the estimate satisfies the sonar sensing constraints. The selected estimation, say , is used as the final target estimation.
We present the VMS approach in detail. Each VMS has
K virtual measurements in total. Let
denote the virtual measurements of the
q-th VMS (
). Using (
1), we generate the virtual measurements of the
q-th VMS (
) as follows. In the case where
q is one, we use
for all
.
Note that the 1-st VMS uses
, true bearing measurement, for target estimation. In the case where
q is not one (i.e.,
), we use
for all
. Here,
is a virtual noise, generated as a random value in the interval
in (
2). (
15) implies that
is generated by adding a random value in the interval
to the bearing measurement of the
k-th sensor (
).
Here, is the true bearing associated to true relative position.
As
Q increases, it is more probable that we have a VMS
(
) satisfying
for all
. (
17) implies that by adding a virtual noise
, we can have a VMS
, which is disturbed with smaller noise compared to true noise
.
By applying the LS solution (
8) to the virtual measurements (
13) of the
q-th VMS (
), we generate a feasible target estimate, say
. Using (
14),
is the LS solution (
8) associated to true bearing measurements.
Among all feasible target estimates (), we select an estimate which has the minimum residual, such that the estimate satisfies the sonar sensing constraints. The selected estimation, say , is used as the final target estimation.
Let
denote the
n-th element in
. The bearing angle associate to
is
Recall that
is the true bearing measurement of the
k-th sensor. We define the angle deviation
as
Here,
returns the phase angle of a complex number
.
exists between
and
. In the case where
exists between
and
,
is identical to
. Using (
19),
when
.
Using
in (
19), the residual associated to the
q-th (
) VMS is set as
Here,
is used for normalizing the angle deviation term
. In addition,
is set as a constant
in the case where
or the straight line segment connecting
and
crosses an obstacle. Here, we use the fact that underwater sound cannot penetrate an obstacle, such as coastline. This fact can be considered as the obstacle constraint.
In (
20), we set
if
and the straight line segment connecting
and
does not cross an obstacle.
implies that
satisfies the range constraint as well as the obstacle constraint.
In (
20),
is a tuning parameter, presenting the penalty for sonar sensing constraints. In the case where
does not satisfy the constraints,
in (
20) is set as
.
Since provides the minimum residual among all (), is utilized as the target estimate.
However, there is an exception case. If
does not satisfy
for any
, then
is not a viable solution. Thus,
is utilized as the final target estimate. Recall that
is the LS-based solution associated to true bearing measurements. If
satisfies (
24) for every
, then
is utilized as the final target estimate.
We next analyze the computational complexity of the VMS method. Since one derives for all , the VMS method has complexity .
3. MATLAB Simulations
Using MATLAB simulations, we demonstrate the superiority of the VMS-based localization method in
Section 2.2. We consider the marine scenario where multiple sonobuoys are located to measure the bearing of an underwater target’s sound. All sonobuoys are deployed to measure the bearing of the target sound in real time. By processing the bearing measurements of each sonobuoy, our goal is to estimate the target’s 2D position.
We compare among the following methods:
VMS-LS method (VMS using LS solution (
8)).
For generating a virtual noise in (
15), we use
. We use the penalty constant as
in (
20). The maximum sensing range
is 15 km for all
. This sensing range information appeared in [
9].
The performance of a cylindrical array sonar may be different from that of another sonar, depending on the geometry of array segments [
14]. Recall that
denotes the
k-th sonar sensor position (
). We consider heterogeneous sonobuoy sensors as follows. In the case where
k mod 2 is zero, the bearing noise is set as
degrees. In the case where
k mod 2 is not zero, the bearing noise is set as
degrees. This noise statistic is typical for sonobuoy systems [
8].
We use 1000 Monte Carlo (MC) simulations of the scenario. The error in the target estimation at the
i-th MC simulation is
Here, is the target estimation at the i-th MC simulation. In addition, is the true position of the target. Since is not accessible by sonar sensors, is not accessible in practice.
Let
(in meters) denote the mean of
for all
. In other words, we use
Let (in seconds) represent the computation time of running MATLAB for each method in one MC simulation. Large implies that the computational burden of the associated method is large.
3.1. Scenario 1
One MC simulation of Scenario 1 is depicted in
Figure 1. At each MC simulation, we randomly deploy 5 sonobuoy sensors inside the circle with radius 2000 m, whose center is at the origin. This random deployment simulates the case where sonobuoys are deployed by airplanes.
Recall that
denotes the
k-th sonar sensor position (
). For random sensor deployment, we use the following equation.
Here, R is a random number selected in the interval . In addition, A is a random number selected in the interval .
In
Figure 1, each sensor with bearing noise
degrees is marked with a blue asterisk. Furthermore, each sensor with bearing noise
degrees is marked with a green asterisk. In addition, the true target located at (−100, 500) is depicted as a black circle. See that sensors are deployed close to the target, since sensors are deployed inside the circle with radius 2000 m, whose center is at the origin.
Considering Scenario 1,
Figure 2 presents the simulation results of various algorithms, mentioned at the beginning of
Section 3. This figure shows the variation of
with respect to
Q.
and
are plotted at
. In
Figure 2,
implies that
Q VMSs are used in the
method. Note that using (
14),
solution in (
8) is identical to
.
Figure 2 shows that increasing
Q improves the estimation accuracy of
algorithm. The algorithm with the lowest
is
. This verifies the performance of using
proposed in our paper.
Figure 3 shows that as
Q increases, computational burden increases in general.
3.2. Scenario 2
One MC simulation of Scenario 2 is depicted in
Figure 4. In this scenario, the target exists far from the sensors. At each MC simulation, we randomly deploy 5 sonobuoy sensors inside the circle with radius 2000 m, whose center is at the origin. For random sensor deployment, we use (
27). Each sensor with bearing noise
degrees is marked with a blue asterisk. Furthermore, each sensor with bearing noise
degrees is marked with a green asterisk. In addition, the true target position is depicted as a black circle.
Considering Scenario 2,
Figure 5 presents the simulation results of various algorithms, mentioned at the beginning of
Section 3. This figure shows the variation of
with respect to
Q.
and
are plotted at
. Note that using (
14),
solution in (
8) is identical to
. The estimation error of
increased, compared to
Figure 2. This shows that as the distance between the sensors and the target increases, the estimation error increases.
Figure 5 shows that increasing
Q enhances the estimation accuracy of
. The algorithm with the lowest
is
. This verifies the performance of
proposed in this paper.
Figure 6 shows that as
Q increases, computational burden increases in general.
3.3. Scenario 3
One MC simulation of Scenario 3 is depicted in
Figure 7. In this scenario, we deploy many sonobuoy sensors. At each MC simulation, we randomly deploy 10 sonobuoy sensors inside the circle with radius 2000 m, whose center is at the origin. For sensor deployment, we use (
27). Each sensor with bearing noise
degrees is marked with a blue asterisk. Furthermore, each sensor with bearing noise
degrees is marked with a green asterisk. Moreover, the true target position is depicted as a black circle.
Considering Scenario 3,
Figure 8 presents the simulation results of various algorithms, stated at the beginning of
Section 3. This figure indicates the variation of
with respect to
Q. The estimation error of
decreased, compared to
Figure 5. This shows that as the number of sensors increases, the estimation error decreases.
Figure 8 indicates that increasing
Q improves the estimation accuracy of
. The algorithm with the lowest
is
. This verifies the performance of using
proposed in our paper.
Figure 9 shows that as
Q increases, computational burden increases in general.
3.4. Scenario 4
One MC simulation of Scenario 4 is depicted in
Figure 10. At each MC simulation, we randomly deploy 10 sensors inside the circle with radius 2000 m, whose center is at the origin. For sensor deployment, we use (
27). Each sensor with bearing noise
degrees is marked with a blue asterisk. Each sensor with bearing noise
degrees is marked with a blue asterisk. See that we use sensors with small bearing noise, compared to Scenario 3. Moreover, the true target position is depicted as a black circle.
Figure 11 presents the simulation results of various algorithms, stated at the beginning of
Section 3. This figure indicates the variation of
with respect to
Q. The estimation error of
decreased, compared to
Figure 8. This shows that as we use sensors with small bearing noise, the estimation error decreases.
Figure 11 indicates that increasing
Q improves the estimation accuracy of
. The algorithm with the lowest
is
. This verifies the performance of using
in our paper.
Figure 12 shows that as
Q increases, computational burden increases in general.
4. Conclusions
This paper addresses localization of an underwater target based on bearings-only measurements of heterogeneous sonobuoys. Sonobuoys have very low bearing accuracy, such as 10 degrees [
8]. We tackle bearings-only localization using heterogeneous sonobuoy sensors, considering the sensing constraints of a sonar sensor. Using MATLAB simulations, we demonstrate the outperformance (considering both location accuracy and time efficiency) of the proposed
, by comparing it with other localization methods.
Using MATLAB simulations, with large Q shows the best performance among all algorithms, considering both time efficiency and localization accuracy. Increasing Q in enhances the localization accuracy, while increasing the computational load. In the future, we will do experiments using real sonar sensors, in order to verify the performance of the proposed localization more rigorously.
The authors of [
36] presented the Gauss Newton (GN) algorithm applied to bearings-only target estimation. The GN method is an iterative method to derive the target estimation. The convergence of the GN estimator is sensitive to the initial guess and the step size [
38]. To further decrease the localization error, we can run the GN method which is initialized with the VMS-based solution. In the case where the distance between a sensor and the target estimation from the GN method is bigger than the maximum sensing range, we can get out of the iteration of the GN method. In this way, we avoid a target solution which is too far from a sensor.
This paper addresses a VMS localization method, considers 2D bearings-only localization. In the future, the proposed VMS approach will be applied in various localization schemes, such as 3D angle-of-arrival localization [
41], time-of-arrival localization [
42,
43,
44], time-difference-of-arrival localization [
45,
46,
47,
48], or doppler-bearings localization [
49,
50].