Low complexity and high accuracy angle of arrival estimation using eigenvalue decomposition with extension to 2D AOA and power estimation

In this paper, an angle of arrival (AOA) estimator is presented. Accurate AOA estimation is very crucial for many applications such as wireless positioning and signal enhancement using space processing techniques. The proposed estimator is based upon applying the eigenvalue decomposition (EVD) method on the crosscorrelation matrix of the received signals at two sides of the antenna array doublets. The proposed method is named the eigenvalue-decomposition-based AOA (EDBA) estimator. In comparison with the ESPRIT algorithm, the EDBA has less complexity because the decomposition in the EDBA method is performed only once and on a smaller matrix dimension than that in the ESPRIT algorithm where the decomposition is performed twice. The other advantage is that the EDBA method has better performance than the ESPRIT algorithm. The EDBA is also extended to twodimensional (2D) AOA estimation with automatic pairing in two ways. The first one performs the 2D AOA estimation by considering the eigenvalues of the crosscorrelation matrix to estimate the azimuth angles and their corresponding eigenvectors to estimate their corresponding elevation angles. Thus, the 2D AOA estimation is performed with automatic pairing and without the need for any pairing or searching techniques. The first 2D extension of the EDBA is named the EDBA-2D estimator. Another 2D AOA estimator that is based upon the EDBA method is also presented and named the Two-EDBA estimator. This second 2D estimator performs the pairing between the azimuth and elevation angles using the alignment of the eigenvalues’ magnitudes. An additional advantage for the EDBA estimator is the fact that it provides an estimate for the received signal power paired automatically with its corresponding AOA estimate. Simulations of the proposed EDBA method and its 2D extensions with the signals’ power estimation are shown to assess their performance.


Introduction
Many applications require accurate estimation of the angle of arrival (AOA) of the received signal. Examples of these applications are signal reception enhancement [1,2] and accurate wireless positioning [3,4]. Many effects cause errors in the AOA estimation. An error in the AOA estimate can cause huge error when locating wireless devices or an increase in the bit error rate when detecting the received signal. Thus, it is of interest to develop techniques to estimate the AOA accurately.
In some cases, it is required to estimate two-dimensional (2D) angles, i.e., the azimuth and elevation angles. Besides the accuracy requirement for estimating both the azimuth and elevation angles, it is of interest to pair them for each of their corresponding source. The complexity of the estimator is another issue that requires the designers of the estimator to take care of. Also, in 2D AOA estimation, the antenna array geometry plays an important role in the estimation procedure. In this paper, we will assume two antenna array geometries and we will provide a 2D estimator, based upon the proposed estimator, for each of these geometries. The first one is the rectangular antenna array shape formulated from two-parallel uniform linear arrays (ULAs). The second antenna array geometry considered in this paper is the L-shaped antenna array formulated from two perpendicular ULAs.
ESPRIT or ESPRIT-based algorithms are proposed in the literature to perform the 2D AOA estimation, such as in [5][6][7][8]. Also, another method proposed to perform the pairing procedure is the modified propagator method (PM) [9,10]. Although these methods do not require complex searching techniques as in [11][12][13][14], they still require pairing techniques in matching the angles estimated. The 2D AOA estimator proposed in [15] uses two-parallel ULAs, whereas the method in [16] is based upon matrix enhancement and the matrix pencil algorithm. The proposed method in [17] constructs a second-order statistic based upon the Schur-Hadamard product steering vector to perform the 2D AOA estimation. In [18], the authors propose an estimator that requires two ULAs to pair the azimuth and elevation angles. In [19][20][21][22], the authors propose 2D AOA techniques that apply singular value decomposition (SVD) on specially constructed covariance matrices.
In this paper, we propose an AOA estimator that enhances the AOA estimation with low computational complexity. The proposed AOA estimator utilizes the eigenvalue decomposition (EVD) method and is named the eigenvalue decomposition-based AOA (EDBA) estimator. The EDBA applies the EVD on a crosscorrelation matrix constructed from both sides of the antenna array doublets.
The first advantage of the proposed EDBA method is its low complexity. In fact, it has lower computational complexity than that of the well-known ESPRIT algorithm due to the fact that the decomposition used in the EDBA method is applied only once and on a smaller matrix dimension than that in the ESPRIT algorithm where the decomposition is applied twice [23]. The second advantage of the proposed methods is its better performance in estimating the AOA than that of the ESPRIT algorithm. Another benefit for using the EDBA method is its ability to estimate the received signals' power with automatic pairing with their corresponding AOA estimates.
The EDBA is also extended for 2D AOA estimation in two ways. The first one is named the EDBA-2D method. In this extension, the azimuth and elevation angles are estimated and paired automatically. The EDBA-2D method starts by formulating a crosscorrelation matrix from the received signals at both sides of two-parallel ULAs. Then, the azimuth angles are estimated from the eigenvalues of the crosscorrelation matrix, and the elevation angles are estimated from the corresponding eigenvectors of the same crosscorrelation matrix. Thus, automatic pairing between the estimated azimuth and elevation angles is provided by the correspondence between the eigenvalues and the eigenvectors of the crosscorrelation matrix. The second way presented to extend the EDBA method for 2D AOA estimation uses the EDBA method twice and separately at the two ULAs of the L-shaped antenna array to estimate the azimuth and elevation angles. Then, we pair them with their corresponding sources using the alignment of the eigenvalues' magnitudes. This method is named the Two-EDBA estimator. The EDBA-2D and Two-EDBA methods have the advantage of not requiring any searching process as well as possessing the advantages of the EDBA method over the ESPRIT algorithm as discussed above.
The paper is organized as follows: Section 2 introduces the system model that forms the foundation for the EDBA estimator. The proposed EDBA estimator is explained in Section 3. Section 4 presents the 2D AOA system models. The proposed 2D extensions for the EDBA method are presented in Section 5. Estimating the received signals' power is presented in Section 6. The simulated performance of the EDBA estimator and its 2D extensions with the signals' power estimation are presented in Section 7. Finally, conclusions are shown in Section 8.

System model for the EDBA estimator
In this section, we present the narrowband received signal model that will be utilized for the AOA estimation. The antenna array is formed from M uniform antenna doublets (i.e., 2M total antenna elements) as shown in Figure 1. Each antenna doublet is formed of two antenna elements spaced by a distance d. We assume a K BPSK sources signals, {š k (t)} K k=1 , impinging upon the antenna array. The signalš k (t) is represented aš s k (t) = s k (t) cos(2π f c t) = {s k (t) exp(j2π f c t)} where f c is the carrier frequency and s k (t) = i α k (i)g(t − iT) , i is an integer which represents the time index and s k (t) is the complex low-pass equivalent ofš k (t), T is the symbol period, a k (i) = A k b k (i) where b k (i) (-1, + 1) which represents the symbol parity and A k is a positive constant which represents the amplitude of s k (t) and g(t) is the (raised cosine) pulse shaping function with where n is an integer so that s k (nT) = a k (n).
In this paper, we use the superscript notation {ˆ} to denote the estimated value of a variable (for example,ĥ is the estimated value of h). So, to start developing the received signal model, first consider the sampled received signal at the antenna elements (after the matched filter stage) that are located on the upper part of the antenna doublets from doublet #1 to doublet #M, shown in Figure 1, and denote it as r 1 (i): where and with l is the signal wavelength, d d is the distance between consecutive antenna doublets, and θ k is the AOA for the kth signal. In addition, n 1 (i) is the noise vector added to the received signal at the antenna elements located at the upper part of the antenna doublets, which is additive white Gaussian noise (AWGN) and has a covariance matrix of s 2 I M × M , where I M × M is the M × M identity matrix. Finally, (.) T represents the transpose operation.
The (M × 1) received signal vector can be written in a matrix form as and From now on, we will drop the term (i) from all terms for simplicity.
The received signal vector at the lower part of the antenna doublets set, shown in Figure 1, will be given the notation r 2 . Thus, with and n 2 is an AWGN at the lower side of the antenna array doublets. Now we showed the received signal model at the antenna arrays, and we will next propose the EDBA estimator.

Proposed eigenvalue-decomposition-based AOA (EDBA) estimator
To implement the EDBA method, we start by formulating a crosscorrelation matrix (R 21 ) between r 2 and r 1 as follows where where p k is the power of the kth received signal.
We use here the assumption that n 1 and n 2 are independent realizations of AWGN.
Considering the eigenvalue decomposition of R 21 as follows where Γ Z is a diagonal matrix with its diagonal containing the eigenvalues of R 21 , i.e., where the largest K values of the diagonal elements of Γ Z correspond to the K sources which we will call g 1 K ≡ g 1 g K . Each element of g 1 K corresponds to one of the K sources. Now, from the EVD of R 21 shown in (10), and from the definition of R 21 in (8), we have and so we can say where T is the appropriate matrix to change the basis vectors.
By substituting (12) into (11), we can write The matrix Γ Z is an M × M matrix, and Z and P s are K × K diagonal matrices. Also, the nonzero diagonal elements of Γ Z are only the first K diagonal elements and then (13) can be written as Looking at (14), we can deduce that the matrix T should be M × K. Also, the lower (M -K) × K part of T should be all zero elements, i.e., and From (14) and the definition of ϒ, P s , and Z, we can deduce ⎡ ⎢ ⎣ But, the left-hand side of (18) is a diagonal matrix, and the matrices P s and Z are diagonal as well. So, from (18), we can deduce that the matrix ϒ is diagonal too, i. e., Now, the eigenvector matrix U is defined as where u k is the kth column eigenvector of R 21 , and the K columns of U that correspond to the K largest (which are also nonzero) eigenvalues of R 21 will be given the notationŪ, i.e., Since the rest (M -K) columns of U in (12) are multiplied by the zero matrix in (17), then the eigenvector matrix (Ū) and the matrix A(θ) span the same signal subspace, and we have where Ω = ϒ -1 and is a diagonal matrix (since ϒ is a diagonal matrix). Thus, each column ofŪ is a rotated scaled version of a corresponding column in A(θ ).
Looking at the vector a(θ k ) which is one of the columns in A(θ ) (as shown in (5)) and defined in (3) and (4), we can deduce θ k as follows where Although A(θ) (from which a(θ k ) and consequently a μ (θ k ) can be deduced) is not available at the receiver, still we could deduce the matrixŪ from U. Thus, looking back at (25) in which ζ k is expressed in terms the elements of a(θ k ), and because u k is a rotated scaled version of a(θ k ), then ζ k that is the average relative phase between two consecutive elements of a(θ k ) can be expressed in terms of two corresponding consecutive elements of u k as follows where u k , μ is the μth element of the vector u k . Thus, θ k can be estimated as followŝ The performance of the proposed method will be shown in the simulation section (Section 7). In comparison with the ESPRIT method, the EDBA uses the EVD on an M × M matrix as shown in (10), whereas in the ESPRIT method, the EVD will be applied on a 2M × 2M for the same system model assumed [23]. Also, the EVD in the EDBA is applied only once, whereas in the ESPRIT method, the EVD is applied twice for estimating the AOA as shown in [23]. Next, we will present the system models that will be used for the extension of the EDBA method for estimating the 2D AOA with automatic pairing.

System model for 2D AOA estimation
For the 2D estimation, we will present two ways to extend the EDBA method to estimate the azimuth and elevation angles. There are two antenna alignments over which the extension will be presented. In the first antenna alignment, we will assume the parallel ULAs shown in Figure 2. The second antenna alignment that we will consider is the L-shaped antenna array illustrated in Figure 3. In the next subsections, we will provide more insight into the system models for the parallel ULA and the L-shaped antenna array.

Parallel ULA
The parallel ULA is formed from M antenna parallel doublets (i.e., 2M total antenna elements) as shown in Figure 2. Each antenna doublet is formed of two antenna elements spaced by a distance d. To start developing the received signal model, first consider the sampled received signal after the matched filter stage (r A ) at the antenna elements that are located on the right-hand side of each antenna doublets, shown in Figure 2, which is represented as and with l is the signal wavelength and d d is the distance between consecutive antenna elements at each ULA. In addition, n A is the noise vector added to the received signal at the antenna elements (from m = 1 M) located along the right-hand side of each antenna doublet, which is AWGN and has a covariance matrix of s 2 I M × M .
The (M × 1) received signal vector can be written in matrix form as The received signal vector at the left-hand side of the antenna doublets set, shown in Figure 2, will be given the notation r B . Thus, with and n B is an AWGN vector at the left-hand side of the antenna array doublets.
The next subsection will provide some insight into the L-shaped antenna array.

The L-shaped antenna array
The L-shaped antenna array is formed of two perpendicular antenna arrays as shown in Figure 3. The sampled received signal (r z ) at the antenna elements (after the matched filter stage) that are located in the z-axis direction is represented as where n z is the AWGN noise vector with a covariance matrix of s 2 I M × M .
Similarly, the sampled received signal at the antenna elements, from m = 1 to m = M, that are located in the x-axis direction is given the notation r x and is represented as Likewise, n x is the AWGN noise vector, and it also has a covariance matrix of s 2 I M × M .
So both (M × 1) received signal vectors can be written in matrix form as Thus, the received signal models for the different 2D antenna arrays structures are presented.

Extension of the EDBA for 2D AOA estimation
In this section, we will propose two methods for which the EDBA can be used to estimate the azimuth and elevation angles with automatic pairing. The first method is named the EDBA-2D method that is based upon the parallel antenna array structure illustrated in Figure 2. The second method is named the Two-EDBA method that is based upon the L-shaped antenna structure in Figure 3.

Proposed EDBA-2D estimator
The EDBA-2D method is based upon the system model introduced in Section 4.1 and illustrated in Figure 2. The EDBA-2D estimation starts by formulating a crosscorrelation matrix (R BA ) between r B and r A as follows Considering the eigenvalue decomposition of R BA as follows where˜ Y is a diagonal matrix with its diagonal containing the eigenvalues of R BA , i.e., where the largest K values of the diagonal elements of correspond to the K sources which we will call γ 1→K ≡γ 1 →γ K . Each element ofγ 1→K corresponds to one of the K sources. Also,Ū is the eigenvector matrix for R BA . Now, from the EVD of R BA shown in (38), and from the definition of R BA in (37), we see that So, we can saỹ whereT is the appropriate matrix to change the basis vectors.
By substituting (40) into (39), we can writẽ The matrix˜ Y is an M × M matrix and P s and Y are K × K diagonal matrices. Also, the nonzero diagonal elements of˜ are only the first K diagonal elements. Thus, (41) can be written as Looking at (42), we can deduce that the matrix˜ should be M × K. Also, the lower (M -K) × K part of should be all zero elements, i.e., From (42) and the definition ofΥ , P s , and Y, we can deduce ⎡ ⎢ ⎣γ But, the left-hand side of (45) is a diagonal matrix, and the matrices P s and Y are diagonal as well. So, from (45), we can deduce that the matrixT is diagonal too, i. e.,Υ Substituting (46) in (42) and looking at the nonzero diagonal elements in both sides of (42), i.e.,γ 1 →γ K , the following equation can be deuced: where (.)* is the complex conjugate operator. Thus, Since |t k,k | 2 is a magnitude square and p k is a power variable, then both terms are not complex variables. Thus, looking at both sides of (47), then But the azimuth angle (j k ) is implicated in ∠y k as can be deuced from (33). Thus, the azimuth angle can be estimated from the eigenvalue of R BA as followŝ For the elevation angle (θ k ), it is estimated from the eigenvectors of R BA in a similar way as discussed in the EDBA estimator explained in Section 3. We will not repeat the discussion and derivation again in this section but rather we will give the final procedure on how to estimate the elevation angle θ k .
We can calculateζ k which is the average relative phase between two consecutive elements ofã(θ k ) in terms of two corresponding consecutive elements ofũ k as follows whereũk,μ is the μth element of the vectorũ k . Thus, θ k can be obtained as followŝ Because of eigenvalue-eigenvector correspondence, estimation of the azimuth angle (j k ) for the kth signal from the kth eigenvalue (γ k ) and its corresponding elevation angle (θ k ) from the corresponding kth eigenvector (ũ k ) is performed with automatic pairing and without the need of any pairing procedure.

Proposed two-EDBA estimator
The Two-EDBA method is based upon the L-shaped antenna structure in Figure 3. The Two-EDBA starts by estimating the azimuth and elevation angles separately. This is a straightforward procedure using the EDBA method on each of the received signals on both ULAs at each axis. Then, after estimating the azimuth (j k ) and elevation (θ k ) angles, the pairing is achieved by the alignment of the eigenvalues' magnitudes. From Section 3, it is shown that the EDBA method estimates the angles from the eigenvectors of the crosscorrelation matrix R 21 . Since the eigenvalues' magnitudes of the crosscorrelation matrix R 21 depend on the received signals' power (which are assumed to be the same on both ULAs on both axis), then the pairing between the estimated azimuth and elevation angles is achieved by aligning the estimated azimuth and elevation angles depending on their corresponding eigenvalues' magnitudes.

Estimating the received signals' power
This section shows how the received signals' power can be estimated after estimating the received signals' AOAs. This is achieved with an automatic pairing between the received signals' power with their corresponding AOAs. Pairing can be performed for the EDBA estimator or its 2D estimator extensions. For sake of simplicity, we will show how the received signals' power is estimated for the simple EDBA estimator explained in Section 3, and the extension for the 2D EDBA estimators can be easily deduced.
To explain how the received signals' power is estimated, let us assume first that the angles of the received signals using the EDBA method are estimated as explained in Section 3. Then, the vector a(θ k ) in (3) and (4) can be reformulated using the estimated angles (θ k ) for k = 1 K. This leads to the ability of reformulating the matrix A(θ) in (5), and we will give it the notation A(θ ). Since the eigenvector matrix U is available from (10), and the matrixÂ(θ ) is known, then the appropriate matrix (T) in (12) can be estimated as followŝ Similarly, the matrix Z depends on θ k in its formulation as shown (6) and (7). Thus, fromθ k , we can reformulate the matrix Z, and we will give it the notationẐ .
Since the eigenvalue matrix Γ Z is available from (10) and then from (13), the received signals' power matrix (P s ) can be estimated as followŝ where (.) + is the pseudo-inverse operator. Thus, the received signal power can be deduced from the diagonal elements ofP s , and we will give it the notationp k . The pairing betweenθ k andp k is automatically provided from the fact thatθ k is estimated from the eigenvectors in U andp k is estimated from their corresponding eigenvalues in Γ Z . Since each eigenvector in U corresponds automatically to a specific eigenvalue in Γ Z , thenθ k is automatically paired with its correspondingp k . Thus, the estimation of the received signal power is estimated with  Similarly, the estimation of the received signal power can be easily applied to the 2D AOA estimation extensions for the EDBA explained in Section 5, resulting in estimating the azimuth/elevation angles jointly with the received signals' power with automatic pairing.
The performance of the proposed methods will be shown in the simulation section (Section 7).

Simulation results
Simulations of the proposed EDBA estimator with its extensions were completed to assess their performance. The elements of each antenna array were separated by a half-wavelength (i.e., d = λ 2 ). Unless mentioned otherwise, the number of antenna doublets was set to M = 10. The number of snap shots was set to W = 256 and W = 1,000 (where W is the number of snap shots over which the correlation matrices were estimated). The number of sources was set to 2. The near far ratio (NFR) was set to -10 dB, where the NFR is defined as 10 log ) . The proposed EDBA and its extensions were compared with the ESPRIT method. Figure 4 shows the root-mean-square error (RMSE) of the angle estimation in degrees for the proposed EDBA method compared with the ESPRIT method for different signal to noise ratios (SNR)s. Since the power for both sources was not equal, then we assumed that the SNR was taken for the first source. The angles for the two sources were set to θ 1 = 75°and θ 2 = 45°. The results in Figure 4 show that the proposed EDBA method gave better performance than the ESPRIT method. Figure 5 shows the RMSE of the angle estimation in degrees for the proposed EDBA-2D method compared with the ESPRIT method for different SNRs. The elevation angles for the two sources were set to θ 1 = 75°a nd θ 2 = 45°. The azimuth angles for the two sources were set to j 1 = 80°and j 2 = 55°. The results in Figure 5 indicate that the proposed EDBA-2D method gave very close performance to the ESPRIT method with the less complexity advantage for the proposed EDBA-2D method. Also, the EDBA-2D performs the estimation of the azimuth and elevation angles with automatic pairing.
The simulation presented in Figure 6 shows the RMSE for the proposed EDBA-2D method compared with the ESPRIT method for different source angle deviation. Table 4 Table of successful pairing rate using the proposed EDBA-2D and Two-EDBA algorithm compared to the algorithm in [14] and the SVD-based algorithm in [22] for different number of antenna elements and W = 1, 000 M Algorithm in [14](%) SVD-based algorithm in [22](%) EDBA-2D (%) Two-EDBA (%)   Table 2 Table of successful pairing rate using the proposed EDBA-2D and Two-EDBA algorithm compared to the algorithm in [14] and the SVD-based algorithm in [22] for different SNRs and W = 1,000 SNR (dB) Algorithm in [14](%) SVD-based algorithm in [22] The angular deviation is performed by setting j 1 to be taken from the following equation j 1 = j 2 + δ where δ is the angular deviation with θ 1 = 75°, θ 2 = 45°, and j 2 = 55°. The SNR for the first source was set to 10 dB. The results in Figure 6 indicate that the proposed EDBA-2D method gave better performance than the ESPRIT method for different AOA deviations. In Figure 7, the simulation shows the power RMSE for the proposed EDBA for different SNRs. The results in Figure 7 show that the EDBA method managed to estimate the received signal power with high accuracy. Tables 1 and 2 show the pairing success rate for the proposed EDBA-2D and Two-EDBA algorithms compared to the algorithm in [14] and the SVD-based algorithm in [22] for different SNRs in dB with W = 256 and W = 1, 000, respectively. The number of antenna elements M was set to 10. Both tables indicate that the proposed EDBA-2D and Two-EDBA algorithms managed to pair the azimuth and elevation angles successfully. Tables 3 and 4 show the pairing success rate for the proposed EDBA-2D and Two-EDBA algorithms compared to the algorithm in [14] and the SVD-based algorithm in [22] for different number of antenna elements with W = 256 and W = 1, 000, respectively. The SNR was set to 15 dB for the first source. Both tables show that the proposed EDBA-2D and Two-EDBA algorithms managed to pair the azimuth and elevation angles successfully.
Thus, the proposed EDBA-2D and Two-EDBA methods have high capability of pairing the estimated azimuth and elevation angles with good accuracy and in an automatic procedure.

Conclusion
In this paper, we propose an AOA estimator that we named the EDBA method. The EDBA method is applied by taking the EVD of the received signal crosscorrelation matrix. The AOA of the received signals is contained in the eigenvectors of the crosscorrelation matrix. So, the AOA of the received signals is deduced from these eigenvectors. Two extension methods for 2D estimation are introduced. The first extension method estimates the elevation angles from the eigenvectors of the crosscorrelation matrix, and the corresponding azimuth angles are estimated from the corresponding eigenvalues. This method is named the EDBA-2D. The second extension method is named the Two-EDBA estimator and uses the alignment of the eigenvalues' magnitudes to pair the azimuth and elevation angles. Also, the EDBA method is extended to estimate the received signals' power with automatic pairing with their corresponding AOAs. Numerical simulation indicated that the EDBA method and its extensions outperformed the ESPRIT AOA estimator and other pairing methods. Also, the EDBA method and its extensions are low complex and do not require any searching or pairing procedure to perform the estimation.
Author details 1 Electrical Engineering Department, University of Ha'il, Ha'il, Saudi Arabia