Angle Estimation andMutual Coupling Self-Calibration in Bistatic MIMO System with Arbitrary Geometry

Ideal array responses are often desirable to a multiple-input multiple-output (MIMO) system. Unfortunately, it may not be guaranteed in practice as the mutual coupling (MC) effects always exist. Current works concerning MC in the MIMO system only account for the uniform array geometry scenario. In this paper, we generalize the issue of angle estimation and MC selfcalibration in a bistatic MIMO system in the case of arbitrary sensor geometry. The MC effects corresponding to the transmit array and the receive array are modeled by two MC matrices with several distinct entities. Angle estimation is then recast to a linear constrained quadratic problem. Inspired by the MC transformation property, a multiple signal classification(MUSIC-) like strategy is proposed, which can estimate the direction-of-departure (DOD) and direction-of-arrival (DOA) via two individual spectrum searches. Thereafter, the MC coefficients are obtained by exploiting the orthogonality between the signal subspace and the noise subspace. The proposed method is suitable for arbitrary sensor geometry. Detailed analyses with respect to computational complexity, identifiability, and Cramer-Rao bounds (CRBs) are provided. Simulation results validate the effectiveness of the proposed method.


Introduction
Multiple-input multiple-output (MIMO) is the technique with the most potential for the next-generation array radar system [1][2][3]. A MIMO system is characterized by multiple transmitting antennas and multiple receiving antennas. Unlike the traditional phased array radar, the MIMO system emits mutual orthogonal waveforms. The spatial diversity and waveform diversity enables the MIMO system to achieve a virtual aperture, which is much larger than the physical aperture of the radar system. Benefiting from the virtual aperture, the MIMO system outperforms the phased array radar with respect to clutter suppression, direction finding, and target tracking. Usually, the MIMO system can be divided into two categories in terms of its antenna configuration, namely, a distributed MIMO system and a colocated MIMO system [4,5]. A distributed MIMO system is equipped with widely displaced transmitting and receiving ðTx/RxÞ antennas, and it can illuminate the area of interest from various perspectives. Therefore, it is capable of extinguishing the radar cross-section (RCS) fluctuation issue. A colocated MIMO system consists of closely spaced Tx/Rx antennas, and it is able to provide super-resolution direction estimation. Herein, we focus on the bistatic MIMO system, which is a typical representative of the colocated MIMO system.
As a canonical issue in bistatic MIMO systems, joint direction-of-departure (DOD) and direction-of-arrival (DOA) estimation has been studied extensively in the past decades. Various methodologies have been proposed, e.g., estimation method of signal parameters via rotational invariance technique (ESPRIT) [6], multiple signal classification (MUSIC) [7], maximum likelihood [8], tensor approaches [9,10], and sparsity-aware methods [11,12]. A common assumption in [6][7][8][9][10][11] is that the sensors are well-calibrated. Nevertheless, owing to the radiation effects [13], sensors may be susceptible to unknown mutual coupling (MC), which would bring perturbation to the signal model and thus lead to degraded estimation performance. In general, the MC effect between sensors can be described by a MC matrix, in which the entities reveal the strength between two sensors. Although the MC matrix can be directly obtained offline via the measurement techniques such as the electromotive force [14], it may be invalid in the scenario where MC varies with time. To pursue online MC calibration, the usage of instrumental sensors or auxiliary sources have been reported [15,16], which require additional resource. As is known to us, the MC intensity of two sensors is inversely proportional to their distance. Therefore, the MC effects between sensors that are far apart can be neglected, i.e., the associated MC coefficients can be approximated by zeros. As a result, sparse arrays, e.g., coprime arrays [17,18], may be free from MC, but it is impossible for applications with limited space.
As mentioned earlier, the MC intensity between two sensors is inversely proportional to their interelement distance. Accordingly, the MC matrix associated with a uniform linear array (ULA) exhibits a special structure, which can be captured by a banded symmetric Toeplitz matrix. For a ULAbased MIMO system, several strategies have been proposed for MC self-calibration. In [19], a MUSIC-like estimator was presented, in which the MC coefficients are extracted into a single vector, and the DODs and DOAs are obtained via two one-dimensional (1D) peak searches. By interpreting the sensors at both ends of the arrays as instrumental sensors, the ESPRIT-like algorithms were introduced in [20,21], where the MC effect is suppressed by using a selection matrix, and then the ESPRIT technique is followed to obtain closeform solutions for DOD and DOA estimation. In order to exploit the tensor structure of the array data, the tensor approaches were subsequently discussed in [22][23][24], which offer more accurate angle estimation performance than the matrix-based ones. Besides, some frameworks have investigated into the two-dimensional (2D) angle estimation issue using a polarization-sensitive array [25][26][27][28].
It should be emphasized that the above algorithms are only effective for MIMO systems with ULA configurations. In practice, the ULA-based MIMO system may occasionally encounter sensor failure [29,30]. Besides, a more complex array manifold would be adopted in the MIMO system [31]. Therefore, the MIMO system with arbitrary Tx/Rx geometries may suffer from unknown MC, and current estimators in [19][20][21][22][23][24] will fail to work. To the best of our knowledge, only a few efforts have been devoted to robust DOA estimation for arbitrary geometry sensor array with unknown MC. In such a case, the MC matrix exhibits no special structure except Hermitian. In [32], an optimization-based framework was investigated. It formulates the DOA estimation as a block sparse inverse problem, which is solved via a joint-sparse recovery algorithm. Like the traditional counterpart, it has high complexity due to the high-dimension nonconvex problem, making it unacceptable for the MIMO system. In [33], an iterative approach was considered, in which DOA and MC coefficients are updated alternately. More recently, the MUSIC-like method was applied to robust DOA estimation in [34]. Although the approaches in [33,34] can be directly extended to angle estimation for a bistatic MIMO system, they are inefficient due to an exhaustive grid search.
In this paper, we generalize the MC issue in a bistatic MIMO system. We consider a general scenario where the Tx/Rx arrays are in arbitrary geometries and MC effects occur in both the Tx/Rx arrays. In such a scenario, the MC matrices exhibit no special structure except symmetry. To tackle the DOD and DOA estimation problem, the eigendecomposition is firstly performed to obtain the noise subspace. Then, by exploiting the orthogonality between the noise subspace and the virtual steering vector, angle estimation is recast to a general quadratic problem. Therefore, the MC transformation property is adopted to transform the optimization problem into two different forms, in which only DOD and DOA are included. Consequently, DOD and DOA can be estimated via two individual spectrum searches. Finally, the MC coefficients can be easily obtained with the estimated DOD and DOA. The proposed method is analyzed in terms of computational complexity, identifiability, and Cramer-Rao bounds (CRBs). Numerical simulations are designed to show the effectiveness of the proposed method.
Notation, bold capital letters, e.g., X, and bold lowercase letters, e.g., x, denote matrices and vectors, respectively. The superscripts ð·Þ T , ð·Þ H , and ð·Þ † account for transpose, Hermitian transpose, and inverse, respectively. The identity matrix is denoted by I, the M × N full one matrix is denoted by 1 M×N , the full zero matrix is denoted by 0. ⊗ , ⊙ , and ⊕ represent, respectively, the Kronecker product, the Khatri-Rao product, and the Hadamard product; diag f·g denotes the diagonalization operation; rank ð·Þ denotes rank operator; and Re ð·Þ and Im ð·Þ return the real part and the image part of a vector, respectively.

Problem Formulation
We consider a bistatic MIMO system with M transmit sensors and N receive sensors, and there are K far-field targets appearing in the same range bins of the radar system. The 2D-DOA pair and the 2D-DOA of the kth ðk = 1, 2,⋯,KÞ target are denoted by Θ t,k = ½θ t,k , ϕ t,k T and Θ r,k = ½θ t,k , ϕ t,k T , where θ t,k and θ r,k are elevation angles and ϕ t,k and ϕ r,k are azimuth angles. Without loss of generality, we assume that both the transmit sensors and the receive sensors are distributed in arbitrary geometries. The coordinate of the mth ðm = 1, 2,⋯,MÞ transmit sensor is p t,m = ½x t,m , y t,m , z t,m T , and the coordinate of the nth ðn = 1, 2,⋯,MÞ receive sensor is p r,n = ½x r,n , y r,n , z r,n T . The response vectors corresponding to the transmit array and the receive array are given by where λ is the carrier wavelength, τ m,k andΔ n,k are the associ- Wireless Communications and Mobile Computing ½cos ðϕ t,k Þ sin ðθ t,k Þ, sin ðϕ t,k Þ sin ðθ t,k Þ, cos ðθ t,k Þ T , and f r,k = ½cos ðϕ r,k Þ sin ðθ r,k Þ, sin ðϕ r,k Þ sin ðθ r,k Þ, cos ðθ r,k Þ T , respectively. DefineA t = ½aðΘ t,1 Þ, aðΘ t,2 Þ,⋯,aðΘ t,k Þ ∈ C M×K and Ar = ½aðΘ r,1 Þ, aðΘ r,2 Þ,⋯,aðΘ r,k Þ ∈ C N×K are the direction matrices corresponding to the transmit array and receive array, respectively. In the absence of the MC effect, the output from the matched filters can be written as where t is the slow time index; s k ðtÞ is the RCS coefficient of the kth targets; nðtÞ is the array noise, which is assumed to be Gaussian white with variance σ 2 ; and sðtÞ = ½s 1 ðtÞ, s 2 ðtÞ,⋯,s k ðtÞ T . Now, we consider a scenario where both the transit sensors and the receive sensors suffer from unknown MC and the MC matrices are denoted by C t ∈ C M×M and C r ∈ C N×N , respectively, which are Hermitian. In the presence of an unknown MC, the matched output in (2) becomes Suppose the RCS coefficients are uncorrelated, and sðtÞ is irrelevant with nðtÞ, and then, the covariance matrix of yðtÞ is where R s = diag fγ 1 , γ 2 ,⋯,γ k g is the covariance matrix of sðtÞ , γ k is the power of s k ðtÞ. With L available snapshots, i.e., t = 1, 2, ⋯, L, R can be estimated viâ 3. The Proposed Algorithm 3.1. DOD and DOA Estimation. Performing eigendecomposition on R, one can obtain the signal subspace and the noise subspace as where Λ s ∈ C K×K and Λ n ∈ C ðMN−KÞ×ðMN−KÞ are two diagonal matrices that contain the K dominate eigenvalues and the remind eigenvalues, respectively. U s ∈ C MN×K and U n ∈ C MN×ðMN−KÞ consist of the corresponding eigenvectors, which are called the signal subspace and the noise subspace, respectively. It is known that U s is orthogonal to U n , i.e., Since U s span the same subspace with the columns of A, we get The above property is similar to that in the 2D-MUSIC algorithm. However, the traditional MUSIC estimator is invalid as both C t and C r are unknown. To further process, the following conclusion in [33] is introduced.
where T ∈ C M×Q with the qth ðq = 1, 2,⋯,QÞ column given by and J q is defined as In this paper, we build on top of Theorem 1 to address the MC self-calibration challenge in the arbitrary geometrybased MIMO system. Suppose that there are P distinct entities c t = ½c t,1 , c t,2 ,⋯,c t,p T in C t , and Q distinct entities c r = ½c r,1 , c r,2 ,⋯,c r,Q T in C r . From Theorem 1, we have where TðΘ t Þ ∈ C M×P and TðΘ r Þ ∈ C N×Q are constructed according to Theorem 1. In combination with the property ðA ⊗ BÞðC ⊗ DÞ = ðACÞ ⊗ ðBDÞ, we can obtain where eðΘ t,k Þ = c t ⊗ ½C t aðΘ t,k Þ and fðΘ r,k Þ = ½C r aðΘ r,k Þ ⊗ c r .
To obtain the DOA, one needs to optimize Notably, (15) is a quadratic optimization issue. To avoid the trivial solution e = 0 MP×1 , the above optimization problem is constrained with d H e = ρ, i.e., where ρ is a constant and d is a vector with the first entity being one and zeros elsewhere. To solve the above issue, we can construct the following Lagrange function: where T is a Lagrange multiplier. Enforcing ∂LðΘÞ/c to zeros yields Plug (18) into (15) and remove the constant item gives By setting a grid to contain all possible DODs, the indexes corresponding to the maximum values in (19) reveal the DODs. Similarly, we can obtain the DOA estimation via computing where g is a NQ × 1 vector with the first entity one and zeros elsewhere.

MC Self-Calibration.
Although DOD Θ t,k and DOA Θ r,k can be separately estimated from (19) and (20), they require further pairing, since there are KðK − 1Þ/2 possible DOD-DOA pairs. Recall the properties in (12) and we can rewrite (8) as where c tr = c r ⊗ c t . Similar to (15), to find the correct DOD-DOA pairs, we need to calculate which can be accomplished via where h is a PQ × 1 vector with the first entity one and zeros elsewhere. List all the possible DOD-DOA pairs, and the real DOD-DOA pairs can be picked out by finding the K maximum values in (23).
It should be pointed out that the optimal solution of (22) is achieved if where α is a constant. After the DODs and DOAs have been paired, c tr can be estimated viâ where b Θ t,k and b Θ r,k denote the estimated Θ t,k and Θ r,k , respectively. Usually, we have c t,1 = 1 and c r,1 = 1. After normalizingĉ t,r , the MC vectors can be estimated viâ

Computational Complexity.
The computational burden (the number of complex multiplication) of the proposed is summarized as follows: Step 1. Estimate the covariance matrixR ⋯ ⋯ ⋯ M 2 N 2 L.
Step 3. Calculate the spectrum function in (19) to get the DOD estimation ⋯⋯ ⋯ Of M 3 N 2 Qg.

Wireless Communications and Mobile Computing
Step 4. Search the maximums of (20) to achieve the DOA estimation ⋯⋯ ⋯ OfM 2 N 3 Pg.

Identifiability. The proposed algorithm is effective if
where y t ≜ yðtÞ, s t = sðtÞ, and n t = nðtÞ: Then, all the measurements are arranged into a "column" vector as y = ½y T t , y T 2 ⋯,y T L T ∈ C MNL×1 . Suppose that the s t is deterministic but unknown. Thereafter, the mean μ ∈ C MNL×1 of y is given by where Define the unknown parameter vectors θ t ≜ ½θ t,1 , θ t,2 ,⋯,θ t,k T , ϕ t ≜ ½ϕ t,1 , ϕ t,1 ⋯,ϕ t,1 T , θr, ½θ r,1 , θ r,1 ,⋯,θ r,k T , ϕ r ≜ ½ϕ r,1 , ϕ r,2 ⋯,ϕ r,k T , α ≜ ½θ T t , θ T r , ϕ T t , ϕ T r T ∈ R 4K×1 , β ≜ ½Re fc T r g, Re fc T t g, Im fc T r g, Im fc T t g T ∈ R 2ðP+QÞ×1 , and γ ≜ ½Re fF T g, Im fF T g T ∈ R 2LK×1 . The whole estimation parame- It is well known that the CRB matrix for ζ is given by where Ψ ≜ ½∂μ/∂α T , ∂μ/∂β T , ∂μ/∂γ T . Next, we will concentrate on each part of Ψ. Before the detailed derivations, the following variables are defined: 1 , 1 , It is easy to find that∂μ/∂θ T T = S ⊙Ã 1 . Stepping further, we can get As ∂μ/∂ Re fc t,p g = S ⊙Ã t,p , and ∂μ/∂ Im fc t,p g = jðS ⊙ A t,p Þ, thus, we have Besides, it is straightforward to get Consequently, it is established that Since we are interested in CRBs on angle and MC estimation only, we will remove the uninterested parts from J. Define Both P −1 ∇ and P −1 ∇ are valid as H H H is nonsingular. Furthermore, the following transform matrix is defined:

Simulation Results
In this section, numerical simulations are carried out to show the effectiveness of the proposed estimator. Specifically, we assume that the MIMO system is configured with M transmit sensors and N receive sensors. Suppose that there are K targets, and each RCS fulfills the Swerling II model and L = 500 snapshots are collected. Search interval in the simulation is set to Δ. The following two scenarios are considered.
In the first example, we test the spatial spectrum of the proposed estimator in Scenario I, where Δ = 0:01°, signalto-noise ratio (SNR) is set to 15 dB, and 10 independent trials are recorded. The results are presented in Figure 2, from which we observe that the proposed estimator can accurately identify all the targets. For comparison purposes, the spectrum results of the traditional reduced-dimension MUSIC (RD-MUSIC) algorithm in [7] are added, which display great disparities for some angles. This is caused by the fact that the traditional RD-MUSIC cannot eliminate the MC effect.
In the second example, we give the 2D spectrum results of the proposed estimator in Scenario II, where Δ = 0:01°and SNR = 15 dB are considered. Figure 3 illustrates the 2D spectrum results of the proposed estimator. Clearly, the proposed estimator can correctly recover the 2D-DODs and the 2D-DOAs, which is evident that the proposed estimator is also suitable for the MIMO system with MC and 3D sensor geometry.
In the third example, we plot the 2D spectrum results of the proposed estimator in Scenario III. Figure 4 gives the results. Similar to the previous simulation, the proposed esti-mator can correctly estimate the 2D-DODs and the 2D-DOAs in such a condition.
In the fourth example, we examine the average running time of the proposed estimator against different search interval Δ in Scenario I, where SNR is fixed at 20 dB. The result is displayed in Figure 5. As expected, a refined search interval would result in more computational burden. For comparison purposes, the performances associated with the reduceddimension MUSIC in [32] (marked with "RD-MUSIC") and the iterative method in [33] (marked with "Iterative method") as well as the MUSIC-like method in [34] (marked with "MUSIC-like") are added. In contrast with RD-MUSIC, the proposed estimator requires a slightly more running time. However, it is shown that the running time of the proposed estimator is two orders of magnitude lower than that of the Iterative method and MUSIC-like. This is because both RD-MUSIC and the proposed estimator require two 1D spectrum searches, while the 2D spectrum search is essential in both the Iterative method and MUSIC-like.
Finally, we plot the average root mean square error (RMSE) curves of the proposed estimator versus SNR in Scenario I, where Δ = 0:01°and 200 Monte Carlo trials are carried out. As shown in Figure 6, the proposed estimator     Wireless Communications and Mobile Computing deteriorates in the low SNR regime (e.g., SNR < 4 dB). With the increasing SNR, the proposed estimator would achieve improved RMSE performance, and it may attain the CRB at high SNR regions. An interesting observation is that the proposed estimator offers better performance than MUSIC-like when SNR is smaller than 6 dB, and it provides very close RMSE performance to the latter when SNR is larger than 6 dB, especially in MC estimation. Nevertheless, both RD-MUSIC and Iterative method provide the unacceptable RMSE performance over the entire SNR regions.

Conclusion
In this paper, we presented a MUSIC-like estimator for bistatic MIMO radar with arbitrary geometries and MC. For the 1D sensor manifold scenario, the proposed estimator only needs two one-dimensional spectrum peak searches to obtain 1D-DODs and 1D-DOAs, while for 2D or 3D sensor geometries, the proposed estimator allows two 2D spectrum searches to achieve 2D-DODs and 2D-DOAs. The proposed estimator is efficient from the viewpoint of complexity. It should have a brighter prospect in numerous applications, e.g., military explorations and internet-of-vehicles.

Conflicts of Interest
The authors declare no conflict of interest.