Height Measurement for Meter-Wave MIMO Radar Based on Sparse Array Under Multipath Interference

Abstract

For meter-wave multiple-input multiple-output (MIMO) radar, the multipath of target echoes may cause severe errors in height measurement, especially in the case of complex terrain where terrain fluctuation, ground inclination, and multiple reflection points exist. Inspired by a sparse array with greater degrees of freedom and low mutual coupling, a height measurement method based on a sparse array is proposed. First, a practical signal model of MIMO radar based on a sparse array is established. Then, the modified multiple signal classification (MUSIC) and maximum likelihood (ML) estimation algorithms based on two classical sparse arrays (coprime array and nested array) are proposed. To reduce the complexity of the algorithm, a real-valued processing algorithm for generalized MUSIC (GMUSIC) and maximum likelihood is proposed, and a reduced dimension matrix is introduced into the real-valued processing algorithm to further reduce computation complexity. Finally, sufficient simulation results are provided to illustrate the effectiveness and superiority of the proposed technique. The simulation results show that the height measurement accuracy can be efficiently improved by using our proposed technique for both simple and complex terrain.

1. Introduction

Height measurement is a hard problem in low-elevation regions due to the reflected multipath signals [1,2,3,4]. The problem can be regarded as how to decorrelate the coherent signals, and then develop many corresponding height measurement methods, such as the spatial smoothing technique [5,6], reconstructing matrix method [7,8], maximum likelihood estimation [9], alternation projection technique [10], compression sensing and sparse reconstruction techniques [11,12,13,14,15,16,17,18]. In a practical topographic environment, multiple reflection paths exist simultaneously. Over the past few decades, a considerable amount of research studies [19,20,21,22,23,24] have been developed to solve this kind of problem. Based on the classical two-rays propagation model, the literature [19] discusses the improvement of low-angle radar by using a sampled aperture radar and maximum likelihood algorithm. A neural network model is proposed in the literature [20], which enhances the phase characteristics of direct signals, reduces the phase distortion caused by multipath signals, and effectively improves the accuracy of the direction of arrival (DOA) estimation. Because the received target echo is seriously damaged by multipath reflection in complex terrain, the above methods may fail. Therefore, a target localization algorithm based on iterative implementation is proposed in the literature [21] to reduce the adverse effects of convolutional coherent interference caused by multipath propagation on target location parameter estimation. The authors of [22] establish a perturbation multipath propagation model for low-angle tracking and localization, and they propose a novel method to resolve the model and estimate the DOA of low-elevation targets. The authors of [23] regard the influence of irregular reflection in complex terrain as the random perturbation of a multipath signal, and they propose a weighted sparse Bayesian learning method. The authors of [24] consider an actual multipath propagation scenario and propose a target localization algorithm based on low-rank decomposition and the sparse representation framework.

It is worth noting that MIMO radar [25,26,27] has shown significant advantages in the detection of low-elevation targets for its waveform diversity, and there have been considerable research efforts [28,29,30,31,32,33,34,35,36]. The authors of [29,30] studied the height measurement of meter-wave polarimetric MIMO radar based on the characteristics of polarization diversity and waveform diversity. The authors of [31] studied the angle estimation of low-elevation targets in meter-wave FDA-MIMO radar by combining the advantages of frequency diversity and waveform diversity. In literature [32], aiming at the beam-splitting phenomenon caused by low-angle multipath reflection, reflected wave and direct wave coupling, a meter-wave MIMO radar virtual subarray-level beam-splitting low-elevation height measurement method was studied. However, in the actual situation, there are undulating ground and multiple reflection paths [37], which cause the above algorithms to fail. The authors of [33] propose a low-altitude measurement method of MIMO radar based on adaptive beamforming, which is suitable for undulating terrain and inclined ground. Reference [38] studies a height measurement algorithm based on a matrix pencil, which is suitable for complex terrain. In addition, references [39,40] provide new insights into MIMO radar systems and have great potential applications in low elevation fields.

The above-mentioned algorithms of low-angle target localization all adopt uniform linear arrays, which have the problems of low estimation accuracy and a large fluctuation of measurement error with an elevation angle change. Considering sparse arrays [41,42] with greater degrees of freedom and low mutual coupling, its equivalent aperture is larger than a uniform linear array, so the performance of direction of arrival estimation will be better. In addition, for a target with multipath reflection, the phase difference caused by the difference between the direct wave and the reflected wave will be different due to the non-uniformity of the array spacing. To sum up, the application of a sparse array in height measurement in a low elevation region, especially in complex terrain, has great advantages.

In this paper, a height measurement method based on a sparse array is proposed for the first time, which is suitable for both flat terrain and complex terrain. In complex terrain, the sparse array is introduced into the MIMO radar system as a transceiver antenna. In flat terrain, the sparse array is introduced into the MIMO radar system as transmitting and receiving antenna, respectively. Firstly, a multipath signal model based on a sparse array is established. Secondly, combined with the GMUSIC algorithm and ML algorithm, a height measurement method based on a coprime array and nested array is proposed. Furthermore, in simple terrain, to reduce the computation cost, a real-valued algorithm and a reduced dimension real-valued algorithm are proposed. Thirdly, the complexity of the three algorithms is analyzed. Finally, the simulations and some comparative experiments are conducted to verify the effectiveness and superiority of the proposed algorithm.

The contributions of our paper can be summarized as follows. The main innovation of this paper is to establish a multipath signal model based on a sparse array, which is different from the previous signal models. Our model is not only closer to reality, but it also has a larger array aperture. On this basis, the real-valued processing algorithm and reduced dimension real-valued processing algorithm are proposed. The proposed algorithm is not only suitable for simple flat terrain but also for complex terrain.

Notations: Superscript ${()}^T$ and ${()}^H$ denote transpose and conjugate transpose, respectively. $\otimes$ denotes the Kronecker product. $I_M$ is an identity matrix with dimension $M$. $\Pi$ is the exchange matrix with ones on its anti-diagonal and zeros elsewhere. $vec(\cdot )$ denotes vectorization processing. $diag(\cdot )$ denotes diagonalization. $\det (\cdot )$ denotes the matrix determinant. $E[\cdot ]$ denotes expectation.

2. Signal Model

Consider a narrowband monostatic MIMO radar system with M-element sensors, which adopts a coprime array or a nested array as transceiver antennas. The structure of the coprime array is shown in Figure 1, and the structure of the nested array is shown in Figure 2. The coprime array consists of two sparse uniform subarrays, the number of elements is $N_1$ and $N_2$ ($N_1<N_2$), the spacing of elements is $N_{2}d$ and $N_{1}d$, $N_1$ and $N_2$ is a pair of coprime integers, and $d$ is half of the wavelength. By a linear superposition combination of two sparse uniform subarrays in the form of first element overlapping, a coprime array containing $N_1+N_2-1$ physical elements is obtained, and the position of the last element is $(N_2-1)N_{1}d$. The nested array is composed of two uniform linear arrays, the total elements of array is $N=N_1+N_2$, the first-level subarray has $N_1$ elements with spacing $d$, and the second-level subarray has $N_2$ elements with spacing $d_2=(N_1+1)d$. Figure 3 shows a multipath reflection signal model based on a sparse array, in which the array antenna is sparse array, $h_a$ is the reference element height (set as the lowest element height), $h_t$ is the target height, ${\theta }_{q,i}$ is the ground inclination angle of the $i\mathrm{th}$ reflection point, ${\theta }_{s,i}$ is the direction of the reflection wave of the $i\mathrm{th}$ reflection point, and B and C stand for the reflection points of the two reflected paths.

Suppose that the transmitted signal of the $m\mathrm{th}$ array element at time $t$ is

$$g_m(t)=s_m(t)e^{j(2\pi f_{0}t)}$$

(1)

where $f_0$ is the carrier frequency, and $s_m(t)$ is the transmitted signal complex envelope of the $m\mathrm{th}$ array element. The transmitted signal of the $m\mathrm{th}$ array element reaching the target through the direct wave is

$$x_{d,m}(t)=s_m(t-{\tau }_{d,m})e^{j2\pi f_0(t-{\tau }_{d,m})}$$

(2)

where ${\tau }_{d,m}=R_{d,m}/c$ is the time delay for the $m\mathrm{th}$ array element of the direct wave, and $R_{d,m}$ is the direct wave distance of the $m\mathrm{th}$ array element. According to the geometric relationship, it can be expressed as

$$R_{d,m}=\sqrt{R_r^2+{(h_t-h_m)}^2}\approx R-h_m\sin {\theta }_d$$

(3)

where $h_m$ is the height of the $m\mathrm{th}$ array element, $R$ is the distance between the target and the radar, and $R_r$ is the horizontal distance between the radar and the target.

Similarly, the transmitted signal of the $m\mathrm{th}$ array element reaching the target through the reflected wave is

$$x_{s,m}(t)={\displaystyle \sum _{i=1}^I{\rho }_{s,i}{s}_m(t-{\tau }_{si,m})e^{j2\pi f_0(t-{\tau }_{si,m})}}$$

(4)

where $I$ is the number of reflected wave paths, and ${\rho }_{si}$ is the reflection coefficient of the $i\mathrm{th}$ reflection wave. ${\tau }_{si,m}=R_{si,m}/c$ is the time delay for the $m\mathrm{th}$ array element of the reflected wave, and $R_{si,m}$ is the reflected wave distance of the $m\mathrm{th}$ array element. Thus, the transmitted signal of the $m\mathrm{th}$ array element reaching the target under the multipath reflection condition is

$$x_m(t)=x_{d,m}(t)+x_{s,m}(t)=s_m(t-{\tau }_{d,m})e^{j2\pi f_0(t-{\tau }_{d,m})}+{\displaystyle \sum _{i=1}^I{\rho }_{s,i}{s}_m(t-{\tau }_{si,m})e^{j2\pi f_0(t-{\tau }_{si,m})}}$$

(5)

Considering the transmitted signal of MIMO radar is a narrowband signal, for a narrowband signal with carrier frequency $f$, $s_m(t+\tau )=s_m(t)\exp (j2\pi f\tau )$. By choosing a small $\tau$, $s_m(t+\tau )=s_m(t)\exp (j2\pi f\tau )\approx s_m(t)$; that is, the time delay will not affect the signal complex envelope. Thus,

$$s_m(t-{\tau }_{d,m})\approx s_m(t-{\tau }_{si,m})=s_m(t)$$

(6)

Using (6) in (5), one obtains

$$\begin{array}{ll}{x}_m(t)& =s_m(t)e^{j2\pi f_0(t-{\tau }_{d,m})}+{\displaystyle \sum _{i=1}^I{\rho }_{s,i}{s}_m(t)e^{j2\pi f_0(t-{\tau }_{si,m})}}\\ & =s_m(t)e^{j2\pi f_0(t-(R-h_a\sin {\theta }_d)/c)}\ast \left(e^{j2\pi f_0{d}_m\sin {\theta }_d/c}+{\displaystyle \sum _{i=1}^I{\rho }_{s,i}{e}^{-j2\pi f_0\Delta R_i/c}{e}^{-j2\pi f_0{d}_m\sin {\theta }_d/c}}\right)\end{array}$$

(7)

where we have exploited the fact that $h_m=h_a+d_m$ and $R_{si,m}=R_{d,m}+\Delta R_i$.

Using $\lambda =c/f$ into (7), one obtains

$$x_m(t)=s_m(t)\stackrel{-}{x_0}(t)\ast \left(e^{j2\pi d_m\sin {\theta }_d/\lambda }+{\displaystyle \sum _{i=1}^I{\rho }_{s,i}{e}^{-j2\pi \Delta R_i/\lambda }{e}^{j2\pi d_m\sin {\theta }_{s,i}/\lambda }}\right)$$

(8)

where $\overline{x_0}(t)=e^{j2\pi f_0(t-(R-h_a\sin {\theta }_d)/c)}$. Therefore, the transmitted signal reaching the target in the multipath condition can be expressed as

$$x(t)={\left[a({\theta }_d)+{\displaystyle \sum _{i=1}^I{\mathrm{e}}^{-\mathrm{j}{\delta }_i}{\rho }_{s,i}a({\theta }_{s,i})}\right]}^{T}s(t)\overline{x_0}(t)$$

(9)

where ${\delta }_i=2\pi \Delta R_i/\lambda$, $s(t)={[s_1(t),\dots ,s_M(t)]}^T$ is the transmitted signal vector of the MIMO radar. $a({\theta }_d),a({\theta }_{s,i})$ are the guiding vectors for the direct and reflected wave, respectively, which can be expressed as follows.

$$\begin{array}{l}a({\theta }_d)={\left[e^{j2\pi d_1\sin {\theta }_d/\lambda },\dots ,e^{j2\pi d_m\sin {\theta }_d/\lambda },\dots ,e^{j2\pi d_M\sin {\theta }_d/\lambda }\right]}^T\\ a({\theta }_{s,i})={\left[e^{j2\pi d_1\sin {\theta }_{s,i}/\lambda },\dots ,e^{j2\pi d_m\sin {\theta }_{s,i}/\lambda },\dots ,e^{j2\pi d_M\sin {\theta }_{s,i}/\lambda }\right]}^T\end{array}$$

(10)

where $d={[d_1,\dots ,d_M]}^T$ is the sparse array sensor location. Considering the multipath effect for the receiving array, the expression of the data received by the receiving array is

$$X(t)=\left[a({\theta }_d)+{\displaystyle \sum _{i=1}^I{\mathrm{e}}^{-\mathrm{j}{\delta }_i}{\rho }_{s,i}a({\theta }_{s,i})}\right]\beta (t)\ast {\left[a({\theta }_d)+{\displaystyle \sum _{i=1}^I{\mathrm{e}}^{-\mathrm{j}{\delta }_i}{\rho }_{s,i}a({\theta }_{s,i})}\right]}^{T}s(t)\overline{x_0}(t)+n(t)$$

(11)

where $\beta (t)$ is the complex reflection coefficient of the target, and $n(t)$ is the Gaussian white noise with zero mean and variance ${\sigma }_n^2$. By matched filtering on Equation (11), we obtain

$$X(t)=\left[a({\theta }_d)+{\displaystyle \sum _{i=1}^I{\mathrm{e}}^{-\mathrm{j}{\delta }_i}{\rho }_{s,i}a({\theta }_{s,i})}\right]\beta (t)\ast {\left[a({\theta }_d)+{\displaystyle \sum _{i=1}^I{\mathrm{e}}^{-\mathrm{j}{\delta }_i}{\rho }_{s,i}a({\theta }_{s,i})}\right]}^T\stackrel{-}{x_0}(t)+N(t)$$

(12)

After vectorization of the data, it can be expressed as

$$Y=\left[a({\theta }_d)+{\displaystyle \sum _{i=1}^I{\mathrm{e}}^{-\mathrm{j}{\delta }_i}{\rho }_{s,i}a({\theta }_{s,i})}\right]\otimes \left[a({\theta }_d)+{\displaystyle \sum _{i=1}^I{\mathrm{e}}^{-\mathrm{j}{\delta }_i}{\rho }_{s,i}a({\theta }_{s,i})}\right]\beta (t)\overline{x_0}(t)+vec(N(t))$$

(13)

3. Height Measurement Method Based on Coprime Array and Nested Array

The gap between the signal model of MIMO radar and the signal model of conventional array radar mainly lies the relationship between the direct wave and reflected wave. In conventional array radar, the direct wave and reflected wave can be equivalent to a coherent signal. While in MIMO radar, especially the sparse array adopted in transmitting and receiving antennas, traditional decoherent algorithms such as the spatial smoothing and matrix reconstruction are not suitable. As such, the generalized MUSIC (GMUSIC) or maximum likelihood (ML) estimation algorithm is studied to obtain a low-elevation angle, and three algorithms based on GMUSIC and the maximum likelihood estimation algorithm for the height measurement of meter-wave MIMO radar are proposed as follows.

3.1. Basic Algorithm

The conventional MUSIC algorithm uses the orthogonality of signal space and noise space to search spectral peaks. Inspired by the conventional MUSIC algorithm, the GMUSIC algorithm and ML algorithm are proposed to match the complex terrain. Furthermore, the sparse array-based meter-wave MIMO radar low-altitude target height measurement method is studied to promote the accuracy of the low elevation angle in case of complex terrain.

To simplify the analysis without losing generality, we assume that there are two reflection paths, and Equation (13) can be simplified as

$$Y=A({\theta }_d,{\theta }_{s,1},{\theta }_{s,2})\Gamma \beta (t)\overline{x_0}(t)+vec(N(t))$$

(14)

$$A({\theta }_d,{\theta }_{s,1},{\theta }_{s,2})=[a({\theta }_d)a({\theta }_{s,1})a({\theta }_{s,2})]\otimes [a({\theta }_d)a({\theta }_{s,1})a({\theta }_{s,2})]$$

(15)

$$\Gamma =[1,{\mathrm{e}}^{-\mathrm{j}{\delta }_1}{\rho }_{s,1}{,\mathrm{e}}^{-\mathrm{j}{\delta }_2}{\rho }_{s,2}{,\mathrm{e}}^{-\mathrm{j}{\delta }_1}{\rho }_{s,1},{\mathrm{e}}^{-\mathrm{j}2{\delta }_1}{\rho }_{s,1}^2,{\mathrm{e}}^{-\mathrm{j}({\delta }_1+{\delta }_2)}{\rho }_{s,1}{\rho }_{s,2}{,\mathrm{e}}^{-\mathrm{j}{\delta }_2}{\rho }_{s,2},{\mathrm{e}}^{-\mathrm{j}({\delta }_1+{\delta }_2)}{\rho }_{s,1}{\rho }_{s,2},{\mathrm{e}}^{-\mathrm{j}2{\delta }_2}{\rho }_{s,2}^2]$$

(16)

where $A({\theta }_d,{\theta }_{s,1},{\theta }_{s,2})$ is the compound guiding vector of the received data and $a({\theta }_d)$, $a({\theta }_{s,1})$, and $a({\theta }_{s,2})$ are guiding vectors based on a sparse array. For example, if the two-level nested array is adopted, and the first-level subarray has four elements, the second-level subarray has three elements, and the two subarrays share the first element; then, the total number of nested array is six. The nested array sensor location and steering vector $a({\theta }_d)$, $a({\theta }_{s,1})$, $a({\theta }_{s,2})$ are, respectively, as follows.

$$d={[d1,d2,d3,d4,d5,d6]}^T={[1d,2d,3d,4d,8d,12d]}^T$$

(17)

$$\begin{array}{ll}a({\theta }_d)& ={[e^{j2\pi d_1\sin {\theta }_d/\lambda },e^{j2\pi d_2\sin {\theta }_d/\lambda },e^{j2\pi d_3\sin {\theta }_d/\lambda },e^{j2\pi d_4\sin {\theta }_d/\lambda },e^{j2\pi d_5\sin {\theta }_d/\lambda },e^{j2\pi d_6\sin {\theta }_d/\lambda }]}^T\\ & ={[e^{j2\pi d\sin {\theta }_d/\lambda },e^{j2\pi 2d\sin {\theta }_d/\lambda },e^{j2\pi 3d\sin {\theta }_d/\lambda },e^{j2\pi 4d\sin {\theta }_d/\lambda },e^{j2\pi 8d\sin {\theta }_d/\lambda },e^{j2\pi 12d\sin {\theta }_d/\lambda }]}^T\end{array}$$

(18)

$$\begin{array}{ll}a({\theta }_d)& ={[e^{j2\pi d_1\sin {\theta }_d/\lambda },e^{j2\pi d_2\sin {\theta }_d/\lambda },e^{j2\pi d_3\sin {\theta }_d/\lambda },e^{j2\pi d_4\sin {\theta }_d/\lambda },e^{j2\pi d_5\sin {\theta }_d/\lambda },e^{j2\pi d_6\sin {\theta }_d/\lambda }]}^T\\ & ={[e^{j2\pi d\sin {\theta }_d/\lambda },e^{j2\pi 2d\sin {\theta }_d/\lambda },e^{j2\pi 3d\sin {\theta }_d/\lambda },e^{j2\pi 4d\sin {\theta }_d/\lambda },e^{j2\pi 8d\sin {\theta }_d/\lambda },e^{j2\pi 12d\sin {\theta }_d/\lambda }]}^T\end{array}$$

(19)

$$\begin{array}{ll}a({\theta }_{s,2})& ={[e^{j2\pi d_1\sin {\theta }_{s,2}/\lambda },e^{j2\pi d_2\sin {\theta }_{s,2}/\lambda },e^{j2\pi d_3\sin {\theta }_{s,2}/\lambda },e^{j2\pi d_4\sin {\theta }_{s,2}/\lambda },e^{j2\pi d_5\sin {\theta }_{s,2}/\lambda },e^{j2\pi d_6\sin {\theta }_{s,2}/\lambda }]}^T\\ & ={[e^{j2\pi d\sin {\theta }_{s,2}/\lambda },e^{j2\pi 2d\sin {\theta }_{s,2}/\lambda },e^{j2\pi 3d\sin {\theta }_{s,2}/\lambda },e^{j2\pi 4d\sin {\theta }_{s,2}/\lambda },e^{j2\pi 8d\sin {\theta }_{s,2}/\lambda },e^{j2\pi 12d\sin {\theta }_{s,2}/\lambda }]}^T\end{array}$$

(20)

Note that $A({\theta }_d,{\theta }_{s,1},{\theta }_{s,2})$ contains three unknown parameters. ${\theta }_{s,2}$, ${\theta }_{s,1}$ and ${\theta }_d$ have a relationship according to the geometric relationship:

$${\theta }_{s,i}\approx (\arctan (\tan {\theta }_d+2h_a/R)+{\theta }_{q,i})$$

(21)

Furthermore, in the actual terrain environment, the inclination angle can be known by topographic reconnaissance in advance. Thus, $A({\theta }_d,{\theta }_{s,1},{\theta }_{s,2})$ can be written as $A({\theta }_d)$, and $A({\theta }_d)$ is marked with $A(\theta )$.

According to the maximum likelihood estimation criterion, the received data covariance matrix can be obtained from the following equation

$$\widehat{R}={\displaystyle \frac{1}{L}}{\displaystyle \sum _{l=1}^{L}y(l)y^H(l)}$$

(22)

where $L$ represents the number of snapshots; that is, $\widehat{R}$ is directly calculated by the snapshots acquired by the actual array. The noise subspace $E_n$ can be obtained by eigenvalue decomposition of the covariance matrix. Different from the steering vector in the conventional MUSIC algorithm, the steering vector $A(\theta )$ is orthogonal to the noise subspace. Therefore, the spectral peak search formula of the generalized MUSIC algorithm is

$$f_{GMUSIC}(\theta )={\displaystyle \frac{\det [A^H(\theta )A(\theta )]}{\det [A^H(\theta )E_n^H{E}_{n}A(\theta )]}}$$

(23)

Correspondingly, the spatial projection matrix of the ML algorithm is constructed by the steering vector matrix proposed in this paper as follows.

$$P_a(\theta )=A(\theta ){(A^H(\theta )A(\theta ))}^{-1}{A}^H(\theta )$$

(24)

Therefore, the spectral peak search formula of ML algorithm is

$$f_{ML}(\theta )={\displaystyle \frac{1}{\det [trace(I_p-P_a(\theta ))\widehat{R}]}}$$

(25)

3.2. Real-Valued Processing Algorithm

Considering MIMO radar enhances the system performance while greatly increasing the amount of computation, a real-valued processing height measurement method is studied to reduce the computation complexity, and they are marked with the UGMUSIC and UML algorithms. It is known that the unitary matrix can transfer the Centro-Hermitian matrix into a real matrix by unitary transformation. The received data covariance matrix $\widehat{R}$ is not a central Hermitian matrix. The covariance matrix is $\widehat{R}$ transformed into the central Hermitian matrix ${\widehat{R}}_{fb}$ by forward and backward spatial smoothing, and the real valued covariance matrix ${\widehat{R}}_U$ can be obtained by unitary transformation as follows.

$${\widehat{R}}_{fb}={\displaystyle \frac{1}{2}}(\widehat{R}+{\Pi }^H{\widehat{R}}^{*}\Pi )$$

(26)

$${\widehat{R}}_U=U^H{\widehat{R}}_{fb}U$$

(27)

where the unitary matrix $\mathbf{U}$ is defined as follows.

$${\mathbf{U}}_{2K}={\displaystyle \frac{1}{\sqrt{2}}}\left[\begin{array}{cc}{{\rm I}}_K& j{{\rm I}}_K\\ {\Pi }_K& -j{\Pi }_K\end{array}\right]$$

(28)

If $P$ is an even number and $K=P/2$, and $P$ is the number of array elements, or

$$U_{2K+1}={\displaystyle \frac{1}{\sqrt{2}}}\left[\begin{array}{ccc}{{\rm I}}_K& 0& j{{\rm I}}_K\\ 0^T& \sqrt{2}& 0^T\\ {\Pi }_K& 0& -j{\Pi }_K\end{array}\right]$$

(29)

$P$ is an odd number and $K=(P-1)/2$. $I_K$ and ${\Pi }_K$ are defined as follows.

$$I_K=\left[\begin{array}{cccc}1& & & \\ & 1& & \\ & & \ddots & \\ & & & 1\end{array}\right]\in R^{K\times K}$$

(30)

$${\Pi }_K=\left[\begin{array}{cccc}& & & 1\\ & & ⋰& \\ & 1& & \\ 1& & & \end{array}\right]\in R^{K\times K}$$

(31)

Similar to Equation (23), the following formula can be obtained

$$f_{UGMUSIC}(\theta )={\displaystyle \frac{\det [{\mathit{A}}_U^H(\theta )A_U(\theta )]}{\det [{\mathit{A}}_U^H(\theta )U_n^H{U}_n{A}_U(\theta )]}}$$

(32)

where $U_n$ is the real-valued noise subspace obtained by the eigenvalue decomposition of real-valued covariance matrix ${\hat{R}}_U$, and $A_U(\theta )$ is a real-valued compound steering vector.

$$A_U(\theta )=U^{H}A(\theta )$$

(33)

Similarly, the UML algorithm search spectrum is as follows.

$$f_{UML}(\theta )={\displaystyle \frac{1}{\det [trace(I_p-P_U(\theta ))\widehat{R_U}]}}$$

(34)

where $P_U(\theta )=A_U(\theta ){({\mathit{A}}_U^H(\theta )A_U(\theta ))}^{-1}{\mathit{A}}_U^H(\theta )$.

3.3. Reduced Dimension Real-Valued Processing Algorithm

Although the real value processing algorithm reduces the computational complexity, it is only one order of magnitude lower than the basic algorithm. To further reduce the computation complexity of the algorithm, the reduced dimension real-valued processing algorithm is proposed.

From the expression of the compound guiding vector, we can obtain

$$\begin{array}{c}A(\theta )=[a({\theta }_d)a({\theta }_{s,1})a({\theta }_{s,2})]\otimes [a({\theta }_d)a({\theta }_{s,1})a({\theta }_{s,2})]\\ =a(\theta )\otimes a(\theta )=Dd(\theta )\end{array}$$

(35)

where $a(\theta )=[a({\theta }_d)a({\theta }_{s,1})a({\theta }_{s,2})]$.

$d(\theta ){=[1,\mathrm{e}}^{{\displaystyle \frac{-j2\pi d\sin {\theta }_d}{\lambda }}},\cdots {,\mathrm{e}}^{{\displaystyle \frac{-j2\pi (2M-2)d\sin {\theta }_d}{\lambda }}}{]}^T$ is the virtual linear array steering vector, and

$$\mathit{D}={\left[\begin{array}{l}\begin{array}{cccccccc}1& 0& \cdots & 0& 0& 0& \cdots & 0\\ 0& 1& \cdots & 0& 0& 0& \cdots & 0\\ ⋮& ⋮& \ddots & ⋮& ⋮& ⋮& \ddots & ⋮\\ 0& 0& \cdots & 1& 0& 0& \cdots & 0\end{array}\}\begin{array}{c}M\end{array}\\ \begin{array}{cccccccc}0& 1& 0& \cdots & 0& 0& \cdots & 0\\ 0& 0& 1& \cdots & 0& 0& \cdots & 0\\ ⋮& ⋮& ⋮& \ddots & ⋮& ⋮& \ddots & ⋮\\ 0& 0& \cdots & \cdots & 1& 0& \cdots & 0\end{array}\}\begin{array}{c}M\end{array}\\ ⋮\\ \begin{array}{cccccccc}0& 0& \cdots & 0& 1& 0& \cdots & 0\\ 0& 0& \cdots & 0& 0& 1& \cdots & 0\\ ⋮& ⋮& \ddots & ⋮& ⋮& ⋮& \ddots & ⋮\\ 0& 0& \cdots & 0& 0& 0& \cdots & 1\end{array}\}\begin{array}{c}M\end{array}\end{array}\}M\right]}_{M^2\times (2M-1)}$$

(36)

Using (35) in (14), one obtains

$$Y=Dd(\theta )\Gamma \beta (t)\overline{x_0}(t)+vec(N(t))$$

(37)

From Equation (37), the received data matrix is a high-dimensional space spanned by a virtual linear array steering vector $d(\theta )$, which can be converted to a low-dimensional space by reduced dimension processing. Let the reduced dimension matrix be the $(2M-1)\ast M^2$ matrix $Q$; then, the received data matrix is

$$Y_{RD}=QDd(\theta )\Gamma \beta (t)\overline{x_0}(t)+N_{RD}$$

(38)

where $N_{RD}=Q\ast vec(N(t))$.

Given that the noise vector after dimension reduction should be Gaussian white noise, it follows that $Q$ satisfies the relation $Q{Q}^H=I_{2M-1}$, and the reduced dimension matrix may be expressed as

$$Q={(D^{H}D)}^{-1/2}{D}^H$$

(39)

Substituting $Q={(D^{H}D)}^{-1/2}{D}^H$ into (37) yields

$$\begin{array}{ll}{Y}_{RD}& ={(D^{H}D)}^{-1/2}{D}^{H}Dd(\theta )\Gamma \beta (t)\overline{x_0}(t)+N_{RD}\\ & ={(D^{H}D)}^{1/2}d(\theta )\Gamma \beta (t)\overline{x_0}(t)+N_{RD}\end{array}$$

(40)

We define $T=D^{H}D$ and $T=diag(1,2,\cdots ,M,\cdots 2,1)$.

The received data covariance matrix after reduced dimension processing is then

$$\begin{array}{ll}{\mathit{R}}_{RD}& =E[{\mathit{Y}}_{RD}{\mathit{Y}}_{RD}^H]\\ & ={\beta }^2(t){\overline{x_0}}^2(t)T^{1/2}d(\theta )\Gamma {\mathit{\Gamma }}^H{\mathit{d}}^H(\theta ){({\mathit{T}}^{1/2})}^H+{\sigma }_n^2{\mathit{I}}_{RD}\end{array}$$

(41)

Notice that after reduced dimension processing, the dimension of the received data covariance matrix becomes $(2M-1)\ast (2M-1)$, and the dimension becomes sharply smaller. Note that the reduced dimension matrix is sparse; thus, the computation increase is not obvious.

The received data covariance matrix after reduced dimension processing can be obtained as follows:

$${\hat{R}}_{RD}={\displaystyle \frac{1}{L}}QY{Y}^H{Q}^H$$

(42)

Note that the received data covariance matrix ${\hat{R}}_{RD}$ is still a complex matrix; thus, its computation complexity can be further reduced by real-valued processing. The reduced dimension real-valued GMUSIC algorithm and reduced dimension real-valued ML algorithm are marked with RD-UGMUSIC and RD-UML.

Thus, the RD-UGMUSIC spectrum formula can be expressed as follows.

$$f_{RD-GMUSIC}(\theta )={\displaystyle \frac{\det [A_{URD}^H(\theta )A_{URD}(\theta )]}{\det [A_{URD}^H(\theta )U_{RDn}^H{U}_{RDn}{A}_{URD}(\theta )]}}$$

(43)

Here, $A_{URD}(\theta )$ is the reduced dimension real valued steering vector, and $U_{RDn}$ is the noise subspace obtained by eigenvalue decomposition of the reduced dimension real-valued received data covariance matrix ${\hat{R}}_{URD}$. The expressions of $A_{URD}({\theta }_d,{\theta }_{s,1},{\theta }_{s,2})$ and ${\hat{R}}_{URD}$ are as follows.

$$A_{URD}({\theta }_d,{\theta }_{s,1},{\theta }_{s,2})=UQA({\theta }_d,{\theta }_{s,1},{\theta }_{s,2})$$

(44)

$${\widehat{R}}_{URD}={\displaystyle \frac{1}{2}}{U}^H({\widehat{R}}_{URD}+{\Pi }^H{\widehat{R}}_{URD}^{*}\Pi )U$$

(45)

Similarly, it is not difficult to obtain the RD-UML spectrum formula:

$$f_{RD-UML}(\theta )={\displaystyle \frac{1}{\det [trace(I_{p_{RD}}-P_{URD}(\theta )){\hat{R}}_{URD}}}$$

(46)

where $P_{URD}(\theta )=A_{URD}(\theta ){({\mathit{A}}_{URD}^H(\theta )A_{URD}(\theta ))}^{-1}{\mathit{A}}_{URD}^H(\theta )$.

The steps of the height measurement method for MIMO radar based on a sparse array under complex terrain are summarized in

Table 1. Steps of algorithm.

Step 1	$\mathrm{Compute}\mathrm{the}\mathrm{steering}\mathrm{vectors}a({\theta }_d)$$a({\theta }_{s,i})$ according to Equation (10).
Step 2	Basic algorithm	$\mathrm{Compute}\mathrm{the}\mathrm{covariance}\mathrm{matrix}\widehat{R}$ according to Equation (22), $\mathrm{and}\mathrm{the}\mathrm{noise}\mathrm{subspace}{E}_n$$\mathrm{is}\mathrm{obtained}\mathrm{by}\mathrm{the}\mathrm{eigenvalue}\mathrm{decomposition}\mathrm{of}\widehat{R}$.
	Real-valued processing algorithm	$\mathrm{Compute}\mathrm{the}\mathrm{real}\mathrm{valued}\mathrm{covariance}\mathrm{matrix}{\widehat{R}}_U$$\mathrm{according}\mathrm{to}\mathrm{Equation}(27),\mathrm{and}\mathrm{the}\mathrm{noise}\mathrm{subspace}{U}_n$$\mathrm{is}\mathrm{obtained}\mathrm{by}\mathrm{the}\mathrm{eigenvalue}\mathrm{decomposition}\mathrm{of}{\widehat{R}}_U$.
	Reduced dimension real-valued processing algorithm	$\mathrm{Compute}\mathrm{the}\mathrm{reduced}\mathrm{dimension}\mathrm{real-valued}\mathrm{covariance}\mathrm{matrix}{\widehat{R}}_{URD}$$\mathrm{according}\mathrm{to}\mathrm{Equation}(45),\mathrm{and}\mathrm{the}\mathrm{noise}\mathrm{subspace}{U}_{RDn}$$\mathrm{is}\mathrm{obtained}\mathrm{by}\mathrm{the}\mathrm{eigenvalue}\mathrm{decomposition}\mathrm{of}{\widehat{R}}_{URD}$.
Step 3	Basic algorithm	$\mathrm{Compute}\mathrm{the}\mathrm{compound}\mathrm{steering}\mathrm{vector}A(\theta )$ according to Equation (15).
	Real-valued processing algorithm	$\mathrm{Compute}\mathrm{the}\mathrm{real-valued}\mathrm{compound}\mathrm{steering}\mathrm{vector}{A}_U(\theta )$ according to Equation (33).
	Reduced dimension real-valued processing algorithm	$\mathrm{Compute}\mathrm{the}\mathrm{reduced}\mathrm{dimension}\mathrm{real-valued}\mathrm{compound}\mathrm{steering}\mathrm{vector}{A}_{URD}(\theta )$ according to Equation (44).
Step 4	Basic algorithm	Compute the target’s elevation estimates via Equations (23) and (25).
	Real-valued processing algorithm	Compute the target’s elevation estimates via Equations (32)and (34).
	Reduced dimension real-valued processing algorithm	Compute the target’s elevation estimates via Equations (43)and (46).
Step 5	Calculated the target’s height by the target’s elevation according to the geometric relationship.

.

4. Complexity Analysis

The complexity of the proposed algorithms is analyzed. Firstly, the complexity of three algorithms mainly includes three parts: the construction of the received data covariance matrix, the eigenvalue decomposition of the covariance matrix and the search for the spectrum peak. In addition, the real-valued processing algorithm and reduced dimension real-valued processing algorithm also involve the complexity of real value processing and reduced dimension real value processing. In real-valued processing algorithms, the real-valued compound steering vector $A_U(\theta )$ and real-valued covariance matrix ${\widehat{R}}_U$ are additionally calculated. In the reduced dimension real-valued processing algorithms, the reduced dimension real-valued compound steering vector $A_{URD}(\theta )$ and reduced dimension real-valued covariance matrix ${\widehat{R}}_{URD}$ are additionally calculated. Considering that the unitary matrix and reduced dimension matrix are all sparse, the complexity of real value processing and reduced dimension real value processing is very small. So, it is ignored. Secondly, the computational complexity of addition is much less than that of multiplication. So, we only consider the multiplication. Thirdly, complex multiplication requires four times as much computation as real multiplication. Above all, the complexity of three algorithms is shown as follows. The complexity of the basic GMUSIC and ML algorithms is marked with $C_{GMUSIC}$ and $C_{GMUSIC}$. The complexity of the real-valued GMUSIC and ML algorithm is marked with $C_{UGMUSIC}$ and $C_{UML}$. The complexity of the reduced dimension real-valued GMUSIC and ML algorithms is marked with $C_{RD-UGMUSIC}$ and $C_{RD-UML}$.

$$C_{GMUSIC}=4M^{4}L+4M^6+4\Theta (8M^2+2M^4)$$

(47)

$$C_{ML}=4M^4(L+M^2)+4\Theta (8M^2+2M^4+M^6)$$

(48)

$$C_{UGMUSIC}=M^{4}L+M^6+\Theta (8M^2+2M^4)$$

(49)

$$C_{UML}=M^4(L+M^2)+\Theta (8M^2+2M^4+M^6)$$

(50)

$$C_{RD-UGMUSIC}=M_{RD}^{2}L+M_{RD}^3+\Theta (8M_{RD}+2M_{RD}^2)$$

(51)

$$C_{RD-UML}=M_{RD}^2(L+M_{RD}^2)+\Theta (8M_{RD}+2M_{RD}^2+M_{RD}^3)$$

(52)

where $\Theta$ is the number of peak searches.

5. Simulation

In this section, the height measurement algorithm based on a sparse array is simulated and compared with some other algorithms based on a uniform linear array. The GMUSIC, UGMUSIC, RD-UGMUSIC, DBF [32] and matrix pencil [38] algorithms based on a uniform linear array are marked with ULA-GMUSIC, ULA-UGMUSIC, ULA-RD-UGMUSIC, ULA-DBF and ULA-jzs in all figures. The GMUSIC, UGMUSIC and RD-UGMUSIC algorithms based on the nested array are marked with NA-GMUSIC, NA-UGMUSIC and NA-RD-UGMUSIC in all figures. The GMUSIC, UGMUSIC and RD-UGMUSIC algorithms based on the coprime array are marked with CPA-GMUSIC, CPA-UGMUSIC and CPA-RD-UGMUSIC in all figures.

In all the following simulations, the RMSE of angle estimation and the RMSE of height measurement are defined as, respectively,

$$\mathrm{RMSE}(\theta )={\displaystyle \frac{1}{P}}\sqrt{{\displaystyle \frac{1}{K}}\left({{\displaystyle \sum _{k=1}^K{\displaystyle \sum _{p=1}^P\left|{\hat{\theta }}_{p,k}-{\theta }_p\right|}}}^2\right)}$$

(53)

$$\mathrm{RMSE}(H)={\displaystyle \frac{1}{P}}\sqrt{{\displaystyle \frac{1}{K}}\left({{\displaystyle \sum _{k=1}^K{\displaystyle \sum _{p=1}^P\left|{\hat{H}}_{p,k}-H_p\right|}}}^2\right)}$$

(54)

where $P$ is the number of targets. $K$ is the number of Monte Carlo trials. ${\hat{\theta }}_{p,k}$ is the angle estimation of the $k\mathrm{th}$ experiment for the $p\mathrm{th}$ target, and ${\theta }_p$ is the angle value for the $p\mathrm{th}$ target. ${\hat{H}}_{p,k}$ is the height estimation of the $k\mathrm{th}$ experiment for the $p\mathrm{th}$ target height, and $H_p$ is the height for the $p\mathrm{th}$ target.

Example 1: Simulation 1 considers a simple terrain, that is, supposing there is only one reflection point. Assume that there are five transmitted array elements and five received array elements. The operating frequency of the radar is $f=300 \mathrm{MHz}$. The transmitted and received sensor location of the coprime array are $d_t=[5d,10d,15d,20d,25d]$ and $d_r=[0,6d,12d,18d,24d]$, respectively. The transmitted and received sensor location of the nested array are $d_t=[0,d,2d,3d,4d]$ and $d_r=[5d,11d,16d,23d,29d]$, respectively. This gives the following direct angle ${\theta }_d=3°$ and indirect angle ${\theta }_{s1}=-3°$. The height of the array (set as the lowest element height) $h_a=5\mathrm{m}$. The distance between the target and the radar $R=\mathrm{200,000}\mathrm{m}$.

Figure 4 shows the RMSE of three GMUSIC algorithms versus the SNR under different array structures in simple terrain and compared with CRB. The derivation of the Cramer–Rao bound can be found in reference [43]. This experiment mainly verifies the effectiveness of the height measurement method based on the sparse array in a simple terrain. From Figure 4, we can see that the RMSE of all the algorithms decreases with the increasing SNR, which agrees with the theory analysis. Compared with the GMUSIC based on ULA, the RMSE of the angle and height estimation based on a coprime array and nested array is much smaller, indicating that the three GMUSIC algorithms based on the coprime array and nested array are all effective. In addition, the three different algorithms based on GMUSIC have their own advantages and disadvantages. Compared with the GMUSIC algorithm, the real-valued GMUSIC and the reduced dimension real-valued GMUSIC algorithms have poor estimation accuracy but lower computation complexity under the case of the same array structure. Figure 5 shows the RMSE versus SNR of different algorithms. From Figure 5, we can see that the ULA-DBF algorithm proposed in [32] outperforms the ULA-GMUSIC algorithm and ULA-MDL algorithm from the aspect of estimation accuracy and compared with CRB. More importantly, the performance of the proposed NA-GMUSIC algorithm is equivalent to the ULA-DBF. Figure 6 shows the RMSE versus SNR of different schemes. From Figure 6, we can see that the MIMO scheme outperforms the SIMO scheme from the aspect of estimation accuracy under the same array structure.

Example 2: Simulation 2 considers a complex terrain; that is, supposing there are two reflection points. Assume that there are ten transmitted array elements and ten received array elements. The operating frequency of the radar is $f=300 \mathrm{MHz}$. The transmitted and received sensor location of the coprime array are $d=[0,5d,6d,10d,12d,15d,18d,20d,24d,25d]$. The transmitted and received sensor location of the nested array are $d=[0,d,2d,3d,4d,5d,11d,17d,23d,29d]$. The direction of the direct wave is ${\theta }_d=3°$; the direction of two reflection waves is ${\theta }_{s1}=-5°$,${\theta }_{s2}=-6°$. The height of the array (set as the lowest element height) $h_a=5\mathrm{m}$. The distance between the target and the radar $R=\mathrm{200,000}\mathrm{m}$.

Figure 7 shows the RMSE of the angle and height estimation versus SNR under different array structures in complex terrain and compared with the CRB. This experiment mainly verifies the effectiveness of the height measurement method based on a sparse array in a complex terrain. As can be seen from Figure 7, the RMSE of the three array structures all decrease with the increase in the SNR, which is consistent with the theoretical analysis. Compared with the GMUSIC based on ULA, the RMSE values of the angle and height estimation based on a coprime array and nested array are much smaller, indicating that the three GMUSIC algorithms based on the coprime array and nested array are all effective. Figure 8 shows the RMSE versus SNR of different algorithms and compared with the CRB. It can be seen from Figure 8 that the CPA-GMUSIC algorithm outperforms the other algorithms proposed in [38]. Figure 9 shows the RMSE versus SNR of different schemes. From Figure 6, we can see that the MIMO scheme outperforms the SIMO scheme on the aspect of estimation accuracy under the same array structure.

Example 3: The computational complexity of the proposed method is studied in Simulation 3. The simulations are implemented by running the Matlab codes on a PC with Intel (R) Core(TM) i7-1165G7 2.8 GHz CPU and 16 GB memory. The average elapsed times against different search intervals used by different methods are recorded in

Table 2. Average computation time (seconds).

Search Interval (°)	GMUSIC		UGMUSIC		RD-UGMUSIC
0.02	NA	0.0463	NA	0.0477	NA	0.0319
	CPA	0.0465	CPA	0.0481	CPA	0.0315
0.005	NA	0.1800	NA	0.1859	NA	0.1241
	CPA	0.1813	CPA	0.1888	CPA	0.1232

. It can be found that the RD-UGMUSIC algorithm outperform other algorithms from the aspect of estimation efficiency. The GMUSIC algorithm and UGMUSIC algorithm are time consuming, since the RD-UGMUSIC algorithm is processed by reduced dimensions during the algorithm procedure.

Example 4: Simulation 4 considers the performance of the proposed method in cases where the number of antennas is smaller or larger. The proposed algorithms with three antennas, five antennas and six antennas based on a coprime array and nested array are marked with CPA-N3, CPA-N5, CPA-N6 and NA-N3, NA-N5, NA-N6, respectively. Figure 10 shows the RMSE versus SNR of different array elements. From Figure 10, we can see that under the same array structure, the more array elements there are, the higher the angle estimation accuracy.

6. Conclusions

In this paper, a novel sparse array-based algorithm, combined with the generalized MUSIC and maximum likelihood estimation algorithm, is proposed in height measurement in the low-elevation region. The proposed algorithm is suitable for both flat terrain and complex terrain. In complex terrain, the sparse array is introduced into the MIMO radar system as a transmitting and receiving antenna. The multipath reflected signal model of sparse array MIMO radar under complex terrain is established. Combined with the generalized multiple signal classification algorithm and the maximum likelihood algorithm, a height measurement method suitable for this model is proposed. In simple terrain, to reduce the computation cost, a real-valued algorithm and a reduced dimension real-valued algorithm are proposed. In addition, the complexity of the three algorithms based on GMUSIC is analyzed. The simulation results show that the proposed algorithm has a good angular measurement performance for a low-altitude target in both simple flat terrain and complex terrain.