Multiple-input multiple-output (MIMO) radar [8, 3, 15, 10, 11, 9, 12, 16, 17, 1, 13], which generally splits into colocated MIMO radar  and widely separated MIMO radar , has received a lot of attention over the past decade due to the advantages in multiple targets detection , parameters estimation [9, 12] and many other applications . MIMO radar simultaneously emits several orthogonal waveforms via colocated or widely separated antennas to achieve waveform/special diversity. For the case of colocated MIMO radar, the waveform diversity can also be achieved in Doppler domain. The corresponding MIMO radar is named as slow-time MIMO radar [16, 17, 13], while the associated waveform design approach is called Doppler division multiple access (DDMA) [1, 13]. The main idea of DDMA is to apply diverse phase modulation values at each transmitter from pulse-to-pulse so that every transmit waveform possesses independent Doppler frequency. Slow-time MIMO radar is approximately equivalent to its conventional MIMO counterpart with reduced Doppler estimation range, but simple waveforme design.
In bistatic colocated MIMO radar, many algorithms for joint direction-of-departure (DOD) and direction-of-arrival (DOA) estimation have been proposed [24, 7, 26, 27, 19, 4, 20, 14, 18, 2]. For example, joint DOD and DOA estimation can be conducted by multiple signal classification (MUSIC) which generally demands two dimensional (2D) spectrum search . By exploiting rotational invariance property (RIP) of signal subspace, estimation of signal parameters via rotational invariance technique (ESPRIT) can be applied to estimate angle information without spectrum search whereas DOD and DOA pairing is still required . In , a generalized algorithm called unitary-ESPRIT (U-ESPRIT) has been introduced to reduce computational complexity. Propagator method (PM) has been also proposed in 
to avoid singular value decomposition (SVD). The aforementioned algorithms can be regarded as signal subspace-based methods, which normally ignore the multi-linear structure of received data and have poor performance at low signal-to-noise ratio (SNR). A possible solution to overcome these disadvantages is to store the received signals in a tensor. In[19, 4, 20, 14, 18, 23], parallel factor (PARAFAC) analysis has been applied to address the problem of poor estimation performance at low SNR. However, conventional tensor model is improper for slow-time MIMO radar because the significantly reduced number of samples in Doppler domain causes a performance loss . Thus, the tensor model-based approach for angle estimation in slow-time MIMO radar needs further investigation.
In this letter, we develop a new tensor model for bistatic slow-time MIMO radar and apply it to joint DOD and DOA estimation. In our model, the received signals are organized in a 3-order tensor. Then modified alternating least squares (ALS) algorithm with extra constant terms is introduced to estimate factor matrices. Finally, the Vandermonde structure of the transmit and receive steering matrices is fully exploited to obtain the angle information from factor matrices. Interestingly, the new tensor model can be regarded as element-wise product of two tensors of different ranks. The first one presumes the same multi-linear structure as that of in the conventional MIMO radar, and the other one contains all phase modulation values for DDMA technique. This enables the proposed method to take advantages of tensor-based algorithms while maintaining the number of samples in Doppler domain. Angle estimation performance for slow-time MIMO radar can hence be significantly improved, which is verified by simulation results.
Ii Slow-Time MIMO Radar Signal Model
Consider a bistatic MIMO radar with collocated transmit and
collocated receive antenna elements. Both transmit and receive arrays are uniform linear arrays (ULAs) whose element spacing are half the working wavelength. The steering vectors of the transmit and receive arrays are then denoted byand , where and are the DOD and DOA, respectively, and denotes the transpose of a matrix/vector. Assuming that there are totally targets in a range cell of interest, the transmit and receive steering matrices are given, respectively, as and .
The matrix of transmit waveforms is denoted by where is the number of snapshots per pulse. To achieve waveform diversity in slow-time MIMO radar, a phase modulation matrix is used at the transmitter during single coherent processing interval (CPI) with pulses. The waveform envelopes at all transmit elements are identical. Typically, a linear frequency modulated (LFM) signal is used. In -th pulse, , the transmitted signal after applying DDMA technique is . According to [16, 17, 1, 13], the phase modulation matrix from pulse to pulse is
where is the radar pulse duration and is the pulse repetition frequency (PRF). For single scatterer at location with Doppler frequency and complex value (defined as radar cross section (RCS) fading coefficient), the received signal of -th pulse at the output of -th receive element can be written as
where is the -th element of and
is the Gaussian white noise. Then the received signal after pulse compression, range gating, and lowpass filtering in slow-time MIMO radar can be expressed as
where and are the -th and -th elements of the corresponding vectors, , is assumed to be an integer, , and is the noise residue after processing. See Appendix for details.
Therefore, all received signals from targets can be collected into the following 3-order tensor of size :
where denotes the outer product and is the noise tensor.
Iii Joint DOD and DOA Estimation for Slow-Time MIMO Radar
Iii-a Conventional Methods
Estimators of and based on signal subspace algorithms have been conventionally conducted on a per-pulse basis. Using these methods, results can be updated from pulse to pulse. Specifically, we can arrive to the conventional signal model just from the mode- unfolding (frontal slices) of (4), given by [20, 14]
where , , , stands for the Hadamard product, denotes the Khatri-Rao product (column-wise Kronecker product), and is the matricized form of of dimension . Inspecting a single pulse of (5), e.g., the -th pulse, the received signal is
which coincides with the signal model used in the conventional signal subspace-based angle estimation algorithms.
Let us reshape (6) into the following matrix:
where , , and denotes the operator that builds a diagonal matrix from a column vector. Model (7) is identical to that in  when Doppler effect is added, except for the reduced number of pulses. For signal subspace-based algorithms, this difference has slight influence. However, it may cause serious performance degradation for tensor-based algorithms.
Iii-B Modified Tensor Decomposition-Based Joint DOD and DOA Estimation in Slow-Time MIMO Radar
To overcome the performance loss caused by the reduced number of pulses, a new tensor model for slow-time MIMO radar is designed.
Recall (24), it can be regarded as an times downsampling sequence after lowpass filtering in Doppler domain applied in order to avoid Doppler ambiguity. In  and , it is shown that PARAFAC decomposition can often be computed by means of a simultaneous matrix decomposition when the tensor is tall in one mode. If the condition is satisfied in a 3-order tensor , the convergence of the ALS algorithm can be improved by applying the singular value decomposition (SVD) of mode-3 matrix . In radar, is a common case. However, the number of efficient pulses for each transmit and receive channel in slow-time MIMO radar is reduced from to that challenges the condition. In the following, we show that the number of pulses can be restored to
Specifically, let the received signal for a single receive element in
pulses be collected. Pulse compression is applied in fast time for each pulse and fast Fourier transform (FFT) is utilized in slow-time.111An example of the output of these operation will be shown in Fig. 3. Phase demodulation and lowpass filtering in Doppler domain are used to separate received signal for each transmit and receive channel. By exploiting times interpolation column-wisely and applying IFFT, the received signal can be reformulated as . Then, is applied again to shift each component with a unique Doppler frequency, which finally gives . It is worth noting that the loss of Doppler information with relatively high frequency means the Doppler frequency in after interpolation is different with in . However, due to the independence of spatial information and Doppler information, this interpolation has no influence on angle estimation.
Therefore, the reformulation of (3) for targets with samples can be approximately written as
where is the -th element of and . Note that there is no index in , meaning that is repeated from one receive element to another without any changes. This property will be used further to separate the phase modulation component. Expression (8) is approximate because of the above mentioned interpolation, which however has no influence on angle estimation due to the independence of spatial and Doppler information.
First, reshape the received signal of the -th pulse in (8) into matrix form as
where . Note that where . The vectorized form of (9) is
where , contains the targets RCS and Doppler information, and is all-one vector.
Let denote the received signal for pulses. It can be expressed as
where , , , and denotes identity matrix. Comparing (11) with (5)–(7), it is important to stress that the matrix (from the DDMA technique) is applied on each transmit-receive channel at every pulse. Moreover, the matrices and in (11) can be regarded as the model-3 unfoldings [20, 14] of the following two tensors that have different ranks
where is identity tensor and symbol , stands for the mode- product of a tensor and a matrix .
The important observation from (12) is that the new signal tensor model for slow-time MIMO radar can be regarded as the Hadamard product of two tensors. One of them is identical to the conventional MIMO radar tensor model, while the other one stands for the phase modulation values applied on the transmit elements. It is also worth noting that the mode-2 unfolding (lateral slices) of are identical to , which exactly explains the essence of DDMA technique.
Consequently, , can be estimated by fully exploiting the Vandermonde structure of the factor matrices . Take for example. Then two subarrays, one without the last element and the other without the first element of the transmit array, can be formed. Using (12), the received signals for these two subarrays can be expressed as
where , with and standing for submatrices of without the last and first row, respectively. Similarly, and are submatrices of without the last and first row, respectively. Since and are Vandermonde matrices, , , where , , and .
Let where stands for the concatenation of two tensors along the -th mode. Particularly, it can be
where and . Note that is fixed. Therefore, the tensor can be used to conduct PARAFAC decomposition via the following modified ALS algorithm
where , and are the mode- unfoldings of the tensors and , respectively, denoted the Frobenius norm of a matrix, and stands for an estimate of . When are fixed, the objective function in (15) is quadratic in . This property remains while the optimization parameter alternates between , and . Thus, at each alternating step, the objective similar to the one in (15) is quadratic with respect to the optimized matrix parameter, and the corresponding PARAFAC decomposition of can be found.
Finally, the matrices and can be extracted from for which the property
should hold. Since has full rank, the least squares method can be used to estimate it, that is, where
stands for the pseudo-inverse of a rectangular matrix. Using eigenvalue decomposition of, we find the eigenvalues representing the estimates of the diagonal elements of . These estimates are then used to compute , e.g., . Here represents natural logarithm.
The parameters can be estimated similarly using matrix instead of in (13)–(16). Specifically, two subarrays without the last and first elements of the receive array are applied, where the receive steering matrices , are obtained from by removing the last and first rows, respectively. The 2-th mode concatenation of two tensors from (12) for subarrays at the receive side is now reformulated as
where is fixed since it is independent on an index of a receive element, and is the concatenated noise residue. Using the modified ALS algorithm above, the second factor matrix can be decomposed. Note that , where . Then the diagonal elements of are estimated by computing the eigenvalues of matrix , and are estimated as .
which can also be written as . If , which is a common case in radar, the maximum number of targets that can be resolved is almost surely .
Iii-C Computational Complexity Analysis
The proposed algorithm for joint DOD and DOA estimation in slow-time MIMO radar requires PARAFAC decomposition and ESPRIT-aided method to compute phase rotations of targets. During each iteration of the ALS algorithm, the number of flops is . The number of flops for matrix computation, i.e., for computing and is . In total, the number of flops needed in our algorithm is , where is the number of iterations. Note that the proposed algorithm requires to perform ALS two times, and it is useful to apply the modified ALS algorithm with better convergence. We refer to [19, 20] and the references therein for more details. Table I summarizes the computational complexity of the proposed algorithm as well as the complexity of most the existing MIMO techniques.
Iv Simulation Results
We demonstrate here the angle estimation performance of the proposed method in comparison to conventional algorithms including PM , ESPRIT , U-ESPRIT, and PARAFAC. The Cramer-Rao lower bound (CRLB) for bistatic MIMO radar is also provided . Throughout our simulations, a slow-time MIMO radar with antenna elements is considered. A chirp signal with MHz and us is used as a waveform envelope . In each CPI, there are pulses with KHz. We assume targets that follow Swerling I model , and thus,
is chosen from a Gaussian distribution as a complex value, which remains fixed from pulse to pulse. The normalized Doppler frequencies of the targets are. The number of Monte Carlo trials is .
In our first example, the DODs and DOAs of the targets are and , respectively. The root mean square errors (RMSEs) of and are computed separately and then combined. As can be seen in Fig. 1, where RMSE is shown versus SNR, the proposed method achieves better performance, especially, at low SNR. Thus, higher estimation accuracy can be achieved by the proposed method. The PM, ESPRIT, and U-ESPRIT methods can be regarded as generalized signal subspace-based approaches, which are sensitive to low SNR. The PARAFAC method exploits the multilinear structure of the received data and avoids degradation at low SNR, but the number of pulses in this method is reduced. By our method, we recover the number of pulses by times.
In our second example, the probability of resolution of two closely spaced targets with and is analyzed. Two targets are considered to be resolved when and . It can be seen from Fig. 2 that all methods exhibit resolution with probability 1 at high SNR, but the proposed method surpasses other methods as it has the lowest SNR threshold. The improved ability of resolving two closely spaced targets can be regarded as the advantage resulted from combining PARAFAC and ESPRIT. Moreover, the concatenation of tensors in (14) for different subarrays approximately doubles the number of samples. Thus, the proposed method achieves higher accuracy and better resolution for angle estimation.
A new tensor model for slow-time MIMO radar that enables improved joint DOD and DOA estimation for multiple targets has been proposed. This tensor can be regarded as an element-wise product of two tensors, where only one of them contains the angular parameters of interest. The model enables us to use PARAFAC decomposition with ESPRIT, and to address the problem of shortage of samples in Doppler domain for slow-time MIMO radar. As a result, the angle estimation performance has been improved as compared to the existing techniques.
Applying matched-filtering to in fast time, we have
where and is noise residue after matched-filtering. Since is an LFM signal, the integral term is known to be approximately a function, and its peak indicates the range cell of the target.
Then, the concatenation of in a single CPI with pulses forms a matrix, i.e., . If fast Fourier transform (FFT) is applied to this matrix column-wisely, peaks with different Doppler frequencies can be found at the slice of target range cell. Each of the peaks corresponds to a unique transmit element as shown in the range-Doppler map shown in Fig. 3. The distance between two adjacent peaks in Doppler domain is determined by . Owing to this Doppler frequency shifts for different transmitted waveforms, it is possible to distinguish each transmit channel at the receiver via filtering in Doppler domain.
For any -th transmit element, the phase demodulation in Doppler domain is . Equivalently
where , , and is the noise residue. After demodulation, each of transmit channels is shifted to baseband in Doppler domain (with a frequency of ). A reduced efficient Doppler range of is distributed to every transmit element (see Fig. 4 for more details). By applying lowpass filtering with cutoff frequency , the second term in (21) can be omitted. Hence, the received signal from -th transmit element to -th receive element at -th pulse can be expressed as
By range gating, the received signal is further expressed as
where is the slice of at the target range cell. Note that in order to avoid ambiguous Doppler returns, it is necessary to ensure that the highest Doppler frequency of interest is smaller than . This implies that the DDMA technique achieves waveform diversity for slow-time MIMO radar at the cost of reduced Doppler frequency estimation range.
Another disadvantage is the decrease of the number of samples in Doppler domain. Recall , we have . Clearly, this is a discrete signal with sampling rate . The uniformly divided Doppler space then leads to the decline of sampling rate to , i.e., . Therefore,
where the number of efficient pulses is reduced to only . Here is assumed to be an integer. Considering the noise term, the result in (3) is obtained.
-  (2009-Feb.) GMTI MIMO radar. In Proc. Int. Waveform Diversity Des. Conf., Vol. , pp. 118–122. External Links: Cited by: §I, §II.
-  (2017-Oct.) Transmit array interpolation for DOA estimation via tensor decomposition in 2-D MIMO radar. ieee_j_sp 65 (19), pp. 5225–5239. External Links: Cited by: §I.
-  (2010-05) On the concept of MIMO radar. In Proc. IEEE Radar Conf., Vol. , pp. 327–332. External Links: Cited by: §I.
-  (2000) A multilinear singular value decomposition. SIAM J. Matrix Anal. Appl. 21 (4), pp. 1253–1278. Cited by: §I.
-  (2008) Decompositions of a higher-order tensor in block terms—part III: alternating least squares algorithms. SIAM Rev. 30 (3), pp. 1067–1083. Cited by: §III-B.
-  (2006) A link between the canonical decomposition in multilinear algebra and simultaneous matrix diagonalization. SIAM Rev. 28 (3), pp. 642–666. Cited by: §III-B.
-  (2008-Jun.) Angle estimation using ESPRIT in MIMO radar. Electron. Lett. 44 (12), pp. 770–771. External Links: Cited by: §I, TABLE I, §IV.
-  (2004-Apr.) MIMO radar: an idea whose time has come. In Proc. IEEE Radar Conf., Vol. , pp. 71–78. External Links: Cited by: §I.
-  (2004-Nov.) Performance of MIMO radar systems: advantages of angular diversity. In Proc. 38th Asilomar Conf. Signals, Syst. Comput., Pacific Grove, CA, pp. 305–309. External Links: Cited by: §I.
-  (2008-Jan.) MIMO radar with widely separated antennas. ieee_m_sp 25 (1), pp. 116–129. External Links: Cited by: §I.
-  (2010-Jun.) Phased-MIMO radar: a tradeoff between phased-array and MIMO radars. ieee_j_sp 58 (6), pp. 3137–3151. External Links: Cited by: §I.
-  (2011-Jun.) Transmit energy focusing for DOA estimation in MIMO radar with colocated antennas. ieee_j_sp 59 (6), pp. 2669–2682. External Links: Cited by: §I.
-  (2010-Nov.) Clutter covariance matrices for GMTI MIMO radar. In Proc. 44th Asilomar Conf. Signals, Syst., Comput., Vol. , pp. 1821–1826. External Links: Cited by: §I, §II.
-  (2009) Tensor decompositions and applications. SIAM Rev. 51 (3), pp. 455–500. Cited by: §I, §III-A, §III-B, §III-B.
-  (2007-Sep.) MIMO radar with colocated antennas. ieee_m_sp 24 (5), pp. 106–114. External Links: Cited by: §I.
-  (2009) MIMO radar signal processing. Wiley, New York. Cited by: §I, §II.
-  (2006-Jul.) MIMO radar space-time adaptive processing for multipath clutter mitigation. In Proc. IEEE Workshop Sensor Array Multichannel Signal Processing, Vol. , pp. 249–253. External Links: Cited by: §I, §II.
-  (2010-Nov.) Tensor algebra and multidimensional harmonic retrieval in signal processing for MIMO radar. ieee_j_sp 58 (11), pp. 5693–5705. External Links: Cited by: §I, §III-A, §IV.
-  (2000-Aug.) Parallel factor analysis in sensor array processing. ieee_j_sp 48 (8), pp. 2377–2388. External Links: Cited by: §I, §III-C.
Tensor decomposition for signal processing and machine learning. ieee_j_sp 65 (13), pp. 3551–3582. External Links: Cited by: §I, §III-A, §III-B, §III-B, §III-C.
-  (1990) Radar handbook. Second edition, McGraw-Hil, New York. Cited by: §IV.
-  (2013) Optimization-based algorithms for tensor decompositions: canonical polyadic decomposition, decomposition in rank-(,,1) terms, and a new generalization. SIAM Rev. 23 (2), pp. 695–720. Cited by: §III-C.
-  (2017) A novel unitary PARAFAC method for DOD and DOA estimation in bistatic MIMO radar. Signal Processing 138, pp. 273 – 279. External Links: Cited by: §I, TABLE I.
-  (2010-Dec.) Direction of departure (dod) and direction of arrival (doa) estimation in MIMO radar with reduced-dimension MUSIC. ieee_j_coml 14 (12), pp. 1161–1163. Cited by: §I, TABLE I.
-  (2011-Nov.) Trilinear decomposition-based transmit angle and receive angle estimation for multiple-input multiple-output radar. IET Radar, Sonar Navig. 5 (6), pp. 626–631. Cited by: §IV.
-  (2012-Feb.) Unitary ESPRIT algorithm for bistatic MIMO radar. Electron. Lett. 48 (3), pp. 179–181. Cited by: §I, TABLE I, §IV.
-  (2011-Jan.) Fast method for multi-target localisation in bistatic MIMO radar. Electron. Lett. 47 (2), pp. 138–139. Cited by: §I, TABLE I, §IV.