Super-Resolution Processing for Multiple Aperture Antenna to Suppress Multipath

Abstract

Angle estimation for low-altitude targets above the sea surface is a challenging problem due to multipath interference from surface reflection signals, and various approaches have been proposed. This paper proposes a matrix pencil method with multiple apertures. The matrix pencil method effectively responds to dynamic scenarios because it performs better when using a single snapshot than other methods. Also, employing multiple apertures is more economical than using one large aperture. Therefore, we propose a computationally efficient approach using this method and structures. The proposed two-stage MP method incrementally improves the resolution in two stages: in stage 1, we extract the denoised signals at each aperture level, and in stage 2, we further improve the resolution with those signals. In comparison with the angular resolution defined by the half-power beamwidth (HPBW) of a uniform linear array (ULA) antenna with an equivalent number of arrays, the proposed method demonstrated a superior resolution of less than 0.087 of the HPBW at a high signal-to-noise ratio (SNR) of 40 dB, and less than 0.31 of it even at a relatively low SNR of 15 dB, based on 90% of the resolving probability. For the multipath problem, the proposed scheme has the advantage of not requiring prior geometric information, and its performance is demonstrated through simulations to be better than the adaptive beamforming method and the composite monopulse method.

1. Introduction

Tracking low-altitude targets with maritime radar is challenging due to the phenomenon known as multipath interference, which arises from the specular reflection and diffractions in addition to the direct signal from the target as described in Figure 1.

In particular, the specular reflection of the dominant component is similar to the target signal in angle, time (distance), and Doppler domains, making it difficult to distinguish using conventional methods. The obvious way to avoid interference is to use a narrow beamwidth to exclude the specular reflection signal from the mainlobe beam. However, reducing beamwidth requires a larger antenna or higher frequency, both of which are generally impractical. Hence, when tracking low-altitude targets, the performance of maritime radars degrades due to multipath interference, and various research efforts have been made to address this issue [1,2]. Studies have focused on the angle estimation method and can be summarized in three categories: the complex monopulse ratio approach, super-resolution techniques, and the nulling method of the surface reflection signal. First, the complex monopulse method was introduced by Sherman [3]. In conventional monopulse systems, the sum and the difference signal for a single target are in-phase so that the ratio of the difference signal to the sum signal is a real number and is directly proportional to the angle [4]. However, a multipath signal introduces a quadrature component into the monopulse ratio, causing erratic angle measurement. Thus, complex monopulse methods have been proposed using the ratio of a complex number including the quadrature component. Sherman defined the complex indicated angle (CA) by the phase of the monopulse ratio and provided CA as a function of the target angle [1,3]. Based on this CA function, Chung et al. proposed a specular reflection error finding function (SREFF) indicating the affected region of the specular reflection in CA data and compensated errors [5]. Blair et al. proposed the probability density function (pdf) of the measured amplitude of the sum signal and amplitude-conditioned pdf of the in-phase and quadrature monopulse ratios for Rayleigh targets in the presence of sea-surface-induced multipath [6]. These methods require a priori geometrical information and estimated depression angle. The second method employs super-resolution techniques of subspace-based estimation to measure the angle in array antennae. These techniques include multiple signal classification (MUSIC), Root–MUSIC, and estimation of signal parameters via rotational invariant techniques (ESPRIT) [7,8,9]. However, subspace-based estimation techniques suffer from rank deficiency of the covariance matrix in the presence of correlated incoming signals, such as in a multipath situation, which makes angle estimation difficult and computationally expensive. Finally, adaptive beamforming is the primary technique used to null the specular reflection [10,11]. After removing the specular signal, they compensate or correct the distorted monopulse slope. Another method is the double-null (DN) technique, using an additional auxiliary null steered in the direction of the specular reflection in the difference beam pattern [12,13]. This method considers only the reflected signal arriving in the difference beam. While the adaptive beamforming method also requires a priori geometrical information to null the specular signal, the DN method is relatively simple and does not demand prior information. However, the angle estimation error in the DN method is on the order of one-tenth beamwidth.

In this paper, we investigate the matrix pencil method (MPM), another high-resolution technique that can be used in array antennas [14,15,16,17], and present a method for combining signals from multiple apertures. The MPM is originally a fast algorithm that uses uniformly sampled data to find the frequency of a superimposed sinusoid and can be employed to find the spatial frequencies, i.e., angles in the array antennae. Unlike subspace methods, it does not calculate the covariance matrix and is, hence, functional even when only a single snapshot is available, such as in dynamic scenarios [18,19]. The authors of [19] reported that MPM outperforms ESPRIT when only a few snapshots are available. They showed that MPM has a lower root-mean-square (RMS) error than the ESPRIT technique for single snapshots, 8 snapshots, and 16 snapshots. This difference increases as the number of snapshots decreases. For example, RMSE is about 1.7 dB lower when using a single snapshot and 0.6 dB lower when using 16 snapshots, according to Figure 1.

The aperture size influences angular resolution across nearly all techniques, albeit differently than in conventional beamforming. It is especially relevant when employing a single snapshot. If high resolution is only required in certain situations, such as multipath conditions in maritime radar, it is not economical or even impractical to increase the size of the antenna. Therefore, recent radar systems consist of multiple antennas that can be operated individually or collectively depending on the situation, e.g., innovative multiple-object tracking radar (IMOTR) [20]. The proposed two-stage MPM is developed for such a multiple aperture system and gradually improves the angular resolutions by stages. In stage 1, MPM is performed to extract the denoised signals from each aperture, and in stage 2, further MPM is performed with these signals to enhance the resolution. The proposed method has shown better performance than adaptive beamforming or the complex monopulse (SREFF) method despite not requiring a priori geometric information compared with other methods.

The remainder of this manuscript is organized as follows. Section 2 provides a summary of the MPM and describes the proposed method. Section 3 presents the signal model of multipath interference and analyzes the performance of the proposed method through simulations. The conclusions are provided in Section 4.

2. Matrix Pencil Method (MPM)

2.1. Summary of MP Method

Consider a uniform linear array (ULA) consisting of N identical sensors and receiving $M_t$ narrow-band coherent signals arriving from directions ${\theta }_1,{\theta }_2,\dots ,{\theta }_{M_t}$ shown in Figure 2.

At any instant, these $M_t$ signals $r_1\left(t\right),r_2\left(t\right),\dots ,r_{M_t}\left(t\right)$ are phase-delayed, amplitude-weighted replicas of one of the signals—say, the first one [21]; hence,

$$r_m\left(t\right)={\alpha }_m{r}_1\left(t\right),m=\mathrm{1,2},\dots ,M_t$$

(1)

where ${\alpha }_m$ is a complex attenuation of the m-th signal with respect to the first signal $r_1$(t), and ${\alpha }_1=1$. The measurement vector obtained by the arrays is expressed as

$$\overrightarrow{y}\left(t\right)=\mathbf{A}\overrightarrow{r}\left(t\right)+\overrightarrow{n}\left(t\right),$$

(2)

where $\overrightarrow{r}\left(t\right)={\left[r_1\left(t\right)r_2\left(t\right)\dots r_{m_t}\left(t\right)\right]}^T$ is the $M_t\times 1$ received signal vector, $\overrightarrow{y}\left(t\right)={\left[y_1\left(t\right) y_2\left(t\right)\dots y_N\left(t\right)\right]}^T$ is the $N\times 1$ array output vector, and $\overrightarrow{n}\left(t\right)$ is the $N\times 1$ noise signal vector generated through the antenna, transmitter, and receiver. A is an $N\times M_t$ matrix with Vandermonde-structured distinct columns $(N>M_t)$, and, hence, is of rank $M_t$. It is expressed as

$$\mathbf{A}=[\overrightarrow{a}\left({\theta }_1\right),\overrightarrow{a}\left({\theta }_2\right),\dots ,\overrightarrow{a}({\theta }_{M_t}\left)\right]$$

(3)

with $\overrightarrow{a}\left({\theta }_m\right)$ being the direction vector associated with the arrival angle ${\theta }_m$, i.e.,

$$\overrightarrow{a}\left({\theta }_m\right)={\left[1,\exp \left(-j\kappa d\sin {\theta }_m\right),\exp \left(-j2\kappa d\sin {\theta }_m\right),\dots \exp (-j\left(N-1\right)\kappa d\sin {\theta }_m)\right]}^T$$

(4)

where $\kappa \left(=2\pi /\lambda \right)$ is the wavenumber and T denotes the matrix transpose. The n-th element of the measurement vector $\overrightarrow{y}$at a specific time $t_s$ is represented by

$$y_n=y_n\left(t_s\right)=\sum _{m=1}^{M_t}{\alpha }_m{r}_1\left(t_s\right)\exp \left(-j(n-1)\kappa d\sin {\theta }_m\right)+n_n\left(t_s\right)=\sum _{m=1}^{M_t}{\overline{r}}_m z_m^{n-1}+{\overline{n}}_n$$

(5)

where ${\overline{n}}_n=n_n\left(t_s\right),$ ${\overline{r}}_m={\alpha }_m{r}_1\left(t_s\right),$ and

$$z_m=\exp \left[j(-\kappa d\sin {\theta }_m)\right], for m=1,\dots ,M_t.$$

(6)

Thus, the matrix $\mathbf{A}$ in (3) and $\overrightarrow{y}\left(t\right)$ in (2) can be rewritten as

$$\mathbf{A}=\left[\begin{array}{cccc}1& 1& \dots & 1\\ z_1& z_2& \dots & z_{M_t}\\ z_1^2& z_2^2& & z_{M_t}^2\\ ⋮& ⋮& \ddots & ⋮\\ z_1^{N-1}& & & z_{M_t}^{N-1}\end{array}\right].$$

(7)

$$\overrightarrow{y}=\overrightarrow{y}\left(t_s\right)=\left[\begin{array}{c}{y}_1\\ y_2\\ ⋮\\ y_N\end{array}\right]=\mathbf{A}\left[\begin{array}{c}\begin{array}{c}{\overline{r}}_1\\ {\overline{r}}_2\end{array}\\ ⋮\\ {\overline{r}}_N\end{array}\right]+\left[\begin{array}{c}{\overline{n}}_1\\ {\overline{n}}_2\end{array}\right]=\mathbf{A}\overline{r}+\overline{n}..$$

(8)

To apply the MPM, we define a Hankel-like matrix Y and two $(N-L)\times L$ matrices ${\mathbf{Y}}_1$ and ${\mathbf{Y}}_2$, as follows [15]:

$$\mathbf{Y}=\left[\begin{array}{ccccc}{y}_0& y_1& \cdots & y_{L-1}& y_L\\ y_1& y_2& \cdots & y_L& y_{L+1}\\ ⋮& ⋮& \ddots & ⋮& ⋮\\ y_{N-L-1}& y_{N-L}& \cdots & y_{N-2}& y_{N-1}\end{array}\right]\in {\mathbf{C}}^{\left(N-L\right)\times \left(L+1\right)},$$

(9)

$${\mathbf{Y}}_1=\left[\begin{array}{cccc}{y}_0& y_1& \cdots & y_{L-1}\\ y_1& y_2& \cdots & y_L\\ ⋮& ⋮& \ddots & ⋮\\ y_{N-L-1}& y_{N-L}& \cdots & y_{N-2}\end{array}\right]\in {\mathbf{C}}^{\left(N-L\right)\times L},$$

(10)

$${\mathbf{Y}}_2=\left[\begin{array}{cccc}{y}_1& \cdots & y_{L-1}& y_L\\ y_2& \cdots & y_L& y_{L+1}\\ ⋮& \ddots & ⋮& ⋮\\ y_{N-L}& \cdots & y_{N-2}& y_{N-1}\end{array}\right]\in {\mathbf{C}}^{\left(N-L\right)\times L}$$

(11)

where L denotes the pencil parameter. For efficient noise filtering, L is set between N/3 and N/2 [15,22]. For these values of L, the noise-induced variance in z is found to be the minimum.

Now consider the matrix pencil

$${\mathbf{Y}}_2-\lambda {\mathbf{Y}}_1$$

(12)

In general, the rank of this matrix pencil will be $M_t$, and z is found as the generalized eigenvalues of the matrix pair $\{{\mathbf{Y}}_2;{\mathbf{Y}}_1\}$. The problem of solving for z can be cast as an ordinary eigenvalue problem for $L\times L$ matrix,

$$det[{\mathbf{Y}}_1^{†}{\mathbf{Y}}_2-\lambda \mathbf{I}]=0,$$

(13)

where I is an identity matrix, ${\mathbf{Y}}_1^{†}$ is the Moore–Penrose pseudoinverse of ${\mathbf{Y}}_1$, and det means the determinant. After finding the eigenvalues $z_m$, the angle can be calculated based on (6), as follows:

$${\theta }_m={\sin }^{-1}\left[-{\displaystyle \frac{\arg \left(z_m\right)}{\kappa d}}\right], for m=1,\dots ,M_t.$$

(14)

Due to the noise, the number of predicted eigenvalues, ${\overline{M}}_t$, may exceed $M_t$, and predicting the number of targets accurately becomes difficult. However, the input signal can be selected using the least-squares method solution of (8), as follows:

$${\overrightarrow{r}}_e={\tilde{\mathbf{A}}}^{†}\overrightarrow{y}, \tilde{\mathbf{A}}=\left[\begin{array}{cccc}1& 1& \dots & 1\\ z_1& z_2& \dots & z_{{\overline{M}}_t}\\ z_1^2& z_2^2& & z_{{\overline{M}}_t}^2\\ ⋮& ⋮& \ddots & ⋮\\ z_1^{N-1}& & & z_{{\overline{M}}_t}^{N-1}\end{array}\right]$$

(15)

where ${\overrightarrow{r}}_e$ is the estimated vector of $\overrightarrow{r}$, and ${\tilde{\mathbf{A}}}^{†}$ is the Moore–Penrose pseudoinverse of $\tilde{\mathbf{A}}$. The m-th element of the estimated vector${\overrightarrow{r}}_e$ represents the signal component in the direction ${\theta }_m$. The number of targets can be determined at this point based on the amplitude of the elements if the SNR is adequately large.

2.2. Forward-and-Backward Matrix Pencil

The forward-and-backward (FB) MPM is a variation of the MPM that can be used for undamped sinusoidal sequences [23], where it improves the accuracy of estimation. This method is applied by modifying (12) as follows:

$${\overline{\mathbf{Y}}}_2-\lambda {\overline{\mathbf{Y}}}_1\triangleq \left[\begin{array}{c}{\mathbf{Y}}_{2f}\\ {\mathbf{Y}}_{2b}\end{array}\right]-\lambda \left[\begin{array}{c}{\mathbf{Y}}_{1f}\\ {\mathbf{Y}}_{1b}\end{array}\right],$$

(16)

where ${\mathbf{Y}}_{1f}={\mathbf{Y}}_1, {\mathbf{Y}}_{2f}={\mathbf{Y}}_2$ in (10) and (11),

$${\mathbf{Y}}_{1\mathit{b}}=\left[\begin{array}{cccc}{y}_{N-1}^{\mathrm{*}}& y_{N-2}^{\mathrm{*}}& \cdots & y_{N-L}^{\mathrm{*}}\\ y_{N-2}^{\mathrm{*}}& y_{N-3}^{\mathrm{*}}& \cdots & y_{N-L-1}^{\mathrm{*}}\\ ⋮& ⋮& \ddots & ⋮\\ y_L^{\mathrm{*}}& y_{L-1}^{\mathrm{*}}& \cdots & y_1^{\mathrm{*}}\end{array}\right] \mathrm{and} {\mathbf{Y}}_{2\mathit{b}}=\left[\begin{array}{cccc}{y}_{N-2}^{\mathrm{*}}& y_{N-3}^{\mathrm{*}}& \cdots & y_{N-L-1}^{\mathrm{*}}\\ y_{N-3}^{\mathrm{*}}& y_{N-4}^{\mathrm{*}}& \cdots & y_{N-L-2}^{\mathrm{*}}\\ ⋮& ⋮& \ddots & ⋮\\ y_{L-1}^{\mathrm{*}}& y_{L-2}^{\mathrm{*}}& \cdots & y_0^{\mathrm{*}}\end{array}\right].$$

(17)

The FB MPM also yields the desired eigenvalues at $z_m\prime s$. When $M/2\le L<M$, the two matrices ${\overline{\mathbf{Y}}}_2$ and ${\overline{\mathbf{Y}}}_1$ have a rank $M_t$ [24]. The FB MPM is more robust to noise than the forward-only or backward-only method if the signal poles are known to be on the unit circle, as in this case. In addition, when this method is applied, $\overrightarrow{y}$ and $\tilde{\mathbf{A}}$ in (15) are replaced by

$${\overrightarrow{r}}_e={\tilde{\mathbf{A}}}^{†}\overrightarrow{y}, \overrightarrow{y}=\left[\begin{array}{c}{y}_{\left(N-1\right)}\\ y_{\left(N-2\right)}\\ ⋮\\ y_1\\ y_0\\ y_1\\ ⋮\\ y_{N-1}\end{array}\right], \mathrm{and} \tilde{\mathbf{A}}=\left[\begin{array}{c}\begin{array}{cccc}{\left(z_1^{\mathrm{*}}\right)}^{N-1}& {\left(z_2^{\mathrm{*}}\right)}^{N-1}& \cdots & {\left(z_L^{\mathrm{*}}\right)}^{N-1}\\ ⋮& ⋮& & ⋮\\ {\left(z_1^{\mathrm{*}}\right)}^2& {\left(z_2^{\mathrm{*}}\right)}^2& \ddots & {\left(z_L^{\mathrm{*}}\right)}^2\\ z_1^{\mathrm{*}}& z_2^{\mathrm{*}}& \cdots & z_L^{\mathrm{*}}\\ 1& 1& \dots & 1\\ z_1& z_2& \dots & z_L\\ z_1^2& z_2^2& & z_L^2\\ ⋮& ⋮& \ddots & ⋮\\ z_1^{N-1}& & & z_L^{N-1}\end{array}\end{array}\right].$$

(18)

The simulations were performed for ULA shown in Figure 2 with 20 arrays ($N=20$) with half-wavelength spacing. The two targets were separated by $\mathsf{\Delta }\theta$, shown in Figure 3, centered on the boresight, and had the same SNR of 20 dB, and the reflective coefficients of the targets had the same amplitude but differed in phase by uniform distribution over [$0, 2\pi ]$. Figure 4 and Figure 5 show the resolving probability and RMS errors from a Monte Carlo simulation with 6000 trials.

The resolving probability is calculated based on the following criteria:

(1)

Power level: there should be only two angles of power levels above 0 dB, based on the estimated amplitude of the signals using (15).

(2)

Accuracy: the estimated angle difference should be less than 40% of the true value, i.e., between $0.6\mathsf{\Delta }\theta$ and $1.4\mathsf{\Delta }\theta$.

Figure 4 shows that the resolving probability of the FB MPM surpasses that of the standard MPM. For example, in the case of targets separated by 3°, while the MPM resolves them with a probability of 0.7, the FB MPM can resolve with a probability of 0.93. The figure also shows that both the standard and FB MPMs outperform the conventional beamforming method with HPBW of approximately 5.1°.

Figure 5 shows that for targets separated by 3°, the rms error of the FB MPM is around 0.34°, which is significantly less than that of the standard MPM (0.56°). The errors were calculated only for the resolved cases.

Figure 6 presents the data histograms of the angle measurements: the left figure shows the case of $\mathsf{\Delta }\theta ={2.4}^{\mathrm{o}}$where the target was present at $\pm {1.2}^{\mathrm{o}}$, corresponding to the probability of 0.7 in Figure 4; the right figure shows the case $\mathsf{\Delta }\theta ={3.0}^{\mathrm{o}}$ where the target was located at $\pm {1.5}^{\mathrm{o}}$, corresponding to the probability of 0.93.

2.3. Proposed Method of Two Stage FB MP

The resolution capability is affected by the aperture size, i.e., the number of arrays in ULA. Increasing the number of arrays faces significant costs, making it uneconomical when high resolution is conditionally required. In this work, we consider the configuration of multiple ULA spaced at an interval D to increase the overall aperture, as shown in Figure 7a. Suppose the spacing D between the apertures is an integer multiple of the array spacing d. In that case, the configuration in (a) can be made by missing some arrays from the configuration in (b), reducing the number of arrays. Although FB MPM for all arrays can be applied by modifying matrix Y in (9), such as inserting zeros at the missing position, there are some issues that will be discussed later.

We propose a two-stage FB MPM to gradually improve the angular resolution through stages. In Stage 1, the FB MPM is used to obtain the eigenvalues with measurements of a single aperture and estimate the angle of the target. Because the apertures are in the same location and the target is oriented in the same direction, the eigenvalues processed by each antenna signal are the same. Therefore, the matrices are substituted in the following manner and the eigenvalue problem of (16) is solved once.

$${\overline{\mathbf{Y}}}_2=\left[\begin{array}{c}{\left[{\mathbf{Y}}_{2f}\right]}_1\\ {\left[{\mathbf{Y}}_{2f}\right]}_2\\ ⋮\\ {\left[{\mathbf{Y}}_{2f}\right]}_{N_r}\\ ----\\ {\left[{\mathbf{Y}}_{2b}\right]}_1\\ {\left[{\mathbf{Y}}_{2b}\right]}_2\\ ⋮\\ {\left[{\mathbf{Y}}_{2b}\right]}_{N_r}\end{array}\right],{\overline{\mathbf{Y}}}_1=\left[\begin{array}{c}{\left[{\mathbf{Y}}_{1f}\right]}_1\\ {\left[{\mathbf{Y}}_{1f}\right]}_2\\ ⋮\\ {\left[{\mathbf{Y}}_{1f}\right]}_{N_r}\\ ----\\ {\left[{\mathbf{Y}}_{1b}\right]}_1\\ {\left[{\mathbf{Y}}_{1b}\right]}_2\\ ⋮\\ {\left[{\mathbf{Y}}_{1b}\right]}_{N_r}\end{array}\right].$$

(19)

Then, with the obtained eigenvalues $z_m(m=1\dots {\overline{M}}_t)$, we estimate the signal vector ${\left[ {\overrightarrow{r}}_e\right]}_k (k=1,\dots ,N_r)$ at the k-th antenna using (18):

$${\left[{\overrightarrow{r}}_e\right]}_k={\tilde{\mathbf{A}}}^{†}{\overrightarrow{y}}_k$$

(20)

The computed ${\left[{\overrightarrow{r}}_e\right]}_k$ becomes the denoised signal in particular direction, allowing the enhancement of resolution in the subsequent stage.

In Stage 2, the FB MPM is executed for signals with the same eigenvalue or the same angle. Due to the separation between the antennas (D) and the arrival angle $\left({\theta }_m\right)$, the signals from different antennas have a phase difference of $\kappa D\sin {\theta }_m$. These signals are sinusoidal signals again and can be estimated using the FB MPM. If the number of signals estimated in Stage 1 is M, then FB MPM in Stage 2 should be performed M times. The measurement vector used in the m-th FB MPM is composed as follows:

$$\overrightarrow{y}=\left[{\left[{\overrightarrow{r}}_e\right]}_{1,m} {\left[{\overrightarrow{r}}_e\right]}_{2,m} \cdots {\left[{\overrightarrow{r}}_e\right]}_{N_r,m}\right] (m=\mathrm{1,2},\dots M)$$

(21)

where ${\left[{\overrightarrow{r}}_e \right]}_{k,m}$ is the m-th component of ${\left[{\overrightarrow{r}}_e\right]}_k$. If the signals have been sufficiently decomposed in Stage 1, ${\left[{\overrightarrow{r}}_e\right]}_{k,m}$ will not be resolved further in Stage 2, because it contains no other signal components. The example below illustrates this further. The measured angle in Stage 2 is determined using (22):

$${\theta }_k={\sin }^{-1}\left[-{\displaystyle \frac{\arg \left(z_k\right)}{\kappa D}}\right]={\sin }^{-1}\left[-{\displaystyle \frac{\arg \left(z_k\right)}{D}}\cdot {\displaystyle \frac{\lambda }{2\pi }}\right].$$

(22)

Here, since $D>0.5\lambda$, the angle may be ambiguous. However, as the second FB MPM is performed for a signal decomposed once in Stage 1, the unambiguous region or the ambiguity number can be selected. In an array antenna, the ambiguity angle—the angle at which grating lobes appear in conventional beamforming—is determined by the following equation [25]:

$$\sin {\theta }_p=\sin {\theta }_o+{\displaystyle \frac{p\lambda }{D}}, p=\pm 1,\pm 2,\pm 3,\dots$$

(23)

where p is defined as the ambiguity number.

For example, if the spacing D between the antennas is $15\lambda$, then the target at ${\theta }_0={1.91}^{\mathrm{o}}$, and targets at $-{5.72}^{\mathrm{o}} \left(p=-2\right),-{1.91}^{\mathrm{o}} \left(p=-1\right)$, ${5.72}^{\mathrm{o}} \left(p=1\right)$, ${7.63}^{\mathrm{o}} \left(p=2\right),$ etc., have the same frequency; thus, the true angle is ambiguous. However, if the angle is known to be within a region based on the result of Stage 1, we can determine the ambiguity number. Alternatively, we can determine the necessary antenna spacing D by calculating the resolving capability of Stage 1 with the minimum SNR.

A total of 6000 Monte Carlo simulations were conducted using the antenna parameters specified in

Table 1. Simulation parameters.

Parameter	Value
$\mathrm{Number} \mathrm{of} \mathrm{arrays} \mathrm{in} \mathrm{one} \mathrm{aperture}, N$	20
$\mathrm{Array} \mathrm{spacing} \mathrm{in} \mathrm{one} \mathrm{aperture}, d$	$0.5 \lambda$
$\mathrm{Number} \mathrm{of} \mathrm{apertures}, N_r$	4
$\mathrm{Distance} \mathrm{between} \mathrm{apertures}, D$	$15 \lambda$

. We evaluated the efficiency of the proposed method on two front targets positioned at $\pm \theta /2,$ while changing θ. The reflective coefficients of the targets have the same amplitude but differ in phase by uniform distribution over [$0, 2\pi ]$.

Figure 8 displays the resolving probability for the two targets. In Stage 1, two targets can be resolved by more than a probability of 0.9 if the angular difference is greater than 1.88°. If Stage 2 is executed subsequently, a resolution of 0.3° can be achieved. Because the FB MPM in Stage 2 operates on one specified component of ${\overrightarrow{r}}_e$ from Stage 1, the target signal need not be resolved further if it has been sufficiently decomposed in Stage 1, as previously mentioned. Therefore, the resolving probability in Stage 2 drops sharply for angular differences above 3°, as observed in Figure 8. In addition, because D is equal to $15\lambda$, the angular ambiguity must be addressed when the angular difference between the targets exceeds 3.86°, according to (23). However, if the targets are decomposed in Stage 1, the ambiguity need not be considered. The final resolution depends on the SNR and is much less than the HPBW, which is 0.85° due to the total aperture size of four antennas spaced 15 wavelengths apart. In addition, the total number of arrays is 80, which is only 67% of the 120 arrays that would be required for filling the same size with arrays uniformly spaced 0.5 wavelengths apart. Hence, the proposed approach is highly cost-effective. Note that the HPBW of the ULA with 80 arrays spaced at 0.5 wavelength is 1.27 degrees.

Figure 9 shows the RMS errors of the estimated angles. The error at $\mathsf{\Delta }\theta =2^o$ is 0.22° in Stage 1, decreasing to 0.03° after Stage 2. The final rms errors after Stage 2 are mostly less than 0.03 for [0.5, 2.0]. As the SNR increases, the resolution probability increases, and the rms error decreases, as shown in Figure 10 and Figure 11. Even at SNR = 15 dB, more than 90% of targets apart 0.03° are resolved.

2.4. Summary and Discussion

The proposed two-stage MPM is summarized with the related equation in

Table 2. Summary of the proposed two-stage FB MPM.

Stage	Description	Related Eqs.
1	$\mathrm{a}. \mathrm{From} \mathrm{the} \mathrm{array} \mathrm{output} y_k(k=0\dots N-1)$$, \mathrm{make} \mathrm{four} \mathrm{matrices} {\left[{\mathbf{Y}}_{1f}\right]}_i,{\left[{\mathbf{Y}}_{2f}\right]}_i,{\left[{\mathbf{Y}}_{1b}\right]}_i,{\left[{\mathbf{Y}}_{2b}\right]}_i$	(17)
	$\mathrm{b}. \mathrm{Make} \mathrm{the} \mathrm{concatenate} \mathrm{matrices} {\overline{\mathbf{Y}}}_1\mathrm{a}\mathrm{n}\mathrm{d}{\overline{\mathbf{Y}}}_2$	(19)
	$\mathrm{c}. \mathrm{Solve} \mathrm{the} \mathrm{eigenvalue} \mathrm{problem} \mathrm{and} \mathrm{find} \mathrm{eigenvalue} z_m$	(16)
	d. Estimate angles with relatively low resolution	(14)
	$\mathrm{e}. \mathrm{Calculate} \mathrm{the} \mathrm{target} \mathrm{signals} {\left[{\overline{r}}_e\right]}_i$ for each aperture	(18)
2	$\mathrm{a}. \mathrm{Make} \mathrm{the} \mathrm{measurement} \mathrm{vector} \overrightarrow{y}$	(21)
	$\mathrm{b}. \mathrm{Make} \mathrm{four} \mathrm{matrices} {\mathbf{Y}}_{1f},{\mathbf{Y}}_{2f},{\mathbf{Y}}_{1b},{\mathbf{Y}}_{2f}$$\mathrm{with} \overrightarrow{y}$	(17)
	$\mathrm{c}. \mathrm{Solve} \mathrm{the} \mathrm{eigenvalue} \mathrm{problem} \mathrm{and} \mathrm{find} \mathrm{eigenvalue} z_m$	(16)
	d. Estimate angles with enhanced resolution	(22)
	e. Calculate the target signals for checking the number of targets	(18)

.

In the configuration of Figure 7a, if we perform a standard FB MPM by inserting zeros into the empty spaces, the size of ${\overline{Y}}_1 \mathrm{o}\mathrm{r} {\overline{Y}}_2$in (16) will be $\left[2\left(N{N}_r+N_z-L\right)\right]\times L$, where $N_r$ is the number of zeros inserted. For $N_r>2$, without loss of generality, the pencil parameter L must be larger than N to achieve optimal MPM. A large L increases the computational complexity of the eigenvalue problem of (16). The computational cost of the eigenvalue problem is known to be practically $O\left(L^3\right)$. Furthermore, due to the periodic presence of zeros in the row vectors of $\overline{Y}$, a spurious signal modulated with this frequency might be incorporated into the measurement, resulting in ambiguous angles.

If L is selected to be smaller than N, the computational cost of eigenvalue problem can be reduced. However, the performance of the MPM degrades. For example, the ${\overline{Y}}_1$ and ${\overline{Y}}_2$ matrices become nearly identical to those in (19) if inefficient row vectors including zeros are eliminated. Then, the resolution also approaches the result of stage 1 in the proposed method. In other words, computational cost is reduced at the expense of performance. Stage 2 of the proposed method can enhance resolution with low computational cost. The size of the measurement vector $\overrightarrow{y}$ in (21) depends on the number of apertures; therefore, the size of the matrix in the eigenvalue problem is small. In addition, this measurement vector is composed of denoised signals from Stage 1, sufficient to improve the resolution.

The proposed method has additional advantages. The aperture spacing D is not necessarily an integer multiple of the array d, and each aperture can be operated independently with different missions, such as tracking targets in different directions if high resolution is not required. The four matrices ${\left[{\mathbf{Y}}_{1f}\right]}_i,{\left[{\mathbf{Y}}_{2f}\right]}_i,{\left[{\mathbf{Y}}_{1b}\right] }_i,{\left[{\mathbf{Y}}_{2b}\right]}_i$ in (19) do not necessarily need to concatenate the matrices of all the apertures, and only when high resolution is required, such as in a multipath situation, can they be concatenated appropriately and can we proceed to Stage 2 to enhance the resolution.

3. Multipath Model and Simulation Result

Figure 1 shows that the multipath interference signal consists of a specular signal and a diffuse signal in addition to the reflected signal from the target. The specular signal is more dominant at low angles and smooth surfaces than the diffuse signal.

In Figure 12, $a_E$ is an Earth radius, and $h_t$, $h_r$ is the height of the target and the antenna, respectively. The distances from the antenna and target to the specular point A are $r_{ar}$ and $r_{at}$. Also, $r_t$ is the distance to the target, and $r_1$ and $r_2$ are the ground distances. Given $h_t$ and $h_r$, all the values in the figure can be obtained from $r\left(=r_1+r_2\right)$ using the equations in

Table 3. List of equations.

Equation	Output
$2r_1^3-3r_1^{2}r+\left[r^2-2a_E\left(h_r+h_t\right)\right]r_1+2a_E{h}_{r}r=0$ Alternatively, $r_1={\displaystyle \frac{r}{2}}-p\sin \left({\displaystyle \frac{\phi }{3}}\right),p={\displaystyle \frac{1}{\sqrt{3}}}\sqrt{4a_E\left(h_r+h_t\right)+r^2} \mathrm{a}\mathrm{n}\mathrm{d} \phi ={\sin }^{-1}\left({\displaystyle \frac{2a_E\left(h_r-h_t\right)r}{p^3}}\right)$	$r_1$
$r_2=r-r_1$	$r_2$
$r_{ar}^2={\left(h_r+a_E\right)}^2+a_E^2-2\left(h_r+a_E\right)a_E\cos \left({\displaystyle \frac{r_1}{a_E}}\right)$	$r_{ar}$
$r_{at}^2={\left(h_t+a_E\right)}^2+a_E^2-2\left(h_t+a_E\right)a_E\cos \left({\displaystyle \frac{r_2}{a_E}}\right)$	$r_{at}$
$r_t^2={\left(h_r+a_E\right)}^2+{\left(h_t+a_E\right)}^2-2\left(h_r+a_E\right)\left(h_t+a_E\right)a_E\cos \left({\displaystyle \frac{r}{a_E}}\right)$	$r_t$
$\psi ={\sin }^{-1}\left({\displaystyle \frac{2a_E{h}_r+h_r^2-r_{ar}^2}{2a_E{r}_{ar}}}\right)$	$\psi$
${\theta }_t={\sin }^{-1}\left({\displaystyle \frac{{\left(h_t+a_E\right)}^2-r_t^2-{\left(h_r+a_E\right)}^2}{2r_t(h_r+a_E)}}\right),{\theta }_r={\sin }^{-1}\left({\displaystyle \frac{2a_E{h}_r+h_r^2+r_{ar}^2}{2\left(a_E+h_r\right)r_{ar}}}\right)$	${\theta }_t, {\theta }_r$

below [10,26,27].

The reflection signal vector, $\overrightarrow{\mathit{r}}\left(t\right)=[r_t\left(t\right),r_i\left(t\right),r_d(t\left)\right]$ consists of the target reflection $r_t\left(t\right)$, the image or specular signal $r_i\left(t\right)$, and the diffuse signal $r_d\left(t\right)$, with the latter two being represented as follows.

$$r_i\left(t\right)=\mathsf{\Gamma }D{\rho }_s{r}_t\left(t\right)e^{-j\alpha }, r_d\left(t\right)={\rho }_d{r}_t\left(t\right)$$

(24)

where $\mathsf{\Gamma }$ is the Fresnel reflection coefficient, ${\rho }_s$ is the scattering coefficient for a rough surface, D is the divergence factor, $\alpha$ is the phase lag given by $2\pi (r_{ar}+r_{at}-r_t)/\lambda$, and ${\rho }_d$ is the diffuse reflection coefficient.

The reflection coefficient $\mathsf{\Gamma }$ is defined as the ratio of the amplitude of the reflected wave to the amplitude of the incident wave and depends on the electrical properties of the surface, grazing angle $\psi$, frequency, and the polarization [6].

$$\mathsf{\Gamma }={\rho }_0{e}^{-j\varphi }=\left\{\begin{array}{c}{\displaystyle \frac{\epsilon \sin \psi -\sqrt{\epsilon -{\cos }^2\psi } }{\epsilon \sin \psi +\sqrt{\epsilon -{\cos }^2\psi }}}, vertical polarization\\ {\displaystyle \frac{\sin \psi -\sqrt{\epsilon -{\cos }^2\psi }}{\sin \psi -\sqrt{\epsilon -{\cos }^2\psi }}}, horizontal polarization\end{array}\right.$$

(25)

where the complex dielectric constant $\epsilon$ is given by

$$\epsilon ={\displaystyle \frac{{\epsilon }_s-{\epsilon }_0}{1+{\omega }_c^2{\tau }^2}}-j\left[{\displaystyle \frac{\left({\epsilon }_s-{\epsilon }_0\right){\omega }_c\tau }{1+{\omega }_c^2{\tau }^2}}-4\pi {\displaystyle \frac{{\sigma }_i}{{\omega }_c}}\right]$$

(26)

where ${\epsilon }_0=4.9$for seawater, ${\omega }_c/2\pi$ is the carrier frequency, and the other parameters are dependent on the seawater temperature, as listed in

Table 4. Parameters of seawater.

Parameter	${10}^{\mathbf{o}}\mathbf{C}$	${20}^{\mathbf{o}}\mathbf{C}$	Units
$\mathrm{Static} \mathrm{dielectric} \mathrm{parameter} {\epsilon }_s$	72.2	69.1	-
$\mathrm{Relaxation} \mathrm{time} \mathrm{of} \mathrm{seawater}, \tau$	$1.21\times {10}^{-11}$	$9.21\times {10}^{-12}$	$\mathrm{s}$
$\mathrm{Ionic} \mathrm{conductivity} \mathrm{of} \mathrm{sea} \mathrm{water}, {\sigma }_i$	$3.6\times {10}^{10}$	$4.7\times {10}^{10}$	${\mathrm{s}}^{-1}$

.

The scattering coefficient ${\rho }_s$ can be obtained by

$${\rho }_s=\left\{\begin{array}{l}{e}^{-8{\pi }^2{g_0}^2}, 0\le g_0<0.1\\ {\displaystyle \frac{0.812537}{1+8{\pi }^2{g_0}^2}}, 0.1\le g_0\end{array}\right.$$

(27)

where $g_0={\sigma }_h/\lambda \sin \psi$ and ${\sigma }_h$ is the RMS sea-surface elevation above the mean level. The divergence factor D can be approximated as the following equation [28].

$$D\approx {\left(1+{\displaystyle \frac{2r_1{r}_2}{a_{e}r\sin \psi }}\right)}^{-1/2}$$

(28)

The diffuse component can be described statistically in nature and is modeled by using a complex Gaussian zero-mean random variable for each direction distributed independently at many points. The simplified one described in [2] considers that the total diffuse reflection is concentrated in one point corresponding to the specular reflection point. The Rayleigh parameter of the diffuse component is described as

$${\sigma }_d=\left\{\begin{array}{l}\sqrt{2}\left|\mathsf{\Gamma }\right|3.68 g_0, 0\le g_0<0.1\\ \sqrt{2}\left|\mathsf{\Gamma }\right|\left(0.454-0.858\right)g_0, 0.1\le g_0<0.5\\ \sqrt{2}\left|\mathsf{\Gamma }\right|0.025, 0.5\le g_0\end{array}\right.$$

(29)

Then, the diffuse reflection coefficient becomes

$${\rho }_d={\sigma }_d\sqrt{{\displaystyle \frac{\left(x_1^2+x_2^2\right)}{2}}} e^{j{\varphi }_d}$$

(30)

where $x_1$ and $x_2$ are independent Gaussian random variables with zero mean and unit variance, and ${\varphi }_d$ is uniformly distributed over $\left[0, 2\pi \right]$.

To compare the proposed method with previous related works, we performed the simulation under the same conditions as Ref. [11].

Table 5. Simulation parameters for the radar and environment.

Parameters	Values
$\mathrm{Antenna} \mathrm{height}, h_r$	15 m
$\mathrm{Number} \mathrm{of} \mathrm{antenna} \mathrm{elements}, N$	$80$
Antenna boresight angle	$0^{\mathrm{o}}$
Radar frequency	10 GHz
Polarization	Vertical
$\mathrm{RMS} \mathrm{sea}-\mathrm{surface} \mathrm{elevation}, {\sigma }_h$	0.4
$\mathrm{Target} \mathrm{height}, h_t$	35 m
Signal-to-noise ratio	40 dB
$\mathrm{Temperature}, T$	${20 }^{\mathrm{o}}\mathrm{C}$

summarizes the radar and environmental parameters. We use the same configuration of multiple antennae described in Table 1, aligned in elevation direction. The target is moving toward the radar with constant velocity and height, and its range varies from 20 to 2 km. The elevation angle is measured at every 200 m interval, and the Monte Carlo simulations are performed for 10,000 runs.

In the proposed method, the two angles with the highest power over the threshold of 0 dB after two-stage FB MP were chosen as the target and the specular direction. The proposed method is compared with the results of three methods in Ref. [11]: The ML estimator, the SREFF method, and adaptive beamforming. Figure 13a,b show the estimated angle and RMS errors, respectively. The ML estimator exhibits superior performance, particularly in the range with a slow varying angle, but requires significant computational power. The proposed method outperforms the technique using SREFF and the adaptive beamforming method. Moreover, the proposed method is computationally simple, and does not require prior geometric information.

Another advantage of the proposed method is that the performance can be estimated according to the SNR. The performance according to the SNR is as presented in the previous section. Unlike the monopulse technique, where the angle error exhibits cyclic behavior [1], the accuracy improves as the SNR increases. In other words, the estimation accuracy can be improved as the target approaches the radar. Figure 14 shows the performance of a simulating target signal with an SNR of 20 dB at 20 km. The SNR is assumed to be inversely proportional to the fourth power of the range. Figure 14a displays two estimated angles in the range where the multipath phenomenon mainly appears, while Figure 14b shows the RMS error according to the SNR. As the SNR decreases, the RMS error tends to increase.

4. Conclusions

This paper proposed a novel two-stage MPM to enhance the angular resolution and accuracy using multiple apertures radar. The proposed method consists of two stages, with denoised signal extraction from each antenna in Stage 1 and combined processing of these signals in Stage 2. The resolving probability of two targets was analyzed with respect to the SNR, and the spacing between antennas to avoid angular ambiguity was provided. The proposed method is beneficial for tracking a sea-surface target in multipath conditions because it requires no prior geometrical information and uses relatively low computation power. We evaluated the performance through simulation, where it demonstrated its superiority over other methods.

It is worth mentioning that multiple aperture antennas can be operated individually or collectively depending on the situation or the required accuracy. Furthermore, the proposed two-stage MPM can be applied to radars using MIMO structures such as automotive radars, because the array signal of the receiving arrays in the MIMO antenna is spatially convolved at each transmit array location and has a configuration similar to that of the multiple aperture structure.