Introduction

Precise point positioning (PPP) is an absolute positioning technique which has advantages in terms of high computational efficiency and low cost (Malys and Jensen 1990; Zumberge et al. 1997; Kouba and Héroux 2001). A higher positioning performance can be obtained with the carrier phase ambiguity fixed to its correct integer, compared with the ambiguity-float PPP (Collins et al. 2008; Ge et al. 2008; Laurichesse et al. 2009). With the rapid development of the new generations of GNSS, most of the satellites are transmitting signals on three or even more frequencies. More attention has been paid to the benefit of the third frequency on ambiguity-fixed PPP.

Since 2009, GPS has started to replace the old satellites with the new Block IIF, which are transmitting the third civil signal L5 (1176.45 MHz) in addition to the existing L1 (1575.42 MHz) and L2 (1227.60 MHz) signals. Currently, there are 12 Block IIF satellites in operation. As a latecomer navigation satellite system, BeiDou and Galileo are both the constellations of which all operating satellites transmit triple-frequency signals. BeiDou and Galileo are currently at the rapid deployment stage. The frequencies are 1561.098 MHz (B1), 1207.14 MHz (B2), and 1268.52 MHz (B3) for BeiDou-2 signals, while 1575.42 MHz (E1), 1176.45 MHz (E5a), and 1207.14 MHz (E5b) for Galileo signals. The B1 band is close to the GPS L1 frequency, and the B2 frequency is identical to Galileo E5b. The BeiDou first launched its regional navigation service in 2012. There are 5 GEO, 7 IGSO, and 3 MEO satellites in operation for BeiDou-2. Galileo provides an open signal in the E1 band and a wideband signal covering the E5 a&b band. Galileo began to launch the Full Operational Capability (FOC) satellites since 2014 (Montenbruck et al. 2017; Sonica et al. 2017). As of June 2019, the operating Galileo satellite number has been increased to 24. The full constellation with 30 satellites is expected to be realized by 2020.

However, the inter-frequency clock bias (IFCB), that is, the difference between the current clock products computed with L1/L2 and the satellite clocks computed with L1/L5, was noticed for the new GPS L5 signals. The IFCB can be larger than 1 dm so that the L1/L2 precise clock products cannot be directly used for L1/L5 PPP (Montenbruck et al. 2012). Consequently, the triple-frequency GPS PPP has been severely affected. Until now, only one GPS triple-frequency PPP AR study has been reported by Geng and Bock (2013) using the simulated observations not affected by IFCB. An extra-wide-lane (EWL), wide-lane (WL), and narrow-lane (NL) sequential fixing strategy is employed in their study. Some studies successfully estimated the GPS IFCB and applied them in triple-frequency ambiguity-float PPP (Pan et al. 2017b, 2018). These IFCB should be further employed for the triple-frequency ambiguity-fixed PPP with real data. Based on real tracking BeiDou data, Pan et al. (2017a) considered IFCB corrections for BeiDou-2 ambiguity-float PPP and achieved an improvement of position estimation by about 10–20%. Gu et al. (2015) conducted triple-frequency ambiguity-fixed PPP with sequential EWL-WL-L1 fixing strategy, while Li et al. (2018) realized a triple-frequency UPD estimation and PPP AR method based on a raw PPP model. Similarly, the efficiency of Galileo and BeiDou triple-frequency PPP AR based on raw PPP model was verified in Xiao et al. (2019), while that based on EWL-WL-L1 strategy were investigated in Li et al. (2019). In these studies, the BeiDou IFCBs are ignored in PPP for simplicity.

To simultaneously fix the triple-frequency PPP ambiguities with the combined GPS, BeiDou-2, and Galileo for a better positioning performance, the time-varying IFCB should be estimated in advance to correct the third-frequency GPS and BeiDou-2 carrier phase. Then, the high-precision triple-frequency ambiguity parameter can be obtained for UPD estimation and PPP AR. To achieve this goal, the raw PPP model is further improved in this study to process the multi-frequency and multi-GNSS observations with the IFCB considered. At the server, the IFCB would be provided together with the UPD at each frequency, while at the user, a modified reference satellite selection strategy is employed to involve the GPS L5 ambiguity in PPP AR. New results on the three systems UPD estimation and PPP AR are provided.

We first formulate the observation equations and stochastic model for PPP based on raw measurements, address the methods applied to IFCB and triple-frequency UPD estimation, and the modified strategy for reference satellite selection in GPS triple-frequency PPP AR. Then, the PPP AR processing strategy and experiment cases are described, followed by an evaluation of the performance of GNSS triple-frequency PPP AR with the comparison with the GPS and GNSS dual-frequency PPP AR. The final section presents conclusions and perspectives.

Method for triple-frequency PPP processing

GPS, Galileo, and BDS satellites transmit signals in the L band using Code Division Multiple Access (CDMA) technology. Their raw observations of pseudorange \(P\) and carrier phase \(L\) can be expressed as:

$$P_{r,i}^{s} = \rho_{r}^{s} + c\left( {t_{r} - t^{s} } \right) + \gamma_{i} \cdot I_{r,1}^{s} + m_{r}^{s} \cdot zwd_{r}^{{}} + b_{r,i} - b_{i}^{s} + e_{r,i}^{s}$$
(1)
$$L_{r,i}^{s} = \rho_{r}^{s} + c\left( {t_{r} - t^{s} } \right) - \gamma_{i} \cdot I_{r,1}^{s} + m_{r}^{s} \cdot zwd_{r}^{{}} + \lambda_{i} \cdot \left( {N_{r,i}^{s} + B_{r,i} - B_{i}^{s} } \right) + IFCB^{s} \cdot \delta_{i3} + \varepsilon_{r,i}^{s}$$
(2)

where \(s\) is the satellite and \(r\) is the receiver; \(i = 1,2,3\) means the carrier frequency, \(\rho\) is the geometric distance between the phase centers of the satellite and receiver antennas, \(c\) is the speed of light in vacuum, \(t_{r}^{{}}\) and \(t^{s}\) are the clock errors for the receiver and satellite clock, respectively. Further, \(\gamma_{i} = {\raise0.7ex\hbox{${f_{1}^{2} }$} \!\mathord{\left/ {\vphantom {{f_{1}^{2} } {f_{i}^{2} }}}\right.\kern-0pt} \!\lower0.7ex\hbox{${f_{i}^{2} }$}}\) with \(f\) is the frequency, \(I_{r,1}^{s}\) is the slant ionospheric delay on the first frequency, \(zwd_{r}^{{}}\) is the wet troposphere delay (ZWD) at zenith, \(b_{r,i}^{{}}\) and \(b_{i}^{s}\) are the code hardware delay from receiver and satellite, \(e\) is the pseudorange measurement noise, \(\lambda_{i}^{{}}\) is the carrier phase wavelength in meters, \(N_{r,i}^{s}\) is the integer carrier phase ambiguity, and \(B_{r,i}^{{}}\) and \(B_{i}^{s}\) are the receiver-dependent and satellite-dependent carrier phase hardware delay, respectively. \(\delta_{i3}^{{}}\) is the Kronecker delta, which is equal to 1 only if \(i = 3\), and to 0 otherwise, and \({\text{IFCB}}^{s}\) is the inter-frequency clock bias in the third-frequency carrier phase. The \({\text{IFCB}}^{s}\) is absent for the L1/L2 but exists in the third-frequency carrier phase. This is because the precise satellite orbit and clock products from IGS are generated with the ionospheric-free combined code and carrier phase observations on the L1/L2 signals. The time-varying phase biases will be grouped with satellite clock offset parameters. Therefore, the satellite clock estimates are different when using different ionospheric-free observations formed by the first + second frequency and the first + third frequency, respectively. This disagreement is the IFCB item in our study. It is demonstrated to be independent on GNSS receiver and antenna types (Li et al. 2012; Pan et al. 2017b). \(\varepsilon\) is the measurement noise of the carrier phase. Other error items whose magnitude is over 1 cm, such as the phase center offset (PCO) and variation (PCV), dry slant troposphere delay, phase wind-up (Wu et al. 1993), and relativistic effect, should be precisely corrected to obtain high-accuracy PPP solution. For simplicity, they are assumed to be precisely corrected with their corresponding models (Petit and Luzum 2010) and not listed in the equations. It should be noted that in case of lacking the precise PCO/PCV information of the third frequency, we use the PCO/PCV corrections of the second frequency instead due to the adjacent frequency.

According to the current International GNSS Service (IGS) analysis convention, after applying precise satellite clock corrections, eventually, the observation equations can be written as:

$$\left\{ {\begin{array}{*{20}l} {P_{r,1}^{s} = \rho_{r}^{s} + c \cdot \bar{t}_{r} + \gamma_{1} \cdot \bar{I}_{r,1}^{s} + m_{r}^{s} \cdot zwd_{r}^{{}} + e_{r,1}^{s} } \hfill \\ {P_{r,2}^{s} = \rho_{r}^{s} + c \cdot \bar{t}_{r} + \gamma_{2} \cdot \bar{I}_{r,1}^{s} + m_{r}^{s} \cdot zwd_{r}^{{}} + e_{r,2}^{s} } \hfill \\ {P_{r,3}^{s} = \rho_{r}^{s} + c \cdot \bar{t}_{r} + \gamma_{3} \cdot \bar{I}_{r,1}^{s} + m_{r}^{s} \cdot zwd_{r}^{{}} + {\text{IFB}}_{r}^{s} + e_{r,3}^{s} } \hfill \\ {L_{r,1}^{s} = \rho_{r}^{s} + c \cdot \bar{t}_{r} - \gamma_{1} \cdot \bar{I}_{r,1}^{s} + m_{r}^{s} \cdot zwd_{r}^{{}} + \lambda_{1} \cdot \bar{N}_{r,1}^{s} + \varepsilon_{r,1}^{s} } \hfill \\ {L_{r,2}^{s} = \rho_{r}^{s} + c \cdot \bar{t}_{r} - \gamma_{2} \cdot \bar{I}_{r,1}^{s} + m_{r}^{s} \cdot zwd_{r}^{{}} + \lambda_{2} \cdot \bar{N}_{r,2}^{s} + \varepsilon_{r,2}^{s} } \hfill \\ {L_{r,3}^{s} = \rho_{r}^{s} + c \cdot \bar{t}_{r} - \gamma_{3} \cdot \bar{I}_{r,1}^{s} + m_{r}^{s} \cdot zwd_{r}^{{}} + \lambda_{3} \cdot \bar{N}_{r,3}^{s} + {\text{IFCB}}^{s} + \varepsilon_{r,3}^{s} } \hfill \\ \end{array} } \right.$$
(3)

In this set of equations, the code and phase hardware delay from receiver and satellite were fully absorbed by the receiver clock, ionosphere, and ambiguity parameter (Li et al. 2018). For the third-frequency pseudorange observables, an extra code bias parameter for each satellite and receiver was added into the estimation. As an alternative treatment, one can fix the satellite inter-frequency pseudorange bias using the corrections provided by Wang et al. (2016) and then only estimate the receiver code bias parameter for the third-frequency pseudorange observation. All the estimated parameters in our PPP models include:

$$X = [x,\;y,\;z,\;\bar{t}_{r} ,\;zwd,\;\bar{I}_{r,1}^{s} ,\;{\text{IFB}}_{r}^{s} ,\;\bar{N}_{r,1}^{s} ,\;\bar{N}_{r,2}^{s} ,\;\bar{N}_{r,3}^{s} ]$$
(4)

The ambiguity parameters and static coordinates of stations are considered as constant. In the kinematic mode, the coordinates of stations are modeled as white noise. Normally, the clock parameter is treated as an epoch-wise parameter for a single-system PPP. For tropospheric processing, the global mapping functions (Boehm et al. 2006) are used to project the slant dry and wet delays to the zenith delays. A spectral density value of 10−8 m2/s is set for the ZWD parameter. Besides, the ionosphere parameters are estimated as a random walk, and the corresponding spectral density value is set to 10−3 m2/s.

The proper weighting of the carrier phase and pseudorange observations is also an important factor for precise parameter estimation. A realistic stochastic model should consider the cross and time-correlation between observations. However, it is numerically challenging to estimate the related coefficients. In practice, the simplified stochastic modeling which is elevation dependent and ignores the correlations between the undifferenced observations is easy to implement, has a small computation burden, and also can produce a good estimation of the coordinate and other parameters (Cetin et al. 2019; Paziewski et al. 2019). Hence, it is assumed that there are no correlations among the triple-frequency observables. The study on the cross and time-correlation of GNSS observations could be investigated as a subject for further research.

The elevation-dependent weighting of observations based on the function 1 + 1/\(\sin^{2} \left( {el} \right)\), where \(el\) is the elevation angle, is used to mitigate the effects of multipath, as well as of residual atmospheric errors (Li et al. 2018). The relative weighting of code and carrier phase observations at the same frequency is chosen as 1:104. It is demonstrated that the precise orbit and clock of BDS-2 is less precise than that of GPS and Galileo (Guo et al. 2016). Besides, the BDS-2 pseudoranges suffer from satellite-induced code bias (Wanninger and Beer 2015). Although a correction model is proposed, the residual code bias still exists. For these reasons, we down weight the BDS observations, especially BDS-GEO. For the observation of the first two frequencies, the system-related weighting of GPS/Galileo/BeiDou MEO + IGSO/BeiDou GEO observations is assumed to be 16:16:4:1. Considering that no precise PCO/PCV information is currently available for the third frequency, we set the weighting for the first two frequencies and the third one to 4:1.

IFCB estimation

As for the third-frequency observations, the IFCB was first noticed in the L5 signal from the GPS Block IIF satellites (Montenbruck et al. 2012; Li et al. 2013). The investigation by Zhang et al. (2017) indicated that not only GPS but also a few BeiDou-2 satellites had significant IFCB in the third-frequency carrier phase. To make full use of the third-frequency observations for AR, the IFCB correction is suggested to be estimated in advance.

The observation equation of IFCB can be expressed as:

$$\begin{aligned} {\text{DIF}} & \equiv {\text{IF}}_{{r,L_{1} L_{2} }}^{s} (t) - {\text{IF}}_{{r,L_{1} L_{3} }}^{s} (t) \\ & = \left( {\frac{{cf_{1}^{{}} }}{{f_{1}^{2} - f_{2}^{2} }} - \frac{{cf_{1}^{{}} }}{{f_{1}^{2} - f_{3}^{2} }}} \right) \cdot \bar{N}_{1}^{s} - \frac{{cf_{2}^{{}} }}{{f_{1}^{2} - f_{2}^{2} }} \cdot \bar{N}_{2}^{s} + \frac{{cf_{3}^{{}} }}{{f_{1}^{2} - f_{3}^{2} }} \cdot \bar{N}_{3}^{s} - \;\frac{{f_{3}^{2} }}{{f_{1}^{2} - f_{3}^{2} }} \cdot {\text{IFCB}}_{{}}^{s} (t) \\ \end{aligned}$$
(5)

where \({\text{IF}}_{{r,L_{i} L_{j} }}^{s} (t) = \frac{{f_{i}^{2} }}{{f_{i}^{2} - f_{j}^{2} }}L_{r,i}^{s} (t) - \frac{{f_{j}^{2} }}{{f_{i}^{2} - f_{j}^{2} }}L_{r,j}^{s} (t)\) denotes the ionospheric-free combination of \(L_{i}\) and \(L_{j}\) carrier phases and \({\text{DIF}}\) means the difference between L1/L2- and L1/L3-based ionospheric-free combination. One common method for IFCB estimation is averaging the epoch-differenced IFCB and then recovering the undifferenced IFCB by constraining the IFCB at the first epoch to zero (Li et al. 2013; Zhang et al. 2017). Instead of this epoch-differenced processing, the least-squares adjustment is employed to estimate all IFCBs in this study simultaneously. After cycle-slip detection, the observations are divided into different arcs. Then, the \({\text{DIF}}\) combinations of each continuous arc are derived. The IFCB is estimated satellite by satellite and a zero-mean a priori constraint over 1 day is added to eliminate the rank defect of 1. In our IFCB estimation, outliers due to remaining cycle slips with small size, multipath, and so on, will degrade the parameter estimation. Therefore, the iteratively reweighted least-squares algorithm is used to mitigate the influence of outliers. The observations whose absolute standard residuals are over 2 would be down weighted to decrease its contribution to the estimation. The iteration would end when no outliers found.

Considering the daily variations of the IFCBs are not over 3 dm, and the IFCB has high short-time stability during a few minutes, it is possible to model the IFCB with a piece-wise linear function model rather than every epoch. We have estimated the IFCB for a period of 30 s and 60 s. The new 30 s IFCB corrections were obtained by linearly interpolating the 60 s corrections and then compared with the original 30 s corrections. The comparison showed that the maximum difference is 3.2 mm, and the root mean square (RMS) of the differences is 1.1 mm. The results indicated that the piece-wise linear model with an interval of 1 min is precise enough for the IFCB estimation. Compared with the 30 s IFCB estimation, the precision loss can be safely ignored, while the computation burden can be reduced significantly. Therefore, we choose to estimate the IFCB using the piece-wise linear model with a time interval of 1 min in this study.

The estimated IFCB corrections on day of year (DOY) 184, 2019 are displayed in Fig. 1. From this figure, we find that the peak–peak amplitude is about 1–2 dm for GPS satellites, 2–3 cm for BeiDou satellites, and less than 1 cm for Galileo ones. Over 99% of single-epoch Galileo IFCB is smaller than 2 mm with a standard deviation of about 1 mm. The Galileo observation is almost not impacted by IFCB errors. Therefore, we are only concerned the IFCB estimation for GPS and BeiDou hereafter. Also, it is noted that the GPS IFCB has a periodic behavior. As discussed in Li et al. (2013), Montenbruck et al. (2012), and Pan et al. (2018), at least the periods of 12, 8, and 6 h are identified. The 12-h period can be explained by that the satellite takes in the same sun illumination, while the 6-h period can be understood in that the two positions of the satellite, before and after 6 h, have the equal heat. The 8-h period cannot be well explained now. More information about the internal and external satellite properties is needed to identify the influence factors.

Fig. 1
figure 1

IFCB estimation of GPS, Galileo, and BeiDou IGSO and MEO (top), and the IFCB amplitude (bottom) of each satellite on DOY 184, 2019

Full size image

UPD estimation and PPP AR

We estimate the undifferenced and uncombined UPD with the high-precision and low-noise integer linear combination of float ambiguities. Three combinations were derived by ambiguity decorrelation from the Z-transformation of LAMBDA method (Teunissen 1995). Noted that, not the same linear combinations would be derived for each receiver-satellite ambiguity arc. Nevertheless, the first two combinations are generally the classic EWL and WL with a long wavelength and low ionospheric delay. A set of equations can be formed with the triple-frequency float PPP ambiguities from a network. Then, the single difference between satellite UPDs is estimated by fixing the UPD of one healthy satellite to zero.

Different from Li et al. (2018), we furthermore improved the float ambiguity estimation by correcting the IFCB in PPP processing. We also compared the GPS and BeiDou UPD estimators with and without IFCB correction to analyze the influence of IFCB on UPD.

Figure 2 shows the RMS of the EWL ambiguity residuals every 30-min session. On average, correcting IFCB can reduce the RMS from 0.023 to 0.022 cycles for BDS, and from 0.039 to 0.012 cycles for GPS. The improvement of EWL residuals RMS is marginal for BDS; however, it reaches to 69.2% for GPS.

Fig. 2
figure 2

RMS of GPS and BDS EWL ambiguities residuals, with and without IFCB correction on DOY 182, 2019

Full size image

The two sub-figures of Fig. 3 show the GPS EWL UPD time series without (top) and with (bottom) IFCB correction. The daily variation of GPS EWL UPD reaches 0.34 cycles without IFCB correction, while it is less than 0.03 cycles with IFCB correction. The stability of EWL UPD is greatly impacted by the GPS IFCB at the decimeter level. Here, the BDS EWL UPD time series without and with IFCB correction are not shown because they nearly overlap.

Fig. 3
figure 3

GPS EWL UPD time series on DOY 182, 2019 without (top) and with (bottom) IFCB correction

Full size image

It indicated that the influence of IFCB on the GPS UPD estimation is significant, while it is marginal on the BeiDou UPD estimation. The IFCB is suggested to be corrected to produce high-quality GPS and BeiDou UPD estimates.

In this study, the three systems UPDs were estimated with an interval of 30 min. The time series of UPD on DOY 182, 2019 is shown in Fig. 4. For the UPD on the L1 frequency, it is found that the daily variation of BeiDou UPD is much larger than that of GPS and Galileo UPD. This is because the number of BeiDou stations is fewer than that of GPS and Galileo. Especially for BeiDou MEO satellites (C11, C12, C14), the full-time tracking by MGEX cannot be guaranteed. Hence, the UPDs for some BeiDou satellites could not be provided for the full 24-h interval. The daily variation of GPS and Galileo L1 UPDs is less than 0.25 cycles. For all satellites, the daily variation of EWL, WL linear combination UPDs is less than 0.03 and 0.05 cycles, respectively. The variations of EWL and WL UPD between adjacent sessions are generally less than 0.01 cycles. The daily variation of WL UPD is more stable than that of L1 one, and the EWL UPD had the highest stability.

Fig. 4
figure 4

UPD estimation of GPS, Galileo, and BeiDou IGSO and MEO on DOY 182, 2019. For each system, the raw L1, WL, and EWL UPD are given in the sub-figures from left to right. The horizontal axis covers 24 h (48 30-min sessions), and the vertical axis is in cycles

Full size image

It is interesting to see that almost all Galileo satellites have a zero EWL UPD value except E24 whose EWL UPD is about -0.063 cycles. We also calculated the Galileo EWL UPD for the whole year of 2019 to confirm these results. The anomaly for E24 was still observed, as shown in the histograms of E24 and all other Galileo EWL UPDs in Fig. 5. The zero-value phenomenon has also been reported by (Li et al. 2019). They formed the EWL ambiguities from the Hatch–Melbourne–Wübbena (HMW) (Hatch 1982; Melbourne 1985; Wübbena 1985) combination measurements based on code and phase observations on E5a and E5b frequencies, and the EWL UPD of all Galileo satellites is zero. In fact, the equivalent EWL/WL ambiguities can be obtained by means of parameter estimation and HMW linear combination. The equivalency for WL UPD from these two methods has been demonstrated in Cheng et al. (2017). Therefore, the EWL/WL UPD is determined by code and phase observations itself and is not related to the precise clock and orbit products. This explains the reason for the same zero EWL UPD for Galileo except for E24 in our study and in Li et al. (2019). However, further studies are required to understand the essential reason for this phenomenon.

Fig. 5
figure 5

Histogram of EWL UPD in the whole year of 2019 for all Galileo satellites except E24 (left) and E24 (right)

Full size image

For the AR processing at the user, after the ambiguity-float PPP processing, the undifferenced and uncombined ambiguity estimation was used to form the single-differenced ambiguities between satellites on each frequency. Then, the fractional parts of the single-differenced ambiguities were corrected with the UPDs. The LAMBDA (Teunissen 1995) algorithm was employed to search for the most likely integer solution. The partial ambiguity resolution method proposed by Li and Zhang (2015) is employed to increase the ambiguity fixing probability.

Different from Galileo and BeiDou systems, only parts of GPS satellites have been updated to provide triple-frequency signals. With the traditional reference satellite selection strategy, it is possible to select one satellite that has only dual-frequency signals. The traditional reference satellite selection strategy should be extended to involve the GPS L5 ambiguity in AR. In this study, we propose to select the reference satellite for each frequency based on the product of elevation angle and continuous epoch count. For the GPS, the reference satellite on L1&L2 frequencies is the same one and may be different from the one on the L5 frequency.

Experiment cases and data processing strategies

The GNSS observations from IGS MGEX, with 30 s sampling intervals, were used to evaluate the performance of the GNSS triple-frequency PPP AR. The information on the available triple-frequency GPS, Galileo, and BeiDou satellites on July 1, 2019 is summarized and given in Table 1.

Table 1 List of the triple-frequency GPS/Galileo/BeiDou satellites
Full size table

Figure 6 shows the distribution of the eight GNSS stations used for the PPP tests. The daily observations from DOY 182 to 196, 2019 were used in the experiments. At the user, daily observables were separated into 24 hourly observation sessions for experiments.

Fig. 6
figure 6

Distribution of the user stations for investigating the performance of triple-frequency PPP AR

Full size image

The ‘GBM’ precise satellite orbit and clock products from the GFZ MGEX analysis center were used (Deng et al. 2014). The elevation-dependent BeiDou satellite-induced code biases were corrected following Wanninger and Beer (2015). The elevation cutoff angle was set to 10 degrees. The geometry-free and HMW linear combinations were used in the preprocessing step to detect the large cycle slips. During the PPP calculation, the detection, identification, and adaptation (DIA) procedure was further applied to find out the remaining cycle slips during the residual-check process (Teunissen 2018). To lower the risk of a wrong fixing, we simultaneously used the bootstrapping success rate and the ratio test (Ji et al. 2010) to validate the integer ambiguities. The success rate and ratio test threshold are set as 0.95 and 3.0, respectively.

Four groups of PPP solutions were performed and compared: dual- and triple-frequency GPS PPP AR, dual- and triple-frequency GPS + Galileo + BeiDou PPP AR. They are named ‘G-2,’ ‘G-3,’ ‘GEC-2,’ and ‘GEC-3’ for short, respectively. Not only the positioning accuracy but also the convergence time of four PPP solutions is analyzed. In this study, the convergence means that the time needed to attain a 3D positioning error less than 10 cm for at least 10 continuous epochs.

PPP AR results and convergence analysis

First, the kinematic results between 02:00 and 05:00 on DOY 182, 2019 from station CUT0 were taken as an example for comparison. The 3D error of four kinematic PPP solutions is shown in Fig. 7. During this period, the number of satellites observed was 10.8, 7.7, and 6.4 for GPS, Galileo, and BeiDou, respectively. There were only five epochs when a cycle-slip was detected in the carrier phase of one satellite among all 360 epochs. This indicates that the observations of this session had good data quality.

Fig. 7
figure 7

3D error time series for four groups of kinematic PPP solutions with observations between 02:00 and 05:00 at station CUT0, on DOY 182, 2019

Full size image

On average, the dual-frequency PPP AR solutions ‘G-2’ and ‘GEC-2’ succeeded in a first convergence with 28.1 and 12.8 min, respectively. With the triple-frequency observations, the solutions ‘G-3’ and ‘GEC-3’ have a convergence time of 26.2 and 12.1 min, respectively. Considering the positioning error RMS after an initialization time of 0.5 h, solutions ‘G-2’ and ‘G-3’ nearly achieved similar 3D positioning accuracy (5.3 and 5.0 cm, respectively). Similarly, the 3D positioning error of GNSS PPP solutions, ‘GEC-2’ and ‘GEC-3,’ is very close to each other (2.7 and 2.5 cm, respectively). Generally, the ‘GEC-3’ solution has the best performance, ‘GEC-2’ the second, ‘G-3’ the third, while ‘G-2’ has the worst performance.

For this example, the convergence time of ‘G-3’ was 2–4 min less than that of ‘G-2.’ The average usable dual-frequency and triple-frequency GPS satellite numbers were 10.8 and 4.9, respectively. Therefore, the effect of the third observation on GPS PPP was greatly limited. A slightly better convergence performance can be reasonably expected in ‘G-3’ solution if all GPS satellites can transmit L5 signals. As for the combined GNSS, the triple-frequency GNSS PPP AR can accelerate the convergence by about 1–2 min. Considering the positioning error after convergence, we saw that the 3D position difference between solution ‘G-2’ and ‘G-3’ was insignificant. This was also the case for the comparison of ‘GEC-2’ and ‘GEC-3.’ A similar phenomenon has been reported in Li et al. (2018) with triple-frequency BDS PPP AR. These results indicate that adding the third-frequency observations had a visible impact on the convergence time but only a marginal impact on the positioning error.

Furthermore, these results provide an analysis of the contribution of multi-GNSS to PPP AR performance. The great benefit and contribution of multi-GNSS to dual-frequency ambiguity-float and -fix PPP has also been confirmed in studies such as Li et al. (2015), Pan et al. (2017c) and Guo et al. (2018).

Moreover, we have calculated the RMS of the positioning errors for each epoch, as shown in Fig. 8, to give overall kinematic positioning performance information, including positioning accuracy and convergence for all test samples. Specifically, the bar plots in Fig. 9 are used to show the 3D errors for different initialization times from 5 to 25 min. It can be seen that in the single GPS solution, adding the L5 observation can improve the position estimate within the initialization time of about 40 min. For multi-GNSS combined PPP, adding the third frequency can accelerate the convergence within the initialization time of about 20 min. The dual-frequency and triple-frequency solutions have almost the same 3D positioning error after an initialization time of about half an hour. Among the four groups of solutions, the fastest convergence and the highest positioning accuracy are achieved by the GPS + Galileo + BeiDou triple-frequency PPP AR. The statistical results of convergence time were 25.3, 23.4, 12.8, and 10.8 min for the four groups of solutions, respectively. Triple-frequency GPS/GNSS PPP AR improved the convergence time by 7.5% and 15.6%, compared with the traditional dual-frequency GPS and GNSS PPP AR solutions, respectively. The statistical results of 3D positioning errors at the last epoch were 4.1 cm for GPS solutions while 2.2 cm for GNSS solutions. It indicated that for ambiguity-fixed PPP, the effect of the third frequency on positioning accuracy is negligible.

Fig. 8
figure 8

RMS Series of 3D positioning error of four groups of kinematic PPP solutions for DOY 182–196, 2019 over all 1-h passes

Full size image
Fig. 9
figure 9

RMS Series of 3D kinematic positioning error for initialization times of 5, 10, 15, 20, and 25 min

Full size image

The performance comparison in static PPP mode was also conducted, as shown in Figs. 10 and 11. A similar phenomenon can be found in kinematic PPP. The conventional dual-frequency GPS PPP AR has an average convergence time of about 18.6 min and an hourly 3D positioning error of 3.1 cm. For the triple-frequency GNSS PPP, the averaged convergence time is 9.3 min and the hourly 3D positioning error is 1.9 cm. For static PPP AR, triple-frequency PPP AR improved the convergence time by about 7.0% and 17.1% compared with ‘G-2’ and ‘GEC-2’ solutions, respectively. As to the hourly positioning accuracy, the difference of the 3D error is less than 0.1 cm between dual-frequency and triple-frequency PPP AR.

Fig. 10
figure 10

RMS Series of 3D positioning error of four groups of static PPP solutions for DOY 182–196, 2019 over a 1-h pass

Full size image
Fig. 11
figure 11

RMS Series of 3D static positioning error for initialization times of 5, 10, 15, 20, and 25 min

Full size image

Conclusions and remarks

This study realized triple-frequency PPP AR for GPS, Galileo, and BeiDou-2 combined system for the first time using the raw observation model. We proposed to estimate the GPS and BeiDou inter-frequency clock bias with a least-squares method using a piece-wise linear model. The IFCB was first corrected to make the L1/L2 based satellite clock bias corrections applicable for L5 observations. With the IFCB, the triple-frequency UPD estimation was realized for the real GPS observation for the first time. Also, the GPS EWL UPD quality has been significantly improved with the IFCB correction. The influence of IFCB on the BeiDou EWL UPD estimation is marginal. Because of not being contaminated by the IFCB, the raw UPD estimation method was directly employed for Galileo, which had 24 satellites in operation. An interesting phenomenon was found that all Galileo satellites except E24 had a zero-value EWL UPD. For the reference satellite selection, considering a special case that not all GPS satellites can transmit three-frequency signals, we proposed to select one reference satellite for the first two frequencies and the other for the third frequency, respectively, to form single difference between satellites.

With the multi-GNSS observations provided by MGEX covering 15 days, the positioning solutions of GPS + Galileo + BeiDou combined triple-frequency PPP AR have been conducted and analyzed. The experimental results using the GNSS observations from MGEX showed that the successful triple-frequency GNSS PPP AR has a significant improvement on the convergence time, compared with dual-frequency ambiguity-fixed PPP. Triple-frequency GNSS kinematic PPP AR produced the fastest convergence time of 10.8 min, and the highest averaged 3D positioning accuracy of 2.2 cm after convergence. Triple-frequency GPS and GNSS PPP AR reduced the average convergence time by 7.5% and 15.6% compared with dual-frequency GPS and GNSS PPP AR, respectively. It is found that adding the third frequency has a smaller contribution to improving the position solutions compared with the integrated multi-GNSS. The additional third frequency has only a marginal contribution to positioning accuracy after convergence.

In our PPP model, the ionospheric parameter is estimated as unknown. In the future, we aim to consider the regional augmentation concept to realize the rapid convergence with the support of the ionospheric and tropospheric corrections from a small regional reference network. This method is also expected to be validated and applied in real-time PPP service.