Angle-Range-Polarization Estimation for Polarization Sensitive Bistatic FDA-MIMO Radar via PARAFAC Algorithm

Abstract

In this paper, we study the estimation of angle, range and polarization parameters of a bistatic polarization sensitive frequency diverse array multiple-input multiple-output (PSFDA-MIMO) radar system. The application of polarization sensitive array in receiver is explored. A signal model of bistatic PSFDA-MIMO radar system is established. In order to utilize the multi-dimensional structure of array signals, the matched filtering radar data can be represented by a third-order tensor model. A joint estimation of the direction-of-departure (DOD), direction-of-arrival (DOA), range and polarization parameters based on parallel factor (PARAFAC) algorithm is proposed. The proposed algorithm does not need to search spectral peaks and singular value decomposition, and can obtain automatic pairing estimation. The method was compared with the existing methods, and the results show that the performance of the method is better. Therefore, the accuracy of the parameter estimation is further improved.

1. Introduction

Multiple-Input Multiple-Output (MIMO) radar is a new system radar that introduces multiple input and multiple output technologies in wireless communication systems into the radar field and combines them with digital array technology. MIMO radar is a new type of radar system that uses multiple transmitting antennas to simultaneously transmit diversity waveforms, and uses multiple receiving antennas to receive echo signals, and processes them in a centralized manner [1]. Unlike phased-array radar, which only uses the degrees of freedom at the receiving end, MIMO radars simultaneously use the degrees of freedom at the receiving and transmitting ends to improve the radar's diversity gain.

However, whether it is a phased array radar or a MIMO radar, when performing beam scanning, the beam pointing is only related to the spatial orientation, and the correlation of the spatial distance is ignored. To this end, researchers introduced the concept of Frequency diversity Array (FDA) [2]. FDA achieves the function of automatic beam scanning by using a small frequency difference between the array elements [3]. FDA form range-dependent beams,so the transmitting waveform of FDA radar is range dependent, and the echo signal of FDA radar contains the range information of the target. Thus, scholars can estimate the angle and range of the target simultaneously [4,5]. Because of the coupling problem of angle-range information. In order to solve the coupling problem of angle-range information, Waseem Khan proposed logarithmic increase frequency offset in [6], and [7] proposed non-linear increase frequency offset, which was solved by increasing square and cubic increase frequency offset. In [8], the angle-range decoupling is achieved by dividing the array into several non-overlapping sub-arrays with different frequency shifts. several non-overlapping sub-arrays with different frequency shifts.

The FDA was introduced into MIMO radar system, and a FDA-MIMO radar technology[9-11] was proposed to solve the coupling problem of angle and range. In [12], a clear range angle estimation method is proposed by utilizing the range angle correlation characteristics of FDA-MIMO radar. In [13], the waveform design based on subarray and decoupling algorithm was adopted to realize the decoupling of the range of the transmission vector range, and the subspace algorithm is used for parameter estimation. In [14], in order to overcome the coupling of direction-of-departure (DOD) and range parameters of FDA transmitting array guidance vector, Can Cai et al. proposed a nonlinear frequency increment method. They use estimation signal parameter via rotational invariance techniques (ESPRIT) method and PARAFAC method to estimate parameters, which reduces the computational complexity.

In addition, a polarimetric MIMO radar system [15] has been proposed by scholars. In [16], a new method for joint estimation of DOD, direction-of-arrival (DOA) and polarization parameters of bistatic MIMO radar is proposed. The joint spectrum in polarization and spatial domains is defined by using multiple signal classification (MUSIC) method. The estimation accuracy and resolving power under different signal-to-noise ratio (SNR) and snapshots are investigated to measure the performance of the method based on this structure [17]. In the polarimetric bistatic MIMO radar system, an algorithm of one-dimensional MUSIC continuous search angle estimation with automatic matching is proposed [18]. In [19], based on the invariance of the virtual array and polarization sensitive array generated by matched filter, a new ESPRIT algorithm is proposed to estimate the angle and polarization parameters of the target. A new radar system is proposed in [20], namely PSFDA-MIMO radar. In [21], anon-search algorithm is proposed to estimate 4-D PSFDA-MIMO radar parameters.

In recent years, tensor method has been proposed as a general tool for storing and manipulating multidimensional data, which can make full use of multidimensional information [22]. Tensor decomposition is the core of signal processing. Tensor decomposition, such as PARAFAC decomposition and Higher-order Singular Value Decomposition (HOSVD), can be used to construct three-dimensional received signal model in radar system. The PARAFAC algorithm has been applied to parameter estimation of MIMO radar [23-25].

In this paper, we study the joint estimation of the DOD, DOA, range and polarization parameters of the polarimetric bistatic FDA-MIMO radar. Tensor decomposition model is used to construct three-dimensional received signal model, and trilinear alternating least squares (TALS) are used to estimate the angle, range and polarization parameters. The algorithm does not need spectral peak search and singular value decomposition, and can realize automatic matching. The results show that the PARAFAC algorithm can effectively estimate the parameters of PSFDA-MIMO radar

The rest of the structure of the paper is given below. The tensor-based signal model of bistatic PSFDA-MIMO radar is introduced, in Section 2. In Section 3, we propose angle-range-polarization estimation for bistatic PSFDA-MIMO radar via PARAFACalgorithm. The theoretical results of Section 3 are proved by simulation results, In Section 4. Finally, the paper ends with a brief conclusion in Section 5.

Notations: Scalars are represented in lowercase letters, e.g., a , vectors are represented in bold lowercase letters, e.g.,a , matrix is represented by capital letters, e.g., Α , tensor by bold Euler script letters, e.g., Χ . The transpose operator is represented as (·)^T·(·)^† represents the pseudo-inverse. The ||·||_F  denote the Frobenius norm, the identity matrix of n× n is expressed as I_n . The tensor [Χ]_(i) stands for the i-mode unfolding of Χ , and its inverse operation, which coverts [Χ]_(i) to Χ , is represented as folding ([Χ]_(i)) . ⊗ denotes the kronecker product, • denotes the Khatri-Rao product (column-wise kronecker product), and 。represents the outer product.

2. Signal Model of Polarimetic Bistatic FDA-MIMO Radar

Considering a bistatic PSFDA-MIMO radar system, which the transmitting array is an uniform linear array (ULA) composed of M transmitting antennas and the receiving array consisting of N crossed dipoles, which is a polarization sensitive ULA. The spacing of the transmitting antennas is denoted by t d and the spacing of the receiving antennas is represented by r d . θ and ϕ represent DOD and DOA, respectively, as you can see Fig. 1.

Fig. 1. Bistatic PSFDA-MIMO radar.

The carrier frequency of the mth element of the transmitting elements can be defined as:

\(f_{m}=f_{1}+(m-1) \Delta f, \quad m=1,2, \ldots, M\)       (1)

where f₁ denotes the reference carrier frequency, and ∆f denotes the frequency increment between adjacent array elements. Assume K targets, the received signal of PSFDA-MIMO radar in the qth pulse is expressed as:

\(\mathbf{X}_{q}=\sum_{k=1}^{K} \mathbf{b}\left(\varphi_{k}, \xi_{k}, \eta_{k}\right) \operatorname{diag}\left(\mathbf{c}_{q}\right) \mathbf{a}_{t}^{T}\left(\theta_{k}, r_{k}\right) \mathbf{S}+\mathbf{W}_{q}, q=1,2, \ldots, Q\)       (2)

where c_q = [\(\beta_{1q}e^{j2\pi{f}_{d1}T_{r}} \),\(\beta_{2q}e^{j2{\pi}f_{d2}T_{r}}\),...,\(\beta_{kq}e^{j2{\pi}f_{dk}T{r}}\)]^T ·β_kq is the complex amplitude and f_dk is the Doppler frequency. Tr is the pulse repetition time. S= [s₁ ,s₂ ,...s_m ]^T is the M × L transmit signal matrix, L represents the pulse code length within a single pulse period. W_q denotes Gaussian white noise vector of N × L with a mean value of zero. Then a_t (θ_k ,r_k) is the transmit steering vectors of the kth target written as:

\(\mathbf{a}_{t}\left(\theta_{k}, r_{k}\right)=\left[1, e^{j \frac{2 \pi}{c}\left(d_{t} f_{1} \sin \left(\theta_{k}\right)-\Delta fr_{k}\right)}, \ldots, e^{j(M-1) \frac{2 \pi}{c}\left(d_{t} f_{1} \sin \left(\theta_{k}\right)-\Delta fr_{k}\right)}\right]^{T}\)       (3)

The receive steering vector for the kth target is expressed as:

b(φ_k,ξ_k,η_k)= a_r(φ_k) ⊗ v(φ_k,ξ_k,η_k)       (4)

\(\mathbf{a}_{r}\left(\varphi_{k}\right)=\left[1, e^{j 2 \pi \frac{f 1}{c} d_{r} \sin \left(\varphi_{k}\right)}, \ldots, e^{j 2 \pi(N-1) \frac{f_{1}}{c} d_{r} \sin \left(\varphi_{k}\right)}\right]^{T}\)       (5)

The kth polarization vector of arrival signal is given by:

\(\mathbf{v}\left(\varphi_{k}, \xi_{k}, \eta_{k}\right)=\left[-\cos \xi_{k}, \cos \varphi_{k} \sin \xi_{k} e^{j \eta_{k}}\right]^{T}\)        (6)

where ξ_k represents the auxiliary polarization angle and 0 ≤ ξ_k ≤\(\frac{\pi}{2}\) , η_k represents polarization phase difference and − \(\pi\) <η_k < \(\pi\) .

Processing the received signal through a matched filter, which output can be expressed as:

\(\mathbf{Y}_{q}=\mathbf{B} \operatorname{diag}\left(\mathbf{c}_{q}\right) \mathbf{A}^{T}+\mathbf{Z}_{q}, q=1,2, \cdots, Q\)       (7)

where\(A = [a_t(\theta_1,r_1),a_t(\theta_2,r_2),...,a_t(\theta_k,r_k) ]\) is the M × K transmit steering matrix, K is the number of targets and B = B = [b(φ₁ξ₁η₁),b(φ₂ξ₂η₂),...,b(φ_kξ_kη_k)] is the N × K receive steering matrix. Ζ_q is the noise matrix after pulse compression.

The traditional bistatic PSFDA-MIMO radar received signal based on two-dimensional signal model is processed by matched filtering, and the data is vectorized, so that the matrix whose dimension is N × M becomes a vector of NM ×1. Then the received signal model under the pulse is given below:

\(\begin{aligned} \mathbf{Y} &=\left[\operatorname{vec}\left(\mathbf{Y}_{1}\right), \operatorname{vec}\left(\mathbf{Y}_{2}\right), \ldots, \operatorname{vec}\left(\mathbf{Y}_{Q}\right)\right] \\ &=(\mathbf{B} \bullet \mathbf{A}) \mathbf{C}^{T}+\mathbf{Z} \end{aligned}\)       (8)

According to the original data structure, the array data can be written as a trilinear model:

\(\mathcal{Y}=\sum_{k=1}^{K} b_{k} \circ a_{k} \circ c_{k}+\mathcal{Z}\)       (9)

Formula (9) can also be expressed as:

\(\mathcal{Y}=\mathcal{I}_{K} \times_{1} \mathbf{B} \times_{2} \mathbf{A} \times_{3} \mathbf{C}+\mathcal{Z}\)       (10)

where \(\mathcal{I}_{K}\) is the K × K × K identity tensor, \(C = [c_1,c_2,...,c_Q]^T \), Z is the 2N × M ×Q noise tensor.

3. Proposed Algorithm

3.1 Trilinear decomposition

Formula (8) can be considered as the 1-mode matrix unfolding of Υ , which is a slice of three-dimensional data along the spatial direction. Trilinear decomposition can also obtain two equivalent matrices as follows:

Y₍₂₎ = (C•B)A^T + Z₂       (11)

Y₍₃₎ = (C•B)A^T + Z3       (12)

TALS algorithm is an efficient trilinear decomposition modulus algorithm. It uses least squares cost function to fit three matrices alternately. When the fitting error reaches the expected range, the algorithm terminates. The steps of the trilinear model are given below:

Step 1: Assuming that the two matrices in Υ₍₁₎ , Υ₍₂₎ and Υ₍₃₎ are known, the LS method is used to fit any of them;

Step 2: Continue to update other matrices with LS method;

Step 3: Repeat step 1 and step 2 the fitting error reaches a certain fitting threshold.

The specific process of TALS algorithm is given below:

Using equation (8), the LS of Υ₍₁₎ is fitted as:

\(f_{\mathbf{Y}_{(1)}}=\min _{\mathbf{B}, \mathbf{A}, \mathbf{C}}\left\|\mathbf{Y}_{(1)}-(\mathbf{B} \bullet \mathbf{A}) \mathbf{C}^{T}\right\|_{F}\)       (13)

Then according to formula (13), the LS estimation of C is

\(\hat{\mathbf{C}}^{T}=(\hat{\mathbf{B}} \cdot \hat{\mathbf{A}})^{\dagger} \mathbf{Y}_{(1)}\)       (14)

where \(\rm \hat{B}\) and \(\rm \hat{A}\) are the estimates ofΒ and Α obtained after the last iteration, respectively. Similarly, using equation (11), the LS of Υ₍₂₎ is fitted as:

\(f_{\mathbf{Y}_{(2)}}=\min _{\mathbf{C}, \mathbf{B}, \mathbf{A}}\left\|\mathbf{Y}_{(2)}-(\mathbf{C} \bullet \mathbf{B}) \mathbf{A}^{T}\right\|_{F}\)       (15)

The LS estimation of Α is

\(\hat{\mathbf{A}}^{T}=(\hat{\mathbf{C}} \bullet \hat{\mathbf{B}})^{\dagger} \mathbf{Y}_{(2)}\)       (16)

where \(\rm \hat{C}\) and \(\rm \hat{B}\) are respectively the estimates of Cand Β obtained in this iteration process. Similarly, the LS fitting of Υ₍₃₎ is

\(f_{\mathbf{Y}_{(3)}}=\min _{\mathbf{A}, \mathbf{C}, \mathbf{B}}\left\|\mathbf{Y}_{(3)}-(\mathbf{A} \bullet \mathbf{C}) \mathbf{B}^{T}\right\|_{F}\)       (17)

The LS estimation of Β is

\(\hat{\mathbf{B}}^{T}=(\hat{\mathbf{A}} \bullet \hat{\mathbf{C}})^{\dagger} \mathbf{Y}_{(3)}\)       (18)

where \(\rm \hat{A}\) and \(\rm \hat{C}\) are estimates of Α and C obtained after previous iterations, respectively. Using the above method, each time a matrix is updated, the matrix updated this time immediately participates in the solution of the next matrix until the algorithm converges.

3.2 angle-range-polarization estimation

Theorem 1. [22] For the trilinear model in equation (8), Α_k , Β_k and C_k respectively are the Kruskal-rank ( k -rank) of Α , Β and C. If the k-rank of the matrix satisfies:

k_A + k_B + k_C \(\geq\) 2K +2       (19)

then Α , Β and Care unique up to permutation and scaling of columns.

The estimates of Α , Β and C are \(\hat{\mathbf{A}} \) , \(\hat{\mathbf{B }}\) and \(\hat{\mathbf{C}}\) respectively, satisfying the following:

\(\begin{array}{l} \hat{\mathbf{A}}=\mathbf{A} \prod \Delta_{1}+\varpi_{1} \\ \hat{\mathbf{B}}=\mathbf{B} \prod \Delta_{2}+\varpi_{2} \\ \hat{\mathbf{C}}=\mathbf{C} \prod \Delta_{3}+\varpi_{3} \end{array}\)       (20)

where ∏ is a permutation matrix, \(\bar{\omega}_1\), \(\bar{\omega}_2\) and \(\bar{\omega}_3\) are the corresponding estimation error matrices. ∆1 , ∆2 and ∆3 are three diagonal matrices, and they satisfy Δ₁Δ₂Δ₃ = I_k .

After TALS process, the estimates of Α and Β can be obtained. \(\hat{a_t}(\theta_k,r_k) \) are the estimation of the transmitting steering vectors of the kth target. First the DOD and distance of the kth target are estimated, the phase angle can be obtained as follows:

\(\hat{h}_{t}=\operatorname{angle}\left(\hat{\mathbf{a}}_{t}\left(\theta_{k}, r_{k}\right)\right)\)       (21)

We construct two matrices as follows:

\(\mathrm{P}=\left[\begin{array}{ccccc} 1 & 1 & \ldots & 1 \\ 0 & 2 \pi \frac{f_{1}}{c} d_{t} & \ldots & (M-1) 2 \pi \frac{f_{1}}{c} d_{t} \\ 0 & -2 \pi \frac{\Delta f}{c} & \ldots & -2 \pi \frac{\Delta f}{c} \end{array}\right]\)       (22)

\(\hat{r}_{k}=\frac{2 \pi d_{t} \frac{f_{1}}{c} \sin \left(\hat{\theta}_{k}\right)-\operatorname{angle}\left(\hat{\mathbf{a}}_{t}\left(\theta_{k}, r_{k}\right)\right)}{2 \pi \frac{\Delta f}{c}}\)        (23)

The least square solution of is u_t

\(\mathbf{u}_{t}=\left[u_{t}, \sin \left(\theta_{k}\right), r_{k}\right]^{T}\)       (24)

Range estimation can be calculated as follow:

\(\hat{r}_{k}=\frac{2 \pi d_{t} \frac{f_{1}}{c} \sin \left(\hat{\theta}_{k}\right)-\operatorname{angle}\left(\hat{\mathbf{a}}_{t}\left(\theta_{k}, r_{k}\right)\right)}{2 \pi \frac{\Delta f}{c}}\)       (25)

Then the kth DOD and range are estimated to be

\(\hat{\theta}_{k}=\arcsin \left(\mathbf{u}_{t}(2)\right) \frac{180^{\circ}}{\pi} \)       (26)

\(\hat{r}_{k}=\mathbf{u}_{t}(3)\)      (27)

Secondly, the DOA is estimated. The estimated value of matrix B includes the estimation of DOA and polarization parameters. According to formula (4), B is the Kronecker product of a_r and υ, which can be expressed as follows:

\(\mathbf{B}\left(\varphi_{k}, \xi_{k}, \eta_{k}\right)=\left[\begin{array}{ccc} \mathbf{v} & \cdots & \mathbf{v} \\ \mathbf{a}_{r 1} \mathbf{v} & \cdots & \mathbf{a}_{r k} \mathbf{v} \\ \vdots & \vdots & \vdots \\ \mathbf{a}_{r 1}^{N-1} \mathbf{v} & \cdots & \mathbf{a}_{r k}^{N-1} \mathbf{v} \end{array}\right]\)       (28)

According to formula (28), the ratio between the n +1 and n rows of matrix B is \([a_{r1},a_{r2},...,a_{r3}] \) , where n N =1,2,.., N . So the average value of the ratio of all rows of matrix B is:

\(\hat{\mathbf{a}}_{r k}=\frac{1}{N-1} \sum_{i=1}^{N-1} \frac{\left[\hat{\mathbf{B}}\left(\varphi_{k}, \xi_{k}, \eta_{k}\right)\right]_{N+1, k}}{\left[\hat{\mathbf{B}}\left(\varphi_{k}, \xi_{k}, \eta_{k}\right)\right]_{N, k}}\)       (29)

So, DOA is estimated to be

\(\hat{\mathbf{a}}_{r k}=\frac{1}{N-1} \sum_{i=1}^{N-1} \frac{\left[\hat{\mathbf{B}}\left(\varphi_{k}, \xi_{k}, \eta_{k}\right)\right]_{N+1, k}}{\left[\hat{\mathbf{B}}\left(\varphi_{k}, \xi_{k}, \eta_{k}\right)\right]_{N, k}}\)       (30)

ς_k is defined as the ratio between the second and first elements of the polarization vector υ(ϕ_kξ_kη_k , , ) in formula (6)

\(\varsigma_{k}=\frac{\sin \xi_{k} \cos \varphi_{k} e^{j \eta_{k}}}{-\cos \xi_{k}}\)       (31)

Therefore, the polarization parameters ξ k and ηk of the kth target are calculated by

\(\xi_{k}=\arctan \left(\left|\frac{\varsigma_{k}}{\cos \varphi_{k}}\right|\right)\)       (32)

\(\eta_{k}=\operatorname{angle}\left(\varsigma_{k}\right)\)       (33)

4. Experimental Results and Analysis

We evaluate the performance of this algorithm in estimating angle-range-polarization parameters of bistatic PSFDA-MIMO radar. Assuming that the number of transmitting antennas at the transmitter of PSFDA-MIMO radar system is M , the number of receiving antennas at the receiver is N , there are K targets, the number of pulses isQ , and the noise considered is white noise, the signal-to-noise ratio (SNR) can be expressed as:

\(\mathrm{SNR}=10 \log _{10}\left(\sum_{q=1}^{Q}\left\|\mathbf{B} \sum_{q} \mathbf{A}^{T} \mathbf{S}\right\|_{F}^{2} / \sum_{q=1}^{Q}\left\|\mathbf{W}_{q}\right\|_{F}^{2}\right)\)       (34)

To assess the performance of ρ parameters estimation, we use the root mean square of error(RMSE) define as:

\(\mathrm{RMSE}=\frac{1}{K} \sum_{k=1}^{K} \sqrt{\frac{1}{200} \sum_{i=1}^{200}\left(\hat{\rho}_{k}-\rho_{k}\right)^{2}}, \rho=\{\theta, r, \phi, \xi, \eta\}\)       (35)

where \(\hat{\rho}_k\) is the estimates of DOD, DOA, range, polarization angle and polarization phase for 200 Monte Carlo trials.

In the first simulation, we use the successive ESPRIT algorithm [20] and ESPRIT-like algorithm [26] to compare with the proposed method. We consider K = 3 targets, the DOD,DOA, range, polarization angle parameters and polarization phase are respectively θ = (10° ,25° ,50° ), ϕ (-5° ,20° ,40° ), r = (10km ,20km ,30km ), ξ = ( π/10, π/4,9 π/20), η = (π/5,4π/5,2π/5) . The radar system parameters for this experiment are set to: M = 20 transmit antennas, N = 20 receive antennas, the number of snapshots is Q =128 , the value of SNR is from −10dB to 20dB with step size 5dB . As can be seen from Fig. 2, RMSE of the three algorithms decreases with the increase of SNR. The algorithm has better performance and estimation accuracy.

Fig. 2. RMSE versus SNR. (a) DOD estimation. (b) DOA estimation. (c) range estimation. (d) polarization angle ξ estimation. (e) polarization phase η estimation.

In the second simulation, we can see Fig. 3, the variation curve of RMSE with the number of snapshots. The radar system parameters for this experiment are set to: M = 25 transmit antennas, N = 25 receive antennas, SNR 10 = − , the number of snapshots is from 100 to 700 every 100 counts. Suppose there are K = 3 targets, the angle, range and polarization parameters information of the target and other parameter information are the same as Experiment 1. It can be seen from the figure that the RMSE curve of the three algorithms decreases with the increase of the number of snapshots, which means that the accuracy of the angle, range and polarization parameters estimation of the algorithm increases with the number of snapshots. Moreover, the proposed in this paper has better effect than the successive ESPRIT.

Fig. 3. RMSE versus number of snapshots. (a) DOD estimation. (b) DOA estimation. (c) range estimation. (d) polarization angle ξ estimation. (e) polarization phaseη estimation.

5. Conclusion

This paper investigates the estimation of the angle, range and polarization parameters of the bistatic PSFDA-MIMO radar system. The signal model of bistatic PSFDA-MIMO radar isestablished, and a joint estimation algorithm of DOD, DOA, range and polarization parameters based on PARAFAC is proposed. The algorithm does not need to search for spectral peaks, and can automatically pair and improve the estimation accuracy. The experimental results demonstrate that the algorithm is effective. In future work, we introduce matrix completion into PSFDA-MIMO radar system to recover the lost data of received signals.