Robust TOA-Based UAS Navigation under Model Mismatch in GNSS-Denied Harsh Environments

: Global Navigation Satellite Systems (GNSS) is the technology of choice for outdoor positioning purposes but has many limitations when used in safety-critical applications such Intelligent Transportation Systems (ITS) and Unmanned Autonomous Systems (UAS). Namely, its performance clearly degrades in harsh propagation conditions and is not reliable due to possible attacks or interference. Moreover, GNSS signals may not be available in the so-called GNSS-denied environments, such as deep urban canyons or indoors, and standard GNSS architectures do not provide the precision needed in ITS. Among the different alternatives, cellular signals (LTE/5G) may provide coverage in constrained urban environments and Ultra-Wideband (UWB) ranging is a promising solution to achieve high positioning accuracy. The key points impacting any time-of-arrival (TOA)-based navigation system are (i) the transmitters’ geometry, (ii) a perfectly known transmitters’ position, and (iii) the environment. In this contribution, we analyze the performance loss of alternative TOA-based navigation systems in real-life applications where we may have both transmitters’ position mismatch, harsh propagation environments, and GNSS-denied conditions. In addition, we propose new robust ﬁltering methods able to cope with both effects up to a certain extent. Illustrative results in realistic scenarios are provided to support the discussion and show the performance improvement brought by the new methodologies with respect to the state-of-the-art.


Introduction
Reliable and precise position, navigation, and timing (PNT) information is fundamental in safety-critical applications such as Intelligent Transportation Systems (ITS), automated aircraft landing, or Unmanned Autonomous Systems (UAS), i.e., unmanned autonomous ground/air vehicles (robots/drones). In addition, in the context of ITS, this is not only of capital importance for the autonomously navigating system, but also for vulnerable road users, such as cyclists and pedestrians, and systems collaterally using this information, i.e., in traffic control and for emergency services. Even if the main source of positioning information are still Global Navigation Satellite Systems (GNSS), they lack of reliability and accuracy in constrained environments such as highly populated areas, which limit their adoption as standalone PNT systems. For instance, GNSS may be affected by attacks such as jamming and spoofing [1], or be severely affected in harsh propagation conditions [2]. Moreover, standard GNSS may not provide the precision needed in the ITS context, i.e., submeter lane-level precision. Several GNSS carrier phase-based precise PNT solutions exist to improve the latter, for instance, Real-Time Kinematics (RTK) ( [3] Chapter 26) or Precise Point Positioning (PPP) techniques ( [3] Chapter 25), but they need either a reference station or precise corrections which do not allow real-time processing, and are even more affected by harsh propagation conditions than standard are given in Section 3. The derivation of new robust filtering techniques under model mismatch is detailed in Section 4. The performance of such methods, compared to state-of-the-art solutions, is shown and discussed in Sections 5 and 6. Finally, conclusions and remarks are drawn in Section 7.

Nominal State-Space Model
Consider the positioning of a mobile agent (i.e., car, bicycle, pedestrian, or UAS), where at time t both 3D position and velocity are to be inferred, p t = [x t , y t , z t ] T and v t = [v x,t , v y,t , v z,t ] T . If L transmitters (Tx) (i.e., UWB anchor nodes or LTE cellular base stations), within the communication range of the target, are located at fixed 3D positions p i = [x i , y i , z i ] T , with i = 1, . . . , L, then the measured Tx to agent distance/range is given by where || · || is the L2-norm, and n i a measurement noise. The full set of available observations (measurement equation) is written as In the nominal case, p i are the true Tx positions and the measurement noise is zero-mean Gaussian, n t ∼ N (0 L , R t ), with known diagonal covariance R t and the individual measurement variance σ 2 r,i driven by the underlying ranging technology. Notice that in real-life applications the number of available observations may be time-varying (i.e., the number of Tx in view or within the agent communication range), then the system must be aware that the size of the measurement vector may be changing over time.
Regarding the dynamics of the agent, it is common practice, and without loss of generality, to consider a constant velocity model (process or state equation), (2) and (3) define the nominal state-space model (SSM) formulation of the TOA-based navigation problem, that is, at every time t estimate x t from available measurements up to time t, z 1:t . In the sequel we introduce the concepts of model mismatch and harsh environments.

Non-Nominal Mismatched State-Space Model
To further introduce the main problem to be considered in this contribution, we define the mismatched SSM and how to account for non-nominal propagation conditions.

Mismatched SSM
In real-life applications, we may have a partial knowledge or uncertainty in a subset of Tx positions. In this case, a Tx position mismatch, on a subset U e of L e ≤ L, is written as where ∆p i = [∆x i ∆y i ∆z i ] T is the position error on the i-th Tx, which can be viewed as a bias in its position. The mismatched measured distance (for i ∈ U e ) is and p t and ∆p i (or equivalentlyp i ) being unknown. We have in this case a new measurement vectorz T t = [z 1,t , . . . ,z L e ,t , z L e+1 ,t , . . . , z L,t ] (the first L e distances are affected by the mismatch and the remaining ones are considered to be obtained under nominal known Tx position conditions). The corresponding measurement equation is theñ wherex t is the augmented state vector of (3) which accounts for the mismatch (refer to Section 3.2).

Harsh Propagation Conditions
It is known that, for instance, in urban environments or indoors, there are several propagation effects which deviate from Gaussianity, then being crucial to design robust methods in order to obtain a reliable solutions. A possible way to account for these conditions is to consider parametric non-Gaussian measurement noise distributions, i.e., Student-t, Laplace, or skew-t [27,28], but then the robust filtering solution is constrained to the imposed parametric noise modeling and thus not flexible to time-varying propagation conditions. Another approach is to consider the typical contamination model arising in robust statistics [2,26]: consider a proportion 1 − ε of observations under nominal Gaussian noise, and another proportion 0 ≤ ε ≤ 1 of observations contaminated by an unknown distribution, where H is an arbitrary contaminating distribution accounting for possible outliers (i.e., corrupted observations), for instance, a non-zero mean Gaussian distribution with variance σ H σ r,i . The main question is how do we estimate the agent's position and velocity under both Tx position mismatch and non-nominal propagation conditions?

Standard Extended KF Solution
Considering the SSM in (2) and (3), it is easy to design an EKF [29], which needs a linearized (approximated) measurement equation, given by the following Jacobian matrix (evaluated at a given point x 0 and with u i,t = ||p 0 − p i ||), The standard TOA-based navigation solution is obtained by applying the EKF to the assumed models (2) and (3): Notice that, for t ≥ 1, the linear KF recursion is only valid if Σ x,0|0 = C x 0 and x 0|0 = E[x 0 ] (mean and covariance of the initial state) [30], which in practice are unknown. An appropriate way to initialize the filter is thus needed.
A practical solution is to use as initial estimate the linear minimum variance distortionless response (LMVDR) estimator, which coincides with the weighted least squares estimator (WLSE) of x 1 , and is given by In this case, due to the nonlinearity of the measurement equation, the WLSE must be replaced by an iterative WLSE.
In case of model mismatch, the standard EKF-based solution is blind because it directly uses the assumed (mismatched) measurementsd i,t = ||p t −p i ||. This solution can always be used and may provide acceptable results (depending on the application requirements) if d 2 i,t >> |ε t (p t , ∆p i )|.

Augmented-State EKF Solution
An alternative to the previous standard EKF, which does not take into account the possible model mismatch, was recently proposed in [24]. The idea is to consider an augmented state, x t = [p t v t p 1,t . . . p L e ,t ] T , that is, to include the (partially) unknown p i into the state to be estimated. The process equation in this case is given by and the measurement equation as in (2). The corresponding Jacobian matrix (dimension L × (6 + 3L e )) isH with ||p−p i || being the i-th LOS unit vector. It has been shown in [24] that this provides better results w.r.t. the standard EKF if not all the Tx positions are partially unknown (i.e., L e < L). In order to take into account the Tx position uncertainty, the measurement noise covarianceR t is modified w.r.t. R t , where the variance of the possibly mismatched measurements is set to A possible choice of σ 2 p,i is to consider a maximum bias, for instance, a uniformly distributed random position bias However, this solution only tackles the Tx mismatch problem and is not designed, as for the standard EKF, to cope with possible harsh propagation conditions, i.e., non-nominal perturbations on the noise n t . Therefore, as it will be shown in Sections 5 and 6, it is not a robust solution and new methodologies must be designed for real-life applicability.

Robust Filtering under Model Mismatch
How to deal with measurements under the contamination model in (7) has been studied for several decades in the field of robust statistics. A recent publication provides an excellent introduction with several applications to practical signal processing problems [26]. In the sequel, we first give a brief introduction on the idea behind the most basic robust methods, and then explain how these concepts can be exploited to design robust filtering methods under model mismatch.

The Idea behind Robust Estimation
We consider as an example the linear regression problem, y = Hx + n, where we want to estimate x and the noise components are i.i.d. If we define the residuals r = y − Hx the solution is given by the least squares (LS) estimator, where we introduced the normalization of the residuals. Notice that a single outlier can destroy our estimate; thus, it is not a robust solution. Instead of using a quadratic function of the residuals, we can consider a general loss function ρ (·) (and its derivative for the minimization, ψ(x) = ∂ρ(x) ∂x ) aŝ and, for instance, ρ LS (x) = x 2 and ρ 1 (x) = |x| correspond to the LS and 1 estimators. The most common ρ (·) is the so-called Huber function: where the value a is chosen to obtain a given efficiency (deviation from the optimal under nominal conditions), for instance, a = 1.345 for a 95% efficiency. The idea is that residuals with large errors are downweighted. To solve the problem in (15) we need the residuals standard deviation σ, or a robust estimate of it, for instance the normalized median absolute deviation (MAD),σ M , defined aŝ where Med(x) is the median and c m a normalizing constant (typically c m = 1.4815). Considering the weight function in (16c) it is easy to see that (15) can be rewritten as which can be solved by an iterative reweighted LS and is the so-called regression M-estimator. Notice that other cost functions ρ (·) can be considered, and more evolved estimators such as the MM-estimator exist [26]. As discussed in [2], for TOA-based navigation problems an important point is the measurement redundancy, that is, having enough observations to be able to find the solution.

Standard Robust Regression EKF
The robust regression M-estimator introduced in Section 4.1 can be used within the EKF framework in order to obtain a robust EKF [26]. Notice that the problem of interest is the one where we may have outliers in the observations z t and the residuals are related to the innovation sequence z t − h t (x t|t−1 ). A first approach is to use a robust score function ψ(·) in order to directly downweight the innovation vector, then not modifying the EKF recursion. Another more general approach is to reformulate the EKF as a regression problem and then use at every t a robust regression M-estimate of x t , being the so-called robust regression EKF (RREKF) ([26] Chapter 7). Notice that both approaches were not developed to cope with a possible model mismatch and they are not directly suited for the navigation problem in this contribution (i.e., as it will be shown in Section 6). Obviously, one can combine the original RREKF with the augmented state in Section 3.2, but the underlying problem to apply these techniques to the mismatched SSM is that the filter is not able to distinguish between true measurement outliers and measurements which deviate from the nominal due to the mismatch. Therefore, new robust solutions under mismatch must be sought.

Robust Regression EKF for Mismatched Models
Carefully analyzing the robust M-estimator downweighting process, that is, ρ σ M , we see that the normalizing termσ M is a key parameter for the filter. Indeed, ifσ M is larger than the nominal variance, large outliers will not be downweighted and degrade the solution. In contrast, if this term is lower than the nominal then meaningful information will be disregarded by the filter. Recall that the problem is to deal with possibly corrupted measurements where the noise is distributed according to the contamination model in (7). In addition, some of the measurements may be biased due to the model mismatch. Then, in practice, we can modify (7) to include the effect of the mismatch as With the nominal distribution being N (0, σ 2 r,i ), the main problem of the RREKF is that it is not able to correctly distinguish among N (0, σ 2 r,i + σ 2 p,i ) and H. Moreover, within the augmented state EKF framework, measurements arising from N (0, σ 2 r,i + σ 2 p,i ) should not be understood as faulty measurements, because these observations are the ones that allow the filter (which includes [p 1,t . . . p L e ,t ] T in the state) to converge to the true Tx positions, and recursively correct the mismatch.
Then, a possible robust EKF solution which is able to cope with model mismatch (named MRREKF), and thus avoid to downweight measurement errors induced by model mismatch is first, consider the augmented-state formulation in Section 3.2 to cope with the Tx mismatch; (ii) second, consider a RREKF-based solution; and (iii) within the M-estimator, instead of the MAD, use a mismatched three sigma rule to normalize the residuals, that is,σ This may be seen as a compromise between the 3σ-rejection or score function type robust EKF ([26] 7.3.1) and the robust regression KF ([26] 7.3.5), in order to control the filter behavior under model mismatch. Surprisingly, to the best of our knowledge, this simple solution has not been explored in the literature. As it will be shown in Section 6, this method behaves extremely well for the problem under study. The price to be paid w.r.t. the standard and augmented state EKFs, is that this solution includes an iterative procedure to find the robust estimate at every Kalman iteration. The iteration procedure terminates either after a fixed number of iterations or when the normalized difference between the actual and previous state estimate is smaller than a fixed threshold. Thus, the computational time will mainly depend on these two parameters. This may be critical in very low cost platforms and real-time applications. This is the reason why we propose an low cost alternative in the sequel.

Computationally Efficient Solution: Robust Weighting Uncertainty Indicator for Robust EKF
One of the limitations of the previous RREKF and MRREKF is their computational complexity, that is, within every EKF time step there is an iterative reweighted LS algorithm to solve the regression M-estimation problem, which, depending on the application at hand, may be a limiting factor. Therefore, we consider another computationally light alternative to the RREKF, which relies on the use of a robust weighting function in order to update the statistical characteristic of the measurement noise variance. Indeed, the performance of the EKF is highly affected by the measurement covariance matrix, R t ; therefore, one or more contaminated measurements will have a significant impact on the filter performance.
The proposed method is based on the values taken by the normalized residuals (i.e., innovations) as introduced in Section 4.1. In order to reduce the impact of the contaminated measurements on the estimator, the idea is to increase the variance of each contaminated measurement depending on the value of |r i /σ M |. To do so, we propose in this study to use the square inverse of Huber's weighting function, as shown in Figure 1, to adjust the measurement noise variances. Thus, in the case of large residuals, i.e., |r i /σ M | > a, and in the framework of Tx position mismatch, the variance can be formulated as Replacing w (r i /σ M ) by (16c), the latter can be rewritten as where we can see that the variance associated with a contaminated measurement is proportional to the square of the normalized residual. In other words, the larger the outlier is, the larger is the associated variance, which in turn implies a direct downweight within the EKF through the Kalman gain. Therefore, with a minor EKF modification, the filter is able to use the robust weight as an outlier uncertainty indicator, which at the end does not impact the augmented-state EKF behavior to correctly deal with model mismatch. This new robust covariance weighting method is subsequently denoted MRCEKF.

Illustrative UWB-Based Indoor Navigation Example
The performance of the new robust filtering methods in Sections 4.3 and 4.4 is first assessed in an indoor environment where only UWB measurements are available. The set of methods considered is as follows.

•
Standard EKF (SEKF) in Section 3. The Huber's score function (a = 1.345), as defined in Section 4, is used to compute the RREKF, MRCEKF, and MRREKF state estimates. The results are obtained from 100 Monte Carlo runs.

Simulation Scenario
We consider L = 8 anchors and a realistic mobile agent trajectory, shown in Figure 2. The position of five anchors (A4-A8) is considered to have a 3-dimensional bias which is drawn randomly in [−0.5, 0.5] m, for each Monte Carlo run. The unknown state vector to be estimated is then x t = [p t , v t ,p 4 , . . . ,p 8 ] T . The initial position and velocity of the agent were set to p m,0|0 = (0, −2, 0.5) T and v m,0|0 = (0, 0, 0) T plus a random initial bias of ±0.5 m and ±0.01 m · s −1 , respectively. The noise of the inlier observations is described by a zero-mean Gaussian distribution whose variance was set to σ 2 r i ,LOS = 0.01 m 2 . Similarly, the multipath/NLOS effect (outliers) in the observations is modeled as a zero-mean Gaussian distribution with a variance α 2 times larger than that of the inlier observations. The mixed LOS/NLOS conditions are modeled as a Markovian process that switches between the LOS and multipath/NLOS conditions. Different percentages, ε, for the the probability of being in multipath/NLOS conditions, and outlier magnitudes, α, are considered and indicated in Table 1.
The number of anchors in multipath/NLOS conditions is drawn randomly within the 8 anchors, and the number of regions and length of each region in mixed LOS/NLOS conditions are also randomly generated to correctly analyze the statistical behavior.

Results
The MATLAB simulations were run on a Linux Ubuntu 18.04 OS-based laptop with an Intel R Core TM i5 processor at 1.7 GHz and 16 Go RAM. First, we show in Figure 2a single realization, and the corresponding five estimates, to clearly illustrate both the impact of the mismatched anchors (blue circles) and harsh propagation conditions (outliers), which are evident in such results. The former can be seen in the 3D plot where the standard SEKF and RREKF are not able to deal with the position mismatch and are affected by a misalignment of the estimated trajectory. Even if the RREKF correctly discards outliers (see 2D plot) it trusts the mismatched anchor positions as being the real ones. The SEKF is not able to cope with neither the mismatch nor outliers. Regarding the augmented state techniques, the MEKF is able to deal with the mismatch but is severely affected by outliers which are understood as Tx position errors. Notice that only the new MRCEKF and MRREKF provide a smooth trajectory, properly dealing with both position mismatch and corrupted measurements.
To give an idea of the order of magnitude in the computational complexity increase of the MRREKF with respect to the MRCEKF (which has no nested loops), the computational time was evaluated using MATLAB's profiler, for the whole simulation. Setting the iteration break to 25 and the threshold to 10 −4 for the MRREKF (and RREKF), it turns out that the total computational time for the MRREKF is of 14.688 s against 4.640 s for the MRCEKF. Further analysis is needed if one wishes to implement the algorithms on a real-time embedded system, which we leave for future work. Figure 3 shows the average horizontal and vertical root mean square errors (RMSE) obtained with the five estimators as a function of the multipath/NLOS standard deviation (i.e., σ r i ,NLOS = ασ r i ,LOS ). We provide the results obtained for two, four, six, and eight anchors under multipath/NLOS conditions, and ε = 10, 25 and 50%.

•
As expected, for a given percentage of contamination, if the outliers are not very strong (α ≤ 1 m), the MEKF performs much better than the SEKF, but it breaks down rapidly as the amplitude of the outliers increases. Indeed, the MEKF tries to compensate the modeling error due to the mismatch, but the measurements, which are corrupted by outliers, are understood by the filter as a bias in the anchors' position. In comparison, the SEKF is less affected by this behavior, with no clear performance breakdown, but in general this method provides a worse estimate when compared to the two robust estimators. • W.r.t. the MEKF behavior, the standard RREKF estimator is relatively stable, thus correctly dealing with outliers, but its nominal performance is degraded because of the mismatch, therefore not bringing a suitable solution for the problem at hand. For low outlier magnitudes, α, it is clear, when comparing the outcome of both MEKF and RREKF, that a method being able to cope with both mismatch and outliers should retain the best qualities of both estimators. • Indeed, the compromise between the MEKF and RREKF is brought by the new MRCEKF and MRREKF, which are relatively stable and cope perfectly with contaminated observations and anchor position mismatch, up to a certain level of contamination. Notice that the interest of the MRREKF is clear with a large outlier amplitude. The MRCEKF outperforms the rest when the number of contaminated measurements is less or equal than 6 and ε ≤ 50%, and the MRREKF provides a good solution up to eight contaminated anchors and for ε = 50%. It can be concluded, from the analysis of these results, that for ε < 50%, α ≤ 3 m and the number of anchors in multipath/NLOS conditions ≤ 6, the MRCEKF approach guarantees horizontal and vertical performances below the tens of centimeters, unlike the other standard approaches. For more extreme propagation conditions, the MRREKF is needed and we can not avoid the iterative M-estimate procedure within every Kalman filter iteration.
These indoor navigation results (i.e., standalone UWB positioning) confirm that both new MRCEKF and MRREKF are promising approaches to provide a high-precision solution in challenging propagation conditions for real-life safety-critical applications. To further support this statement, in the sequel we analyze a realistic outdoor urban environment.

Realistic UWB/LTE TOA-Based Urban Navigation Example
In this section, we want to confront the performance of the SEKF with the new MRCEKF and MRREKF, considering a realistic urban environment, with both UWB and LTE measurements available. In order to correctly assess the behavior of the proposed methods, the choice was to use realistic synthetic data instead of a real measurement campaign. It is known that the latter may lead to unexpected behaviors in the hardware or the data acquisition process, out of the scope of this contribution. Therefore, in order to characterize the filtering methods we resort to a randomized Monte Carlo analysis, which is not environment-dependent, and which may account for all possible situations encountered in real-life environments (i.e., tunnels, tree-lined streets, or tall buildings) by inducing different types of noise, outliers, filter initialization, etc.
The scenario implies the use of a mobile agent, i.e., a low-speed ground vehicle that will follow a simulated path composed of four narrow streets. The algorithms operate in a 3D state-space, so they are also valid for low-speed aerial vehicles such as drones. High-speed vehicles are out of the scope of this contribution. The trajectory was generated using a realistic trajectory generator in which the model is based on a constant velocity model. Smooth transitions occur between segments, which use a continuous acceleration profile.

Realistic Urban Scenario Definition
We have chosen a realistic urban scenario in downtown Toulouse, France, where the buildings are tall enough and the streets narrow enough to be in severe urban canyon propagation conditions. We consider that a total of 32 mismatched UWB anchors are placed on four-meter-high street lamps (alternatively on both sides of the street) all along a 4-street trajectory. This scenario is chosen because of the presence of 3 LTE antennas located in the surroundings, installed on the roof-top of buildings, which are 3 fixed and perfectly known anchors (LTE eNodeB type in which characteristics are defined in [17]). In a real-life scenario, an RTK-based solution could be used to determine the positions of the LTE base stations with a centimeter-level accuracy. The origin, i.e., (0,0,0), of the East-North-Up (ENU) coordinate frame is defined by the LTE base station number 1. The other two base station coordinates, as well as the UWB anchors, are then expressed in this ENU frame. Refer to Figures 4-6 for details.   We consider 32 identical UWB anchors, that is, they have the same statistical properties in terms of Gaussian noise standard deviation, measurement range bias, maximum position bias, and multipath/NLOS noise amplitude. Their true position is precisely known in order to obtain the ranges and trajectory ground truths, which will not be available to the user (indeed, we recall that all 32 UWB anchors have a position mismatch). The mobile agent follows the 4-street trajectory counterclockwise and spends 50 s in each street. A street is defined by the trajectory between two consecutive corners. We arbitrarily choose the mobile agent's starting point, with a random initial position bias of ±3 m and zero initial velocity with a random initial bias of ±0.01 m·s −1 . Because of the urban canyon nature of this scenario, the mobile agent does not have a communication link available to all UWB anchors, but only to the 7 nearest ones: the reference anchor, i.e., the closest to the agent; 3 upstream anchors; and 3 downstream anchors. Considering the 3 LTE base stations, the mobile agent has a total of 10 signals available at each time step.
An example of this situation is depicted on Figure 7, where we can see the real trajectory of the mobile agent (in blue), the three LTE antennas with their labels, and the 32 UWB anchors represented by the black spots near the trajectory. As we use 3 LTE antennas, the UWB anchors will have labels from 4 to 35. The starting point of the mobile agent is shown by the blue cross near anchor 4. In this situation, the agent is the red spot located near anchor 12, which is its current reference anchor, thus the 7 UWB anchors having a communication link with the mobile agent will be anchors 9 to 15. Moreover, we also consider that the percentage of outliers affecting the different signals available is not the same for all visible anchors, and changes depending on the mobile agent's location. We therefore make the following simulation hypothesis.

•
Anchors located in the agent's reference anchor street have an outlier percentage of ε = 10, that is, mild multipath probability due to the narrow streets set-up.

•
Anchors located in a corner of the trajectory have an outlier percentage of ε = 25, that is, a larger multipath conditions probability.

•
The rest of the anchors not on the same agent's reference anchor street have an outlier percentage of ε = 50, induced by the lack of direct visibility, i.e., NLOS conditions.
Because the percentage of outliers changes during the agent's trajectory every time it switches reference anchor, we will neither study the influence of the number of mismatched anchors, nor the number of anchors under multipath/NLOS conditions, but rather the influence of the anchor position bias and the outlier magnitude α. We consider anchors' maximum position bias varying from 1 m to 10 m, and α changing from 30 to 300, which means a multipath/NLOS standard deviation from 3 m to 30 m. The set-up is summarized in Table 2.

Results
The average horizontal and vertical RMSE have been evaluated for each estimator over 1000 Monte Carlo runs, as shown in Figures 8 and 9. To obtain a representation of any urban environment, random anchors' position bias, noise, and different types of outliers (Table 2) were generated at each Monte Carlo run.
At first sight, the SEKF and RREKF provide similar results as in Section 5, i.e., they are not able to cope neither with the position bias nor with outliers, whereas the new robust algorithms better deal with both effects, providing a horizontal RMSE which is almost always sub-metric, no matter the position bias or the outlier magnitude. However, when we focus on the curves' trend, we can notice that increasing the position bias tends to flatten the curves, as if the algorithms were more sensitive to the anchors' position bias than to the outlier magnitude, which means that the robustness to harsh propagation conditions is always effective. For instance, the MRREKF is more responsive to the multipath/NLOS standard deviation when ∆p i ∈ [−1, +1] (the curve is gently increasing with α). In comparison, when the maximum position bias increases to 10 m, the curve is almost flat and the horizontal RMSE is sub-metric whatever the multipath/NLOS magnitude. In that case, the algorithm struggles to deal with the high position bias and is not impacted by the fact that the magnitude of outliers is more and more important.
Concerning the vertical RMSE, all algorithms show meter-level performances in most of the cases. Indeed, the vertical position estimate suffers from large uncertainties due to the poor vertical diversity of the UWB anchors, and the LTE base stations as well. As pointed out in the introduction, the network geometry has an impact on the minimum achievable performance, but such analysis is out of the scope of this contribution. Nevertheless, these results once again depict the performance improvement provided by the two new robust algorithms, which are the MRCEKF and the MRREKF.

Conclusions and Outlook
In this contribution, we derived two new robust filtering techniques able to deal with both model mismatch (i.e., transmitters' position uncertainty) and harsh propagation conditions (i.e., measurements corrupted by outliers). Such methodologies rely on an augmented state, which includes the uncertain anchors' position, and tools arising from robust statistics. It was also shown that standard techniques suffer from a non-negligible performance loss if mismatch and outliers are not properly taken into consideration; therefore, the new methodologies are promising solutions for real-life safety-critical applications. The first methodology is based on the standard robust regression M-estimation EKF, with the penalty of needing an iterative procedure within every KF step, which may not be a good option for very low cost platforms. The second methodology focused on the latter, to reduce the computational complexity, by leveraging a robust weighting function to adapt the measurement noise statistics. Overall, both methods provide a significant performance improvement under several configurations with respect to the state-of-the-art.
Note that the vehicle model used in this contribution is a low-speed ground vehicle, which could correspond to a car moving through small and narrow streets in a city center. Notice that the algorithms operate in a 3D state-space, so they are also valid for low-speed aerial vehicles such as drones. High-speed vehicles are out of the scope of the contribution, and indeed, even if the KFs can also be used in this case in such environments, they are not likely to have strong outliers (open-sky environments or rural areas). Therefore, even if the new filtering approaches are valid for any type of vehicle, dynamics, and propagation conditions, their clear benefit is in situations where GNSS is not available and the propagation conditions are severe. Notice that the focus of this contribution was on GNSS-denied harsh environments such as indoors and deep urban canyons, but the new methodologies proposed in the article can also be applied in the context of GNSS and multisensor data fusion (where other types of mismatch may arise), setting the basis for a robust multisensor navigation framework to be further explored.
Another issue worth mentioning is time synchronization. First, notice that the standard way to exploit UWB for positioning is through a two-way ranging protocol within the UWB network. This avoids the need to synchronize the UWB Tx/Rx. In addition, there are several ways to synchronize or time-stamp LTE measurements together with UWB measurements, for instance, using the NTP networking protocol for clock synchronization. One could also imagine that a GNSS receiver is also available but not used for positioning due to its very degraded performance. Notice that even in very harsh propagation conditions, a GNSS receiver provides a clock synchronization in the order of 100 ns, which is more than enough to avoid worrying about the LTE data synchronization. Therefore, the proposed approach without taking into account the LTE clock errors is valid to statistically characterize the new filtering methods, which are intended to deal with model mismatch and outliers.
As a short-term perspective, these techniques must be assessed with real-data experiments, which is currently under study because other problems related to the hardware and the specific platform appear. In the long-term, such techniques could be embedded in a complete GNSS/UWB/LTE/IMU/Lidar/Vision platform to deal with a full set of heterogeneous measurements and model mismatches.