## 2 Methodology  
The linearization of the pseudorange and carrier phase observation equations containing the parameters to be sought can be expressed as [24]:  
\[p_{r,1}^{s}=\mu_{r}^{s}\mathbf{x}+dt_{r}-dt^{s}+V_{r,trop}^{s}+I_{r,1}^{s}+(d_{r,1}-d _{1}^{s}) \tag{1}\]  
\[I_{r,1}^{s}=\mu_{r}^{s}\mathbf{x}+dt_{r}-dt^{s}+V_{r,trop}^{s}-I_{r,1}^{s}+\lambda_{1 }^{s}(N_{r,1}^{s}+b_{r,1}-b_{1}^{s}) \tag{2}\]  
where \(s\) indicates satellite; \(r\) represents the receiver; \(p_{r,1}^{s}\) and \(I_{r,1}^{s}\) denote observed minus computed (OMC) values of the pseudorange and carrier phase observables on the first frequency, respectively; \(\mu_{r}^{s}\) represents the unit vector from the satellite to the receiver; \(x\) is the vector of the receiver position increments relative to the a priori position; \(dt_{r}\) and \(dt^{s}\) are the clock bias of the receiver and satellite, respectively; \(V_{r,trop}^{s}\) is the troposphere delay; \(I_{r,1}^{s}\) is the slant ionospheric delay on the first frequency \(t_{1};d_{r,1}\) is the frequency-dependent receiver uncalibrated code delay (UCD) with respect to satellite \(s\); \(d_{1}^{s}\) is the frequency-dependent satellite UCD; \(N_{r,1}^{s}\) is the integer phase ambiguity; \(b_{r,1}\) and \(b_{1}^{s}\) are the uncalibrated phase delays (UPDs) for the receiver and satellite, frequency-dependent. In order to speed up the convergence, the GIMs product provided by IGS can be used as a pseudo-observation to constrain the ionospheric delay. The mathematical model can be written as:  
\[I_{r,GIM}^{s}=I_{r,1,GIM}^{s}+\varepsilon_{GIM},\ \ \ \ \ \varepsilon_{GIM}\in N(0, \sigma_{GIM}^{2}) \tag{3}\]  
where \(I_{r,1,GIM}^{s}\) is calculated from the GIMs product value and the ionosphere pierce point (IPP) angle [25]; \(\sigma_{GIM}^{2}\) is a priori variance for the noise \(\varepsilon_{GIM}\) of the GIMs product. Then the ionospheric delay correction (\(v_{I}\)) equation is:  
\[v_{I}=I_{r,GIM}^{s}-I_{r}^{s} \tag{4}\]  
IGS also provides tropospheric delay products (i.e., zenith tropospheric delay, ZTD), which can be used as a constraint, and expressed as:  
\[V_{r,ZTD}^{s}=V_{r,1,ZTD}^{s}+\varepsilon_{ZTD},\ \ \ \ \ \varepsilon_{ZTD}\in N(0, \sigma_{ZTD}^{2}) \tag{5}\]  
where \(V_{r,1,ZTD}^{s}\) can be expressed as a ZTD value and mapping function [26], \(\sigma_{ZTD}^{2}\) is a priori variance for its noise \(\varepsilon_{ZTD}\) of the ZTD product. Then the tropospheric delay correction (\(v_{T}\)) equation is:  
\[v_{T}=V_{r,ZTD}^{s}-V_{r,trop}^{s} \tag{6}\]  
The Kalman filtering method [27] is a parameter estimation method proposed by [PERSON] around 1960. The method predicts the state of the next moment by the noise of the measurement and state vector at the priori moment. The true value of the current moment is estimated more accurately by observing the noise and the estimated noise. For the GNSS raw observations (pseudorange and carrier phase), the observation model can be expressed as:  
\[\mathbf{Z}_{k+1}=\mathbf{H}_{k+1}\mathbf{X}_{k+1}+v_{k+1},\ \ \ \ \ v_{k+1}\in N(0,\mathbf{R}_{k+1}) \tag{7}\]  
where \(\mathbf{Z}_{k+1}\) is the vector of measurements at epoch \(k+1\) of the GNSS raw observations; \(\mathbf{H}_{k+1}\) is the design matrix; vector \(\mathbf{X}_{k+1}\) represents the parameters; \(v_{k+1}\) represents the observation noise with the priori covariance matrix \(\mathbf{R}_{k+1}\). For single-frequency observations, the \(\mathbf{Z}_{k+1}\) can be written as:  
\[\mathbf{Z}_{k+1}=\left[ \begin{array}{c}p_{r,1}^{s}\\ I_{r,1}^{s}\end{array} \right]+\left[ \begin{array}{c}v_{p_{1}}\\ v_{l_{1}}\end{array} \right] \tag{8}\]where \(v_{p_{1}}\) represents the pseudorange observation noise; \(v_{l_{1}}\) represents the carrier phase observation noise. According to the mathematical functions above, the state parameter vector \(\mathbf{X}_{k+1}\) can be written as:  
\[\mathbf{X}_{k+1}=\left[ \begin{array}{cc}\mathbf{x}&dt_{r}&V^{s}_{r,trop}&\mathbf{I}^{s}_{r }&\mathbf{N}^{s}_{r}\end{array}