## 3 Multi-GNSS PPP/INS Mathematical Model
### MEMS-INS Model  
Commonly used coordinate systems in inertial navigation mainly include inertial frame (\(i\)-frame), Earth-centered Earth-fixed frame (\(e\)-frame), body frame (\(b\)-frame), and navigation frame (\(n\)-frame). With the dead reckoning algorithm, the dynamic equation of INS in \(e\)-frame can be summarized as [37]:  
\[\begin{pmatrix}\dot{r}_{\text{INS}}^{e}\\ \dot{v}_{\text{INS}}^{e}\\ C_{b}^{e}\end{pmatrix}=\begin{pmatrix}\dot{v}_{\text{INS}}^{e}\\ C_{b}^{e}f_{ib}^{b}-2\Omega_{ie}^{e}v_{\text{INS}}^{e}+g^{e}\\ C_{b}^{e}\Omega_{eb}^{b}\end{pmatrix} \tag{12}\]  
where \(r_{\text{INS}}^{e}\) and \(v_{\text{INS}}^{e}\) denote the position and velocity vectors of the IMU in the \(e\)-frame, respectively; \(C_{b}^{e}\) denotes the direction cosine matrix from \(b\)-frame to \(e\)-frame, and \(f_{ib}^{b}\) denotes the specific force vector of \(e\)-frame with respect to \(i\)-frame in the \(e\)-frame; \(\Omega_{ie}^{e}\) denotes the skew symmetric matrix of the earth rotation vector of \(\omega_{ie}^{e}\) relative to \(i\)-frame, expressed in \(e\)-frame; \(g^{e}\) denotes the gravity vector in the \(e\)-frame; \(\Omega_{eb}^{b}\) denotes the skew symmetric matrix of \(\omega_{eb}^{b}\), and \(\omega_{eb}^{b}\) can be written as:  
\[\omega_{eb}^{b}=\omega_{ei}^{b}+\omega_{ib}^{b}=\omega_{ib}^{b}-\omega_{ie}^{b }=\omega_{ib}^{b}-C_{e}^{e}\omega_{ie}^{e} \tag{13}\]  
where \(\omega_{ib}^{b}\) denotes gyroscope angular velocity, \(C_{e}^{b}\) denotes the inverse matrix of \(C_{b}^{e}\).  
By using the dynamic equation, the error state model of INS can be expressed as [38]:  
\[\begin{pmatrix}\dot{\delta r}_{\text{INS}}^{e}\\ \dot{\delta v}_{\text{INS}}^{e}\\ \dot{\vartheta}^{e}\end{pmatrix}=\begin{pmatrix}\delta v_{\text{INS}}^{e}\\ -2\Omega_{ie}^{e}\delta v_{\text{INS}}^{e}+\left[\left(C_{b}^{e}f_{ib}^{b} \right)\times\right]q^{e}-C_{b}^{e}\delta f_{ib}^{b}\\ -\Omega_{ie}^{e}\varphi^{e}+C_{b}^{e}\delta\omega_{ib}^{b}\end{pmatrix} \tag{14}\]  
where \(\left[\left(C_{b}^{e}f_{ib}^{b}\right)\times\right]\) denotes the skew symmetric matrix of \(\left(C_{b}^{e}f_{ib}^{b}\right)\). The errors of specific force and gyroscope angular velocity can be written as follows:  
\[\delta f_{ib}^{b}=\delta b_{a}+\varepsilon_{a} \tag{15}\]  
\[\delta\omega_{ib}^{b}=\delta b_{g}+\varepsilon_{g} \tag{16}\]  
Then, the linearized state equation of INS can be expressed as:  
\[\begin{bmatrix}\dot{\delta r}_{\text{INS}}^{e}\\ \dot{\delta v}_{\text{INS}}^{e}\\ \dot{\vartheta}^{e}\\ \dot{\delta b}_{a}\\ \dot{\delta b}_{g}\end{bmatrix}=\begin{bmatrix}0&I&0&0&0\\ 0&-2\Omega_{ie}^{e}&\left[\left(C_{b}^{e}f_{ib}^{b}\right)\times\right]&-C_{b}^ {e}&0\\ 0&0&-\Omega_{ie}^{e}&0&C_{b}^{e}\\ 0&0&0&0&0\\ 0&0&0&0&0\end{bmatrix}\begin{bmatrix}\delta r_{\text{INS}}^{e}\\ \delta v_{\text{INS}}^{e}\\ \delta v_{\text{INS}}^{e}\\ \delta v_{\text{INS}}^{e}\\ \delta v_{\text{INS}}^{e}\end{bmatrix}+\begin{bmatrix}0&0\\ -C_{b}^{e}&0\\ 0&C_{b}^{e}\\ 0&0\\ 0&0\end{bmatrix}\begin{bmatrix}\varepsilon_{a}\\ \varepsilon_{b}\end{bmatrix} \tag{17}\]  
where \(\delta r_{\text{INS}}^{e}\), \(\delta v_{\text{INS}}^{e}\), and \(\varphi^{e}\) denote the position, velocity, and attitude error vector in \(e\)-frame; \(\delta b_{a}\) and \(\delta b_{g}\) denote the bias errors of accelerometer and gyroscope, respectively; \(\varepsilon_{a}\) and \(\varepsilon_{g}\) denote the processing noise vectors of accelerometer and gyroscope biases, respectively.