IEEE INTERNET OF THINGS JOURNAL, VOL. 7, NO.

4, APRIL 2020 3133

Kalman-Filter-Based Integration of IMU

and UWB for High-Accuracy Indoor
Positioning and Navigation
Daquan Feng , Chunqi Wang, Chunlong He , Member, IEEE, Yuan Zhuang , Member, IEEE,
and Xiang-Gen Xia, Fellow, IEEE

Abstract—The emerging Internet of Things (IoT) applications, Index Terms—Extended Kalman filter (EKF), indoor
such as smart manufacturing and smart home, lead to a huge positioning system (IPS), inertial measurement unit (IMU),
demand on the provisioning of low-cost and high-accuracy Internet of Things (IoT), ultrawideband (UWB), unscented
positioning and navigation solutions. Inertial measurement unit Kalman filter (UKF).
(IMU) can provide an accurate inertial navigation solution in a
short time but its positioning error increases fast with time due to
the cumulative error of accelerometer measurement. On the other
hand, ultrawideband (UWB) positioning and navigation accuracy I. I NTRODUCTION
will be affected by the actual environment and may lead to uncer- OW-COST and high-accuracy positioning and navigation
tain jumps even under line-of-sight (LOS) conditions. Therefore,
it is hard to use a standalone positioning and navigation system
to achieve high accuracy in indoor environments. In this article,
L solutions for indoor mobile robots have become critical
in the Internet of Things (IoT) applications, such as smart
we propose an integrated indoor positioning system (IPS) com- manufacturing and smart home [1]. The inertial navigation
bining IMU and UWB through the extended Kalman filter (EKF) system (INS) is based on kinematics and Newton classical
and unscented Kalman filter (UKF) to improve the robustness mechanics [2]. The core of the INS is the inertial measure-
and accuracy. We also discuss the relationship between the geo-
metric distribution of the base stations (BSs) and the dilution of ment unit (IMU), which consists of a three-axis accelerometer
precision (DOP) to reasonably deploy the BSs. The simulation and a three-axis gyroscope [3], [4]. The IMU can obtain the
results show that the prior information provided by IMU can attitude information and motion characteristics of the car-
significantly suppress the observation error of UWB. It is also rier, such as acceleration, angular velocity, and angle [5].
shown that the integrated positioning and navigation accuracy Without using any reference base stations (BSs), the posi-
of IPS significantly improves that of the least squares (LSs) algo-
rithm, which only depends on UWB measurements. Moreover, tion of the carrier can be directly calculated by mathematical
the proposed algorithm has high computational efficiency and integrations of acceleration. Because of its low cost, low envi-
can realize real-time computation on general embedded devices. ronmental impact and high accuracy in a short time period,
In addition, two random motion approximation model algorithms INS has been widely used in mobile object positioning and
are proposed and evaluated in the real environment. The exper- navigation scenarios, such as aircrafts, vehicles, and pedestri-
imental results show that the two algorithms can achieve certain
robustness and continuous tracking ability in the actual IPS. ans, but errors increase rapidly with time [6]. On the other
hand, many researchers have considered to adopt the ultraw-
Manuscript received August 2, 2019; revised November 11, 2019 and ideband (UWB) technology in the indoor positioning system
December 9, 2019; accepted December 25, 2019. Date of publication (IPS) [7] and lots of work has been done, including channel
January 9, 2020; date of current version April 14, 2020. This work was
supported in part by the National Natural Science Foundation of China model [8], multipath component estimation [9], and theoretical
under Grant 61701317, in part by the Young Elite Scientists Sponsorship lower band of positioning errors [10]. UWB is a communica-
Program by CAST under Grant 2018QNRC001, in part by the Guangdong tion technology that uses nanosecond nonsinusoidal narrow
Natural Science Foundation under Grant 2017A030310371, in part by the
Shenzhen Overseas High-Level Talents Innovation and Entrepreneurship under pulse signal to transmit data, it has become an effective trans-
Grant KQJSCX20180328093835762, in part by the Tencent Rhinoceros Birds- mission technology in location-aware sensor networks [11].
Scientific Research Foundation for Young Teachers of Shenzhen University, Inherently, the UWB-based ranging technology has the advan-
in part by the Natural Science Foundation of SZU, and in part by the Start-up
Fund of Peacock Project. (Corresponding author: Chunlong He.) tages of short pulse interval and high time resolution and can
Daquan Feng, Chunqi Wang, and Chunlong He are with the Guangdong achieve centimeter-level ranging accuracy [12]. In addition,
Province Engineering Laboratory for Digital Creative Technology and it has good robustness to against the multipath effect [13].
Guangdong Key Laboratory of Intelligent Information Processing, Shenzhen
University, Shenzhen 518060, China (e-mail: fdquan@gmail.com; chun- However, due to the high frequency band of UWB, it is only
longhe@163.com). suitable for line-of-sight (LOS) conditions. When there exit
Yuan Zhuang is with the State Key Laboratory of Information Engineering opaque objects, the raging accuracy will be greatly reduced.
in Surveying, Mapping and Remote Sensing, Wuhan University, Wuhan
430079, China. Therefore, if only IMU or UWB is used, it is difficult to
Xiang-Gen Xia is with the College of Electronics and Information achieve high accuracy in complex indoor environments. To this
Engineering, Shenzhen University, Shenzhen 518060, China, and also with the end, the researcher has considered to take advantages of their
Department of Electrical and Computer Engineering, University of Delaware,
Newark, DE 19716 USA. complementary characteristics to improve the positioning and
Digital Object Identifier 10.1109/JIOT.2020.2965115 navigation accuracy [14]. A multisensor fusion architecture
2327-4662 
c 2020 IEEE. Personal use is permitted, but republication/redistribution requires IEEE permission.
See https://www.ieee.org/publications/rights/index.html for more information.
Authorized licensed use limited to: SOUTH CHINA UNIVERSITY OF TECHNOLOGY. Downloaded on May 13,2025 at 07:21:27 UTC from IEEE Xplore. Restrictions apply.
 3134 IEEE INTERNET OF THINGS JOURNAL, VOL. 7, NO. 4, APRIL 2020

is proposed for the IPS in [15]. When UWB is available

and reliable, the long-term drift error of IMU is calibrated
by UWB. Otherwise, the system will switch from UWB to
IMU immediately to get the positioning and navigation status.
In [16], an integrated positioning solution of IMU and UWB
is proposed, which can provide reliable and continuous posi-
tion, especially in the case of non-LOS (NLOS) conditions.
However, in [16], to address the time synchronization issue, Fig. 1. System architecture of IMU and UWB integration.
all the UWB modules are connected to a central controller
via the fiber line. This will increase the complexity of the
deployment and the cost of the system. Both [15] and [16] and angle measurement. Compared with the traditional
use the time of arrival (TOA) ranging algorithm and have three BSs least squares (LSs) positioning algorithm, the
not considered the clock drift effect. Moreover, although the proposed algorithm can greatly reduce the complex-
acceleration noise from IMU is considered in the accelera- ity and cost of BS deployment when the positioning
tion measurements, the effect of acceleration noise on the accuracy requirement is not high.
velocity and displacement due to the integration is not con- 2) UKF fusion positioning algorithm based on single BS
sidered in [15] and [16]. In [17], a drift-free and real-time and three BSs are proposed to reduce the computational
localization and tracking algorithm combining IMU, UWB and complexity. The simulation results show that the UKF
region-specific sensor is proposed. It is shown that the fusion fusion algorithm can achieve better positioning accuracy
algorithm can obtain more accurate 3-D velocity and height than the EKF fusion positioning algorithm.
information by using the lower body biomechanical model. 3) AUAM and AUM approximate motion models are
Similarly, a magnetometer-free lower body motion capture proposed, which can greatly reduce the jitter of posi-
(MOCAP) algorithm is introduced in [18], where it combines tioning data. The experimental results show that the
IMU with the UWB positioning system and the human lower positioning trajectories of the two models are smoother,
body biomechanical model for 3-D positioning and attitude and the performance of AUM is better than that of
tracking. Note that most of the existing indoor positioning AUAM.
and navigation systems are based on the human lower body In Section II, the basic principles of IPS based on IMU
MOCAP system [19], which increases the cost and deploy- and UWB are described. In Section III, the framework and
ment difficulty. Moreover, although the drift error in yaw angle algorithm of the integrated positioning system based on IMU
estimation can be eliminated by the usage of the magnetome- and UWB are proposed. In Section IV, simulation and exper-
ter in [20], it limits the applications under the long-standing imental results are presented and positioning performances of
magnetic interference in the indoor environment. different algorithms are compared and discussed. In Section V,
Therefore, to achieve low cost and high accuracy of the some conclusions are given.
indoor positioning and navigation system, this article focuses
on the integration of IMU and UWB based on Kalman fil- II. P OSITIONING AND NAVIGATION S YSTEM BASED ON
ter (KF) algorithms. The data obtained by IMU are used for THE I NTEGRATION OF IMU AND UWB
the state equation while the data obtained by UWB are used
In this section, we first describe the architecture of the
for the observation equation of KF. In UWB ranging, to bet-
indoor positioning and navigation system based on the inte-
ter solve the clock synchronization and clock drift problems,
gration of IMU and UWB, and then introduce the IMU-based
we propose an enhanced asymmetric double-sided two-way
and UWB-based positioning algorithms.
ranging (EADS-TWR) algorithm. Specifically, in order to fur-
Fig. 1 describes the architecture of an indoor positioning
ther improve the indoor positioning accuracy and reduce the
and navigation system based on the integration of IMU and
deployment cost of BSs, the extended KF (EKF) algorithm
UWB. The IMU sensor consists of a three-axis accelerometer
and unscented KF (UKF) algorithm based on three BSs and
and a three-axis gyroscope. The UWB sensors consist of an
single BS are adopted in our scheme, respectively. In addi-
unknown position tag and three BSs with known positions.
tion, to reduce the jitter of positioning data, two approximate
motion models, namely, approximate uniform motion (AUM)
and approximate uniform acceleration motion (AUAM), are A. Positioning Processing Algorithms Based on IMU
proposed to make the positioning results smoother and more Much work has been done on the transformations of dif-
stable. Considering that the positioning accuracy is not only ferent coordinate systems and the update calculation of the
affected by the distance measurement accuracy but also by the attitude matrix in strapdown INS (SINS) [21]. This article
geometric distribution of the BS, we discuss the influence of mainly introduces the basic principles, common coordinate
the geometric distribution of BSs and the dilution of precision systems, and related coordinate transformation matrix theory.
(DOP) on the positioning accuracy. The main principle of the SINS is shown in Fig. 2.
The main contributions of this article can be summarized To describe the space motion state of the carrier well,
as follows. it is necessary to choose an appropriate coordinate system.
1) We propose a direct positioning algorithm (DPA) based Assume that the body coordinate system (the B system) is
on EKF fusion for single BS with single distance represented by Oxb yb zb and the navigation coordinate system

Authorized licensed use limited to: SOUTH CHINA UNIVERSITY OF TECHNOLOGY. Downloaded on May 13,2025 at 07:21:27 UTC from IEEE Xplore. Restrictions apply.
 FENG et al.: KF-BASED INTEGRATION OF IMU AND UWB FOR HIGH-ACCURACY INDOOR POSITIONING AND NAVIGATION 3135

Let Cn1 , C12 , and C2b denote the basic rotations from the
N system to the first rotation system, from the first rotation
system to the second rotation system and from the second
rotational system to the B system, respectively. Then, the
coordinate transformation matrix from the N system to the
B system, Cnb , is expressed as follows:

Cnb = C2b C12 Cn1 (1)

⎡ ⎤
1 0 0
Fig. 2. Principle of SINS. C2b = ⎣ 0 cos γ sin γ ⎦ (2)
0 − sin γ cos γ
⎡ ⎤
cos θ 0 − sin θ
C12 = ⎣ 0 1 0 ⎦ (3)
sin θ 0 cos θ
⎡ ⎤
cos ψ sin ψ 0
Cn1 = ⎣ − sin ψ cos ψ 0 ⎦. (4)
0 0 1
On the other hand, the transformation matrix between two
Cartesian coordinate systems is an orthogonal matrix. Thus
 −1  T
Fig. 3. Determination of carrier space angle position. Cbn = Cnb = Cnb (5)
⎡ ⎤
cos ψ cos θ C12 C13
(the N system) is represented by Oxn yn zn as shown in Fig. 3. Cbn = ⎣ cos θ sin ψ C22 C23 ⎦ (6)
The origin coincides with the center of mass of the carrier. − sin θ sin γ cos θ cos θ cos γ
In the B system, the xb -axis and yb -axis point to the front
along the longitudinal direction and the right along the trans- where Cbn is the coordinate transformation matrix from the B
verse direction, respectively, and the zb -axis is vertical along system to the N system and
the carrier and forms the right-hand coordinate system with
the xb -axis and yb -axis. In the N system, the east–north–down C12 = sin γ sin θ cos ψ − cos γ sin ψ
(END) coordinate system is selected. Specifically, the xn -axis C13 = sin γ sin ψ + cos γ sin θ cos ψ
and the yn -axis refer to the east and north directions in the C22 = cos ψ sin γ + sin γ sin θ sin ψ
local horizontal plane, respectively, and the zn -axis points to C23 = sin ψ sin θ cos γ − sin γ cos ψ.
the down along the vertical line of the earth to form the right-
hand system. The coordinate transformation matrix between In the system, attitude updating refers to the real-time cal-
the B system and the N system is the attitude matrix of the culation of matrix Cbn based on the output of IMU. Let ωnb
carrier. The coordinates of the accelerometer and gyroscope in denote the angular velocity of the B system relative to the N
the system belong to the B system. The final output of accel- system. Then, the components of ωnb in the B system, ωnb b , is

eration, velocity, and position belong to the N system. Note given as follows:
that both Euler angle and quaternion methods can be used for ⎡ b ⎤ ⎡ ⎤ ⎡ ⎤ ⎡ ⎤
the attitude updating in our IPS. However, considering that the ωnbx 0 0 γ
b ⎦
⎣ ωnby = C2b C12 ⎣ 0 ⎦ + C2b ⎣θ ⎦ + ⎣ 0 ⎦. (7)
Euler angle method is more intuitive and easy to understand
than the quaternion method, we adopt the Euler angle method ωnbz
b −ψ 0 0
for the attitude transformation between the N system and the Then, by expanding and merging, the Euler angle differential
B system. Taking the N system as the reference coordinate equation is obtained as follows:
system. The heading angle of the carrier is Yaw (expressed
⎡ ⎤ ⎡ ⎤−1 ⎡ b ⎤
in ψ), the pitch angle is Pitch (expressed in θ ), and the roll γ 1 0 − sin θ ωnbx
angle is Roll (expressed in γ ). The parameters ψ, θ , and γ ⎣ θ ⎦ = ⎣0 cos γ sin γ cos θ ⎦ ⎣ ωnby b ⎦
. (8)
are a set of Euler angles, which describe the carrier space ψ 0 − sin γ cos γ cos θ ωnbz
b
angular position relationship between the N system and the
B system as shown in Fig. 3. When coordinates are rotated The acceleration in the B system, ab , is measured by the
with Euler angles, the product of matrices cannot be exchanged three-axis accelerometer
since different products represent different rotation orders. The  T
transformation matrix is the multiplication of the transforma- ab = abx aby abz . (9)
tion matrices determined by the basic rotations which will be Thus, the acceleration in the N system, an1 , is obtained by the
mathematically shown in detail below. The sequence of the coordinate transformation
multiplication is arranged from right to left in the order of the  T
basic rotations. an1 = an1 x an1
y an1
z = Cbn ab . (10)

Authorized licensed use limited to: SOUTH CHINA UNIVERSITY OF TECHNOLOGY. Downloaded on May 13,2025 at 07:21:27 UTC from IEEE Xplore. Restrictions apply.
 3136 IEEE INTERNET OF THINGS JOURNAL, VOL. 7, NO. 4, APRIL 2020

If the gravity vector g is removed from the acceleration of the

N system, then the acceleration in the N system, an , is obtained
⎡ n ⎤ ⎡ n1 ⎤ ⎡ ⎤
ax ax 0
an = ⎣ any ⎦ = ⎣ an1
y
⎦ − ⎣ 0 ⎦. (11)
anz an1
z g
When the sample interval is short, the carrier is approxi-
mately subjected to constant force and with uniform acceler-
ation of linear motion. Let vn denote the velocity variation
in the N system. Then, by substituting Newton’s second law
for momentum conservation equation, we have
⎡ n⎤ ⎡ n⎤
vx ax
⎣ vny ⎦ = ⎣ any ⎦t. (12)
vnz anz
Let vn (t) denote the velocity in the N system at time t. Then,
the velocity in the N system at time t + 1, vn (t + 1), can be
obtained by the acceleration integral as follows:
⎡ n ⎤ ⎡ n ⎤ ⎡ n⎤
vx (t + 1) vx (t) vx
⎣ vny (t + 1) ⎦ = ⎣ vny (t) ⎦ + ⎣ vny ⎦. (13)
vnz (t + 1) vnz (t) vnz
Let χ n denote the displacement variation in the N system. Fig. 4. Double-sided two way ranging protocol. (a) SDS-TWR protocol.
Then (b) EADS-TWR protocol.
⎡ ⎤ ⎡ n⎤ ⎡ n⎤
χxn vx ax
1
⎣ χyn ⎦ = ⎣ vny ⎦t + ⎣ any ⎦t2 . (14)
2 n TreplyA and TreplyB as the time delay of devices A and B,
χzn vnz az
respectively. Then, according to Fig. 4(b), it is easy to see
Let χ n (t) denote the position in the N system at time t, and that TroundA and TroundB can be expressed as follows:
then the position in the N system at time t + 1, χ n (t + 1), is
described as follows: TroundA = 2Tprop + TreplyB
⎡ n ⎤ ⎡ n ⎤ ⎡ ⎤
χx (t + 1) χx (t) χxn TroundB = 2Tprop + TreplyA . (16)
⎣ χyn (t + 1) ⎦ = ⎣ χyn (t) ⎦ + ⎣ χyn ⎦. (15) Thus, we can get the TOF between devices A and B as follows:
χzn (t + 1) χzn (t) χzn
TroundA × TroundB − TreplyA × TreplyB
Tprop = . (17)
B. Positioning Processing Algorithm Based on UWB TroundA + TroundB + TreplyA + TreplyB
The wireless signal-based positioning methods are generally Then, distance d can be expressed as
divided into the ranging-based and nonranging-based methods.
The UWB positioning method usually adopts the ranging- d = cTprop (18)
based algorithm and has two steps. The first step is to measure where c is the propagation velocity of the electromagnetic
the distance and angle information, and the second step is wave.
to calculate the position using the measurement distance and After completing the ranging process, we adopt the multi-
angle information. lateration to determine the position of the tag based on the
Time of flight (TOF) ranging is the common method to distance measurements. Specifically, Fig. 5 shows an example
measure the distance between two nodes. Most existing works of trilateration, where the distances from the tag to three BSs
in TOF adopt the symmetric double-sided two-way ranging are measured. Obviously, in the 2-D plane, the tag should be
(SDS-TWR) technology [22] and as shown in Fig. 4(a) to located at the intersection of three circles centered on three
mitigate the influence of clock synchronization between nodes. BSs. As long as the three BSs are not in a straight line, the
However, the frequency drift caused by the crystal clock drift result of the trilateration is unique.
cannot be solved. Besides, the SDS-TWR technology requires Assuming that the unknown tag is located at (x, y) and the
a long processing time [23]. To overcome the above problems, ith BS is located at (xi , yi ). Then, the true distance between
we propose an improved EADS-TWR optimization algorithm the unknown tag and the ith BS, di , can be expressed as
in our model as shown in Fig. 4(b).
Define Tprop as the TOF between device A and device di = (xi − x)2 + (yi − y)2 . (19)
B, TroundA as the time duration from device A sending the
polling message to receiving the response message from device Let di denote the measured distance between the unknown tag
B, TroundB as the time from device B sending the response and the ith BS. Then, the difference between the true distance
message to receiving the final message from device A, and and the measured distance can be expressed as ρi = di −di . To

Authorized licensed use limited to: SOUTH CHINA UNIVERSITY OF TECHNOLOGY. Downloaded on May 13,2025 at 07:21:27 UTC from IEEE Xplore. Restrictions apply.
 FENG et al.: KF-BASED INTEGRATION OF IMU AND UWB FOR HIGH-ACCURACY INDOOR POSITIONING AND NAVIGATION 3137

III. M ULTISENSOR DATA F USION A LGORITHM BASED

ON IMU AND UWB
In this section, we first introduce the EKF fusion position-
ing algorithm according to three distance measurements. Then,
we propose the UKF fusion positioning algorithm by only sin-
gle distance and angle measurements. Finally, we develop two
approximate motion models, namely, AUM and AUAM, to
further reduce the jitter of positioning data.
In the SINS, the positioning error may increase rapidly with
Fig. 5. LS position algorithm of three BSs. time due to the cumulative error of accelerometer measure-
ment [24]. On the other hand, UWB can provide centimeter-
level ranging and positioning accuracy. However, the UWB
deal with the ranging noise, we adopt the classic LS method signal is vulnerable to the obstacles and only suitable for LOS
to minimize the value of ni=1 ρi2 . Specifically, each distance conditions. To this end, we consider to improve the position-
determines an equation for the position of an unknown tag ing accuracy of the carrier by using KF to merge the data from
⎧ 2 IMU and UWB sensors.
⎪ d12 = (x1 − x)2 + (y1 − y)2
2 2
⎪ It is known that linear KF can derive the optimal carrier state
⎪
⎨ d = (x2 − x) + (y2 − y)
2 under the condition of the linear Gaussian model, in which
⎪ .. (20)
the noises of IMU and UWB sensors are independent with
⎪
⎪ .
⎩ 2 each other and both obey Gaussian distribution with zero mean
dn = (xn − x)2 + (yn − y)2 .
and variances σ 2 . However, there always exist some degrees
Let all the equations subtract the first equation, then we get of nonlinearity in the actual system, such as the square rela-
tion and trigonometric function relation in the state equation
Aε = b (21) or observation equation. To more accurately approximate the
nonlinear systems, various filtering algorithms are adopted to
where deal with the nonlinear factors. In [25], a dual filter integra-
⎡ ⎤ tion for microelectromechanical system (MEMS) sensors and
x2 − x1 y2 − y1
⎢x3 − x1 y3 − y1 ⎥   WiFi fingerprints are proposed, where the EKF algorithm is
⎢ ⎥ x
A=⎢ . .. ⎥ , ε = used to obtain the smooth constrained fingerprints.
⎣ .. . ⎦ y
xn − x1 yn − y1
⎡ 2  ⎤ A. EKF Algorithm Based on Distance Measurements
x2 + y22 − d22 − x12 + y21 − d12 
In the 2-D plane, assume that the tag is moving in a
1⎢
2 ⎥
⎢x3 + y3 − d3 − x1 + y1 − d1 ⎥
2 2 2 2 2
b= ⎢ .. ⎥. (22) straight line with uniform acceleration, define X(k) as the state
2⎣ . ⎦ vector at time k, contains position, x(k), velocity, v(k), and

xn2 + y2n − dn2 − x12 + y21 − d12 acceleration, a(k), which can be expressed as follows:

Thus, the LS solution of ε is X(k) = xx (k) xy (k) vx (k) vy (k) ax (k) ay (k) .
T

 −1 (26)
ε = AT A AT b. (23)
When the sample time period is T, Tω(k) denotes the pro-
In the LS method, each ranging value has adopted the same
cess noise of acceleration, (T 2 /2)ω(k) denotes the process
weight. Obviously, in the process of ranging, when the carrier
noise of velocity and (T 3 /6)ω(k) denotes the process noise
is closer to the BS, the ranging error is smaller. Therefore,
of position due to the double integral of acceleration, accord-
we choose a larger weight for the smaller ranging value, the
ing to the equation of uniform acceleration motion at time
positioning accuracy will be further improved. To this end,
k + 1. Thus, the state equation can be expressed as follows:
we propose the weighted LSs (WLSs) algorithm to solve the
problem, in which the weighting coefficient η is expressed by ⎧ T3
⎪ xx (k + 1) = xx (k) + vx (k)T + 2 ax (k)T + 63 ωx (k)
⎪ 1 2
the reciprocal of the ranging value d as ⎪
⎪
⎪
⎪ xy (k + 1) = xy (k) + vy (k)T + 12 ay (k)T 2 + T6 ωy (k)
⎡ 1 ⎤ ⎪
⎨ 2
d2 0 0 0 vx (k + 1) = vx (k) + ax (k)T + T2 ωx (k) (27)
⎢0 1
0⎥ ⎪ T2
⎢ 0 ⎥ ⎪
⎪ v (k + 1) = v (k) + a (k)T + ω (k)
η=⎢
d3 ⎥. (24) ⎪
⎪
y y y 2 y
⎢ .. ⎥ ⎪
⎣0 0 . 0⎦ ⎩ ax (k + 1) = ax (k) + Tωx (k)
⎪
1 ay (k + 1) = ay (k) + Tωy (k).
0 0 0 dn
Then, the state equation in matrix form is expressed as
Then, the WLS solution of 
ε is expressed as follows:
 −1

ε = AT ηA AT ηb. (25) X(k + 1) = FX(k) + GW(k) (28)

Authorized licensed use limited to: SOUTH CHINA UNIVERSITY OF TECHNOLOGY. Downloaded on May 13,2025 at 07:21:27 UTC from IEEE Xplore. Restrictions apply.
 3138 IEEE INTERNET OF THINGS JOURNAL, VOL. 7, NO. 4, APRIL 2020

where F denotes the state transition matrix, G denotes the Algorithm 1 EKF Algorithm
 T
noise driving matrix, W(k) = ωx (k) ωy (k) denotes pro- Initialize: State mean U(0) = E[X(0)], state covariance matrix
cess noise vector with zero mean and covariance matrix P(0) = var[X(0)], In is n × n unit matrix.
Q = diag(σax2 , σ 2 ) at time k 1: Predict state.
ay 
X (k|k − 1) = FX(k − 1)
⎡ 2 ⎤ ⎡ T3 ⎤
1 0 T 0 T2 0 6 0 2: Predict state covariance matrix.
⎢0 1 0 T T2 ⎥
⎢ T3 ⎥
⎢ 0 ⎥ ⎢0 6 ⎥ P(k|k − 1) = FP(k − 1|k − 1)F T + GQGT
⎢0 0 1 0 2 ⎥ ⎢ ⎥
⎢ T2 ⎥
F=⎢ ⎢0 0 0 1
T 0 ⎥, G = ⎢ 2
⎥
0 ⎥. 3: Calculate Kalman filter gain.
⎢ 0 T ⎥ ⎢ T 2⎥
⎣0 0 0 0 ⎢ 0 2 ⎥  −1
1 0⎦ ⎣T 0⎦ K = P(k|k − 1)H T (k) H(k)P(k|k − 1)H T (k) + R
0 0 0 0 0 1 0 T
4: Update state.
(29)
X(k|k) = 
X (k|k − 1) + K[Z(k) − Z(k|k − 1)]
Let Z(k) denote the observation vector, containing the true
distance di (k) with the observation noise vi (k) at time k. Then, 5: Update state covariance matrix.
the observation equation can be expressed as follows: P(k|k) = [In − KH(k)]P(k|k − 1)
⎡ ⎤
d1 (k) + v1 (k)
⎢d2 (k) + v2 (k)⎥
⎢ ⎥
Z(k) = ⎢ .. ⎥ = H(k)X(k) + V(k) (30)
⎣ . ⎦
dn (k) + vn (k)
where
 H(k) represents the observation matrix, and V(k) =
T
v1 (k) v2 (k) . . . vn (k) represents the observation
noise vector with zero mean and covariance matrix R =
2 , σ 2 , . . . , σ 2 ) at time k. In addition, the detailed
diag(σd1 d2 dn
distance equation is shown as follows:
⎡  2 ⎤
⎡ ⎤ (xx (k) − x1 ) + xy (k) − y1
2
d1 (k) ⎢
⎢d2 (k)⎥ ⎢  2 ⎥
⎥
Fig. 6. DPA of single BS.
⎢ ⎥ ⎢ ⎢ (x (k) − x )2
+ x (k) − y ⎥
⎢ .. ⎥ = ⎢
x 2 y 2 ⎥. (31)
⎣ . ⎦ ⎢ .
.. ⎥
⎥
⎣ ⎦ B. UKF Algorithm Based on Distance and Angle
dn (k)  2
(xx (k) − xn )2 + xy (k) − yn Measurements
Because of the EKF algorithm only uses the first-order
Since (31) is nonlinear, it needs to be linearized and the EKF
approximation in the Taylor expansion, it inevitably intro-
algorithm can be adopted. At each time step, by taking the
duces the linearization error. In the UKF algorithm, it does
first-order Taylor expansion, the Jacobian matrix, H(k), is
not use the linearization process of the nonlinear function.
obtained as
⎡ ∂d (k) ∂d1 (k)
⎤ Particularly, for the one-step prediction, the mean and vari-
1
∂xx (k) ∂xy (k) 0 0 0 0 ance of the equation undergoing nonlinear transformation are
⎢ ∂d2 (k) ∂d2 (k) ⎥
⎢ 0 0 0 0⎥ captured by the unscented transformation (UT) [26]. The UKF
⎢ ∂xx (k) ∂xy (k) ⎥
H(k)  ⎢ . .. .. .. .. .. ⎥ (32) algorithm approximates the probability density distribution of
⎢ .. . . . . .⎥
⎣ ⎦ the nonlinear function. It obtains the posterior probability den-
∂dn (k) ∂dn (k)
∂xx (k) ∂xy (k) 0 0 0 0 sity of the state through a set of deterministic samples rather
than approximating the nonlinear function by the derivative of
where the Jacobian matrix. In this way, it can effectively overcome
⎧ ∂di (k) xx (k)−xi
⎨ ∂xx (k) =
⎪
(xx (k)−xi )2 +(xy (k)−yi )
2
the limitations of low accuracy and poor stability of the EKF
algorithm.
∂di (k) xy (k)−yi
⎪
⎩ ∂xy (k) = . In this article, to further reduce the deployment cost of the
(xx (k)−xi )2 +(xy (k)−yi )
2
BSs, we propose an efficient positioning algorithm based on
The main idea of the EKF algorithm is to linearize the non- single BS with known position (x0 , y0 ). In particular, UWB
linear state or observation equations by Taylor’s expansion detects distance d and angle ϕ between the tag and the BS to
and retains its first-order approximation term. However, this establish the observation equation, where the angle is obtained
procedure inevitably introduces the linearization errors. If the by applying the arctangent operation to the position difference.
linearization assumption is not true, the performance of the IMU detects the acceleration to establish the state equation.
EKF algorithm will degrade and diverge. In addition, it is not The single BS positioning process is shown in Fig. 6.
easy to calculate the Jacobian matrix, which increases the com- Assuming that the tag is moving in the 2-D plane with
putational complexity of the algorithm. The detailed process uniform acceleration, the state vector containing position,
of the EKF algorithm is shown in Algorithm 1. velocity, and acceleration is the same as that in the EKF

Authorized licensed use limited to: SOUTH CHINA UNIVERSITY OF TECHNOLOGY. Downloaded on May 13,2025 at 07:21:27 UTC from IEEE Xplore. Restrictions apply.
 FENG et al.: KF-BASED INTEGRATION OF IMU AND UWB FOR HIGH-ACCURACY INDOOR POSITIONING AND NAVIGATION 3139

Algorithm 2 UKF Algorithm observation equation can be expressed as follows:

Initialize: State mean Uu (0) = E[Xu (0)], state covariance matrix  
Pu (0) = var[Xu (0)]. d(k) + vd (k)
Zu (k) = = h(X(k)) + Vu (k) (34)
1: Calculate sigma points. ϕ(k) + vϕ (k)
⎡ ⎤T
X√u (k|k)
(i) where h(X(k)) is a nonlinear
 function of the observation
Xu (k|k) = ⎣ X u (k|k) + √(n + λ)Pu (k|k) ⎦ , T
X u (k|k) − (n + λ)Pu (k|k) equation, and Vu (k) = vd (k) vϕ (k) is the observation
noise vector with zero mean and covariance matrix Ru =
2: Predict sigma points. diag(σd2 , σϕ2 )
(i) (i)
Xu (k + 1|k) = FXu (k|k), i = 1, . . . , 2n + 1, ⎡
   2 ⎤
3: Calculate state mean and state covariance matrix. d(k) (xx (k) − x0 )2 + xy (k) − y0
=⎣   ⎦. (35)
2n ϕ(k) arctan y
x (k)−y0
(i) (i) xx (k)−x0
X u (k + 1|k) = ωm Xu (k + 1|k),
i=0
2n   Unlike EKF, UKF does not need to compute the Jacobian
(i) (i)
Pu (k + 1|k) = ωc  Xu (k + 1|k) − X u (k + 1|k)
i=0
matrix at each time and adopts UT to make the statistics
 T keep consistent for the random variables undergoing nonlin-
 (i)
Xu (k + 1|k) − X u (k + 1|k) + Qu ,
ear transformation [27]. For the Gaussian distribution, through
4: Update sigma points. carefully selecting the sample points, UT can capture the mean
⎡ ⎤T and covariance accurately to the third order. The selection
(k + 1|k)
X u√
(i) of sample points is based on the prior mean and covariance
Xu (k + 1|k) = ⎣X u (k + 1|k) + √(n + λ)Pu (k + 1|k)⎦ ,
X u (k + 1|k) − (n + λ)Pu (k + 1|k) matrix. Define a nonlinear transformation y = f (x), where the
random variable x is n-dimensional with mean x and covari-
5: Predict observation. ance matrix P. Then, the statistical characteristics of y can be
 
(i) (i)
Zu (k + 1|k) = h Xu (k + 1|k) , i = 1, . . . , 2n + 1, calculated by forming a matrix χ of 2n + 1 sigma points χ (i)
with the corresponding weight ω as follows:
6: Calculate observation mean and observation covariance matrix.
⎧ (0)
2n
(i) (i) ⎨ χ = x √ 
Z u (k + 1|k) = ωm 
Zu (k + 1|k),
χ (i) = x + √(n + λ)Pi , i = 1, . . . , n (36)
i=0
  ⎩ (i)
2n
(i) (i) χ = x − (n + λ)P i , i = n + 1, . . . , 2n
Pzk zk = ωc Zu (k + 1|k) − Z u (k + 1|k)
 i=0 T √
 (i)
Zu (k + 1|k) − Z u (k + 1|k) + Ru ,
where ( P)i represents the ith column of the square root of a
2n   matrix. The corresponding weights of these sample points are
(i) (i)
Pxk zk = ωc  Xu (k + 1|k) − X u (k + 1|k) calculated as follows:
 i=0 T ⎧ (0)
(i) ⎪ λ
Zu (k + 1|k) − Z u (k + 1|k) , ⎨ ωm = n+λ  
(0) λ
ωc = n+λ + 1 − α2 + β (37)
7: Calculate Kalman gain matrix. ⎪
⎩ (i) (i)
ωm = ωc = 2(n+λ) , i = 1, 2, . . . , 2n
1
K = Pxk zk P−1
zk zk ,
(i)
8: Update state and state covariance matrix. where ωm is the weight of mean and ωc(i) is the weight of
 covariance of the sigma points, the superscript i is the index
Xu (k + 1|k + 1) = X u (k + 1|k) + K Zu (k + 1) − Z u (k + 1|k) ,
Pu (k + 1|k + 1) = Pu (k + 1|k) − KPzk zk K T . of sample points, and λ  α 2 (n+κ)−n is a scaling parameter
to reduce the total prediction error. In practice, α is usually set
to a small positive value to keep the mean of the sigma points
around x, κ is set to zero to ensure that matrix (n + λ)P is
algorithm as shown in (26) and (27). Then, the state equation a semipositive-definite matrix, and β is a nonnegative weight
in matrix form can be expressed as follows: coefficient to incorporate the prior distribution of x. According
(i)
to (36) and (37), we can get a set of sigma points Xu (k|k),
Xu (k + 1) = FXu (k) + Wu (k) (33) and then further predict them to get the state mean and state
covariance matrix as X u (k + 1|k) and Pu (k + 1|k), respectively.
where F denotes the state transition matrix and is the same Similarly, we can get the observation mean and observation
as that in (28) and (29) in the EKF algorithm, Wu (k) T
= covariance matrix as Z u (k + 1|k) and Pzk zk , respectively, and
ωxx (k) ωxy (k) ωvx (k) ωvy (k) ωax (k) ωay (k) also the cross-covariance matrix as Pxk zk between the state
denotes the process noise vector with zero mean and vector and the observation vector. Then, we can calculate the
covariance matrix Qu = diag(σx2x , σx2y , σv2x , σv2y , σa2x , σa2y ). Kalman gain, K, and finally obtain the state and state covari-
Let Zu (k) denote the observation vector. It contains the true ance matrix updates as Xu (k + 1|k + 1) and Pu (k + 1|k + 1),
distance d(k), true angle ϕ(k) with observation distance noise respectively. The detailed process of UKF algorithm is shown
vd (k), and observation angle noise vϕ (k) at time k. Then, the in Algorithm 2.

Authorized licensed use limited to: SOUTH CHINA UNIVERSITY OF TECHNOLOGY. Downloaded on May 13,2025 at 07:21:27 UTC from IEEE Xplore. Restrictions apply.
 3140 IEEE INTERNET OF THINGS JOURNAL, VOL. 7, NO. 4, APRIL 2020

C. Dilution of Precision Models Algorithm 3 AUAM Filtering Algorithm

Initialize:UWB sample time period t=0.02s (the rate of updates by
The positioning accuracy is not only affected by the distance
the UWB is 50Hz), IMU sample time period T=0.01s (the rate of
measurement accuracy but also by the geometric distribution updates by the IMU is 100Hz)
of the BS. The DOP represents the geometric distribution 1: Measure the distance.
of tag and BSs, which plays an important role in the posi- 
tioning accuracy. The DOP mainly includes the following d = d1 d2 . . . dn
parameters: geometric DOP (GDOP), position DOP (PDOP), 2: Calculate the LS position.
horizontal DOP (HDOP), vertical DOP (VDOP), and time  −1
DOP (TDOP). This article focuses on the positioning of the ε(k) = AT A AT b
tag in the 2-D plane. Under the same distance measurement 3: Calculate the differential velocity.
accuracy, the smaller the HDOP, the higher the positioning
ε(k) − ε(k − 1)
accuracy. Therefore, it is important to arrange the BS distri- v(k) =
bution reasonably to reduce the HDOP value. In detail, the t
HDOP can be defined as 4: Convert the acceleration coordinate system.

σx2 + σy2 an = anx any = Cbn ab
HDOP = (38)
σd 5: Update the position state.

where σd is the standard deviation of distance measurement ε(k + 1) = ε(k) + v(k)t + 12 an T 2

error, σx2 and σy2 are the variances of x-axis and y-axis, respec-
tively. According to (19), the derivatives at the approximate
position (x , y ) can be given as follows:
⎡ ⎤ ⎡ ∂d1 ∂d1 ⎤ D. AUAM and AUM Approximate Motion Models
d1 − d1 ∂x ∂y
⎢d2 − d ⎥ ⎢ ⎢
∂d2 ∂d2 ⎥  Since the motion of the tag is usually random, it is difficult
⎢ 2⎥ ⎢ ∂x ∂y ⎥⎥ x − x to establish an accurate motion model. Therefore, we pro-
⎢ .. ⎥ ⎢ . = .. ⎥ y − y (39)
⎣ . ⎦ ⎣ .. . ⎦ pose two approximate motion models in a short time period:
dn − dn ∂dn ∂dn 1) one is the AUAM model and 2) the other is the AUM model.
∂x ∂y
Specifically, the AUAM model combines the position and the
where differential velocity information from the UWB as well as the
⎧ ∂d xi −x
⎨ ∂xi = √ acceleration information from the IMU. The whole procedure
(xi −x )2 +(yi −y )2
yi −y of the AUAM filtering algorithm is shown in Algorithm 3.
⎩ ∂di = √ .
∂y (xi −x )2 +(yi −y )2 Without considering the acceleration from the IMU sensor, the
When n ≥ 3, the LS method can be used to solve (39) and AUAM model will degenerate into the AUM model which only
we get integrates the position and differential velocity information.
 −1
C = BT B BT D (40) IV. S IMULATION AND E XPERIMENTAL R ESULTS
where In this section, we evaluate the performance of the proposed
⎡ ⎤ ⎡ ∂d1 ∂d1 ⎤
d1 − d1 ∂x ∂y
algorithm by simulation and experiment. First, we compare
⎢d2 − d ⎥ ⎢ ∂d2 ∂d2 ⎥   the EKF fusion positioning algorithm with the LS and WLS
⎢ 2⎥ ⎢ ∂x ∂y ⎥ x − x
D = ⎢ . ⎥, B = ⎢ ⎢ ..
⎥
.. ⎥, C = y − y . positioning algorithms. Then, we evaluate the performances
⎣ .. ⎦ ⎣ . . ⎦ of the UKF fusion algorithm and the DPA algorithm based
dn − dn ∂dn ∂dn on a single BS. In addition, we compare the performances of
∂x ∂y
the proposed EKF and UKF algorithms. Finally, we carry out
The covariance matrix of distance error R is the same as that
experiments in the laboratory to verify the actual effect of the
in the EKF algorithm. The positioning error covariance matrix
proposed algorithms.
Qc caused by distance error can be expressed as follows:
 −1
Qc = BT R−1 B BT R−1 σu2 A. EKF Simulation Results
 −1 In the simulations, it is assumed that the tag moves in the xy
= BT B σu2 plane with the initial position at (0, 0), the initial horizontal and
= Eσu2 (41) vertical velocities are both at 0.15 m/s, the initial horizontal
where E = (BT B)−1is a symmetric matrix representing the and vertical accelerations are both at 0.002 m/s2 , the sample
influence of distance error on positioning error, σu2 is the vari- time period T = 1 s, and the total running time N = 50 s.
ance of user equivalent ranging error (UERE), and the smaller Fig. 7 demonstrates the performances of the LS and the
UERE means better positioning accuracy [28]. Let proposed EKF algorithms in terms of the positioning trajectory
  and error, where Q and R denote the process and obser-
E E12
E = 11 . (42) vation noise variances, respectively. As shown in Fig. 7(a),
E21 E22 the LS trajectory has oscillation and the larger the obser-
√
Finally, we can get the HDOP = E11 + E22 . vation noise variance is, the larger the oscillation appears.

Authorized licensed use limited to: SOUTH CHINA UNIVERSITY OF TECHNOLOGY. Downloaded on May 13,2025 at 07:21:27 UTC from IEEE Xplore. Restrictions apply.
 FENG et al.: KF-BASED INTEGRATION OF IMU AND UWB FOR HIGH-ACCURACY INDOOR POSITIONING AND NAVIGATION 3141

(a) (b) (c)

(d) (e) (f)

Fig. 7. Performances of LS and EKF algorithms. (a) LS trajectory (Q = 10−8 ). (b) EKF trajectory (Q = 10−8 ). (c) Position error CDF (Q = 10−8 ).
(d) LS trajectory (R = 0.01). (e) EKF trajectory (R = 0.01). (f) Position error CDF (R = 0.01).

(a) (b) (c)

Fig. 8. Performances of LS and WLS algorithms. (a) LS and WLS position. (b) Position error. (c) Position error CDF.

It indicates that the observation noise has great influence of the LS algorithm is over 31 cm and the average position
on the measurements. However, after using the proposed error is about 14 cm. For the proposed EKF algorithm, the
EKF algorithm, the trajectory accuracy can be significantly maximum position error is about only 20 cm and the average
improved and the results are much closer to the true trajectory position error is less than 9 cm. The positioning accuracy has
as shown in Fig. 7(b). Moreover, it can be seen that when been improved by about 50%. When R = 0.01 and Q = 10−4 ,
the observation noise variance increases, the EKF trajectory the maximum position error of the LS algorithm is over 57 cm
keeps consistent. As shown in Fig. 7(d), when the process and the average position error is about 18 cm. However, the
noise variance becomes large, the carrier’s motion will be maximum position error is about only 22 cm and the aver-
curvilinear. However, after using the proposed EKF algorithm, age position error is reduced to less than 10 cm when the
the trajectory accuracy can be significantly improved and the proposed EKF fusion algorithm is used. In addition, to eval-
results are much closer to the true trajectory as shown in uate the performance of the WLS algorithm, we randomly
Fig. 7(e). In Fig. 7(c) and (f), the cumulative distribution func- generate four BSs with known positions and two hundred ref-
tions (CDF) of positioning errors are illustrated. Specifically, erence points on a 10-m by 10-m area. We also assume that the
when R = 0.01 and Q = 10−8 , the maximum position error observation distance noise is white Gaussian noise with zero

Authorized licensed use limited to: SOUTH CHINA UNIVERSITY OF TECHNOLOGY. Downloaded on May 13,2025 at 07:21:27 UTC from IEEE Xplore. Restrictions apply.
 3142 IEEE INTERNET OF THINGS JOURNAL, VOL. 7, NO. 4, APRIL 2020

(a) (b) (c)

(d) (e) (f)

Fig. 9. Performances of DPA and UKF algorithms. (a) DPA trajectory (Q = 10−8 ). (b) UKF trajectory (Q = 10−8 ). (c) Position error CDF (Q = 10−8 ).
(d) DPA trajectory (Q = 10−4 ). (e) UKF trajectory (Q = 10−4 ). (f) Position error CDF (Q = 10−4 ).

(a) (b) (c)

Fig. 10. Performances of EKF and UKF algorithms. (a) Position trajectory. (b) Position error. (c) Position error CDF.

mean and variance 0.01 m. The simulation results in Fig. 8 Fig. 9 compares the performances of the DPA without fil-
show that the positioning accuracy of the WLS algorithm is tering and UKF algorithms. As shown in Fig. 9(a), at the
improved compared with the LS algorithm. initial stage, the DPA trajectory can continuously track the
true trajectory. However, after several iterations, the oscillation
appears. This is mainly due to the accumulation of observa-
B. UKF Simulation Results tion angle noise. On the other hand, as shown in Fig. 9(b),
In the UKF simulations, it is also assumed that the tag the UKF trajectory can still track the true trajectory, which
moves in the xy plane with the initial position at (0, 0), the indicates that the proposed UKF algorithm can improve the
initial horizontal velocity and the vertical velocity are both at positioning accuracy based on single observation distance and
0.15 m/s, the initial horizontal acceleration and the vertical angle. Moreover, it is shown that when the observation noise
acceleration are both at 0.002 m/s2 , the sample time period variance increases, the UKF trajectory is still consistent with
T = 1 s, and the total running time N = 50 s. The correlation the true trajectory. Comparing Fig. 9(a) with Fig. 9(d) and
coefficients α = 0.01, β = 2, κ = 0, and state dimension Fig. 9(b) with Fig. 9(e), as the process noise becomes larger,
n = 6, the observation dimension m = 2 in UT. In the follow- the carrier’s motion is to the curvilinear motion. In addition,
ing, Rd and Ra denote the observation distance noise variance the CDF of the positioning error is shown in Fig. 9(c) and (f).
and angle noise variance, respectively. Specifically, when Rd = 0.1, Ra = 0.001, and Q = 10−8 , the

Authorized licensed use limited to: SOUTH CHINA UNIVERSITY OF TECHNOLOGY. Downloaded on May 13,2025 at 07:21:27 UTC from IEEE Xplore. Restrictions apply.
 FENG et al.: KF-BASED INTEGRATION OF IMU AND UWB FOR HIGH-ACCURACY INDOOR POSITIONING AND NAVIGATION 3143

(a) (b) (c)

Fig. 11. DOP. (a) LS trajectory. (b) DOP results. (c) Position error CDF.

Fig. 12. Hardware architecture of the positioning module.

maximum position error of DPA is over 78 cm and the average

position error is about 33 cm. However, for the proposed UKF
algorithm, the maximum position error is about only 18 cm
and the average position error is less than 7 cm. The position-
ing accuracy has been improved by about 78%. In order to
better evaluate the performances of two positioning algorithms
and two filtering algorithms, we compare their performances
with the same assumption R = 0.01, Rd = 0.01, Ra = 0.01,
and Q = 10−8 .
Fig. 10 shows the results of different positioning algorithms
Fig. 13. Indoor positioning and navigation experiment. (a) Experimental
based on single BS and three BSs. The average positioning environment. (b) Ground grid.
error of the LS algorithm based on three BSs is 15 cm, while
that of the DPA algorithm based on single BS is 53 cm. The
average positioning error of the EKF fusion algorithm based BSs can reduce the HDOP value. However, the effect becomes
on three BSs is 6 cm, and that of a single BS is 18 cm. The less and less and it will increase the system complexity with
average positioning error of the UKF fusion algorithm based the increase of BSs. The CDF result of the system positioning
on three BSs is 5 cm, and that of single BS is 13 cm. In error is shown in Fig. 11(c).
summary, when the positioning accuracy requirement is high,
multiple BS positioning algorithms can be selected, otherwise, D. Experiments and Results
single BS positioning algorithm can be selected to reduce the
complexity and cost of deploying BSs. The hardware used in this experiment is shown in Fig. 12,
which is with the “baseplate and module” structure. The base-
plate adopts STM32F103T8U6 microcomputer as the main
C. DOP Simulation Results control chip. The peripheral circuit includes power, LED indi-
In Fig. 11, we evaluate the influence of the distribution cator, lithium battery management, LIS3DH, accelerometer,
of BS. In the triangle model, the three BSs are located at etc. The modules are configured with DWM1000 of Decawave
(−5, −5), (5, −5), and (5, 5); in square model, the four BSs and MAX2000 of YCHIOT. Note that the hardware can be
are located at (−5, −4.5), (5, −4.5), (5, 4.5), and (−5, 4.5); configured as BS or tag through the USB instruction.
in X model, the five BSs are located at (−5, −5.5), (5, −5.5), The experimental environment is shown in Fig. 13(a), where
(5, 5.5), (−5, 5.5), and (0, 0) as shown in Fig. 11(a). In each smaller grid on the ground is a square with a side length
Fig. 11(b), it can be seen that the HDOP value is relatively of 0.5 m. To accurately measure the true position of the tag,
small in the center area and gradually increases outwards the all the brackets and the BSs are deployed in the same way.
edge of BSs. It is also shown that increasing the number of Then, the true position of the tag is determined by the vertical

Authorized licensed use limited to: SOUTH CHINA UNIVERSITY OF TECHNOLOGY. Downloaded on May 13,2025 at 07:21:27 UTC from IEEE Xplore. Restrictions apply.
 3144 IEEE INTERNET OF THINGS JOURNAL, VOL. 7, NO. 4, APRIL 2020

(a)
Fig. 15. Experiments of AUAM and AUM models.

(b)

Fig. 14. Experiment of ranging calibration. (a) Calibration experiment. (a)

(b) Ranging error CDF.

projection of the center points of small red, green, blue, and

yellow squares from the top of the module antenna to the
ground grid, as shown in Fig. 13(b).
Since the actual ranging measurement is affected by the
environment, including latitude, longitude, air quality, obstacle
height, and other factors, we must first calibrate the module.
Particularly, the calibration coefficients can be obtained by
measuring the distance between BS and tag. There are many
fitting formulas to fit the data and the simplest linear equation
is used in the experiments. Obviously, the ranging accuracy (b)
will become better after calibrating as shown in Fig. 14.
In the first experiment scenario, the tag starts from (1, 1)
and moves two rounds along a square with a side length of
4 m as shown in Fig. 15. The experimental results show that
the LS trajectory oscillates relative to the true trajectory due
to the influence of noise. The AUAM filtering trajectory has
a good tracking effect in the initial stage, but the error caused
by the acceleration becomes larger with the passage of time.
The AUM filtering trajectory keeps continuous tracking. It can
also be seen that the positioning trajectories of the AUAM and
AUM algorithms are smoother than the LS algorithm.
In the second experiment scenario, the tag moves along a (c)
square with a side length of 6 m as shown in Fig. 16. In
Fig. 16. Experiments of spiral and circuitous. (a) Experiments of spiral.
the spiral experiment, the results show that the root mean (b) Experiments of circuitous. (c) Position error CDF.
square error (RMSE) of x-axis is 7.79 cm, and that of y-axis is
6.32 cm. In the circuitous experiment, the results show that the
RMSE of x-axis is 5.50 cm, and that of y-axis is 5.03 cm. It is
because there are fewer sample points on the edge of the BSs expanding the deployment of the BSs and ensuring that the
in the circuitous experimental. Thus, when the tag is closer to tag always moves in the center of the coverage can effectively
the center, the positioning accuracy is better. It implies that improve the positioning accuracy.

Authorized licensed use limited to: SOUTH CHINA UNIVERSITY OF TECHNOLOGY. Downloaded on May 13,2025 at 07:21:27 UTC from IEEE Xplore. Restrictions apply.
 FENG et al.: KF-BASED INTEGRATION OF IMU AND UWB FOR HIGH-ACCURACY INDOOR POSITIONING AND NAVIGATION 3145

V. C ONCLUSION [14] Y. Zhuang, Y. Cao, N. El-Sheimy, and J. Yang, “Guest editorial: Special
issue on toward positioning, navigation, and location-based services
In this article, a fusion positioning system based on IMU (PNLBS) for Internet of Things,” IEEE Internet Things J., vol. 5, no. 6,
and UWB was studied. The data obtained by IMU were used pp. 4613–4615, Dec. 2018.
for the state equation of KF while the data obtained by UWB [15] L. Yao, Y.-W. A. Wu, L. Yao, and Z. Liao, “An integrated IMU and
UWB sensor based indoor positioning system,” in Proc. IEEE Int. Conf.
was used for the observation equation of KF. First, the EKF Indoor Position. Indoor Navig. (IPIN), Sep. 2017, pp. 1–8.
algorithm based on multiple observation BSs was introduced to [16] Y. Zhong, T. Liu, B. Li, L. Yang, and L. Lou, “Integration of UWB
improve the positioning accuracy. Then, an efficient UKF algo- and IMU for precise and continuous indoor positioning,” in Proc. IEEE
Ubiquitous Position. Indoor Navig. Location Based Services (UPINLBS),
rithm based on single observation BS was proposed to reduce Mar. 2018, pp. 1–5.
the deployment complexity of BS and lower the cost. In addi- [17] P. K. Yoon, S. Zihajehzadeh, B.-S. Kang, and E. J. Park, “Robust
tion, to reduce the jitter of positioning data, two approximate biomechanical model-based 3-D indoor localization and tracking method
using UWB and IMU,” IEEE Sensors J., vol. 17, no. 4, pp. 1084–1096,
motion models, namely, AUM and AUAM, were proposed to Feb. 2017.
make the positioning results smoother and more stable. In [18] S. Zihajehzadeh and E. J. Park, “A novel biomechanical model-aided
addition, the experimental results showed that the proposed IMU/UWB fusion for magnetometer-free lower body motion capture,”
IEEE Trans. Syst., Man, Cybern., Syst., vol. 47, no. 6, pp. 927–938,
EKF algorithm improves the positioning accuracy by 54% Jun. 2017.
compared to the LS algorithm and the proposed UKF algo- [19] S. Zihajehzadeh, P. K. Yoon, B. S. Kang, and E. J. Park, “UWB-
rithm improves the positioning accuracy by 50% compared to aided inertial motion capture for lower body 3-D dynamic activity
and trajectory tracking,” IEEE Trans. Instrum. Meas., vol. 64, no. 12,
the DPA algorithm. Note that in our system, we have only pp. 3577–3587, Dec. 2015.
considered the Gaussian noise distribution. Thus, for reducing [20] T. Seel and S. Ruppin, “Eliminating the effect of magnetic disturbances
the complexity, we have only studied the EKF and UKF meth- on the inclination estimates of inertial sensors,” in Proc. 20th IFAC World
Congr., vol. 50, Jul. 2017, pp. 8798–8803.
ods. In the future work, we will consider a generalized noise [21] L. Chang, J. Li, and S. Chen, “Initial alignment by attitude estimation
distribution, and focus on the application and the optimization for strapdown inertial navigation systems,” IEEE Trans. Instrum. Meas.,
of other filter methods, e.g., particle filter. vol. 64, no. 3, pp. 784–994, Mar. 2015.
[22] G. Wang, S. Qian, Q. Lv, H. Wei, H. Lin, and B. Liang, “UWB and
IMU system fusion for indoor navigation,” in Proc. IEEE 37th Chin.
R EFERENCES Control Conf. (CCC), Jul. 2018, pp. 4946–4950.
[23] A. Wang and Y. Song, “Improved SDS-TWR ranging technology in
[1] Y. Zhuang, Q. Wang, M. Shi, P. Cao, L. Qi, and J. Yang, “Low-power UWB positioning,” in Proc. IEEE Int. Conf. Sensor Netw. Signal
centimeter-level localization for indoor mobile robots based on ensemble Process. (SNSP), Oct. 2018, pp. 222–225.
Kalman smoother using received signal strength,” IEEE Internet Things [24] A. Noureldin, T. B. Karamat, M. D. Eberts, and A. El-Shafie,
J., vol. 6, no. 4, pp. 6513–6522, Aug. 2019. “Performance enhancement of MEMS-based INS/GPS integration for
[2] Y. Zhuang and N. El-Sheimy, “Tightly-coupled integration of WiFi and low-cost navigation applications,” IEEE Trans. Veh. Technol., vol. 58,
MEMS sensors on handheld devices for indoor pedestrian navigation,” no. 3, pp. 1077–1096, Mar. 2009.
IEEE Sensors J., vol. 16, no. 1, pp. 224–234, Jan. 2016. [25] Y. Zhuang, Y. Li, L. Qi, H. Lan, J. Yang, and N. El-Sheimy, “A two-
[3] Y. Li, Y. Zhuang, P. Zhang, H. Lan, X. Niu, and N. El-Sheimy, “An filter integration of MEMS sensors and WiFi fingerprinting for indoor
improved inertial/WiFi/magnetic fusion structure for indoor navigation,” positioning,” IEEE Sensors J., vol. 16, no. 13, pp. 5125–5126, Jul. 2016.
Inf. Fusion, vol. 34, pp. 101–119, Mar. 2017. [26] S. Wang, Y. Lyu, and W. Ren, “Unscented transformation based
[4] Y. Li, Z. He, Z. Gao, Y. Zhuang, C. Shi, and N. El-Sheimy, “Toward distributed nonlinear state estimation: Algorithm, analysis, and exper-
robust crowdsourcing-based localization: A fingerprinting accuracy indi- iments,” IEEE Trans. Control Syst. Technol., vol. 27, no. 5,
cator enhanced wireless/magnetic/inertial integration approach,” IEEE pp. 2016–2029, Sep. 2019.
Internet Things J., vol. 6, no. 2, pp. 3585–3600, Apr. 2019. [27] E. A. Wan and R. Van Der Merwe, “The unscented Kalman filter for non-
[5] Y. Zhuang, J. Yang, L. Qi, Y. Li, Y. Cao, and N. El-Sheimy, “A per- linear estimation,” in Proc. IEEE Adapt. Syst. Signal Process. Commun.
vasive integration platform of low-cost MEMS sensors and wireless Control Symp. (ASSPCC), Oct. 2000, pp. 153–158.
signals for indoor localization,” IEEE Internet Things J., vol. 5, no. 6, [28] C. Tang et al., “Improvement of orbit determination accuracy for
pp. 4616–4631, Dec. 2018. Beidou navigation satellite system with two-way satellite time frequency
[6] C. Ren, Q. Liu, and T. Fu, “A novel self-calibration method for MIMU,” transfer,” Adv. Space Res., vol. 58, no. 7, pp. 1390–1400, May 2016.
IEEE Sensors J., vol. 15, no. 10, pp. 5416–5422, Oct. 2015.
[7] K. Guo et al., “Ultra-wideband-based localization for quadcopter navi-
gation,” Unmanned Syst., vol. 4, no. 1, pp. 23–24, Feb. 2016.
[8] A. Chandra et al., “Frequency-domain in-vehicle UWB channel
modeling,” IEEE Trans. Veh. Technol., vol. 65, no. 6, pp. 3929–3940,
Jun. 2016. Daquan Feng received the Ph.D. degree in
[9] S. Wang, G. Mao, and J. A. Zhang, “Joint time of arrival esti- information engineering from the National
mation for coherent UWB ranging in multipath environment with Key Laboratory of Science and Technology on
multi user interference,” IEEE Trans. Signal Process., vol. 67, no. 14, Communications, University of Electronic Science
pp. 3743–3755, Jul. 2019. and Technology of China, Chengdu, China, in
[10] A. Alarifi et al., “Ultra wideband indoor positioning technolo- 2015.
gies: Analysis and recent advances,” Sensors, vol. 16, no. 5, p. 707, He was a Research Staff with State Radio
May 2016. Monitoring Center, Beijing, China, and then a
[11] S. Maranò, W. M. Gifford, H. Wymeersch, and M. Z. Win, “NLOS Postdoctoral Research Fellow with the Singapore
identification and mitigation for localization based on UWB experimen- University of Technology and Design, Singapore.
tal data,” IEEE J. Sel. Areas Commun., vol. 28, no. 7, pp. 1026–1035. He was a visiting student with the School of
Sep. 2010. Electrical and Computer Engineering, Georgia Institute of Technology,
[12] Y. Lu, J. Yi, L. He, X. Zhu, and P. Liu, “A hybrid fusion algorithm for Atlanta, GA, USA, from 2011 to 2014. He is currently an Assistant Professor
integrated INS/UWB navigation and its application in vehicle platoon with the Guangdong Province Engineering Laboratory for Digital Creative
formation control,” in Proc. Int. Conf. Comput. Sci. Electron. Commun. Technology and Guangdong Key Laboratory of Intelligent Information
Eng. (CSECE), Feb. 2018, pp. 157–161. Processing, College of Electronics and Information Engineering, Shenzhen
[13] M. Gunia, F. Protze, N. Joram, and F. Ellinger, “Setting up an ultra- University, Shenzhen, China. His research interests include URLLC commu-
wideband positioning system using off-the-shelf components,” in Proc. nications, LTE-U, and massive IoT networks.
IEEE 13th Workshop Position. Navig. Commun. (WPNC), Oct. 2016, Dr. Feng is an Associate Editor of IEEE C OMMUNICATIONS L ETTERS
pp. 1–6. and IEEE ACCESS.

Authorized licensed use limited to: SOUTH CHINA UNIVERSITY OF TECHNOLOGY. Downloaded on May 13,2025 at 07:21:27 UTC from IEEE Xplore. Restrictions apply.
 3146 IEEE INTERNET OF THINGS JOURNAL, VOL. 7, NO. 4, APRIL 2020

Chunqi Wang received the B.S. degree in opto- Yuan Zhuang (Member, IEEE) received the
electronics information science and engineering bachelor’s degree in information engineering and
from the College of Testing and Opto-Electronic the master’s degree in microelectronics and solid-
Engineering, Nanchang Hangkong University, state electronics from Southeast University, Nanjing,
Nanchang, China, in 2017. He is currently pursuing China, in 2008 and 2011, respectively, and the Ph.D.
the master’s degree in information and communi- degree in geomatics engineering from the University
cation engineering with the Guangdong Province of Calgary, Calgary, AB, Canada, in 2015.
Engineering Laboratory for Digital Creative He was an Algorithm Designer with Trusted
Technology and Guangdong Key Laboratory of Positioning Inc., Calgary, and the Lead Scientist with
Intelligent Information Processing, College of Bluvision Inc., Fort Lauderdale, FL, USA. He is cur-
Electronics and Information Engineering, Shenzhen rently a Professor with the State Key Laboratory of
University, Shenzhen, China. Information Engineering in Surveying, Mapping and Remote Sensing, Wuhan
His current research interests include pedestrians and robots positioning University, Wuhan, China. He has coauthored over 50 academic papers and
and navigation. 11 patents. His current research interests include multisensors integration,
real-time location system, personal navigation system, wireless positioning,
Internet of Things, and machine learning for navigation applications.
Prof. Zhuang has received over ten academic awards. He is an Associate
Editor of IEEE ACCESS, the Guest Editor of the IEEE I NTERNET OF T HINGS
J OURNAL and IEEE ACCESS, and a reviewer of over ten IEEE journals.

Xiang-Gen Xia (Fellow, IEEE) received the

B.S. degree in mathematics from Nanjing Normal
University, Nanjing, China, in 1983, the M.S. degree
in mathematics from Nankai University, Tianjin,
China, in 1986, and the Ph.D. degree in electri-
cal engineering from the University of Southern
California, Los Angeles, CA, USA, in 1992.
He was a Senior/Research Staff Member with
Hughes Research Laboratories, Malibu, CA, USA,
from 1995 to 1996. In September 1996, he joined the
Department of Electrical and Computer Engineering,
Chunlong He (Member, IEEE) received the M.S. University of Delaware, Newark, DE, USA, where he is the Charles Black
degree in communication and information sci- Evans Professor. He has authored the book Modulated Coding for Intersymbol
ence from Southwest Jiaotong University, Chengdu, Interference Channels (New York: Marcel Dekker, 2000). His current research
China, in 2010, and the Ph.D. degree from Southeast interests include space-time coding, MIMO and OFDM systems, digital signal
University, Nanjing, China, in 2014. processing, and SAR and ISAR imaging.
From September 2012 to September 2014, he Dr. Xia received the National Science Foundation Faculty Early Career
was a visiting student with the School of Electrical Development (CAREER) Program Award in 1997, the Office of Naval
and Computer Engineering, Georgia Institute of Research Young Investigator Award in 1998, and the Outstanding Overseas
Technology, Atlanta, GA, USA. Since 2015, he has Young Investigator Award from the National Nature Science Foundation of
been with the College of Electronics and Information China in 2001. He was the Technical Program Chair of the Signal Processing
Engineering, Shenzhen University, Shenzhen, China, Symposium, Globecom 2007 in Washington, DC, USA, and the General Co-
where he is currently an Associate Professor. His research interests include Chair of ICASSP 2005 in Philadelphia. He is currently serving/has served
communication and signal processing, green communication systems, channel as an Associate Editor for numerous international journals, including the
estimation algorithms, and limited feedback techniques. IEEE W IRELESS C OMMUNICATIONS L ETTERS, the IEEE T RANSACTIONS
Dr. He is a member of the Institute of Electronics, Information, and ON S IGNAL P ROCESSING , the IEEE T RANSACTIONS ON W IRELESS
Communication Engineering. He is currently an Associate Editor of IEEE C OMMUNICATIONS, the IEEE T RANSACTIONS ON M OBILE C OMPUTING,
ACCESS. and the IEEE T RANSACTIONS ON V EHICULAR T ECHNOLOGY.

Authorized licensed use limited to: SOUTH CHINA UNIVERSITY OF TECHNOLOGY. Downloaded on May 13,2025 at 07:21:27 UTC from IEEE Xplore. Restrictions apply.