Two-Dimensional DOA Estimation for Coprime Planar Arrays Based on Self-Correlation Tensor

In the coprime planar array (CPA), the existing tensor DOA estimation has the problem that the statistics are not fully utilized.We propose a two-dimensional DOA estimation method based on tensor self-correlation, which realizes the high-resolution and high-precision joint estimation of elevation angle and azimuth angle. Firstly, we represent the received signals of two subarrays with tensors and then obtain the self-correlation covariance tensor of the subarrays themselves and the cross-correlation covariance tensor of the two subarrays. en, we extract the covariance tensor corresponding to the maximum continuous virtual array and prove the expression of the maximum continuous virtual array aperture of the proposed method. Compared with the existing methods, the proposed method eectively improves the maximum aperture of the continuous virtual array. Finally, the signal subspace is solved by tensor expansion and tensor decomposition. Simulation results show that under the same conditions, the proposed method has higher estimation accuracy and degree of freedom than the cross-correlation tensor method, and the resolution is also improved signicantly.


Introduction
Direction of arrival (DOA) is a technology that processes and analyzes the data received from sensor array to obtain target location information, which is an important research topic in array signal processing [1]. It has been widely applied in radar, sonar, wireless communication [2][3][4][5][6], and other areas. Most subspace DOA methods were initially proposed based on non-sparse uniform linear arrays [7][8][9][10], which can complete DOA estimation in the case of overdetermination, that is, M array sensors can complete M − 1 signal estimation at most [11][12][13]. However, in the actual situation, the number of located target sources is greater than the number of array sensors, that is, underdetermined estimate.
To improve the degree of freedom of the array and achieve a practical estimation of multiple signal sources under unknown conditions, experts and scholars have proposed sparse array [14][15][16][17], which can e ectively expand the aperture of the array and improve the degrees of freedom and estimation accuracy, and typical sparse arrays include nested arrays [18,19], coprime arrays [20], and minimum redundant arrays (MRAs) [21]. Also, many algorithms are proposed based on these sparse arrays, such as spatial smoothing MUSIC (SSM) [22], compressive sensing (CS) [23], and discrete Fourier transform (DFT) [24]. However, the above arrays are all one-dimensional arrays. Usually, one-dimensional arrays can only estimate one-dimensional DOA. To achieve two-dimensional DOA estimation, many experts and scholars have generalized traditional subspace methods based on one-dimensional arrays to two-dimensional arrays. Aiming to realize DOA estimation with high degrees of freedom as much as possible, two-dimensional MUSIC methods (2-D MUSIC), two-dimensional ESPRIT methods (2-D ESPRIT) [25], and PM methods [26] based on coprime planar arrays and coprime L-arrays have also been proposed, which combine the advantages of sparse formations and realize two-dimensional DOA estimation. However, these methods are based on matrix operation processing, which will blur the data structure and bring loss to the accuracy of the DOA estimation. is operation will blur the data structure and bring loss to the accuracy of the DOA estimation.
In this regard, experts and scholars have introduced tensors into DOA estimation [27][28][29][30][31][32]. Tensors [33][34][35] are high-dimensional data, which can save information of different physical meanings in other dimensions. By comparing the performance of tensor-based DOA estimation method with matrix-based DOA estimation method, the authors in [36][37][38] showed that the former has better estimation accuracy. However, these methods are mainly based on the cross-correlation tensor (CCT) for processing the subarray received signal, which fails to make full use of the position information of the array elements, which will lead to the narrowing of the aperture of the virtual array, resulting in loss of degrees of freedom and estimation accuracy.
In this paper, we propose a self-correlation tensor-based method for estimating 2-D DOA in the coprime planar array. Compared with the cross-correlation tensor-based method, the proposed method can achieve a higher degree of freedom and higher accuracy of joint estimation of elevation angle and azimuth angle. e proposed method uses tensors to represent the received signal. rough the self-correlation processing of the signal, we can obtain the covariance tensor of the signal and then extract the covariance tensor to get the covariance tensor of the corresponding maximum continuous virtual uniform rectangular array. Finally, we can estimate the azimuth angle and elevation angle of the signal by CP decomposition (CANDECOMP/PARAFAC decomposition (CPD)) [33,39]. Compared with CCT and SSM, simulation results show that the proposed method improves estimation accuracy and degree of freedom. e main contribution of this paper is summarized as follows: (1) We propose a new method for estimating the twodimensional (2-D) DOA in the coprime planar array based on SCT, which achieves the high-resolution and high-precision joint estimation of elevation angle and azimuth angle. It uses tensors to model signals, taking into account the dimensional structure of the data. Compared with the CCT and SSM methods, the estimation accuracy of the proposed method is better than that of CCT and SSM. Compared with the PSS methods, when the number of sources is less than the number of minimum subarray elements, the performance of the proposed method is close to that of the PSS method, and when the number of sources is more than the number of small subarray elements, the proposed method is much better than the PSS method. (2) Aiming at the problem that the virtual array is not constructed sufficiently in the existing tensor DOA estimation method, the self-covariance information and the cross-covariance information are combined to increase the maximum degree of freedom and estimation accuracy of the DOA estimate. However, in order to apply this method to coprime planar arrays of different size, the maximum continuous virtual planar array aperture expression must be obtained, so we derive the expression of the largest continuous virtual planar array in the appendix, which makes this method more practical. By proof derivation, the aperture expression of the maximum e rest of the paper is organized as follows. We give the signal model in Section 2. In Section 3, we introduce the twodimensional DOA estimation based on the proposed method and analyze the performance of the method. e simulation results are provided in Section 4, and we conclude this paper in Section 5.
Notations. Scalars are denoted by italic letters (M, m), matrices and column vectors by bold letters (C, c), tensors by (R, N), and collections by hollow letters (N, W). Let (•) T , (•) * , (•) H , and E • { } correspond to transpose, conjugate, Hermitian transpose, and statistical expectation. e Kronecker product of two matrices A, B is denoted by A ⊗ B, and the outer product of two tensors A, B is denoted by A°B.

Signal Model.
e model of the coprime planar array shown in Figure 1 comprises two sparse uniform planar arrays, named as D 1 and D 2 , and the spacing between adjacent elements is Nλ/2, Mλ/2, where N and M are coprime integers and λ is the wavelength of impinging signals. D 1 and D 2 only coincide at the origin, the total number of physical elements is M 2 + N 2 − 1, and the reference position of the ith element to the origin is (x i , y i ). e mathematical expression of the position of the whole array element is Ω � (mNd 0 , mNd 0 ) ∪ (nMd 0 , nMd 0 )|m ∈ [0, M − 1], n ∈ [0, N − 1]}. Suppose there are k far-field narrowband independent signals impinging on the array with the signals from where θ i ∈ (− π, π), φ i ∈ (0, π/2) are the azimuth angle and elevation angle of the ith signal.
According to the geometric relation of the element and the tensor theory [27], the received signals of D 1 can be represented as a 3-D tensor: where are steering vectors of D 1 and u k � cos θ k cos φ k , v k � sin θ k cos φ k are directional factors in the x and y directions. e signal data vector where L denotes number of snapshots and N 1 is additive white Gaussian noise of D 1 . Similarly, the expression of D 2 received signal can be represent as where are steering vectors of D 2 and N 2 is additive white Gaussian noise. en, we can slice the tensors X i along the snap dimension [40] and get X i (1), X i (2), . . . , X i (l), . . . , X i (L), i � 1, 2. e 4-D covariance tensor R X i X j can be expressed as However, the number of snapshots is limited. Here, we use statistical average instead of time average, and 4-D tensor of the received signal covariance can be represented as and obtain P X i X j .

Tensor DOA Estimation
In this section, we propose a new tensor-based two-dimensional DOA estimation method that effectively incorporates auto-related information of signals. Compared with the methods in [38], the proposed method in this paper effectively increases the aperture of the virtual array and improves the estimation accuracy.

Construction of Virtual Array.
We define the dimensional sets Q 1 � 1, 3 { }, Q 2 � 2, 4 { } and then use tensor expansion [33,41] to merge dimensions with the same directional factor: Figure 2 shows the virtual array position of a coprime planar array with M � 3 and N � 2. Consider that coarray tensor U i is the single snapshot, which cannot meet the uniqueness principle of tensor decomposition, and the number of sources that can be estimated is limited. To solve this problem, we generate new structural coarray tensors H i by structural recombination of coarray tensors U i . However, only continuous virtual planar arrays can be used for structural recombination of coarray tensors U i . us, we divide the continuous sensors in the virtual array into four parts, which are represented as In Figure 2, we show four continuous virtual area arrays wrapped in rectangular boxes, respectively. e elements in U i,j are reorganized to map the augmented virtual URA U 1 , U 2 , U 3 , and U 4 , where the equivalent 2-D coarray tensor U 1 , U 2 , U 3 , U 4 can be represented as

Tensor Decomposition of Structure Coarray.
Coarray tensor U i is the single snapshot, which cannot meet the uniqueness principle of tensor decomposition, the identification of tensor decomposition is limited, and only a single Mathematical Problems in Engineering source can be detected. In order to estimate more sources, the idea of spatial smoothing algorithm can be introduced. We e position of the first sensors of each subarray can be represented as , so the coarray tensor of the first subarray can be represented as where E � [1, e jπ(Mu k − Nv k ) , e jπ(Nu k − Mv k ) , e jπ(Nu k − Nv k ) ] is used to compensate for the phase difference caused by the different positions of G i x,y in the four virtual URA, are the steering vectors of H i 1 . erefore, by space smoothing for virtual URA u i , smooth coarray tensor H i can be represented as where d a ix (u k ) � e jaπu k , a ∈ 1, 2, . . . , J i , and d b iy (u k ) � e − jbπv k , b ∈ 1, 2, . . . , J i serve as the shifting factor vectors along the direction of x-axis and y-axis, respectively. We can concatenate H 1 , H 2 , H 3 , and H 4 along the third dimension as 3-D tensor T ∈ C P×P×4 , which is represented as and then we do a tensor decomposition of T; we define the rank(a) as Kruskal's rank of a. According to the tensor decomposition uniqueness criterion, we have where . In order to make full use of the DOFs of the structured coarray tensor T, we should maximize f(P, M, N). Take M � 2, N � 3 for example; when P � 4, degree of freedom (DOF) can reach the maximum 15, and the signal subspace U s can be represented as en, the noise subspace can be represented as where I is the identity matrix. With the steering function of the virtual N for the scanning bin (u, v) defined as , the coarray tensor MUSIC spectrum of MUSIC spectrum based on SCT can be calculated as As such, the 2-D DOA estimates (θ k , φ k ) can be obtained via spectrum searching on P(u, v). e detailed steps of the method are shown in Table 1.

Complexity Analysis.
We analyze the computational complexity of the proposed method in the coprime planar array and compare it with the CCT method, SSM method, and PSS method. Let us first analyze the computational complexity of the proposed method. For a coprime planar array, it can be divided into two uniform planar subarrays of sizes M × M and N × N, the proposed method constructs covariance tensors by received signals, and two subarrays can construct four covariance tensors, so the computational complexity of this step is O(L(M 4 + N 4 + 2M 2 N 2 )). en, by tensoring the covariance tensor, the structure covariance matrix T is obtained, and the computational complexity of this step is O(2KP 2 (2M + N − P) 2 ). en, we do a tensor decomposition of T, the computational complexity of this step is (16KP 2 (2M + N − P) 2 + 3K 3 ), and the complexity required for the spectral peak search to calculate the elevation angle and azimuth angle is O(n g n l P 4 ), where n g � 360°/0.01°, n l � 90°/0.01°, and n g × n l represents the number of meshes searched by the spectral peak. erefore, the total complexity of the SCT method is O(L(M 4 + N 4 + 2M 2 N 2 )+ 2KP 2 (2M + N − P) 2 + 16KP 2 (2M + N − P) 2 +3K 3 + n g n l P 4 ).
e complexity of the PSS method is O(n g n l (N 4 /M 2 + M 4 /N 2 )).
For a more intuitive representation, Table 2 summarizes the complexity of the different methods. As can be seen from Table 2, the computational complexity of the four methods is mainly derived from the spectral peak search, and the   Figure 3, where the number of snapshots is 5000, SNR � 15 dB, and the search step is 0.1°.

Simulation Results
It can be seen that the DOA estimated by the proposed method matches the true ones. erefore, the proposed method also has good estimation performance under underdetermined conditions.

Comparison of Angle Resolution Capabilities.
We compare the proposed method with the method based on CCT and SSM, respectively, and the number of sources is 2. First, we compare the case where the direction of the signals is close to each other, and the azimuth and elevation angles of the two sources are (20°, 20°) and (21°, 21°), respectively. To avoid the loss of accuracy of the algorithm estimate due to the small aperture of the array, we set M � 7, N � 8, the signal-to-noise ratio (SNR) is 10 dB, the number of snapshots is 2000, and the search step is 0.01°. It can be seen from Figure 4 that the proposed method in this paper can distinguish the signal direction very well, while the other two methods cannot distinguish two sources. Compared with SSM, the proposed method has the same virtual array aperture, but the proposed method uses tensors to model the signal and estimate the DOA from multiple dimensions. Compared with the CCT, both the proposed method and CCT are based on tensor modeling of signals, but the proposed method has a larger virtual array aperture in DOA estimation, thereby improving the resolution and accuracy, so the proposed method has higher estimation accuracy and angle resolution.

Comparison of RMSE Performance under the Condition of Different Angle Intervals.
We compare the RMSE performance of the proposed method with the other three methods under the condition of different angle intervals. Consider M � 3, N � 2, the search step is 0.01°, the number of snapshots � 2000, and SNR � 10 dB. Firstly, we set the DOAs of the two sources as (20°, 20°), (23°, 23°), and each time, the angle interval change is 1°, so the next DOAs of the sources are (20°, 20°), (24°, 24°). By analogy, the last DOAs are (20°, 20°), (35°, 35°), respectively. e simulation results are shown in Figure 5. When the angle interval is close, the advantage of high tensor dimension makes the RMSE performance of CCT method better than that of SSM, and the proposed method not only has the advantage of high tensor dimension but also has a larger array aperture, making the RMSE performance better than CCT and SSM. e PSS method exchanges a smaller degree of freedom for a larger array aperture, and the RMSE of proposed method is close to that of the PSS method when the angle interval is not less than 5°.

Comparison of RMSE Performance under Two Sources.
We compare the estimation accuracy of the proposed method with the subspace-based SSM, PSS [42], and CCT and define the root mean square error as the performance metric, and the RMSE of azimuth can be calculated as where R stands for the number of Monte Carlo trails, R � 100, and φ r,k is the estimated azimuth angle of the kth source for the rth trial, respectively. Consider M � 3, N � 2, the search step is 0.01°, and the azimuth and elevation angles of the two sources are (20°, 20°) and (30°, 30°), respectively. Figure 6 shows the RMSE performance varying with the Table 1: e steps of the proposed method.
Step 1 Construct 3-D tensor of the received signal X i Step 2 Construct 4-D self-correlation coarray tensor R X i X j from the received signal X i of the subarray Step 3 Transform the dimension order of the 4-D tensor R X i X j and then merge the dimensions with the same directional factor to get U i,j Step 4 Take the coarray tensor U 1 , U 2 , U 3 , U 4 corresponding to the largest continuous virtual array from U i,j Step 5 Transform the dimension of the coarray tensor U i by tensor expansion to obtain structured coarray tensor T Step 6 Perform CPD on the structured coarray tensor T to obtain matrix factors c x (u k ), c y (v k ) containing direction information and calculate (θ k , φ k ) by spectral peak search SCT    Figure 7 shows the RMSE performance varying with SNR, where the simulation conditions of SNR are -10 dB, -5 dB, 0 dB, 5 dB, 10 dB, and 15 dB, and the number of snapshots � 2000. As shown in Figures 6 and 7, the proposed method has better estimation performance than CCT and SSM. Under the condition of low SNR and low number of snapshots, the proposed method has more advantages and better robustness. PSS is slightly better than the proposed method in precision. e specific reasons are as follows: (1) Compared with CCT, the proposed method increases the self-correlation of the subarrays, which can increase the aperture of the maximum continuous virtual array, the degrees of freedom, and the estimation accuracy. (2) Compared with SSM, the proposed method effectively incorporates the multidimensional information, which is not utilized in the matrix-based methods. (3) Under the condition of the same physical sensor, although PSS has higher estimation accuracy and its estimation accuracy is related to the array aperture of the large subarray, its degree of freedom is limited by the number of array elements of the small subarray, that is, DOF ≤ min(M 2 , N 2 ) − 1. e DOA estimation method of PSS is to estimate the signal direction with two sparse uniform subarrays, and the interval between the elements of the subarrays is greater than half a wavelength, so ambiguous peaks will appear to interfere with the DOA estimation, but the corresponding MUSIC spectrums of the two subarrays uniquely generate a common peak at the direction of true DOA, so the estimated accuracy of the PSS depends on the array aperture of the large subarray. When estimating two sources, the main reason why the PSS estimation accuracy is higher than the proposed algorithm when estimating the two sources is that the array aperture of the PSS in the DOA estimate is ((M − 1)Nλ) 2 . In the proposed method, according to (11), the virtual array aperture is ((P − 1)λ) 2 . In the RMSE performance simulation, P � 4, M � 3, N � 2, and the array aperture of the proposed method is smaller than that of the PSS method, but the theoretical degree of freedom of the PSS algorithm is only 3, which is much lower than the theoretical degrees of freedom of the proposed method.

Comparison of RMSE Performance under Four Sources.
We compare the RMSE performance of the proposed method with the other three methods under the condition of estimating four sources. Consider M � 3, N � 2, the search step is 0.01°, and the azimuth and elevation angles of the two sources are (15°, 15°), (30°, 30°), (45°, 45°), and (60°, 60°), respectively. Figure 8 shows the RMSE performance varying with the number of snapshots where the simulation conditions of the number of snapshots are 100, 500, 1000, 2000, and 5000 and the SNR � 15 dB. Figure 9 shows the RMSE performance varying with SNR, where the simulation conditions of SNR are -5 dB, 0 dB, 5 dB, 10 dB, and 15 dB and the number of snapshots � 5000. e PSS method reduces the degrees of freedom in exchange for a higher estimation accuracy, and the maximum number of predictable sources is DOF ≤ min(M 2 , N 2 ) − 1. As shown in Figures 8 and 9, we can see that the PSS method fails in the DOA estimation    when the number of sources is not less than the number of array elements of the small array. e proposed method has high degrees of freedom after adding the self-covariance information, and the estimation accuracy is better than the other three methods.

Conclusion
In this paper, we propose a two-dimensional DOA estimation method based on self-correlation tensor for the problem of statistics of signals underutilized in tensor DOA estimation and solve the problem that tensors are insufficient in virtual array aperture construction in DOA estimation. First, we use tensor to represent the received signal, then construct the self-correlation coarray tensor of the received signal, then merge the dimensions with the same direction factor and extract the coarray tensor corresponding to the largest continuous virtual array, obtain the structural coarray tensor with multirank by tensor expansion, and finally approach the CPD to the structural coarray tensor to obtain DOA estimation with super-resolution.
rough simulation and proof, the proposed method effectively expands the maximum virtual array aperture. e proposed method can effectively enlarge the aperture of the maximum continuous virtual planar array and increase the degree of freedom. Compared with CCT and SSM, simulation results show that under the same conditions, the proposed method has better performance. Where In conclusion, since Q � − Q, for any mutual prime number M, N, there is always k 1 , p 1 , such that en, we arbitrarily take two elements from the sets H, H, H 6 , and H 7 to form a plane C, C � C 1 ⋃ C 2 ⋃ C 3 ⋃ C 4 , where C 1 � x 1 , y 1 |x 1 ∈ H, y 1 ∈ H C 2 � x 2 , y 2 |x 2 ∈H, y 2 ∈H C 3 � x 3 , y 3 |x 3 ∈ H 6 , y 3  We define Ω � (x, y)|1 − M − N ≤ x, y ≤ M + N − 1 . e planar region consisting of C 1 , C 2 is called Ω 1 , and the planar region consisting of C 3 , C 4 is called Ω 2 . Since there are some breakpoints in H and H, the plane Ω 1 formed by the element in H and H also has breakpoints, and a collection of all its breakpoints can be represented as erefore, the largest continuous part in Ω is Data Availability e data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest
e authors declare that they have no conflicts of interest.