2.1.2. 3 − *σ* and 1−*σ* Dual w-Test

Different from the traditional w-test, this paper firstly adopts 3 − *σ* w-test. The first 3 − *σ* w-test is to prevent the observation noise variance from being too small, as well as too strict, for the corresponding w-test threshold. At the same time, in the 3 − *σ* w-test, the method of excluding satellites is not to eliminate the |*wi*|max corresponding satellite, but it is to eliminate the corresponding satellite with the largest absolute value of the predicted pseudorange residual when the |*wi*|max exceeds the threshold. The predicted pseudorange residual is calculated as (8):

$$
\Delta \rho = \rho\_{IMII} - \rho\_{GNSS} = r\_{IMII} + c \left( dt\_R - dt^S \right) + I\_\rho + T\_\rho - \rho\_{GNSS} \tag{8}
$$

where, *ρGNSS* is the observed pseudorange, *ρIMU* is the pseudorange predicted by *IMU*, *rIMU* is the geometric range between the observed satellite and the user position estimated by *IMU*. *dtR* and *dt<sup>S</sup>* are the receiver and satellite clock errors, respectively, *I<sup>ρ</sup>* and *T<sup>ρ</sup>* are tropospheric and ionospheric corrections, respectively. The *rIMU* can be calculated as (9):

$$r\_{IMU} = \sqrt{\left(X\_k^G - X\_k^{IMU}\right)^2 + \left(Y\_k^G - Y\_k^{IMU}\right)^2 + \left(Z\_k^G - Z\_k^{IMU}\right)^2} \tag{9}$$

 *X<sup>G</sup> <sup>k</sup>* ,*Y<sup>G</sup> <sup>k</sup>* , *<sup>Z</sup><sup>G</sup> k* is the satellite position at epoch *k*, *XIMU <sup>k</sup>* ,*YIMU <sup>k</sup>* , *<sup>Z</sup>IMU k* is user positions estimated using IMU data at epoch *k*. However, due to the complexity of urban environments, it is impossible to ensure correct detection using the 3 − *σ* w-test. Therefore, the positions calculated before and after each 3 − *σ* w-test are saved until either no faulty satellite measurements are detected or the number of remaining observed satellites is insufficient. Then, in order to ensure that multiple faults can be detected, this paper takes each subset obtained after the 3 − *σ* w-test, removing a satellite each time, and performing the 1−*<sup>σ</sup>* w-test on *<sup>C</sup>*<sup>1</sup> *<sup>m</sup>* each subset. The results can be one of four cases:


The fault conditions at a given epoch can be determined by considering the test results in these four cases. In case 1, we consider that there is no faulty satellite at this epoch, as the universal set and all subsets have passed the w-test. In case 2, the high correlation of each satellite will result in the universal set passing the test, while the low correlation of the faulty satellite in the subset, after one satellite exclusion, can result in the subset not passing the test. Therefore, in this case, we consider that there are multiple faults. In case 3, as a single satellite fault can lead to the universal set not passing the w-test, the subset can only pass the w-test in the case that the faulty satellite is excluded. Therefore, a single fault case is considered in this case. In case 4, faulty satellites in the universal set and subsets can lead to the failure to pass the w-test for a part or all of the subsets. Therefore, the existence of multiple faults is considered in this case. Satellite selection is then made according to the fault conditions. In case 1, all of the satellites at this epoch are selected for a further GNSS/IMU integration. In case 3, the satellites in the subset, which passed the w-test, are selected for further fusion. Considering cases 2 and 4 with multiple faults, the *C*2 *<sup>m</sup>* subsets are further generated, which are then subjected to range detection. The range is calculated by the difference between the predicted position estimated by the IMU data and

the position calculated by the selected subset in the proposed algorithm. The expression for range detection is:

$$
\begin{bmatrix} c \\ n \\ u \end{bmatrix} = \begin{bmatrix} -\sin\lambda\_0 & \cos\lambda\_0 & 0 \\ -\sin\varphi\_0 \cos\lambda\_0 & -\sin\varphi\_0 \sin\lambda\_0 & \cos\varphi\_0 \\ \cos\varphi\_0 \cos\lambda\_0 & \cos\varphi\_0 \sin\lambda\_0 & \sin\varphi\_0 \end{bmatrix} \begin{bmatrix} x\_s - x\_0 \\ y\_s - y\_0 \\ z\_s - z\_0 \end{bmatrix} \tag{10}
$$

where *λ*<sup>0</sup> and *ϕ*<sup>0</sup> are, respectively, the latitude and longitude corresponding to the predicted position. *xs*, *ys*, *zs* and *x*0, *y*0, *z*0, respectively, are the coordinates of the calculated position and the predicted position in the WGS-84 coordinate system. Then,|*e*|,|*n*|, and |*u*| are compared with the range threshold. Here, the threshold of the range value is set as 17 m, as the city speed limit is around 60 km/h (i.e., 17 m/s). Only the subsets that pass the range detection test are used further for the optimal subset selection.

#### *2.2. Scoring Strategy Based Optimal Subset Selection*

After range detection, the subsets that pass the test are selected. The optimal subset within these selected subsets is chosen, and the corresponding measurements in the optimal subset are used to integrate with the IMU data to calculate position. The strategy uses a scoring mechanism to subtract the positions calculated using the selected subsets from the predicted position at the current epoch. The predicted position can be obtained from that of the previous epoch combined with inertial navigation information. The difference in position is then scored according to the following formula, based on a weighting method, in which the smaller the *JointCost* the higher the score. Finally, the satellites corresponding to the position difference with the highest score are selected to be combined with the inertial navigation. The *JointCost* is calculated as:

$$\text{J/JointCost} = \frac{\text{Cost}(1) - \text{Cost}1min}{\text{Cost}1max - \text{Cost}1min} + \frac{\text{Cost}(2) - \text{Cost}2min}{\text{Cost}2max - \text{Cost}2min} + \frac{\text{Cost}(3) - \text{Cost}3min}{\text{Cost}3max - \text{Cost}3min}}{\text{Cost}3max - \text{Cost}3min} \tag{11}$$

Here, *Cost*1*max*, *Cost*2*max*, *Cost*3*max* are the maximum values of longitude, latitude, and height difference among all the position differences. *Cost*1*min*, *Cost*2*min*, and *Cost*3*min* are the minimum values of longitude, latitude, and height difference among all the position differences. *Cost*(1),*Cost*(2),*Cost*(3) are all the longitude, latitude, and height difference among all the position differences.
