Battery powered asset trackers with a long device autonomy are a deep-rooted ambition for localisation solutions in retail, logistics and healthcare. This leaves a limited energy budget available for sensing, computing and communicating on the mobile sensor node. In this paper, we explore the attainable precision and accuracy in hybrid radio-acoustic, indoor 3D positioning for such energy-constrained devices. A state-of-the-art low-power ultrasonic-RF ranging technique is extended into a 3D indoor positioning system. The system is evaluated in a real-life measurement setup including over 30 000 measurements. 2D advanced positioning estimation algorithms are adapted for the 3D application, and an in-depth comparison is performed. Positioning in favourable conditions shows that the Euclidean distance error can be narrowed down to 0.0102 m. A road map with improvements on the hardware, algorithms and energy consumption is proposed to obtain this centimetre accuracy in the non-ideal situations. All data and algorithms are shared in open source, encouraging ongoing research to test and improve their innovative positioning methods.
Positioning solutions in retail, logistics and healthcare aim for low-power asset tracking, as
neglecting the power consumption results in a limited battery lifetime, higher maintenance
and hardware costs. In the context of these emerging applications, the constrained energy
budget should not counter the need for precise and accurate positioning. The approach for 3D
positioning in this paper comes from this application-driven perspective and tries to answer
the following question: what accuracy can be obtained with energy-constrained sensor nodes?
Active ultrasonics is considered an attractive solution. Thanks to the relatively low propagation
speed of sound, centimetre accurate positioning can be achieved without high-speed and energy
consuming processing. The low system complexity and emerging microelectromechanical
systems (MEMS) microphone technology enable the design of dedicated low-power and
lowLGOBE
cost acoustic sensor nodes. The first series of mature, ultrasound positioning systems were
presented as early as the year 2000. The leading-edge systems were the Active Bat system [
A short introduction to the system architecture, ranging method and measurement setup is presented in the next section. Section 3 explores the measurement results and describes the used positioning methods. Next, a comparison between these algorithms is performed, to be ifnalised with a road map towards centimetre accuracy and a conclusion. 2. Set-up
This study adopts the hybrid ranging system that is described in [
B1 Beacon
Mz r1
B2
RF Wake-Up Signal CAhciorupsStiicgnal
M2 r2
M1 r3 M3
M
Mobile Node
B3 The ranging information is thus comprised in the received audio signal and depends on the distance to the beacon. At the sensor node, high frequency, bit based sampling eliminates the usage of a power hungry ADC and limits the amount of data that needs to be transmitted over the energy consuming RF connection. Distances are calculated by using binary template cross correlation, based on the matching procedure of Sokal [18] and by performing a peakseeking procedure on the correlated data. The accuracy of this ranging method depends on the used chirp frequency bandwidth (Δ ), the mobile sensor node wake-up time ( ) and the binary sampling frequency ( , ). At the mobile sensor nodes, acoustic acquisition, amplification, binary sampling, packet handling and RF-communication are all done with 90.62 µJ. In this paper, we extend the aforementioned ranging method into a 3D positioning system and evaluate the performance in a challenging real-life environment. Note that this method difers from typical propagation delay or round-trip time measurements in respectively synchronised and unsynchronised systems. It calculates the distance out of the pulse compression data instead of using an internal timer. Nevertheless, it provides the necessary ranges to perform true-range multilateration algorithms. It is generally known that at least four ranges are necessary for 3D-indoor positioning [19], leading to a setup expansion with a minimum of four beacons. A Time Division Multiple Access (TDMA) scheme is chosen to obtain the corresponding distances, as this prevents acoustic signal collision eficiently at the cost of a lower positioning update rate. With this, the main objective of this paper is to find a positioning method with the smallest mean 3D Euclidean distance at all test positions in the acoustically harsh set-up.
All experiments are performed in the Techtile test environment [20] depicted in Fig. 2. This 8 x 6 x 2.4 m room, comprising 140 modular wooden tiles, represents a real-life, acoustically harsh environment, with high reverberation time (Table 1). Four acoustic beacons are spread across the test bed, with three of them close to the outer edges of the room, maximising the spatial resolution on all three axes. Each beacon consists of an ultrasonic speaker (Kemo L010) directed to the centre of the room and an of-the-shelf amplifier with a frequency range over 45 kHz. These four beacons are connected to a multi-channel DAC [21], where the chirp signals are generated and TDMA beacon selection is performed. The sensor node consists of a single MEMS microphone with distinctive amplifier and band filter combined with a sub-GHz wireless MCU launchpad [22]. A central computer connected to another launchpad handles the received binary sampled audio chirps and calculates the measured distances.
Measurements are performed on 150 fixed positions. Three layers with an uniform planar array of 50 points on fixed heights are created. On each position at least 50 distances are measured to each of the four beacons, resulting in a total of over 30 000 ranging measurements. These time expanded ranging measurement are necessary to average out errors introduced by the dynamic Techtile environment.
3.20 2.60 2.00 1.40 0.80 3.20 2.60 2.00 1.40 0.80
Beacon 1
Beacon 3
To understand the potential positioning inaccuracies and single out advanced positioning algorithms, we first explore the ranging measurements. A helpful data visualisation technique for errors between the actual and measured distance is a heat map, displaying the error range’s magnitude in a 2D projection of the room. Fig. 3 displays four such heat maps, one for each beacon, showing the Mean Absolute Error (MAE) of 50 positions on a fixed height (1.4 m). The speaker position and orientation are illustrated by a symbol. A couple of observations can be made. Firstly, the error increases with the distance. This can be attributed to signal attenuation and can be easily described by the inverse square law. In this system lower signal-to-noise ratios are detrimental for the zero-crossing based, binary template cross correlation. Secondly, larger errors can be found closer to the edge of the room. A potential hypothesis can be found in the high reverberation time of the room. And thirdly, in all four heat maps, a directional speaker pattern can be noticed. To clarify the latter two, a more elaborate statistical analysis and extra measurements are performed.
Table 2 gives more insight in the accuracy and precision of the distance measurements between a single beacon and a sensor node close to a wall. To measure the influence of reverberation, two sets of 140 ranging measurements were conducted on the same position, one with and another without acoustic absorbers added on the nearest wall. The mean ( ), diference between minimum and maximum value (range), median (P50), P95 and standard deviation ( ) values are all spectacularly improved when absorbers are added, cancelling out reflections and decreasing the reverberation locally. In this case, the low median and mean distance error values result in a high accuracy whilst the standard deviation of only 5.2 cm and range of 22.8 cm reveal a small error variation with low outliers and thus a high precision. To evaluate the influence of the speaker, a far-field, 2D speaker directivity measurement was performed. The normalised polar plot with an angular resolution of 1° for a sine wave of 25 kHz can be found in Fig. 4. A half power beam width of 40° confirms the directionality hypothesis.
180° 135° 225° 90° 270°
These ranging results are promising, with centimetre-accuracy given favourable circumstances. Yet it will be interesting how well the 3D positioning algorithms perform for mobile sensor nodes close to walls or out of the main acoustic lobes.
Range-based positioning optimisation is intrinsically finding a position estimate with the smallest error term fitting all the distance measurements. This is often achieved by solving a set of nonlinear equations in a least square sense. A fair simulation-based comparison of diferent iterative and closed-form nonlinear algorithms for 2D positioning can be found in [23]. We have adapted these algorithms for a 3D environment and extended the simulations with experimental data from ranging measurements. This section gives an overview of the implemented algorithms accompanied with a short explanation. For a full derivation and description we refer to the original papers. Notations are the same in all algorithms with: • is the position vector of beacon i containing its x, y and z coordinates • is the number of beacons (4) • is the ranging measurement to beacon i • ̂ is the position estimate The pseudocode of the algorithms can be found in Appendix A. The Python implementation of the distance measurement, positioning algorithms and the data sets can be found in [17].
Classical multilateration interpret the measured ranges as circles centred around the beacon. The combination of several of these circular ‘lines of positions’ (LOPs), generate an intersection. The latter is considered to represent the best estimate for the position of the mobile sensor node. However, this geometrical approach is susceptible to noise. As errors are introduced, the circular LOPs do not intersect in a specific point, but create an area with potential estimates. In [24], linear LOPs are created by drawing a straight line through the intersecting points of two ranging circles. The intersection of these newly generated lines indicates the position estimation. The advantage of this method is that it is still possible to generate linear LOPs when, due to errors, two or more circles do not intersect. The Moore-Penrose inverse, least squares estimate can be used to find an algorithmic rather than a geometric solution with a minimum number of linear LOPs of − 1 .
Bancroft [25] uses an algebraic, non-iterative technique to solve the GPS equations that is computationally eficient and numerically stable. In contrast to iterative techniques, higher order terms of the Taylor series expansions are not discarded, resulting in a better performance when beacons are poorly arranged. Although there is a lot of similarity with the Simple Intersection method, the Moore-Penrose inverse is not used for minimising the diference between the exact position and measurement residual, but as a projection operator in the n-dimensional space of these residuals. 3.2.3. Beck 3.2.4. Chueng Instead of a range-based least squares (R-LS) minimisation, a squared-range-based least squares method (SR-LS or least quartic) can be used to find a position estimate [ 26]. An exact, optimal solution can be found by using a bisection algorithm to find a single root of a univariate strictly monotonous function. This method shows a better mean squared position error in comparison with the R-LS method.
An alternative way to compute the SR-LS estimate is presented in [27]. For a 3-dimensional space, all roots of a seven root equation are calculated by rewriting it in a seventh degree polynomial. In contrast with the original paper, the absolute of the roots are taken instead of eliminating complex roots. The root closest to zero is used to compute the estimate of the globally optimal constrained weighted least squares.
Also known as Weighted Iterative Least Squares (WILS) is a straightforward method using only algebraic operations, in the form of a Jacobian matrix or the derivative of the measurement equation, to find an estimation. The disadvantages of this method are the need for prior information in the form of a initial positioning guess and the fact that it is computationally harder. Sirola [23] suggests the simple intersection method as a starting point for the position estimation. When used in time series, the algorithm extends naturally to a Kalman filter.
To compare the diferent 3D position estimation algorithms, we have generated the Cumulative Density Function (CDF) of the mean Euclidean error distance of 50 measurement in each of the 150 positions for the diferent methods, as plotted in Fig. 5. At first glance, it is clear that implementing more advanced estimation algorithms greatly improves the position estimation. When taking a closer look, we can divide the plots into three groups. The first one is the Chueng method, showing the steepest incline and outperforming all the other algorithms. With half the measurements below a distance of 0.73 m and 95% just above the 1.54 m, it demonstrates to be accurate even with erroneous ranging measurements. The three methods belonging to the second group are Bancroft, Beck and Gauss-Newton. They have a similar curvature and a P90 value below 3.5 m, with the Gauss Newton method having a single outlier at 69.8 m. The last curve consists of the basic intersection method, with a low accuracy and high susceptibility to noise (P90 of 10.34 m).
1 0.8 0.6 F D C 0.4 0.2 0 0 Simple Intersection Bancroft Beck Chueng
Gauss-Newton
To investigate the influence of the ranging error on the precision, the position estimates of the two best performing position estimation techniques are plotted in Fig. 6. Two scenario’s are depicted, with the yellow dot representing the actual position. The sensor node position in the first scenario is closer to the wall, with larger average distance measurement errors to two beacons: 0.923 m and 0.626 m to respectively Beacon 1 in the upper left corner, and Beacon 4, centrally at the far right. The second sensor node position is more centred, and has a maximum mean ranging error below 5 cm. The Chueng algorithm behaves similarly in both scenarios: forming a dense cloud in the neighbourhood of the actual position. The mean Euclidean distance error is 0.490 m for the first scenario and 0.554 m for the second. This indicates a lower accuracy, while at the same time having a high precision for both scenarios, hence the dense cloud. The Bancroft forms a more scattered arrangement, performing worse on both accuracy and precision. Even with precise distance measurements, the estimates spread out in the z-direction, where the beacon range resolution is lowest.
Position 8
Room Setup Range Bancroft
Speaker Chueng
Actual Position
Position 34
Room Setup Range Bancroft
Speaker Chueng
Actual Position (a) (b)
In the most favourable conditions, the Euclidean distance is 0.0102 m, indicating that the intended
accuracy can be achieved. Undeniable, there is still a road to go to reach this centimetre accuracy
in non-ideal conditions, and the error to be bridged is at least an order of magnitude. For more
consistent exact location estimates of energy-constraint nodes, we propose the following
promising upgrades:
• Hardware and system based improvements. To prevent signal attenuation due to
speaker directivity, the system could be extended on a hardware level by adding extra
speakers to the beacons. Care should be taken in the design, as constructive interference
and spatial aliasing occur at higher frequencies when multiple speakers are used, resulting
in non-uniform and, again, directional sound pressure fields [ 28]. Beacon placement
should be optimised as well, potentially increasing the number of beacons. Even with
precise ranging measurements, positioning errors are largest in the z-dimension of the
room. The resolution on this axis can be increased by placing a beacon in the ceiling
or embed it in the floor. When such beacons are added to the system, adaptive beacon
selection algorithms could be implemented to minimise the amount of active anchors and
increase the position reliability.
• Improve, combine or implement new positioning algorithms. The chosen
optimisation algorithms can be considered as safe and credible. They are an excellent way to test
the feasibility of the system at the cost of accuracy and precision. We strongly believe that
progress can be made if this research is extended with more recent positioning methods.
Advanced position estimation techniques could be combined. Algorithms like non-linear
least square estimation and Gauss-Newton require an initial positioning guess. In this
paper, the poorly performing simple intersection method is used as a starting point as
suggested by [23]. Future research could explore the Chueng algorithm as a starting point
or investigate the potential accuracy gain of time-series Kalman filters.
• Reduce the distance measurement error. Current processing uses the index of
maximum correlation peak for the distance calculations. [
In the present study, we extended a state-of-the-art RF-acoustic ranging strategy into a lowpower, 3D positioning system. Four existing advanced 2D position estimation algorithms have been adapted to operate in 3D. Their performance has been assessed and compared based on the actual distance measurement data in an acoustically challenging environment. With these constrains, a median Euclidean distance error below 0.75 m was obtained. The current TDMA implementation results in a sensor node consumption of only 362.45 µJ for a single position estimate. In further research we aim to improve the energy consumption and positioning accuracy by adapting the infrastructure, combining or implementing (new) optimisation algorithms and researching multiple access methods. The data of over 30 000 ranging measurements and applied position estimation algorithms are shared openly [17], encouraging ongoing research to test and improve their innovative positioning methods. Letters 4 (2020) 1–4. doi:1 0 . 1 1 0 9 / L S E N S . 2 0 2 0 . 2 9 9 0 2 1 3 . [17] B. Cox, Hybrid 3D-positioning, 2021. URL: https://github.com/BertCox/
Hybrid3DPositioning. doi:1 0 . 5 2 8 1 / z e n o d o . 3 2 6 9 7 4 7 . [18] R. R. Sokal, The Principles and Practice of Numerical Taxonomy, Taxon 12 (1963) 190 – 199. [19] A. Kupper, Location-based services : fundamentals and applications., Wiley, Chichester, 2005. [20] G. Callebaut, J. V. Mulders, G. Ottoy, L. V. der Perre, A Primer on Techtile: An R&D Testbed for Distributed Communication, Sensing and Positioning, 2021. a r X i v : a r X i v : 2 1 0 5 . 0 6 7 4 0 . [21] N. Instruments, USB-6212: 16 AI (16-Bit, 400 kS/s), 2 AO (250 kS/s), Up to 32 DIO USB
Multifunction I/O Device, National Instruments, 2017. Rev. [22] T. Instruments, SimpleLink Sub-1 GHz CC1310 wireless MCU LaunchPad development kit,
Texas Instruments, 2018. Rev. A. [23] N. Sirola, Closed-form algorithms in mobile positioning: Myths and misconceptions, in: 2010 7th Workshop on Positioning, Navigation and Communication, 2010, pp. 38–44. doi:1 0 . 1 1 0 9 / W P N C . 2 0 1 0 . 5 6 5 3 7 8 9 . [24] J. Cafery, A new approach to the geometry of TOA location, in: Vehicular Technology Conference Fall 2000. IEEE VTS Fall VTC2000. 52nd Vehicular Technology Conference (Cat.
No.00CH37152), volume 4, 2000, pp. 1943–1949 vol.4. doi:1 0 . 1 1 0 9 / V E T E C F . 2 0 0 0 . 8 8 6 1 5 3 . [25] S. Bancroft, An Algebraic Solution of the GPS Equations, IEEE Transactions on Aerospace and Electronic Systems AES-21 (1985) 56–59. doi:1 0 . 1 1 0 9 / T A E S . 1 9 8 5 . 3 1 0 5 3 8 . [26] A. Beck, P. Stoica, J. Li, Exact and Approximate Solutions of Source Localization Problems,
IEEE Transactions on Signal Processing 56 (2008) 1770–1778. doi:1 0 . 1 1 0 9 / T S P . 2 0 0 7 . 9 0 9 3 4 2 . [27] K. Cheung, H. So, W.-K. Ma, Y. Chan, Least squares algorithms for time-of-arrival-based mobile location, IEEE Transactions on Signal Processing 52 (2004) 1121–1130. doi:1 0 . 1 1 0 9 / T S P . 2 0 0 4 . 8 2 3 4 6 5 . [28] C. Quested, A. Moorhouse, B. Piper, B. Hu, An analytical model for a dodecahedron loudspeaker applied to the design of omni-directional loudspeaker arrays, Applied Acoustics 85 (2014) 161–171. doi:h t t p s : / / d o i . o r g / 1 0 . 1 0 1 6 / j . a p a c o u s t . 2 0 1 4 . 0 3 . 0 2 3 . [29] F. Yuan, Q. Wei, E. Cheng, Multiuser chirp modulation for underwater acoustic channel based on VTRM, International Journal of Naval Architecture and Ocean Engineering 9 (2016). doi:1 0 . 1 0 1 6 / j . i j n a o e . 2 0 1 6 . 0 9 . 0 0 4 . [30] C. Bernard, P.-J. Bouvet, A. Pottier, P. Forjonel, Multiuser Chirp Spread Spectrum Transmission in an Underwater Acoustic Channel Applied to an AUV Fleet, Sensors (Basel, Switzerland) 20 (2020) 1527. doi:1 0 . 3 3 9 0 / s 2 0 0 5 1 5 2 7 .
Algorithm 1 Simple Intersection
with −1 the identity matrix gives: Algorithm 2 Bancroft ifnd the roots of: = ( )
−1 = ( )
−1 1 ‖‖ 2 2 + (2 − 1 ) + ‖‖ 2 = 0,
̂ = + , = 1, 2 and from the two position candidates: pick the one with the smallest residual as the final estimate.
Algorithm 3 Beck Use a bisection method to obtain a solution for in the declining interval [−1/ 1, ∞], where
The estimate ̂ is given by the first 3 elements of .̂ Algorithm 4 Chueng
Compute and Λ from the eigenvalue decomposition with the eigenvalues of the matrix ( Ψ−1 ) = ⎢⎢⎢ 00 ⎥⎥⎥ Find the absolute of the complex root ∗ closest to zero of the following seven-root equation 0 = 4 4 − 2 4 4 + ∑ 3 + ∑ + 2 3
∑ =1 (1 + ) 4 =1 (1 + )2 .
=1 1 +
2 =1 1 + − 3 ∑
3 2 −
∑ 2 =1 3 ( + ) (1 + ) 2 and compute the estimate of the constrained weighted least squares with the first three elements of ̂ as the position estimate .̂ Algorithm 5 Gauss-Newton . As a 3: Set +1 = + Δ , where Δ the least square solution is of:
1 1 (Σ− 2 (ℎ ( )− )) Δ = − (Σ− 2 ) 4: If stopping condition ‖Δ ‖ < is not satisfied and ≤ max, increment and repeat from Step 2.