DOA Estimation for Sources with Large Power Differences

Sources with large power differences are very common, especially in complex electromagnetic environments. Classical DOA estimation methods suffer from performance degradation in terms of resolution when dealing with sources that have large power differences. In this paper, we propose an improved DOA algorithm to increase the resolution performance in resolving such sources.*e proposedmethod takes advantage of diagonal loading and demonstrates that the invariant property of noise subspace still holds after diagonal loading is performed. We also find that the Cramer–Rao bound of the weak source can be affected by the power of the strong source, and this has not been noted before. *e Cramer–Rao bound of the weak source deteriorates as the power of the strong source increases. Numerical results indicate that the improved algorithm increases the probability of resolution while maintaining the estimation accuracy and computational complexity.


Introduction
Direction-of-arrival (DOA) estimation has been extensively investigated due to its wide applications in radar, sonar, wireless communication, and navigation [1][2][3][4]; e.g., DOA estimation is of significance for source location in MIMO radar [5], and the issue of DOA estimation for noncircular sources is valuable in practical communications [6]. Subspace-based methods are preferred among numerous DOA estimation methods due to their high superresolution and ease of implementation [7][8][9]. For passive radar, the echo from the target may be several dB lower than the direct signal or intentional interference. e estimation performance deteriorates dramatically if the weak target is masked by the strong target [10]. However, the classical MUSIC, Root-MUSIC, and Capon algorithms suffer from severe performance degradation when resolving closely spaced sources with large power differences [11,12], especially in the threshold region [13].
ere are two main types of methods for DOA estimation with sources that have power differences. One type involves estimating the weak sources after mitigating the strong interference. Generally, this type of method requires prior knowledge of the strong interference. With prior DOA knowledge of the strong interference, the interference jamming method (IJM) is proposed in [14], where strong interferences are eliminated by transforming the M × Ndimensional array manifold matrix into M × (N − j) dimensions, in which M is the number of array elements, N is the number of incident sources, and j is the number of strong interferences. However, the IJM method suffers from the loss of array aperture because the steering vector belonging to the strong interferences is removed from the array manifold matrix. Another DOA estimation algorithm based on the elimination of interference is proposed in [15], and it shows better performance than the jamming jam method (JJM). Its performance is limited when the power difference between two sources is not large. Provided with accurate prior DOA knowledge, the constrained MUSIC approach [16], the weighted prior-MUSIC approach [17], and the improved L1-SVD method [18] achieve accurate DOA estimation for weak sources, whereas the performance degrades if the prior DOA knowledge of strong interference is inaccurate. Given the number of strong interference, the DOA estimation of weak sources in the presence of strong jamming has been studied based on the extended noise subspace [19]. Prior knowledge of strong interference is needed; otherwise, the extended noise subspace is unavailable. Additionally, a modified spectrum using eigenvector projections has proven valid in dealing with sources exhibiting large power differences based on the number of known sources [20]. However, prior knowledge regarding the input sources is difficult to acquire, especially in noncooperative situations, and limits the applications of methods that need prior knowledge. e other type of DOA estimation method for sources with power differences does not need prior knowledge. However, this type mostly employs numerous iterations at the cost of increased computational complexity. An oblique projection filter can be used to constrain sources from undesired directions without any prior knowledge of strong interference [21]. Multiple singular value decomposition is needed in this method, and this increases the complexity of the algorithm. A hybrid method without any prior knowledge utilizes compressive sensing (CS) to estimate the DOAs of strong sources and then removes the effects of strong sources via orthogonal complements [22]. However, it suffers from resolution loss if the DOAs of strong sources are not accurately estimated. e number of sources and the DOAs of weak sources can be jointly obtained from the Eigenbeam mCapon method [23]. Repeated iterations are employed at the cost of considerable computation complexity, and this is also true for the RELAX algorithm [24]. Recently, sparse theory and sparse arrays have attracted much attention and have provided a new method for DOA estimation [6]. In [25], the authors analyze the loss of weak targets in the sparse asymptotic minimum variance (SAMV) method and propose the robust sparse asymptotic minimum variance (RSAMV) algorithm. A robust orthogonal projection method based on an optimized sparse array is proposed without prior knowledge in which the DOAs of strong and weak targets are precisely estimated [26], but many iterations are still needed.
Both iterations and prior knowledge are not required for the DOA estimation algorithm, which utilizes the invariant property of noise subspace (IPNS) [27,28]. In this paper, we propose an improved DOA estimation algorithm to increase the resolution performance based on the IPNS. e improved algorithm takes advantage of diagonal loading, and we demonstrate that the IPNS still holds after diagonal loading. e Cramer-Rao bound (CRB) of the weak source is studied in the presence of a strong source to examine the influence caused by the power of the strong source. e numerical results indicate that our approach is highly effective in solving two closely spaced sources with large power differences. e rest of this paper is organized as follows. In Section 2, the signal model is introduced. In Section 3, an improved algorithm for DOA estimation for sources with large power differences is proposed. We demonstrate that the IPNS is still valid after diagonal loading and analyze the complexity of the improved algorithm. In Section 4, the CRBs of weak sources under different powers of strong sources are studied. Numerical simulations and discussions are presented in Section 5. Finally, we conclude our work in Section 6.

Signal Model
Consider P independent narrowband sources impinging on a uniform linear antenna array from the far-field, as shown in Figure 1. e distance between any two antenna elements is d � λ 0 /2, where λ 0 is the carrier wavelength. e number of receiver antenna elements is M (M > P). Suppose that the ith source impinges on the uniform linear antenna array from direction θ i . e DOAs of all sources are defined as e output of the uniform linear antenna array at time t can be expressed as follows: where T is the received signal vector and (·) T denotes the matrix transpose operation. s(t) � [s 1 (t), s 2 (t), . . . , s P (t)] T is the vector of incident signals. n(t) � [n 1 (t), n 2 (t), . . . , n M (t)] T is the additive and independent Gaussian noise vector whose covariance and mean are σ 2 n I (I denotes the identity matrix) and zeros, respectively. A � [a(θ 1 ), a(θ 2 ), . . . , a(θ P )] is the array e covariance matrix of x(t) is given by the following: where E[·] and (·) H denote the expectation and the conjugate transpose operators, respectively. R s � E s(t) s(t) H ] is the covariance matrix of the input signal. Actually, the covariance matrix R is unavailable in real applications and is typically replaced by a finite sample covariance matrix: where N is the number of snapshots.

Improved Algorithm Based on the Invariant Property of
Noise Subspace. After introducing a virtual source, the noise subspace remains the same as before only when the DOA of the virtual source overlaps with the actual sources. is property, called the invariant property of noise subspace, is demonstrated in [27] and can be used to estimate the DOAs of actual sources. Diagonal loading is a widespread approach to improve robustness against mismatch errors [29,30]. We find that diagonal loading is also effective in improving the resolution of the IPNS method. us, we propose an improved algorithm based on the IPNS. Diagonal loading effectively reduces the power differences between the strong and weak sources, thereby increasing the resolution of the improved algorithm. Next, we certify that the IPNS still holds after diagonal loading so we can use it for DOA estimation.

2
International Journal of Antennas and Propagation e eigenvalue decomposition (EVD) of the covariance matrix R is given by the following: where U � U S U N and Σ � diag(Σ S , Σ N ); diag(·) denotes the diagonal matrix; and U S , U N , Σ S , and Σ N are the eigenvectors and eigenvalues of the signal subspace and noise subspace, respectively. e space spanned by the impinging sources and the noise are orthogonal; that is, are the eigenvalues of the signal subspace. Σ N is also a di- are the eigenvalues of the noise subspace. Suppose that a virtual source impinges on the receiving antenna array from the direction θ ′ . e virtual source does not physically exist but is mathematically needed. en, the new covariance matrix including the virtual source can be defined as follows: where k > 0 is the power of the virtual source. e EVD of the new covariance matrix R ′ can be given by the following: where U ′ � U S ′ U N ′ and Σ ′ � diag(Σ S ′ , Σ N ′ ), U S ′ and U N ′ are the eigenvectors of the signal subspace and the noise subspace, respectively, and are the signal and noise eigenvalues corresponding to the new covariance matrix R ′ , respectively. Suppose that the direction of the virtual source θ ′ overlaps with the ith actual source, that is, θ ′ � θ i , i ∈ [1, 2, 3, . . . , P]. Define a selective unit vector e i as follows: where e i is a unit vector consisting of P elements and all the elements are 0 except the ith element, which is 1. en, the new covariance matrix R ′ can be written as follows: where It is evident in equation (10) that the virtual source does not affect the noise power but only increases the power of the ith source when the direction θ ′ overlaps with that of the ith actual source. us, we can obtain that the noise subspace eigenvalues Σ N ′ remain consistent with Σ N when θ ′ overlaps with the ith actual source, that is, When the direction of the virtual source is different from those of any of the actual sources, that is, In this case, the noise subspace eigenvalues Σ N ′ are different from Σ N . e eigenvalues of R ′ in descending order are as follows: e eigenvalues of the noise subspace in equation (13) Diagonal loading has been demonstrated to be an effective approach for improving robustness in beamforming. We employ diagonal loading to improve the resolution of the IPNS method. After employing diagonal loading, the covariance matrix can be expressed as follows: where σ 2 is the coefficient of the diagonal loading. After the EVD of R, we can obtain the following: where Σ � diag(Σ S , Σ N ). It is apparent in equation (15) that the diagonal matrix Σ N of the noise subspace can be written as follows: After introducing the virtual source, we define the covariance matrix R ′ based on R as follows: International Journal of Antennas and Propagation Similarly, it can also be certified that the virtual source does not affect the noise subspace of R ′ compared to that of R if the direction θ ′ overlaps with the actual sources. Assume that the virtual source impinges on the receiving antenna array from the direction θ ′ . en, R ′ can be rewritten as follows: where R ′ � AR S ′ A H + σ 2 n I. e EVD of R ′ is given by the following: where When the direction θ ′ overlaps with the actual sources, that is, θ ′ � θ i , i ∈ [1, 2, 3, . . . , P], we can obtain from equation (16) that In contrast, as the direction of the virtual source is different from those of any of the actual sources . Now, we have demonstrated that the invariant property of the noise subspace still holds after employing diagonal loading.
e noise subspace eigenvalues of R ′ remain consistent with those of R when θ ′ � θ i , i ∈ [1, 2, 3, . . . , P]. Hence, we can still use the IPNS for DOA estimation after diagonal loading.
According to the above discussion, we can use the IPNS after diagonal loading to estimate the DOAs of the incident sources. e spatial zero spectrum f(θ ′ ) can be expressed as follows: e DOAs of the impinging sources are those θ for which 1/f(θ ′ ) reaches their maxima, which can be expressed as follows: where θ s and θ e denote the start angle and end angle of the spectral search, respectively. To implement the improved algorithm, the detailed steps are listed in Algorithm 1. Generally, the covariance matrix R is not available in practical applications, and the finite sample covariance matrix R is used instead of R.

Complexity Analysis of the Improved Algorithm.
e computational complexity of the algorithm has a direct relationship with its time consumption. Let J denote the number of spectral points. We must compute ‖a H (θ)U N ‖ 2 for each spectral point in the MUSIC method, where ‖·‖ 2 is the Frobenius norm. e dimensionality of U N is M-P; then, the computational complexity of the MUSIC method in the spectral search step is J(M − P)(M + 1) flops [31,32]. To obtain the noise subspace eigenvector U N , the EVD of the covariance matrix R is needed. If the FSD method in [33] is employed for the EVD of R, the complexity is given by M 2 (P + 2) flops. us, the total computational complexity of the classical MUSIC method is given by the following: However, the EVD of R ′ is required in each spectral point for the IPNS method.
us, the computational complexity of the IPNS method in the spectral search step is JM 2 (P + 2) flops. Before the spectral search step, JM 2 flops are needed for computing ka(θ ′ )a(θ ′ ) H to obtain R ′ . e total computational complexity of the IPNS method is as follows: For the proposed algorithm, diagonal loading requires another M flop. e total computational complexity of the improved algorithm is as follows:

CRB of a Weak Source in the Presence of a
Strong Source e CRB provides the lower bound for the unbiased estimators and provides fundamental lower limits for the variance of the parameter estimator [34][35][36]. us, studies on the CRB of weak sources in the presence of strong sources with different powers can be used as references for unbiased estimators. For multiparameter estimation, the CRB inequality has the following form [37]: where θ i is the DOA estimation of the ith source, J is the Fisher information matrix, and (·) −1 denotes the matrix inverse operation. e ijth element of the Fisher information matrix J is given by the following [34]: 4 International Journal of Antennas and Propagation For independent incident sources, the input signal covariance matrix is diagonal, that is, R s � diag(ρ 1 , ρ 2 , . . . , ρ P ), where ρ i represents the power of the ith source. Hence, the covariance matrix R can be simplified as follows: e partial derivative of the array covariance matrix R with respect to θ i is as follows: where the partial derivative of the steering vector a(θ i ) with respect to θ i is and q(θ i ) � e −jk 0 dsin θ i . dq(θ i ) m− 1 /dθ i is given by the following: (33) When the distance between the receiving antenna elements is λ 0 /2, we can obtain that dq(θ i ) m−1 /dθ i � −j(m − 1)π cos θ i e −j(m− 1)π sin θ i , m � 1, 2, . . . , M. Substituting equation (33) into equation (32) and then into equation (31), we obtain zR/zθ i . en, we can obtain the ijth element of the Fisher matrix J ij from equation (29). e CRB of θ � [θ 1 , θ 2 , . . . , θ P ] can be obtained, that is,

Numerical Simulation
Computer simulations regarding the CRB of the weak source, the spatial spectrum of the proposed method, the resolution probability, the root-mean-square error (RMSE), and the computational complexity are conducted to demonstrate the performance of our improved algorithm. We run 1000 independent Monte Carlo simulations for each experiment to study the probability of resolution and the RMSE of the proposed method, while the other two algorithms, the MUSIC algorithm and IPNS method, are also discussed as comparative simulations. e input signal-tonoise ratio (SNR) is defined as SNR � 10 log 10 (ρ 2 s /σ 2 n ), where ρ 2 s is the power of the source. In the following simulations, we set source 1 with SNR 1 as the weak source and source 2 with SNR 2 as the strong source. e interference-to-signal ratio (ISR) is defined as SNR 2 -SNR 1 . e ISR variation is the result of increasing SNR 2 , while SNR 1 does not change. Figure 1, suppose two independent sources impinge on the receiving antenna array with M � 9 element numbers, and the power of the noise is σ 2 n � 1. e study of the CRB of the weak source calculated by ���� θ CRB is conducted in the following simulations. e CRB increase percentage is defined as

CRB Analysis. As shown in
is the CRB of the source with ISR � 0 dB, which means that the incident two sources have the same power.
In the first simulation, we are interested in the change in the CRB of the weak source when the ISR increases, but SNR 1 remains the same. e numerical simulation results of the CRB as a function of the ISR under different SNR 1 conditions are shown in Figure 2. Figure 2 shows that the CRB of the weak source rises as the ISR increases and eventually reaches a fixed value, and this indicates that the CRB of the weak source deteriorates as the power of the strong interference increases. We can also see that CRB deteriorations exist for each SNR 1 , meaning that increasing the power of the strong source has a negative effect on the weak source, especially when the SNR 1 of the weak source is low. e increase percentages of the CRB of the weak source against the ISR under different SNR 1 are shown in Figure 3. From Figure 3, we can see that the CRB increases percentage (i) Step 1. Obtain the output of the uniform linear antenna array x(t) from equation (1) and compute the sample covariance matrix R as in equation (3). (ii) Step 2. Perform the EVD of R to obtain the noise subspace eigenvalue matrix Σ Ν � diag(λ P+1 , λ P+2 , . . . , λ M ). (iii) Step 3. Apply diagonal loading on R to obtain R � R + σ 2 I through equation (14). σ 2 is equal to λ 1 , which is the largest eigenvalue of the signal subspace. e noise subspace eigenvalues of R are Σ Ν � diag(λ P+1 , λ P+2 , . . . , λ M ). (iv) Step 4. Suppose that the virtual source impinges on the antenna array from direction θ′; then, obtain the matrix R′ � R + ka(θ′)a(θ′) H as in equation (18), where k � tr(R)/M, tr(·) denotes the trace of R and M is the number of array elements. (v) Step 5. Perform the EVD of R′ to obtain the eigenvalues in descending order as λ 1 ′ ≥ · · · ≥ λ P ′ ≥ λ P+1 ′ ≥ · · · ≥ λ M ′ . (vi) Step 6. For each θ′, compute the spatial zero spectrum f(θ′) � M j�P+1 (λ j ′ − λ j ) through equation (23). e DOAs of the sources are P sets of θ′ that enable 1/f(θ′) to reach P maximum values, as in equation (24). International Journal of Antennas and Propagation for the weak source increases sharply from 0 to 11% when SNR 1 � 0 dB, but it only increases from 0 to 5.5% when SNR 1 � 8 dB. is demonstrates that the CRB of the weak source is more susceptible to the power of the strong source when its SNR is low. Furthermore, it can be seen from Figure 3 that the CRB increase percentage increases from 11% to 12% as the ISR increases from 10 dB to 20 dB when SNR 1 � 0 dB. is means that the influence of the power of the strong source is limited when its power exceeds that of the weak source by too much.
In the second simulation, we further study the influence caused by the power increase of the strong source on the CRB of the weak source. e CRBs of the weak source as a function of the SNR 1 under different ISRs are shown in Figure 4. It can be seen in Figure 4 that the CRB of the weak source increases as the strong source power increases at each SNR 1 , and this once again demonstrates that the power of the strong source leads to CRB degradation for the weak source. Moreover, this reveals that the estimation performances of weak sources are affected by the power of strong sources.
e CRB increase percentage for the weak source against SNR 1 under different ISRs is shown in Figure 5. It is clearly seen in Figure 5 that the CRB increase percentage changes from 18.5% to 4.4% as SNR 1 increases from -10 dB to 10 dB when ISR � 20 dB. e CRB degradation is most apparent when SNR 1 is low, and this can also be found when ISR � 2 dB, ISR � 5 dB, and ISR � 10 dB; therefore, the weak source is more susceptible when its SNR is lower and the power of the strong interference is higher.
In the third simulation, the CRB of the weak source against the angular separation between the two sources is studied. We suppose that two independent sources impinge on the antenna array from ±θ; then, the angular separation between the two sources is 2θ. When θ varies from 1°to 7.5°, the CRB of the weak source versus the angular separation is shown in Figure 6. Figure 6 indicates that the CRB of source 1 (the weak source) is lowest when the two sources have the same power. e CRB of the weak source rises when the strong interference power increases under each angular separation.
From the above discussion, we can see that the CRB of the weak source is susceptible and tends to be influenced by the power of the strong source, especially when its SNR is low; this has not been studied before. However, sources with large power differences are very common in complex electromagnetic environments. e DOA estimation for sources with different powers should receive special attention.

Spatial Spectra.
In this simulation, we compare the spatial spectra of different algorithms when resolving two and three sources. Suppose that two independent sources impinge on the receiving antenna array (M � 9) from 2°a nd −2°, respectively. Source 1 is the weak source with SNR 1 � 2 dB, and source 2 is the strong source with SNR 2 � 12 dB. e spatial spectra are shown in Figure 7(a).
Suppose that three independent sources impinge on the receiving antenna array (M � 9) from −2°, 2°and 6°with SNR 1 � 17 dB, SNR 2 � 7 dB, and SNR 3 � 12 dB, respectively. e spatial spectra of the different algorithms are shown in Figure 7(b). From Figure 7, we can see that only the strong source is readily estimated by the MUSIC algorithm when two sources exist. ree sources cannot even be resolved using the MUSIC algorithm. e spectra of the improved algorithm and IPNS method are similar. As shown in Figure 7(b), three peaks are more distinct for the improved algorithm than those of the IPNS method.     International Journal of Antennas and Propagation

Probability of Resolution.
In this section, we focus on the performance of the improved algorithm in terms of resolving two impinging sources. For the first simulation, we suppose that two independent sources (source 1 and source 2) impinge on the receiving antenna array (M � 9) from 1.5°a nd −1.5°, respectively. e probability of resolution as a function of SNR 1 is shown in Figure 8. It can be seen in Figure 8 that the improved algorithm has a higher probability of resolution than those of the other two methods. When SNR 1 � 1 dB and ISR � 10 dB, the probability of resolution for the improved algorithm exceeds those of the IPNS method and the MUSIC algorithm by 17.8% and 62.0%, respectively. Furthermore, it is 22.4% and 55.1% higher than those of the IPNS method and the MUSIC algorithm when ISR � 0 dB and SNR 1 � 5 dB.
For the second simulation, we suppose that two independent sources impinge on the receiving antenna array (M � 9) with SNR 1 � 10 dB. e probability of resolution as a function of the angular separation between the input sources under different ISRs is shown in Figure 9, which indicates that the improved algorithm has better resolution performance than the other methods. When the angular separation is 1.4°and ISR � 10 dB, the probability of resolution for the improved algorithm exceeds those of the IPNS method of the MUSIC algorithm by 18.4% and 42.3%, respectively. Moreover, it is 25.4% and 56.7% higher than those of the IPNS method and the MUSIC algorithm when the angular separation is 2.15°and ISR � 0 dB, respectively. e probability of resolution for each method is low when SNR 1 is low or the angular separation is small, and the proposed method is no exception. en, the increase in the probability of resolution yielded by the proposed method is not so apparent when compared with the other two methods. Similarly, the increase in the probability of resolution for the proposed method is not so significant when SNR 1 is high or the angular separation is large. e probability of resolution for each method is high under such conditions.
From the abovementioned discussion, we can see that the improved algorithm has better resolution performance for sources with larger power differences; even for sources with the same power, the proposed algorithm still outperforms the other two methods. From the abovementioned discussion, we can see that diagonal loading is effective in improving the probability of resolution of the DOA estimation method based on the IPNS. Diagonal loading reduces the ratio of λ P+1 /λ M , thus increasing the probability of resolution. However, the increase in the probability of      International Journal of Antennas and Propagation 7 resolution generated by diagonal loading is limited. Diagonal loading also decreases the signal-to-noise ratio (SNR).

Root-Mean-Square Error.
In this section, we study the RMSE of the proposed method. Suppose that two independent sources impinge on the receiving antenna array (M � 9) from 1.5°and −1.5°, respectively. e RMSE of the weak source as a function of SNR 1 is shown in Figure 10. As a comparison, the CRB is also provided in Figure 10. We can find in Figures 10(a)-10(c) that the improved algorithm achieves better RMSE performance for the weak source than those of the other methods under each ISR. More specifically, the RMSE of the improved algorithm is similar but slightly better than that of the IPNS method when ISR � 5 dB and ISR � 10 dB, and this is clearest when ISR � 0 dB. e RMSE of the strong source and the CRB as a function of SNR 1 are shown in Figure 11. We can see in     Figure 11 that the RMSE of the improved algorithm is nearly the same as those of the other two methods when ISR � 5 dB. e RMSE of the improved algorithm is not as good as those of the other two methods when ISR � 10 dB, which means that the MUSIC algorithm has better estimation accuracy for strong sources. However, when the power of the strong source decreases to ISR � 0 dB, the RMSE of the improved algorithm is the lowest. However, the challenge in DOA estimation with respect to sources that have large power differences is the estimation performance of the weak source. Hence, it is objective to conclude that the RMSE of the proposed algorithm is similar to that of the IPNS method, which is even better than that of the MUSIC algorithm.

Computational Complexity.
In this section, we study the computational complexity of the proposed method. If a rough grid of 3°and a fine grid of 1°are employed, the spectral points 60(J 1 � 180°/3°) and 180(J 2 � 180°/1°) should be computed in the spectral search, respectively. We compare the computational complexities of different algorithms and display the results in Table 1. It is evident in Table 1 that the IPNS method has more computational complexity than the MUSIC algorithm; moreover, the computational complexity of the IPNS method grows faster than that of the MUSIC algorithm as M and P increase. However, the improved algorithm that exploits diagonal loading has nearly the same computational complexity as that of the IPNS method.

Conclusions
In this paper, we propose an improved algorithm to increase the resolution performance with regard to resolving sources with large power differences and further study the Cramer-Rao bound of the weak source in the presence of the strong source. We demonstrate that the IPNS still holds after diagonal loading. Simulation results show that the proposed algorithm has a better probability of resolution than the IPNS and MUSIC algorithms. Moreover, the proposed method does not result in increased computational complexity and maintains the estimation accuracy. We also find that the Cramer-Rao bound of the weak source suffers from degradation as the power of the strong source increases, especially when the SNR of the weak source is low. Our work provides a new way to estimate sources with different powers, and this is valuable in noncooperative situations, such as passive radar applications.
Data Availability e data that support the findings of this study are available from the corresponding author upon reasonable request.

Conflicts of Interest
e authors declare that there are no conflicts of interest regarding the publication of this paper.