Keywords

1 Introduction

1.1 Motivation

Locating a moving object in real-time has become an attractive topic in the academic and research community worldwide recently, owing to its great applicability potential in both military and commercial fields [1,2,3,4,5]. Taking advantage of the existing technologies, e.g., terrestrial radio frequency sources, to provide a solution to the real-time object localization problem is strongly encouraged, due to cost minimization. These include for example, extracting time of arrival, received signal strength (RSS), angle of arrival (AoA) information from the received radio signal, or a combination of them [5,6,7,8,9,10,11]. Which approach to use depends mainly on the available hardware. In this work, we employ combined RSS and AoA measurements, because nowadays, practically all devices can measure the RSS information, and the AoA information can be extracted from RSS measurements by using rotating a directional antenna and choosing the direction from which the highest RSS measurement is obtained [12, 13].

1.2 Research Question

Designing real-time localization algorithms is a very difficult task, because of numerous challenges that have to be taken into consideration, such as accuracy, execution time and limited energy resources to name a few. Therefore, the research question of this work can be formulated as:

How to design an efficient (highly accurate and computationally low complex) localization algorithm, robust to network topology and channel characteristics, applicable in real-time?

In order to address the research question, we have begun with the following hypothesis:

An efficient real-time localization algorithm can be developed by linearizing the considered non-linear measurement model and by following Bayesian approach, where prior knowledge is integrated with observations to enhance the estimation accuracy of an estimator.

1.3 Related Work

Works [6,7,8,9,10,11] investigated the classical target localization problem, where the proposed estimators were based on radio measurements only, and prior knowledge was disregarded. The authors in [2, 3, 5] considered a problem of tracking of a moving target, and they combined the available radio measurements with prior knowledge given from the target state model. However, in [2, 3, 5] only pure RSS-based target tracking problem was investigated. In [4], the target tracking problem which employs hybrid, RSS and AoA, measurements was considered. The authors in [4] proposed a Kalman filter (KF) and a particle filter (PF) to solve the tracking problem. Furthermore, they proposed a generalized pattern search method for estimating the path loss exponent (PLE) for each link in every time step.

1.4 Contributions

In this paper, we investigate the RSS-AoA-based target tracking problem for unknown target transmit power. This setting is of practical interest for low-cost systems in which testing and calibration is not the priority. Also, due to sensor’s battery drain during time, the true value transmit power becomes not perfectly known over time. By integrating prior knowledge given by the target transition model and by employing the maximum a posteriori (MAP) approach, we propose a tight approximation of the MAP estimator. Based on a well-known recursive approach, typical for Bayesian methods, we develop an iterative algorithm which updates the mean and covariance of the target state in each time step.

2 Relationship to Smart Systems

Smart systems embrace functions of sensing, actuation, and control with a goal of describing and analyzing an environment, and making decisions based on the available data in a predictive or adaptive manner. A good example of such a system is a wireless sensor network, composed of a number of scattered sensors that cooperate between themselves in order to respond in an attentive, adaptive and active way to the changes in the environment registered by sensors. However, in many applications, the information gathered by the sensors is meaningless if not correlated to accurate location of where the changes are occurring (e.g., a system might be set up to respond locally to changes in sensor data).

Therefore, determining accurate location of sensors is a very important task in forming a smart system. Moreover, accurate localization of objects and people in both indoor and outdoor environments enables new applications in emergency and commercial services that can improve safety and efficiency in everyday life (e.g., smart parking, monitoring of storage conditions and goods, assistance for elderly or people with disabilities) [14].

Another important example where target localization helps building a smart system is video surveillance. It has been used to monitor security sensitive areas such as banks, highways, borders, etc. Owing to rapid advances high-speed network infrastructure and computing power, as well as large capacity storage devices, multi sensor video surveillance systems have been developed recently. However, traditional video outputs that were controlled by human operators became overwhelming both for operators and storage devices, due to the increased number of cameras. Therefore, in order to filter out redundant information and increase the response time to forensic events motivated the development of smart video surveillance systems. Such systems require fast, reliable and robust algorithms for moving object detection, classification, tracking and activity analysis [15].

3 Problem Formulation

Let \( \varvec{x}_{t} = \left[ {x_{\text{x}} , x_{\text{y}} } \right]^{T} \) and \( \varvec{a}_{i} = \left[ {a_{{i{\text{x}}}} ,a_{{i{\text{y}}}} } \right]^{T} \), for \( i = 1, \ldots , N \), denote the unknown location of a moving target at time instant \( t \) and known location of the \( i \)-th static anchor, respectivelyFootnote 1. For simplicity, we assume a constant velocity target state transition model [1,2,3,4,5, 16, 17], i.e.,

$$ \varvec{\theta}_{t} = \varvec{S\theta }_{t - 1} + \varvec{r}_{t} , $$
(1)

where \( \varvec{\theta}_{t} = \left[ {\varvec{x}_{t} , \varvec{v}_{t} } \right]^{T} \) represents the target state at time t (described by its location and velocity, \( \varvec{v}_{t} \), in a 2-dimensional plane) and \( \varvec{r}_{t} \) is the state process noise [1,2,3,4,5]. This process noise is assumed to be zero-mean Gaussian with a covariance matrix Q, i.e., \( \varvec{r}_{t} \sim {\mathcal{N}}(0,\varvec{Q}) \). Covariance \( \varvec{Q} \) and state transition matrix, \( \varvec{S} \), are defined as

$$ \varvec{Q} = q\left[ {\begin{array}{cccc} {\Delta ^{3} /3} & 0 & {\Delta ^{2} /2} & 0 \\ 0 & {\Delta ^{3} /3} & 0 & {\Delta ^{2} /2} \\ {\Delta ^{2} /2} & 0 & \Delta & 0 \\ 0 & {\Delta ^{2} /2} & 0 & \Delta \\ \end{array} } \right], $$
$$ \varvec{ S} = \left[ {\begin{array}{*{20}c} 1 & 0 & \Delta & 0 \\ 0 & 1 & 0 & \Delta \\ 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 1\\ \end{array} } \right], $$

with \( q \) and \( \Delta \) denoting the state process noise intensity and the sampling interval between two consecutive time steps [1, 3, 16, 17].

At each time step, \( t \), the target emits a signal to the anchors which then withdraw the RSS and AoA information from it. Thus, the observation model can be written as

$$ \varvec{z}_{t} = \varvec{h}(\varvec{x}_{t} ) + \varvec{n}_{t} , $$
(2)

where \( \varvec{z}_{t} = \left[ {\varvec{P}_{t}^{T} ,\phi_{t}^{T} } \right] (\varvec{z}_{t} \in {\mathbb{R}}^{2N} ) \) is the observation vector composed of RSS, \( \varvec{P}_{t} = \left[ {P_{i}^{\left( t \right)} } \right]^{T} \), and AoA, \( \phi_{t} = \left[ {\phi_{i}^{\left( t \right)} } \right]^{T} , \) measurements at time instant \( t \). The function \( \varvec{h}(\varvec{x}_{t} ) \) in (2) is defined as \( h_{i} \left( {\varvec{x}_{t} } \right) = P_{0} - 10\gamma \log_{10} \frac{{\left\| {\varvec{x}_{t} - \varvec{a}_{i} } \right\|}}{{d_{0} }} \) for \( i = 1, \ldots , N \) [18], and \( h_{i} \left( {\varvec{x}_{t} } \right) = \tan^{ - 1} \left( {\frac{{x_{\text{y}} - a_{{i{\text{y}}}} }}{{x_{\text{x}} - a_{{i{\text{x}}}} }}} \right) \) for \( i = N + 1, \ldots , 2N \), [19], where \( P_{0} \) is the reference power at a distance \( d_{0} \), and \( \gamma \) is the PLE. The measurement noise is modeled as \( \varvec{n}_{t} \in {\mathcal{N}}(0,\varvec{C}) \), where the noise covariance is defined as \( \varvec{C} = {\text{diag}}\left( {\left\lfloor {\sigma_{{n_{i} }}^{2} , \sigma_{{m_{i} }}^{2} } \right\rfloor } \right) \otimes \varvec{I}_{4} \), with \( \sigma_{{n_{i} }}^{2} \) (dB) and \( \sigma_{{m_{i} }}^{2} \) (rad) being the variances of the RSS and AoA measurement noise, respectively, \( \varvec{I}_{M} \) representing the identity matrix of size \( M \) and symbol \( \otimes \). noting the Kronecker product. Often in practice, testing and calibration are not the priority in order to keep low network implementation costs; hence, the target transmit power, \( P_{T} \), is often not calibrated, i.e., not known. Not knowing \( P_{T} \) corresponds to not knowing \( P_{0} \) in (2) [18, 20].

We employ (1) and (2) to build the marginal posterior probability distribution function (PDF), \( p(\varvec{\theta}_{t} |\varvec{z}_{1:t} ) \), from which we can quantify the confidence we have in the values of the state \( \varvec{\theta}_{t} \) given all the past measurements \( \varvec{z}_{1:t} \). From \( p(\varvec{\theta}_{t} |\varvec{z}_{1:t} ) \), we can obtain an estimate at any desired time step.

Below we show a recursive procedure, typical for Bayesian methods [1,2,3,4,5], for the evaluation of \( p(\varvec{\theta}_{t} |\varvec{z}_{1:t} ) \) at any time instant.

  • Initialization: The marginal posterior PDF at \( t = 0 \) is set to the prior PDF \( p(\varvec{\theta}_{0} ) \) of \( \varvec{\theta}_{0} \).

  • Prediction: According to the state transition model (3), the predictive PDF of the state at \( t \) is found as

    $$ p\left( {\varvec{\theta}_{t} |\varvec{z}_{1:t - 1} } \right) = \mathop \smallint \nolimits p\left( {\varvec{\theta}_{t} |\varvec{\theta}_{t - 1} } \right)p\left( {\varvec{\theta}_{t - 1} |\varvec{z}_{1:t.1} } \right)d\varvec{\theta}_{t - 1} . $$
    (3)
  • Update: By following the Bayes’ rule [1, 13], we have

    $$ p\left( {\varvec{\theta}_{t} |\varvec{z}_{1:t} } \right) = \frac{{p(\varvec{z}_{t} |\varvec{\theta}_{t} )p(\varvec{\theta}_{t} |\varvec{z}_{1:t - 1} )}}{{\mathop \smallint \nolimits p\left( {\varvec{z}_{t} |\varvec{\theta}_{t} } \right)p\left( {\varvec{\theta}_{t} |\varvec{z}_{1:t - 1} } \right)d\varvec{\theta}_{t} }}. $$
    (4)

It is worth noting that the denominator in (4) is just a normalizing constant. In general, the marginal PDF at \( t - 1 \) cannot be calculated analytically, and the integral in (3) cannot be obtained analytically if the state model is non-linear. Therefore, some approximations are required in order to obtain \( p(\varvec{\theta}_{t} |\varvec{z}_{1:t} ) \).

4 The Proposed MAP Algorithm

A state estimate, \( \hat{\varvec{\theta }}_{t|t} \), of \( \varvec{\theta}_{t} \) can be obtained from \( p(\varvec{\theta}_{t} |\varvec{z}_{1:t} ) \) by employing the MAP criteria [13], i.e., by maximizing the marginal posterior PDF

$$ \hat{\varvec{\theta }}_{t|t} = \mathop {\text{arg max}}\limits_{{\varvec{\theta}_{t} }} \,p\left( {\varvec{\theta}_{t} |\varvec{z}_{1:t} } \right) \approx \mathop {\text{arg max}}\limits_{{\varvec{\theta}_{t} }} \,p(\varvec{z}_{t} |\varvec{\theta}_{t} )p(\varvec{\theta}_{t} |\varvec{z}_{1:t - 1} ). $$
(5)

The problem in (5) resembles the maximum likelihood estimator, apart from the existence of the prior PDF. This problem is highly non-convex and its analytical solution cannot be obtained in general. Therefore, we approximate (5) by another estimator whose solution is readily available.

First, for small noise power we can write from the RSS and AoA model respectively

$$ \lambda_{i} \left\| {\varvec{x}_{t} - \varvec{a}_{i} } \right\| \approx \eta d_{0} , {\text{for}}\, i\, = \, 1, \ldots ,N, $$
(6)
$$ \varvec{c}_{i}^{T} \left( {\varvec{x}_{t} - \varvec{a}_{i} } \right) \approx \,0, {\text{for}}\, i\, = \,N + 1, \ldots ,2N, $$
(7)

where, \( \lambda_{i} = 10^{{\frac{{P_{i}^{\left( t \right)} }}{10\gamma }}} \), \( \eta = 10^{{\frac{{P_{0} }}{10\gamma }}} \), and \( \varvec{c}_{i} = \left[ { - \sin \phi_{i}^{\left( t \right)} , \cos \phi_{i}^{\left( t \right)} } \right]^{T} \). Next, by transforming from Cartesian to polar coordinates we express \( \varvec{x}_{t} - \varvec{a}_{i} = r_{i} \varvec{u}_{i} : r_{i} \ge 0, \left\| {\varvec{u}_{i} } \right\| = 1 \), and make use of the available azimuth angle observations to define \( \varvec{u}_{i} = [\cos \phi_{i}^{\left( t \right)} , \sin \phi_{i}^{\left( t \right)} ] \), and rewrite (6) as

$$ \lambda_{i} r_{i} \approx \eta d_{0} . $$

By multiplying the left hand side of the above equation by \( \varvec{u}_{i}^{T} \varvec{u}_{i} \), we arrive at

$$ \lambda_{i} \varvec{u}_{i}^{T} (\varvec{x}_{t} - \varvec{a}_{i} ) \approx \eta d_{0} . $$
(8)

Hence, we can approximate the RSS and AoA measurement models by (8) and (7), respectively, which written in a vector form and applying the weighted least squares (WLS) criterion leads to

$$ \hat{\varvec{y}}_{t|t} = \mathop {{ \arg }\,{ \hbox{min} }}\limits_{{\varvec{y}_{t} = \left[ {\varvec{x}_{t} ,\eta } \right]^{T} }} \left\| {\varvec{W}(\varvec{Ay}_{t} - \varvec{b})} \right\|^{2} , $$
(9)

where \( \varvec{W} = \varvec{I}_{2} \otimes {\text{diag}}(\varvec{w}) \), \( \varvec{w} = \left[ {\sqrt {w_{i} } } \right]^{T} \) and \( w_{i} = 1 - \frac{{P_{i}^{\left( t \right)} }}{{\sum\nolimits_{i = 1}^{N} {P_{i}^{\left( t \right)} } }} \)

$$ \varvec{A} = \left[ {\begin{array}{*{20}c} \vdots & \vdots \\ {\lambda_{i} \varvec{u}_{i}^{T} } & {d_{0} } \\ \vdots & \vdots \\ {\varvec{c}_{i}^{T} } & 0 \\ \vdots & \vdots \\ \end{array} } \right],\varvec{ b} = \left[ {\begin{array}{*{20}c} \vdots \\ {\lambda_{i} \varvec{u}_{i}^{T} \varvec{a}_{i} } \\ \vdots \\ {\varvec{c}_{i}^{T} \varvec{a}_{i} } \\ \vdots \\ \end{array} } \right]. $$

The solution of (9) is given by \( \hat{\varvec{y}}_{t|t} = \left( {\varvec{A}^{T} \varvec{W}^{T} \varvec{WA}} \right)^{ - 1} (\varvec{A}^{T} \varvec{W}^{T} \varvec{b}) \). Once we have an estimate of the target location, we can use it to find the maximum likelihood estimate of \( P_{0} \), \( \hat{P}_{0} \), from the RSS measurement model as

$$ \hat{P}_{0} = \frac{{\sum\nolimits_{i = 1}^{N} {P_{i}^{\left( t \right)} + 10\gamma \,\log_{10} \frac{{\left\| {\varvec{x}_{t} - \varvec{a}_{i} } \right\|}}{{d_{0} }}} }}{N}. $$
(10)

At this point, we can benefit from this estimated value by using it to resume the estimation process as if \( P_{0} \) is known, i.e., we can compute \( \hat{\eta } = 10^{{\frac{{\hat{P}_{0} }}{10\gamma }}} \) in order to get

$$ \widehat{\varvec{A}} = \left[ {\begin{array}{*{20}c} \vdots & \vdots & \vdots \\ {\lambda_{i} \varvec{u}_{i}^{T} } & 0 & 0 \\ \vdots & \vdots & \vdots \\ {\varvec{c}_{i}^{T} } & 0 & 0 \\ \vdots & \vdots & \vdots \\ \end{array} } \right],\varvec{ \hat{b}} = \left[ {\begin{array}{*{20}c} \vdots \\ {\lambda_{i} \varvec{u}_{i}^{T} \varvec{a}_{i}^{T} + \hat{\eta }d_{0} } \\ \vdots \\ {\varvec{c}_{i}^{T} \varvec{a}_{i}^{T} } \\ \vdots \\ \end{array} } \right]. $$

As far as the prior knowledge part in (5) is concern, we can assume that \( p(\varvec{\theta}_{t - 1} |\varvec{z}_{1:t - 1} ) \) has Gaussian distribution with mean \( \hat{\varvec{\theta }}_{t - 1|t - 1} \) and covariance \( \widehat{\varvec{M}}_{t - 1|t - 1} \) [17]. Then, according to (3), we obtain

$$ p\left( {\varvec{\theta}_{t} |\varvec{z}_{1:t - 1} } \right) \approx \frac{1}{k}\exp \left\{ { - \frac{1}{2}\left( {\varvec{\theta}_{t} - \hat{\varvec{\theta }}_{t|t - 1} } \right)^{T} \widehat{\varvec{M}}_{t| - 1}^{ - 1} \left( {\varvec{\theta}_{t} - \hat{\varvec{\theta }}_{t|t - 1} } \right)} \right\}, $$
(11)

where \( k \) is a constant, and \( \hat{\varvec{\theta }}_{t|t - 1} \) and \( \hat{\varvec{M}}_{t|t - 1} \) are respectively the mean and covariance of the one-step predicted state, obtained through (1) as

$$ \begin{aligned} \hat{\varvec{\theta }}_{t|t - 1} & = \varvec{S\hat{\theta }}_{t - 1|t - 1} , \\ \widehat{\varvec{M}}_{t|t - 1} & = \varvec{SM}_{t - 1|t - 1} \varvec{S}^{T} + \varvec{Q}. \\ \end{aligned} $$
(12)

Therefore, (5) can be written as

$$ \begin{aligned} \hat{\varvec{\theta }}_{t|t} & = \mathop {{ \arg }\,{ \hbox{min} }}\limits_{{\varvec{\theta}_{t} }} \left( {\varvec{z}_{t} - \varvec{h}\left( {\varvec{x}_{t} } \right)} \right)^{T} \varvec{C}^{ - 1} \left( {\varvec{z}_{t} - \varvec{h}\left( {\varvec{x}_{t} } \right)} \right) \\ & \quad + \left( {\varvec{\theta}_{t} - \hat{\varvec{\theta }}_{t|t - 1} } \right)^{T} \widehat{\varvec{M}}_{t| - 1}^{ - 1} \left( {\varvec{\theta}_{t} - \hat{\varvec{\theta }}_{t|t - 1} } \right) . \\ \end{aligned} $$
(13)

Similarly as in (9), the problem in (13) can be rewritten as

$$ \hat{\varvec{\theta }}_{t|t} = \mathop {{ \arg }\,{ \hbox{min} }}\limits_{{\varvec{\theta}_{t} }} \left\| {\widehat{\varvec{W}}(\varvec{H\theta }_{t} - \varvec{f})} \right\|^{2} , $$
(14)

where \( \varvec{H} = \left[ {\widehat{\varvec{A}};\widehat{\varvec{M}}_{t|t - 1}^{{ - \frac{1}{2}}} } \right],\varvec{f} = \left[ {\hat{\varvec{b}};\widehat{\varvec{M}}_{t|t - 1}^{{ - \frac{1}{2}}} \hat{\varvec{\theta }}_{t|t - 1} } \right] \), \( \widehat{\varvec{W}} = \varvec{I}_{2} \otimes {\text{diag}}(\hat{\varvec{w}}) \), \( \hat{\varvec{w}} = \left[ {\hat{w}_{i} } \right]^{T} ,\hat{w}_{i} = 1 - \frac{{\hat{d}_{i} }}{{\sum\nolimits_{i = 1}^{N} {\hat{d}_{i} } }} \), \( \hat{d}_{i} = d_{0} 10^{{\frac{{\hat{P}_{0} - P_{i}^{\left( t \right)} }}{10\gamma }}} \). The solution to (14) is given by \( \hat{\varvec{\theta }}_{t|t} = \left( {\varvec{H}^{T} \varvec{W}^{T} \varvec{WH}} \right)^{ - 1} (\varvec{H}^{T} \varvec{W}^{T} \varvec{f}) \).

The step by step proposed MAP algorithm for unknown \( P_{T} \), is outlined in Algorithm 1.

5 Simulation Results

This section validates the performance of the proposed algorithm through computer simulation. All simulations in this work were performed in MATLAB.

5.1 Simulation Set-Up

Two scenarios, in which the target takes sharp maneuvers and a smoother one, are investigated. The state model (1) is used for the target state transition, and the radio measurements are acquired according to (2) at each time instant. Three anchors are fixed at \( [\left[ {70,10} \right]^{T} ,\left[ {40,70} \right]^{T} , \left[ {10,40} \right]^{T} ] \), and a sample is taken every \( \Delta = 1 \) s during \( T = 150 \) s trajectory duration in each Monte Carlo, \( M_{c} = 1000 \), run. The reference distance is set to \( d_{0} = 1 \) m, the reference power to \( P_{0} = - 10 \) dBm, and the PLE to \( \gamma = 3 \), \( \sigma_{{n_{i} }} = 9 \) dB, \( \sigma_{{m_{i} }} = 4\pi /180 \) rad, and \( q = 2.5 \times 10^{ - 3} {\text{m}}^{2} /{\text{s}}^{3} \). The performance metric used here is the root mean square error (RMSE), defined as \( {\text{RMSE}} = \sqrt {\sum\nolimits_{i = 1}^{{M_{c} }} {\frac{{\left\| {\varvec{x}_{i,t} - \hat{\varvec{x}}_{i,t} } \right\|^{2} }}{{M_{c} }}} } \), where \( \hat{\varvec{x}}_{i,t} \) represents the estimate of the true target location, \( \varvec{x}_{i,t} \), in the \( i \)-th Mc run at time instant \( t \).

The performance of the proposed MAP algorithm is compared with a classical approach described by (9), which makes use of the radio measurements only and disregards the prior knowledge, labeled here as WLS.

5.2 Results

Figure 1 illustrates the true target trajectory in the two considered scenarios as well as the estimated trajectories given by the considered approaches. In the first scenario, the initial target’s location was set to \( \left[ {21,20} \right]^{T} \), whereas in the second scenario the starting point was set to \( \left[ {35,15} \right]^{T} \). We can see that the proposed MAP estimator performs considerably better in both scenarios than the WLS one. However, it seems like there is still room for further improvement of the proposed approach, since the results are not particularly smooth. This might be done by better utilization of the prior knowledge, i.e., in the update step of the state covariance matrix, and is left for future work. Figure 2 illustrates the RMSE (m) versus \( t \) (s) performance comparison of the considered approaches. From it, one can see that the proposed algorithm performs well in both scenarios, obtaining an average RMSE below 3 m. Furthermore, the superiority of the Bayesian approach, which combines the observations with prior knowledge, over the classical one, which utilizes the observations exclusively, is clearly observed in every time step.

Fig. 1.
figure 1

The true target trajectory in: (a) the first and (b) the second considered scenario and the estimated ones.

Fig. 2.
figure 2

RMSE versus \( {\text{t}} \)(s) comparison for (a) the first and (b) the second considered scenario.

6 Conclusions

In this work, we have investigated the target tracking problem which makes use of the combined RSS and AoA measurements, in which the target transmit power is considered not known. By resorting to Bayesian approach and relying on the MAP criterion, we proposed a tracking algorithm, which efficiently solves the target tracking in all considered scenarios and significantly outperforms the considered classical approach in every time instant. Although the proposed algorithm offers excellent estimation accuracy, it seems that there is still room for further improvement, which will be the topic of our future work.

This work is a part of our ongoing research, and we have validated the performance of our algorithms by means of simulations only. As a part of our future work, we plan to validate their performance using real measurements from [12, 13].