Recurrent Algorithm for TDOA Localization in Sensor Networks

ABSTRACT:

Using the mathematical apparatus of the extended Kalman Filter, the recurrent algorithm of the passive location in sensor networks - based on the Time Difference of Arrival method in case of correlated errors of measurements - is developed. The initial estimates of Radio Frequency Sources coordinates and the correlation matrix of the vector estimation are determined based on the method of the least squares in case of 3 difference measurement distances. Efficiency analysis of recurrent adaptive algorithm and its comparison with the quadratic correction one are performed by statistical modeling. A comparison of them with the lower limit of the Cramer-Rao was carried out. The implementation of the recurrent adaptive algorithm requires 2.7 times less computational cost than the quadratic correction one.

KEYWORDS:
Passive location; Time Difference of Arrival method; Extended Kalman Filter; Recurrent adaptive algorithm; Sensor network

INTRODUCTION

The problem of passive position determination of Radio Frequency Sources (RFS) is widely met in the monitoring of the surrounding space, disaster management, in intelligent transport and security systems. Currently, sensor networks are used for its solution (Rullan-Lara et al. 2013; Amar and Leus 2010).

One of the main approaches of passive position determina- tion of RFS is based on the application of the Time Difference of Arrival method (TDOA), which uses the time difference of reception of signals received by the various sensors and the network reference sensor. This method has a significant advantage in the ease of implementation, being widely used in practice ITU (2014).

The accuracy of determining the RFS coordinates based on the TDOA depends on the errors from measuring time signal reception sensors of the sensor network. The errors from the determination of the difference in signal reception times are correlated, because they contain the error in the reference sensor measurement.

In the known methods for determining the RFS coordinates based on the TDOA (Buzuverov 2008), the coordinate calculation is performed after receiving the measurements from all sensors. In this study, based on the mathematical apparatus of the extended Kalman Filter, the algorithm is developed, which - after the formulation of the initial conditions based on the measurement of time for receiving signals from 4 sensors - allows to recurrently specify the location of RFS as the measurement proceeds from the other sensors. The developed algorithm also evaluates the time error from receiving the signal by the reference sensor measurement that allows considering the synthesized algorithm as adaptive.

FORMULATION OF THE PROBLEM

The sensors network transducers have the coordinates (x_i^S, y_i^S), i = $\stackrel{─}{0,\mathit{n}}$. The position of the RFS is characterized by a point with coordinates (x,y). When determining the coordinates of RFS on the x-y plane, the sensor network should consist of no less than 4 sensors. Figure 1 shows block diagrams of sensor networks on the x-y plane, consisting of 9 (n = 8) and 4 (n = 3) sensors.

Figure 1
Configuration of the sensor network with (a) 9 and (b) 4 sensors

When using TDOA, the time difference in the reception of signals between the sensors i = $\stackrel{─}{1,\mathit{n}}$ and the reference probes is measured:

(1) ${\Delta }_{\mathit{i}0}={\mathit{t}}_{\mathit{i}}-{\mathit{t}}_0+{\mathit{n}}_{\mathit{i}}-{\mathit{n}}_0,\mathit{i}=\stackrel{─}{1,\mathit{n},}$

where: t_i is the time of signal reception of the i − m sensor; t₀ is the time of signal reception by the reference sensor; n_i is the uncorrelated error of time measurements of signal by the reception of the i − m sensor (Buzuverov 2008) with dispersion σ_n², i = $\stackrel{─}{0,\mathit{n}}$.

Because of the transformation in Eq. 1, the TDOA equations for network can be represented as:

(2) $\begin{array}{c}{\mathit{d}}_{\mathit{i}0}=\mathit{c}\left({\mathit{t}}_{\mathit{i}}+{\mathit{n}}_{\mathit{i}}-\mathit{t}\right)-\mathit{c}\left({\mathit{t}}_0+{\mathit{n}}_0-\mathit{t}\right)={\mathit{R}}_{\mathit{i}}-{\mathit{R}}_0+\hfill \\ {\mathit{v}}_{\mathit{i}}-{\mathit{v}}_0={\mathit{R}}_{\mathit{i}}-{\mathit{R}}_0+{\mathit{v}}_{\mathit{i}0},\mathit{i}=\stackrel{─}{1,\mathit{n}},\hfill \end{array}$

where: R_i is the distance between the i − m sensor and the RFS, i = $\stackrel{─}{1,\mathit{n}}$; R₀ is the distance between the reference sensor and the RFS; d_i0 is the measured difference of distances i = $\stackrel{─}{1,\mathit{n}}$; c is the spreading speed of electromagnetic waves; t is the moment of RFS signal emission; v_i0 is the error of difference measurement between the distances:

(3) ${\mathit{v}}_{\mathit{i}0}={\mathit{v}}_{\mathit{i}}-{\mathit{v}}_0,\mathit{i}=\stackrel{─}{1,\mathit{n}};$

where: v_i are the uncorrelated random variables having the meaning of the distance measurement errors between the RFS and the network sensors with a dispersion σ_v² = c² σ_n², i = $\stackrel{─}{0,\mathit{n}}$.

The difference in distance R_i - R₀ is determined by the formula

(4) ${\mathit{R}}_{\mathit{i}}-{\mathit{R}}_0=\sqrt{{\left(\mathit{x}-{\mathit{x}}_{\mathit{i}}^{\mathit{S}}\right)}^2+{\left(\mathit{y}-{\mathit{y}}_{\mathit{i}}^{\mathit{S}}\right)}^2}-\sqrt{{\mathit{x}}^2+{\mathit{y}}^2}$

In Eq. 4, the coordinates of the reference sensor were relied to be 0. The errors of difference measurement between the distances v_i0, i = $\stackrel{─}{1,\mathit{n}}$, are correlated because they contain the reference sensor measurement error v₀.

The presence of correlated errors makes it difficult to use traditional recurrent target coordinates evaluation algorithms. This difficulty can be avoided through the introduction of v₀ into the state vector of the estimated parameters.

It is believed that, during the difference measurement of the distances between the RFS and the network sensors, its coordinates do not change. Synthesizing a recurrent algorithm is required, and, after the formation of the initial conditions based on the measurement of signals receiving time from 4 sensors, it allows to specify recurrently the location RFS in process of receipt of measurements from the other sensors and to estimate the error of the reference sensor measurement.

DEVELOPMENT

The coordinates of the RFS position (x_k, y_k) should be assessed; the measurement time of distance difference between the RFS and the network sensors, as well as its coordinates, do not change. The error in measurements of the reference sensor from Eq. 2 does not change as well in 1 measurement cycle, so the equation describing the dynamics of the estimated parameters has the form

(5) ${\mathit{u}}_{\mathit{k}}={\mathit{u}}_{\mathit{k}-1},$

where: u_k = (x_k, y_k, v_0k)^T is the state vector including position coordinates of the RFS and the measurement error of the reference sensor v_0k on the current time step k. The index k characterizes the incomings sequence of measured differences in distances.

The measurement equation describing the measured k ^th distance difference considering Eqs. 2, 4 and 5 can be viewed as

(6) ${\mathit{r}}_{\mathit{k}}=\mathit{h}\left({\mathit{u}}_{\mathit{k}}\right)+{\mathit{v}}_{\mathit{k}},$

where: h(u_k) is a non-linear function described by the non-linear measurement expression:

(7) $\mathit{h}\left({\mathit{u}}_{\mathit{k}}\right)=\sqrt{{\left({\mathit{x}}_{\mathit{k}}-{\mathit{x}}_{\mathit{i}}^{\mathit{S}}\right)}^2+{\left({\mathit{y}}_{\mathit{k}}-{\mathit{y}}_{\mathit{i}}^{\mathit{S}}\right)}^2}-\sqrt{{\mathit{x}}_{\mathit{k}}^2+{\mathit{y}}_{\mathit{k}}^2}-{\mathit{v}}_{0\mathit{k}}$

where: v_k is the measurement error of the sensor with coordinates x_i, y_i.

Using the model (Eqs. 6 and 7), a recurrent estimation algorithm of the state vector can be obtained based on the extended Kalman Filter (Welch and Bishop 2006), being described by

(8) ${\mathit{K}}_{\mathit{k}}={\hat{\mathit{P}}}_{\mathit{k}-1}·\frac{\partial {\mathit{h}}^{\mathit{T}}\left({\hat{\mathit{u}}}_{\mathit{k}-1}\right)}{\partial {\mathit{u}}_{\mathit{k}}}{\left[\frac{\partial \mathit{h}\left({\hat{\mathit{u}}}_{\mathit{k}-1}\right)}{\partial {\mathit{u}}_{\mathit{k}}}{\hat{\mathit{P}}}_{\mathit{k}-1}\frac{\partial {\mathit{h}}^{\mathit{T}}\left({\hat{\mathit{u}}}_{\mathit{k}-1}\right)}{\partial {\mathit{u}}_{\mathit{k}}}+{\mathit{\sigma }}_{\mathit{v}}^2\right]}^{-1}$

(9) ${\hat{\mathit{u}}}_{\mathit{k}}=\hat{\mathit{u}}{}_{\mathit{k}-1}+{\mathit{K}}_{\mathit{k}}\left[{\mathit{r}}_{\mathit{k}}-\mathit{h}\left({\hat{\mathit{u}}}_{\mathit{k}-1}\right)\right]$

(10) ${\hat{\mathit{P}}}_{\mathit{k}}={\hat{\mathit{P}}}_{\mathit{k}-1}-{\mathit{K}}_{\mathit{k}}\frac{\partial \mathit{h}\left({\hat{\mathit{u}}}_{\mathit{k}-1}\right)}{\partial {\mathit{u}}_{\mathit{k}}}{\hat{\mathit{P}}}_{\mathit{k}-1},$

where: û_k is the estimation of the state vector u_k at the k^th step; P̂_k is the correlation matrix of the state vector of the estimation error u_k at the k^th step; K_k is the coefficient of the filter gain:

(11) $\begin{array}{c}\frac{\partial \mathit{h}\left({\hat{\mathit{u}}}_{\mathit{k}-1}\right)}{\partial {\mathit{u}}_{\mathit{k}}}=\left[\frac{{\hat{\mathit{x}}}_{\mathit{k}-1}-{\mathit{x}}_{\mathit{i}}^{\mathit{S}}}{\sqrt{{\left({\hat{\mathit{x}}}_{\mathit{k}-1}-{\mathit{x}}_{\mathit{i}}^{\mathit{S}}\right)}^2+{\left({\hat{\mathit{y}}}_{\mathit{k}-1}-{\mathit{y}}_{\mathit{i}}^{\mathit{S}}\right)}^2}}\right.-\hfill \\ \frac{{\hat{\mathit{x}}}_{\mathit{k}-1}}{\sqrt{{\hat{\mathit{x}}}_{\mathit{k}-1}^2+{\hat{\mathit{y}}}_{\mathit{k}-1}^2}};\frac{{\hat{\mathit{y}}}_{\mathit{k}-1}-{\mathit{y}}_{\mathit{i}}^{\mathit{S}}}{\sqrt{{\left({\hat{\mathit{x}}}_{\mathit{k}-1}^2-{\mathit{x}}_{\mathit{i}}^{\mathit{S}}\right)}^2+{\left({\hat{\mathit{y}}}_{\mathit{k}-1}-{\mathit{y}}_{\mathit{i}}^{\mathit{S}}\right)}^2}}-\hfill \\ \left.\frac{{\hat{\mathit{y}}}_{\mathit{k}-1}}{\sqrt{{\hat{\mathit{x}}}_{\mathit{k}-1}^2+{\hat{\mathit{y}}}_{\mathit{k}-1}^2}};-1\right]\hfill \end{array}$

The resulting algorithm (Eqs. 8 - 10) is non-linear and belongs to the class of the adaptive ones because it along with the estimation of the RFS the coordinates, an estimation of the unknown error v₀ is determined. The initial conditions of û₀ and P̂₀ must be set for the implementation of adaptive filtering. The initial evaluation of the vector û₀^T is û₀^T = (x̂₀, ŷ₀, 0) The initial estimates of RFS coordinates x̂₀, ŷ₀ are determined based on the method of least squares (LS) in case of 3 difference measurement distances (Amar and Leus 2010) using the Eq. 12:

(12) $\mathit{\omega }=0.5{\left({\mathit{A}}^{\mathit{T}}{\Sigma }^{-1}\mathit{A}\right)}^{-1}{\mathit{A}}^{\mathit{T}}{\Sigma }^{-1}\mathit{b}$

where: ω^T = (x̂₀, ŷ₀, R̂₀) is the vector consisting of the RFS assessment;

(13) $\begin{array}{c}\mathit{A}=\left|\begin{array}{ccc}{\mathit{x}}_1^{\mathit{S}}& {\mathit{y}}_1^{\mathit{S}}& {\mathit{d}}_{10}\\ {\mathit{x}}_2^{\mathit{S}}& {\mathit{y}}_2^{\mathit{S}}& {\mathit{d}}_{20}\\ {\mathit{x}}_3^{\mathit{S}}& {\mathit{y}}_3^{\mathit{S}}& {\mathit{d}}_{30}\end{array}\right|;\mathit{b}=\left|\begin{array}{ccc}{\mathit{x}}_1^{\mathit{S}}+& {\mathit{y}}_1^{\mathit{S}}-& {\mathit{d}}_{10}^2\\ {\mathit{x}}_2^{\mathit{S}}+& {\mathit{y}}_2^{\mathit{S}}-& {\mathit{d}}_{20}^2\\ {\mathit{x}}_3^{\mathit{S}}+& {\mathit{y}}_3^{\mathit{S}}-& {\mathit{d}}_{30}^2\end{array}\right|;\hfill \\ \Sigma =\left[\begin{array}{ccc}2{\mathit{\sigma }}_{\mathit{v}}^2& {\mathit{\sigma }}_{\mathit{v}}^2& {\mathit{\sigma }}_{\mathit{v}}^2\\ {\mathit{\sigma }}_{\mathit{v}}^2& 2{\mathit{\sigma }}_{\mathit{v}}^2& {\mathit{\sigma }}_{\mathit{v}}^2\\ {\mathit{\sigma }}_{\mathit{v}}^2& {\mathit{\sigma }}_{\mathit{v}}^2& 2{\mathit{\sigma }}_{\mathit{v}}^2\end{array}\right].\hfill \end{array}$

The correlation matrix Ω of the vector estimation error ω is defined (Amar and Leus 2010) by

(14) $\Omega ={\left({\left(\stackrel{-0}{\mathit{A}}\right)}^{\mathit{T}}\mathit{H}^{-1}{}_{\mathit{A}}{}^{-0}\right)}^{-1}$

where: H is the matrix ​​determined by the formula H = B∑B; B = diag {R̂₁, R̂₂, R̂₃}; R̂_j = √(x̂₀ - x̂₀^S)² + (ŷ₀ - ŷ₀^S)²; j = $\stackrel{─}{1,3}$.

The initial correlation matrix P̂₀ has the block form

(15) ${\hat{\mathit{P}}}_0=\left[\begin{array}{cc}{\Omega }_{2\times 2}& 0\\ 0& {\mathit{\sigma }}_{\mathit{v}}^2\end{array}\right]$

where Ω_2×2 is determi 56'[p-;;;;;;n 56734r5 ned based on Ω, by deleting the third row and the third column.

After the formation of the initial conditions based on the time measurements to receive signals from the 4 sensors, the synthesized algorithm (Eqs. 8 - 10) allows to recurrently specify, at each step k, the location of the RFS as the measurement proceeds from the other sensors k = $\stackrel{─}{1,\mathit{n}-3}$.

THE EFFECTIVENESS OF THE ALGORITHM ANALYSIS

The efficiency analysis of recurrent adaptive algorithm (Eqs. 8 - 10) and its comparison with the quadratic correction one is performed by statistical modeling. The quadratic correction algorithm provides the highest accuracy among those considered in Amar and Leus (2010). It consists of 2 stages: (i) evaluate the RFS coordinates, which depends on the distance R_0, substitute in the initial functionality then re-solve the linear optimization problem; (ii)-adjust the solution for a quadratic relation.

The modeling of algorithms is performed for the sensor network configuration (Fig. 1a), which consists of 9 sensors to determine the RFS position at the coordinates: S0(0; 0), S1(0; 20), S2(20$\sqrt{2}$; 20$\sqrt{2}$), S3(20; 0), S4(20$\sqrt{2}$; -20$\sqrt{2}$), S5(0; −20), S6(-20$\sqrt{2}$; -20$\sqrt{2}$), S7(−20; 0) and S8(-20$\sqrt{2}$; 20$\sqrt{2}$). Figure 1b shows a sensor network with a minimal number of sensors to form the initial adaptive filtering conditions with the coordinates: S0(0; 0), S1(20$\sqrt{2}$; 20$\sqrt{2}$), S2(20$\sqrt{2}$; -20$\sqrt{2}$) and S3(−20; 0). The RFS is placed on a circle with a radius of 100 km relative to the reference sensor D₀. The root mean square (RMS) of the measurement error is s_v = 30 m. As an indicator of the efficiency, the circular standard deviation was used, $\hat{\mathit{\sigma }}=\sqrt{{\hat{\mathit{\sigma }}}_{\mathit{x}}^2+{\hat{\mathit{\sigma }}}_{\mathit{y}}^2}$.

Figure 2 shows the dependence of the actual ${\hat{\mathit{\sigma }}}_{\mathit{LS}}^{\mathit{MC}}$ (curve 1) circular RMS of the RFS position estimation error with 4 sensors (Fig. 1b) obtained by Monte Carlo simulation using Eq. 12, which corresponds to the initial conditions of the adaptive filter. This figure also shows the dependence of the theoretical (curve 2) circular RMS of the RFS position estimation error, which is calculated based on relevant elements of the correlation matrix of estimation errors determined using Eq. 14. RMS values ${\hat{\mathit{\sigma }}}_{\mathit{LS}}^{\mathit{MC}}$ ​range from 2 to 8 km.

Figure 2
Circular RMS of the RFS position estimation error for the initial conditions

Figure 3 shows the dependence of the actual ${\hat{\mathit{\sigma }}}_{\mathit{FK}}^{\mathit{MC}}$ (curve 1) and theoretical $\hat{\mathit{\sigma }}$_FK (curve 2) circular RMS of the RFS position estimation error with 9 sensors (Fig. 1) for the recurrent algorithm. The received actual and theoretical RMS values correspond well, indicating the proper operation of the algorithm. MSE values ${\hat{\mathit{\sigma }}}_{\mathit{FK}}^{\mathit{MC}}$ range from 1.3 to 1.9 km. The application of recurrent adaptive algorithm reduces the circular RMS of the RFS location estimation error in 1.5 - 4.2 times (Fig. 3).

Figure 3
Circular RMS of the RFS position estimation error for recurrent adaptive algorithm

Figure 4 shows the dependence of the actual ${\hat{\mathit{\sigma }}}_{\mathit{QC}}^{\mathit{MC}}$ (curve 1; where QC is the quadratic correction) and theoretical $\hat{\mathit{\sigma }}$_FK (curve 2) circular RMS of the RFS location estimation error for the quadratic correction algorithm.

Figure 4
Circular RMS of the RFS position estimation error for the quadratic correction algorithm

Figure 5 shows a circular RMS of the RFS location estimation error, which corresponds to the lower limit of the Cramér-Rao bound and characterizes the potential accuracy of the possible RFS coordinates. The values of circular RMS error ​of the RFS positioning of the recurrent adaptive and quadratic correction algorithms, positioned close to the corresponding values ​of circular RMS of the Cramér-Rao bound lower limit, indicate their high efficiency.

Figure 5
Circular RMS of the RFS position estimation error of the Cramér-Rao bound lower limit

Figure 6 shows the change of the actual (curve 1) and theoretical (curve 2) circular RMS error estimation of the RFS locations using recurrent adaptive algorithm for fixed RFS coordinates (x, y). The use of recurrent adaptive algorithm reduces the circular RMS of the RFS location estimation error in 2.5 times relative to the initial conditions.

Figure 6
Circular RMS dynamics of the RFS locations estimation error under sequential data arrival

It is interesting to compare the calculation of the costs required in the implementation of the recurrent adaptive and quadratic correction algorithms. They can be assessed by determining the number of required multiplications (divisions), since they are performed from 100 to 150 times slower than addition and subtraction. For the example with the implementation of the recurrent adaptive algorithm, 461 multiplications are required and, for the quadratic correction one, this value is 1,246. Thus, the application of the developed algorithm reduces the computational cost by 2.7 times.

CONCLUSIONS

Using the mathematical apparatus of Kalman filtering, the algorithm is developed, which, after the formation of the initial conditions, is based on measurements of the receiving time of the signals from 4 sensors. This allows specifying the location of the recurrent RFS as the measurement proceeds from the other sensors. It belongs to the class of adaptive algorithms, because, along with the estimation of the RFS coordinates, it determines the estimation of unknown error of reference sensor measurement. According to the simulation results, the use of recurrent algorithm can reduce the circular RMS of the RFS location estimation error by 1.5 - 4.2 times, providing characteristics similar to the potentially achievable.

The implementation of the recurrent adaptive algorithm requires 2.7 times less computational cost than the quadratic correction one. The resulting algorithm can also be easily extended to the case of RFS filtering trajectory at which its motion parameters are estimated (Chiang et al. 2012).

REFERENCES

• International Telecommunication Union (2015) Comparison of methods of determining the geographic location of the signal source based on the time difference of arrival and angle of arrival of the signal Report ITU-R SM.2211-1. Geneva: ITU.
• Amar A, Leus G (2010) A reference-free time difference of arrival source localization using a passive sensor array. Proceedings of the IEEE Sensor Array and Multichannel Signal Processing Workshop; Jerusalem, Israel.
• Buzuverov GV (2008) Algorithms for passive location in a distributed network of sensors for the time difference of arrival method. In: Buzuverov GV, Gerasimov OI, editors. Information-measuring and control systems. Moscow: Radio Engineering.
• Chiang CT, Tseng PH, Feng KT (2012) Hybrid Unified Kalman Tracking Algorithms for heterogeneous wireless location systems. IEEE Trans Veh Technol 61(2):702-715. doi: 10.1109/TVT.2011.2180939.
  » https://doi.org/10.1109/TVT.2011.2180939
• Rullan-Lara JL, Sanahuja G, Lozano R, Salazar S, Garcia-Hernandez R, Ruz-Hernandez JA (2013) Indoor localization of a quadrotor based on WSN: a real-time application. Int J Adv Robotic Sy 10(1):1-9. doi: 10.5772/53748.
  » https://doi.org/10.5772/53748
• Welch G, Bishop G (2006) An introduction to the Kalman Filter. University of North Carolina at Chapel Hill, Department of Computer Science; [accessed May 24, 2016]. https://www.cs.unc.edu/~welch/media/pdf/kalman_intro.pdf
  » https://www.cs.unc.edu/~welch/media/pdf/kalman_intro.pdf

Publication Dates

• Publication in this collection
  02 Oct 2017
• Date of issue
  Oct-Dec 2017

History

• Received
  21 June 2016
• Accepted
  30 Mar 2017