# TOA and DOA Estimation for IR-UWB Signals: An Efficient Algorithm by Improved Root-MUSIC

Authors: Yang Li; Jia Yu; Liyuan He
Journal: Wireless Communications and Mobile Computing
 (2021)

---

## Abstract

The automatically paired time of arrival (TOA) and direction of arrival (DOA) can be jointly estimated via a high-precision multidimensional spectral peak search- (SPS-) based multiple signal classification (MUSIC) algorithm in the impulse radio ultrawideband (IR-UWB) positioning system, while heavy computational burden is required. To tackle this issue, we propose an improved root-MUSIC algorithm for joint TOA and DOA estimation. After modelling the frequency domain form of the received signal, the algorithm first uses the signal subspace to establish the relationship between the two antennas. Then, the MUSIC spatial spectrum function is reconstructed with this relation, which enables it to offer a spectrum function in regard to the one-dimensional (1D) parameter of time delay. For further reducing the complexity, the TOA estimates of one antenna are obtained via 1D polynomial root finding instead of SPS, and the TOA estimates of the other antenna can be calculated by the established relationship. Finally, the DOA estimation can be achieved with the estimated TOAs. Due to the relationship between two antennas with signal subspace, the parameters estimated by the proposed algorithm are autopaired. Numerical simulations substantiate the superiority of the proposed algorithm.


---

## Body

## 1. Introduction

Impulse radio ultrawideband (IR-UWB) technique, which yields the advantages of a high data rate, low power consumption, low computational complexity, low power spectrum, and great immunity to multipath fading, has been widely applied in various fields, e.g., radar, location, and wireless communications [1–5]. An IR-UWB signal has great potential in precise positioning applications due to its extremely narrow pulse width [6]. Time of arrival (TOA) estimation-based ranging technology plays an important role in the IR-UWB positioning system, as it can make full use of the extremely high time resolution of the IR-UWB signal. And the acquisition of direction of arrival (DOA) is of benefit to reduce the number of nodes in the positioning system. Therefore, the joint TOA and DOA estimation has attracted extensive interest in recent years.For TOA estimation, the parameter estimation method based on time domain signal processing mainly includes a coherent detection method using pulse template matched filter [7], incoherent TOA estimation algorithm based on energy detection [8–11], and two-step TOA estimation algorithm based on energy detection and coherent detection [12]. These algorithms obtain the TOA estimation with the arrival time of the direct path (DP) component. Unfortunately, due to the multipath effect of the UWB signal, especially in non-line-of-sight (NLOS) conditions, the resolution of the aforementioned methods declines because DP is no longer the strongest path (SP) in general. To address this issue, various TOA estimation approaches with superresolution based on frequency domain processing were proposed, such as the minimum norm spectrum estimation algorithm [13], multiple signal classification algorithm [14], and propagation operator method (PM) [15]. These algorithms first model the frequency domain channel impulse response and subsequently employ spectral peak search (SPS) to obtain TOA estimation, which leads to the increase of system complexity and the deterioration of real-time performance.Numerous DOA estimation approaches have been developed for DOA estimation for UWB signals. A kind of data model similar to the data model of narrowband DOA estimation was presented in [16], which processes the received signal in the frequency domain and then conducts the MUSIC algorithm to achieve DOA estimation. In [17], an iterative quadratic maximum likelihood DOA estimation method is proposed, and a beamforming DOA estimation algorithm with frequency invariant in the frequency domain based on beamspace was studied in [18]. Considering the fact that the UWB signal has high time resolution, extracting DOA information from TOA can lead to high DOA estimation accuracy; many researchers have turned to the study of joint TOA and DOA estimation algorithms [19–22]. A low-complexity matrix pencil algorithm was proposed in [19, 20], whereas the estimation performance is susceptible to noise since only one snapshot data is used, and extra parameter pairing of the estimated TOAs is required. In [21], the coarse TOA estimation is obtained by energy estimation and minimum distance criterion; then, the refined TOAs are estimated by the threshold comparison method based on the signal power delay spectrum, and finally, the DOA estimation is achieved by minimum variance unbiased estimation with the estimated TOAs of each antenna. A three-step joint TOA and DOA estimation algorithm was studied in [22], where the conventional threshold correlation approach is employed to obtain the time delay for the signal to reach each antenna, and then, the least mean square algorithm is employed to estimate the TOA and the time difference of arrival (TDOA) jointly. Finally, more accurate DOAs can be obtained by improving the accuracy of the TDOA. In [23], a successive PM algorithm was proposed to reduce the computational burden of the conventional 2D-PM algorithm, which first obtains coarse initial TOA estimation with low complexity and then performs two one-dimensional (1D) local SPS with the initial TOA estimation. The two-dimensional (2D) MUSIC algorithm is able to deal with the TOA pairing problem in [19, 20], which can estimate automatically paired TOAs with high resolution, whereas the 2D total SPS leads to a tremendous amount of computation. A root-MUSIC-based algorithm was presented in [24] to tackle this issue, where the closed-form solutions of TOA estimation can be obtained based on the polynomial root principle, and then, the DOAs can be calculated with the estimated TOA. However, since the received signals of the two antennas are processed separately, the cross-correlation information between the two antennas is neglected and extra pairing of TOA estimates is involved.In this paper, the issue of joint TOA and DOA estimation in the IR-UWB system is investigated and a subspace-based improved MUSIC algorithm is proposed. Specifically, the proposed algorithm first models the received signals in the frequency domain and uses the signal subspace to establish the relationship between the two antennas. Furthermore, the MUSIC spatial spectrum function is reconstructed with this relationship, in which the two-dimensional parameter estimation is transformed into one-dimensional parameter estimation, which can be solved via polynomial root finding with remarkably reduced complexity. The proposed algorithm requires no spectral peak search and can directly give the closed-form solutions of the estimated parameters. Owing to the established relationship with signal subspace, the algorithm can estimate autopaired TOAs. Simulation results demonstrate that the proposed algorithm has significantly better parameter estimation performance than the matrix pencil, PM, and ESPRIT algorithms, and its estimation performance is close to that of the 2D-MUSIC and root-MUSIC algorithms.The main contributions of our work are summarized as follows:(1)
We propose an improved root-MUSIC algorithm for joint TOA and DOA estimation in the IR-UWB positioning system. By constructing and exploiting the relationship between the received signal of the two antennas, the 2D SPS is transformed into 1D SPS and the polynomial root finding is performed to further reduce complexity(2)
The proposed algorithm can fully use the signal subspace and noise subspace to obtain autopaired TOA estimates(3)
The proposed algorithm possesses close estimation performance to the MUSIC and root-MUSIC algorithms and outperforms the matrix pencil, PM, and ESPRIT algorithmsThe rest of this paper is organized as follows. The data model is introduced in Section2, and Section 3 presents the proposed algorithm. In Section 4, we analyze the complexity of the proposed algorithm and list the advantages of the algorithm. Section 5 exhibits the simulation results, and the conclusions are drawn in Section 6.Notation.
Vectors (matrices) are represented by bold lower (upper) case letters.⋅T,⋅H, and ⋅−1 denote the transpose, conjugate-transpose, and inverse operation, respectively. ⋅+ represents the pseudoinverse operation. IN×N is an N×N identity matrix. ⋅F is the Frobenius norm.


## 2. Data Model

Assume that the UWB positioning system uses the second derivative of the Gaussian pulse as the UWB transmission signal. The transmission signal adopts Direct Sequence Binary Phase Shift Keying (DS-BPSK) modulation and can be expressed by [25, 26]
(1)st=∑i=−∞+∞∑n=0Nc−1bicnpt−iTs−nTc,where Nc represents the pulse repetition number of a single binary data symbol, bi∈−1,+1 is the modulated binary data symbol sequence, and cn∈−1,+1 represents a pseudorandom sequence used to implement multiple access communication. Ts is the period of a binary data symbol, and Tc is the pulse repetition period. pt denotes the second derivative of the Gaussian pulse and can be represented by
(2)pt=e−2πt2/Γ21−4πt2Γ2,where Γ is the pulse forming factor related to the pulse width.UWB signals generally propagate in a multipath environment. According to the well-known Saleh-Valenzuela (S-V) model [27], we assume that a pulse of the transmitted signal generates multiple multipath components after passing through the channel, and these multipath components arrive at the receiver in the form of clusters. Suppose that the signal generates K clusters through the channel, and each cluster has L multipaths; the channel impulse response model of the kth cluster of the UWB channel is
(3)hkt=∑l=1Lαlkejϕlkδt−τlk,where αlk stands for the channel attenuation coefficient in the kth cluster of the lth path and obeys Rayleigh distribution; the phase ϕlk is a random variable uniformly distributed on 0,2π. δt denotes the Dirac function. τlk represents the channel delay of the lth path in the kth cluster. In other words, it can be assumed that the multipath delay of the k-cluster received signal is the same, that is, τlk=τl. Let αlkejϕlk=βlk denote the amplitude of random complex fading; then, (3) can be rewritten as
(4)hkt=∑l=1Lβlkδt−τl.The time domain form in thekth cluster signal received by the UWB system can be also written as
(5)ykt=st∗hkt+wkt,where wkt is the additive white Gaussian noise and ∗ represents convolution.Converting the received signal into frequency domain form, we can get [23, 25]
(6)Ykω=SωHkω+Wkω=∑l=1LβlkSωe−jωτl+Wkω,where Ykω,Sω,Hkω, and Wkω stand for the Fourier transform of ykt,st,hkt, andwkt, respectively.Sample the signal in (6) with sampling interval Δω=2π/NN>L, the discrete frequency domain form of the sampled received signal can be obtained from
(7)Ykωn=SωnHkωn+Wkωn=∑l=1LβlkSωne−jωnτl+Wkωn,where ωn=nΔω and n=0,1,⋯,N−1. (7) can also be represented by
(8)yk=SEτβk+wk,where yk=Ykω0,⋯,YkωN−1T denotes equal interval sampling of N points in the frequency domain of ykt. S=diagSω0,⋯,SωN−1 is a diagonal matrix, and the diagonal elements are the N-point frequency domain equal interval sampling values of st. Eτ=e1,⋯,eL∈ℂN×L represents the delay matrix, and el=1,e−jΔωτl,⋯,e−jN−1ΔωτlT. The coefficients of kth cluster channel complex fading are included in βk=β1k,⋯,βLkT. wk=Wkω0,⋯,WkωN−1T is the sampling vector in the frequency domain of noise.


## 3. The Proposed Joint TOA and DOA Estimation Algorithm

The UWB system model is shown in Figure1, where only two antennas are required. Assuming that the UWB signals are in the far field of the antenna array, the signals incident to the array can be regarded as a beam of parallel waves.Figure 1
Antenna array structure.LetY1 and Y2 be the received signals in the frequency domain of the two antennas, respectively. Let τ=τ1,τ2,⋯,τL and ς=ς1,ς2,⋯,ςL denote the vectors of TOAs corresponding to the two antennas, τl and ςl stands for the TOAs associated with the lth l=1,2,⋯,L path. Then, according to the signal model in Section 2, Y1 and Y2 can be written as
(9)Y1=SEτΒ+W1,(10)Y2=SEςΒ+W2,where Β=β1,β2,⋯,βK∈ℂL×K is a complex fading coefficient of the channel. W1=w11,⋯,w1K∈ℂN×K and W2=w21,⋯,w2K∈ℂN×K is the noise samples received by antenna 1 and antenna 2. Eτ and Eς are the delay matrix which can be represented by
(11)Eτ=11⋯1e−jΔωτ1e−jΔωτ2⋯e−jΔωτL⋮⋮⋱⋮e−jN−1Δωτ1e−jN−1Δωτ2⋯e−jN−1ΔωτL,(12)Eς=11⋯1e−jΔως1e−jΔως2⋯e−jΔωςL⋮⋮⋱⋮e−jN−1Δως1e−jN−1Δως2⋯e−jN−1ΔωςL.Denote as the difference of the TOAs of thelth path. As can be concluded from Figure 1, Δτ^l can be obtained from
(13)Δτ^l=dsinθlc,where d represents the distance between the two antennas, θl is the DOA of the lth path, and c denotes the speed of light. From (12), the DOA estimates can be calculated by
(14)θ^l=arcsinΔτ^lcd,l=1,2,⋯,L.It can be seen from (13) that TOAs of two antennas should be estimated first to achieve DOA estimation.


### 3.1. Relationship Conductiond

Before estimating TOA, we first use the signal subspace to establish the relationship between the two antennas. Construct the overall output matrixZ∈ℂ2N×K with the received signal matrix in the frequency domain of two antennas Y1∈ℂN×K and Y2∈ℂN×K as
(15)Z=Y1Y2=SEτSEςΒ+W1W2.Let,(16)Aτ,ς=A1A2=SEτSEς∈ℂ2N×L,W=W1W2∈ℂ2N×K,and Aτ,ς can be specified by [23]
(17)Aτ,ς=Sω0Sω0⋯Sω0Sω1e−jΔωτ1Sω1e−jΔωτ2⋯Sω1e−jΔωτL⋮⋮⋱⋮SωN−1e−jN−1Δωτ1SωN−1e−jN−1Δωτ2⋯SωN−1e−jN−1ΔωτLSω0Sω0⋯Sω0Sω1e−jΔως1Sω1e−jΔως2⋯Sω1e−jΔωςL⋮⋮⋱⋮SωN−1e−jN−1Δως1SωN−1e−jN−1Δως2⋯SωN−1e−jN−1ΔωςL=aτ1,ς1,aτ2,ς2,⋯,aτL,ςL,where aτl,ςl=Sω0,Sω1e−jΔωτl,⋯,SωN−1e−jN−1Δωτl,Sω0,Sω1e−jΔωςl,⋯,SωN−1e−jN−1ΔωςlT.Then, (15) can be represented by
(18)Z=Aτ,ςΒ+W.The covariance matrix ofZ can be approximately calculated by R^=ZZH/K. By performing eigendecomposition, R^ can be decomposed as
(19)R^=U^sD^sU^sH+U^ND^NU^NH,where U^s∈ℂ2N×L is the signal subspace [28] spanned by the eigenvectors corresponding to the largest L eigenvalues and U^N∈ℂ2N×2N−L is the noise subspace [29] consisting of the eigenvectors corresponding to the remaining N−L eigenvalues. According to subspace theory, in noise-free condition, the signal subspace Us can be represented by
(20)Us=A1A2T=SEτSEςT,where T∈ℂL×L is a nonsingular matrix. Then, we partition Us into two parts as follows:
(21)Us=Us1Us2=A1TA2T,where Us1∈ℂN×L,Us2∈ℂN×L denote the signal subspace of antenna 1 and antenna 2, respectively.AsUs1=A1T and Us2=A2T, define a transformation matrix H as
(22)H=Us2Us1+,where ⋅+ represents the pseudoinverse operation. As T is full rank, then H=A2TT−1A1+=A2A1+ and the relationship between A1 and A2 can be represented by
(23)A2=HA1.It can be deduced from (22) and (23) that the relation between A1 and A2 can be obtained with the signal subspace Us of the overall output, which can be utilized to simplify the two-dimensional parameter estimation to one-dimensional estimation.


### 3.2. TOA Estimation Based on Root-MUSIC Algorithm

According to the orthogonal relationship between each column vector ofA and noise subspace UN, the 2D-MUSIC spatial spectrum function can be constructed as
(24)P2D‐MUSICτ,ς=1aτ,ςHU^NU^NHaτ,ς=1a1τa2ςHU^NU^NHa1τa2ς,where aτ,ς, a1τ, and a2ς are the column vectors of Aτ,ς, A1, and A2, respectively. By locating the largest L peaks of P2D‐MUSICτ,ς, τ and ς that correspond to the peaks are exactly to be estimated, whereas expensive computational cost is required due to the two-dimensional global spectral peak search. In addition, the polynomial root-finding technique is difficult to be implemented directly since a two-dimensional parameter to be estimated is involved. To tackle this issue, we use the relationship established between A1 and A2 to transform the two-dimensional SPS into one-dimensional SPS, and then, the polynomial root-finding technique is exploited to replace SPS, which can reduce the computational complexity significantly.According to (23), we can derive the relation a2ς=Ha1τ; then, (24) can be reconstructed as
(25)P2D‐MUSICτ,ς=1a1τHa1τHU^NU^NHa1τHa1τ=1a1τHQQHa1τ,where Q=IN,HHU^N, a1τ=Sbτ, and bτ=1,e−jΔωτ,⋯,e−jN−1ΔωτT represent the column vector of delay matrix Eτ. As can be seen from (25), only 1D SPS for τ is involved to obtain the TOA of the first antenna, where the complexity can be reduced substantially compared with 2D SPS.Subsequently, the polynomial root-finding technique is employed to further reduce complexity. Denotez=e−jΔωτ, then we have bτ=bz=1,z,⋯,zN−1T as well as a1τ=a1z=Sbz, the delay estimation of (25) can be thereby transformed into the following polynomial rooting problem:
(26)a1T1zQQHa1z=0.TheL roots z^1,z^2,⋯,z^L of (26) nearest to the unit circle are exactly required, and zl=e−jΔωτl. Accordingly, the TOA estimates of antenna 1 can be calculated by
(27)τ^l=−anglez^lΔω,l=1,2,⋯,L.Once the rootsz^l=e−jΔωτ∧l are obtained, a^1τl=Sb^zl=S1,z∧l,⋯,z∧lN−1T can be calculated and subsequently a^2ςl can be estimated by
(28)a^2ςl=Ha^1τl.By usingA^2ς=a^2ς1,a^2ς2,⋯,a^2ςL, E^ς can be estimated via
(29)E^ς=S−1A^2ς.According to (12), define Ea and Eb as matrices composed of the first and the last N−1 rows of Eς, respectively, and we can get Eb=EaΦ, where Φ=diage−jΔως1,e−jΔως2,⋯,e−jΔωςL. Correspondingly, define E^a and E^b as the matrices composed of the first and the last N−1 rows of E^ς in (29), the TOAs of antenna 2 can be obtained by
(30)ς^l=−angleEa+EbΔω,l=1,2,⋯,L,where ⋅+ represents pseudoinverse operation. It is worth noting that the TOA estimates τ^l and ς^l are automatically paired, as the relationship between A1 and A2 is exploited in (28) with H.


### 3.3. The Procedure of the Proposed Algorithm

The major steps of the proposed improved root-MUSIC algorithm in the IR-UWB system are as follows:(1)
CalculateR^ of the total received frequency domain signal Z in (15), then obtain U^s and U^N with eigenvalue decomposition of R^(2)
ComputeH and reconstruct the MUSIC spatial spectral function according to (22)–(25)(3)
Construct the polynomial in (26) and conduct polynomial rooting to find the L largest roots z^1,z^2,⋯,z^L to estimate TOA estimates corresponding to antenna 1 τ^ll=1,2,⋯,L via (27)(4)
Calculatea^2ςl and E^ς via (28) and (29)(5)
Obtain the TOA estimates corresponding to antenna 2τ^l through (30)(6)
CalculateΔτ^l=ς^l−τ^l and estimate the DOAs θ^l according to (14)Remark 1.
It is assumed that the multipath raysL are known in this paper. If L is unknown, it can be estimated through the method in [30].Remark 2.
We assume that the number of rays in each cluster is the same, since it is generally believed that the channel remains stable and the rate of change of the channel is very slow during the observation time.


## 3.1. Relationship Conductiond

Before estimating TOA, we first use the signal subspace to establish the relationship between the two antennas. Construct the overall output matrixZ∈ℂ2N×K with the received signal matrix in the frequency domain of two antennas Y1∈ℂN×K and Y2∈ℂN×K as
(15)Z=Y1Y2=SEτSEςΒ+W1W2.Let,(16)Aτ,ς=A1A2=SEτSEς∈ℂ2N×L,W=W1W2∈ℂ2N×K,and Aτ,ς can be specified by [23]
(17)Aτ,ς=Sω0Sω0⋯Sω0Sω1e−jΔωτ1Sω1e−jΔωτ2⋯Sω1e−jΔωτL⋮⋮⋱⋮SωN−1e−jN−1Δωτ1SωN−1e−jN−1Δωτ2⋯SωN−1e−jN−1ΔωτLSω0Sω0⋯Sω0Sω1e−jΔως1Sω1e−jΔως2⋯Sω1e−jΔωςL⋮⋮⋱⋮SωN−1e−jN−1Δως1SωN−1e−jN−1Δως2⋯SωN−1e−jN−1ΔωςL=aτ1,ς1,aτ2,ς2,⋯,aτL,ςL,where aτl,ςl=Sω0,Sω1e−jΔωτl,⋯,SωN−1e−jN−1Δωτl,Sω0,Sω1e−jΔωςl,⋯,SωN−1e−jN−1ΔωςlT.Then, (15) can be represented by
(18)Z=Aτ,ςΒ+W.The covariance matrix ofZ can be approximately calculated by R^=ZZH/K. By performing eigendecomposition, R^ can be decomposed as
(19)R^=U^sD^sU^sH+U^ND^NU^NH,where U^s∈ℂ2N×L is the signal subspace [28] spanned by the eigenvectors corresponding to the largest L eigenvalues and U^N∈ℂ2N×2N−L is the noise subspace [29] consisting of the eigenvectors corresponding to the remaining N−L eigenvalues. According to subspace theory, in noise-free condition, the signal subspace Us can be represented by
(20)Us=A1A2T=SEτSEςT,where T∈ℂL×L is a nonsingular matrix. Then, we partition Us into two parts as follows:
(21)Us=Us1Us2=A1TA2T,where Us1∈ℂN×L,Us2∈ℂN×L denote the signal subspace of antenna 1 and antenna 2, respectively.AsUs1=A1T and Us2=A2T, define a transformation matrix H as
(22)H=Us2Us1+,where ⋅+ represents the pseudoinverse operation. As T is full rank, then H=A2TT−1A1+=A2A1+ and the relationship between A1 and A2 can be represented by
(23)A2=HA1.It can be deduced from (22) and (23) that the relation between A1 and A2 can be obtained with the signal subspace Us of the overall output, which can be utilized to simplify the two-dimensional parameter estimation to one-dimensional estimation.


## 3.2. TOA Estimation Based on Root-MUSIC Algorithm

According to the orthogonal relationship between each column vector ofA and noise subspace UN, the 2D-MUSIC spatial spectrum function can be constructed as
(24)P2D‐MUSICτ,ς=1aτ,ςHU^NU^NHaτ,ς=1a1τa2ςHU^NU^NHa1τa2ς,where aτ,ς, a1τ, and a2ς are the column vectors of Aτ,ς, A1, and A2, respectively. By locating the largest L peaks of P2D‐MUSICτ,ς, τ and ς that correspond to the peaks are exactly to be estimated, whereas expensive computational cost is required due to the two-dimensional global spectral peak search. In addition, the polynomial root-finding technique is difficult to be implemented directly since a two-dimensional parameter to be estimated is involved. To tackle this issue, we use the relationship established between A1 and A2 to transform the two-dimensional SPS into one-dimensional SPS, and then, the polynomial root-finding technique is exploited to replace SPS, which can reduce the computational complexity significantly.According to (23), we can derive the relation a2ς=Ha1τ; then, (24) can be reconstructed as
(25)P2D‐MUSICτ,ς=1a1τHa1τHU^NU^NHa1τHa1τ=1a1τHQQHa1τ,where Q=IN,HHU^N, a1τ=Sbτ, and bτ=1,e−jΔωτ,⋯,e−jN−1ΔωτT represent the column vector of delay matrix Eτ. As can be seen from (25), only 1D SPS for τ is involved to obtain the TOA of the first antenna, where the complexity can be reduced substantially compared with 2D SPS.Subsequently, the polynomial root-finding technique is employed to further reduce complexity. Denotez=e−jΔωτ, then we have bτ=bz=1,z,⋯,zN−1T as well as a1τ=a1z=Sbz, the delay estimation of (25) can be thereby transformed into the following polynomial rooting problem:
(26)a1T1zQQHa1z=0.TheL roots z^1,z^2,⋯,z^L of (26) nearest to the unit circle are exactly required, and zl=e−jΔωτl. Accordingly, the TOA estimates of antenna 1 can be calculated by
(27)τ^l=−anglez^lΔω,l=1,2,⋯,L.Once the rootsz^l=e−jΔωτ∧l are obtained, a^1τl=Sb^zl=S1,z∧l,⋯,z∧lN−1T can be calculated and subsequently a^2ςl can be estimated by
(28)a^2ςl=Ha^1τl.By usingA^2ς=a^2ς1,a^2ς2,⋯,a^2ςL, E^ς can be estimated via
(29)E^ς=S−1A^2ς.According to (12), define Ea and Eb as matrices composed of the first and the last N−1 rows of Eς, respectively, and we can get Eb=EaΦ, where Φ=diage−jΔως1,e−jΔως2,⋯,e−jΔωςL. Correspondingly, define E^a and E^b as the matrices composed of the first and the last N−1 rows of E^ς in (29), the TOAs of antenna 2 can be obtained by
(30)ς^l=−angleEa+EbΔω,l=1,2,⋯,L,where ⋅+ represents pseudoinverse operation. It is worth noting that the TOA estimates τ^l and ς^l are automatically paired, as the relationship between A1 and A2 is exploited in (28) with H.


## 3.3. The Procedure of the Proposed Algorithm

The major steps of the proposed improved root-MUSIC algorithm in the IR-UWB system are as follows:(1)
CalculateR^ of the total received frequency domain signal Z in (15), then obtain U^s and U^N with eigenvalue decomposition of R^(2)
ComputeH and reconstruct the MUSIC spatial spectral function according to (22)–(25)(3)
Construct the polynomial in (26) and conduct polynomial rooting to find the L largest roots z^1,z^2,⋯,z^L to estimate TOA estimates corresponding to antenna 1 τ^ll=1,2,⋯,L via (27)(4)
Calculatea^2ςl and E^ς via (28) and (29)(5)
Obtain the TOA estimates corresponding to antenna 2τ^l through (30)(6)
CalculateΔτ^l=ς^l−τ^l and estimate the DOAs θ^l according to (14)Remark 1.
It is assumed that the multipath raysL are known in this paper. If L is unknown, it can be estimated through the method in [30].Remark 2.
We assume that the number of rays in each cluster is the same, since it is generally believed that the channel remains stable and the rate of change of the channel is very slow during the observation time.


## 4. Performance Analysis

### 4.1. Complexity Analysis

The complexity of the proposed algorithm is mainly caused by the transformation of the signals into the frequency domain, calculation of covariance, eigenvalue decomposition of covariance, and polynomial rooting. Specifically, the complexity of transforming the transmitted and received signals into the frequency domain model isON2+2KN2. Calculating covariance R^ needs O4KN2, and the complexity of eigenvalue decomposition is O8N3. Calculating the transformation matrix H=Us2Us1+ requires ON3+2LN2+L2N, and construct Q needs O2N22N−L. The procedure of polynomial rooting costs ON−13. Consequently, the total complexity of the proposed algorithm is ON−13+13N3+6K+1N2+L2N, while the 2D-MUSIC algorithm costs O8N3+6K+1N2+m22N+12N−L, where m is the peak search times. And the root-MUSIC algorithm [24] requires O10N3+3+4L!+8KN2+4L!−2LL!−2L+6N−LL!−2L−4.Figure2 shows the complexity of the above methods, where m=2000, K=100, and L=3. It can be concluded from Figure 2 that the proposed algorithm requires the lowest computational complexity. And the complexity of the polynomial rooting-based algorithm is much lower than that of the SPS-based algorithm.Figure 2
Computational complexity of different algorithms.


### 4.2. Advantages of the Proposed Algorithm

The proposed algorithm has the following advantages:(1)
The proposed algorithm can achieve automatically paired two-dimensional TOA estimation(2)
The proposed algorithm has an estimation performance close to that of the 2D-MUSIC algorithm, but remarkably lower computational burden is required, since the 2D SPS is transformed into 1D polynomial rooting(3)
The proposed algorithm possesses approximate estimation performance compared with the root-MUSIC algorithm and outperforms the matrix pencil [19, 20], ESPRIT, and PM [15] algorithms


## 4.1. Complexity Analysis

The complexity of the proposed algorithm is mainly caused by the transformation of the signals into the frequency domain, calculation of covariance, eigenvalue decomposition of covariance, and polynomial rooting. Specifically, the complexity of transforming the transmitted and received signals into the frequency domain model isON2+2KN2. Calculating covariance R^ needs O4KN2, and the complexity of eigenvalue decomposition is O8N3. Calculating the transformation matrix H=Us2Us1+ requires ON3+2LN2+L2N, and construct Q needs O2N22N−L. The procedure of polynomial rooting costs ON−13. Consequently, the total complexity of the proposed algorithm is ON−13+13N3+6K+1N2+L2N, while the 2D-MUSIC algorithm costs O8N3+6K+1N2+m22N+12N−L, where m is the peak search times. And the root-MUSIC algorithm [24] requires O10N3+3+4L!+8KN2+4L!−2LL!−2L+6N−LL!−2L−4.Figure2 shows the complexity of the above methods, where m=2000, K=100, and L=3. It can be concluded from Figure 2 that the proposed algorithm requires the lowest computational complexity. And the complexity of the polynomial rooting-based algorithm is much lower than that of the SPS-based algorithm.Figure 2
Computational complexity of different algorithms.


## 4.2. Advantages of the Proposed Algorithm

The proposed algorithm has the following advantages:(1)
The proposed algorithm can achieve automatically paired two-dimensional TOA estimation(2)
The proposed algorithm has an estimation performance close to that of the 2D-MUSIC algorithm, but remarkably lower computational burden is required, since the 2D SPS is transformed into 1D polynomial rooting(3)
The proposed algorithm possesses approximate estimation performance compared with the root-MUSIC algorithm and outperforms the matrix pencil [19, 20], ESPRIT, and PM [15] algorithms


## 5. Simulation Results

In this section, simulations are provided to evaluate the estimation performance of the proposed algorithm. The signal to noise ratio (SNR) is defined bySNR=10lgytF2/wtF2, where yt represents the received time domain signal, wt is the received noise, and ⋅F stands for Frobenius norm. The root mean squared error (RMSE) is presented for evaluation which is defined as
(31)RMSE=1L∑l=1L1M∑m=1Mχ∧l,m−χl2,where χl and χ^l,m denote the theoretical parameter, i.e., TOA or DOA, and the estimate of the lth path in the mth trial, respectively. M denotes the number of Monte Carlo trials.In the simulations, we assume that the UWB pulse wave function ispt=exp−2πt2/Γ21−4πt2/Γ2, where Γ=0.25ns represents the shaping factor for the pulse. The repetition of every symbol is Nc=5, the chip duration is Tc=2ns, and the symbol duration is Ts=NcTc=10ns. We plot the UWB pulse wave function pt and the transmitted BPSK-UWB signal st in Figures 3 and 4, respectively. We suppose that the UWB channel parameters βlk are known, which are random complex fading amplitudes. There are N=64 frequency samples, K=100 observed clusters, and L=3 rays of BPSK-UWB arriving signals whose TOAs and DOAs are 0.2 ns, 0.3 ns, and 0.4 ns and 10°,30°,and45°, respectively.Remark 3.
Besides the Gaussian pulse, some other types of pulse wave functions for use in IR-UWB communications have been proposed such as Laplacian monopulse waveshape [31] and Hermite pulse shapes [32]. In practice, however, the second derivative of a Gaussian pulse is widely used [33]. The pulse waveform can be changed by regulating the shaping factor Γ.Remark 4.
Because of the motion of people and equipment in and around the building, the UWB channel parametersτlk and βlk are random time-varying functions. However, the rate of the variations is very slow compared to any useful signaling rate that is likely considered, e.g., higher than tens of kbit/s. Thus, these parameters can be treated as virtually time-invariant random variables.Remark 5.
We suppose that each cluster has the same number of rays. Normally, we consider that the variation rate of the channel is slow and the channel remains stable in the observing time, so we can assume that the corresponding delays in each cluster are the same.Remark 6.
We set the number of frequency samples as an integer power of 2. As a result, fast Fourier transform (FFT) can be used, which is more efficient than DFT.Figure 3
UWB pulse waveformpt.Figure 4
Transmitted signalst.We provide the parameter estimation results of the proposed algorithm in Figures5 and 6, where SNR=20dB and M=50. As illustrated in Figures 5 and 6, the proposed algorithm can estimate TOAs and DOAs accurately.Figure 5
TOA estimation of the proposed algorithm (SNR=20dB).Figure 6
DOA estimation of the proposed algorithm (SNR=20dB).Figures7 and 8 compare the parameter estimation performance with different methods, including the proposed algorithm, 2D-MUSIC, root-MUSIC [24], matrix pencil [19, 20], PM [15], and ESPRIT algorithms, where N=64, L=3, K=100, and M=500. It is observed that the proposed algorithm, root-MUSIC algorithm, and 2D-MUSIC algorithm outperform the other algorithms in parameter estimation. The matrix pencil algorithm performs worst, as this algorithm is proposed for single snapshot data, and its estimation performance is very sensitive to noise. In the case of multiple snapshots, we average the estimated values under each snapshot, while its performance is still poor. The PM algorithm has poor performance at low SNR; this is because the frequency domain division of the received signal and the transmitted signal is used to obtain the frequency domain channel impulse response in [15], which enhances the noise at low SNR. Furthermore, the proposed algorithm possesses the same parameter estimation as the root-MUSIC algorithm which is close to that of the 2D-MUSIC algorithm, while the proposed algorithm has lower computational complexity.Figure 7
RMSE of TOA estimation with different algorithms.Figure 8
RMSE of DOA estimation with different algorithms.The RMSE results with differentL of the proposed algorithm are provided in Figures 9 and 10, where N=64, K=100, and M=500. As Figures 9 and 10 demonstrate, the parameter estimation performance of the proposed algorithm improves with the decrease of L.Figure 9
RMSE of TOA estimation with differentL.Figure 10
RMSE of DOA estimation with differentL.


## 6. Conclusions

In this paper, we have presented a subspace-based improved root-MUSIC algorithm with an IR-UWB system for joint TOA and DOA estimation. The relation between the two received antennas is first established by partitioning the signal subspace, with which the 2D spatial spectral function can be converted into a 1D one. Then, the polynomial root-finding technique is exploited to obtain the TOA of one antenna instead of SPS for further reducing the computational burden, and the TOA of the other antenna can be estimated by utilizing the relation between the two received antennas. Finally, the DOAs can be calculated with the difference of the TOAs. Closed-form solutions to the parameter estimation can be attained, and no extra parameter pairing is required. Simulations demonstrate the validity of the proposed algorithm for the IR-UWB system with regard to computational complexity and parameter estimation performance.



---

*Source: 9960273-2021-05-03.xml*