Underwater Ambiguity Elimination Method Based on Co-Prime Sensor Array

Recently, the direction of arrival estimation with co-prime arrays has gradually been applied in underwater scenarios because of its significant advantages over traditional uniform linear arrays. Despite the advantages of co-prime arrays, the spatial spectra obtained directly from conventional beamforming can be degraded by grating lobes due to the sparse spatial sampling in passive sensing applications, which will seriously deteriorate the estimation performance. In this paper, capon beamforming is applied to a co-prime sensor array as a pretreatment before high-resolution direction of arrival (DOA) estimation methods. The amplitudes extracted from the beam-domain outputs of two subarrays and the phases extracted from the cross-spectrum of the spatial spectrum are exploited to suppress the spurious peaks in beam patterns and eliminate ambiguities. Consequently, interference can be further mitigated, and the performance of high-resolution DOA methods will be guaranteed. Simulations show that the method proposed can improve the reliability and accuracy of DOA estimation with great value in practice.


Introduction
Direction of arrival (DOA) estimation is an important array signal processing technique that finds broad applications in underwater passive sonar systems. An important goal for DOA estimation is to be able to locate closely spaced sources in the presence of considerable noise. Utilizing this technique, passive sonar systems can detect, localize, and track underwater acoustic sources by monitoring with an array with spatially separated sensors. DOA estimation has been an active research area during the last two decades. Conventional techniques for DOA estimation, such as beamforming, have been widely applied in the field of passive sonar. However, they can only provide limited angular resolution. Well-known methods to improve angular resolution are implementation of longer arrays or super-resolution subspace-based methods. Unfortunately, the applicability of both techniques is constrained when applied in a complex underwater environment. Firstly, the implementation of long towed linear arrays (TLA) dedicated to platforms such as unmanned underwater vehicles (UUVs) is frequently impossible because of structure constraints (e.g., heavy cables). Secondly, in realistic underwater environments, the received signals are temporally correlated as a result of multi-path arrivals. Under these conditions, the performance of subspace-based methods, such as multiple signal classification (MUSIC), will significantly degrade. Therefore, how to approach the problems mentioned above has attracted much interest in recent years [1,2].

Array Structure, Signal Model, and Beamforming
In this paper, the CA structure is composed of two co-prime subarrays. The positions of the 1st sensor in the two subarrays coincide, and the sensors in each subarray are arranged evenly and equidistantly. The inner spacing of sensors in the two subarrays is co-prime. As shown in Figure 1, the co-prime pair of the two subarrays is (M, N), where M and N are both positive integers (M, N ≥ 2), and the array expansion factor is α. In subarray one, the number of sensors is αN + 1 and the inner spacing is Md. In subarray two, the number of sensors is αM + 1 and the inner spacing is Nd. d = λ/2 is the half-wavelength of the impinging signal. Conventional beamforming (CBF) is employed as the signal processing method in this paper for conciseness, which can be further upgraded to MVDR methods or other robust beamforming algorithms. The scanning range is θ ∈ (−π/2, π/2). Assume that all the targets are in the far-field and the sources are all narrowband signals. As is shown in Figure 2, the directions of k = 1, ...K targets are θ k , k = 1, ..., K. Suppose s k (t) is the signal model of the kth target, then s k (t) = A k (t)e −j(2πct/λ k (t)+ϕ k (t)) k = 1, 2, ..., K In Equation (1), A k (t) is the amplitude of s k (t), ϕ k (t) is the phase of s k (t), and c is the velocity of sound. In an underwater sensing environment, all the parameters of the signals s k (t) are assumed to be random in this signal model. As a matter of fact, the original source signal parameters are basically different from each other. The amplitudes at the receiving point are related to the propagation distance from the sources, so A k (t) is random. When narrowband coherent sources exist or the multi-path effects appear, the wavelength of received signals is influenced, and then λ k (t) is random. For the influence of the environment, the phase ϕ k (t) of the received signals is also random, which is in the range of (0, 2π) [29,30]. For the amplitude and phase of the far field signal at time t and time t + τ, we have A k (t + τ) ≈ A k (t) and ϕ k (t + τ) ≈ ϕ k (t). Then the relationship between s k (t + τ) and s k (t) is In Equation (2), s k (t) is the kth signal received at time t, and the the impinging signal angle is θ k , k = 1, ...K. Then the received signal x n (t) and x m (t) can be presented as Equation (3), In Equation (3), n n (t) is the additive white Gaussian noise with unknown variance at the nth receiver sensor of subarray one, n = 1, ..., αN + 1 and n m (t) is the noise at the mth receiver sensor of subarray two, m = 1, ..., αM + 1. τ nk is the relative time lag for the kth signal at the nth sensor in subarray one, and τ mk is the relative time lag for the kth signal at the mth sensor in subarray two. With Equation (2), Equation (3) can be expressed as In Equation (4), θ k is the direction of the kth signal, k = 1, ..., K. Then, the received signal vectors of two subarrays can be expressed as In Equation (5), x 1 (t) = [x 11 (t), ..., x 1(αN+1) (t)] T and x 2 (t) = [x 21 (t), ..., x 2(αM+1) (t)] T are the vector form of signal received at time t. n 1 (t) = [n 11 (t), ..., n 1(αN+1) (t)] T and n 2 (t) = [n 21 (t), ..., n 2(αM+1) (t)] T are the corresponding noise at time t. s(t) = [s 1 (t), ..., s K (t)] T is the signal vector at time t. A 1 (θ) = [a 11 (θ) ... a 1K (θ)] and A 2 (θ) = [a 21 (θ) ... a 2K (θ)] are the steering matrices of two subarrays. a 1k (θ) = [1 e −j2πMdsinθ k /λ ... e −j2παNMdsinθ k /λ ] T and a 2k (θ) = [1 e −j2πNdsinθ k /λ ... e −j2παMNdsinθ k /λ ] T are the steering vector, k = 1, ..., K.
The first sensor of the two subarrays is shared at the coordinate origin as reference sensor. The relationship between the sensors in two subarrays and the reference sensor is shown in Figure 3. There is a wave-path difference between the two sensors in each subarray. The wave-path difference between the nth sensor and the reference sensor of subarray one is τ n = nMdsinθ k /c, n = 1, ..., αN, and the wave-path difference between the mth sensor and the reference sensor of subarray two is τ m = mNdsinθ k /c, m = 1, ..., αM. Different wave-path differences lead to different phase lags In Equation (6), f 0 is the center frequency of the signal, and because of its narrowband property, f ≈ f 0 . Then the phase lags ϕ n and ϕ m can be expressed as From Equation (7), when the phase lags between two sensors are known, the DOA of the signals can be estimated. In beamforming processing flow, in order to receive the signals from direction θ k , k = 1, ...K and suppress signals from other directions, a weight vector w is needed to make the main beam pointing to the direction of θ k . Then the weighted outputs of the two subarrays are In Equation (8), y 1 (t) is the weighted outputs of subarray one, and y 2 (t) is the weighted outputs of subarray two. When the ideal weight vectors w 1 and w 2 are fixed, the power P 1 (w 1 ) and P 2 (w 2 ) of the outputs in direction of θ k can be calculated: In Equation (9), R 1x is the covariance matrix of x 1 (t) in size of M × M, and R 2x is the covariance matrix of x 2 (t) in size of N × N. Because w 1 and w 2 are the weight vectors pointed to direction θ k , then Equation (9) can be expressed as For Equation (10), in fact we do not know where the real θ k is, and we wish to estimate it. We can search within the range θ ∈ (−π/2, π/2) in order to get the power spectrum: In Equation (11), the CBF takes w 1 (θ) = a 1 (θ k ) and w 2 (θ) = a 2 (θ k ). Then the spatial power spectra of the two subarrays are obtained: For Equation (12), in real-world scenarios, the finite sample estimates of the array covariance matrix R 1x and R 2x are formed: In Equation (13), p is the snapshot index, x 1 (z) and x 2 (z) are the finite discrete samples, z = 1...p. Therefore, the spatial spectra of the two subarrays, P 1 (θ) and P 2 (θ), can be obtained. However, due to the sparse spatial sampling, the spectra of the two subarrays are severely affected by grating lobes, which will influence the results of DOA estimation.

Single-Target Case
Let us analyze the co-prime properties of subarrays and the phase difference of the adjacent sensors, ∆ϕ, in each subarray: In Equation (14), mod is the modulus operator, ∆D = Qd, and Q is the co-prime factor M or N. k is an integer, and k ∈ [−∆ϕ/2π, D/λ − ∆ϕ/2π] since θ ∈ (−π/2, π/2) and 0 ≤ sin θ ≤ 1. For a true signal direction θ t and a grating lobe direction θ g , the relationship between their directions can be derived from Equation (14). Equation (15) shows the possible intervals between true signal directions and grating lobe directions of two subarrays, In Equation (15), R 1 and R 2 are both integers, 1]. Therefore, we have R 1 /N = R 2 /M. Because of the co-prime properties of (M, N), there will be only one pair of (R 1 , R 2 ) to get the true peak of the spatial spectrum, except when R 1 = R 2 = 0, [24]. Hence, the beam-domain outputs of the two subarrays can be processed by the product processor to get a common peak, as is shown in Figure 4 [21]. In each bearing, the beam outputs of different subarrays are correlated and conjugated, i.e., In Equation (16), y 1 (θ) and y 2 (θ) are the beam outputs of subarrays. Take the minimum of |P 1,2 (θ)| and |P 2,1 (θ)|, i.e., P(θ) = min(|P 1,2 (θ)| , |P 2,1 (θ)|) (17) In Equation (17), P 1,2 (θ) and P 2,1 (θ) are the conjugate correlation spatial pseudo spectrum. The highest peaks are selected on the spectra that are associated with the target directions. Then the directions are obtained, which can be represented as (θ 1 , θ 2 , ..., θ m ), where m is a positive integer. When there is only one target present, i.e., m = 1, there will be only one common spectral peak in the true signal direction θ t according to the deduction after Equation (14). The grating lobes can be completely suppressed in single-target scenarios with the methods in this subsection.

Multiple-Targets Case
For multi-targets, in addition to the true target direction, different targets in the spatial spectra of the two subarrays may generate grating lobes in the same direction, so the method for a single-target case will not be valid in some special cases. Assume that there are two targets, with target one in direction θ 1 and target two in direction θ 2 . Two subarrays are applied to carry out beamforming simultaneously, and the grating lobes generated are related to the target directions. Therefore, the spectrum peaks obtained by subarray one include the true target directions θ 1 , θ 2 , and a group of grating lobe directions θ 1g1 , ..., θ 1gu (u is a positive integer). Similarly, the spectrum peaks obtained by subarray two include the true target directions θ 1 , θ 2 , and another group of grating lobe directions θ 2g1 , θ 2g2 , ..., θ 2gv (v is a positive integer). When the two groups of grating lobes overlap, i.e., θ 1gu = θ 2gv = θ s , there will be a spurious peak in the direction θ s , which cannot be eliminated. For a spurious peak in the direction θ s , the following conditions are satisfied: In Equation (18), R 3 and R 4 are also integers, 1]. When the equations in Equation (18) are valid, the spurious peak θ s caused by the overlapping of the grating lobes cannot be eliminated, and the directions of spurious peaks satisfy Equation (19): When there are multiple targets, this problem will still exist, and with an increase in the number of targets, the number of spurious peaks will continue increasing. The amplitudes of spurious peaks vary in spatial spectrum, and higher ones may even mask true signals. In underwater scenarios, the amplitudes and phases of different signals in different directions to the reference sensor are generally not the same as is mentioned in Section 2. Therefore, the amplitude information of beam-domain outputs and the phase information of the spatial cross-spectrum of the two subarrays can be utilized as the theoretical basis for evaluating whether there is a true signal or not. The theories will be elaborated in the next several paragraphs.
For the beam-domain output of each subarray, it contains the waveform information of the signals. In ideal scenarios, for the direction of true target θ t , the amplitudes of normalized outputs, A 1 (θ t ) and A 2 (θ t ), of the two subarrays are the same, i.e., A 1 (θ t ) = A 2 (θ t ). However, for the direction of spurious peak θ s , the normalized outputs, A 1 (θ s ) and A 2 (θ s ), of the two subarrays are generally different, i.e., A 1 (θ s ) = A 2 (θ s ). In practical scenarios, due to the ambient noise and interference, the actual beam-domain outputs,Â 1 (θ t ) andÂ 2 (θ t ), of subarrays deviate from the ideal values A 1 (θ t ) and A 2 (θ t ). The amplitude information obtained in the true target direction is unequal, but the deviation is small, i.e.,Â 1 (θ t ) ≈Â 2 (θ t ). Therefore, the spurious peaks can be further distinguished by exploiting this feature. Fourier Transform is applied to the beam-domain output of the two subarrays: (20) In Equation (20), it is assumed that the noise is white and not correlated with signals. Then the cross-spectrum of the two subarrays G 12 (ω) is obtained.
In Equation (21), due to the properties of noise, N 1,2 (ω) = 0. The phase information of the spatial cross-spectrum is extracted, i.e., ϕ = ϕ 2 − ϕ 1 . The element spacings are different in the two co-prime subarrays; hence, the delays from element to element should be different as well if a plane wave is observed. Ideally, in the true target direction, the phase differences extracted from the beam-domain outputs of the two subarrays are the same, while in the direction of the spurious peaks, the phase information is generally different. Thus, to deal with the spurious peak problem, the phase feature can also be exploited. Therefore, there are two processing stages to get the true peaks. The inputs and outputs in each stage are shown in Figure 5.

Amplitude Selection
By deduction of the above sections, the signal amplitude characteristics are introduced to identify the spurious peaks. In consideration of the amplitude characteristic of signals in real-world scenarios, we make an amplitude hypothesis. For the true target directions, the amplitude differences of the normalized outputs between the two subarrays are small; while for the spurious peak directions, the amplitude differences of the normalized outputs between the two subarrays are big. This is true for ideal sensors, but in practice each sensor has its own gain, which is direction dependent. In this paper, the sensors are assumed to have ideal omnidirectional gain patterns. The outputs of subarrays are normalized to the reference sensor, which will only eliminate the noise. As a matter of fact, the possibility that two sources with the same amplitude, phase, and distance appear in different directions is low [29]. Therefore, this case will not be included in this paper. Because of this, there will be an amplitude threshold G 1 to evaluate the peaks (θ 1 , θ 2 , ..., θ m ) from the product selection. The directions that cannot meet the threshold will be eliminated, as is shown in Figure 6. The normalized beam-domain outputs in the directions of (θ 1 , θ 2 , ..., θ m ) are Then the amplitude information A 1 (θ j ) and A 2 (θ j ) can be obtained from the outputs y 1 and y 2 by The amplitude ratio ρ(θ j ) of the beam-domain output of the two subarrays can be calculated, i.e., In Equation (24), ρ(θ j ) is the amplitude ratio with ideal value ρ(θ j ) ideal = 1. In the real world, there can be an amplitude threshold G 1 that is a little bit less than 1. The spurious peaks will be eliminated when the amplitude ratio ρ(θ j ) is lower than G 1 . With the amplitude selection, true directions (θ 1 , θ 2 , ..., θ n ), n ≤ m can be found.

Phase Selection
After the amplitude selection, there may still be some spurious peaks that satisfy the amplitude threshold, and the phase information of the spatial cross-spectrum can be exploited. In consideration of the phase characteristics of signals in real-world scenarios, we make a phase hypothesis. For the true target directions, the phase differences extracted from the beam-domain outputs of the two subarrays are small; while for the spurious peak directions, the phase differences extracted from the beam-domain outputs of the two subarrays are big. The phase differences extracted above correspond to the spatial cross-spectra of the two subarrays. Since the sensor displacements are different, the delays will be different in the two co-prime array components. Take two targets in directions of 0 • and 10 • for example, the phase is extracted in Table 1, and Figure 7 also illustrates the case. Because of this, there can be a phase threshold G 2 to evaluate the peaks (θ 1 , θ 2 , ..., θ n ) from the amplitude selection. The directions that cannot meet the threshold will be eliminated, as is shown in Figure 8. The normalized beam-domain outputs in the directions of (θ 1 , θ 2 , ..., θ n ) are   The phase information φ 1 (θ j ) and φ 2 (θ j ) can be obtained from the outputs y 1 (θ j ) and y 2 (θ j ) by The phase difference ∆φ(θ j ) of the beam-domain outputs of the two subarrays is In Equation (27), ∆φ(θ j ) is a phase difference corresponding to the phase information of the spatial cross-spectrum of the two subarrays with ideal value ∆φ(θ j ) ideal = 0, which can be derived from Equation (21). In the real world, there will be a phase threshold G 2 a little bit more than 0. The spurious peaks can be further eliminated when ∆φ(θ j ) is higher than G 2 . In this way, true target directions (θ 1 , θ 2 , ..., θ k ), k ≤ n can be obtained in phase selection.

Processing Flow
In this paper, co-prime array geometry is employed. The processing flow is shown in Figure 9. To start with, the subarray CBF is utilized to obtain the initial directions. Then, the beam-domain outputs of subarrays in each bearing are conjugated and multiplied to obtain the peaks of the correlation spectrum with the ambiguity being eliminated. Based on this, the amplitude information from two subarrays is exploited to further discriminate the spurious peaks with the amplitude threshold. After that, the phase information is exploited to pick out spurious peaks that have not been distinguished by the amplitude selection stage. Finally, all spurious peaks in the spatial spectrum can be identified, and the ambiguity will be eliminated. Following the steps in the processing flow, the true directions of targets are more deterministic, which is beneficial to interesting fields such as UUV sonar arrays [31]. This paper focuses on the ambiguity elimination problem rather than an underwater detection problem, so it is assumed that the SNR is enough to execute the processing flow. The amplitude and phase information can be estimated in a reasonable range.

Simulation and Discussion
The parameters of CA are chosen first: the co-prime pair is (M, N) = (3, 4), and the extension factor is α = 3. Thus, the sensor number of the two subarrays is (M 1 , M 2 ) = (13, 10). The total number of sensors generated by the subarrays is M SUM = 19, and the sensor number of ULA with the same aperture is M ULA = 37. The resolution of a certain direction is directly related to the changing rate of the steering vector near the direction, which can be given in Equation (28): The Rayleigh Criterion for the ULA is 3.1 • . Assume that all the sensors are identical, reciprocal, and omnidirectional. The center frequency of signals is set at f 0 = 1.5 kHz, and the sampling frequency is f s = 15 kHz. The underwater sound velocity is c = 1500 m/s, and the sensor spacings in subarrays are ∆d 1 = 1.5 m (1.5 λ) and ∆d 2 = 2.0 m (2.0 λ), respectively. CBF is adopted, and the bearing range is θ ∈ (−π/2, π/2).

Processing Method Simulation
Assume a scenario that K = 2 far-field point targets are simulated. The rationality and correctness of the algorithm are verified by changing the parameters of targets. The amplitude threshold and the phase threshold are set at (G 1 , G 2 ) = (0.95, 4). The solid lines represent the true target directions, and the dotted lines represent the spurious peak directions.

Product Processing
The directions of two targets are (θ 1 , θ 2 ) = (0 • , 12.5 • ). The SNR level is (SNR 1 , SNR 2 ) = (0, 0), and the phases of targets are set randomly. The number of snapshots is 1000, and the simulation result is shown in Figure 10a. In this situation, the spurious peaks can be eliminated by product processing, but when the direction of targets changes to (θ 1 , θ 2 ) = (0 • , 10 • ), the spurious peaks cannot be eliminated, as is shown in Figure 10b. Since the main beamwidth of the beamformer is about 3.1 • , the second target is located out of the beam range.

Amplitude Selection
The directions of two targets are (θ 1 , θ 2 ) = (0 • , 10 • ). The SNR level is (SNR 1 , SNR 2 ) = (0, 5), and the phases of targets are set randomly. The number of snapshots is 1000, and the simulation result is shown in Figure 11a. In this situation, the spurious peaks cannot be distinguished by product processing. With the amplitude selection, the spurious peaks can be identified. However, when two targets have similar SNRs, for example when (SNR 1 , SNR 2 ) = (0, 0), the amplitude selection is not as satisfactory, as is shown in Figure 11b. Additionally, individual gain patterns make this even more difficult when the sensors are not directionless.

Phase Selection
The directions of two targets are (θ 1 , θ 2 ) = (0 • , 10 • ). The SNR level is (SNR 1 , SNR 2 ) = (0, 0), and the phases of targets are set randomly. The number of snapshots is 1000, and the simulation result is shown in Figure 12a. In this situation, the spurious peaks cannot be distinguished by amplitude selection. With phase selection, the spurious peaks can be identified, as is shown in Figure 12b. However, when a true signal overlaps with a spurious peak, or two signals overlap with each other, the processing does not work well enough. The resolution limit of the array is 3.1 • , and the sources are in the same direction (within the main beam). The mainlobe interference elimination, blind source separation, and independent component analysis can be exploited to solve this problem in further work.

Comparison
After the processing flow is simulated, the spatial spectrum obtained by CBF of whole array, the spatial spectrum obtained by MUSIC of whole array, the proposed method with co-prime array in this paper, and the half-wavelength ULA with the same aperture are compared.
As is shown in Figure 13, comparing with the CBF method applied to the whole array, the proposed method in this paper can effectively distinguish spurious peaks and accurately estimate the true direction of the targets. The MUISC algorithm is also applied to the whole array under the condition that the number of signals is known, and the spurious problem still exists, which is similar to the situation when CBF is applied to the whole array. Comparing with the half-wavelength ULA with same aperture, the spatial spectral peaks estimated by the proposed method are close to the width of the ULA results, which can meet the estimation requirements.

The Selection of Thresholds
The selection thresholds G 1 and G 2 are highly related to false alarms and underreporting, which has considerable influence on the PoS of the proposed method in this paper. Before discussing the PoS for the two targets case and the multiple targets case, the selection of thresholds should be simulated. As is stated in Section 3, the amplitude threshold G 1 is a little bit less than 1, and the phase threshold G 2 is a little bit more than 0, theoretically.
The conditions of simulation remain unchanged, and the parameters of targets remain the same as in the last section. In the processing flow, the amplitude selection is achieved before the phase selection, and the relationship between PoS and G 1 is simulated. The range of amplitude threshold G 1 is adjusted from 0.8 to 1 with the interval 0.02. For each G 1 , 500 Monte Carlo experiments were conducted, and the PoS versus G 1 is obtained. As is shown in Figure 14a, the PoS goes steadily up with the increase of G 1 . The amplitude threshold G 1 gradually eliminates the angle which is not qualified to it. When G 1 = 0.94, the PoS reaches the max with PoS = 81.2% and the value of G 1 near 0.95 is recommended to be the amplitude threshold. In this range, the amplitude selection can guarantee at least PoS = 70%. Then, the PoS goes down until PoS=0 when G 1 = 1. In this situation, the threshold is too high and approaches the ideal circumstance, which is not a good range for the selection of G 1 . After selection of the amplitude threshold G 1 , the selection of phase threshold is simulated. The conditions of simulation remain unchanged, and the parameters of targets remain the same. The amplitude G 1 is set at 0.95, as a reliable value from the simulation result. The range of phase threshold G 2 is adjusted from 0 to 10 with the interval 0.5. For each G 2 , 500 Monte Carlo experiments were conducted, and the PoS with the corresponding G 2 is obtained. As is shown in Figure 14b, the PoS increases sharply when G 2 ≤ 3. When G 2 is in the range of 3.5 to 4.5, the PoS reaches the top with PoS = 98%, and the values near here are recommended to be the phase threshold. Then, the PoS decreases slowly, for G 1 has provided enough guarantee for the phase selection. The simulations are based on two targets here and in practice, and the two thresholds can be set artificially according to empirical experience or adaptive methods.

PoS for Two Targets
In order to verify the ambiguity elimination effect when the directions of targets change, the SNR level is fixed at (SNR 1 , SNR 2 ) = (0, 0). The direction of target one remains unchanged θ 1 = 0 • , and the direction of target two changes within the range θ 2 ∈ (−90 • : 90 • ). In each bearing, 1000 Monte Carlo experiments are conducted, and the PoS versus corresponding θ 2 is obtained.
As is shown in Figure 15, when target two is located in the direction of θ 2 ∈ (−3 • : 3 • ), the PoS decreases significantly. This is caused by the overlap of the two signals, which has been discussed in Section 3. When θ 2 ∈ (−89 • : 90 • ), the end-fire direction of the array, the DOA ability of the array decreases obviously. While in most directions, the PoS remains stable with little fluctuation.

PoS for Multiple Targets
As is shown in Section 3, when the number of targets increases, the number of spurious peaks will continue increase. The simulation of PoS versus the number of targets in a certain range is carried out. The SNR of targets are fixed between 0 dB and 3 dB. The direction of target one remains unchanged at 0 • , and the direction of other targets changes with different angle intervals from 3 • to 15 • , as is shown in Figure 16a. In each case, 1000 Monte Carlo experiments were conducted, and the PoS was obtained. As is shown in Figure 16b, when the interval of targets is 3 • , the PoS remains stable near 100% with the number of targets under 8. However, as the number of targets increases, the PoS decreases, and when the number of targets reaches 11 the method becomes invalid. This is because when more targets appear, the component of spurious peaks becomes more complex. When the amplitude and phase information extracted from this spurious peak are similar to that of true targets, the existence of targets cannot be evaluated, i.e., the spurious peak will be identified as a true peak. When the spurious peak and the true peak overlap in the same direction, the true peak would also be identified as a spurious peak. When the interval of targets increases, the number of identifiable targets decreases. Therefore, the method proposed in this paper has better performance in scenarios with fewer targets. When the scenarios become complicated, it is recommended to apply and combine more methods like independent component analysis to accomplish the estimation.

Statistical Performance Analysis
The root-mean-square error (RMSE) of the method proposed in this paper versus SNR and the number of snapshots are next investigated over 1000 independent Monte Carlo simulations. The RMSE of DOA estimation is calculated as In Equation (29), ψ is the number of Monte Carlo simulations, κ = 1, ..., ψ, K is the number of sources, θ k is the kth real direction of the signal, andθ kκ is the estimated direction of the kth source in the κth Monte Carlo simulation, k = 1, ..., K.
The directions of targets are fixed at (θ 1 , θ 2 ) = (0 • , 10 • ). The number of snapshots is 1000, and the SNR level is adjusted from −15 dB to 0 dB. The thresholds are set at (G 1 , G 2 ) = (0.95, 4). The RMSE versus SNR with the proposed method is shown in Figure 17a. The RMSE appears to decrease gradually as the SNR of the signals increases. When the SNR level is low, the presence of the noise seriously affects the accuracy of DOA estimation. When the SNR level is sufficiently high, the method proposed in this paper presents the similar accuracy with little RMSE. When the SNR level reaches −9 dB, the RMSE remains at 0. The RMSE versus snapshots is then simulated. The SNR level of targets is fixed at 0 dB, i.e., (SNR 1 , SNR 2 ) = (0, 0), and the number of snapshots is adjusted from 100 to 1000. As is shown in Figure 17, when the number of snapshots increases, the RMSE descends gradually. When there are sufficient snapshots, the RMSE remains stable near 0. However, compared to the RMSE versus SNR level, RMSE versus snapshots has less influence on the proposed method. From the statistical performance analysis simulation, the proposed method can guarantee the DOA estimation performance under sufficient SNR and achieve similar performance with insufficient snapshots.

Prospects For Future Work
The method in this paper is proposed to solve the problem of grating lobes because of the sparse spatial sampling in passive sensing applications. In practical sonar system applications, how to quickly eliminate ambiguity is a practical problem. Strong demand for wide-band signal processing and vector hydrophone array processing with sparse arrays under random interferences still exists. In this paper, the reduction of mutual coupling and relevant material effect are qualitatively discussed, which can be verified with more experiments. The threshold selection in different cases depends on the parameter setting, and the rules will be further studied.

Conclusions
CAs are gradually applied in underwater scenarios because they have significant advantages over traditional ULAs. For the same sensor number, CAs provide better accuracy and higher resolution. While for the same aperture, CAs can reduce the physical cost. Because of the expansion of inter-spacing, the mutual coupling can be reduced. Despite the advantages of CAs, the spatial spectra obtained directly from CBF can be degraded by grating lobes due to the sparse spatial sampling in passive sensing applications, which will seriously deteriorate the estimation performance. Therefore, the ambiguity elimination based on CAs is worthy of study. The ambiguity elimination in this paper is solved taking advantage of the co-prime property of subarrays in CAs. The amplitude and phase information of beam-domain outputs are exploited based on this idea.
The method proposed in this paper is applied to the CA as a pretreatment before the high-resolution DOA estimation algorithms like MUSIC, etc. In the processing flow, the subarray CBF is firstly utilized to obtain the initial direction set. Then, the beam-domain outputs of subarrays in each bearing are conjugated and multiplied to obtain the peaks of the correlation spectrum. On the basis of this, the amplitude information from two subarrays is exploited to further discriminate the spurious peaks with the amplitude threshold. After the amplitude selection, the phase information is exploited to identify spurious peaks that are not distinguished by the amplitude selection. Finally, all the spurious peaks in the spatial spectrum can be identified, and the ambiguity is eliminated.
Simulation results validated the effectiveness of the proposed method in this paper. The selection of the threshold is discussed, and the recommended values are given. The statistical performance analysis of the proposed method is studied and proofed to guarantee the performance. To summarize, the method with CA utilizes fewer sensors to estimate the DOA accurately and can eliminate the ambiguity caused by at least three targets. When the angle range of interest is given, this method performs better. In practical scenarios like underwater applications, the proposed method has great practical significance. It can effectively promote the related array signal processing algorithm in reality, save costs, and is easy to implement.

Conflicts of Interest:
The authors declare no conflict of interest.