# Improved ESPRIT Method for Joint Direction-of-Arrival and Frequency Estimation Using Multiple-Delay Output

Authors: Wang Xudong; Xiaofei Zhang; Jianfeng Li; Jinchao Bai
Journal: International Journal of Antennas and Propagation
 (2012)

---

## Abstract

An automatic pairing joint direction-of-arrival (DOA) and frequency estimation is presented to overcome the unsatisfactory performances of estimation of signal parameter via rotational invariance techniques- (ESPRIT-) like algorithm of Wang (2010), which requires an additional pairing. By using multiple-delay output of a uniform linear antenna arrays (ULA), the proposed algorithm can estimate joint angles and frequencies with an improved ESPRIT. Compared with Wang’s ESPRIT algorithm, the angle estimation performance of the proposed algorithm is greatly improved. The frequency estimation performance of the proposed algorithm is same with that of Wang’s ESPRIT algorithm. Furthermore, the proposed algorithm can obtain automatic pairing DOA and frequency parameters, and it has a comparative computational complexity in contrast to Wang’s ESPRIT algorithm. By the way, this proposed algorithm can also work well for nonuniform linear arrays. The useful behavior of this proposed algorithm is verified by simulations.


---

## Body

## 1. Introduction

Uniform linear antenna (ULA) arrays  have been used in radar, sonar, electron reconnaissance, seismic data processing, and so on [1–4]. The direction-of-arrival (DOA) and frequency estimation of signals impinging on the ULAs are two fundamental problems in array processing [5–10]. The problem of joint DOA and frequency estimation arises in the applications of fields mentioned above. For example, these parameters can be applied to locate the targets for radars and to locate pilot tones in electron reconnaissance systems [11]. Furthermore, a precise estimation of these parameters is helpful to attain a better pulse descriptor word (PDW) and thus enhances the system performance. Optimal techniques based on maximum likelihood [12] are often applicable but might be computationally prohibitive. Some ESPRIT-based joint angle and frequency estimation methods have been proposed in [13–17]. Mathews [13] discuss this problem in the context of radar applications. Pro-ESPRIT is proposed to estimate angle and frequency. Haardt and Nossek [14] discuss the problem in the context of mobile communications for space division multiple access applications. Their method is based on unitary-ESPRIT, which involves a certain transformation of the data to real valued matrices. Multiresolution ESPRIT is used for joint angle frequency estimation in [15]. ESPRIT method is used for frequency and angle estimation under uniform circular array in [16, 17]. These ESPRIT-based joint angle and frequency estimation methods can give satisfactory estimation performances, but require additional pairing.The authors in [18] proposed a joint angle and frequency estimation using multiple-delay output based on ESPRIT, which outperform the conventional ESPRIT algorithm, but it requires an additional pairing. The additional pairing fails to work for low signal-to-noise ratio (SNR), which can be shown from the simulation results in Figure 1. Furthermore, the angle estimation of Wang’s ESPRIT algorithm can be improved further.Figure 1
Angle-frequency scatter of the Wang’s ESPRIT algorithm at SNR = 0 dB.An improved ESPRIT for joint direction-of-arrival (DOA) and frequency estimation is presented in this paper. Compared with Wang’s ESPRIT algorithm, the algorithm has the improved performance of DOA estimation. Moreover, the proposed algorithm can obtain automatically paired DOA and frequency estimation, while Wang’s ESPRIT algorithm requires the additional pairing. The proposed algorithm has a comparative computational complexity in contrast to Wang’s ESPRIT algorithm. Moreover, this proposed algorithm can also work well for nonuniform linear arrays.Notation 1.
(
·
)
T, (·)H, (·)-1, (·)+, and ⊙ denote transpose, conjugate transpose, pseudoinverse operations, and Hadamard product, respectively. Angle(·) refers to the phase of a complex number. diag(v) stands for diagonal matrix, the diagonal of which is a vector v.


## 2. The Data Model

We consider a uniform linear array with spacingd. There is an array of M sensors on which K incident waves impinge (M>K), if the signals are all within the assumed narrow (with respect to the center frequency) band of the receiver, the signal received at the mth antenna is
(1)ym(t)=∑k=1Kej2πmdfksin(θk)/csk(t)+nm(t),m=0,1,2,…,M-1,
where sk(t) is the narrow-band signal of the kth source. c is velocity of light. θk,fk are the DOA and the frequency of the kth signal, respectively. nm(t) represents an additive noise term which is assumed to be zero mean and stationary.The outputs of the uniform linear antenna arrays (ULA) is(2)Y0=[y0(n),y1(n),y2(n),…,yM-1(n)]T,n=1,2,3,…,N.Based on this available samples, the problem is to estimate the angles and the frequencies of all sources from theN times snapshots. Suppose the number of signals K is assumed known. Denote the state vector, we get
(3)Y0=AS+N0,
where
(4)S=[s1,s2,…,sK]T∈CK×N,N0=[n0,n1,…nM-1]T∈CM×N,A=[11⋯1e-jα1e-jα2⋯e-jαK⋮⋮⋱⋮e-j(M-1)α1e-j(M-1)α2⋯e-j(M-1)αK],
where αk=2πdfksin(θk)/c, k=1,2,…,K.In order to joint estimate DOA and frequency, we addP-1 delayed outputs for the received signal of array antenna. We suppose that0<(P-1)τ<1/max(fi).The delayed signal for (1) with delay τ is
(5)ym(t-τ)=∑k=1Kej2πmdfksin(θk)/csk(t-τ)+nm′(t)=∑k=1Kej2πmdfksin(θk)/csk(t)e-j2πfkτ+nm′(t).The delayed signal for (5) with τ can be denoted as
(6)Y1=AΦS+N1,
where
(7)Φ=diag{e-jβ1,e-jβ2,…,e-jβK},
where βk=2πfkτ, k=1,2,…,K.The delayed signal for (5) with pτ can be denoted as
(8)Yp=AΦpS+Np,p=0,1,…,P-1.According to (3), (6), and (8), we define
(9)Y=[Y0Y1⋮YP-1]=[AAΦ⋮AΦP-1]S+[N0N1⋮NP-1].


## 3. Joint Angle and Frequency Estimation

We can use received signal to attain the direction matrixA and the delay matrix Φ, and then estimate angle and frequency. The covariance matrix of the received signal can be reconstructed via RY=YYH. Using eigenvalue decomposition of RY, we can get the signal subspace Es. In the free-noise case, Es can be denoted as
(10)Es=[AAΦ⋮AΦP-1]T,
where T is a K×Kfull-rank matrix.


### 3.1. Frequency Estimation

According to (10), we define E1 and E2(11)E1=[AAΦ⋮AΦP-2]T,E2=[AΦAΦ2⋮AΦP-1]T.
According to (11),
(12)E2=[AΦAΦ2⋮AΦP-1]T=[AAΦ⋮AΦP-2]TT-1ΦT=E1T-1ΦT.LetΨ=T-1ΦT, so Ψ=E1+E2. We use eigenvalue decomposition of Ψ to get the estimate of Φ. Also we obtain an estimation of the matrix T-1 via the eigenvectors of Φ. In the no-noise case,
(13)T^-1=T-1Π(14)Φ^=ΠTΦΠ,
where Π is a permutation matrix, and Π-1=ΠT.BecauseΨ has the same eigenvalues as Φ, we can get βk, k=1,2,…,K from the eigenvalues λk, k=1,2,…,K of matrix Ψ. And then estimate frequency fk, k=1,2,…,K(15)f^k=12πτangle(λk).


### 3.2. Angle Estimation

There exists a transformation matrixH corresponding to the finite number of row interchanged operations for subspace Es such that
(16)Ec=HEsT^-1=HEsT-1Π=FΠ,
where
(17)F=HEsT-1=H[AAΦ⋮AΦP-1]=[GGΦ⋮GΦM-1],
where
(18)G=[11⋯1e-jβ1e-jβ2⋯e-jβK⋮⋮⋱⋮e-j(P-1)β1e-j(P-1)β2⋯e-j(P-1)βK],Φ=diag{e-jα1,e-jα2,…,e-jαK}∈CK×K,
where
(19)αk=2πdfksin(θk)/c,k=1,2,…,K.
According to (16) and (17), take the first (M-1)P rows of Ec, we define matrix Ec1 as follows:
(20)Ec1=[GGΦ⋮GΦM-2]Π.
Taking the last (M-1)P rows of Ec, we define matrix Ec2(21)Ec2=[GΦGΦ2⋮GΦM-1]Π.
We define
(22)C=Ec1+Ec2.
According to (20) and (21), in free-noise case, C can be denoted as
(23)C=ΠTΦΠ.
Then, the DOA estimation is obtained via
(24)θ^k=arcsin(c2πfkdangle(ξk)),k=1,2,…,K,
where ξk is the kth diagonal element of matrix C.The following summarizes the major steps of this proposed algorithm.(1)
Compute the covariance matrix from the received signal viaRY=YYH.(2)
Get the signal subspaceEs by using EVD of RY. Construct the matrices E1 and E2 from Es via (11) and (12).(3)
Employ EVD onE1+E2 to get T^-1 and Φ^. Estimate frequency via (15).(4)
ComputeEc by (16) and construct matrices Ec1 and Ec2 via (20) and (21).(5)
Compute the matrixC=Ec1+Ec2 and estimate DOA by (24).Remark 1.
The proposed algorithm can obtain automatically parameter estimation.From (14) and (23), DOA and frequency estimations suffer from the same permutation ambiguity, so we can get automatically paired DOA and frequency.


### 3.3. Complexity Analysis

Note thatM, K, N and P are the number of the antennas, the sources, the snapshots, and delays, respectively. Computational complexity of Wang’s ESPRIT algorithm and this proposed algorithm is a function of M, K, N, and P.In contrast to Wang’s ESPRIT algorithm, the proposed algorithm has a comparative computational load. For the proposed algorithm, the covariance matrix estimation costsO(M2P2N); eigen decomposition of RY needs O(M3P3); the computation of E1+E2 requires O(2K^2M(P-1)+2K^3); eigen decomposition of E1+E2 needs O(K3); compute the matrix C=Ec1+Ec2 needs O(2K^2(M-1)P+2K^3). The major computational complexity of this algorithm is
(25)O(M2P2N+M3P3+(2K2M(P-1)+2K3+K3)+(2K2(M-1)P+2K3)),
while Wang’s ESPRIT algorithm requires
(26)O(M2P2N+M3P3+(2K2M(P-1)+2K3+K3)+(2K2(M-1)+2K3+3K3)).
According to [19], we derive the CRB for the frequency estimation and angle estimation as
(27)CRB=σ22N{Re[DHΠA⊥D⊙P^sT]},-1
where  D=[(∂a1/∂f1), (∂a2/∂f2)⋯(∂aK/∂fK), (∂a1/∂θ1), (∂a2/∂θ2)×⋯(∂aK/∂θK)],   with ak being the kth column of A(28)P^s=1N∑t=1NststH,ΠA⊥=IPM-A(AHA)-1AH.


### 3.4. Discussion

The performance of this proposed algorithm with nonuniform linear array will be discussed in this section. The distance between themth antenna and the left reference element is dm(m=0,…,M-1), where d0=0.There is an array ofM sensors on which K incident waves impinge (M>K). The direction matrix A′ is
(29)A′=[11⋯1e-jd1γ1e-jd1γ2⋯e-jd1γK⋮⋮⋱⋮e-jdM-1γ1e-jdM-1γ2⋯e-jdM-1γK],
where γk=2πfksin(θk)/c, k=1,2,…,K.Replacing matrixAwith A′ in (3)–(12), we can get the frequency estimation by (15).For a nonuniform linear array, (20) and (21) can be rewritten as follows:
(30)E′c1=[GGΦ1⋮GΦM-2]Π=[E′c11E′c12⋮E′c1(M-1)],E′c2=[GΦ1GΦ2⋮GΦM-1]Π=[E′c21E′c22⋮E′c2(M-1)],
where
(31)Φj=diag{e-jdjγ1,e-jdjγ2,…,e-jdjγK}∈CK×K,j=1,2,…,M-1,E′c1j and E′c2j represent the (j-1)P+1 to jPth rows of matrix E′c1 and E′c2, respectively.We define(32)Cj=Ec1j+Ec2j,j=1,2,…,M-1.
Then, the DOA estimation can be obtained via
(33)θ^k=1M-1∑j=1M-1arcsin(c2πfk(dj-dj-1)angle(ξkj)),k=1,2,…,K,
where ξkj is the kth diagonal element of matrix Cj.


### 3.5. The Advantages of the Proposed Algorithm

The proposed algorithm has the following advantages.(1)
The proposed algorithm has a better angle estimation performance than Wang’s ESPRIT algorithm.(2)
The proposed algorithm can obtain automatically paired parameter estimation, while Wang’s ESPRIT algorithm requires additional pairing.(3)
The proposed algorithm has a comparative computational complexity in contrast to Wang’s ESPRIT algorithm.(4)
This proposed algorithm also suit for nonuniform linear arrays.


## 3.1. Frequency Estimation

According to (10), we define E1 and E2(11)E1=[AAΦ⋮AΦP-2]T,E2=[AΦAΦ2⋮AΦP-1]T.
According to (11),
(12)E2=[AΦAΦ2⋮AΦP-1]T=[AAΦ⋮AΦP-2]TT-1ΦT=E1T-1ΦT.LetΨ=T-1ΦT, so Ψ=E1+E2. We use eigenvalue decomposition of Ψ to get the estimate of Φ. Also we obtain an estimation of the matrix T-1 via the eigenvectors of Φ. In the no-noise case,
(13)T^-1=T-1Π(14)Φ^=ΠTΦΠ,
where Π is a permutation matrix, and Π-1=ΠT.BecauseΨ has the same eigenvalues as Φ, we can get βk, k=1,2,…,K from the eigenvalues λk, k=1,2,…,K of matrix Ψ. And then estimate frequency fk, k=1,2,…,K(15)f^k=12πτangle(λk).


## 3.2. Angle Estimation

There exists a transformation matrixH corresponding to the finite number of row interchanged operations for subspace Es such that
(16)Ec=HEsT^-1=HEsT-1Π=FΠ,
where
(17)F=HEsT-1=H[AAΦ⋮AΦP-1]=[GGΦ⋮GΦM-1],
where
(18)G=[11⋯1e-jβ1e-jβ2⋯e-jβK⋮⋮⋱⋮e-j(P-1)β1e-j(P-1)β2⋯e-j(P-1)βK],Φ=diag{e-jα1,e-jα2,…,e-jαK}∈CK×K,
where
(19)αk=2πdfksin(θk)/c,k=1,2,…,K.
According to (16) and (17), take the first (M-1)P rows of Ec, we define matrix Ec1 as follows:
(20)Ec1=[GGΦ⋮GΦM-2]Π.
Taking the last (M-1)P rows of Ec, we define matrix Ec2(21)Ec2=[GΦGΦ2⋮GΦM-1]Π.
We define
(22)C=Ec1+Ec2.
According to (20) and (21), in free-noise case, C can be denoted as
(23)C=ΠTΦΠ.
Then, the DOA estimation is obtained via
(24)θ^k=arcsin(c2πfkdangle(ξk)),k=1,2,…,K,
where ξk is the kth diagonal element of matrix C.The following summarizes the major steps of this proposed algorithm.(1)
Compute the covariance matrix from the received signal viaRY=YYH.(2)
Get the signal subspaceEs by using EVD of RY. Construct the matrices E1 and E2 from Es via (11) and (12).(3)
Employ EVD onE1+E2 to get T^-1 and Φ^. Estimate frequency via (15).(4)
ComputeEc by (16) and construct matrices Ec1 and Ec2 via (20) and (21).(5)
Compute the matrixC=Ec1+Ec2 and estimate DOA by (24).Remark 1.
The proposed algorithm can obtain automatically parameter estimation.From (14) and (23), DOA and frequency estimations suffer from the same permutation ambiguity, so we can get automatically paired DOA and frequency.


## 3.3. Complexity Analysis

Note thatM, K, N and P are the number of the antennas, the sources, the snapshots, and delays, respectively. Computational complexity of Wang’s ESPRIT algorithm and this proposed algorithm is a function of M, K, N, and P.In contrast to Wang’s ESPRIT algorithm, the proposed algorithm has a comparative computational load. For the proposed algorithm, the covariance matrix estimation costsO(M2P2N); eigen decomposition of RY needs O(M3P3); the computation of E1+E2 requires O(2K^2M(P-1)+2K^3); eigen decomposition of E1+E2 needs O(K3); compute the matrix C=Ec1+Ec2 needs O(2K^2(M-1)P+2K^3). The major computational complexity of this algorithm is
(25)O(M2P2N+M3P3+(2K2M(P-1)+2K3+K3)+(2K2(M-1)P+2K3)),
while Wang’s ESPRIT algorithm requires
(26)O(M2P2N+M3P3+(2K2M(P-1)+2K3+K3)+(2K2(M-1)+2K3+3K3)).
According to [19], we derive the CRB for the frequency estimation and angle estimation as
(27)CRB=σ22N{Re[DHΠA⊥D⊙P^sT]},-1
where  D=[(∂a1/∂f1), (∂a2/∂f2)⋯(∂aK/∂fK), (∂a1/∂θ1), (∂a2/∂θ2)×⋯(∂aK/∂θK)],   with ak being the kth column of A(28)P^s=1N∑t=1NststH,ΠA⊥=IPM-A(AHA)-1AH.


## 3.4. Discussion

The performance of this proposed algorithm with nonuniform linear array will be discussed in this section. The distance between themth antenna and the left reference element is dm(m=0,…,M-1), where d0=0.There is an array ofM sensors on which K incident waves impinge (M>K). The direction matrix A′ is
(29)A′=[11⋯1e-jd1γ1e-jd1γ2⋯e-jd1γK⋮⋮⋱⋮e-jdM-1γ1e-jdM-1γ2⋯e-jdM-1γK],
where γk=2πfksin(θk)/c, k=1,2,…,K.Replacing matrixAwith A′ in (3)–(12), we can get the frequency estimation by (15).For a nonuniform linear array, (20) and (21) can be rewritten as follows:
(30)E′c1=[GGΦ1⋮GΦM-2]Π=[E′c11E′c12⋮E′c1(M-1)],E′c2=[GΦ1GΦ2⋮GΦM-1]Π=[E′c21E′c22⋮E′c2(M-1)],
where
(31)Φj=diag{e-jdjγ1,e-jdjγ2,…,e-jdjγK}∈CK×K,j=1,2,…,M-1,E′c1j and E′c2j represent the (j-1)P+1 to jPth rows of matrix E′c1 and E′c2, respectively.We define(32)Cj=Ec1j+Ec2j,j=1,2,…,M-1.
Then, the DOA estimation can be obtained via
(33)θ^k=1M-1∑j=1M-1arcsin(c2πfk(dj-dj-1)angle(ξkj)),k=1,2,…,K,
where ξkj is the kth diagonal element of matrix Cj.


## 3.5. The Advantages of the Proposed Algorithm

The proposed algorithm has the following advantages.(1)
The proposed algorithm has a better angle estimation performance than Wang’s ESPRIT algorithm.(2)
The proposed algorithm can obtain automatically paired parameter estimation, while Wang’s ESPRIT algorithm requires additional pairing.(3)
The proposed algorithm has a comparative computational complexity in contrast to Wang’s ESPRIT algorithm.(4)
This proposed algorithm also suit for nonuniform linear arrays.


## 4. Simulation Results

We present Monte Carlo simulations that are used to assess the angle and frequency estimation performance of this algorithm. The number of Monte Carlo trials is 1000. Note thatM is the number of antennas; P is the number of the delays; N is the number of snapshots; K is the number of the sources.Define(34)RMSE=(1K)*∑k=1K(11000)∑m=11000[amk-a0k]2,
where amk is the kth estimated angle/frequency, and a0k is the kth perfect angle/frequency.Simulation 1
The performance of Wang’s ESPRIT algorithm and this proposed algorithm is investigated.M=12, K=3, P=3, and N=400 in this simulation. Their DOAs are 10°, 20°, and 30°, and their carrier frequencies are 500kHz, 700kHz, and 900kHz. Figures 1 and 2 show the performance of Wang’s ESPRIT algorithm and this proposed algorithm with SNR=0dB. From Figures 1 and 2, we find the proposed algorithm can obtain automatically paired parameter estimation, while Wang’s ESPRIT algorithm cannot, so it works well.Figure 2
Angle-frequency scatter of proposed method, SNR = 0 dB.Simulation 2
We compare this proposed algorithm with Wang’s ESPRIT algorithm and CRB. From Figures3 and 4, we find that this proposed algorithm has better angle estimation performance than Wang’s ESPRIT algorithm and has the same frequency estimation accuracy. Figure 5 to Figure 7 give more angle estimation performance comparison results between these two methods at different M, K, P, and N. Since the frequency estimation performance of the proposed is same with that of Wang’s ESPRIT, it is not necessary to plot it from Figures 5, 6, and 7.Figure 3
Angle estimation performance comparison atM = 12, K = 3, P = 3, and N=100.Figure 4
Frequency estimation performance comparison atM = 12, K = 3, P = 3, and N=100.Figure 5
Angle estimation performance comparison atM = 12, K = 4, P = 3, and N=400.Figure 6
Angle estimation performance comparison atM = 16, K = 3, P = 3, and N=400.Figure 7
Angle estimation performance comparison atM = 12, K = 3, P = 4, and N=400.Simulation 3
This proposed algorithm performance under different snapshotsN is investigated in this simulation. M=12, K=3, and P=3are used in this simulation. Figure 8 shows the angle-frequency estimation performance under different N. We find that the angle-frequency estimation performance of this algorithm is improved with N increasing.Angle-frequency estimation with different snapshotN.
(a)
Angle estimation
(b)
Frequency estimationSimulation 4
The performance of this algorithm under different source numberK is investigated in the simulation. M=12, P=3, and N=400 are used in this simulation. The source number K is set to 2, 3, and 4. This proposed algorithm has different performance under different source numbers, as shown in Figure 9. From Figure 9, we find that angle and frequency estimation performance of this algorithm degrades with the increasing of the source number K.Angle-frequency estimation with different sources.
(a)
Angle estimation
(b)
Frequency estimationSimulation 5
The performance of this algorithm under different antenna numberM is investigated in the simulation. K=3, P=3, and N=400 are used in this simulation. The antenna number M is set to 8, 12, and 16. This proposed algorithm has different performance under different antenna numbers, as shown in Figure 10. From Figure 10, we find that angle and frequency estimation performance of this algorithm is improved with M increasing.Angle-frequency estimation with different antennas.
(a)
Angle estimation
(b)
Frequency estimationSimulation 6
The performance of this algorithm under different delay numberP is investigated in the simulation. M=12, K=3, and N=400 are used in this simulation. The delay number P is set to 2, 3, 4, and 5. This proposed algorithm has different performance under different delay numbers, as shown in Figure 11. From Figure 11, we find that angle and frequency estimation performance of this algorithm is improved with P increasing.Angle-frequency estimation with different delay number.
(a)
Angle estimation
(b)
Frequency estimationSimulation 7
The performance of this proposed algorithm with nonuniform linear arrays is investigated.M=8, K=3, P=3, and N=400 are used in this simulation. Their DOAs are 10∘, 20∘, and 30∘, and their carrier frequencies are 500kHz, 700kHz, and 900kHz. Figure 12, shows this proposed algorithm with nonuniform linear arrays with SNR=0dB. From Figure 12, we find the proposed algorithm can also work well in the case of nonuniform linear arrays.Figure 12
Angle-frequency scatter for nonuniform linear arrays, SNR = 10 dB.


## 5. Conclusion

This paper has presented an improved joint angle-frequency estimation method, which has better angle estimation performance than Wang’s ESPRIT algorithm and has the same frequency estimation accuracy. The computational complexity of this proposed algorithm is comparative in contrast to Wang’s ESPRIT algorithm. Since the DOA and frequency estimations suffer from the same permutation ambiguity, this novel method can obtain automatically paired DOA and frequency. This advantage is more obvious when the input SNR is below 0 dB. Furthermore, the proposed algorithm can also work well in the case of nonuniform linear arrays.



---

*Source: 309269-2012-08-27.xml*