## 2 Multi-Constellation PPP Model
### Multi-Constellation PPP  
Ionosphere-free combination observations are normally used in PPP to remove the first-order ionospheric delay. Their code and carrier phase observation can be expressed as:  
\[\left\{ \begin{array}{l}\frac{f_{1}^{2}}{f_{1}^{2}-f_{2}^{2}}\rho_{1}-\frac{f _{2}^{2}}{f_{1}^{2}-f_{2}^{2}}\rho_{2}=R+c(dt-dT)+T_{trop}+\varepsilon_{P}\\ \frac{f_{1}^{2}}{f_{1}^{2}-f_{2}^{2}}\lambda_{1}\phi_{1}-\frac{f_{2}^{2}}{f_{1 }^{2}-f_{2}^{2}}\lambda_{2}\phi_{2}=R+c(dt-dT)+\lambda_{IF}N_{IF}+T_{trop}+ \varepsilon_{\phi}\end{array} \right. \tag{1}\]  
where \(f_{i}\) (\(i=1,2\)) are carrier-phase frequencies in Hertz; \(\rho_{i}\) is the code observation at the \(i\)th frequency in metres; \(\phi_{i}\) is the carrier phase observation in cycles; \(dt\) is the receiver clock offset in seconds; \(dT\) is the satellite clock offset in seconds; \(c\) is the speed of light in metres per second; \(T_{trop}\) is the tropospheric delay in metres; \(N_{IF}\) are the parameters of the float ambiguity after being redefined; \(R\) is the geometric range in metres; \(\varepsilon_{\rho}\) and \(\varepsilon_{\phi}\) are the code and carrier phase observation noise, respectively, including the multipath in metres; and \(\lambda_{IF}\) is the wavelength after being redefined. The linearization formula of Equation (1) is:  
\[ \begin{array}{l}\rho_{IF}=R+c(\delta t-\delta T)+T_{trop}+\varepsilon_{p}\\ \qquad\qquad=\frac{x^{i}-x_{p}}{r}dx_{p}+\frac{v^{i}-y_{p}}{r}dy_{p}+\frac{z^{ i}-z_{p}}{r}dz_{p}+cdt+MFT_{ztd}-cdT+R+\varepsilon_{p}\end{array} \]  
\[ \begin{array}{l}\lambda_{IF}\Phi_{IF}=R+c(\delta t-\delta T)+\lambda_{IF}N_{ IF}-T_{trop}+\varepsilon_{\phi}\\ \qquad\qquad=\frac{x^{i}-x_{p}}{r}dx_{p}+\frac{y^{i}-y_{p}}{r}dy_{p}+\frac{z^{ i}-z_{p}}{r}dz_{p}+cdt+MFT_{ztd}\ \ +\lambda_{IF}N_{IF}-cdT+R+\varepsilon_{\phi}\end{array}  \tag{2}\]  
In Equation (2), \(r\) is the geometric range computed using the linearization point. \(M\)F is the mapping function; \(T_{ztd}\) is the zenith tropospheric wet delay; \(\rho_{IF}\) are the code observations; and \(\Phi_{IF}\) are carrier observations. \(x^{i}\), \(y^{i}\), \(z^{i}\) are the position of satellites, and \(x_{p}\), \(y_{p}\), \(z_{p}\) are the coordinates of the stations. Then, the least squares residuals can be written as:  
\[\left( \begin{array}{c}V_{\rho IF}\\ V_{\phi IF}\end{array} \right)=\left( \begin{array}{ccc}B&A&M&O\\ B&A&M&\lambda\end{array} \right)\left( \begin{array}{c}dX_{p}\\ dt\\ T_{ztd}\\ N_{IF}\end{array} \right)-\left( \begin{array}{c}L_{\rho IF}\\ L_{\phi IF}\end{array} \right) \tag{3}\]  
where the parameters \(dX_{p}\), \(dt\), \(T_{ztd}\) and \(N_{IF}\) are the receiver position, receiver clock offset, zenith tropospheric wet delay and ambiguity of the ionosphere-free combination, respectively; \(B\), \(A\), \(M\), and\(\lambda\) are the corresponding coefficient matrices, respectively. \(V_{\rho IF}\) and \(V_{\phi IF}\) are the observation residuals. For multi-constellation kinematic PPP, the code hardware delay biases should be absorbed into the receiver clock and system time difference, and the carrier-phase biases related to the frequencyshould be absorbed into the ambiguity. In addition, the system time difference and receiver clock offset were linearly correlated. Generally, instead of estimating a receiver clock parameter for each satellite system observation, we always introduce a system time difference parameter to reflect the difference between the different system times. The multi-constellation kinematic PPP observation equations can be expressed as follows:  
\[\left( \begin{array}{c}V_{\rho IF}\\ V_{\phi IF}\\ V_{\rho IF_{1}}\\ V_{\phi IF_{1}}\\ \end{array} \right)=\left( \begin{array}{cccccc}B&A&O&M&O&O\\ B&A&O&M&\lambda&O\\ B_{1}&O&A&M&O&O\\ B_{1}&O&A&M&O&\lambda_{other}\\ \end{array} \right)\left( \begin{array}{c}dX_{p}\\ dt_{gps}\\ dt_{other}\\ T_{ztd}\\ N_{IF}\\ \end{array} \right)-\left( \begin{array}{c}L_{\rho IF}\\ L_{\phi IF}\\ L_{\phi IF_{1}}\\ \end{array}