Source Localization and Tracking for Dynamic Radio Cartography using Directional Antennas

TL;DR

This paper introduces a unified compressive sensing framework combining RSS and DoA for source localization using directional antennas, improving accuracy and efficiency in dynamic radio cartography.

Contribution

It proposes a novel formulation that unifies RSS and DoA in a CS framework, optimizing antenna beams for better localization and spectrum sensing.

Findings

Significant improvement in source localization accuracy.

Enhanced spectrum sensing efficiency with directional antennas.

Effective dynamic RF map estimation over time.

Abstract

Utilization of directional antennas is a promising solution for efficient spectrum sensing and accurate source localization and tracking. Spectrum sensors equipped with directional antennas should constantly scan the space in order to track emitting sources and discover new activities in the area of interest. In this paper, we propose a new formulation that unifies received-signal-strength (RSS) and direction of arrival (DoA) in a compressive sensing (CS) framework. The underlying CS measurement matrix is a function of beamforming vectors of sensors and is referred to as the propagation matrix. Comparing to the omni-directional antenna case, our employed propagation matrix provides more incoherent projections, an essential factor in the compressive sensing theory. Based on the new formulation, we optimize the antenna beams, enhance spectrum sensing efficiency, track active primary users accurately and monitor spectrum activities in an area of interest. In many practical scenarios there is no fusion center to integrate received data from spectrum sensors. We propose the distributed version of our algorithm for such cases. Experimental results show a significant improvement in source localization accuracy, compared with the scenario when sensors are equipped with omni-directional antennas. Applicability of the proposed framework for dynamic radio cartography is shown. Moreover, comparing the estimated dynamic RF map over time with the ground truth demonstrates the effectiveness of our proposed method for accurate signal estimation and recovery.

Figures

Equations

$$a_{n,m,p}=\frac{1}{R_{n,p}^{\eta }}{e}^{-j\frac{2\pi d(m-1)}{{\lambda }_c}sin({\theta }_{n,p})},$$

$$r_n={\mathbf{w}}_n^H{\mathbf{A}}_n\mathbf{s}+{\mathbf{w}}_n^H{\mathbf{z}}_n,$$

$$y_n(t)=$$

$$+$$

$$v_n(t)={\epsilon }_n(t)+m_{v_n},$$

$$\mathbf{y}(t)=\mathbf{\Gamma }(t)\mathbf{x}(t)+{\mathbf{m}}_v+\mathbf{\epsilon }(t),$$

$$<\hat{\mathbf{x}}(t),\hat{\mathbf{W}}(t)>$$

$$\text{s.t.}\mathbf{\Gamma }(t)=h(\mathbf{W}(t))$$

$$\hat{\mathbf{x}}(t)=\underset{\mathbf{x}}{\mathrm{argmin}}\Vert \mathbf{y}(t)-\mathbf{\Gamma }\mathbf{x}-{\mathbf{m}}_v{\Vert }_2^2+\lambda \Vert \mathbf{x}{\Vert }_1.$$

$$\hat{\mathbf{x}}(t)$$

$$\hat{\mathbf{W}}(t+1)$$

$${\mathbf{\Gamma }}^B(t+1)=h(\mathbf{W}(t)).$$

$$\boldsymbol{\Gamma}=\underset{\boldsymbol{\Gamma}}{\text{argmin}}\;\mathbb{E}_{x}\big{\{}\|\boldsymbol{y}(t)-\boldsymbol{m}_{v}^{B}\;-\boldsymbol{\Gamma}{\boldsymbol{x}}\|_{2}^{2}\big{\}}.$$

$$\boldsymbol{\Gamma}=\underset{\boldsymbol{\Gamma}}{\text{argmin}}\;\text{trace}\Big{(}\sum_{n}\boldsymbol{\gamma}_{n}\boldsymbol{\gamma}_{n}^{T}\Big{)}^{-1}.$$

$${\mathbf{u}}_n(t+1)$$

$${\tilde{\mathbf{w}}}_n(t+1)$$

$${\hat{\mathbf{w}}}_n(t+1)$$

$${\tilde{\mathbf{w}}}_n(t+1)$$

$${\hat{\mathbf{w}}}_n(t+1)$$

$${\hat{\mathbf{x}}}^{(n)}(t)=\underset{\mathbf{x}}{\mathrm{argmin}}$$

$$+\lambda \Vert \mathbf{x}{\Vert }_1+\alpha \sum _{i\in {\mathcal{N}}_n}\Vert {\hat{\mathbf{x}}}^{(i)}(t-1)-\mathbf{x}{\Vert }_2^2.$$

$$\{{\mathbf{u}}_j(\tau )\}\text{where,}$$

$$\forall j\in {\mathcal{N}}_n\text{and}\tau =t-B+1,\cdots ,t-1.$$

$$RSS(\mathbf{g})=\sum _{p=1}^P\frac{{\hat{x}}_p(t)}{R_{g,p}^{\eta }}.$$

$$\text{normalized error}(t)=\frac{\Vert \mathbf{x}(t)-\hat{\mathbf{x}}(t){\Vert }_2}{\Vert \mathbf{x}(t){\Vert }_2},$$

Full text

I Introduction

The cognitive radio (CR) paradigm is a promising solution to alleviate today’s spectrum deficiency caused by an increasing demand for ubiquitous wireless access [1, 2, 3, 4]. In the static spectrum allocation strategies, the licensed holders of the spectrum (a.k.a. primary users or PUs) often under-utilize this valuable resource [1]. The CR paradigm allows the unlicensed or secondary users (SUs) to coexist with PUs and to access the spectrum as long as they do not interfere with PUs. The under-utilized spectrum bands that can be used by the SUs are called spectrum holes [5]. One approach to monitoring (sensing) spectrum and identifying the spectrum holes is through the localization of emission sources (PUs). Given the locations of PUs and the signal propagation parameters one can estimate spectrum power of an area of interest. Collaborative source localization techniques are divided into two main categories, i) centralized techniques, and ii) distributed techniques. Centralized techniques are based on the assumption that signals from all SUs are collected at a fusion center (FC), that is tasked with fusing the received signals and localizing the sources. On the other hand, distributed techniques are based on the assumption that there is no FC. In this case, each SU receives signals from its neighbors and performs source localization locally [6].

Source localization can generally be performed by estimating time-difference of arrival (TDoA), received-signal-strength (RSS), or direction of arrival (DoA). The first approach implies synchronization of SUs, which is infeasible in many scenarios. Most of the previous works on spectrum sensing exploit omni-directional antennas for SUs [7, 8]. Thus, DoA cannot be used for source localization in these works. The most popular DoA estimation methods include MUSIC, Capon and ESPRIT [9]. However, the performance of these methods are limited in low SNRs or when the sources are placed close to each other. In [10, 11] the RSS and DoA of a PU is estimated using energy measurements from a sectorized antenna and the performance and theoretical bounds of DoA/RSS estimation are studied. The orientation of a PU with respect to SU’s location is determined in [12, 13] based on the RSS where an SU is equipped with a reconfigurable antenna. The signal received at an antenna array is sparse in the spatial domain and compressive sensing can be employed for DoA estimation using less number of radio frequency (RF) chains [14]. In [15] the problem of DoA estimation using electronically steerable parasitic array radiator (ESPAR) antenna based on compressive sensing is exploited for a sensor equipped with a single RF chain. However, in the present paper directional antennas are employed and their scanning pattern are updated adaptively. Moreover, compressive sensing formulation is employed for improving the accuracy of source parameters estimation.

In the present paper, we assume some sensors (or SUs) are randomly deployed in the area of interest to sense spectrum and the goal is to estimate locations of PUs and generating power spectrum map for the whole area of interest. Each sensor is equipped with an antenna with uniform linear array (ULA). The weighting (beamforming) vector corresponding to each ULA can be adapted/updated to steer/direct the antenna beam and thus localize PUs within the area and track them accurately (if they are mobile). We develop both centralized and distributed collaborative source localization algorithms, in which RSS and DoA information are integrated in a unified formulation. In our formulation, we assume an area of interest is divided into $P$ grid points. PUs (emitting sources) may only reside in these grid point (Figure 1). PUs locations and powers are modeled by a sparse vector, in which a non-zero element indicates the presence of a PU in the corresponding grid point. The received signal at each SU is a superposition of the emitted signals by sources (whose locations and signal powers are unknown). In the centralized setup, we establish a compressive sensing problem at the FC, where the compression operator is modeled in terms of the signal propagation parameters and the sensing patterns of the SUs directional antennas (encompassing ULA structure and DoA information as well as weighting/beamforming vectors). In the distributed setup, each SU shares its received signal with its neighbors (determined by the connectivity graphs), and hence each SU can establish a compressive sensing problem, based on the collected received signals from its neighbors. For both centralized and distributed setups, we provide source localization and tracking algorithms that exploit the inherent sparsity and adaptive beamforming to accomplish accurate source localization and tracking within the area. The main contributions of our paper are summarized as follows:

• •

RSS and DoA information are integrated in a unified formulation based on compressive sensing.

• •

Sensing patterns (or weighting vectors) of directional antennas at SUs are optimized, in order to quickly and accurately discover the PUs’ activities and efficiently track them if they move within the area.

• •

In the distributed setup, neighboring SUs localize sources collaboratively and track them via adapting their weighting vectors and steering/directing their beams.

• •

As a by-product we can estimate the power spectrum map of the area (so-called radio cartography).

Notation: Throughout this paper, vectors and matrices are written as bold lowercase and uppercase letters, respectively. Norm $\ell_{2}$ of a vector is denoted by $\|.\|_{2}$ and $t$ mod $B$ is equal to the remainder of the division $t/B$.

The rest of the paper is organized as follows. Section II states the problem and explains the employed system model. Sections III and IV present our proposed centralized and distributed source localization and tracking algorithms. Section V elaborates how we can leverage on our source estimation findings in the previous sections to obtain the power map of the entire area. Section VI exhibits our experimental results.

II System Model and Problem Statement

Our system model can be viewed as the extension of the system model in [7], for the case when SUs are equipped with directional antennas (instead of omni-directional antennas). Suppose an area of interest is divided into $P$ grid points as shown in Fig. 1 and PU transmitters are assumed to be located in a subset of these grid points, unknown to SUs (i.e., the number and the locations of PU transmitters as well as their transmission power are unknown to SUs). Moreover, there are $N$ SUs (or spectrum sensors) with known locations, where each sensor receives a superposition of the PUs’ signals subject to a zero mean measurement noise. We suppose antenna of each SU is a ULA consisting of $M$ elements with adjacent element spacing of $d$. Assuming $\theta_{n,p}$ is the orientation of grid point $p$ with respect to sensor $n$, the manifold matrix for sensor $n$ is given by $\boldsymbol{A}_{n}=[a_{n,m,p}]$, where

[TABLE]

in which $\lambda_{c}$ is the wave-length corresponding to the operating frequency of each SU. Moreover, $R_{n,p}$ is the distance between grid point $p$ and sensor $n$, and $\eta$ is the path-loss exponent [9]. We note that in general, the received signal at sensor $n$ would be affected by channel fading gains between PU transmitters and sensor $n$. These gains can be obtained by training sequences, which is beyond the scope of this paper. For simplicity, similar to [7, 8], we only consider path loss effect in our signal modeling in (1), to indicate that the received signal will be attenuated inversely proportional to $R_{n,p}$111The propagation model can be extended for fading channels by assuming an uncertainty on the channel gains.. We write the received signal at sensor $n$ as

[TABLE]

where $\boldsymbol{w}_{n}=[w_{n,1}\;\ldots\;w_{n,M}]^{T}$ is the weighting (or beamforming) vector at sensor $n$, and we assume that $\|\boldsymbol{w}_{n}\|_{2}^{2}=C_{1},\forall n$ where $C_{1}$ is the antenna gain. In (2) $\boldsymbol{s}=[s_{1}\;\ldots\;s_{P}]^{T}$ is the transmitted signal vector from PUs, and $\boldsymbol{z}_{n}=[z_{n,1}\;\ldots\;z_{n,M}]^{T}$ is the received noise vector at sensor $n$ corresponding to its ULA, where $z_{n,m}\sim{\cal C}{\cal N}(0,\sigma_{z_{n}}^{2})$ are zero-mean complex Gaussian independent and identically distributed (i.i.d). Fig. 1 illustrates an example of antenna spatial gain (or sensing pattern) of sensors at a single time slot. However, the sensing pattern can be a different one at the next time slot. The beam of antennas in a set of consecutive time slots should (i) cover the whole area in order to discover a new appearing PU and (ii) scan the location of previously estimated active PUs to track their changes.

We assume that the signals of PUs’ transmitters are independent and let $y_{n}(t)=|r_{n}(t)|^{2}=r_{n}(t)r^{*}_{n}(t)$ denote the RSS of sensor $n$ at time $t$. We can write $y_{n}(t)$ as the following:

[TABLE]

where the vector $\small{\boldsymbol{\gamma}_{n}(t)=[\boldsymbol{\gamma}_{1n}(t)\;\;\boldsymbol{\gamma}_{2n}(t)\;\;\ldots\;\boldsymbol{\gamma}_{Pn}(t)]^{T}}$, and $x_{p}(t)=|s_{p}(t)|^{2}$ represents the propagation power at the grid point $p$ and time instant $t$, $\boldsymbol{x}(t)=[x_{1}(t)\;x_{2}(t)\;\ldots\;x_{P}(t)]^{T}$ is the propagation power vector, $\boldsymbol{w}_{n}(t)$ is the beamforming vector that sensor $n$ uses at time $t$, and $v_{n}(t)$ is distributed as an exponential random variable with mean $m_{v_{n}}=\|\boldsymbol{w}_{n}\|_{2}^{2}\sigma^{2}_{z_{n}}=C_{1}\sigma^{2}_{z_{n}}$. Matrix $\boldsymbol{C}_{p}$ is a $P\times P$ matrix with all entries equal to [math] except the diagonal entry of index $(p,p)$. Note that $\gamma_{n,p}$ includes the pathloss effect as well as the beamforming vector $\boldsymbol{w}_{n}$. We can write $v_{n}(t)$ in (II) as

[TABLE]

where $\varepsilon_{n}(t)$ is a zero-mean random variable with variance $m_{v_{n}}^{2}$. We assume that sensor $n$ can adjust its antenna sensing pattern and steer its beam, via optimizing its beamforming vector $\boldsymbol{w}_{n}$. Examining (II) we note that vector $\gamma_{n}$ depends on the optimization variable $\boldsymbol{w}_{n}$. Moreover, vector $\boldsymbol{x}(t)$ is unknown and sparse. Assuming the RSS $y_{n}(t)$ for some $n$ values (in the centralized setup, $n\in\{1,\ldots,N\}$, and for the distributed setup, $n\in{\cal N}_{n}$, where ${\cal N}_{n}$ is the set of neighboring sensors of sensor $n$) are available, our main goal is to exploit the steering capability of SUs’ antennas in order to improve the accuracy and efficiency of spectrum sensing, for the general scenario when the number and the locations of PU transmitters are unknown and time varying. This means that we aim at finding an accurate estimate of the sparse vector $\boldsymbol{x}(t)$ as time changes, as we allow SU sensors to adjust their antenna sensing patterns. By allowing SUs to steer their beams and optimize their beamforming vectors, we enable discovering of new PU signals and tracking those that have been discovered previously. Having an accurate estimate of $\boldsymbol{x}(t)$ we can build the dynamic radio power map of the field, via estimating the RSS at each grid point and time instant $t$. We note that as we solely use RSS for estimating $\boldsymbol{x}(t)$ and the power map, no synchronization is required among SU sensors.

Considering centralized estimation of $\boldsymbol{x}(t)$, let $\boldsymbol{y}(t)=[y_{1}(t)\;y_{2}(t)\;\ldots\;y_{N}(t)]^{T}$ denote the vector containing the RSS of all $N$ sensors. Moreover, let $\boldsymbol{\Gamma}$ denote an $N\times P$ matrix, whose $n$-th row is $\boldsymbol{\gamma}_{n}^{T}$. Matrix $\boldsymbol{\Gamma}$ is referred to as the measurement matrix in the literature of compressive sensing and we refer to it as propagation matrix here. The vector of RSS at the FC can be written as

[TABLE]

where $\boldsymbol{m}_{v}=[m_{v_{1}}\;m_{v_{2}}\;\dots\;m_{v_{N}}]^{T}$ is the vector of averaged RSS and vector $\boldsymbol{\varepsilon}(t)$ is the zero-mean residual of all sensors at time instant $t$ . The centralized estimation of the sparse signal $\boldsymbol{x}(t)$ can be viewed as the solution to the following constrained optimization problem

[TABLE]

where $\boldsymbol{W}=[\boldsymbol{w}_{1}\;\ldots\;\boldsymbol{w}_{N}]$ is an $M\times N$ matrix containing beamforming vectors and $h(\cdot)$ is a function of beamforming matrix $\boldsymbol{W}$ which can be derived using (II). Parameter $C_{0}$ tunes the impact of $\ell_{1}$ regularization and parameter $C_{1}$ restricts antenna gain of SUs. Fig. 1 shows a sensor network with $3$ SUs and $100$ candidate points for PUs where only $3$ of them are active and propagate signals. SUs continuously scan the space to discover signals of new PUs and track the previously detected ones.

Previous works have focused on estimation of variable $\boldsymbol{x}$ subject to a given propagation matrix. The following minimization problem has been proposed for estimating $\boldsymbol{x}(t)$ for each time independently [7, 16]:

[TABLE]

In Problem (7), propagation matrix is a function of distance of SUs and PUs and it can not be optimized. However, in the present work we employ directional antennas and beamforming is needed to track existing PUs and discover any appearing PUs. Next, we will discuss the centralized and distributed solutions of our constrained optimization problem.

III Centralized Source Localization and Tracking

In this section, we address Problem (II). To achieve this, we need to link between the entries of matrix $\boldsymbol{\Gamma}$ and dynamic of network through beamforming matrix $\boldsymbol{W}$, and the last estimate of the propagation power vector $\boldsymbol{x}(t\!-\!1)$. Each SU sweeps a different pattern at each time slot. Assume at time slot $t$ we jointly estimate the propagation power vector $\boldsymbol{x}(t)$, and optimize the beamforming vector $\boldsymbol{w}_{n}(t)$. We exploit $B$ recent measurements from $B$ distinguished sensing patterns. Each $B$ consecutive time slots make a time block of size $B$. Let us break Problem (II) into two alternating subproblems as the following

[TABLE]

The first subproblem estimates $\boldsymbol{x}(t)$ using $B$ recent measurements of $N$ sensors and the second subproblem updates the next beamforming vector for the next time slot. In (8a), \boldsymbol{y}^{B}(t)\!=\!\big{[}\boldsymbol{y}^{T}\!(t\!-\!B\!+\!1)\;\boldsymbol{y}^{T}\!(t\!-\!B\!+\!2)\;\ldots\;\boldsymbol{y}^{T}\!(t)\big{]}^{T}, is a vector containing $NB$ sensed RSS measurements and \boldsymbol{\Gamma}^{B}\!(t)=[\boldsymbol{\Gamma}^{T}\!(t\!-\!B\!+\!1)\;\boldsymbol{\Gamma}^{T}\!(t\!-\!B\!+\!2)\;\ldots\;\boldsymbol{\Gamma}^{T}\!(t)\big{]}^{T} is an $NB\times P$ matrix. Note that in each time slot $N$ new beamforming vectors are updated for $N$ sensors. In other words, $N$ new rows are added to matrix $\boldsymbol{\Gamma}^{B}$ and $N$ obsolete rows are removed. Moreover, in $\boldsymbol{\Gamma}^{B}(t\!+\!1)$ the first $(B\!-\!1)\times N$ rows are shared with $\boldsymbol{\Gamma}^{B}(t)$. However, the last $N$ rows of $\boldsymbol{\Gamma}^{B}(t\!+\!1)$ are subject to optimization. In the rest of this section we explain how to solve this problem.

Solving Problem (8a) w.r.t. $\boldsymbol{x}$ is a straightforward problem which is discussed in the context of sparse regression and compressive sensing. However, Problem (8b) is a challenging problem and the main contribution of this paper is solving (8b) efficiently, i.e., updating beamforming vectors efficiently.

III-A Recovery of sparse sources

Subproblem (8a) is a classic regression problem. Usually number of sensors ($N$) are much smaller than the number of grid points ($P$), which results in a compressive sensing problem. $\ell_{1}$ regularization promotes sparsity that plays a key role for sparse recovery. There have been extensive efforts for solving this problem efficiently. Some of the state-of-the-art algorithms are least absolute shrinkage and selection operator (LASSO) [8] and iterative re-weighted least squares (IRLS) [15, 17]. There are two essential factors in compressive sensing theory (i) sparsity of the underlying process and (ii) incoherent projections from the process [18]. In the next subsection we study the second factor in order to provide a set of incoherent projections of $\boldsymbol{x}$ in (8).

III-B Successive optimization of beamforming vectors

Problem (8b) can be regarded as updating matrix $\boldsymbol{\Gamma}$ in order to project the most information from $\boldsymbol{x}$ into $\boldsymbol{y}$ in each time slot. There are two phases for updating entries of $\boldsymbol{\Gamma}$: phase i) A part of $\boldsymbol{\Gamma}$ is optimized for sensing a general unknown $\boldsymbol{x}$ (discovering), phase ii) A part of $\boldsymbol{\Gamma}$ is allocated and optimized for tracking a known previously estimated $\boldsymbol{x}$ (tracking). In the first phase unseen spectrum activities are discovered. While in the second phase an estimate for $\boldsymbol{x}$ is available and small changes are tracked via optimizing beamforming vectors of sensors.

$\bullet$ Phase 1: Discovery of unknown sources: First we write the problem in terms of optimization of $\boldsymbol{\Gamma}$.

[TABLE]

In the case of Gaussian and independent noise for each sensor it is easy to show that this problem is equal to solving the following problem,

[TABLE]

This problem enforces rows of $\boldsymbol{\Gamma}$ to be linearly independent. In the case of linearly dependent set of rows of $\boldsymbol{\Gamma}$ the cost function in (10) will be infinity. Intuitively, as we measure less correlated projections of an unknown vector we can recover it more accurately222We say two vectors are less correlated when their inner vector product has a smaller absolute value.. Equation (10) indicates that we should regularize Problem (8b) accordingly. On the other hand, Problem (8b) is expressed in terms of the beamforming matrix $\boldsymbol{W}$. According to (II), row space of $\boldsymbol{\Gamma}$ is spanned by $\boldsymbol{w}_{n}^{H}\boldsymbol{A}_{n}$ for all $n$. Problem (10) should be cast in terms of these bases. Let $\boldsymbol{u}_{n}^{T}(t)=\boldsymbol{w}_{n}(t)^{H}\boldsymbol{A}_{n}$ be an auxiliary row vector and matrix $\boldsymbol{U}(t)$ be a matrix whose $N$ rows are $\boldsymbol{u}_{n}^{T}(t)$ for $n=1,\ldots,N$. The beam update rule, inspired by (10), is formulated as

[TABLE]

where $\boldsymbol{U}_{B}$ is the concatenation of all scanned beams in the recent $B-1$ time slots and a new row $\boldsymbol{u}$. In this problem, $\boldsymbol{u}$ is the last row of $\boldsymbol{U}_{B}$ which is the subject of optimization. This problem will be repeated $N$ times to find a new scanning beam for all $N$ sensors. Problem (10) is referred as A-optimal solution in the literature of optimization [19]. However, (11) is modified to use D-optimality by minimizing determinant [19]. D-optimality provides attractive properties which is helpful for efficient optimization. D-optimal optimization is a sub-modular problem which can be solved optimally in a greedy manner. In other words, at time slot $t$ when the beams of $N$ SUs is updated, we can use greedy approach to update each beam independently and add the new found beam to $\boldsymbol{U}_{B}$. In (11), the collection of $\boldsymbol{u}_{n}$ vectors in $\boldsymbol{U}_{B}$ constructs a polygonal in $\mathbb{R}^{P}$ in which each vertex corresponds to a previously scanned beam. The proposed optimization problem finds a new vertex on the unit sphere such that the volume of the new polygonal is maximized. This criterion has many applications in sensor selection, sensor placement, and data reduction [20, 21, 22]. We refer to the solution of this problem as the discovery-based beamforming vectors. These beamforming vectors are uncorrelated to each other and each one tries to discover a new direction which is not covered by the span of other beams. Finally, $\boldsymbol{w}_{n}(t\!+\!1)$ must be computed using estimated $\boldsymbol{u}_{n}(t\!+\!1)$ in (11) as follows.

[TABLE]

The number of grid points for potential PUs ($P$) should be greater than $M$ to have an invertible $\boldsymbol{A}_{n}\boldsymbol{A}_{n}^{T}$. The solution to (12) provides a set of incoherent vectors as the measurement matrix in the compressive sensing formulation. Incoherency of measurement vectors plays a key role in compressive sensing of a general unknown $\boldsymbol{x}$ [18]. However, the current estimation of $\boldsymbol{x}$ gives us valuable information about eventful parts of the area of interest which is discussed next.

$\bullet$ Phase 2: Tracking of sources with a priori knowledge: We can obtain an initial (rough) estimate of propagation power vector $\boldsymbol{x}$ using any randomly chosen set of beamforming vectors. This rough estimate of $\boldsymbol{x}$ is considered as the initial result of source localization and is used for tracking the discovered sources more efficiently at the next time slot. To achieve this, it is sufficient to update $\boldsymbol{w}_{n}$ according to the equations given by

[TABLE]

This update rule performs similar to a matched filter (matched to the received signal) which maximizes the RSS. In each time block of operation, i.e., $B$ time slots, one beam is dedicated to be matched with the dynamic of network. The task of $B\!-\!1$ beams is to discover unseen areas of network to discover any new activity. Once a new activity is detected by estimation of $\boldsymbol{x}$, it will be refined and tracked via the dedicated beam. Assuming $NB\leq M$, then discovery-based beams can be updated via selection from the null space of beams in the recent $B\!-\!1$ time slots. Null space of these beams associates with a space that is not covered yet and the next beam is restricted to be selected from this space. Fig. 2 illustrates a time block of sensing operation. Assume we have a rough estimation of propagation power vector $\boldsymbol{x}$ and at the first time slot a sensing beam is optimized such that the previously detected source can be tracked. The second beam (at the second time slot) is a vector which is orthogonal to $\boldsymbol{w}_{1}(1)$ and at the third time slot $\boldsymbol{w}_{1}(3)$ is optimized to be orthogonal to the two-dimensional subspace spanned by $\boldsymbol{w}_{1}(1)$ and $\boldsymbol{w}_{1}(2)$. In general, we use the recent $B$ beams for estimation in each time. In any $B$ consecutive time slots there is one beam for tracking and $B\!-\!1$ orthogonal/uncorrelated beams333In the case of $NB\leq M$ the optimum beams will be orthogonal and otherwise we can just optimize them to be uncorrelated.. Each of these beams is orthogonal to the recent tracking beam and orthogonal to the rest of beams.

In the case that $NB>M$ after $M$ time slots the null space will be empty and $\boldsymbol{w}(M\!+\!1)$ can not be selected from the null space. In this case we can not employ the solution suggested by null space selection and Problem (11) must be solved for each time. Alg. 1 summarizes the steps of the proposed algorithm. First beamforming vectors are initialized randomly and the corresponding propagation matrix $\boldsymbol{\Gamma}$ is constructed. Then, for each time slot measurements are sensed and the joint $\boldsymbol{x}$ estimation and beamforming problem is solved.

IV Distributed Implementation

In some practical scenarios there is no FC to aggregate measurements and to perform the centralized algorithm. Still, we can leverage on limited communication between neighboring SUs to enable collaborative spectrum sensing in a distributed fashion, using the beam steering capabilities of the antennas. Each sensor processes its received signal separately. As neighboring SUs communicate over inter-sensor links (governed by the connectivity graphs), each sensor learns and collects signals of its neighbors. Consequently SUs can reach a consensus over the result of spectrum sensing and their estimates of propagation power vector $\boldsymbol{x}$. Fig. 3 shows a simple distributed network in which $10$ SUs estimate the locations and the signal powers of $4$ active PUs.

Sensor $n$ is connected to a subset of SUs which is determined by set $\mathcal{N}_{n}$. Distributed version of Problem (8a) is proposed as the following:

[TABLE]

In this problem $\hat{\boldsymbol{x}}^{(n)}$ is the obtained estimation of propagation power vector at the $n^{\text{th}}$ SU. Vector $\boldsymbol{y}^{(n)}$ in the collection of sensing measurements from $\mathcal{N}_{n}$ sensors for recent $B$ time slots and $\boldsymbol{\Gamma}^{(n)}$ is the corresponding propagation matrix. The last term of optimization applies consensus among connected SUs for each estimation. Parameter $\alpha$ regularizes the impact of the consensus term in the main cost function. Alg. 2 shows the steps of the distributed implementation of our proposed algorithm.

Beamforming optimization is similar to the centralized solution. However, each SU only has access to the scanned beams of itself and its neighbors to construct matrix $\boldsymbol{U}^{B}$ in problem (11). We define matrix $\boldsymbol{U}^{(n)}(t)$ for each time which corresponds to the already scanned beams in recent $B$ time slots via sensor $n$ and its neighbors. Mathematically, the set of rows of matrix $\boldsymbol{U}^{(n)}(t)$ can be expressed by,

[TABLE]

The goal is to add a new row such that the criterion in Problem (11) is optimized.

V Radio Cartography

An important side product of estimating $\boldsymbol{x}(t)$ in Alg. 2 and Alg. 1 is generating the RF power map for any arbitrary point in the network. This is also referred to as radio cartography. Consider a point $\boldsymbol{g}$ with coordinate $(g_{1},g_{2})$ in the network. To estimate the RSS at point $(g_{1},g_{2})$, it is sufficient to consider the estimated sources at their corresponding locations and apply the underlying propagation model. Propagation pattern of sources is considered to be omni-directional. The estimated RSS at location $\boldsymbol{g}$ and time instant $t$ can be interpolated as follows:

[TABLE]

where $R_{g,p}$ is the distance between location $g$ and the $p$-th PU. Fig. 4 shows the block diagram of dynamic spectrum cartography via our proposed localization and beamforming algorithm.

VI Experimental Results

In this section, experiment results for synthetic data in a cognitive radio network are presented. Simulations are performed using CVX toolbox in Matlab. The area of interest is $100$m $\times 100$m. The locations of PUs on this grid are unknown and SUs are placed randomly on this grid. Fig. 5 shows the impact of directional antennas on the localization accuracy. There are $8$ active PUs in the network and 10 SUs are sensing the environment. In this figure there is no noise on the RSS and our proposed method reaches almost exact reconstruction after few time slots. Fig. 5a shows the ground truth power map and location of $10$ SUs in the network. The estimated power map using omni-directional antennas is depicted in Fig. 5b.

Employing Directional antennas by randomly chosen beamforming weights provides a significant improvement compared to omni-directional antennas which is shown in Fig. 5c. Our proposed method reaches exact reconstruction after $50$ time slots. Fig. 5d shows the reconstructed power map using our algorithm that is identical to the ground truth.

Fig. 6 exhibits the numerical comparison corresponding to the simulation of Fig. 5. This graph shows the performance of localization in terms of normalized error of estimation as follows:

[TABLE]

in which $\hat{\boldsymbol{x}}(t)$ is the estimated vector that contains location and power of each PU at time $t$. As it can be seen, after $20$ time slots our algorithm retrieves the true locations for PUs.

The effect of additive noise is studied in Fig. 7 over time. As it can be seen, the performance of localization when $\text{SNR}=10$ dB is close to the noiseless regime. The noise is added according to (2). A Monte Carlo simulation with $50$ different noise realization is performed for three level of noise. In this figure mean square error (MSE) is evaluated over time which is equal to average of normalized errors.

Fig. 8 indicates the normalized error after $20$ time slots of the algorithm versus different SNRs.

In Fig. 9(a) we investigate the impact of the number of elements $M$ in each antenna array on the normalized error of localization. As this figure suggests, $M=8$ is sufficient for accurate localization within the area. Fig. 9(b) shows performance of localization in presence of different number of active PUs in a network with $100$ grid points and $10$ sensors each one is equipped with $8$ elements.

Fig. 10 considers the situation where there are 8 PUs in the area and the location of one PU is changing every $20$ time slots. Each location change is discovered and tracked. However, the location change occurred at time slot $40$ takes a longer time to be detected and tracked. The performance of our distributed localization algorithm for detecting and tracking $8$ PUs is compared with the case where sensors have omni-directional antennas. Our algorithm tracks the location changes and refines the estimated power map over time.

The connectivity of nodes plays a key role in the performance of our distributed algorithm. Fig. 11a and Fig. 11b compare two graphs of connectivity for $10$ spectrum sensors in which nodes are more connected to each other in Fig. 11b. The normalized error of localization for these two distributed networks are compared with the centralized solution in Fig. 11c. As the network is more connected, sensors reach consensus faster. However, the centralized solution converges to the final solution in much less number of time slots.

Fig. 12 shows performance of our proposed distributed algorithm for localization of PUs versus parameter $\alpha$ which regularizes the impact of consensus. The location of PUs are assumed to be constant over time and the distributed algorithm reaches consensus about location of PUs over time. This plot is drawn in terms of averaged normalized error of all sensors. At the first time slot different sensors may have very different estimation of PUs’ locations. However, after several time slots the accuracy of their estimations increases and their solutions converge to each other. Fig. 13 compares the estimated power map of two distant SUs. The ground truth power map is shown and the estimated power map using SU1 and SU2 are compared at time slots $1$ and $100$. At time slot $100$, all sensors have used information of all other sensors and a consensus has been reached on the PUs’ locations estimates.

VII Conclusion

The impact of directional antennas on the performance of source localization and tracking is studied. Omni-directional antennas provides redundant measurements over time. While, directional antennas are able to capture a set of incoherent projections from the underlying spectrum activity pattern. It improves the accuracy of source localization significantly. As the number of elements in an ULA increase, the quality of power spectrum map increases. However, employing $8$ or $16$ elements in each ULA is sufficient for practical scenarios. The proposed localization is used in a dynamic spectrum cartography system. Moreover, the distributed implementation of our proposed method is presented and evaluated.