# Reduced-Complexity Algorithms for Indoor Map-Aware Localization Systems

Authors: Francesco Montorsi; Fabrizio Pancaldi; Giorgio M. Vitetta
Journal: International Journal of Navigation and Observation
 (2015)

---

## Abstract

The knowledge of environmental maps (i.e., map-awareness)
can appreciably improve the accuracy of optimal
methods for position estimation in indoor scenarios. This improvement,
however, is achieved at the price of a significant complexity
increase with respect to the case of map-unawareness,
specially for large maps. This is mainly due to the fact that
optimal map-aware estimation algorithms require integrating
highly nonlinear functions or solving nonlinear and nonconvex
constrained optimization problems. In this paper, various
techniques for reducing the complexity of such estimators are
developed. In particular, two novel strategies for restricting the
search domain of map-aware position estimators are developed
and the exploitation of state-of-the-art numerical integration and
optimization methods is investigated; this leads to the development
of a new family of suboptimal map-aware localization
algorithms. Our numerical and experimental results evidence that
the accuracy of these algorithms is very close to that offered
by their optimal counterparts, despite their significantly lower
computational complexity.


---

## Body

## 1. Introduction

Indoor localization systems have found widespread application in a number of different areas, both military and civilian [1]. In many cases, the maps of the environments where the users are supposed to lie (e.g., the floor plans of a given building or a spatial road network) are perfectly known and this form of a priori knowledge (dubbedmap-awareness in the following) can be employed for improving localization accuracy. In particular, in the technical literature about localization systems map information is usually exploited to (a) improve the quality of position estimates generated by standard localization technologies; (b) set simple geometric constraints in the search space offingerprinting systems (e.g., see [2–4]); (c) acquire an accurate knowledge of propagation of wireless signals on the basis of two-dimensional or three-dimensionalray-tracing techniques [5, 6]; (d) developmap-based statistical models for the measurements acquired by multiple anchors in a given indoor environment. It is important to point out the following:(i)
Point (a) refers to the wide family ofmap-matching algorithms (see [7] and references therein), which includes both simple search techniques exploiting geometric and topological information provided by maps, and more refined techniques based on signal processing algorithms, such as Kalman filtering and particle filtering [8]. Distinct types of map-matching algorithms may offer a substantially different complexity-performance trade-off, but the improvement in localization accuracy is limited by the fact that map information is not exploited in the first stage of position estimation (i.e., in ranging). Note also that these algorithms have been mainly developed for transport applications (where they are employed to refine the estimated position provided by GPS or GPS integrated withdead-reckoning in outdoor environments). In principle, they can be also employed in indoor scenarios, provided that a proper map (e.g., a schematic representation of the possible user positions such as a Voronoi diagram) is available for the considered environment; for instance, map-matching methods based on particle filtering have been proposed in [9–13] for localization systems operating in indoor scenarios.(ii)
Fingerprinting methods (see point (b)) make a limited use of the geometrical properties of the propagation scenario, since they mainly rely on a training database. In addition, the preliminary measurement campaign necessary to generate such a database can represent a formidable task in large buildings and/or time-varying scenarios.(iii)
Ray-tracing techniques (see point (c)) can certainly outperform the above-mentioned techniques in terms of accuracy; however, real-time ray-tracing is often unfeasible because of its large computational complexity and the detailed knowledge it requires about various physical properties (e.g., electrical permittivity of walls and boundary conditions) of localization scenarios.(iv)
Point (d) refers to the so-calledmap-aware statistical localization techniques, which have been developed fortime of arrival (TOA),time difference of arrival (TDOA), andreceived signal strength- (RSS-) based indoor localization systems [14, Sec. 4.11.5], [15, 16]. Such techniques rely on the availability ofmap-aware statistical signal models and are characterized by the following relevant features: (a) they are able to compensate for thenon-line-of-sight (NLOS)bias (i.e., for the extra delay or extra attenuation due to propagation through obstructions, mainly walls), which represents a major source of error in indoor localization; (b) their use requires measurement campaigns substantially less time consuming than those needed for fingerprinting techniques in the same scenarios. In addition, previous work illustrated in [15, 16] has evidenced that these models can be easily combined with optimal estimation techniques to derive novel localization algorithms.This paper represents a follow-up to [15], where new map-aware statistical models based on experimental measurements have been exploited to develop optimalmap-aware minimum mean square error (MMSE) andmaximum a posteriori (MAP) estimators. Our previous work has evidenced that these estimators may substantially outperform their map-unaware counterparts at the price, however, of a significant complexity increase, specially for large maps. In this paper, the problem of complexity reduction of map-aware MMSE and MAP estimators is tackled and the following novel contributions are developed:(1)
Two novel algorithms for restricting the domain over which the agent position is searched for in map-aware estimation are devised. These algorithms, dubbeddistance-reduced domain (DRD) andprobability-reduced domain (PRD) in the following, can reduce the rate at which the complexity of map-aware estimators increases with map size without substantially affecting localization accuracy.(2)
The application of specific mathematical tools (namely,cubature rules for numerical integration [17, 18] anddirect-search methods [19, 20]), which can be exploited to reduce the implementation complexity of map-aware algorithms, is analysed.(3)
Novel suboptimal map-aware algorithms resulting from the combination of the above-mentioned algorithms and mathematical tools are proposed and are compared, in terms of bothroot median square error (RMSE) and computational complexity (namely, overall number offloating point operations (FLOPs)), tomaximum likelihood (ML) estimators (note that ML estimators are optimal formap-unawarelocalization and have been widely adopted for their limited computational complexity [21]). The following is worth pointing out:(i)
Most of the performance results illustrated in this paper mainly rely on a set of data generated according to the statistical model developed in [15], which, in turn, is based on experimental measurements. This approach, which is commonly adopted in the technical literature (e.g., see [22, 23]), is motivated by the fact that experimental databases referring to localization systems have usually a limited size; this is mainly due to the time consuming tasks required by any measurement campaign in this field. However, some performance results evaluated on the basis of our experimental database are also shown.(ii)
In our work only RSS-based localization algorithms are considered, even if the proposed approach can be easily extended to TOA and TDOA based systems, since all rely on similar statistical signal models [15].(iii)
Our contribution focuses on the role played bymap-awareness in the localization of still targets in static environments (consequently, the potential improvement deriving from the knowledge of mobility models is not taken into consideration).The remaining part of this paper is organized as follows. In Section 2, some models for indoor maps and measurements in map-aware RSS-based localization systems are illustrated. In Section 3, optimal and suboptimal localization algorithms based on the proposed models are developed, whereas in Section 4 various problems arising from their implementation are analysed and some solutions are proposed. In Section 5 various numerical results about the performance and complexity of the devised algorithms are illustrated. Finally, in Section 6, some conclusions are drawn.Notations. Theprobability density function (pdf) of a random vector R evaluated at the point r is denoted f(r); Nr;m,Σ denotes the pdf of a Gaussian random vector R having mean m and covariance matrix Σ, evaluated at the point r; diag{·} denotes a square matrix having the arguments on its main diagonal and zeros elsewhere; S denotes the cardinality of the set S; xi denotes the ith element of the vector x and x denotes its size; x denotes the floor of x; projS,p denotes the Euclidean projection of p over a domain S, that is, the point q∈S which minimises the Euclidean distance with p; centreS denotes the centre of the rectangle S; ex (ey) is the unit vector along the x (y) axis.


## 2. Reference Scenario and Signal Model

In the following, we focus on atwo-dimensional (2D) RSS-based localization system employing Na devices, calledanchors, and whose positions piai=1Na are known, to estimate the position p≜x,yT of anagent (see Figure 1). The mathematical models defined in [15] for the map, the wireless connectivity between the agent and the anchors, and the measurements acquired for position estimation are adopted; their essential features are summarised below.Figure 1
Localization system employingNa=6 anchors (green squares) to localise a single agent (blue circle) and operating in a bounded map (whose support R is identified by the grey area). The dashed line denotes the map bounding box BR. LOS (NLOS) links are evidenced by continuous (dashed) lines (note that Nobs=4 in the considered case, since the agent is not connected to anchors 4 and 5). The dotted line represents the outer edge of a possible rectangular covering, whereas the Ns=6 segments representing the environmental obstructions are evidenced in red.


### 2.1. Map Model

Arectangular uniform map, having support R and area AR, is assumed in the following, since it provides a good approximation of the floor plans of typical buildings (see [15, eq. (1)]). This map consists of the union of Nr nonoverlapping rectangles {Rk⊂R2}k=1Nr having their sides parallel to the axes of the adopted reference frame and representing the so-calledrectangular covering of R, so that R=⋃k=1Nr‍Rk (see Figure 1). Note that, generally speaking, the support R of a rectangular map can be approximated by different coverings {Rk} and there is no one-to-one mapping between the support R of a (rectangular) map and its covering Rk. However, various algorithms are already available in the technical literature to partition a generic polygon into a set of rectangles (e.g., see [24]) and, in particular, to generate aminimal nonoverlapping covering (MNC) [25, 26] of R that is a partition consisting of aminimal number of rectangles for the considered domain. In the following, a MNC of R, which is a partition {RkMNC} consisting of aminimal number of rectangles for this domain, or adense nonoverlapping covering (DNC) {RkDNC}, which originates from splitting the MNC rectangles whose areas or side lengths exceed a given threshold, is considered, so that(1)R=⋃k=1NrMNCRkMNC=⋃k=1NrDNCRkDNC;here NrMNC (NrDNC) is the number of rectangles in the selected MNC (DNC) of R. As it will become clearer later, a MNC of R can be very useful to reduce the computational load of some map-aware localization algorithms; some algorithms, however, perform better if a DNC is selected. In localization systems, R is usually time invariant and such coverings can be computed offline resorting to a number of optimized (but complicated) algorithms (e.g., see [25, 26] and references therein); in our work, however, the simpler algorithm described in [27] is adopted, since it can handle maps containing “holes” (e.g., inaccessible areas of a building floor) and its input data can be easily extracted from the floor plans of typical buildings.Note that our approach to map-aware localization requires not only the knowledge of a MNC (or, equivalently, a DNC) of the map supportR, but also that of the obstructions (e.g., walls) the map contains; this knowledge is exploited to evaluate the function No(p1,p2), which represents the number of obstructions interposed between two arbitrary points p1 and p2 of R. On the contrary,map-unaware localization systems, which are introduced later for comparison, rely on the knowledge of thebounding boxBR of R only.


### 2.2. Connectivity Model

In our work, the coverage regionR(i) of the ith anchor (with i=1,2,…,Na) is assumed to consist of a circle centred at pia whose radius dmax,i depends on the transmitted power and signal propagation conditions; inside this region, if the agent is connected with the associated anchor, it acquires a singleobservationzi (in particular, a RSS measurement) for localization purposes by exchanging wireless signals with the anchor itself [15, Sec. II.B.] (if multiple measurements are acquired by the ith anchor, they are averaged or filtered in order to generate a single observation zi). In practice, in the presence of harsh propagation conditions, the shape of the coverage region is likely to substantially differ from a circular shape, so that the number Nobs of observations available to the agent cannot be predicted theoretically; for this reason, the parameter Nobs≜Z is defined, where Z is the set of indices associated with the anchorstruly connected with the agent. It is easy to show that if i∈Z, then p∈Ri⋂R; consequently, in amap-aware localization system the agent position has to be searched for inside the domain R(Z)≜⋂i∈Z‍Ri⋂R. In our work, the region R(Z) is assumed to be well approximated by the rectangular region that consists of all the covering rectangles intersecting with the exact R(Z); that is,(2)RZ≃⋃k∈KRk,where(3)K≜k=1,…,Nr∣Rk∩Ri≠∅,∀i∈Zand Rk belongs to a MNC or a DNC set (in the following, unless explicitly stated, the symbol Rk (Nr) denotes RkDNC (NrDNC) or RkMNC (NrMNC)). Note that applying this connectivity model to the considered localization problem results in a preliminarydomain reduction; in fact, instead of considering the NrMNC or NrDNC rectangles covering the whole map support, all the proposed localization algorithms restrict their search domain to NrSRC≜K rectangles approximating R(Z) (2). In amap-unaware system, instead, the agent position is expected to belong to(4)PZ≜⋂i∈ZRi∩BR.In this case, it is assumed that for any i the domain R(i) can be approximated by a square region whose centre is pia and whose side is 2dmax,i; this implies that P(Z) (4) consists of a single rectangle.Finally, we note that the setsZ and K in (3) and the domains R(Z) (2) and P(Z) (4) are all available to the agent localization algorithm immediately after all RSS observations have been acquired.


### 2.3. Statistical Modelling of Observations

We consider a wireless localization system inferring the agent positionp from a set of Nobs RSS observations zi,i∈Z, with Nobs≤Na. In [15] it is shown that, given the trial agent position p~∈R(Z), themap-aware likelihood for the observation vector z≜zi,∀i∈ZT∈RNobs can be expressed as (see [15, Eq. (13)])(5)fz∣p~=∏i∈ZNzi;dip~+μb,ip~,σb,i2p~+σn,i2p~.Here di(p~)≜p~-pia is the Euclidean distance between the trial position p~ and the ith anchor and(6)μb,ip~≜μb,0+μb,mNop~,piauNo,σb,ip~≜max⁡σb,0-σb,mNop~,pia,0denote the mean and the standard deviation, respectively, of thebiasbi(p)~N(μb,i(p),σb,i2(p)) affecting the ith observation zi and originating from signal propagation inside the obstructions, and(7)σn,ip~≜σn,0dip~d0βn represents the standard deviation of the bias-unrelated and position-dependentnoiseni(p)~N(0,σn,i2(p)) affecting the same observation. Note that the likelihood function f(z∣p~) (5) is influenced by the structure of the map support R through the parameters μb,i(p~) and σb,i(p~), since both these quantities depend on Nop~,pia.Map-unaware systems, instead, cannot rely on the knowledge of the functionNo·,· in bias modelling. Moreover, state-of-the-art map-unaware systems often employ a LOS/NLOS detection stage before localization, whose output Z^NLOS (Z^LOS) is the subset of Z containing the indices of the anchors detected to experience NLOS (LOS) conditions with the agent. Given Z^NLOS and the trial agent position p~∈P(Z), themap-unaware likelihood(see [15, Eq. (15)])(8)fz;p~=∏i∈ZNzi;dip~+μb,iMUp~,σb,iMUp~2+σn,iMUp~2is adopted for the same observation vector z considered in (5). Here μb,iMU≜κb and σb,iMU≜γb for i∈Z^NLOS (both are zero otherwise) and σn,iMU(p~)=σn,0MUdi(p~)/d0βnMU.


### 2.4. Calibration

The parameters appearing in both the map-aware and map-unaware likelihood functions (see (5) and (8), resp.) need to be extracted from a set of experimental measurements, as illustrated in [15]. In our work, the performance of the developed localization algorithms has been assessed on the floors of 3 different buildings of the University of Modena and Reggio Emilia (these buildings mainly host offices and laboratories). In each scenario, 15–20 measurement sites have been identified in order to acquire RSS data referring to almost 100 independent links. Then, the method illustrated in [15] has been exploited to extract the values of the unknown parameters appearing in (5) and (8); this has produced the following values: μb,0=12.6 m, μb,m=2.53 m, σb,0=7.07 m, σb,m=3.0 m, σn,0=2.47 m, βn=0.21, and d0=1 m for the map-aware model (5) and the parameters κb=21 m, γb=2.81 m, σn,0MU=4.47 m, and βnMU=0.19 for the map-unaware model (8) (please see [15, Table I] for further details). It is important to note the following:(i)
In our experimental campaign, we have also found that the same values of the considered parameters can accurately describe the statistical models to be adopted for RSS measurements in different scenarios; this result is motivated by the fact that actually similar propagation conditions are experienced in distinct buildings (which can be classified aslight commercial in standard terminology), even if the areas of the considered propagation scenarios (floor plans) are significantly different (see Section 5).(ii)
In principle, a preliminary measurement campaign is always required to adapt models (5) and (8) to different propagation scenarios. This need is common to various localization techniques, which may require a time consuming training phase.(iii)
Map-aware solutions require a proper preprocessing of the map to compute its MNC or DNC and to evaluate the functionNo(p1,p2) (see Section 2.1); however, this step can be carried out offline and does not appreciably affect the overall computational complexity (see also Section 4.1).


## 2.1. Map Model

Arectangular uniform map, having support R and area AR, is assumed in the following, since it provides a good approximation of the floor plans of typical buildings (see [15, eq. (1)]). This map consists of the union of Nr nonoverlapping rectangles {Rk⊂R2}k=1Nr having their sides parallel to the axes of the adopted reference frame and representing the so-calledrectangular covering of R, so that R=⋃k=1Nr‍Rk (see Figure 1). Note that, generally speaking, the support R of a rectangular map can be approximated by different coverings {Rk} and there is no one-to-one mapping between the support R of a (rectangular) map and its covering Rk. However, various algorithms are already available in the technical literature to partition a generic polygon into a set of rectangles (e.g., see [24]) and, in particular, to generate aminimal nonoverlapping covering (MNC) [25, 26] of R that is a partition consisting of aminimal number of rectangles for the considered domain. In the following, a MNC of R, which is a partition {RkMNC} consisting of aminimal number of rectangles for this domain, or adense nonoverlapping covering (DNC) {RkDNC}, which originates from splitting the MNC rectangles whose areas or side lengths exceed a given threshold, is considered, so that(1)R=⋃k=1NrMNCRkMNC=⋃k=1NrDNCRkDNC;here NrMNC (NrDNC) is the number of rectangles in the selected MNC (DNC) of R. As it will become clearer later, a MNC of R can be very useful to reduce the computational load of some map-aware localization algorithms; some algorithms, however, perform better if a DNC is selected. In localization systems, R is usually time invariant and such coverings can be computed offline resorting to a number of optimized (but complicated) algorithms (e.g., see [25, 26] and references therein); in our work, however, the simpler algorithm described in [27] is adopted, since it can handle maps containing “holes” (e.g., inaccessible areas of a building floor) and its input data can be easily extracted from the floor plans of typical buildings.Note that our approach to map-aware localization requires not only the knowledge of a MNC (or, equivalently, a DNC) of the map supportR, but also that of the obstructions (e.g., walls) the map contains; this knowledge is exploited to evaluate the function No(p1,p2), which represents the number of obstructions interposed between two arbitrary points p1 and p2 of R. On the contrary,map-unaware localization systems, which are introduced later for comparison, rely on the knowledge of thebounding boxBR of R only.


## 2.2. Connectivity Model

In our work, the coverage regionR(i) of the ith anchor (with i=1,2,…,Na) is assumed to consist of a circle centred at pia whose radius dmax,i depends on the transmitted power and signal propagation conditions; inside this region, if the agent is connected with the associated anchor, it acquires a singleobservationzi (in particular, a RSS measurement) for localization purposes by exchanging wireless signals with the anchor itself [15, Sec. II.B.] (if multiple measurements are acquired by the ith anchor, they are averaged or filtered in order to generate a single observation zi). In practice, in the presence of harsh propagation conditions, the shape of the coverage region is likely to substantially differ from a circular shape, so that the number Nobs of observations available to the agent cannot be predicted theoretically; for this reason, the parameter Nobs≜Z is defined, where Z is the set of indices associated with the anchorstruly connected with the agent. It is easy to show that if i∈Z, then p∈Ri⋂R; consequently, in amap-aware localization system the agent position has to be searched for inside the domain R(Z)≜⋂i∈Z‍Ri⋂R. In our work, the region R(Z) is assumed to be well approximated by the rectangular region that consists of all the covering rectangles intersecting with the exact R(Z); that is,(2)RZ≃⋃k∈KRk,where(3)K≜k=1,…,Nr∣Rk∩Ri≠∅,∀i∈Zand Rk belongs to a MNC or a DNC set (in the following, unless explicitly stated, the symbol Rk (Nr) denotes RkDNC (NrDNC) or RkMNC (NrMNC)). Note that applying this connectivity model to the considered localization problem results in a preliminarydomain reduction; in fact, instead of considering the NrMNC or NrDNC rectangles covering the whole map support, all the proposed localization algorithms restrict their search domain to NrSRC≜K rectangles approximating R(Z) (2). In amap-unaware system, instead, the agent position is expected to belong to(4)PZ≜⋂i∈ZRi∩BR.In this case, it is assumed that for any i the domain R(i) can be approximated by a square region whose centre is pia and whose side is 2dmax,i; this implies that P(Z) (4) consists of a single rectangle.Finally, we note that the setsZ and K in (3) and the domains R(Z) (2) and P(Z) (4) are all available to the agent localization algorithm immediately after all RSS observations have been acquired.


## 2.3. Statistical Modelling of Observations

We consider a wireless localization system inferring the agent positionp from a set of Nobs RSS observations zi,i∈Z, with Nobs≤Na. In [15] it is shown that, given the trial agent position p~∈R(Z), themap-aware likelihood for the observation vector z≜zi,∀i∈ZT∈RNobs can be expressed as (see [15, Eq. (13)])(5)fz∣p~=∏i∈ZNzi;dip~+μb,ip~,σb,i2p~+σn,i2p~.Here di(p~)≜p~-pia is the Euclidean distance between the trial position p~ and the ith anchor and(6)μb,ip~≜μb,0+μb,mNop~,piauNo,σb,ip~≜max⁡σb,0-σb,mNop~,pia,0denote the mean and the standard deviation, respectively, of thebiasbi(p)~N(μb,i(p),σb,i2(p)) affecting the ith observation zi and originating from signal propagation inside the obstructions, and(7)σn,ip~≜σn,0dip~d0βn represents the standard deviation of the bias-unrelated and position-dependentnoiseni(p)~N(0,σn,i2(p)) affecting the same observation. Note that the likelihood function f(z∣p~) (5) is influenced by the structure of the map support R through the parameters μb,i(p~) and σb,i(p~), since both these quantities depend on Nop~,pia.Map-unaware systems, instead, cannot rely on the knowledge of the functionNo·,· in bias modelling. Moreover, state-of-the-art map-unaware systems often employ a LOS/NLOS detection stage before localization, whose output Z^NLOS (Z^LOS) is the subset of Z containing the indices of the anchors detected to experience NLOS (LOS) conditions with the agent. Given Z^NLOS and the trial agent position p~∈P(Z), themap-unaware likelihood(see [15, Eq. (15)])(8)fz;p~=∏i∈ZNzi;dip~+μb,iMUp~,σb,iMUp~2+σn,iMUp~2is adopted for the same observation vector z considered in (5). Here μb,iMU≜κb and σb,iMU≜γb for i∈Z^NLOS (both are zero otherwise) and σn,iMU(p~)=σn,0MUdi(p~)/d0βnMU.


## 2.4. Calibration

The parameters appearing in both the map-aware and map-unaware likelihood functions (see (5) and (8), resp.) need to be extracted from a set of experimental measurements, as illustrated in [15]. In our work, the performance of the developed localization algorithms has been assessed on the floors of 3 different buildings of the University of Modena and Reggio Emilia (these buildings mainly host offices and laboratories). In each scenario, 15–20 measurement sites have been identified in order to acquire RSS data referring to almost 100 independent links. Then, the method illustrated in [15] has been exploited to extract the values of the unknown parameters appearing in (5) and (8); this has produced the following values: μb,0=12.6 m, μb,m=2.53 m, σb,0=7.07 m, σb,m=3.0 m, σn,0=2.47 m, βn=0.21, and d0=1 m for the map-aware model (5) and the parameters κb=21 m, γb=2.81 m, σn,0MU=4.47 m, and βnMU=0.19 for the map-unaware model (8) (please see [15, Table I] for further details). It is important to note the following:(i)
In our experimental campaign, we have also found that the same values of the considered parameters can accurately describe the statistical models to be adopted for RSS measurements in different scenarios; this result is motivated by the fact that actually similar propagation conditions are experienced in distinct buildings (which can be classified aslight commercial in standard terminology), even if the areas of the considered propagation scenarios (floor plans) are significantly different (see Section 5).(ii)
In principle, a preliminary measurement campaign is always required to adapt models (5) and (8) to different propagation scenarios. This need is common to various localization techniques, which may require a time consuming training phase.(iii)
Map-aware solutions require a proper preprocessing of the map to compute its MNC or DNC and to evaluate the functionNo(p1,p2) (see Section 2.1); however, this step can be carried out offline and does not appreciably affect the overall computational complexity (see also Section 4.1).


## 3. Estimation Algorithms

In this Section, optimal map-aware (namely, MMSE and MAP) and map-unaware (namely, ML) localization algorithms are analysed first. Then, two novel heuristic (suboptimal) strategies are proposed to reduce the computational load required by optimal estimation algorithms.


### 3.1. Optimal Map-Unaware Estimator

In a map-unaware context (i.e., in the absence of prior information about the position of the agent), the ML approach can be adopted to generate the optimal estimatep^ML(z)≜argmaxp~lnf(z;p~) of p [28], where f(z;p~) is given by (8). In practice, this estimate can be evaluated as(9)p^MLz=arg⁡minp~∈PZ⁡∑i∈Zzi-dip~-μb,iMU2σn,iMUp~2+σb,iMU2+2ln⁡∏i∈Zσn,iMUp~2+σb,iMU2,where P(Z) is a single rectangle (see Section 2.2). The last expression deserves the following comments:(1)
The cost function appearing in (9) exhibits a discontinuous behaviour, even if it does not contain map-aware functions (e.g., No(·,·)), since μb,iMU and σb,iMU are assigned a constant value or zero depending on Z^NLOS. Consequently, in principle, optimization methods involving gradients or Hessians cannot be used, unless the cost function is, at leastapproximately, continuous.(2)
The evaluation ofp^ML(z) in (9) requires solving a constrained least squares problem with a quadratic regularization term (represented by the product appearing in the right-hand side of (9)) which penalizes trial positions p~ affected by large noise and bias. Unluckily, standard optimization methods do not necessarily achieve theglobal minimum; when this occurs, large localization errors are found (note that some methods developed to mitigate the problem of local minima, likeprojection onto convex sets (POCS) andmultidimensional scaling (MDS), cannot be straightforwardly applied in this case because of the significant complexity of our observation model).(3)
The complexity of the cost function in (9) can be related to Nobs as ONobs.


### 3.2. Optimal Map-Aware Estimators

In a map-aware context, a MMSE or a MAP approach can be employed in position estimation. In particular, it is not difficult to show that the MMSE estimatep^MMSE(z) of p is given by [28](10)p^MMSEz≜Epp∣z=∬Rxfz∣pdp∬Ryfz∣pdp∬Rfz∣pdp.Substituting (5) in (10) yields, after some manipulation,(11)p^MMSEz=∬RZxhz,p/∏i∈Zσn,i2p+σb,i2pdp∬RZyhz,p/∏i∈Zσn,i2p+σb,i2pdp∬RZhz,p/∏i∈Zσn,i2p+σb,i2pdp,where(12)hz,p≜exp⁡-12∑i∈Zzi-dip-μb,ipσn,ip2.The last result deserves the following comments:(1)
The MMSE estimation procedure, unlike its MAP and ML counterparts, does not involve an optimization step and, consequently, does not suffer from problem of local minima. However, it requires numerical multidimensional integration.(2)
Our numerical results have evidenced that (a) the product∏i∈Z‍σn,i2(p)+σb,i2(p) appearing in (11) exhibits small variations inside the considered map if the anchor density is approximately constant and (b) the accuracy of MMSE estimation is negligibly influenced by the presence of this term for any p∈R such that di(p)≫1∀i∈Z (i.e., for any agent position far from the anchors). Then, discarding this term yields the approximate (and computationally simpler) variant(13)p^MMSEz≃∬RZxhz,pdp∬RZyhz,pdp∬RZhz,pdpof the MMSE estimator.(3)
The computational load required by the evaluation ofNo(·,·) is O(Ns) (see Section 2.1), where Ns is thenumber of segment intersection tests accomplished in processing a given observation; consequently, the complexity required by h(z,p) (12) is approximately ON¯s·Nobs, where N¯s is the average of Ns evaluated over the observations collected in z (note that N¯s is nonlinearly related to the shape of obstructions, the trial agent position, and the anchor positions).An alternative to the MMSE estimator is offered by the so-called MAP estimator, which can be expressed as p^MAP(z)≜argmaxp~lnf(z∣p~) [28] or, if (5) and Bayes’ rule are exploited, as(14)p^MAPz=arg⁡minp~∈RZ⁡fcp~,where(15)fcp~≜-2ln⁡fp~∣z=∑i∈Zzi-dip~-μb,ip~2σn,i2p~+σb,i2p~+2ln⁡∏i∈Zσn,i2p~+σb,i2p~is the MAPcost function. This function deserves the following comments:(1)
It is not differentiable, since it depends onNo(·,·). Consequently, its gradient/Hessian matrix cannot be evaluated analytically and, in principle, steepest descent methods cannot be used; in practice, such methods can be exploited if it is assumed that the cost function isapproximately continuous.(2)
For largeNobs, it is characterized by a better numerical stability than that of the function h(z,p) (12). In fact, for large Nobs, the argument of the function exp(·) appearing in h(z,p) becomes small, so that the resulting numerical values of h(z,p) itself may quickly drop below machine precision. On the contrary, this problem is not experienced with fc(p~) (15), because of the presence of a natural logarithm in its expression.(3)
Its computational complexity is similar to that ofh(z,p), which is ON¯s·Nobs.


### 3.3. Distance-Reduced Domain Map-Aware Estimator

The computational complexity of the MMSE and MAP estimators described above becomes unacceptable when the considered map is large. To mitigate this problem, a new technique for restricting the domain over which the agent position is searched for has been developed; it consists of the following three steps:(1)
Map-Unaware Estimation. A raw estimate p^raw(z) is generated by the map-unaware estimator (9), without performing any NLOS correction (this is equivalent to assuming Z^LOS≡Z and Z^NLOS=∅ in (9)).(2)
Domain Reduction. A portion(16)RredZ≜⋃k∈K^RkofR close to p^raw is extracted, where K^ is a subset of the values of the index k associated with the set of rectangles Rk forming R(Z) (2); in practice, K^ is generated according to the heuristic criterion(17)K^≜k∈K:p^raw-projRk,p^raw<dth,wheredth≜mediani∈Zγσb,i2p^raw+σn,i2p^raw is athreshold distance, γ is a fixed parameter, and K is defined by (3). Note that the definition given above for dth entails that only the covering rectangles within γ standard deviations from p^raw are selected.(3)
Map-Aware Estimation. The final estimate p^DRD(z) is generated using either the MMSE (13) or the MAP estimator (14) under the constraint that p∈Rred(Z); intuitively, since the subset Rred(Z) is smaller than R(Z) (2), the computational load required by this step is lower than that required by (13) or (14).Note that the proposed technique mitigates localization complexity by restricting the search domain R(Z) on the basis of a distance criterion, in the sense that only a subset of rectangles, close to the raw estimate p^raw(z), is taken into consideration. For this reason, in the following, this procedure is dubbeddistance-reduced domain MMSE (DRD-MMSE) ordistance-reduced domain MAP (DRD-MAP), if a MMSE or a MAP estimator is employed in the last step, respectively. Note that DRD performance strongly depends on (a) the accuracy of the p^raw(z) estimate and (b) the criterion adopted in generating the set K^ (17) (if the proposed criterion fails to select the map portion where the agent truly lies, then a large error is made in step (3), where it is assumed that p∈Rred(Z)). In principle, other heuristic criteria (e.g., based on selecting a fixed number of rectangles close to p^raw) could be adopted; however, (17) turned out to provide the best accuracy/complexity trade-off among the “distance-based” criteria we tested in our computer simulations.


### 3.4. Probability-Reduced Domain Map-Aware Estimator

An alternative to the DRD approach is based on a sort ofprobabilistic criterion and consists of the following two steps:(1)
Map-Aware Raw Estimation. Similarly to step (1) of the DRD algorithm, a portion Rred(Z) (16) of R is extracted. However, in this case, the integer set K^ is generated by taking the values of the index k associated with the Nr,min largest elements of the set(18)F≜ln⁡fz∣p~=pkc,∀k∈K∩N,whereNr,min is a fixed parameter, pkc≜centreRk, f(z∣p~) is given by (5), K is defined by (3), and N is a proper integer set, whose target is reducing the cardinality of F (and thus the overall computational load of the estimation algorithm). In practice, N≜k∣pkc-pic>dth,∀i∈N∖k has been selected, where dth is a fixed real and positive parameter, so that only rectangles spaced by at least dth meters are considered (instead of all rectangles forming the search domain R(Z) (2)). The effect of the choice expressed by (18) is to generate the reduced domain Rred(Z) (16) by selecting the Nr,min rectangles that, on the basis of their centers, exhibit the largest likelihood values.(2)
Map-Aware Estimation. The final estimate p^PRD(z) is evaluated exploiting either the MMSE (13) or the MAP estimator (14) under the constraint that p∈Rred(Z).The resulting estimation technique is dubbedprobability-reduced domain MMSE (PRD-MMSE) orprobability-reduced domain MAP (PRD-MAP), if a MMSE or a MAP estimator, respectively, is employed in the last step of the procedure illustrated above. It is important to point out that this approach offers some advantages with respect to its DRD counterpart, since (a) the criterion (18) is more closely related to the optimal MAP criterion than (17) and (b) all the PRD steps rely on map-awareness, that is, on the most accurate modelling we propose. However, it should be also taken into account the fact that the first step of the PRD approach may entail a significant complexity, since the evaluation of a map-aware likelihood is required.


## 3.1. Optimal Map-Unaware Estimator

In a map-unaware context (i.e., in the absence of prior information about the position of the agent), the ML approach can be adopted to generate the optimal estimatep^ML(z)≜argmaxp~lnf(z;p~) of p [28], where f(z;p~) is given by (8). In practice, this estimate can be evaluated as(9)p^MLz=arg⁡minp~∈PZ⁡∑i∈Zzi-dip~-μb,iMU2σn,iMUp~2+σb,iMU2+2ln⁡∏i∈Zσn,iMUp~2+σb,iMU2,where P(Z) is a single rectangle (see Section 2.2). The last expression deserves the following comments:(1)
The cost function appearing in (9) exhibits a discontinuous behaviour, even if it does not contain map-aware functions (e.g., No(·,·)), since μb,iMU and σb,iMU are assigned a constant value or zero depending on Z^NLOS. Consequently, in principle, optimization methods involving gradients or Hessians cannot be used, unless the cost function is, at leastapproximately, continuous.(2)
The evaluation ofp^ML(z) in (9) requires solving a constrained least squares problem with a quadratic regularization term (represented by the product appearing in the right-hand side of (9)) which penalizes trial positions p~ affected by large noise and bias. Unluckily, standard optimization methods do not necessarily achieve theglobal minimum; when this occurs, large localization errors are found (note that some methods developed to mitigate the problem of local minima, likeprojection onto convex sets (POCS) andmultidimensional scaling (MDS), cannot be straightforwardly applied in this case because of the significant complexity of our observation model).(3)
The complexity of the cost function in (9) can be related to Nobs as ONobs.


## 3.2. Optimal Map-Aware Estimators

In a map-aware context, a MMSE or a MAP approach can be employed in position estimation. In particular, it is not difficult to show that the MMSE estimatep^MMSE(z) of p is given by [28](10)p^MMSEz≜Epp∣z=∬Rxfz∣pdp∬Ryfz∣pdp∬Rfz∣pdp.Substituting (5) in (10) yields, after some manipulation,(11)p^MMSEz=∬RZxhz,p/∏i∈Zσn,i2p+σb,i2pdp∬RZyhz,p/∏i∈Zσn,i2p+σb,i2pdp∬RZhz,p/∏i∈Zσn,i2p+σb,i2pdp,where(12)hz,p≜exp⁡-12∑i∈Zzi-dip-μb,ipσn,ip2.The last result deserves the following comments:(1)
The MMSE estimation procedure, unlike its MAP and ML counterparts, does not involve an optimization step and, consequently, does not suffer from problem of local minima. However, it requires numerical multidimensional integration.(2)
Our numerical results have evidenced that (a) the product∏i∈Z‍σn,i2(p)+σb,i2(p) appearing in (11) exhibits small variations inside the considered map if the anchor density is approximately constant and (b) the accuracy of MMSE estimation is negligibly influenced by the presence of this term for any p∈R such that di(p)≫1∀i∈Z (i.e., for any agent position far from the anchors). Then, discarding this term yields the approximate (and computationally simpler) variant(13)p^MMSEz≃∬RZxhz,pdp∬RZyhz,pdp∬RZhz,pdpof the MMSE estimator.(3)
The computational load required by the evaluation ofNo(·,·) is O(Ns) (see Section 2.1), where Ns is thenumber of segment intersection tests accomplished in processing a given observation; consequently, the complexity required by h(z,p) (12) is approximately ON¯s·Nobs, where N¯s is the average of Ns evaluated over the observations collected in z (note that N¯s is nonlinearly related to the shape of obstructions, the trial agent position, and the anchor positions).An alternative to the MMSE estimator is offered by the so-called MAP estimator, which can be expressed as p^MAP(z)≜argmaxp~lnf(z∣p~) [28] or, if (5) and Bayes’ rule are exploited, as(14)p^MAPz=arg⁡minp~∈RZ⁡fcp~,where(15)fcp~≜-2ln⁡fp~∣z=∑i∈Zzi-dip~-μb,ip~2σn,i2p~+σb,i2p~+2ln⁡∏i∈Zσn,i2p~+σb,i2p~is the MAPcost function. This function deserves the following comments:(1)
It is not differentiable, since it depends onNo(·,·). Consequently, its gradient/Hessian matrix cannot be evaluated analytically and, in principle, steepest descent methods cannot be used; in practice, such methods can be exploited if it is assumed that the cost function isapproximately continuous.(2)
For largeNobs, it is characterized by a better numerical stability than that of the function h(z,p) (12). In fact, for large Nobs, the argument of the function exp(·) appearing in h(z,p) becomes small, so that the resulting numerical values of h(z,p) itself may quickly drop below machine precision. On the contrary, this problem is not experienced with fc(p~) (15), because of the presence of a natural logarithm in its expression.(3)
Its computational complexity is similar to that ofh(z,p), which is ON¯s·Nobs.


## 3.3. Distance-Reduced Domain Map-Aware Estimator

The computational complexity of the MMSE and MAP estimators described above becomes unacceptable when the considered map is large. To mitigate this problem, a new technique for restricting the domain over which the agent position is searched for has been developed; it consists of the following three steps:(1)
Map-Unaware Estimation. A raw estimate p^raw(z) is generated by the map-unaware estimator (9), without performing any NLOS correction (this is equivalent to assuming Z^LOS≡Z and Z^NLOS=∅ in (9)).(2)
Domain Reduction. A portion(16)RredZ≜⋃k∈K^RkofR close to p^raw is extracted, where K^ is a subset of the values of the index k associated with the set of rectangles Rk forming R(Z) (2); in practice, K^ is generated according to the heuristic criterion(17)K^≜k∈K:p^raw-projRk,p^raw<dth,wheredth≜mediani∈Zγσb,i2p^raw+σn,i2p^raw is athreshold distance, γ is a fixed parameter, and K is defined by (3). Note that the definition given above for dth entails that only the covering rectangles within γ standard deviations from p^raw are selected.(3)
Map-Aware Estimation. The final estimate p^DRD(z) is generated using either the MMSE (13) or the MAP estimator (14) under the constraint that p∈Rred(Z); intuitively, since the subset Rred(Z) is smaller than R(Z) (2), the computational load required by this step is lower than that required by (13) or (14).Note that the proposed technique mitigates localization complexity by restricting the search domain R(Z) on the basis of a distance criterion, in the sense that only a subset of rectangles, close to the raw estimate p^raw(z), is taken into consideration. For this reason, in the following, this procedure is dubbeddistance-reduced domain MMSE (DRD-MMSE) ordistance-reduced domain MAP (DRD-MAP), if a MMSE or a MAP estimator is employed in the last step, respectively. Note that DRD performance strongly depends on (a) the accuracy of the p^raw(z) estimate and (b) the criterion adopted in generating the set K^ (17) (if the proposed criterion fails to select the map portion where the agent truly lies, then a large error is made in step (3), where it is assumed that p∈Rred(Z)). In principle, other heuristic criteria (e.g., based on selecting a fixed number of rectangles close to p^raw) could be adopted; however, (17) turned out to provide the best accuracy/complexity trade-off among the “distance-based” criteria we tested in our computer simulations.


## 3.4. Probability-Reduced Domain Map-Aware Estimator

An alternative to the DRD approach is based on a sort ofprobabilistic criterion and consists of the following two steps:(1)
Map-Aware Raw Estimation. Similarly to step (1) of the DRD algorithm, a portion Rred(Z) (16) of R is extracted. However, in this case, the integer set K^ is generated by taking the values of the index k associated with the Nr,min largest elements of the set(18)F≜ln⁡fz∣p~=pkc,∀k∈K∩N,whereNr,min is a fixed parameter, pkc≜centreRk, f(z∣p~) is given by (5), K is defined by (3), and N is a proper integer set, whose target is reducing the cardinality of F (and thus the overall computational load of the estimation algorithm). In practice, N≜k∣pkc-pic>dth,∀i∈N∖k has been selected, where dth is a fixed real and positive parameter, so that only rectangles spaced by at least dth meters are considered (instead of all rectangles forming the search domain R(Z) (2)). The effect of the choice expressed by (18) is to generate the reduced domain Rred(Z) (16) by selecting the Nr,min rectangles that, on the basis of their centers, exhibit the largest likelihood values.(2)
Map-Aware Estimation. The final estimate p^PRD(z) is evaluated exploiting either the MMSE (13) or the MAP estimator (14) under the constraint that p∈Rred(Z).The resulting estimation technique is dubbedprobability-reduced domain MMSE (PRD-MMSE) orprobability-reduced domain MAP (PRD-MAP), if a MMSE or a MAP estimator, respectively, is employed in the last step of the procedure illustrated above. It is important to point out that this approach offers some advantages with respect to its DRD counterpart, since (a) the criterion (18) is more closely related to the optimal MAP criterion than (17) and (b) all the PRD steps rely on map-awareness, that is, on the most accurate modelling we propose. However, it should be also taken into account the fact that the first step of the PRD approach may entail a significant complexity, since the evaluation of a map-aware likelihood is required.


## 4. Implementation of Estimation Algorithms

All the map-aware and map-unaware estimation algorithms illustrated in Section3 require numerical integration and/or the application of optimization methods and involve nonlinear and nondifferentiable cost functions. For these reasons, their implementation raises various problems; some possible solutions are illustrated below.


### 4.1. Map-Aware Implementations

In this section, two different implementations for the evaluation of the MMSE estimate on the basis of (13) are proposed and analysed in detail. The first implementation, called MMSEc, is based on the use of the so-calledcubature formulas (e.g., see [17, 29–31]), which approximate any multidimensional definite integral as a finite sum of terms depending on its integrand. In particular, in the following, any integration over the search domain R(Z) is expressed as the sum of NrSRC distinct integrals, each referring to a distinct rectangular domain (see (1)), and the cubature formulas illustrated in [30] are exploited for the evaluation of each of the NrSRC integrals. Then, from (13) the approximate expression(19)p^MMSEcz≃∑k=1NrSRC∑j=1Nn,kwk,jαk,jxhz,αk,j∑k=1NrSRC∑j=1Nn,kwk,jαk,jyhz,αk,j∑k=1NrSRC∑j=1Nn,kwk,jhz,αk,jis easily inferred, where α(k,j)=[α(k,j)(x),α(k,j)(y)]T≜ck+Lkα(k,j)(F) is the jthnode referring to the kth rectangle Rk≜[lk;rk]×[bk;tk] in R(Z) (2), wk,j≜wk,j(F)detLk is the correspondingweight, ck≜1/2lk+rk;bk+tkT, and Lk≜1/2diagrk-lk,tk-bk, whereas α(k,j)(F)∈[-1;1]×[-1;1] and wk,j(F) denote the jthnode and the jthweight, respectively, of the cubature formula of order Nn,k. It is worth mentioning that, in principle, the larger the size of Rk is, the stronger the fluctuations of the integrand function to be expected over this domain are; consequently, high cubature orders (i.e., large values of Nn,k) should be selected for large rectangles, at the price, however, of an increase in the computational load. In our computer simulations, the heuristic formula Nn,k=NnmARk+dmin has been adopted for selecting the number of nodes; here, Nn(·) denotes the function defined in [30] to map the degree d of a polynomial into the number of integration nodes, mARk+dmin is the polynomial degree heuristically associated with the (nonpolynomial) integrand functions appearing in (13), ARk=detLk is the area of Rk, whereas m and dmin are real parameters whose values are listed in Table 1. It is also important to point out that (a) all the sums appearing in (19) require the evaluation of the same quantities hz,α(k,j) and (b) MMSEc is noniterative and its approximate complexity is ON¯s·Nobs·Neval, where Neval=Nn¯·NrSRC and Nn¯ denotes the average number of nodes selected for integrating over each of the NrSRC rectangles. Further details about our implementation of MMSEc can be found in [27, 32].Table 1
Covering and specific values selected for the implementations illustrated in Section4.
Algorithm
Cov.
Parameters
M
M
S
E
c
MNC
m
=
1.4
m
-
2, dmin=5
M
M
S
E
r
MNC
Δ
=
0.3
m
M
A
P
f
MNC
Algorithm = active-set, maxiter =103, tolfun = tolx = 10-2
M
A
P
g
MNC
Δ
1
=
2
m, Δfinal=0.3m, Nr,min=3
M
A
P
c
MNC
Δ
1
=
2
m, Δfinal=0.3m
D
R
D
-
E
c, DRD-Pg
MNC
γ
=
0.6, Nr,min=5
P
R
D
-
E
c, PRD-Pg
MNC
d
t
h
=
3.5
m, Nr,min=5
P
R
D
-
P
p
DNC
d
t
h
=
3.5
m, Nr,min=10
M
L
f
—
The same as forMAPf
M
L
m
—
Δ
=
10
m; others the same as for MAPf
M
L
c
—
Δ
1
=
2
m, Δfinal=0.3mThe second implementation, calledMMSEr, is based on the well-knowntrapezoidal integration rule; then, from (13), the approximate formula(20)p^MMSErz≃∑k=1NrSRC∑i=1Nn,kx∑j=1Nn,kywk,i,jαxk,i,jhz,αk,i,j∑k=1NrSRC∑i=1Nn,kx∑j=1Nn,kywk,i,jαyk,i,jhz,αk,i,j∑k=1NrSRC∑i=1Nn,kx∑j=1Nn,kywk,i,jhz,αk,i,jis easily inferred. Here, thenodes{αk,i,j≜lk+i-1Δ,bk+j-1ΔT,i=1,2,…,Nn,kx,j=1,2,…,Nn,ky} of the kth integration formula correspond to the vertices of a regular grid extending over Rk and characterized by a spacing Δ (see Table 1), whereas the corresponding weights wk,i,j are evaluated as(21)
w
k
,
i
,
j
≜
Δ
2
4
for
i
∈
1
,
N
n
,
k
x
,
j
∈
1
,
N
n
,
k
y
Δ
2
2
for
i
∈
1
,
N
n
,
k
x
,
j
∈
2
,
…
,
N
n
,
k
y
-
1
Δ
2
2
for
i
∈
2
,
…
,
N
n
,
k
x
-
1
,
j
∈
1
,
N
n
,
k
y
Δ
2
for
i
∈
2
,
…
,
N
n
,
k
x
-
1
,
j
∈
2
,
…
,
N
n
,
k
y
-
1with Nn,kx≜(rk-lk)/Δ+1 and Nn,ky≜(tk-bk)/Δ+1. It is important to point out the following: (a) similarly to (19), the three sums in (20) require the evaluation of the same quantities hz,α(k,i,j); (b) MMSEr is noniterative and its approximate complexity is ON¯s·Nobs·Neval, where Neval=Nnx¯·Nny¯·NrSRC and Nnx¯ (Nny¯) denotes the average number of integration points selected along the x (y) direction.In implementing the MAP estimator (14), the search over the space R(Z) has been turned into Nr searches over the distinct rectangles Rk,k=1,2,…,Nr covering R(Z). This choice is motivated by the fact that (a) the constraint p~∈R(Z), unlike p~∈Rk, cannot be directly formulated as a set of linear inequalities, as required by standard optimization techniques; (b) solving Nr distinct optimization problems, each involving a search domain much smaller than R(Z), strongly reduces the probability of reaching local minima; (c) when iterative optimization algorithms are adopted over each of the rectangles Rk, fast convergence can be usually achieved if their centres are selected as initial points. Note that in this case direct-search methods [19, 20] should be preferred to standard Newton-based or gradient-based methods, since the first class of methods usually admits simple implementations for constrained problems and does not require smooth cost functions. These considerations have led us to develop 3 different implementations, called MAPf, MAPg, and MAPc, of the MAP strategy (14). The first one (namely, MAPf) is summarised in Algorithm 1 and is based on the use of the MATLABfmincon routine to solve the optimization problem appearing in step (2). This implementation is characterized by various adjustable parameters; their values, which are listed in Table 1, have been selected to minimise the overall computational complexity.Algorithm 1:Implementation of the MAP strategy (14) based on the fmincon routine (MAPf).
Require: covering Rk,k∈K of R(Z) (2); fcp (15).
(1)fork∈Kdo
(2)pmin,k←arg⁡minp~∈Rkfcp~
(3)fmin⁡,k←fcpmin⁡,k
(4)end for
(5)k^←arg⁡mink{fmin,k}, p^MAPf←pmin,k^The second implementation (namely,MAPg) is summarised in Algorithm 2 and is based on discretisation of the search space and on the use of iterative direct-search methods (see [19, Sec. 3]). In practice, in its kth iteration, for each rectangle of the selected covering, the cost function fc(·) (15) is evaluated at the vertices of a rectangular grid characterized by spacing Δk (see lines (4)-(5)); then, on the basis of the resulting values, only the most promising rectangles are saved for the next iteration (see line (8)), in which a halved step size is adopted (see line (9)). Note that both the computational complexity and the accuracy of MAPg depend on the initial and the final values of the step size (denoted by Δ1 and Δfinal, resp.), so that a proper trade-off has to be achieved in the selection of these parameters (the values we adopted are listed in Table 1).Algorithm 2:Implementation of the MAP strategy (14) based on an iterative grid search (MAPg).
Require: covering Rk=[lk;rk]×[bk,tk],k∈K of R(Z) (2); fcp (15); Δ1, Δfinal and Nr,min (Table 1).
(1)Δ←Δ1, K~←K
(2)whileΔ>Δfinaldo
(3)fork∈K~do
(4)I←1,…,(rk-lk)/Δ, J←1,…,(tk-bk)/Δ
(5)(i^,j^)←arg⁡mini∈I,j∈Jfc(lk,bk)+Δ(i,j)
(6)pmin,k←(lk,bk)+Δ(i^,j^), fmin⁡,k←fcpmin⁡,k
(7)end for
(8)K~← indices of the Nr,min smallest elements of {fmin⁡,k}
(9)Δ←Δ/2
(10)end while
(11)k^←arg⁡mink∈K~{fmin⁡,k}, p^MAPg←pmin⁡,k^The third implementation (namely,MAPc) is summarised in Algorithm 3 and is based on a different direct-search method, known as thecompass algorithm [20, Sec. 8.1]. Similarly to MAPg, MAPc exploits a discretisation of the search space and an iterative reduction of the step size; however, the grid it employs is not regular. In fact, in each rectangle, the domain is explored moving from a “current trial position” p~ on the basis of given displacement vectors ±ex or ±ey (see line (6)) in order to identify the direction along which the cost function decreases (see line (8)). Even in this case, the overall computational complexity and the accuracy of this implementation depend on the initial and final values of the step size (see Table 1), which need to be carefully selected.Algorithm 3:Implementation of the MAP strategy (14) based on an iterative compass search (MAPc).
Require: covering Rk,k∈K of R(Z) (2); fcp (15); Δ1 and Δfinal (Table 1).
(1)fork∈Kdo
(2)Δ←Δ1
(3)pmin,k←centreRk, fmin,k←fcpmin,k
(4)whileΔ>Δfinaldo
(5)bImproved←False
(6)for d∈ex,ey,-ex,-eydo
(7)p~←pmin,k+Δ⋅d
(8)iffcp~<fmin,k and p~∈Rkthen
(9)fmin,k←fcp~, pmin,k←p~
(10)     bImproved← True, break
(11)end if
(12)end for
(13)Δ←Δ/2 only if bImproved = False
(14)end while
(15)end for
(16)k^←arg⁡minkfmin,k, p^MAPc←pmin,k^It is not difficult to show that the overall computational complexity of each of the MAP implementations proposed above isON¯s·Nobs·Neval, where Neval denotes the overall number of times the cost function fc(·) (15) is evaluated. Unluckily, Neval cannot be easily related to the other parameters of the proposed MAP implementations because of its nonlinear dependence on such parameters; consequently, a more accurate estimate of the computational burden required by the proposed implementations cannot be provided.Let us focus now on the problem of combining the implementation of the MMSE and MAP estimators illustrated above with the DRD and PRD techniques described in Section3.3. As far as the DRD technique is concerned, we note that all the proposed MAP/MMSE implementations can be employed in its last step. In the following, however, the implementations of DRD-MMSE and DRD-MAP techniques are always based on MMSEc and MAPg, respectively; for this reason, they are denoted by  DRD-Ec and DRD-Pg, respectively. It is important to point out the following:(1)
The algorithmsDRD-Ec and DRD-Pg are not intrinsically iterative. However, they can be easily modified to include a mechanism for iteratively updating the value of the parameter γ in order to improve their robustness against inaccuracies in the initial setup (see Algorithm 4).(2)
The complexity ofDRD-Ec and DRD-Pg is approximately O(NrSRC+Nobs(N¯s·Neval+Neval,raw)), where Neval denotes the overall number of evaluations of the MMSE integrand (for DRD-Ec) or that of the MAP cost function (for DRD-Pg), and Neval,raw represents the overall number of evaluations of the ML cost function.Algorithm 4:Implementation of the iterative DRD estimator (see Section 3.3) employed by DRD-Ec (DRD-Pg).
Require: covering Rk,k∈K of R(Z) (2); p^MMSEc (p^MAPg) estimator; constants γ and Nr,min (Table 1).
(1)p^raw←p^MLm, computed assuming Z^NLOS=Z
(2)dth←mediani∈Zγσb,i2p^raw(z)+σn,i2p^raw(z)
(3)dk←p^raw-projRk,p^raw for k∈K (see (17))
(4)whileK^<Nr,mindo
(5)K^←k∣dk<dth, dth←1.5dth
(6)end while
(7)p^DRD-Ec←p^MMSEc (p^DRD-Pg←p^MAPg) evaluated in the reduced map Rred(Z)=⋃k∈K^‍RkIn our work,MMSEc and MAPg have been also selected for the last step of the PRD technique; the resulting implementations are denoted by PRD-Ec and PRD-Pg, respectively, in the following. It is worth mentioning that(1)
bothPRD-Ec and PRD-Pg are not iterative;(2)
their complexity isO(N¯s·Nobs·(Neval,1+Neval,2)), where Neval,1 and Neval,2 denote the overall number of evaluations of the MMSE integrand (for PRD-Ec) and that of the MAP cost function (for PRD-Pg) in the first step and in the second step, respectively (see Section 3.4);(3)
the values of parametersdth and Nr,min need to be carefully selected.Our numerical results (see Section 5) have evidenced that whatever option for map-aware estimation is selected, the parameter N¯s (i.e., the average number of intersection tests accomplished in processing each observation) appearing in the complexity estimate of all the implementations developed in this paragraph plays a significant role in penalizing them with respect to their map-unaware counterparts. To reduce N¯s, the following two techniques have been devised:(1)
Segment Partitioning (SP). This consists of an offline phase, where (a) R is partitioned in multiple cells using a nonuniform rectangular grid (in a way that each anchor lies on a vertex of the grid) and (b) the obstruction segments falling inside each cell are identified. Then, when No(pia,p2) is evaluated, only the segments of the cells intersecting with the (pia,p2) segment are taken into account, thus reducing N¯s.(2)
Precomputing (PC). This is based on an offline evaluation of the function No(·,·) over Npt points (uniformly) distributed inside each rectangle of the available MNC or DNC; this requires storing Npt·N¯a·NrMNC or Npt·N¯a·NrDNC integers in a table, where N¯a denotes the average number of connected anchors. This technique is “memory hungry” for large maps; however, the precomputed table replaces that memorising information about the position of Ns obstruction segments and can be stored efficiently using a limited number of bits for the integer representation of the precomputed values of No(·,·). Moreover, this memory requirement is often more than counterbalanced by the reduction of computational complexity (since N¯s=0 in this case).In our computer simulations, the first technique has been adopted for all the map-aware estimators. The second technique, instead, has been used in the first step of PRD-MAP estimation only; the resulting implementation is dubbed PRD-Pp. This choice is motivated by the fact that this estimator can be easily adapted to the use of precomputed quantities, since it evaluates the map-aware likelihood (8) at a set of fixed points {pkc} (see Section 3.4).


### 4.2. Map-Unaware Implementations

Before illustrating our implementations of the ML estimator in (9), it is important to mention that if the constraint p~∈P(Z) was neglected, standard algorithms (e.g., theNelder-Mead simplex method [33]) could be employed to solve the optimization problem it involves. Unluckily, such a constraint plays an important role, specially in large maps, since (a) minima different from the global one are substantially less likely to be contained in the domain P(Z) than in the whole R2 space and (b) the map-unaware likelihood (8) is not defined for p~∉P(Z), as already explained in [15, Sec. II.C.]. For these reasons, in the 3 different implementations of the ML estimator (9) proposed below, the constraint p~∈P(Z) is always accounted for. The first implementation, called MLf, relies on the MATLABfmincon routine, which is employed to solve the optimisation problem of (9). The second one, dubbed MLm, is based on partitioning P(Z) into Nsub subrectangles (whose sides have a length not exceeding a threshold Δ) and on running the MATLABfmincon routine over each of the subrectangles; in this case, the final estimate is selected among Nsub candidates according to a minimum-cost criterion (note that a similar approach has been adopted in Algorithm 1). Finally, the third implementation, called MLc, employs thecompass algorithm (see Algorithm 3) under the assumption that Nr=1 and R1≡P(Z). Note that the complexity of each of the proposed implementations is ONobs·Neval, where Neval is the number of times the ML cost function is evaluated. The values of the parameters selected for each of them are listed in Table 1.Finally, a diagram summarising all the implementations developed in this section is provided in Figure2.Figure 2
Implementations of optimal map-aware and map-unaware estimators proposed in this paper.


## 4.1. Map-Aware Implementations

In this section, two different implementations for the evaluation of the MMSE estimate on the basis of (13) are proposed and analysed in detail. The first implementation, called MMSEc, is based on the use of the so-calledcubature formulas (e.g., see [17, 29–31]), which approximate any multidimensional definite integral as a finite sum of terms depending on its integrand. In particular, in the following, any integration over the search domain R(Z) is expressed as the sum of NrSRC distinct integrals, each referring to a distinct rectangular domain (see (1)), and the cubature formulas illustrated in [30] are exploited for the evaluation of each of the NrSRC integrals. Then, from (13) the approximate expression(19)p^MMSEcz≃∑k=1NrSRC∑j=1Nn,kwk,jαk,jxhz,αk,j∑k=1NrSRC∑j=1Nn,kwk,jαk,jyhz,αk,j∑k=1NrSRC∑j=1Nn,kwk,jhz,αk,jis easily inferred, where α(k,j)=[α(k,j)(x),α(k,j)(y)]T≜ck+Lkα(k,j)(F) is the jthnode referring to the kth rectangle Rk≜[lk;rk]×[bk;tk] in R(Z) (2), wk,j≜wk,j(F)detLk is the correspondingweight, ck≜1/2lk+rk;bk+tkT, and Lk≜1/2diagrk-lk,tk-bk, whereas α(k,j)(F)∈[-1;1]×[-1;1] and wk,j(F) denote the jthnode and the jthweight, respectively, of the cubature formula of order Nn,k. It is worth mentioning that, in principle, the larger the size of Rk is, the stronger the fluctuations of the integrand function to be expected over this domain are; consequently, high cubature orders (i.e., large values of Nn,k) should be selected for large rectangles, at the price, however, of an increase in the computational load. In our computer simulations, the heuristic formula Nn,k=NnmARk+dmin has been adopted for selecting the number of nodes; here, Nn(·) denotes the function defined in [30] to map the degree d of a polynomial into the number of integration nodes, mARk+dmin is the polynomial degree heuristically associated with the (nonpolynomial) integrand functions appearing in (13), ARk=detLk is the area of Rk, whereas m and dmin are real parameters whose values are listed in Table 1. It is also important to point out that (a) all the sums appearing in (19) require the evaluation of the same quantities hz,α(k,j) and (b) MMSEc is noniterative and its approximate complexity is ON¯s·Nobs·Neval, where Neval=Nn¯·NrSRC and Nn¯ denotes the average number of nodes selected for integrating over each of the NrSRC rectangles. Further details about our implementation of MMSEc can be found in [27, 32].Table 1
Covering and specific values selected for the implementations illustrated in Section4.
Algorithm
Cov.
Parameters
M
M
S
E
c
MNC
m
=
1.4
m
-
2, dmin=5
M
M
S
E
r
MNC
Δ
=
0.3
m
M
A
P
f
MNC
Algorithm = active-set, maxiter =103, tolfun = tolx = 10-2
M
A
P
g
MNC
Δ
1
=
2
m, Δfinal=0.3m, Nr,min=3
M
A
P
c
MNC
Δ
1
=
2
m, Δfinal=0.3m
D
R
D
-
E
c, DRD-Pg
MNC
γ
=
0.6, Nr,min=5
P
R
D
-
E
c, PRD-Pg
MNC
d
t
h
=
3.5
m, Nr,min=5
P
R
D
-
P
p
DNC
d
t
h
=
3.5
m, Nr,min=10
M
L
f
—
The same as forMAPf
M
L
m
—
Δ
=
10
m; others the same as for MAPf
M
L
c
—
Δ
1
=
2
m, Δfinal=0.3mThe second implementation, calledMMSEr, is based on the well-knowntrapezoidal integration rule; then, from (13), the approximate formula(20)p^MMSErz≃∑k=1NrSRC∑i=1Nn,kx∑j=1Nn,kywk,i,jαxk,i,jhz,αk,i,j∑k=1NrSRC∑i=1Nn,kx∑j=1Nn,kywk,i,jαyk,i,jhz,αk,i,j∑k=1NrSRC∑i=1Nn,kx∑j=1Nn,kywk,i,jhz,αk,i,jis easily inferred. Here, thenodes{αk,i,j≜lk+i-1Δ,bk+j-1ΔT,i=1,2,…,Nn,kx,j=1,2,…,Nn,ky} of the kth integration formula correspond to the vertices of a regular grid extending over Rk and characterized by a spacing Δ (see Table 1), whereas the corresponding weights wk,i,j are evaluated as(21)
w
k
,
i
,
j
≜
Δ
2
4
for
i
∈
1
,
N
n
,
k
x
,
j
∈
1
,
N
n
,
k
y
Δ
2
2
for
i
∈
1
,
N
n
,
k
x
,
j
∈
2
,
…
,
N
n
,
k
y
-
1
Δ
2
2
for
i
∈
2
,
…
,
N
n
,
k
x
-
1
,
j
∈
1
,
N
n
,
k
y
Δ
2
for
i
∈
2
,
…
,
N
n
,
k
x
-
1
,
j
∈
2
,
…
,
N
n
,
k
y
-
1with Nn,kx≜(rk-lk)/Δ+1 and Nn,ky≜(tk-bk)/Δ+1. It is important to point out the following: (a) similarly to (19), the three sums in (20) require the evaluation of the same quantities hz,α(k,i,j); (b) MMSEr is noniterative and its approximate complexity is ON¯s·Nobs·Neval, where Neval=Nnx¯·Nny¯·NrSRC and Nnx¯ (Nny¯) denotes the average number of integration points selected along the x (y) direction.In implementing the MAP estimator (14), the search over the space R(Z) has been turned into Nr searches over the distinct rectangles Rk,k=1,2,…,Nr covering R(Z). This choice is motivated by the fact that (a) the constraint p~∈R(Z), unlike p~∈Rk, cannot be directly formulated as a set of linear inequalities, as required by standard optimization techniques; (b) solving Nr distinct optimization problems, each involving a search domain much smaller than R(Z), strongly reduces the probability of reaching local minima; (c) when iterative optimization algorithms are adopted over each of the rectangles Rk, fast convergence can be usually achieved if their centres are selected as initial points. Note that in this case direct-search methods [19, 20] should be preferred to standard Newton-based or gradient-based methods, since the first class of methods usually admits simple implementations for constrained problems and does not require smooth cost functions. These considerations have led us to develop 3 different implementations, called MAPf, MAPg, and MAPc, of the MAP strategy (14). The first one (namely, MAPf) is summarised in Algorithm 1 and is based on the use of the MATLABfmincon routine to solve the optimization problem appearing in step (2). This implementation is characterized by various adjustable parameters; their values, which are listed in Table 1, have been selected to minimise the overall computational complexity.Algorithm 1:Implementation of the MAP strategy (14) based on the fmincon routine (MAPf).
Require: covering Rk,k∈K of R(Z) (2); fcp (15).
(1)fork∈Kdo
(2)pmin,k←arg⁡minp~∈Rkfcp~
(3)fmin⁡,k←fcpmin⁡,k
(4)end for
(5)k^←arg⁡mink{fmin,k}, p^MAPf←pmin,k^The second implementation (namely,MAPg) is summarised in Algorithm 2 and is based on discretisation of the search space and on the use of iterative direct-search methods (see [19, Sec. 3]). In practice, in its kth iteration, for each rectangle of the selected covering, the cost function fc(·) (15) is evaluated at the vertices of a rectangular grid characterized by spacing Δk (see lines (4)-(5)); then, on the basis of the resulting values, only the most promising rectangles are saved for the next iteration (see line (8)), in which a halved step size is adopted (see line (9)). Note that both the computational complexity and the accuracy of MAPg depend on the initial and the final values of the step size (denoted by Δ1 and Δfinal, resp.), so that a proper trade-off has to be achieved in the selection of these parameters (the values we adopted are listed in Table 1).Algorithm 2:Implementation of the MAP strategy (14) based on an iterative grid search (MAPg).
Require: covering Rk=[lk;rk]×[bk,tk],k∈K of R(Z) (2); fcp (15); Δ1, Δfinal and Nr,min (Table 1).
(1)Δ←Δ1, K~←K
(2)whileΔ>Δfinaldo
(3)fork∈K~do
(4)I←1,…,(rk-lk)/Δ, J←1,…,(tk-bk)/Δ
(5)(i^,j^)←arg⁡mini∈I,j∈Jfc(lk,bk)+Δ(i,j)
(6)pmin,k←(lk,bk)+Δ(i^,j^), fmin⁡,k←fcpmin⁡,k
(7)end for
(8)K~← indices of the Nr,min smallest elements of {fmin⁡,k}
(9)Δ←Δ/2
(10)end while
(11)k^←arg⁡mink∈K~{fmin⁡,k}, p^MAPg←pmin⁡,k^The third implementation (namely,MAPc) is summarised in Algorithm 3 and is based on a different direct-search method, known as thecompass algorithm [20, Sec. 8.1]. Similarly to MAPg, MAPc exploits a discretisation of the search space and an iterative reduction of the step size; however, the grid it employs is not regular. In fact, in each rectangle, the domain is explored moving from a “current trial position” p~ on the basis of given displacement vectors ±ex or ±ey (see line (6)) in order to identify the direction along which the cost function decreases (see line (8)). Even in this case, the overall computational complexity and the accuracy of this implementation depend on the initial and final values of the step size (see Table 1), which need to be carefully selected.Algorithm 3:Implementation of the MAP strategy (14) based on an iterative compass search (MAPc).
Require: covering Rk,k∈K of R(Z) (2); fcp (15); Δ1 and Δfinal (Table 1).
(1)fork∈Kdo
(2)Δ←Δ1
(3)pmin,k←centreRk, fmin,k←fcpmin,k
(4)whileΔ>Δfinaldo
(5)bImproved←False
(6)for d∈ex,ey,-ex,-eydo
(7)p~←pmin,k+Δ⋅d
(8)iffcp~<fmin,k and p~∈Rkthen
(9)fmin,k←fcp~, pmin,k←p~
(10)     bImproved← True, break
(11)end if
(12)end for
(13)Δ←Δ/2 only if bImproved = False
(14)end while
(15)end for
(16)k^←arg⁡minkfmin,k, p^MAPc←pmin,k^It is not difficult to show that the overall computational complexity of each of the MAP implementations proposed above isON¯s·Nobs·Neval, where Neval denotes the overall number of times the cost function fc(·) (15) is evaluated. Unluckily, Neval cannot be easily related to the other parameters of the proposed MAP implementations because of its nonlinear dependence on such parameters; consequently, a more accurate estimate of the computational burden required by the proposed implementations cannot be provided.Let us focus now on the problem of combining the implementation of the MMSE and MAP estimators illustrated above with the DRD and PRD techniques described in Section3.3. As far as the DRD technique is concerned, we note that all the proposed MAP/MMSE implementations can be employed in its last step. In the following, however, the implementations of DRD-MMSE and DRD-MAP techniques are always based on MMSEc and MAPg, respectively; for this reason, they are denoted by  DRD-Ec and DRD-Pg, respectively. It is important to point out the following:(1)
The algorithmsDRD-Ec and DRD-Pg are not intrinsically iterative. However, they can be easily modified to include a mechanism for iteratively updating the value of the parameter γ in order to improve their robustness against inaccuracies in the initial setup (see Algorithm 4).(2)
The complexity ofDRD-Ec and DRD-Pg is approximately O(NrSRC+Nobs(N¯s·Neval+Neval,raw)), where Neval denotes the overall number of evaluations of the MMSE integrand (for DRD-Ec) or that of the MAP cost function (for DRD-Pg), and Neval,raw represents the overall number of evaluations of the ML cost function.Algorithm 4:Implementation of the iterative DRD estimator (see Section 3.3) employed by DRD-Ec (DRD-Pg).
Require: covering Rk,k∈K of R(Z) (2); p^MMSEc (p^MAPg) estimator; constants γ and Nr,min (Table 1).
(1)p^raw←p^MLm, computed assuming Z^NLOS=Z
(2)dth←mediani∈Zγσb,i2p^raw(z)+σn,i2p^raw(z)
(3)dk←p^raw-projRk,p^raw for k∈K (see (17))
(4)whileK^<Nr,mindo
(5)K^←k∣dk<dth, dth←1.5dth
(6)end while
(7)p^DRD-Ec←p^MMSEc (p^DRD-Pg←p^MAPg) evaluated in the reduced map Rred(Z)=⋃k∈K^‍RkIn our work,MMSEc and MAPg have been also selected for the last step of the PRD technique; the resulting implementations are denoted by PRD-Ec and PRD-Pg, respectively, in the following. It is worth mentioning that(1)
bothPRD-Ec and PRD-Pg are not iterative;(2)
their complexity isO(N¯s·Nobs·(Neval,1+Neval,2)), where Neval,1 and Neval,2 denote the overall number of evaluations of the MMSE integrand (for PRD-Ec) and that of the MAP cost function (for PRD-Pg) in the first step and in the second step, respectively (see Section 3.4);(3)
the values of parametersdth and Nr,min need to be carefully selected.Our numerical results (see Section 5) have evidenced that whatever option for map-aware estimation is selected, the parameter N¯s (i.e., the average number of intersection tests accomplished in processing each observation) appearing in the complexity estimate of all the implementations developed in this paragraph plays a significant role in penalizing them with respect to their map-unaware counterparts. To reduce N¯s, the following two techniques have been devised:(1)
Segment Partitioning (SP). This consists of an offline phase, where (a) R is partitioned in multiple cells using a nonuniform rectangular grid (in a way that each anchor lies on a vertex of the grid) and (b) the obstruction segments falling inside each cell are identified. Then, when No(pia,p2) is evaluated, only the segments of the cells intersecting with the (pia,p2) segment are taken into account, thus reducing N¯s.(2)
Precomputing (PC). This is based on an offline evaluation of the function No(·,·) over Npt points (uniformly) distributed inside each rectangle of the available MNC or DNC; this requires storing Npt·N¯a·NrMNC or Npt·N¯a·NrDNC integers in a table, where N¯a denotes the average number of connected anchors. This technique is “memory hungry” for large maps; however, the precomputed table replaces that memorising information about the position of Ns obstruction segments and can be stored efficiently using a limited number of bits for the integer representation of the precomputed values of No(·,·). Moreover, this memory requirement is often more than counterbalanced by the reduction of computational complexity (since N¯s=0 in this case).In our computer simulations, the first technique has been adopted for all the map-aware estimators. The second technique, instead, has been used in the first step of PRD-MAP estimation only; the resulting implementation is dubbed PRD-Pp. This choice is motivated by the fact that this estimator can be easily adapted to the use of precomputed quantities, since it evaluates the map-aware likelihood (8) at a set of fixed points {pkc} (see Section 3.4).


## 4.2. Map-Unaware Implementations

Before illustrating our implementations of the ML estimator in (9), it is important to mention that if the constraint p~∈P(Z) was neglected, standard algorithms (e.g., theNelder-Mead simplex method [33]) could be employed to solve the optimization problem it involves. Unluckily, such a constraint plays an important role, specially in large maps, since (a) minima different from the global one are substantially less likely to be contained in the domain P(Z) than in the whole R2 space and (b) the map-unaware likelihood (8) is not defined for p~∉P(Z), as already explained in [15, Sec. II.C.]. For these reasons, in the 3 different implementations of the ML estimator (9) proposed below, the constraint p~∈P(Z) is always accounted for. The first implementation, called MLf, relies on the MATLABfmincon routine, which is employed to solve the optimisation problem of (9). The second one, dubbed MLm, is based on partitioning P(Z) into Nsub subrectangles (whose sides have a length not exceeding a threshold Δ) and on running the MATLABfmincon routine over each of the subrectangles; in this case, the final estimate is selected among Nsub candidates according to a minimum-cost criterion (note that a similar approach has been adopted in Algorithm 1). Finally, the third implementation, called MLc, employs thecompass algorithm (see Algorithm 3) under the assumption that Nr=1 and R1≡P(Z). Note that the complexity of each of the proposed implementations is ONobs·Neval, where Neval is the number of times the ML cost function is evaluated. The values of the parameters selected for each of them are listed in Table 1.Finally, a diagram summarising all the implementations developed in this section is provided in Figure2.Figure 2
Implementations of optimal map-aware and map-unaware estimators proposed in this paper.


## 5. Simulation Results

In this section, the implementations described in Section4 are compared in terms of accuracy and complexity. In our analysis, the three different maps, shown in Figure 3, are considered. Such maps are characterized by different sizes and obstruction densities (see Table 2), by a specific deployment of anchors (whose positions are identified by green squares in Figure 3), and by a similar number of anchors per unit area (1 anchor each ≃100 m2), so that a fair comparison to the data referring to different maps is ensured. In addition, in all cases, the anchor coverage radius (see Section 2.2) is dmax,i=dmax=30 m (with i=1,2,…,Na).Table 2
Values of some parameters characterizing the maps considered in Section5 and their MNC/DNC coverings. The number of anchors (Na) and intersection tests (Ns) is also provided for each map.
Map
N
r
M
N
C
N
r
D
N
C
N
s
A
R
N
a
1
340
632
574
540
m
2
6
2
529
925
1029
800
m
2
8
3
773
1461
1375
1310
m
2
13Figure 3
Support regions of map #1 (a), map #2 (b), and map #3 (c). Note that different scales have been used in the representation of these maps.
(a)
(b)
(c)The remaining part of this section is organized as follows. In Section5.1, the method employed to estimate the computational complexity of the proposed implementations is described. In Section 5.2, the structure of the software simulator developed to assess both the accuracy and the complexity of our implementations is illustrated. Finally, in Sections 5.3 and 5.4, various numerical results based on statistical models and experimental measurements, respectively, are analysed.


### 5.1. Evaluation of the Computational Complexity of Different Implementations

Map-aware implementations are expected to outperform their map-unaware counterparts at the price of a larger computational complexity. However, in this case, an accurate assessment of computational complexity is a challenging task, mainly because of the adaptive and data-dependent nature of the proposed algorithms; for instance, the termNeval (appearing in all the formulas summarised in Table 3) cannot be easily related to the size of the search domain or other implementation-specific parameters. For this reason, the computational complexity of each implementation has been assessed by estimating the average amount F¯ of FLOPs required in our computer simulations; this task has been accomplished under the following assumptions:(1)
Basic algorithmic steps (e.g., basic checks and assignments in Algorithms1–4 for map-aware estimation and LOS/NLOS detection in map-unaware estimation) play a negligible role.(2)
The computational complexity of each implementation mainly depends on the average numberN¯eval of evaluations of integrands (in MMSE estimators) or of cost functions (in MAP and ML estimators); as already pointed out in Section 4, the complexity of such evaluations depends, in turn, on N¯s and N¯obs. For this reason, in our simulations, counters for Neval, Ns, and Nobs have been used in order to estimate N¯eval, N¯s, and N¯obs, respectively (see Table 3).(3)
Since distinct types of operations are characterized by different computational weights, the expressionF¯=N¯eval∑t‍wtnt has been adopted, where nt is the number of FLOPs of type t performed in each evaluation of integrands or cost function and wt is the associated computational weight. In particular, in our simulations, the weights {3,3,3,2,2,1,0} have been assigned to square roots, exponentials, logarithms, divisions, multiplications, additions, and comparisons, respectively.Our approach is expected to provide statistically meaningful results about F¯ only if a large number (in practice, at least 104) of statistically independent measurements are available. Acquiring such a large number of experimental data would require a burdensome experimental campaign involving at least 150 measurement sites (so that 1502≃104 independent links would be available). To circumvent this problem, the computational complexity of the devised algorithms has been assessed mainly through computer simulations, following the approach always adopted in related work (e.g., see [22, 23]). However, some numerical results based on the available measurements are also shown (see Section 5.4) to support our conclusions.Table 3
Computational complexity of the implementations proposed in Section4 for the ML, MAP, and MMSE estimators.
Algorithm
Complexity
M
M
S
E
c
O
N
-
s
⋅
N
obs
⋅
N
eval
M
M
S
E
r
O
N
-
s
⋅
N
obs
⋅
N
eval
M
A
P
f, MAPg, MAPc
O
N
-
s
⋅
N
obs
⋅
N
eval
D
R
D
-
E
c, DRD-Pg
O
N
r
S
R
C
+
N
obs
N
-
s
⋅
N
eval
+
N
eval
,
r
a
w
P
R
D
-
E
c, PRD-Pg
O
(
N
-
s
⋅
N
obs
⋅
(
N
eval
,
1
+
N
eval
,
2
)
)
M
L
f, MLm, MLc
O
N
obs
⋅
N
eval


### 5.2. Description of the Developed Simulator

A software simulator has been developed (the MATLAB code used to generate the numerical results is publicly available (together with other results not included in this paper for space limitations) athttp://frm.users.sourceforge.net/publications.html, in accordance with the philosophy of reproducible research standard [34]) to test the proposed implementations in the harsh propagation conditions characterizing typical indoor environments [15]. In each simulation run, the following steps are accomplished sequentially:(1)
Generation of the Agent Position. The true agent position p is drawn from the selected map according to the rule p~f(p), where f(p)=1/AR for p∈R and f(p)=0 elsewhere (in other words, theuniform statistical map modelof [15] has been adopted).(2)
Selection of the Connected Anchors. The ith anchor is deemed to be “connected” with the agent only if the following conditions are jointly satisfied: (a) Nop,pia<Nw,max≜8, (b) σn,i2(p)+σb,i2(p)<σmax2≜30 m2, and (c) dip<d~max≜27 m. In addition, it is assumed that each agent-anchor link is characterized by a “failure probability” Pf≜0.1 (due, e.g., to synchronization or connection problems between the agent and the anchor itself).(3)
Generation of the Noisy Observations. The noisy observation vector z is drawn from the map-aware likelihood (5) according to rule z~f(z∣p~=p); the associated sets Z^LOS and Z^NLOS are randomly generated assuming the presence of errors with probability PeNLOS=0.1 (to account for the presence of a nonideal LOS/NLOS detector).(4)
Position Estimation. All the estimator implementations discussed in Section 4 are run sequentially; they are all fed by z and R(Z) (2) (map-aware estimators ) or P(Z) (4) (for map-unaware ones). Their output is either a position estimate p^(z) or a failure indication (occurring, e.g., when optimization routines do not converge).In our simulations, steps (1)–(4) have been repeated for at least Nruns=104 times in order to generate accurate estimates of the RMSE ϵ (which is more robust to outliers than rootmean square error) and of the mean number of FLOPs F¯ for each implementation.


### 5.3. Accuracy and Complexity Comparison

Figure4 shows the accuracy (in terms of RMSE ϵ) of the implementations discussed in Section 4 for the 3 considered maps. These results show the following:(1)
The accuracy achieved by the DRD-MMSE (DRD-MAP) and PRD-MMSE (DRD-MAP) estimator is similar to that of MMSE (MAP) estimators (e.g.,ϵPRD-Ec≃ϵDRD-Ec≃1.06ϵMMSEc for map #3).(2)
The implementationsMMSEc and MMSEr provide a 10%–15% performance improvement with respect to MAPf, MAPg, and MAPc (e.g., ϵMAPg≃1.16ϵMMSEc and ϵMAPg≃1.15ϵMMSEr for map #3).(3)
Map-awareness significantly improves estimation accuracy (up to 250%), as already evidenced in [15] (e.g., ϵMLm≃2.46ϵMMSEc and ϵMLm≃2.12ϵMAPg for map #3).It is important to point out that the estimation accuracy of the proposed estimators is mainly related to (a) the average spatial density of obstruction segments (i.e., to the ratio Ns/AR) and (b) the average spatial density of anchors (i.e., to the ratio Na/AR). Since the propagation scenario is similar in the three maps (light commercial buildings) and a similar anchor density has been adopted for the sake of a fair comparison (see Table 2), the RMSEs provided by the considered estimators are similar over the three maps and do not necessarily increase with AR (e.g., map #3 is larger than map #2, but a lower RMSE is achieved for map #3).Figure 4
RMSE performance provided by various ML, MAP, and MMSE estimation algorithms.Note also that Figure4 does not unveil some important aspects of the considered estimation algorithms, which emerge, instead, from the related box plot of Figure 5 describing the distribution of the RMSE outcomes for map #3. Here, the box represents the 75th percentile, the line within the box denotes the median, the black dashed whiskers are upper adjacents, and the cross markers represent the outliers. From Figure 5, the following can be easily inferred:(1)
The behaviour of all the estimators is characterized by the presence of multiple RMSE outliers associated withϵ>8 m as well as of upper adjacents and values of 75th percentiles extending well beyond their median values. This is due to the statistical model of Section 2, which reproduces harsh propagation conditions degrading RSS measurements in practical systems. Note that nonlinear filtering (e.g., particle filtering) could be employed to reject, at least in part, outliers when tracking moving agents.(2)
The implementationsMMSEc and MMSEr are characterized by the smallest dispersion of estimation errors (in fact, upper adjacents are at 7.3 m for both).Some numerical results about the complexity of the estimators implementations are shown in Figure 6 (note that a logarithmic scale is adopted for F¯ since complexities range roughly from 103 to 109 FLOPs). These results lead to the following conclusions:(1)
The MAP and MMSE estimators based on the DRD and PRD techniques are characterized by a computational complexity 4–11 times lower than that of their counterparts not employing such domain reduction techniques, especially for large maps (e.g.,F¯MAPg≃3.9F¯DRD-Pg and F¯MAPg≃11F¯PRD-Pp for map #3). This result can be related to the mean ratio between the areas Rred(Z) (16) and R(Z) (2), that is, to the amount of “domain reduction.” In fact, our simulations have evidenced that, in the case of map #3, on average, the DRD and PRD-based estimators restrict their search domain to 13% and 4% only, respectively, of the whole domain R(Z). This makes the complexities of DRD and PRD estimators much closer to that of map-unaware estimators (at the price of a negligible sacrifice in accuracy).(2)
Even if all map-aware implementations exhibit similar accuracy, they may require significantly different complexities. In fact, the parameterF¯ ranges from F¯MMSEr≃1.6·108 (for the most demanding map-aware algorithm) to F¯PRD-Pp≃1.1·106 (for the most computationally efficient map-aware algorithm), if map #3 is considered. Moreover, there is a significant gap between the FLOP counts of map-aware estimators and those of their map-unaware counterparts (e.g., F¯MMSEc≃4.1·103F¯MLm for map #3).(3)
Similarly to the RMSE, the parameterF¯ is significantly influenced by the size of R(Z) (2) and gets larger as the average spatial density of obstruction segments increases.(4)
PRD-based estimators are less computationally demanding than their DRD counterparts, specially in their MMSE variants (e.g.,F¯DRD-Ec≃3.1F¯PRD-Ec for map #3). Additional numerical results about the complexity of such algorithms (not shown here for space limitations) can be found in [27, 32].Figure 5
Box plot of the RMSEϵ provided by various ML, MAP, and MMSE estimation algorithms. Only map #3 is considered.Figure 6
FLOP countF¯ referring to various ML, MAP, and MMSE estimation algorithms.


### 5.4. Experimental Results

The conclusions inferred from Figures4–6 have been also supported by additional numerical results based on a set of experimental measurements referring to map #2. In fact, in the last part of our work, the simulator described in Section 5.2 has been modified in order to replace the computer-based generation of the noisy observation vector z (see step (3) of the simulator) with a random selection of observations stored in the experimental database described in [15]. Note that such a database, like other experimental databases publicly available, contains a limited set of statistically independent measurements (referring to 116 distinct wireless links) and this unavoidably reduces the overall number of independent simulations that can be run in processing such measurements; for this reason, from a statistical point of view, the results analysed in the previous Sections 5.2 and 5.3 are characterized by a higher level of confidence than those based on experimental measurements (see our comments at the end of Section 5.1). In practice, in this case, the number of simulation runs has been reduced to Nruns=103. The results shown in Figures 7 and 8 lead to the following conclusions:(1)
The RMSE achieved by map-aware estimators on the basis of experimental measurements is on the order of 2 m and is about 30% smaller than that of the simulated counterpart; these results can be related to the harsh propagation environment reproduced by the statistical model (5).(2)
The RMSE provided by map-unaware estimators on the basis of experimental measurements is close to that of its simulated counterpart; the good agreement between experimental and simulated results can be related to the availability of an ideal LOS/NLOS detector (actually, a small probability of errorPeNLOS=0.1 was assumed in Section 5.2).(3)
The computational complexity involved by the considered algorithms (both map-aware and map-unaware) is substantially the same in the case of simulated and experimental data; this means that the shape of the likelihood functions (5) and (8) does not significantly change substituting the experimental data with the observations generated through the statistical simulator.Finally, the possibility of generalising the proposed approach to other environments deserves some comments. In fact, in this study, three distinct scenarios, characterized by different maps and areas, but by similar propagation conditions, have been taken into consideration. Our measurements have led us to the conclusion that all the measurements acquired in these scenarios can be described by the same statistical model, which is that the same parameters can be exploited by map-aware and map-unaware estimators for each of the three considered environments. This suggests that the methodology adopted in this work can be easily adapted to other buildings characterized by offices, laboratories, and open spaces, like those encountered in light commercial buildings. Nonetheless, it is important to stress that, in principle, indoor environments characterized by different propagation conditions require a dedicated calibration phase for extracting a statistical model for measurements, as already mentioned in Section 2.4. Given such a model, the proposed techniques can be easily applied.Figure 7
RMSE performance provided by various ML, MAP, and MMSE estimation algorithms. These results, referring to map #2, are based on measurements.Figure 8
FLOP countF¯ for various ML, MAP, and MMSE estimation algorithms. These results, referring to map #2, are based on measurements.


## 5.1. Evaluation of the Computational Complexity of Different Implementations

Map-aware implementations are expected to outperform their map-unaware counterparts at the price of a larger computational complexity. However, in this case, an accurate assessment of computational complexity is a challenging task, mainly because of the adaptive and data-dependent nature of the proposed algorithms; for instance, the termNeval (appearing in all the formulas summarised in Table 3) cannot be easily related to the size of the search domain or other implementation-specific parameters. For this reason, the computational complexity of each implementation has been assessed by estimating the average amount F¯ of FLOPs required in our computer simulations; this task has been accomplished under the following assumptions:(1)
Basic algorithmic steps (e.g., basic checks and assignments in Algorithms1–4 for map-aware estimation and LOS/NLOS detection in map-unaware estimation) play a negligible role.(2)
The computational complexity of each implementation mainly depends on the average numberN¯eval of evaluations of integrands (in MMSE estimators) or of cost functions (in MAP and ML estimators); as already pointed out in Section 4, the complexity of such evaluations depends, in turn, on N¯s and N¯obs. For this reason, in our simulations, counters for Neval, Ns, and Nobs have been used in order to estimate N¯eval, N¯s, and N¯obs, respectively (see Table 3).(3)
Since distinct types of operations are characterized by different computational weights, the expressionF¯=N¯eval∑t‍wtnt has been adopted, where nt is the number of FLOPs of type t performed in each evaluation of integrands or cost function and wt is the associated computational weight. In particular, in our simulations, the weights {3,3,3,2,2,1,0} have been assigned to square roots, exponentials, logarithms, divisions, multiplications, additions, and comparisons, respectively.Our approach is expected to provide statistically meaningful results about F¯ only if a large number (in practice, at least 104) of statistically independent measurements are available. Acquiring such a large number of experimental data would require a burdensome experimental campaign involving at least 150 measurement sites (so that 1502≃104 independent links would be available). To circumvent this problem, the computational complexity of the devised algorithms has been assessed mainly through computer simulations, following the approach always adopted in related work (e.g., see [22, 23]). However, some numerical results based on the available measurements are also shown (see Section 5.4) to support our conclusions.Table 3
Computational complexity of the implementations proposed in Section4 for the ML, MAP, and MMSE estimators.
Algorithm
Complexity
M
M
S
E
c
O
N
-
s
⋅
N
obs
⋅
N
eval
M
M
S
E
r
O
N
-
s
⋅
N
obs
⋅
N
eval
M
A
P
f, MAPg, MAPc
O
N
-
s
⋅
N
obs
⋅
N
eval
D
R
D
-
E
c, DRD-Pg
O
N
r
S
R
C
+
N
obs
N
-
s
⋅
N
eval
+
N
eval
,
r
a
w
P
R
D
-
E
c, PRD-Pg
O
(
N
-
s
⋅
N
obs
⋅
(
N
eval
,
1
+
N
eval
,
2
)
)
M
L
f, MLm, MLc
O
N
obs
⋅
N
eval


## 5.2. Description of the Developed Simulator

A software simulator has been developed (the MATLAB code used to generate the numerical results is publicly available (together with other results not included in this paper for space limitations) athttp://frm.users.sourceforge.net/publications.html, in accordance with the philosophy of reproducible research standard [34]) to test the proposed implementations in the harsh propagation conditions characterizing typical indoor environments [15]. In each simulation run, the following steps are accomplished sequentially:(1)
Generation of the Agent Position. The true agent position p is drawn from the selected map according to the rule p~f(p), where f(p)=1/AR for p∈R and f(p)=0 elsewhere (in other words, theuniform statistical map modelof [15] has been adopted).(2)
Selection of the Connected Anchors. The ith anchor is deemed to be “connected” with the agent only if the following conditions are jointly satisfied: (a) Nop,pia<Nw,max≜8, (b) σn,i2(p)+σb,i2(p)<σmax2≜30 m2, and (c) dip<d~max≜27 m. In addition, it is assumed that each agent-anchor link is characterized by a “failure probability” Pf≜0.1 (due, e.g., to synchronization or connection problems between the agent and the anchor itself).(3)
Generation of the Noisy Observations. The noisy observation vector z is drawn from the map-aware likelihood (5) according to rule z~f(z∣p~=p); the associated sets Z^LOS and Z^NLOS are randomly generated assuming the presence of errors with probability PeNLOS=0.1 (to account for the presence of a nonideal LOS/NLOS detector).(4)
Position Estimation. All the estimator implementations discussed in Section 4 are run sequentially; they are all fed by z and R(Z) (2) (map-aware estimators ) or P(Z) (4) (for map-unaware ones). Their output is either a position estimate p^(z) or a failure indication (occurring, e.g., when optimization routines do not converge).In our simulations, steps (1)–(4) have been repeated for at least Nruns=104 times in order to generate accurate estimates of the RMSE ϵ (which is more robust to outliers than rootmean square error) and of the mean number of FLOPs F¯ for each implementation.


## 5.3. Accuracy and Complexity Comparison

Figure4 shows the accuracy (in terms of RMSE ϵ) of the implementations discussed in Section 4 for the 3 considered maps. These results show the following:(1)
The accuracy achieved by the DRD-MMSE (DRD-MAP) and PRD-MMSE (DRD-MAP) estimator is similar to that of MMSE (MAP) estimators (e.g.,ϵPRD-Ec≃ϵDRD-Ec≃1.06ϵMMSEc for map #3).(2)
The implementationsMMSEc and MMSEr provide a 10%–15% performance improvement with respect to MAPf, MAPg, and MAPc (e.g., ϵMAPg≃1.16ϵMMSEc and ϵMAPg≃1.15ϵMMSEr for map #3).(3)
Map-awareness significantly improves estimation accuracy (up to 250%), as already evidenced in [15] (e.g., ϵMLm≃2.46ϵMMSEc and ϵMLm≃2.12ϵMAPg for map #3).It is important to point out that the estimation accuracy of the proposed estimators is mainly related to (a) the average spatial density of obstruction segments (i.e., to the ratio Ns/AR) and (b) the average spatial density of anchors (i.e., to the ratio Na/AR). Since the propagation scenario is similar in the three maps (light commercial buildings) and a similar anchor density has been adopted for the sake of a fair comparison (see Table 2), the RMSEs provided by the considered estimators are similar over the three maps and do not necessarily increase with AR (e.g., map #3 is larger than map #2, but a lower RMSE is achieved for map #3).Figure 4
RMSE performance provided by various ML, MAP, and MMSE estimation algorithms.Note also that Figure4 does not unveil some important aspects of the considered estimation algorithms, which emerge, instead, from the related box plot of Figure 5 describing the distribution of the RMSE outcomes for map #3. Here, the box represents the 75th percentile, the line within the box denotes the median, the black dashed whiskers are upper adjacents, and the cross markers represent the outliers. From Figure 5, the following can be easily inferred:(1)
The behaviour of all the estimators is characterized by the presence of multiple RMSE outliers associated withϵ>8 m as well as of upper adjacents and values of 75th percentiles extending well beyond their median values. This is due to the statistical model of Section 2, which reproduces harsh propagation conditions degrading RSS measurements in practical systems. Note that nonlinear filtering (e.g., particle filtering) could be employed to reject, at least in part, outliers when tracking moving agents.(2)
The implementationsMMSEc and MMSEr are characterized by the smallest dispersion of estimation errors (in fact, upper adjacents are at 7.3 m for both).Some numerical results about the complexity of the estimators implementations are shown in Figure 6 (note that a logarithmic scale is adopted for F¯ since complexities range roughly from 103 to 109 FLOPs). These results lead to the following conclusions:(1)
The MAP and MMSE estimators based on the DRD and PRD techniques are characterized by a computational complexity 4–11 times lower than that of their counterparts not employing such domain reduction techniques, especially for large maps (e.g.,F¯MAPg≃3.9F¯DRD-Pg and F¯MAPg≃11F¯PRD-Pp for map #3). This result can be related to the mean ratio between the areas Rred(Z) (16) and R(Z) (2), that is, to the amount of “domain reduction.” In fact, our simulations have evidenced that, in the case of map #3, on average, the DRD and PRD-based estimators restrict their search domain to 13% and 4% only, respectively, of the whole domain R(Z). This makes the complexities of DRD and PRD estimators much closer to that of map-unaware estimators (at the price of a negligible sacrifice in accuracy).(2)
Even if all map-aware implementations exhibit similar accuracy, they may require significantly different complexities. In fact, the parameterF¯ ranges from F¯MMSEr≃1.6·108 (for the most demanding map-aware algorithm) to F¯PRD-Pp≃1.1·106 (for the most computationally efficient map-aware algorithm), if map #3 is considered. Moreover, there is a significant gap between the FLOP counts of map-aware estimators and those of their map-unaware counterparts (e.g., F¯MMSEc≃4.1·103F¯MLm for map #3).(3)
Similarly to the RMSE, the parameterF¯ is significantly influenced by the size of R(Z) (2) and gets larger as the average spatial density of obstruction segments increases.(4)
PRD-based estimators are less computationally demanding than their DRD counterparts, specially in their MMSE variants (e.g.,F¯DRD-Ec≃3.1F¯PRD-Ec for map #3). Additional numerical results about the complexity of such algorithms (not shown here for space limitations) can be found in [27, 32].Figure 5
Box plot of the RMSEϵ provided by various ML, MAP, and MMSE estimation algorithms. Only map #3 is considered.Figure 6
FLOP countF¯ referring to various ML, MAP, and MMSE estimation algorithms.


## 5.4. Experimental Results

The conclusions inferred from Figures4–6 have been also supported by additional numerical results based on a set of experimental measurements referring to map #2. In fact, in the last part of our work, the simulator described in Section 5.2 has been modified in order to replace the computer-based generation of the noisy observation vector z (see step (3) of the simulator) with a random selection of observations stored in the experimental database described in [15]. Note that such a database, like other experimental databases publicly available, contains a limited set of statistically independent measurements (referring to 116 distinct wireless links) and this unavoidably reduces the overall number of independent simulations that can be run in processing such measurements; for this reason, from a statistical point of view, the results analysed in the previous Sections 5.2 and 5.3 are characterized by a higher level of confidence than those based on experimental measurements (see our comments at the end of Section 5.1). In practice, in this case, the number of simulation runs has been reduced to Nruns=103. The results shown in Figures 7 and 8 lead to the following conclusions:(1)
The RMSE achieved by map-aware estimators on the basis of experimental measurements is on the order of 2 m and is about 30% smaller than that of the simulated counterpart; these results can be related to the harsh propagation environment reproduced by the statistical model (5).(2)
The RMSE provided by map-unaware estimators on the basis of experimental measurements is close to that of its simulated counterpart; the good agreement between experimental and simulated results can be related to the availability of an ideal LOS/NLOS detector (actually, a small probability of errorPeNLOS=0.1 was assumed in Section 5.2).(3)
The computational complexity involved by the considered algorithms (both map-aware and map-unaware) is substantially the same in the case of simulated and experimental data; this means that the shape of the likelihood functions (5) and (8) does not significantly change substituting the experimental data with the observations generated through the statistical simulator.Finally, the possibility of generalising the proposed approach to other environments deserves some comments. In fact, in this study, three distinct scenarios, characterized by different maps and areas, but by similar propagation conditions, have been taken into consideration. Our measurements have led us to the conclusion that all the measurements acquired in these scenarios can be described by the same statistical model, which is that the same parameters can be exploited by map-aware and map-unaware estimators for each of the three considered environments. This suggests that the methodology adopted in this work can be easily adapted to other buildings characterized by offices, laboratories, and open spaces, like those encountered in light commercial buildings. Nonetheless, it is important to stress that, in principle, indoor environments characterized by different propagation conditions require a dedicated calibration phase for extracting a statistical model for measurements, as already mentioned in Section 2.4. Given such a model, the proposed techniques can be easily applied.Figure 7
RMSE performance provided by various ML, MAP, and MMSE estimation algorithms. These results, referring to map #2, are based on measurements.Figure 8
FLOP countF¯ for various ML, MAP, and MMSE estimation algorithms. These results, referring to map #2, are based on measurements.


## 6. Conclusions

In this paper, specific solutions to the problem of implementing optimal map-aware and map-unaware estimators for indoor localization have been developed. These solutions are based on state-of-the-art numerical integration and optimization methods as well as novel techniques for reducing the size of search domains. Their accuracy and complexity in terms of RMSE and FLOP count, respectively, have been analysed through computer simulations relying on both statistical models (validated with measurements) and experimental databases. Numerical results have evidenced that (a) map-awareness can appreciably improve RMSE performance of RSS-based localization algorithms; (b) the proposed domain reduction techniques significantly reduce the computational complexity of map-aware algorithms; (c) the best accuracy-complexity trade-off is achieved by thePRD-Ep implementation, since it achieves nearly optimal RMSE performance and entails a computational burden much smaller than that of all the other map-aware algorithms in the considered scenarios; (d) experimental results substantially confirm the simulation results in terms of computational complexity, whereas the RMSE achieved over the experimental database is about 30% smaller than that expected by simulations because of the harsh propagation environment reproduced by the adopted statistical model.



---

*Source: 562680-2015-10-28.xml*