Effect Analysis of GNSS/INS Processing Strategy for Sufficient Utilization of Urban Environment Observations

Abstract

The occlusion of buildings in urban environments leads to the intermittent reception of satellite signals, which limits the utilization of observations. This subsequently results in a decline of the positioning and attitude accuracy of Global Navigation Satellite System (GNSS)/Inertial Navigation System (INS) integrated system (GNSS/INS). This study implements a smooth post-processing strategy based on a tightly coupled differential GNSS/INS. Specifically, this strategy used the INS-estimated position to reinitialize integer ambiguity. The GNSS raw observations were input into the Kalman filter to update the measurement. The Rauch–Tung–Striebel smoothing (RTSS) algorithm was used to process the observations of the entire period. This study analyzed the performance of loosely coupled and tightly coupled systems in an urban environment and the improvement of the RTSS algorithm on the navigation solution from the perspective of fully mining the observations. The experimental results of the simulation data and real data show that, compared with the traditional tightly coupled processing strategy which does not use INS-aided integer ambiguity resolution and RTSS algorithm, the strategy in this study sufficiently utilized INS observations and GNSS observations to effectively improve the accuracy of positioning and attitude and ensure the continuity of navigation results in an obstructed environment.

1. Introduction

With the development of the geospatial information service industry, the demand for rapid and accurate geospatial information access has increased. As an important means of obtaining spatial position information, Global Navigation Satellite System (GNSS)/Inertial Navigation System (INS) technology has wide application prospects in mobile measurement [1], autonomous driving [2], and intelligent service [3] fields. GNSS has high long-term absolute positioning accuracy, but its anti-interference ability is weak, and its sampling rate is low. INS can navigate autonomously without being affected by the external environment, and it exhibits high short-term accuracy and sampling rate, but it has one significant problem by which errors accumulate over time [4,5]. Therefore, GNSS/INS integration can achieve the complementary advantages of high-precision positioning and attitude determination. GNSS/INS integration is typically applied in urban environments, which are extremely complex, with dense high-rise buildings, viaducts, tunnels, and other infrastructure. Of these, the most typical environment is the urban canyon, which is formed by dense building blocks and tall trees. In this environment, the GNSS signal suffers frequent lock loss, and the existing observation data are not fully utilized and do not meet the demand for high-precision positioning and pose determination. Therefore, it is necessary to analyze the performance of GNSS/INS to improve the accuracy of positioning and attitude determination in this environment.

In recent years, domestic and foreign researchers have investigated partial and complete GNSS signal lock of loss. The improved methods mostly increased the accuracy of the integrated system by adding observations; these can be divided according to whether additional sensors are added or not. One method involves maintaining heading and speed updates using observations from external sensors or equipment, such as odometers and laser scanners. In [6], GNSS, INS, and an odometer are integrated into a single system, and speed information measured by the odometer is used to assist the INS with navigation and calculation when the GNSS signal is interrupted. In [7], GNSS, INS, and light detection and ranging (LiDAR) technology observations are fused using an extended Kalman filter. When there are no GNSS observations, LiDAR observations are used to assist the INS in estimating and compensating for the sensor error. In [8], based on graph optimization theory, the GNSS position, inertial measurement unit (IMU) pre-integration result, and relative pose matched by the LiDAR scan are fused, which limits the horizontal position error during satellite signal interruption. Although adding sensors can partially improve the accuracy of the GNSS/INS combined system in a satellite denied environment, it also increases the cost and complexity of the system. Another method is to construct non-holonomic constraints based on the state and law of carrier motion to suppress the drift error of the INS inertial device. In [9], a velocity constraint equation was constructed under the assumption that the carrier is always in contact with the ground and does not slip; an elevation constraint was proposed based on small changes in the carrier elevation over a short period of time. By increasing the redundancy of the measured values, the accumulation of inertial navigation errors is controlled, and the speed constraint can provide higher precision pose results than those of the height constraint. In [10], a heading angular velocity constraint that can enhance dynamic vehicle navigation accuracy was proposed. Analyzing non-holonomic constraint observations provides a strategy for selecting observations for different constraints. This method avoids the addition of more sensors by constructing a constrained observation vector with fixed observations. The actual observations are not constant due to changes in terrain and carrier motion, and as a result, improved navigation accuracy cannot be guaranteed. In addition, this method can only be used in a carrier movement environment that meets the limited constraints.

Scholars have also studied how to improve the positioning and attitude accuracies of integrated systems in an urban environment using GNSS/INS integration methods and filtering models without any additional observations. Currently, GNSS/INS integration is primarily divided into three methods: loosely coupled, tightly coupled, and deeply coupled. The deeply coupled method combines GNSS and INS observations and GNSS signal tracking into one filter and uses GNSS to track satellite signals through INS. It is an integration on the hardware level, which is difficult to implement and has not been widely used [11,12]. The loosely coupled method integrates GNSS with assisted INS, which has the advantage of a simple structure. However, when the number of observation satellites is less than four, the Kalman filter stops working, and the integrated system degenerates into a single system with only INS [13]. The tightly coupled method is a form of mutual assistance between the GNSS and INS. It inputs the original GNSS observations as measured values into a Kalman filter to correct INS navigation errors and inertial device errors. INS simultaneously aid GNSS for cycle slip detection and repair [14] and integer ambiguity resolution (IAR) [15]. In [16], recognizing the high cost of high-precision inertial navigation, a low-cost sensor micro-electro-mechanical system IMU is used for tight coupling. The improvement of tightly coupled positioning performance is analyzed through experiments in three typical actual urban environments. Compared with the loosely coupled systems, tightly coupled systems exhibit superior robustness and are more suitable for satellite signal multi-occlusion environments, such as urban canyons. In [17], by simulating an observation environment with an insufficient number of satellites, the effect of different satellite numbers on the performance of the Precise Point Positioning (PPP)/INS tightly coupled with non-holonomic constraints is analyzed, and it was concluded that the greater the number of observation satellites, the better the positioning accuracy.

Differential GNSS (DGNSS) based on the carrier phase is a commonly used high-precision positioning method. Under the premise that the integer ambiguity is fixed correctly, the DGNSS/INS can achieve centimeter-level positioning accuracy. In the DGNSS/INS tightly coupled model, the integer ambiguity solution is primarily divided into two methods: a centralized Kalman filter [18] and a decentralized Kalman filter [19]. In [20], a tightly coupled real-time kinematic (RTK)/INS algorithm based on a centralized Kalman filter was designed, and a float solution update was used to avoid difficulties in fixing the ambiguity when the number of satellites was less than four. This algorithm can effectively use less than four carrier phase observations to limit the accumulation of inertial device errors and quickly resolve integer ambiguity. However, a centralized filter expands integer ambiguity into the state vector and estimates it with INS errors. Although this method can constrain the ambiguity of adjacent epochs, it increases the filter order and calculation amount. Additionally, for this filter type, the ambiguity at the previous moment affects the current state estimation, and the accuracy and stability of the filter can be seriously affected when the GNSS signal is frequently blocked and the fault tolerance is poor [21]. Decentralized Kalman filtering involves placing the ambiguity parameters into a sub-filter to be individually fixed and then uses accurate carrier phase observations to measure and update the primary filter. In [22], the ambiguity solution and the estimation of the navigation state and inertial device errors are realized using a dual Kalman filter. However, it is difficult to fix the ambiguity in this model when there are fewer than four satellites, which directly affects the accuracy of the navigation solution. In addition, two filters are more difficult to implement. [23] tackles the issues of large computational burden and poor robustness in tightly coupled multi-constellation, multi-frequency, and multi-type observation models by utilizing a sequential Kalman filter for data fusion between the two sensors. For GNSS channel failure detection and prevention, a robust estimation method based on the Gaussian test was proposed. In addition, optimal smoothing based on the Kalman filter can use measurement data to provide the optimal estimate for a certain measurement period. Therefore, when real-time processing is not required, a post-processing smoothing algorithm can make full use of the period observations to improve the pose accuracy [24]. Currently, commonly used post-processing smoothing algorithms include the Rauch–Tung–Striebel smoothing (RTSS) algorithm and the forward-backward smoothing (FBS) algorithm. FBS requires backward filtering, which increases the amount of data storage required, and it only improves the position accuracy. In contrast, in the tightly coupled DGNSS/INS, the double-difference ambiguity is easier to fix, and the accuracy of the position result meets the requirements of surveying and mapping applications. Therefore, there is no need to further process the position solution through FBS [25]. RTSS only requires two steps, forward filtering and backward smoothing, and its data storage capacity is small and easy to program. However, when the ambiguity parameter is estimated as a part of the state vector, it is difficult to improve the accuracy of the navigation solution during initialization and reconvergence [26].

Currently, research on the tightly coupled GNSS/INS in an urban environment is focused on improving the filtering model and adding redundant observations. There are few studies regarding making full use of available observations in an environment with frequent losses of satellite signals. Therefore, this study aims to implement a smooth post-processing strategy based on a tightly coupled carrier phase DGNSS/INS. This strategy utilizes a decentralized filter to fix integer ambiguity, thereby avoiding the existing problem of the RTSS algorithm in the tightly coupled GNSS/INS. After the GNSS signal is interrupted, the ambiguity is quickly fixed by the INS in the ambiguity initialization stage. Then, the result is further smoothed by the RTSS. This study first analyzes the positioning and attitude determination performance of the loosely coupled and tightly coupled DGNSS/INS when the number of observable satellites is insufficient through a simulated satellite lock-loss experiment, and verifies the effectiveness of the RTSS in improving the pose accuracy in an occlusion environment. Then, field experiments are used to verify method effectiveness. The results show that the proposed strategy can effectively improve the pose accuracy even during the intermittent reception of satellite signals. It can provide a more reliable and high-precision solution for navigation and positioning tasks in challenging environments.

2. Materials and Methods

2.1. Differential Tightly Coupled GNSS/INS

In the tightly coupled DGNSS/INS extended Kalman filter model, Correct IAR is the premise driving high-precision carrier phase positioning. Using decentralized filters, the tightly coupled model is divided into two parts: a tightly coupled filter module and an IAR module, which lowers the dimensionality of the filter matrix and improves fault tolerance. The primary function of the IAR module is to estimate the integer ambiguity parameters. When a new satellite is observed, the position information output by the filter assists with IAR. After the satellite signal interruption ends and the integer ambiguity must be fixed again, the position information calculated by the INS is used to assist the reinitialization of the integer ambiguity. The tightly coupled filter module exhibits the simplest system model because the ambiguity is not extended to the state vector. The state vector, covariance matrix, INS pose, and other sequentially stored information are inputted into the RTSS module to optimally smooth the entire filtering process. Figure 1 displays the tightly coupled data processing strategy used in this study.

2.1.1. System Model.

The system model is based on the INS error an Earth-centered Earth-fixed (ECEF) frame, the error equation in the ECEF is selected to construct the system model more conveniently. The INS psi-angle error model in the ECEF is given as [27]

$$\begin{array}{c}\left(\begin{array}{c}\dot{{\mathit{\varphi }}^e}\\ \begin{array}{c}\delta {\dot{\mathit{v}}}^e\\ \delta {\dot{\mathit{r}}}^e\end{array}\end{array}\right)=\left(\begin{array}{c}-\left[{\mathit{\omega }}_{ie}^e\times \right]{\mathit{\varphi }}^e-{\mathit{C}}_b^e\delta {\mathit{\omega }}_{ib}^b\\ \begin{array}{c}{\mathit{C}}_b^e\delta {\mathit{f}}_{ib}^b+\left[\left({\mathit{C}}_b^e{\mathit{f}}_{ib}^b\right)\times \right]{\mathit{\varphi }}^e-2\left[{\mathit{\omega }}_{ie}^e\times \right]\delta {\mathit{v}}^e+\delta {\mathit{g}}^e\\ \delta {\mathit{v}}^e\end{array}\end{array}\right)\end{array}$$

(1)

where · is the differential symbol; $i$, $e$, and $b$ represent the Earth-centered inertial frame, ECEF frame, and IMU body frame, respectively; $\left(\mathit{a}\right)\times$ represents a skew-symmetrical matrix composed of vector $\mathit{a}$; $\mathit{\varphi }$, $\delta {\mathit{v}}^e$ and $\delta {\mathit{r}}^e$ are the attitude error vector, velocity error vector, and position error vector, respectively; ${\mathit{C}}_b^e$ is the direction cosine matrix of the IMU body frame relative to the ECEF; ${\mathit{\omega }}_{ie}^e$ is the earth rotation;$\delta {\mathit{g}}^e$ is the gravity error vector in $e$-frame; ${\mathit{f}}_{ib}^b$ is specific force vector; and $\delta {\mathit{f}}_{ib}^b$ and $\delta {\mathit{\omega }}_{ib}^b$ are the measurement error for accelerometers and gyroscope, respectively. The specific derivation process of Equation (1) can be found in Appendix A.

The tactical-grade IMU used in this study does not consider the scale factor error or cross-coupling error, but the random errors of inertial devices need to be considered. The IMU measurement error model is as follows

$$\begin{array}{c}\begin{array}{c}\delta {\mathit{\omega }}_{ib}^b=\delta {\mathit{b}}_g+{\mathit{w}}_{\omega }\\ \delta {\mathit{f}}_{ib}^b=\delta {\mathit{b}}_a+{\mathit{w}}_f\end{array}\end{array}$$

(2)

where ${\mathit{w}}_{\omega }$ and ${\mathit{w}}_f$ are white noise of gyroscope and accelerometer, respectively; $\delta {\mathit{b}}_g$ and $\delta {\mathit{b}}_a$ are the gyroscope and accelerometer bias errors respectively, which can be modeled as first-order Gauss–Markov processes:

$$\begin{array}{c}\begin{array}{c}\delta {\dot{\mathit{b}}}_g=-\frac{1}{{\mathit{T}}_g}\delta {\mathit{b}}_g+{\mathit{w}}_b\\ \delta {\dot{\mathit{b}}}_a=-\frac{1}{{\mathit{T}}_a}\delta {\mathit{b}}_a+{\mathit{w}}_a\end{array}\end{array}$$

(3)

where ${\mathit{T}}_g$ and ${\mathit{T}}_a$ are the correlation time of gyroscope and accelerometer, respectively; ${\mathit{w}}_b$ and ${\mathit{w}}_a$ are the driven white noise of gyroscope and accelerometer, respectively.

Based on Equations (2) and (3), the bias error is expanded to the state parameters for estimation. The system equation can be expressed as follows:

$$\begin{array}{c}\dot{\mathit{x}}\left(t\right)=\mathit{F}\left(t\right)\mathit{x}\left(t\right)+\mathit{G}\left(t\right)\mathit{w}\left(t\right)\end{array}$$

(4)

where,

$$\begin{array}{c}\mathit{x}={\left[{\left(\mathit{\varphi }\right)}^T {\left(\delta {\mathit{r}}^e\right)}^T {\left(\delta {\mathit{v}}^e\right)}^T {\left(\delta {\mathit{b}}_g\right)}^T {\left(\delta {\mathit{b}}_a\right)}^T\right]}^T\end{array}$$

(5)

$\mathit{F}$ is the continuous system state transition matrix; $\mathit{G}$ is the continuous system noise distribution matrix; $\mathit{w}$ is the system noise vector. The specific form of the above matrix can be found in [28].

Equation (4) is a continuous dynamics equation, and it must be discretized as follows:

$$\begin{array}{c}{\mathit{x}}_{\mathit{k}}={\mathit{\Phi }}_{\mathit{k}/\mathit{k}-1}{\mathit{x}}_{\mathit{k}-1}+{\mathit{w}}_{\mathit{k}-1}\end{array}$$

(6)

The specific form of each matrix in Equation (6) can be found in [29]. In this study, the power spectral density matrix of the system noise matrix is set up by angle random walk and speed random walk in the inertial navigation manual [30].

2.1.2. Measurement Model

GNSS uses satellites in orbit to transmit radio signals to provide users with position information through passive ranging. The basic observation equation is as follows [31]:

$$\begin{array}{cc}\hfill {\lambda }_i{\phi }_{i,r}^S(t_r)& ={\rho }_r^S(t_r,t^S)-c(d{T}_r-d{T}^S)-{\lambda }_i{N}_{i,r}^S-I_{i,r}^S+T_r^S-(d_{i,r,\phi }-D_{i,\phi }^S)\hfill \\ & +M_{ir,\lambda \phi }^S+R_r^S+{\epsilon }_{i,r,\lambda \phi }^S\hfill \\ \hfill P_{i,r}^S(t_r)& ={\rho }_r^S(t_r,t^s)-c(d{T}_r-d{T}^S)+I_{i,r}^S+T_r^S-(d_{i,r,P}-D_{i,P}^S)\hfill \\ & +M_{i,r,P}^S+R_r^S+{\epsilon }_{i,r,P}^S\hfill \\ \hfill {\dot{P}}_{i,r}^S& =-\frac{c}{f_i}{D}_i\hfill \end{array}$$

(7)

where $i$ is the signal frequency; ${\lambda }_i$ is the wavelength of the carrier; $f_i$ is the frequency of the carrier; $r$ is the receiver number; $s$ is the satellite number; $P_{i,r}^S\left(t_r\right)$ is the code measurement of the pseudorange observations; ${\lambda }_i{\phi }_{i,r}^S\left(t_r\right)$ is the carrier phase observations; ${\dot{P}}_{i,r}^S$ is the pseudorange rate; $N_{i,r}^S$ is the integer ambiguity; ${\rho }_r^S\left(t_r,t^S\right)$ is the geometric distance from the satellite to the receiver; $c$ is the speed of light in vacuum; $d{T}_r$ is the receiver clock error;$d{T}^S$ is the satellite clock error; $I_{i,r}^S$ is the ionospheric delay; $T_r^S$ is the tropospheric delay; $d_{i,r,\phi }$ and $d_{i,r,P}$ are the receiver hardware phase and code delay, respectively;$D_{i,\phi }^S$ and $D_{i,P}^S$ are the satellite hardware phase and code delay, respectively; $M_{i,r,P}^S$ is the multipath effect of the code measurement pseudorange; $R_r^S$ is the multipath effect; ${\epsilon }_{i,r}^S$ is the observation error; $D_i$ is the Doppler frequency shift, where $D_i={\lambda }_i\frac{\left[{\mathit{P}}^S-{\mathit{P}}_r\right]·\left[{\mathit{V}}_r-{\mathit{V}}^S\right]}{\rho }$; $\rho$ is the geometric distance between the satellite and receiver; ${\mathit{P}}^S$ and ${\mathit{P}}_r$ are the satellite and station coordinates, respectively; and ${\mathit{V}}_r$ and ${\mathit{V}}^S$ are the satellite and station velocities, respectively.

There are a variety of error terms in Equation (7), which inhibit integer characteristic maintenance by non-difference ambiguity; thus, the accuracy of the solution directly using the non-difference observation equation is limited. The use of double-difference observations for positioning calculations can eliminate error items such as the receiver clock error, satellite clock error, and receiver and satellite hardware phase delay, and it can greatly reduce the impact of ionospheric and tropospheric delays on the observations to ensure the integer characteristics of double-difference ambiguity [31]. The observation equation of the double difference between the stations is as follows:

$$\begin{array}{cc}\hfill {\lambda }_i\nabla \Delta {\phi }_{i,r_{1,2}}^{S_{1,2}}& =\nabla \Delta {\rho }_{r_{1,2}}^{S_{1,2}}-{\lambda }_i\nabla \Delta N_{i,r_{1,2}}^{S_{1,2}}-\nabla \Delta I_{i,r_{1,2}}^{S_{1,2}}+\nabla \Delta T_{r_{1,2}}^{S_{1,2}}+\nabla \Delta M_{i,r_{1,2},\lambda \phi }^{S_{1,2}}\hfill \\ & +\nabla \Delta d_{r_{1,2},other}^{S_{1,2}}+\nabla \Delta {\epsilon }_{i,r_{1,2},\lambda \phi }^{S_{1,2}}\hfill \\ \hfill \nabla \Delta P_{i,r_{1,2}}^{S_{1,2}}& =\nabla \Delta {\rho }_{r_{1,2}}^{S_{1,2}}-\nabla \Delta I_{i,r_{1,2}}^{S_{1,2}}+\nabla \Delta T_{r_{1,2}}^{S_{1,2}}+\nabla \Delta M_{i,r_{1,2,},P}^{S_{1,2}}\hfill \\ & +\nabla \Delta d_{r_{1,2},others}^{S_{1,2}}+\nabla \Delta {\epsilon }_{i,r_{1,2},P}^{S_{1,2}}\hfill \\ \hfill \nabla \Delta {\dot{\rho }}_{i,r_{1,2}}^{S_{1,2}}& =-{\lambda }_i\nabla \Delta D_{i,r_{1,2}}^{S_{1,2}}+\nabla \Delta {\epsilon }_{\dot{\rho }}\hfill \end{array}$$

(8)

where $\nabla \Delta$ is the double-difference operator; $r_{1,2}$ represents the base station and rover, respectively; and $S_{1,2}$ represents the two satellites simultaneously observed by $r_{1,2}$.

Simplified as:

$$\begin{array}{cc}\hfill {\lambda }_i\nabla \Delta \phi & =\nabla \Delta \rho -{\lambda }_i\nabla \Delta N_i-\nabla \Delta I_i+\nabla \Delta T+\nabla \Delta M_{i,\lambda \phi }+\nabla \Delta d_{others}+\nabla \Delta {\epsilon }_{i,\lambda \phi }\hfill \\ \hfill \nabla \Delta P_i& =\nabla \Delta \rho -\nabla \Delta I_i+\nabla \Delta T+\nabla \Delta M_{i,p}+\nabla \Delta d_{others}+\nabla \Delta {\epsilon }_{i,\lambda \phi }\hfill \\ \hfill \nabla \Delta \dot{\rho }& =-{\lambda }_i\nabla \Delta D_i+\nabla \Delta {\epsilon }_{\dot{\rho }}\hfill \end{array}$$

(9)

when the differential GNSS and INS are integrated, the approximate values of $\nabla \Delta \rho$ and $\nabla \Delta D_i$ are derived from the position and velocity outputs by the INS. Equation (9) is the basis for constructing the tightly coupled extended Kalman filter observation model in this study.

The position and speed corrected by the lever arm value are considered in Equation (9) and linearized to obtain the double-difference observation equation:

$$\left[\begin{array}{c}{\lambda }_1\nabla \Delta {\mathit{\phi }}_1-\nabla \Delta {\mathit{\rho }}_{INS}\\ \nabla \Delta {\mathit{P}}_1-\nabla \Delta {\mathit{\rho }}_{INS}\\ \nabla \Delta \dot{{\mathit{\rho }}_1}-\nabla \Delta {\mathit{v}}_{\rho ,INS}\end{array}\right]=\left[\begin{array}{ccccc}{\mathit{H}}_1& {\mathit{H}}_2& \mathbf{0}& \mathbf{0}& \mathbf{0}\\ {\mathit{H}}_1& {\mathit{H}}_2& \mathbf{0}& \mathbf{0}& \mathbf{0}\\ {\mathit{H}}_3& \mathbf{0}& {\mathit{H}}_4& {\mathit{H}}_5& \mathbf{0}\end{array}\right]\left[\begin{array}{c}\mathit{\varphi }\\ \delta {\mathit{r}}^e\\ \delta {\mathit{v}}^e\\ {\mathit{b}}_g\\ {\mathit{b}}_a\end{array}\right]+\left[\begin{array}{c}\nabla \Delta {\mathit{\epsilon }}_{{\mathsf{\Phi }}_1}\\ \nabla \Delta {\mathit{\epsilon }}_{P_1}\\ \nabla \Delta {\mathit{\epsilon }}_{{\dot{\rho }}_1}\end{array}\right]$$

(10)

where $\nabla \Delta {\mathit{\phi }}_1$, $\nabla \Delta {\mathit{P}}_1$, and $\nabla \Delta \dot{{\mathit{\rho }}_1}$ are the carrier phase, pseudorange, and pseudorange rate double-difference observations, respectively; $\nabla \Delta {\mathit{\rho }}_{INS}$ and $\nabla \Delta {\mathit{v}}_{\rho ,INS}$ are the distance and speed double-difference observation values calculated by INS, respectively; $\nabla \Delta \mathit{\epsilon }$ is the observation error vector; and $\mathit{H}$ is the coefficient matrix of the state parameter given as

$${\mathit{H}}_1=\left[\begin{array}{ccc}\begin{array}{c}\begin{array}{c}\nabla \Delta p_{\mathsf{\Phi },1}\\ \nabla \Delta p_{\mathsf{\Phi },2}\end{array}\\ ⋮\\ \nabla \Delta p_{\mathsf{\Phi },n}\end{array}& \begin{array}{c}\begin{array}{c}\nabla \Delta h_{\mathsf{\Phi },1}\\ \nabla \Delta h_{\mathsf{\Phi },2}\end{array}\\ ⋮\\ \nabla \Delta h_{\mathsf{\Phi },n}\end{array}& \begin{array}{c}\begin{array}{c}\nabla \Delta i_{\mathsf{\Phi },1}\\ \nabla \Delta i_{\mathsf{\Phi },2}\end{array}\\ ⋮\\ \nabla \Delta i_{\mathsf{\Phi },n}\end{array}\end{array}\right] {\mathit{H}}_2=\left[\begin{array}{ccc}\begin{array}{c}\begin{array}{c}\nabla \Delta {\mathit{l}}_{\mathsf{\Phi },1}\\ \nabla \Delta {\mathit{l}}_{\mathsf{\Phi },2}\end{array}\\ ⋮\\ \nabla \Delta {\mathit{l}}_{\mathsf{\Phi },n}\end{array}& \begin{array}{c}\begin{array}{c}\nabla \Delta m_{\mathsf{\Phi },1}\\ \nabla \Delta m_{\mathsf{\Phi },2}\end{array}\\ ⋮\\ \nabla \Delta m_{\mathsf{\Phi },n}\end{array}& \begin{array}{c}\begin{array}{c}\nabla \Delta n_{\mathsf{\Phi },1}\\ \nabla \Delta n_{\mathsf{\Phi },2}\end{array}\\ ⋮\\ \nabla \Delta n_{\mathsf{\Phi },n}\end{array}\end{array}\right] {\mathit{H}}_3=\left[\begin{array}{ccc}\begin{array}{c}\begin{array}{c}\nabla \Delta p_{\dot{\mathit{\rho }},1}\\ \nabla \Delta p_{\dot{\mathit{\rho }},2}\end{array}\\ ⋮\\ \nabla \Delta p_{\dot{\mathit{\rho }},n}\end{array}& \begin{array}{c}\begin{array}{c}\nabla \Delta h_{\dot{\mathit{\rho }},1}\\ \nabla \Delta h_{\dot{\mathit{\rho }},2}\end{array}\\ ⋮\\ \nabla \Delta h_{\dot{\mathit{\rho }},n}\end{array}& \begin{array}{c}\begin{array}{c}\nabla \Delta i_{\dot{\mathit{\rho }},1}\\ \nabla \Delta i_{\dot{\mathit{\rho }},2}\end{array}\\ ⋮\\ \nabla \Delta i_{\dot{\mathit{\rho }},n}\end{array}\end{array}\right]$$

$${\mathit{H}}_4=\left[\begin{array}{ccc}\begin{array}{c}\begin{array}{c}\nabla \Delta {\mathit{l}}_{\dot{\mathit{\rho }},1}\\ \nabla \Delta {\mathit{l}}_{\dot{\mathit{\rho }},2}\end{array}\\ ⋮\\ \nabla \Delta {\mathit{l}}_{\dot{\mathit{\rho }},n}\end{array}& \begin{array}{c}\begin{array}{c}\nabla \Delta m_{\dot{\mathit{\rho }},1}\\ \nabla \Delta m_{\dot{\mathit{\rho }},2}\end{array}\\ ⋮\\ \nabla \Delta m_{\dot{\mathit{\rho }},n}\end{array}& \begin{array}{c}\begin{array}{c}\nabla \Delta n_{\dot{\mathit{\rho }},1}\\ \nabla \Delta n_{\dot{\mathit{\rho }},2}\end{array}\\ ⋮\\ \nabla \Delta n_{\dot{\mathit{\rho }},n}\end{array}\end{array}\right]\hspace{0.17em}{\mathit{H}}_5=\left[\begin{array}{ccc}\begin{array}{c}\begin{array}{c}\nabla \Delta u_{\dot{\mathit{\rho }},1}\\ \nabla \Delta u_{\dot{\mathit{\rho }},2}\end{array}\\ ⋮\\ \nabla \Delta u_{\dot{\mathit{\rho }},n}\end{array}& \begin{array}{c}\begin{array}{c}\nabla \Delta t_{\dot{\mathit{\rho }},1}\\ \nabla \Delta t_{\dot{\mathit{\rho }},2}\end{array}\\ ⋮\\ \nabla \Delta t_{\dot{\mathit{\rho }},n}\end{array}& \begin{array}{c}\begin{array}{c}\nabla \Delta g_{\dot{\mathit{\rho }},1}\\ \nabla \Delta g_{\dot{\mathit{\rho }},2}\end{array}\\ ⋮\\ \nabla \Delta g_{\dot{\mathit{\rho }},n}\end{array}\end{array}\right]$$

The details of above matrix are given as follows:

$$\begin{gathered} p_{\mathsf{\Phi }} =\frac{\Delta y{\left({\tilde{\mathit{C}}}_b^e{\mathit{l}}^b\right)}_z-\Delta z{\left({\tilde{\mathit{C}}}_b^e{\mathit{l}}^b\right)}_y}{{\rho }_{INS}} h_{\mathsf{\Phi }}=\frac{\Delta z{\left({\tilde{\mathit{C}}}_b^e{\mathit{l}}^b\right)}_x-\Delta x{\left({\tilde{\mathit{C}}}_b^e{\mathit{l}}^b\right)}_z}{{\rho }_{INS}} i_{\mathsf{\Phi }}=\frac{\Delta x{\left({\tilde{\mathit{C}}}_b^e{\mathit{l}}^b\right)}_y-\Delta y{\left({\tilde{\mathit{C}}}_b^e{\mathit{l}}^b\right)}_x}{{\rho }_{INS}} \\ {\mathit{l}}_{\mathsf{\Phi }}=\frac{\Delta x}{{\rho }_{INS}} m_{\mathsf{\Phi }}=\frac{\Delta y}{{\rho }_{INS}} n_{\mathsf{\Phi }}=\frac{\Delta z}{{\rho }_{INS}} \\ {p}_{\dot{\rho }}=-\frac{\Delta x{\mathit{A}}_{11}+\Delta y{\mathit{A}}_{21}+\Delta z{\mathit{A}}_{31}}{{\rho }_{INS}} h_{\dot{\rho }}=-\frac{\Delta x{\mathit{A}}_{12}+\Delta y{\mathit{A}}_{22}+\Delta z{\mathit{A}}_{32}}{{\rho }_{INS}} i_{\dot{\rho }}=-\frac{\Delta x{\mathit{A}}_{13}+\Delta y{\mathit{A}}_{23}+\Delta z{\mathit{A}}_{33}}{{\rho }_{INS}} \\ {u}_{\dot{\rho }}=-\frac{\Delta x{\mathit{B}}_{11}+\Delta y{\mathit{B}}_{21}+\Delta z{\mathit{B}}_{31}}{{\rho }_{INS}} t_{\dot{\rho }}=-\frac{\Delta x{\mathit{B}}_{12}+\Delta y{\mathit{B}}_{22}+\Delta z{\mathit{B}}_{32}}{{\rho }_{INS}} g_{\dot{\rho }}=-\frac{\Delta x{\mathit{B}}_{13}+\Delta y{\mathit{B}}_{23}+\Delta z{\mathit{B}}_{33}}{{\rho }_{INS}} \\ {\mathit{l}}_{\dot{\rho }}=\frac{-\Delta x}{{\rho }_{INS}} m_{\dot{\rho }}=\frac{-\Delta y}{{\rho }_{INS}} n_{\dot{\rho }}=\frac{-\Delta z}{{\rho }_{INS}} \end{gathered}$$

(11)

In Equation (11), ${\rho }_{INS}$ is the distance between the satellite and receiver; $\Delta x$,$\Delta y$,$\Delta z$ is the coordinate difference between satellite and receiver; $\left(x^s,y^s,z^s\right)$ are the satellite coordinates; $\left(\tilde{x},\tilde{y},\tilde{z}\right)$ are the receiver coordinates estimated by INS; ${\mathit{l}}^b$ is lever arm vector; where

$$\begin{array}{c}\begin{array}{c}\begin{array}{c}\Delta x=x^s-\tilde{x}\\ \Delta y=y^s-\tilde{y}\\ \Delta z=z^s-\tilde{z}\end{array}\\ {\rho }_{INS}=\sqrt{{\left(\Delta x\right)}^2+{\left(\Delta y\right)}^2+{\left(\Delta z\right)}^2}\end{array}\\ \mathit{A}=\left[{\mathit{\omega }}_{ie}^e\times \right]\left[\left({\tilde{\mathit{C}}}_b^e{\mathit{l}}^b\right)\times \right]+\left[\left({\tilde{\mathit{C}}}_b^e\left({\mathit{l}}^b\times {\widehat{\mathit{\omega }}}_{ib}^b\right)\right)\times \right]\\ \mathit{B}={\tilde{\mathit{C}}}_b^e\left[{\mathit{l}}^b\times \right]\end{array}$$

Further, the observations from Equation (10) can be expressed as

$$\begin{array}{c}{\mathit{z}}_k={\mathit{H}}_k{\mathit{x}}_k+{\mathit{v}}_k\end{array}$$

(12)

For the covariance matrix of the observation noise ${\mathit{v}}_k:$

$$\begin{array}{c}{\mathit{R}}_k=diag\left(\begin{array}{ccc}\mathit{D}{\mathit{R}}_{{\mathsf{\Phi }}_1}{\mathit{D}}^T& \mathit{D}{\mathit{R}}_{P_1}{\mathit{D}}^T& \mathit{D}{\mathit{R}}_{{\dot{\rho }}_1}{\mathit{D}}^T\end{array}\right)\end{array}$$

where $\mathit{D}$ is the single difference matrix,

$$\mathit{D}=\left[\begin{array}{lllll}-1& 1& 0& \cdots & 0\\ -1& 0& 1& \cdots & 0\\ ⋮& ⋮& ⋮& \ddots & ⋮\\ -1& 0& 0& \cdots & 1\end{array}\right]$$

${\mathit{R}}_{{\mathsf{\Phi }}_1}$,${\mathit{R}}_{P_1}$,${\mathit{R}}_{{\dot{\rho }}_1}$ represent the measurement noise covariance of the carrier phase, pseudorange and pseudorange rate, respectively.

$$\begin{gathered} {\mathit{R}}_{{\mathsf{\Phi }}_1}=diag\left(\begin{array}{cccc}2{\sigma }_{{\mathsf{\Phi }}_{1,1}}^2& 2{\sigma }_{{\mathsf{\Phi }}_{1,2}}^2& \cdots & 2{\sigma }_{{\mathsf{\Phi }}_{1,n}}^2\end{array}\right); \\ {\mathit{R}}_{P_1}=diag\left(\begin{array}{cccc}2{\sigma }_{P_{1,1}}^2& 2{\sigma }_{P_{1,2}}^2& \cdots & 2{\sigma }_{P_{1,n}}^2\end{array}\right); \\ {\mathit{R}}_{{\dot{\rho }}_1}=diag\left(\begin{array}{cccc}2{\sigma }_{{\dot{\rho }}_{1,1}}^2& 2{\sigma }_{{\dot{\rho }}_{1,2}}^2& \cdots & 2{\sigma }_{{\dot{\rho }}_{1,n}}^2\end{array}\right); \end{gathered}$$

where

$${\sigma }_{{\mathsf{\Phi }}_1}=\frac{{\sigma }_{{\mathsf{\Phi }}_{1,0}}}{\sin E}; {\sigma }_{P_1}=\frac{{\sigma }_{P_{1,0}}}{\sin E}; {\sigma }_{{\dot{\rho }}_1}=\frac{{\sigma }_{{\dot{\rho }}_{1,0}}}{\sin E};$$

$E$ is the satellite elevation; ${\sigma }_{{\mathsf{\Phi }}_{1,0}}$, ${\sigma }_{P_{1,0}}$, and ${\sigma }_{{\dot{\rho }}_{1,0}}$ are the standard deviations of the carrier phase, pseudorange, and pseudorange rate observation errors, respectively, which can be set according to the corresponding observation accuracy of the GNSS receiver.

2.2. INS-Aided IAR

INS-assisted IAR uses high-precision prior position information provided by the INS to assist the GNSS in obtaining a more precise ambiguity float solution, reducing the search range and improving its accuracy and search efficiency. In particular, when the GNSS signal is partially interrupted, the tightly coupled filter still predicts the parameter error and continuously corrects the navigation result. When the satellite signal is received again, the position and variance matrices obtained by the INS are used to initialize the GNSS navigation module and serve as additional observations. After the integer ambiguity parameters are solved, the GNSS observations and the integer ambiguities are inputted into the tightly coupled filter to update the measurement. In this study, we use the position information predicted by INS as virtual observations. The observation equation for the INS-assisted acquisition of GNSS ambiguity float solution is as follows:

$$\begin{array}{c}\left[\begin{array}{c}{\mathit{L}}_{P_1}\\ {\mathit{L}}_{{\mathit{\varphi }}_1}\\ \mathbf{0}\end{array}\right]=\left[\begin{array}{cc}\mathit{B}& \underset{n-1\times m}{0}\\ \mathit{B}& -{\lambda }_1·\underset{n-1\times m}{\mathit{I}}\\ \underset{3\times 3}{\mathit{I}}& \underset{3\times m}{\mathbf{0}}\end{array}\right]\left[\begin{array}{c}\delta {\mathit{p}}_{r_2}\\ \nabla \Delta \mathit{N}\end{array}\right]+\left[\begin{array}{c}{\epsilon }_{p_1}\\ {\epsilon }_{{\mathsf{\Phi }}_1}\\ {\epsilon }_{INS}\end{array}\right]\end{array}$$

(13)

where $\mathit{I}$ is the identity matrix:

$${\mathit{L}}_P{}_1=\left[\begin{array}{c}\nabla \Delta P_1^{S_1,ref}-\nabla \Delta {\rho }_{INS,0}^{S_1,ref}\\ \nabla \Delta P_1^{S_2,ref}-\nabla \Delta {\rho }_{INS,0}^{S_2,ref}\\ \begin{array}{c}⋮\\ \nabla \Delta P_1^{S_{n-1},ref}-\nabla \Delta {\rho }_{INS,0}^{S_{n-1},ref}\end{array}\end{array}\right], {\mathit{L}}_{{\mathit{\varphi }}_1}=\left[\begin{array}{c}{\lambda }_1\nabla \Delta \phi { }_1^{S_1,ref}-\nabla \Delta {\rho }_{INS,0}^{S_1,ref}\\ {\lambda }_1\nabla \Delta {\phi }_1^{S_2,ref}-\nabla \Delta {\rho }_{INS,0}^{S_2,ref}\\ \begin{array}{c}⋮\\ {\lambda }_1\nabla \Delta {\phi }_1^{S_{n-1},ref}-\nabla \Delta {\rho }_{INS,0}^{S_{n-1},ref}\end{array}\end{array}\right]$$

To improve the ambiguity initialization efficiency and solution accuracy, after obtaining the ambiguity float point solution, we implement different IAR strategies for different navigation stages. In the initialization phase, the float solution of the double-difference wide-lane ambiguity is directly rounded, and the two nearest integer values of each ambiguity are combined with other ambiguity values. Each combination is considered in the error equation for the least-squares solution, and the combination with the smallest error σ is used to calculate the ratio of the second smallest ${\sigma }_1$ to the smallest error $\sigma$.

$$\begin{array}{c}ratio=\frac{{\sigma }_1}{\sigma }\end{array}$$

(14)

when $\sigma$ and the ratio meet certain conditions, the ambiguity is considered to be successfully fixed [32]. In the continuous navigation phase, the Lambda algorithm [33] is used to obtain the integer ambiguity solution. After the satellite signal is interrupted or the innovation vector in the filter measurement update becomes very large, the integer ambiguity must be reinitialized. The GNSS/INS integration module transmits the position estimated by the INS to the GNSS navigation module for IAR.

2.3. RTSS Algorithm Based on Extended Kalman Filter

In the GNSS/INS integration system, the forward Kalman filter can only use historical information to estimate the current state, and it is susceptible to the influence of satellite signal interruption, which leads to great error accumulation. The RTSS algorithm can effectively limit the accumulation of pose errors by fully utilizing forward observations from the initial time to the current time and the backward dynamic constraint information from the end time to the current time [25]. Figure 2 displays the error curve processed by the RTSS when the GNSS signal is not locked. The error curve after the RTSS exhibits a rising and then a falling trend.

Reversing the RTSS algorithm based using forward Kalman filtering is conducted as follows. The forward Kalman filter solution is first used to obtain the state vector ${\widehat{\mathit{X}}}_{f,k}$ of each sampling point, covariance matrix ${\mathit{P}}_{f,k}$, state transition matrix ${\mathit{\Phi }}_{k,k-1}$, prediction covariance matrix ${\mathit{P}}_{f,k/k-1}$, and navigation information obtained updating the INS strap-down algorithm. After forward filtering ends at the $t_e$ epoch, the RTSS algorithm is implemented in reverse order:

$$\begin{array}{cc}\hfill {\mathit{K}}_{s,k}& ={\mathit{P}}_{f,k}{\mathsf{\Phi }}_{k+1,k}^T{\mathit{P}}_{f,k+1/k}^{-1}\hfill \\ \hfill {\widehat{\mathit{X}}}_{s,k}& ={\widehat{\mathit{X}}}_{f,k}+{\mathit{K}}_{s,k}({\widehat{\mathit{X}}}_{s,k+1}-{\widehat{\mathit{X}}}_{f,k+1/k})\hfill \\ \hfill {\mathit{P}}_{s,k}& ={\mathit{P}}_{f,k}-{\mathit{K}}_{s,k}({\mathit{P}}_{f,k+1/k}-{\mathit{P}}_{s,k+1}){\mathit{K}}_{s,k}^T\hfill \\ \hfill k& =t_e-1,t_e-2,\dots ,2,1,0\hfill \end{array}$$

(15)

where ${\mathit{P}}_{s,t_e}={\mathit{P}}_{f,t_e}$, ${\widehat{\mathit{X}}}_{s,t_e}={\widehat{\mathit{X}}}_{f,t_e}$, and ${\mathit{K}}_{s,k}$ is the smooth gain matrix. Because the tightly coupled system uses a closed-loop feedback mechanism, when the pose information of the INS strap-down algorithm is updated and the random constant bias is corrected, the state vector ${\widehat{\mathit{X}}}_{f,k}$ must be set to zero. Therefore, the state vector ${\widehat{\mathit{X}}}_{f,k+1/k}$ of the one-step prediction is always $\mathbf{0}$. Similarly, after each RTSS, a closed-loop correction is also performed. This closed-loop correction minimizes the linearization error of the system model and ensures that the integration system is not influenced by divergence.

3. Results and Discussion

To study the performance of the DGNSS/INS tightly coupled RTSS post-processing strategy in this study under the GNSS signal lock loss environment, we analyze and experimentally verified the following three aspects from the perspective of fully mining observations: (1) the ability of loosely coupled and tightly coupled to use observations in an environment with insufficient satellites, (2) the effectiveness of the RTSS algorithm during a period of GNSS signal interruption, and (3) the influence of INS-assisted IAR on the position and attitude accuracy in an environment with frequent lock loss. Experiments were conducted near Shandong University of Science and Technology (Qingdao, China). The split closed-loop fiber optic integrated navigation system SPAN-LCI manufactured by NovAtel (Figure 3) was used in the experiment, which contained a GNSS receiver and an IMU-LCI tactical fiber IMU. The primary performance indicators of SPAN-LCI [30] are displayed in

Table 1. SPAN-LCI performance indicators.

	Gyroscope	Accelerometer
Bias stability	$<1.0°/\mathrm{h}$	$<1.0 \mathrm{mg}$
Random walk	$<0.05°/\sqrt{\mathrm{h}}$	$50 \mathsf{\mu }\mathrm{g}/\sqrt{\mathrm{Hz}}$
Scale factor	$100 \mathrm{ppm}$	$250 \mathrm{ppm}$
Sampling rate	$200 \mathrm{Hz}$	$200 \mathrm{Hz}$

. Before the experiment, the lever arm value was calibrated.

The GNSS receiver used was NovAtel ProPark6, which can realize the data collection of GPS, BDS, and GLONASS data. The sampling frequency was 20 Hz, and the GNSS antenna type was NOV703GGG. The reference station was set up at a known point in the marine survey comprehensive experimental field of Shandong University of Science and Technology, with an open surrounding environment. The data collection time of the mobile station was approximately 2.5 h, the baseline length typically exceeded 10 km, and the satellite suffered lock loss in several places. The experimental trajectory on Google Earth is displayed in Figure 4d. Combining the number of common-view satellites of the base station with the rover station in Figure 4a, the position dilution of precision (PDOP) of the rover station satellites (Figure 4b) and the sky plot (Figure 4c), further reflects the experimental environment. In the experiment, according to Table 1, the accelerometer random walk was set to $0.05°/\sqrt{\mathrm{h}}$, and the gyroscope random walk was set to $50 \mathsf{\mu }\mathrm{g}/\sqrt{\mathrm{Hz}}$, which is used to set the system noise matrix. According to the observation accuracy of the GNSS receiver,${\sigma }_{{\mathsf{\Phi }}_{1,0}}$, ${\sigma }_{P_{1,0}}$ and ${\sigma }_{{\dot{\rho }}_{1,0}}$ were set to $0.005 \mathrm{m}$, $0.3 \mathrm{m},$ and $0.1 \mathrm{m},$ respectively. In this study, the NovAtel high-precision integrated navigation post-processing software Inertial Explorer 8.90 was used to perform DGNSS/INS tightly coupled smoothing in both directions, and the navigation result was used as a reference value.

3.1. Simulation Experiment

3.1.1. Test 1

To analyze and verify the performance of the tightly coupled and loosely coupled systems in an environment where the number of available satellites is insufficient, experimental data within 6500 to 6560 epochs from the experimental data are intercepted. In a realistic urban environment, satellites with lower altitudes are more likely to be blocked, and satellites with higher altitudes are more suitable for simulating GNSS satellite loss-of-lock experiments. Therefore, G27, G26, and G16 were selected when simulating an observation environment with only three visible satellites, and G27 and G26 were selected when simulating an observation environment with only two visible satellites. What needs to be emphasized here is that the integer ambiguity of these selected satellites has been correctly fixed in the previous observation epoch.

In Figure 5 and

Table 2. RMSE statistics of 3D position and attitude of the tightly coupled and loosely coupled systems under 1 min of different observable satellite numbers.

		Tightly Coupled (TC)			Loosely Coupled (LC)
Number of Satellites		>4	3	2	>4	3
Position Error (m)	3D	0.028	0.033	0.115	0.036	0.714
Attitude Error (°)	Roll	0.0058	0.0058	0.0058	0.0066	0.0076
	Pitch	0.0114	0.0114	0.0114	0.0126	0.0128
	Heading	0.0238	0.0239	0.0250	0.0235	0.0289

, when there are sufficient satellites available, the accuracy of the position calculated by the loosely and tightly coupled systems are both at the centimeter level, and the attitude accuracy is at the same magnitude. Therefore, when the positioning accuracy of the GNSS is high, there is no significant difference between the solution results of the loosely and tightly coupled systems. When only three satellites are observed, the 3D position error calculated by loosely coupled increases exponentially over time and the maximum 3D position error reaches 1.3 m. Under the premise that the integer ambiguity is fixed correctly, the 3D position error calculated by tightly coupled system can be maintained at the centimeter level, and the root mean square error (RMSE) is 0.033 m. This result indicates that when the number of satellites is less than four, the tightly coupled system can still use the existing satellite observations for measurement updates to limit the error drift. However, the GNSS module in the loosely coupled system cannot use the existing three satellite observations to perform differential positioning solutions, which is equivalent to the GNSS signal being completely interrupted. When there are two observation satellites, the 3D position error calculated by tightly coupled reaches the decimeter level, because there is no redundant observation value.

In tightly coupled system, within a short period, the number of observable satellites substantially influences the position error, but the influence on the attitude error is not apparent. The RMSEs of the roll and pitch are 0.0058° and 0.0114°, respectively. The observability of the heading is poor, and when the number of observable satellites is sufficient, the RMSE is 0.0238°. The heading angle also conforms to the law that there are more available observations and sufficient accuracy. In general, when more observations are used, more accurate pose results can be obtained. Under the premise that the IAR is correct, the tightly coupled system can be measured and updated normally when there are fewer than four satellites, this can guarantee navigation solution accuracy. Therefore, the tightly coupled system is more suitable for observation environments where the satellite is partially locked.

3.1.2. Test 2

To analyze the improvement effect of the RTSS algorithm on the navigation solution, this experiment performed loosely coupled forward filtering processing and loosely coupled reverse RTSS processing on simulated data in Test 1.

The position and attitude error curves in Figure 6 and the error statistics in

Table 3. Position and attitude RMSEs under different processing methods.

		Satellite Number ≥ 4		Satellite Number = 3
		LC-EKF	LC-RTS	LC-EKF	LC-RTS
Position Error (m)	3d	0.036	0.032	0.714	0.035
Attitude Error (°)	Roll	0.0066	0.0063	0.0076	0.0066
	Pitch	0.0126	0.0119	0.0128	0.0118
	Heading	0.0235	0.0121	0.0289	0.0143

demonstrate that the 3D position RMSE calculated by the loosely coupled system is 0.035 m after RTSS, which is 95% more accurate than the position solution before smoothing. Compared with real-time processing, RTSS is expected to bring accuracy improvement because more observations are used. The accuracy of the roll and pitch after RTSS improves slightly, by approximately 13% and 8%, respectively, and the error curve is smoother than that obtained by forward filtering alone. The most apparent improvement in attitude is that of the heading. It can be seen from Figure 6 that the heading after forward filtering exhibits a systematic error. After smoothing, the system error of the heading is significantly reduced by approximately 50%. Comparing the position and attitude error curves before and after smoothing, it can be seen that even if the GNSS signal is not interrupted, RTSS can improve the position and attitude accuracy and also effectively suppress the accumulation of INS estimation errors between measurement updates. Therefore, in the case of a partial or complete lock loss, RTSS as a bridging algorithm can sufficiently utilize the observations during the observation epoch before and after the lock loss, ensuring the accuracy of the navigation solution in the short term when the number of available satellites is insufficient.

3.1.3. Test 3

In the partial satellite lock-out simulation experiments of Tests 1 and 2, satellites whose integer ambiguities have been correctly fixed in the previous epoch are selected. In an actual urban canyon environment, the satellite signal is received intermittently after the satellite signal is completely blocked for a short time, and the integer ambiguity must be fixed again. To analyze the performance of the INS-assisted IAR strategy in this study, we intercepted data within 374,470–374,560 epochs and designed four sets of experimental schemes.

Plan 1: The satellite signal is completely interrupted for 10 s and then continues to observe three satellites for 30 s.

Plan 2: The satellite signal is completely interrupted for 20 s and then continues to observe three satellites for 30 s.

Plan 3: The satellite signal is completely interrupted for 30 s and then continues to observe three satellites for 30 s.

Plan 4: First, the satellite signal is completely out of lock, and only three satellites are observed; finally, over four satellites are observed. The above three observation environments last for 10 s each. The specific experimental situation of plan4 is shown in Figure 7.

Based on the pose errors in Figure 8, Figure 9 and Figure 10 and

Table 4. Maximum error (MAX) and RMSE of position statistics for different plans.

	RMSE (m)				MAX (m)
	INS-Aided		Not INS-Aided		INS-Aided		Not INS-Aided
	EKF	RTS	EKF	RTS	EKF	RTS	EKF	RTS
Plan 1	0.036	0.027	0.085	0.057	0.070	0.032	0.218	0.088
Plan 2	0.062	0.029	0.203	0.032	0.146	0.044	0.373	0.043
Plan 3	0.117	0.030	0.250	0.057	0.242	0.044	0.443	0.095
Plan 4	0.037	0.028	0.098	0.029	0.090	0.040	0.202	0.044

and

Table 5. RMSE of attitude for different plans.

		INS-Aided		Not INS-Aided
		EKF	RTS	EKF	RTS
Roll (°)	Plan 1	0.0059	0.0056	0.0059	0.0056
	Plan 2	0.0062	0.0058	0.0062	0.0059
	Plan 3	0.0061	0.0057	0.0061	0.0057
	Plan 4	0.0058	0.0049	0.0058	0.0050
Pitch (°)	Plan 1	0.0070	0.0067	0.0070	0.0068
	Plan 2	0.0082	0.0079	0.0082	0.0079
	Plan 3	0.0087	0.0082	0.0087	0.0082
	Plan 4	0.0114	0.0110	0.0114	0.0110
Heading (°)	Plan 1	0.0247	0.0099	0.0246	0.0096
	Plan 2	0.0247	0.0100	0.0249	0.0099
	Plan 3	0.0267	0.0111	0.0271	0.0114
	Plan 4	0.0238	0.0089	0.0242	0.0090

, it can be seen that when the position information calculated by the INS to assist the IAR is not used, the maximum 3D position error accumulates to 0.443 m, and the 3D position RMSE is 0.25 m in Plan 3. When using the prior position, the maximum position error is 0.242 m, and the RMSE is 0.117 m. When there is no INS-aided IAR, the ambiguity cannot be fixed. Therefore, even if three satellites are continuously observed after the satellite is interrupted, the observation values of these satellites cannot be used to update the measurement, resulting in the continuous accumulation of position errors. After using the prior position provided by the INS to assist the IAR, the integer ambiguities of the three satellites are quickly fixed. The observations of these three satellites are effectively used for measurement updates, and the increase in position error is significantly reduced. The accuracy of the attitude primarily depends on the inertial navigation, and the tactical-level inertial navigation is used in this experiment. Therefore, the influence of the INS-aided IAR on the attitude accuracy is not apparent. After RTSS, the 3D position RMSE within 1 min of the satellite signal interruption can be controlled within 0.06 m, and the attitude is also significantly improved. Specifically, the roll, pitch, and heading are improved by approximately 15%, 5%, and 60%, respectively.

In addition, after the satellite signal is interrupted for 10 s to 20 s, the position accuracy obtained after the INS-aided integer ambiguity is quickly fixed and remains at the centimeter level. After the satellite signal was interrupted for 30 s, the position error estimated by the INS reached the decimeter level. Even if the integer ambiguity was quickly fixed, and the measurement was updated, the position error still accumulated, and the final position error was at the decimeter level. The effectiveness of the INS-aided IAR therefore depends on the position accuracy of the INS estimation.

Based on Figure 11 and Table 4 and Table 5, it can be seen that in an environment with fewer than four satellites, intermittently observed after the satellite signal is interrupted, the position error calculated by the tightly coupled without INS-aided IAR accumulates. In contrast, for INS-assisted IAR, once there are observations, the position error is greatly reduced. Similarly, after smoothing, the accuracy of the position and posture also improved. Comparing the position error curves in Figure 8, Figure 9, Figure 10 and Figure 11 and the position error statistics in Table 4 and Table 5 indicates that the INS-aided IAR and smoothing post-processing tightly coupled strategy in this study does not require additional information, thereby fully utilizing the existing observations to maximize the limitation of the accumulation of errors and improve navigation solution accuracy.

3.2. Field Experiment in Urban Environment

To verify the performance of the tightly coupled post-processing solution strategy in this study in an urban environment, data from the final hour of the Qingdao medium and a short baseline experiment were intercepted. Then, the tightly coupled solution without INS-aided IAR (Method 1), the tightly coupled solution with INS-aided IAR (Method 2), and the tightly coupled solution strategy in this study (Method 3) were used to process the data. Figure 12 displays the experimental observation environment. The total number of common-view satellites of the base and rover stations stabilized at approximately 16 in most periods, and the satellites frequently suffered lock loss during certain periods. The experimental trajectories at jumping points A, B, and C marked in Figure 12, Figure 13, Figure 14 and Figure 15 correspond to the three experimental trajectories in Figure 13. At jumping point A, there are tall buildings, and the satellite changes quickly. At this time, the observation structure is poor, and there are many integer ambiguity reinitializations. At jumping point B, the test vehicle passes a viaduct, causing the GNSS signal to suffer complete lock loss. When relocking the satellite signal, the satellite geometry is poor, and the observation accuracy is low. At jumping point C, the test vehicle passes through the south gate of Shandong University of Science and Technology, and the upper beam of the door is wide, and hence shields the GNSS signal.

In Figure 14 and Figure 15 and

Table 6. Position and attitude error statistics obtained using three different methods.

	RMSE(m)			MAX(m)
	Not INS-Aided	INS-Aided	INS-Aided RTS	Not INS-Aided	INS-Aided	INS-Aided RTS
Position	0.063	0.054	0.049	1.159	0.347	0.253
Roll	0.0085	0.0085	0.0072	0.0546	0.0543	0.0539
Pitch	0.0115	0.0115	0.0113	0.0883	0.0886	0.0897
Heading	0.0510	0.0490	0.0425	0.1684	0.1541	0.1679

, it can be seen that when the observation environment and number of visible satellites are sufficient, the position errors acquired by the three solutions are maintained at the centimeter level. Of these, the position curve obtained using Method 3 is smoother, and the attitude accuracy is improved, the heading accuracy improvement is the greatest. At jumping points A, B, and C, the position and attitude accuracy are affected to varying degrees. The 3D position RMSE calculated using Method 1, 2, and 3, was 0.063 m, 0.054 m, and 0.049 m, respectively; meanwhile, the maximum error is 1.159 m, 0.374m, and 0.253m, respectively. Compared with the other two methods, the accuracy of the attitude calculated using Method 3, especially that of the heading, is significantly improved. Therefore, the solution strategy in this study effectively limits the increase in position and attitude errors when the satellite is out of lock by sufficiently utilizing the observations before and after the satellite signal occlusion period. This is an effective means to improve the positioning and attitude accuracy of the GNSS/INS integrated system in an urban environment.

4. Conclusions

In an actual urban environment, due to the influence of tall buildings, trees, and other occlusion factors, the satellite signal frequently suffers lock loss, causing a reduction in the accuracy of the GNSS/INS integrated system, rendering it unable to meet application requirements. This study implemented an INS-aided IAR DGNSS-INS tightly coupled smoothing post-processing strategy. This strategy uses INS-aided IAR in the ambiguity initialization stage after the GNSS signal is interrupted to achieve fast IAR, followed by RTSS to determine the navigation solution during lock-out period processing. Moreover, this solution strategy does not require additional observations but improves the accuracy of the integrated system by fully utilizing existing observations. The primary results obtained through experimental analysis in this study are as follows.

(1) The performances of the loosely and tightly coupled systems and the improvement effect of the RTSS algorithm on the navigation solution were analyzed through simulation experiments where available satellites is insufficient. The experimental results demonstrate that when the number of observable satellites was sufficient, there is no significant difference in the positioning and attitude determination performance of the loosely and tightly coupled satellites. When the number of observation satellites was less than four, compared to the loosely coupled system, the tightly coupled system fully utilized the observations of satellites to update the measurements and constrain the drift error of the inertial device. In addition, during the lock-out period, the accuracy of the loosely coupled system after RTSS was significantly improved. The position accuracy was improved by approximately 95%. The degree of improvement in the roll and pitch was approximately 10%, and that in the heading was approximately 50%.

(2) The performances of the tightly coupled system without INS-aided IAR, the tightly coupled system with INS-aided IAR, and the solution strategy were compared and analyzed through simulation and field experiments. The results demonstrate that the strategy proposed in this study can use the position provided by the INS to quickly fix integer ambiguity of the three satellites after the satellite signal is completely interrupted, allowing the observations of the three satellites to become available observations. The strategy in this study effectively limited the growth of the pose error, and the position and attitude accuracy of the navigation solution was further improved after smoothing.

In this study, the performance of the tightly coupled DGNSS/INS under the condition of short-term satellite signal multi-occlusion was analyzed, and solutions were provided. However, it was difficult to obtain centimeter-level positioning accuracy with GNSS signals that are completely out of lock for a long period. Future work should focus on special environments where GNSS is completely out of lock for a long period, such as garages and tunnels, to further enrich quality control methods and to obtain more reliable and accurate navigation information.