## Methodology
### Mathematic Relationship of Residuals on Different Frequencies  
Since the coordinates are fixed and the zenith tropospheric delays and ambiguities are treated as constants during a period of time, the residuals, which are primarily the multipath errors, are related to the time-varying parameters, that is, receiver clock and ionospheric delays. Therefore, Equation (2) can be rearranged as follows:  
\[\textbf{V=AX-L}, \tag{3}\]  
\[\textbf{A}=\left[ \begin{array}{cc}\textbf{E}&-\textbf{I}\\ \textbf{E}&-\gamma_{2}\cdot\textbf{I}\end{array} \right]\textbf{X}=\left[ \begin{array}{c}\overline{\textbf{i}_{r}}\\ \overline{\textbf{i}_{1}^{*}}\end{array} \right]\textbf{L}=\left[ \begin{array}[ ]{c}\overline{\textbf{i}_{1}^{*}}\\ \overline{\textbf{i}_{2}^{*}}\end{array} \right]\textbf{V}=\left[ \begin{array}[ ]{c}\textbf{v}_{1}^{*}\\ \textbf{v}_{2}^{*}\end{array} \right]\textbf{,} \tag{4}\]  
where **A** denotes the design matrix; **E** denotes a column vector of n-dimension with value one; **I** is an identity matrix of n-dimension; **X** is the parameter vector of receiver clock and ionospheric delays; **L** refers to the vector of unmodeled errors, and each element of \(\overline{\textbf{i}_{j}}\) can be expressed as \(l_{j}^{s}+\textbf{e}_{r}^{s}\cdot\textbf{r}_{r}-T_{r}^{s}-\lambda_{j}^{s} \overline{N_{j}^{*}}\) and **V** represents the vector of residuals. According to the least-squares criterion, the estimated parameters read as follows:  
\begin{table}
\begin{tabular}{c c} \hline
**Item** & **Processing Information** \\ \hline Estimator & Least-squares estimator for generating phase residuals \\ Observations & Raw pseudo-range and carrier-phase observations from GPS, GLONASS, and BDS \\ Sampling rate & 1 s \\ Elevation cutoff & 7\({}^{\circ}\) \\ \multirow{4}{*}{Weighting scheme} & Elevation-dependent weight; 3 dm and 3 mm for GPS pseudo-range and carrier-phase; 4.5 dm and 3 mm for GLONASS pseudo-range and carrier-phase; 9.0 dm for BDS pseudo-range; and 5 mm and 15 mm for IGSO/MEO and GEO \\ \multirow{4}{*}{Satellite orbit/clock} & \multirow{2}{*}{GBM final precise orbit/clock products generated by GFZ ([PERSON] et al. [29])} \\  & The zenith hydrostatic delay corrected by [PERSON]’s model (30); the zenith wet \\ \multirow{2}{*}{Tropospheric delay} & delay and the horizontal gradients estimated as piecewise constants every hour and \\  & six hours, respectively; Global Mapping Function (GMF) applied \\ \multirow{2}{*}{Ionospheric delay & Estimated epoch by epoch \\ Satellite/Receiver & GPS/GLONASS: Corrected both at satellite and receiver \\ antenna phase center & BDS: PCO and PCV corrected at satellite, while replaced by GPS at receiver \\ \multirow{2}{*}{Phase-windup effect} & \multirow{2}{*}{ISB estimated as white noise, GPS as reference, whereas IFB estimated as constant for \\  & a whole day \\ \multirow{2}{*}{Station displacement & Solid Earth tide, pole tide, ocean tide loading, IERS Convention 2010 \\ Receiver coordinate & Estimated as constants for daily solution \\ \multirow{2}{*}{Receiver clock & Estimated as white noise \\ \multirow{2}{*}{Ambiguity}} & \multirow{2}{*}{Estimated as constant for each arc: float value} \\ \hline GPS, Global Positioning System; GLONA  
\[\overset{\wedge}{\textbf{X}}=(\textbf{A}^{T}\textbf{A})^{-1}\textbf{A}^{T}\textbf{L}. \tag{5}\]  
By substituting Equation (5) into (3), the residuals are derived:  
\[\textbf{V}=(\textbf{A}(\textbf{A}^{\textbf{T}}\textbf{A})^{-1}\textbf{A}^{T}- \textbf{I})\textbf{L}. \tag{6}\]  
The term \(\textbf{A}^{T}\textbf{A}\) can also be simplified as follows:  
\[\textbf{A}^{\textbf{T}}\textbf{A}=\left[ \begin{array}{cc}\textbf{a}_{11}& \textbf{a}_{12}\\ \textbf{a}_{21}&\textbf{a}_{22}\end{array} \right], \tag{7}\]  
where \(\textbf{a}_{11}=2n\), \(\textbf{a}_{12}=\textbf{a}_{21}^{T}=-(1+\gamma_{2})\textbf{E}^{T}\), and \(\textbf{a}_{22}=(1+\gamma_{j\delta}^{2})\textbf{I}\).