A Propagator Method for Bistatic Coprime EMVS-MIMO Radar

In this paper, a novel two-dimensional (2D) direction-of-departure (DOD) and 2D direction-of-arrival (DOA) estimate algorithm is proposed for bistatic multiple-input multiple-output (MIMO) radar system equipped with coprime electromagnetic vector sensors (EMVS) arrays. Firstly, we construct the propagator to obtain the signal subspace. (en, the ambiguous angles are estimated by using rotation invariant technique. Based on the characteristic of coprime array, the unambiguous angles estimation is achieved. Finally, all azimuth angles estimation is followed via vector cross product. Compared to the existing uniform linear array, coprimeMIMO radar is occupying large array aperture, and the proposed algorithm does not need to obtain signal subspace by eigendecomposition. In contrast to the state-of-the-art algorithms, the proposed algorithm shows better estimation performance and simpler computation performance. (e proposed algorithm’s effectiveness is proved by simulation results.


Introduction
In recent ten years, multiple-input multiple-output (MIMO) radar has aroused extensive attention. MIMO radar has many important branches, for example, detection, parameter estimation, waveform design, and synchronization. Among them, the estimation of DOD and DOA is one of the important research fields of bistatic MIMO radar. Up to now, many algorithms have been proposed to solve the above issue, for example, signal parameters estimation method based on rotational invariance technique (ESPRIT) [1], maximum-likelihood estimator [2], spectrum peak search method [3], tensor approach [4,5], sparsity-aware strategy [6], and propagator method (PM) [7]. However, most of the existing algorithms can only solve the one-dimensional angle estimation problem.
ere are just a small number of algorithms focusing on the 2D-DOD and 2D-DOA estimation in MIMO radar [8][9][10]. For the 2D estimation algorithms in [8][9][10], the key of them is to utilize nonlinear geometry of the scalar sensor array [8][9][10]. On the other hand, some works propose the use of electromagnetic vector sensor (EMVS) to solve the 2D-DOD and 2D-DOA estimation in bistatic MIMO radar. Different from the scalar sensor array, a signal EMVS can not only perform the task of two-dimensional direction angle estimation but also provide the polarization status of the incoming signal. e EMVS is firstly proposed to be applied to bistatic MIMO radar in [11], where an EMVS array is mounted at transmitting and an EMVS array is at receiving. In [12], a new method for 2D-DOD and 2D-DOA estimation utilizing transmitting and receiving EMVS arrays is proposed. e framework used the ESPRIT-based method to estimate the 2D-DOD and 2D-DOA. en, vector cross-product strategy is used to estimate the polarization status and azimuth angles. In [13], a method similar to PM is proposed to avoid the step of eigendecomposition in [12]. However, none of the above methods can make full use of the MIMO radar's virtual aperture when we apply the selection matrix to the array measurement. In addition, the mentioned methods in [12][13][14] require extensive pair calculation and offer limited identifiability. For avoiding existing shortcomings in above algorithms, a parallel factor (PARAFAC) scheme was proposed in [15]. However, low computational efficiency and sensitivity to initialization are the drawbacks of PARAFAC decomposition. In [16], an improved PM-like method was introduced, which can realize the rotation invariance of the whole virtual steering vector. In addition, it can automatically pair the elevation angles associated with the receiving and transmitting arrays. More recently, the ESPRIT-like approach was investigated in [17], which is suitable for distributing EMVS sensor.
In order to eliminate the effect of phase ambiguity, less than half wavelength is the limited requirement of distance between sensors, which limits the receive and transmit array's virtual aperture. In addition, the mutual coupling problem may occur when sensors are closely displaced [18], which may reduce the estimation performance. e coprime array has a promising prospect to settle the above problems, so it has attracted extensive attention [19][20][21]. Also, the EMVS with nested array geometry has also been considered in [22,23]. Two sparse and uniform linear arrays constitute a coprime linear array (ULA), in which the internal element spacings are coprime numbers. It is ambiguous to estimate the angle of two subarrays separately, but the coprime property can uniquely determine it. Due to the increase of array aperture, coprime array performs better than ULA in angle estimation. e results show that the coprime EMVS array can effectively improve the performance of parameter estimation [24]. However, there is little research on EMVS-MIMO radar's parameter estimation.
is paper presents a bistatic coprime EMVS-MIMO radar framework different from ULA manifold [12][13][14][15][16], since both receivers and transmitters are composed of coprime EMVS. On this basis, the method of using propagator to construct signal subspace is proposed. en the transmitting elevation and receiving elevation are estimated and paired one-to-one, but the estimation result is ambiguous. Fortunately, the coprime characteristic allows a unique elevation. Finally, the cross-product technique is used to estimate the azimuth angle, and the results are automatically paired to estimate the elevation angle. Based on coprime array, this method increases the aperture, so it has better estimation performance. Numerical examples show that the estimation effect is improved. roughout the paper, lowercase letters represent vectors and uppercase letters denote matrices, respectively. An M × M identity matrix is represented by I M , 0 M stands for the M × M full-zeros matrix, and 1 n represents a row vector whose n − th entity is one and others are zero; ⊗ represents the Kronecker product and ⊙ is the Khatri-Rao product. { } denotes a diagonal matrix whose diagonal entities consist of the n − th row of B; ‖ · ‖ F represents the Frobenius norm; angle(·) means to gain the phase. real(·) returns the real part in the bracket, and ⊕ is Hadamard product. e vector cross product between two column vectors e 1 � [e 1 , e 2 , e 3 ] T and e 2 � [e 4 , e 5 , e 6 ] T is defined as e 1 * e 2 � 0 −e 3 e 2 e 3 0 −e 1 −e 2 e 1 0 ⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦.

The Proposed Algorithm
2.1. Problem Formulation. As shown in Figure 1, suppose that there is a bistatic EMVS-MIMO radar system, whose receiving and transmitting are distributed in the coprime geometries. In terms of the transmit array, let us assume that Subarray 1 is equipped with M t EMVS and Subarray 2 is equipped with N t EMVS, where M t and N t are coprime integers. Subarray 1's adjacent distance is N t λ/2 and Subarray 2's adjacent distance is M t λ/2, where λ represents the carrier wavelength. Similar to transmitting array, two EMVS subarrays constitute the receiving array, which have M r and N r elements, respectively. What is more, it is assumed that there are K far-field targets appearing in same range bin, and k − th target's 2D-DOA pair is (θ r,k , ϕ r,k ), and 2D-DOD pair is (θ t,k , ϕ t,k ), where the elevation angles are denoted by θ t,k and θ r,k , and the azimuth angles are denoted by ϕ t,k and ϕ r,k . e matched outputs can be expressed as where τ represents the slow-time index (pulse index), and s(τ) � [s 1 (τ), s 2 (τ), . . . , s K (τ)] T represents the target reflect coefficient vector; a t,k � [e jN t (M t − 1)π sin θ t,k , . . . , e jN t π sin θ t,k , 1, e − M t π sin θ t,k , . . . , e − jM t (N t − 1)π sin θ t,k ] T is the k − th transmit steering vector, and a r,k � [e jN r (M r − 1)π sin θ r,k , . . . , e jN r π sin θ r,k , 1, e − M r π sin θ r,k , . . . , e − jM r (N r − 1)π sin θ r,k ] T is the k − th receive steering vector; b t,k ∈ C 6×K represents the k − th polarization response vectors of transmit array, and b r,k ∈ C 6×K represents the k − th polarization response vectors of the receive array; stand for the associate polarization response matrices, respectively. e Gaussian noise vector with zero mean and variance of σ 2 is represented by n(τ). Moreover, where c t,k and c r,k represent the k − th transmit and receive auxiliary polarization angles, respectively. e corresponding polarization phase differences are indicated, respectively, by η t,k and η r,k . Let b tp,k and b rq,k be the p-th and It should be emphasized that b t1,k and b r1,k represent the electric field information vectors, while b t2,k and b r2,k account for the magnetic field information vectors. In addition, When the targets are uncorrelated, the y(τ)'s covariance matrix is where . . , r K is the s(τ)'s covariance matrix, and r k represents the power of s k (τ). In practical application, the estimation of R can be obtained via L available snapshots; that is, Herein, estimating the 2D-DOD and 2D-DOA from R is our goal.

Propagator.
Let the first K rows of C be C 1 ∈ C K×K , and let the rest be C 2 ∈ C (36MN− K)×K ; that is, Since C 1 is a nonsingular matrix, taking a unique linear transformation of C 1 can get C 2 ; that is, where P c ∈ C (36MN− K)×K represents the propagator. Let Denote the noiseless item of R as R; that is, where G ∈ C 36MN×K and H ∈ C 36MN× (36MN− K) . Obviously, According to the relationship in (10), we have Consequently, P c can be estimated by fitting where G is the estimation of G and H is the estimation of H. According to equation (15), we can obtain the P c 's leastsquares solution by ereafter, we can obtain the P's estimation via Based on equation (11), there is a full-rank matrix T so that Obviously, P plays a similar role as the signal subspace that is obtained from eigendecomposition.
Let ψ t1 � (J N2 P) † J N1 P and ψ t2 � (J N4 P) † J N3 P. Obviously, Θ t1 and Θ t2 consist of the eigenvalues of ψ t1 and ψ t2 , respectively. Moreover, the T's k − th column is an eigenvector that corresponds to the k − th eigenvalue of both ψ t1 and ψ t2 . Next, we have the eigendecomposition on ψ t1 so that we get the eigenvector matrix T as well as the associate eigenvalue matrix Θ t1 � diag λ t1,1 , λ t1,2 , . . . , λ t1,K ; that is, Multiply ψ t2 left by T and right by T − 1 to get the Θ t2 's estimation, which is denoted by Θ t2 . Let the k − th diagonal entity of Θ t2 be λ t2,k , and then the transmit elevation angle can be obtained from Subarray 1 and Subarray 2 as Due to the fact that the mapping function angle(·) is wrapped within the range [−π, π], but N t π sin θ t,k and M t π sin θ t,k are within [−N t π, N t π] and [−M t π, M t π], the transmit elevation angles and receive elevation angles estimated from (24a24b) are ambiguous.
Similarly, we have

Mathematical Problems in Engineering
Define Obviously, Since the eigenvector matrix T has been estimated previously, we can get the estimation of Θ r1 and Θ r2 (denoted by Θr1 and Θr2, respectively) by left and right multiplying operation on ψ r1 and ψ r2 , respectively. Let λ r1,k and λ r2,k be the k-th diagonal elements of Θr1 and Θr2, respectively. e receive elevation angles can be estimated via Similar to the estimated transmit elevation angles, due to the coprime characteristic of the two subarrays, the receiver elevation angles estimated from (28a 28b) are ambiguous.

Unique Elevation Angle Determination.
It is worth noting that θ (est) t1,k and θ (est) r1,k are obtained via the rotational invariance properties of Subarray 1, while θ (est) t2,k and θ (est) r2,k are achieved via Subarray 2's rotational invariance characteristics. Because of the coprime between M t , N t and M r N r , θ (est) t1,k , θ (est) r1,k and θ (est) t2,k , θ (est) r2,k can be uniquely determined even though they are ambiguous.
Since the transmit Subarray 1's interelement spacing is N t λ/2, there are N t possible answers, which include one from (24a). e connection between the n t -th (n t � 1, 2, . . . , N t ) possible solution θ (n t ) t1,k and the estimation θ (est) t1,k can be expressed as For the transmit Subarray 2, the k-th estimated transmit elevation angle should have M t possible solutions, including one derived from (24b); the m t − th (m t � 1, 2, . . . , M t ) possible solution θ eoretically, the unique solution can be determined by finding the coincident solutions in (29) and (30). As displayed in Figure 2, the results of (29) and (30) are presented, where θ � 30°, M t � 3, and N t � 4 are considered. In practice, the recovered angles may be closely distributed rather than coincide. erefore, it can be estimated by the average of the closest solutions of sin θ (n t ) t1,k and sin θ where θ t1,k and θ t2,k denote the two nearest solutions' associate angles. Similar to the transmit, the unique receive elevation angle estimation θ r,k can be obtained by us.

Azimuth Angle Estimation.
To estimate the azimuth angle, we need to estimate B t and B r first. According to (17), we can estimate C via Define J 1 � 1 p ⊗ I 6N ∈ C 6M×36MN , p ∈ 1, 2, . . . , 6M { }, and we can find that We have that Once the elevation angles have been estimated, the direction matrices A t and A r can be reconstructed (denoted as A t and A r , respectively). en B r can be estimated via where C r � J 1 C � J 1 PT − 1 . ereafter, we can compute where b r1,k represents a vector that is made up of the first three entities of B r 's k − th column and b r2,k is a vector that consists of the last three entities of the k − th column of B r , respectively. Finally, we can estimate ϕ r,k by where f r,k (1) is the f r,k 's first entity and f r,k (2) represents the second entity of f r,k . Similarly, we can define J 2 � I 6M ⊗ 1 q ∈ C 6M×36MN , q ∈ 1, 2, . . . , 6N { }, and construct C t � J 2 P. According to (35)-(37), the receive azimuth angle ϕ t,k can be obtained. Since the azimuth angle is one-to-one mapping with the elevation angle and the perturbation has been compensated, therefore, the estimated angles can be matched automatically.
2.6. 2D-TPA and 2D-RPA Estimation. As long as we complete the 2D-DOD and 2D-DOA estimates, we can get the estimates of D tk and D rk , represented by D tk and D rk , respectively. Based on equations (3a) and (3b), we can get the estimations of v tk and v rk via Finally, 2D-TPA and 2D-RPA estimation is finished via Because the premise of estimating 2D-TPA and 2D-RPA is to get 2D-DOD and 2D-DOA estimations, angles θ t,k , ϕ t,k , c t,k , η t,k and θ r,k , ϕ r,k , c r,k , η r,k are automatically paired.
In order to understand the algorithm proposed in this paper more efficiently, the following are the algorithm steps of the estimator proposed in this paper: Step 1: obtain the estimation of the covariance R firstly; Step 2: calculate the propagator P c by (16) and construct matrix P; Step 3: build the selection matrices J t,n and J r,n (n � 1, 2, 3, 4), and get ψ t via (21a), (21b) and (22a), (22b); Step 4: obtain T and corresponding eigenvalues by having eigendecomposition on ψ t , and get θ est t1/2,k via (24a) and (24b); Step 5: obtain ψ r via (26a) and (26b), and get the ambiguous elevation angles in receiver through (28a) and (28b); Step 6: recover all the estimations based on (29) and (30), and determine the unique elevation from (30); Step 7: construct J 1 , J 2 , A r , and A t , and then calculate C t and C r . ereafter, recover B r and B t . Finally, achieve the azimuth angles via (36) and (37).

Complexity.
In terms of the proposed algorithm, computing P c results in major complexity. Since the proposed algorithm can pair angle estimation automatically, the final complexity is 72MNK 2 . Next, the complexity of PM [16], PARAFAC [15], and ESPRIT [12] algorithms is listed in Table 1 to facilitate comparison with the proposed algorithm. As it can be seen from Table 1, compared with the ESPRIT, PM, and PARAFAC algorithms, the proposed algorithm is more efficient.

Simulation Results
We use the Monte Carlo method to verify that the proposed method is effective in this section. In the simulation, we equip the bistatic coprime EMVS-MIMO radar with M transmit sensors and N receive sensors, respectively. e receive array and transmit array are both distributed in coprime geometries, respectively, and the subarrays' sensor numbers are M t , N t (M � M t + N t − 1) and M r , N r (N � M r + N r − 1), respectively. Assume that there are K � 3 and θ t � (40°, 20°, 30°), ϕ t � (15°, 25°, 35°), c t � (10°, 22°, 35°), η t � (36°, 48°, 56°), θ r � (24°, 38°, 16°), ϕ r � (21°, 32°, 55°), and c r � (42°, 33°, 60°) with η r � (17°, 27°, 39°), respectively. In the simulation, we set p and q to be 1 and 2, respectively. e reflected coefficients are modeled with a complex Gaussian process, and there are L snapshots. Each simulation curve is based on 200 independent trials. At the same time, the root mean square error (RMSE) and average running time (ART) of elevation estimation are used to evaluate the performance. Herein, RMSE is defined as where θ t/r,k (transmit or elevation) is the θ t/r,k 's estimation of the k − th trial. For the purpose of comparison, the performances of PM [16], PARAFAC [15], and ESPRIT [12] (all these algorithms are based on ULA geometry with N receivers and M transmitters) as well as CRB are added. e RMSE performance of different methods versus Signal-to-Noise Ratio (SNR) is presented in Figure 3, where we assume that M t � 3, N t � 4, M r � 4, N r � 5, and L � 200 are considered. It is shown that all algorithms achieve lower RMSE with higher SNR. When the SNR is lower than −4 dB, the proposed algorithm's performance is worse than those of the other algorithms, because it is impossible to accurately estimate the propagator in the low SNR region. However, the estimation result of the proposed algorithm is better than those of the other algorithms when the SNR is larger than −4 dB, which implies that the coprime manifold is helpful to improve estimation accuracy. Figure 4 depicts the ART performance comparison versus SNR. e results show that the proposed algorithm runs less than ESPRIT algorithm in [12], which benefits from the fact that the proposed algorithm can avoid the eigendecomposition of the high-dimensional covariance matrix. Figure 5 shows the RMSE of polarization angles estimation versus SNR, from which we can see that, with the increase of SNR, the RMSE of polarization angles estimation is decreasing. Figure 6 shows the different algorithms' estimation performance versus diverse snapshots, where the value range of snapshots L varies from 50 to 1000 with interval of 50. Meanwhile, we suppose that M t � 3, N t � 4 and M r � 4, N r � 5 with SNR � 10 dB. It is displayed that the performance of angle estimation which belongs to all algorithms improves with increase of L. Notably, it is obvious that the algorithm performs better than ESPRIT, PM, and PAR-AFAC. Figure 7 illustrates the ART performance. e result shows that the proposed algorithm's ART as a whole is much smaller than that of ESPRIT. Besides, the ART performance of the proposed algorithm is quite close to that of PARAFAC and PM. Moreover, with the increased number of L, the ART curve of the proposed algorithm performs basically the same as that of PM. e RMSE performance of polarization angles estimation versus L is shown in Figure 8, and it can be clearly seen that the higher L is, the lower RMSE is. Figure 9 presents the RMSE performance versus N, where SNR � 10 dB, L � 200, M t � 3, N t � 4, and M r � 1, respectively. N r is set to be varying from 3 to 27 with interval of 4, and the reason for operating above is to ensure M r and N r as coprime numbers. As expected, all algorithms provide better RMSE when N increases. Interestingly, when N > 7, the proposed algorithm provides obviously better RMSE performance compared to PM, ESPRIT, and PARAFAC in angle estimation. However, the RMSE performance of the proposed algorithm is worse than the remaining three when N < 5, which is owing to the fact that the propagator cannot be accurately estimated in the presence of smaller sensor number. Figure 10 gives the comparison results of ART with different N. Obviously, the number of N has a significant impact on ART performance. It is clearly manifested that the ART required by all algorithms increases significantly with N getting more. In addition, the ART required by the proposed algorithm is much smaller compared with those of ESPRIT, PM, and PARAFAC. e RMSE performance of polarization angles estimation is shown in Figure 11. It can be seen that RMSE decreases with the increase of N.

Conclusion
is paper proposes a new angle estimation algorithm for bistatic coprime MIMO radar. e core of the algorithm is to complete 2D-DOA and 2D-DOD estimation by constructing propagator. Because it does not need to involve high eigendecomposition and use coprime geometry with larger aperture, it is shown to be more computationally efficient than the state-of-the-art subspace algorithm. Compared with the existing algorithms, the proposed algorithm can achieve better estimation performance due to the larger virtual aperture. e effectiveness and improvement of the proposed algorithm are verified by numerical simulation.

Data Availability
No data were used in this paper.

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