Elektrotehniški vestnik 83(5): 236-242, 2016 Original scientific paper Trajectory estimation of a moving target from Ultra-wideband ranging measurements Ziga Zupanec, Fabio Ricciato, Luka Sajn Faculty of Computer and Information Science, University of Ljubljana, Vecna pot 113, 1000 Ljubljana, Slovenia Abstract. In this paper we propose a method for trajectory estimation of a moving node based on minimizing the residual sum defined as the difference between a reported and actual distance from the anchor nodes. We devise extensive and complex indoor experiments with exploratory data analysis and interpretation of the results. Our findings show a slight improvement over an existing point-based localization system in an indoor environment. Keywords: Indoor positioning system (IPS), Ranging measurements, Non-linear least squares (NLS), Non-line-of-sight (NLOS), Ultrawideband (UWB), blind node, localization error, mobile node, target node, velocity estimation Ocenjevanje trajektorije primikajoče tarče z merjenjem razdalj v UWB V prispevku predstavljamo rekonstrukcijo začetne pozicije ter hitrosti in smeri premikajoče se tarče v realnem okolju. Predlagamo postopek, ki ta dva vektorja izračuna na podlagi ocenjenih razdalj med premikajočo se tarčo in skupkom stačionarnih oddajnikov. Bistvo postopka je optimizačija kriterijske funkčije z metodo najmanjših kvadratov. Opravili smo obsezne eksperimente tako na prostem, kot v zaprtem prostoru. Podrobno smo analizirali dobljene rezultate in jih poko-mentirali. Pri eksperimentih na prostem je napaka med končnim mestom trajektorij in dejansko lokačijo tarče manjša od enega metra pri razdalji 16m. V zaprtem prostoru smo rahlo izboljsšali natančšnost obstoječšega komercialnega sistema za določanje lokačije. 1 Introduction An indoor positioning system (IPS) is a system that determines the position (ločation) and/or movement of an obječt within a čonfined spače. The system čonsists of a set of fixed, statič nodes čalled ančhors and an undetermined position of a blind node roaming freely within the čonfined spače. Perhaps the most widely known positioning system is the Global Positioning System (GPS) [1]. This system employs a network of earth orbiting satellites and terrestrial GPS rečeivers. Eačh GPS rečeiver determines its ločation based on data rečeived from the satellites [2]. If we draw a parallel to Received 12 September 2016 Accepted 11 November 2016 GPS, the anchor nodes of IPS resemble satellites and the blind node resembles a GPS receiver device. Anchor nodes, blind node and processing unit are integral parts of IPS. There are cases where one cannot use GPS since it relies on a direct satellite-signal reception [3]. For instance, indoor environments, such as warehouses and shopping malls, underground facilities, such as basements, tunnels, underground mines, etc. One can also consider our solution as a replacement for GPS where power consumption is of concern. Yet another case of usage are devices that are not equipped with a GPS receiver or have the GPS receiver switched off most of the time. The aforementioned situations are the reason to look closer at an indoor localization system. A multitude of technologies is suitable for obtaining an indoor position. These technologies include and are not limited to WLAN [4], [5], [6], RFID [7], cellular-based [8], UWB [9], BT [10] and others. Liu et al. published a research of various commercially available wireless-based indoor positioning systems and solutions, which suggests that the most widely used technology uses WLAN's RSS to determine the target's position [11]. A common practice of determining the position of a moving target is to condense the information into one single point of the location. This method yields successive points of the position along the path of a target movement. Because of various errors affecting the point-based localization, the tracking path is often nonlinear and does not reflect the actual movement of the target [5]. Our intention is to improve the positioning accuracy TRAJECTORY ESTIMATION OF A MOVING TARGET FROM ULTRA-WIDEBAND RANGING MEASUREMENTS 237 of a blind node by enhancing the Multilux's LIPS — point-based positioning — system with a trajectory estimation of a linearly moving node. Multilux Ltd. is a Slovenian company that develops indoor positioning systems. Table 1. Example of a measurement packet. 2 Theory 2.1 Scenario Imagine a moving object within a room bound by four walls (Fig. 1). The object lies somewhere in the timestamp anchor ID ranging status 1452967702506 1 21.49 OK 1452967702506 2 12.64 OK 1452967702506 3 0.00 ERR 1452967702506 4 2.80 OK 1452967702506 5 13.51 OK 1452967702506 6 11.29 OK collinear anchors and four non-collinear anchors to get the blind node trajectory in a 3D space. 2.2 Technology All necessary hardware used in our experiments was provided by Multilux Ltd*. Their proprietary point-based positioning solution is called LIPS. The LIPS communication devices function at the 2.45 GHz frequencies, which belong to the designated ISM [13]. They operate at very low energy levels and use ranging sensors developed by Nanotron Technologies. A blind node pings the anchor nodes in a round-robin fashion performing two-way ranging measurements (Fig. 2). Nanotron Technologies have developed a time-of-flight Figure 1. Reference scenario showing a moving blind node surrounded by four anchors. room with initial position P and is moving towards a wall with constant velocity v. To determine the object current position, we multiply its velocity v by elapsed time t and add the result to the object's initial position Pt = P + vt. The moving object is called a (moving) blind node because its properties are a priori unknown. We can determine the blind node properties by using ranging measurements (reported by the blind node) obtained from pinging anchor nodes. The anchor nodes are an integral part of the model. Their position is fixed and known to the outside observer. The ranging measurement is the distance between a blind node and an anchor. It is the path of a signal travelling from one end to another. Signal velocity vs is the speed at which the signal wave carries data, vs = -j^, where c is the speed of light in a vacuum and er is the relative dielectric constant of the air. Knowing the signal velocity and TOF, f, which is the time a signal needs to travel from a transmitter to a receiver measured by a common time reference [12], the signal path is computed as follows: d = vsf. To obtain the blind node trajectory in a 2D space, we need the ranging measurements from at least three non- PC Figure 2. Abstract scheme of a blind node 'A' performing ranging measurements from the anchor nodes (magenta color). method that employs a ranging signal sent by a reader and an acknowledgment sent back from the tag to cancel out the requirements for clock synchronization [14]. The ranging measurement is a distance between the blind node and an anchor. The blind node periodically sends data to a processing unit that determines its position and trajectory. The provided ranging information is as close to the "raw data" as we can get without accessing the physical layer or the device firmware (which we have no control of). Being able to access the internal variables (f) *http://multilux.eu 238 ZUPANEC, RICCIATO, SAJN and modify the device to support broadcasting mode, we employ a differential time-difference of the arrival (DTDoA) approach that gives us more data from all nearby devices [6]. Working with (only) the ranging data pose a constraint that is considered in our problem formulation. A blind node creates a measurement packet every 130 ms or about eight times per second. The packet contains a relative timestamp, reflecting the time of its formation and an n-tuple with a ranging information, where n is the number of anchor nodes deployed in a given scenario (Table 1). The ranging information is a distance/status pair for ¿-th anchor node. A stream of packets is then sent to a processing unit via a serial-over-USB connection. The processing unit utilizes the NLS formula to estimate the blind node trajectory. 2.3 Problem formulation The trajectory estimation is a non-linear least-squares problem defined as a set of given ranging measurements over a period of time to determine the blind-node initial position P = [x,y,z]J and its velocity v = [vx,vy ,vz ]T. With an array of ranging measurements D recorded for the duration of a single walk (one linear trajectory), the problem can be expressed as an NLS optimization function: k [Po, V] = argmin^ ^ (||P + vtm - Pk\\ - dk,mf P ,v k= 1 mEDk (1) where: • P is the initial position of the blind node; • v is the velocity vector of the blind node; • K is the number of deployed anchors; • Dk is an array of recorded measurements for k-th anchor; • tm is the relative time of m-th measurement starting from t0 = 0; • Pk is the position of k-th anchor; • dk,m is the distance between the blind node and k-th anchor taken from m-th measurement. This function computes such P0 and v where the squared difference between measured distance dk,m and the computed distance based on time ||P + vtm|| to all the anchors is minimal. This problem can be solved with a MATLAB unconstrained solver fminunc. The solver uses a Quasi-Newton method to find the local minimum. To speed up the process, we resort to FMUG. The described routine can be used in MATLAB as: The least-squares solution represents the maximum likelihood estimate for the Gaussian error with a zero mean [4], [5]. This error is caused by delays in data processing between communication devices. If an error is strictly positive, the cause is attributed to NLOS Code 1 The use of fminunc in MATLAB._ 1:x0 = [0, 0, 0, 0, 0, 0]; 2:% Initial guess [x_0, y_0, z_0, v_x, v_y, v_z]. 3: 4:fun = @(y)objfun_gradient(y, anchor_positions , D); 5:% Gradient of the objective function with additional parameters . 6: 7:options = optimoptions(@fminunc, 'SpecifyObjectiveGradient' , true, 'Algorithm', ' quasi—newton'); 8:% Tell fminunc to use cutom gradient and Quasi—Newton method. 9: 10:[ traj_est , fval] = fminunc(fun, x0, options); 11:% Solve NLS. between the node and an anchor (heavy multipath error). In order to handle this error, we developed several heuristic methods that work in conjunction with our resolution routine. 2.4 Error model In practical experiments, the measured distance is not equal to the ground-truth distance due to various elements. We analyzed two types of the physical error: • small, symmetric distance error (sources are hardware delays and multipath) • large, positive distance error (source is NLOS with a heavy multipath). A multipath is a phenomenon where the signal path between a receiver and transmitter comes from direct and multiple indirect sources (Fig. 3). A multipath error depends on the amplitude, phase and on the received signal itself [16]. In the line-of-sight (LOS) conditions, the strongest path of the signal typically corresponds to the first path, while in the NLOS conditions, weak components typically precede the strongest path [17]. The cases where there is no direct signal and the cases where indirect signals are predominant are labelled a heavy multipath. A heavy multipath is particularly present in indoor environments. There are no formal methods for detecting NLOS from the ranging measurements. To tackle this problem, we employ three different heuristic methods. 2.5 Pruning heuristics In order to trim the erroneous data, a number of different statistical methods was employed in the process. Simpler methods based on the arithmetic mean and the median proved to be inferior to the more sophisticated ones. As a result, the three pruning heuristics described in this chapter served as our filtering methods. The 5th percentile heuristic and the residual pruning method are strictly pre- and postprocessing mehtods, respectively, while the robust mean method functions in both steps. Depending on the quality of the gathered data, each TRAJECTORY ESTIMATION OF A MOVING TARGET FROM ULTRA-WIDEBAND RANGING MEASUREMENTS 239 pruning heuristic displayed various degrees of accuracy. The preliminary research on the pruning methods showed that other studies mostly depend on the Kalman filter [15] for its simplicity. In our case, the Kalman filter is not a viable solution since it is very sensitive to outliers and large errors [5]. _""""--■--..■.B Figure 3. Multipath — the signal path between transmitter A and receiver B taking various routes. 2.5.1 5th percentile of distances: This method was developed as a more robust method compared to the statistical median method due to the distribution of the reported distances from an anchor node. A limitation of this method is that it can only be applied to segments where the target node is stationary. The distance errors are obtained by subtracting the measured data from the ground truth and used to form a cumulative distribution function (CDF) for each anchor. Given the properties of a physical environment, especially indoors, we anticipate a bimodal distribution of the collected data (Fig. 4). The 5th percentile of the CDF distribution was chosen to avoid the negative distance error component between devices on the lower end and to cut off measurements with positive distance error components caused by NLOS (Fig. 5). 10 15 reported distance [metres] Figure 4. Graphical representation of a bimodal distribution. 2.5.2 Robust mean over n-tuple: Clustering and outlier rejection further reduce the error caused by NLOS as well as errors from other sources. This is based on the fact that the direct path signals of multiple anchor nodes will localize the blind node close to the ground-truth location, while the reflection path signal will localize the blind node at random locations [16]. For all experiments conducted in this work a 3-tuple, (g) formation was selected for the six anchors. This provided us with 20 combinations of triplets in the Figure 5. CDF plot of a bimodal distribution from Figure 4. preprocessing step. A point-based localization was then used to form a cloud of points for each triplet. In the postprocess, a centroid was determined among all the points as an arithmetic mean value over all clouds. The furthest 1 % of the points was eliminated and a new centroid was elected [6]. This process was iterated until there were 50 % of the initial number of the points left. A visual representation of the outcome of the described procedure is depicted in Figure 6. 9 10 11 distance [metres] Figure 6. Robust mean over a triplet with a 50 % retention rate (green cluster) and iteration dropout of 1 %. 2.5.3 Residual pruning: Residual pruning is a method that repeats the resolution routine until a stop condition is met. After its initial application, the outcome are two vectors — one is initial position P0, the other one is velocity vector v. The residuals used in the resolution routine are then unrolled, sorted and 0.5% of the packet measurements (0.5% of m from D), representing the largest residuals, are discarded. This process is iteratively repeated until there are 80 % to 90 % packet measurements left in D. The expected outcome of this method is to give a better approximation to the ground truth trajectory by identifying and eliminating the outliers after localization. 10 30 25 20 4 15 12 13 10 0 0 20 240 ZUPANEC, RICCIATO, SAJN 2.5.4 Pruning heuristics analysis: Of the three described heuristic methods, only one is used for a given computation. The best performing method is determined by the quality of the input data which is further described in the following chapters. Figure 7. Localization pipeline. Preprocessing heuristics require distance data as input. Processed distance data DD is then piped into a resolution routine together with anchor positions A and the estimated initial position of the target node. The result of the entire pipeline is target node initial position P0 and its velocity v. The speed of the pipeline depends on the heuristic method properties, the most defining ones being the amount of data and the possible iterative nature of the heuristic. 3 Results and Conclusion Indoor experiments were conducted in an empty faculty cafeteria. This location was chosen because it is spacious and has a high ceiling. Six anchors were placed on the ceiling with an antenna facing downwards. The anchors were at different heights (3.9 m - 4.5 m) because the roof is inclined towards one end. The coverage area of an anchor can be calculated as follows: R = H tan (60), where R is the radius of the area covered on the ground and H is the anchor height. All anchor positions were projected on the floor using a laser measurer and marked with an 'X' in order to get the anchors exact heights. The distance between anchor A5 and anchor A6 served as a referential segment. Other anchors coordinates were obtained by measuring their perpendicular distance onto our referential segment. Coordinates were, as before, submitted to the blind nodes. An IV pole served as an improvised moving object with two blind nodes attached at the height of 1.3 m (the L node) and 2 m (the H node), respectively. The lower blind node was mounted with its antenna facing perpendicular to the ground, while the upper blind node antenna faced upwards. A tour is a sequence of segments. Three different tours were designed and conducted at two different paces, resulting in six independent experiments. The idea was to test the areas covered with many anchors as well as the areas with obstructions, corners and similar which is reflected in the tour patterns. The keypoints of the tours were labeled with markers on the ground. The coordinates of these keypoints were assigned with the aid of the aforementioned referential segment. The IV pole was moved along these segments at a constant speed at a slow and a faster pace. Ten-second stops were made at each keypoint. Simultaneously, the ranging information was captured from the blind nodes. At the stationary periods, the collected data was analyzed using the three pruning heuristics. Figure 8 is a graphical representation of the results. Polar charts map the distance error from the ground-truth location at each stationary point. Each symbol represents a different statistical pruning heuristic. The greatest deviations from the ground-truth were found in the red and green stationary periods. This was attributed to a low anchor coverage. Overall, the best results were obtained by the residual pruning heuristic, represented with an asterisk. While the other two heuristics were inferior to the residual pruning heuristic, they both yielded comparable results. The moving periods for Experiment 2 are represented in Figure 9. When moving toward the keypoints with a poor anchor coverage, the trajectory errors increased. In general, the direction of the trajectories reflected the ground-truth, whereas the trajectory scalar values were overestimated. The best optimization of the full data trajectories was achieved using the residual pruning heuristic with the retention rate of 80 %. This work elaborates the trajectory-based localization. The goal was to determine the blind node trajectory, defined as a function of its initial position and velocity in a given time frame. The resolution routine, the input of which are the ranging measurements from the industrial, scientific and medical radio bands, used to estimate the trajectory parameters, is an optimization of the nonlinear least-square model. Additional statistical pruning heuristics were combined with our resolution routine to filter the biggest outlying measurements. The indoor environment presented a greater set of challenges for the system in determining the location of the blind node. Various elements affected the accuracy of the data regarding the blind node position. Overall, the most resilient heuristic of residual pruning turned out to be the most accurate for the stationary as well as moving segments. With this heuristic, 90% of all distance errors were under 4 m, which is in range with indoor positioning systems that are currently available on the market. The heuristics used during the experiments were able to handle small and symmetric distance error anomalies. Although there was a moderate success with the described pruning heuristics in combating the lack of LOS in a harsh indoor environment, it is up to the manufacturer to filter out the multipath signals as the phenomenon is detected on the physical level. Regarding our initial prediction that our 5th percentile heuristic would perform better in handling a bimodal distribution, we are now able to confirm its designated D TRAJECTORY ESTIMATION OF A MOVING TARGET FROM ULTRA-WIDEBAND RANGING MEASUREMENTS 241 Figure 8. Static position estimation for the H node at the height of 2 m — Experiment 2 (fast pace). Residual pruning retention rate: 90 %, robust mean with a 50 % retention rate. --Full data trajecto .......»«• Residual pruning --»- GT trajectory O Start point x Mid point • End point 3 3 1 2 1 2 0 0 "vrt— "T7" \ I / X ■ * x : * / ; ! / / / i / i » / ■■ ' 1 / ' 3 4 3 4 22 11 »••»-.; 0 • -.......................................................................... ____" -------B ---j2___^ ^ _____^ J>____ _I_I_I_l_ 5 10 15 20 distance [metres] 30i— 25 - 20 15 ^ 10 5 - 0 - 0 25 Figure 9. Experiment 2 (fast pace) trajectories with the H node (height: 2 m). Residual pruning retention rate: 90 %. Opacity indicates the course of the experiment, trajectories fade with time. 242 ZUPANEC, RICCIATO, SAJN function. Still, this heuristic is greatly inferior to the residual pruning heuristic. For further work on this topic we suggest that our localization pipeline be used with other wireless technologies. Since our testing equipment was very expensive and small scale, using more widely implemented technologies, i.e. cellular networks, bluetooth, WLAN, etc., would provide us with a valuable comparison for determining appropriateness of these low-cost solutions as IPS. [17] S. Marano, W. M. Gifford, H. Wymeersch, M. Z. Win, NLOS identification and mitigation for localization based on UWB experimental data, IEEE Journal on Selected Areas in Communications 28.7, IEEE, 2010, pp. 1026-1035. References [1] Official U.S. Government information about the Global Positioning System (GPS) and related topics, http://www.gps.gov [2] B. Hofmann-Wellenhof, H. Lichtenegger, J. Collins, Global positioning system: theory and practice, Springer Science & Business Media, 2012, Ch. Overview of GPS, pp. 11-25. [3] A. Kleusberg, R. B. Langley, The limitations of GPS, GPS World 1.2, 1990. [4] A. Coluccia, F. Ricciato, On ML estimation for automatic RSS-based indoor localization, in: 2010 5th IEEE International Symposium on Wireess Pervasive Computing (ISWPC), IEEE, 2010, pp. 495Z502. [5] A. Coluccia, F. Ricciato, Maximum likelihood trajectory estimation of a mobile node from RSS measurements, in: 2012 9th Annual Conference on Wireless On-demand Network Systems and Services (WONS), IEEE, 2012, pp. 151-158. [6] N. Facchi, F. Gringoli, F. Ricciato, A. Toma, Emitter localisation from reception timestamps in asynchronous networks, Computer Networks 88, Elsevier, 2015, pp. 202-217. [7] H. Zou, H. Wang, L. Xie, Q.-S. Jia, An RFID indoor positioning system by using weighted path loss and extreme learning machine, in: 2013 IEEE 1st International Conference on Cyber-Physical Systems, Networks, and Applications (CPSNA), IEEE, 2013, pp. 66-71. [8] B. S. Lakmali, D. Dias, Database correlation for GSM location in outdoor & indoor environments, in: 2008 4th International Conference on Information and Automation for Sustainability, ICIAFS 2008, IEEE, 2008, pp. 42-47. [9] M. Kuhn, C. Zhang, B. Merkl, D. Yang, Y. Wang, M. Mahfouz, A. Fathy, High accuracy UWB localization in dense indoor environments, in: 2008 IEEE International Conference on UltraWideband, ICUWB 2008, Vol. 2, IEEE, 2008, pp. 129-132. [10] C. Frost, C. S. Jensen, K. S. Luckow, B. Thomsen, R. Hansen, Bluetooth indoor positioning system using fingerprinting, in: Mobile lightweight wireless systems, Springer, 2011, pp. 136150. [11] H. Liu, H. Darabi, P. Banerjee, J. Liu, Survey of wireless indoor positioning techniques and systems, IEEE Transactions on Systems, Man, and Cybernetics, Part C: Applications and Reviews 37.6, IEEE, 2007, pp. 1067-1080. [12] S. Lanzisera, D. Zats, K. S. Pister, Radio frequency time-of-flight dis- tance measurement for low-cost wireless sensor localization, IEEE Sensors Journal 11.3, IEEE, 2011, pp. 837-845. [13] Multilux Ltd., LIPS Local Indoor Positioning System (2015). URL https://www.multilux.eu/docs/MRFUM01.pdf [14] Nanotron Technologies GmbH, Real time location systems (rtls), in: A White Paper from Nanotron Technologies GmbH, Nanotron Technologies GmbH, 2007. URL http://www.nanotron.com/EN/pdf/WP_RTLS.pdf [15] Y. Bar-Shalom, X. R. Li, T. Kirubarajan, Estimation with applications to tracking and navigation: theory algorithms and software, John Wiley & Sons, 2004. [16] J. Xiong, K. Sundaresan, K. Jamieson, Tonetrack: Leveraging frequency-agile radios for time-based indoor wireless localization, in: Proceedings of the 21st Annual International Conference on Mobile Computing and Networking, ACM, 2015, pp. 537549. Ziga Zupanec received his B.Sc. and M.Sc. degrees in computer science from the University of Ljubljana, Slovenia in 2013 and 2016, respectively. He is currently working on an application to help physicians annotate medical histories with the WHO's ICF classification. Fabio Ricciato received a Laurea degree in 1999 and a doctoral degree in 2003 from the University La Sapienza of Rome, Italy. His current research interests are in various topics in the area of Communication Networks & Wireless systems. His recent research activity involves problems related to Network Measurements, Traffic Monitoring, Wireless Localization, Wireless Sensor Networks, RFID, Network Security & Privacy, and others. Luka Sajn received his PhD from the Faculty of Computer and Information Science, University of Ljubljana, Slovenia in 2007. His research interests are in: cell counting, multi-resolution pattern parametrization, ischaemic heart-disease diagnosing from scintigraphic images and whole-body bone scintigraphy segmentation. Currently, he is a researcher and an assistant professor at the same faculty.