GPS / BDS-2 / Galileo Precise Point Positioning Ambiguity Resolution Based on the Uncombined Model

: In this study, an uncombined precise point positioning (PPP) model was established and was used for estimating fractional cycle bias (FCB) products and for achieving ambiguity resolution (AR), using GPS, BDS-2, and Galileo raw observations. The uncombined PPP model is ﬂexible and e ﬃ cient for positioning services and generating FCB. The FCBs for GPS, BDS-2, and Galileo were estimated using the uncombined PPP model with observations from the Multi-GNSS Experiment (MGEX) stations. The root mean squares (RMSs) of the ﬂoat ambiguity a posteriori residuals associated with all of the three GNSS constellations, PPP method, is modeled as a white noise process. The di ﬀ erences of the ISB series between BDS-2 and Galileo indicate that the clock datum bias of the satellite clock o ﬀ set estimation accounts for the variation of the ISB series.

in ionospheric modeling [38,39], differential code bias (DCB) estimation [40], and FCB estimation [8]. Li (2013) used the uncombined PPP method to shorten the convergence time with ambiguity-fixed solutions [41]. Gu (2015) and  estimated FCB products from raw ambiguities using the uncombined PPP method with BDS-2 triple-frequency observations [42,43]. Xiao (2019) proposed a common uncombined model to achieve multi-frequency PPP AR using the triple-frequency signals from Galileo and BDS-2 [44]. The equivalence of the FCB estimation between observations combined PPP model and uncombined PPP models has been confirmed and verified [8,45]. The advantages of flexibly and efficiently estimating FCBs using the uncombined PPP model have also been analyzed [8].
With the rapid development of GNSS, the potential of the uncombined PPP method for shortening convergence time and for improving positioning accuracy with ambiguity resolution should be further delineated.
Furthermore, the inter-system bias (ISB) resulting from system-specific differences, such as different time and coordinate systems, needs to be properly modeled in multi-GNSS PPP. Theoretically, ISB mainly originates from the differences of the receiver/antenna hardware delays between different GNSSs. However, the datum differences between different GNSSs in satellite clock offset estimation are introduced into the multi-GNSS PPP model, which are reflected in the ISBs. Generally, three different strategies are introduced into multi-GNSS PPP models, where the ISB parameter is modeled either with a random constant, a random walk, or a white noise process, respectively [17,[46][47][48][49]. Zeng (2017) analyzed the ISB stability of BDS-2 relative to GPS with six types of receivers, and the great stable variation of the one day ISB series was shown in the experiment [50]. Zang (2020) further investigated the daily stability of the ISBs among GNSSs and proposed to add constraints using the stable property of the ISB to improve the positioning performance [47]. Jiang (2017) also proposed a model for GPS/BDS-2 PPP, constrained with the predicted ISBs [51]. Zhou (2019) compared the positioning differences using three different ISB estimation models, which indicated that the ISB variation was also related to the clock datum differences among the different GNSSs in satellite clock correction estimation [46]. Therefore, the ISB parameters should be carefully considered to precisely deal with satellite clock corrections in multi-GNSS PPP AR.
This work contributes to the literature by demonstrating an uncombined PPP model for estimating the FCB products and by analyzing the positioning performance using GPS, BDS-2, and Galileo observations. In addition, the ISBs developed with the GPS, BDS-2, and Galileo observations are also carefully evaluated. Herein, we firstly present the uncombined PPP model with multi-GNSS observations, and then, FCB products are estimated with raw ambiguities. Then, the data process strategy is elaborated in detail. The performance of the ISB stability and FCB estimations is verified with the GPS, BDS-2, and Galileo dual-frequency raw observations from the globally distributed Multi-GNSS Experiment (MGEX) stations. Furthermore, we evaluate the improved performance of the PPP ambiguity resolution with the combined model using GPS/BDS-2/Galileo observations.

Methodology
Firstly, the uncombined PPP model was established using multi-GNSS raw observations with the presented basic observation equation. The stochastic models of ISB in the proposed multi-GNSS PPP model are given in detail. The FCBs were estimated using the uncombined PPP model and were also used for ambiguity resolution.

Basic Code and Carrier Phase Observation Equation
The code and carrier phase measurement for receiver r and satellite s of the GNSS q are described as [52]: where P q,s r, f and L q,s r, f denote the observed minus the computed values of the code and carrier phase observations at the frequency f ( f = 1, 2) (m); u q,s r is the unit vector of the component from the receiver r to the satellite s; x is the vector of the receiver positioning increments (m) relative to a priori position; c is the light speed; dt q r is the receiver clock; dt q,s is the satellite clock; Mw q,s r is the elevation-dependent wet mapping function of the tropospheric zenith wet delay (ZWD), ZWD r (m); I q,s r,1 is the ionospheric delay along the line-of-sight of a receiver to a satellite at the first frequency and γ q,s f are the respective code and phase biases at the frequency f for the satellite; and ε P, f and ε L, f are the respective code and phase measurement noises, which include the multipath effect.

Uncombined PPP with Ambiguity Resolution
In practice, the final precise satellite clock corrections produced by the ionospheric-free combination are adopted for fixing the satellite clock offset parameters, which assimilate the ionospheric-free combination of the satellite code bias. Hence, the clock parameters for the receiver and satellite, as well as for the ionospheric delays and ambiguities, are reparameterized as [8,53,54]: where DCB r = d r,1 − d r,2 and DCB s = d s 1 − d s 2 are the receiver and satellite differential code biases (DCBs), respectively. After removing the satellite clock offset parameters in the following equations, the uncombined multi-GNSS PPP model can be rewritten as: Considering the differences of the receiver clock offset between the different GNSSs, the ISB parameters are introduced into the multi-GNSS PPP model. Choosing the GPS receiver clock offset parameter as the reference, the linearized observation model of the uncombined multi-GNSS PPP model using the GPS, BDS-2, and Galileo observations can be rewritten as: The estimated parameter vector X for the uncombined multi-GNSS PPP model is presented as: where I r,1 is the vector of the ionospheric delay parameters in Equation (2) for all observed satellites and N r, f is the vector of the ambiguity parameters in Equation (2) of each frequency for all observed satellites. It has been confirmed that the ISBs between the different GNSSs are relatively stable over a short time span. However, the ISB is not only caused by the differences in the receiver/antenna hardware delays between the different GNSSs, but also by the datum differences of the satellite clock offset. The precise corrections for the satellite clock from external products are used to fix the satellite clock parameters, which crucially determine the positioning accuracy in the PPP model. Hence, the ISB estimation cannot be simply estimated as constant, as per previous studies. Therefore, we estimated the ISB parameters with the white noise process, which is an efficient method to avoid the effects of satellite clock datum differences and to analyze their stability in one day. In the proposed multi-GNSS PPP method, the ISB parameters estimated with the white noise process can be described as [46]: where ISB q r (k) SPP refers to the initial values estimated with pseudorange measurements at the epoch k, and σ 2 ISB is the a priori variance of the ISB parameters.

FCB Estimation
The float ambiguity shown in Equation (2) can be denoted as the combination of the integer value and its fractional cycle biases in the receiver and satellite ends. Hence, the estimated ambiguities are rewritten as: where N q,s r, f is the integer ambiguity, and B q r, f and B q,s f are the FCBs for the receiver and the satellite, respectively.
In the combined PPP model, the ionospheric-free ambiguity is conveniently reformed by the WL and NL ambiguities for ambiguity resolution. The WL ambiguity is easily fixed to an integer value for its relatively longer wavelength [44]. In the uncombined PPP model, the created linear combinations of the ambiguities also have higher accuracy and are less correlated than the individual ambiguities of each frequency. These selected combinations of the carrier phase are expected to have relatively long wavelengths and to avoid significant ionosphere errors. After systematic investigation, two integer linear combinations, i.e., (4, −3) and (1, −1), with respective lane widths of 11 and 86 cm, were chosen [43]. Hence, the combined ambiguities can be formed with the original ambiguities as [43]: where N q,s r, (4,−3) and N q,s r,(1,−1) are denoted as NL and WL ambiguities for convenience. Consequently, the receiver or the satellite FCBs are also reformed as: The FCBs can be estimated using the reformed undifferenced ambiguities from a referenced network [29,30]. Both of the above ambiguity combinations in cycle can be written as: where ∆n is the fractional part of the float ambiguity solution N s r (cycle); N s r is the integer value; B r is the FCB for the receiver r, while the FCB of the satellite s is B s . Assuming n satellites are tracked by m stations, the combined float ambiguities in Equation (9) of all station-satellite pairs in an FCB estimation network are expressed as [8]: where R r and S s are the coefficient vectors for the FCB of the respective receiver and satellite. R r has m elements, among which the element r is 1 and the other elements are zero. S s has n elements, among which the element s is −1 and the other elements are zero. To avoid the rank defect in Equation (12), one receiver FCB parameter is fixed to zero. The FCBs are estimated with 15 min intervals. All of the fractional parts of the float ambiguities from each station are inputted into the Kalman filter to estimate the receiver and satellite FCBs and to maintain the consistency of the FCBs between the different receivers for one satellite [41,55]. In the Kalman filter of the FCB estimation, the receiver FCBs are estimated as a white noise model. For the satellite FCBs, they can be characterized as constant parameters in the WL combination and as random walk processes in the NL combination with noise 0.05 cycles. After estimating the WL and NL combination FCBs, the individual FCBs of each frequency can be recovered using the relationship in Equation (10). Hence, users can adopt the individual FCBs of each frequency to correct float ambiguity for achieving PPP AR.

Data and Processing Strategy
In this paper, FCB estimation, ISB estimation, and PPP positioning evaluations were performed with the data from MGEX on 3 September 2019. Dual-frequency observations (GPS: L1/L2; BDS-2: B1/B2; Galileo: E1/E5a) with 30 s intervals were used in the processing. The cutoff angle for elevation was 7 • for the PPP processing and 15 • for the AR, while it was 30 • for the ambiguities of participating in the FCB estimation. The precise ephemeris and clock offset from GeoForschungsZentrum (GFZ) were used to fix the satellite orbit and the clock parameters. Considering the abnormal stability of the ISB estimated from GFZ products [46], the white noise model was used for the ISB parameters. The zenith tropospheric delay was estimated as a random walk process, and the ionospheric delays were estimated as a white noise process. The receiver clock offset was estimated as a white noise process while the ambiguities were estimated as constant parameters. In the processing, the phase center offset (PCO) and variation (PCV) values for the BDS-2 satellite published by the European Space Agency (ESA) were corrected [56], and the values in the IGS14.atx file were used for the other GNSSs. For the satellite-induced code multipath effects of the BDS-2 inclined geosynchronous orbit (IGSO) and medium Earth orbit (MEO) satellites, the recommended corrections in [57] were adopted. For the other terms, such as relativistic delay, Sagnac effect, phase windup effect, and tide displacement, they were corrected with the corresponding models. For ambiguity resolution, the optimal integer solution was fixed using the least-squares ambiguity decorrelation adjustment (LAMBDA) algorithm [58]. The threshold value of the ratio test was 2.0. Partial fixing of the ambiguities was achieved with a minimum of four satellites for each GNSS system.

ISB Stability
We selected 24 stations in the MGEX network, as shown in Figure 1, which can observe the GPS, BDS-2, and Galileo satellite signals for PPP processing and ISB analysis. Firstly, the ISBs were estimated in the uncombined PPP method with float solutions. The one day solutions were used to analyze the stability. In Figure 2, the mean values of the ISB one day solutions for the 24 stations are presented for BDS-2 and Galileo, relative to GPS. The magnitude of the BDS-2 ISB values is dozens of nanoseconds, which is significantly larger than that of Galileo. In addition, the magnitude of the ISB values is strongly correlated to the type of receiver. In Figure 2, the 24 stations are grouped as five receiver types to illustrate the ISB differences. As shown in Table 1, for "A" and "B", they have the same brand receivers, and thus, obtained the same magnitude of ISB, the values of which are 38-45 ns for BDS-2, while there are significant differences for Galileo. For BDS-2, the mean ISB values of "A", "B", "C", "D", and "E" are 43. 16

ISB Stability
We selected 24 stations in the MGEX network, as shown in Figure 1, which can observe the GPS, BDS-2, and Galileo satellite signals for PPP processing and ISB analysis. Firstly, the ISBs were estimated in the uncombined PPP method with float solutions. The one day solutions were used to analyze the stability. In Figure 2, the mean values of the ISB one day solutions for the 24 stations are presented for BDS-2 and Galileo, relative to GPS. The magnitude of the BDS-2 ISB values is dozens of nanoseconds, which is significantly larger than that of Galileo. In addition, the magnitude of the ISB values is strongly correlated to the type of receiver. In Figure 2, the 24 stations are grouped as five receiver types to illustrate the ISB differences. As shown in Table 1, for "A" and "B", they have the same brand receivers, and thus, obtained the same magnitude of ISB, the values of which are 38-45 ns for BDS-2, while there are significant differences for Galileo. For BDS-2, the mean ISB values of "A", "B", "C", "D", and "E" are 43. 16

ISB Stability
We selected 24 stations in the MGEX network, as shown in Figure 1, which can observe the GPS, BDS-2, and Galileo satellite signals for PPP processing and ISB analysis. Firstly, the ISBs were estimated in the uncombined PPP method with float solutions. The one day solutions were used to analyze the stability. In Figure 2, the mean values of the ISB one day solutions for the 24 stations are presented for BDS-2 and Galileo, relative to GPS. The magnitude of the BDS-2 ISB values is dozens of nanoseconds, which is significantly larger than that of Galileo. In addition, the magnitude of the ISB values is strongly correlated to the type of receiver. In Figure 2, the 24 stations are grouped as five receiver types to illustrate the ISB differences. As shown in Table 1, for "A" and "B", they have the same brand receivers, and thus, obtained the same magnitude of ISB, the values of which are 38-45 ns for BDS-2, while there are significant differences for Galileo. For BDS-2, the mean ISB values of "A", "B", "C", "D", and "E" are 43. 16 Figure 3 shows the variation of ISBs with the means removed for BDS-2. A significant gradient trend is presented in the left panel for all ISB series. This changing trend enlarged the average STD of all ISB series to 1.32 ns, which indicates that the clock datum differences in the satellite clock correction estimation between the different GNSSs are assimilated into the ISBs. After removing the mean trend, the variations of the one day ISBs are shown in the right panel, in which the average STD is improved to 0.29 ns. Figure 4 shows the one day ISB variations for Galileo, in which the average STD is 0.25 ns.
Remote Sens. 2020, 12, x FOR PEER REVIEW 8 of 20  Figure 3 shows the variation of ISBs with the means removed for BDS-2. A significant gradient trend is presented in the left panel for all ISB series. This changing trend enlarged the average STD of all ISB series to 1.32 ns, which indicates that the clock datum differences in the satellite clock correction estimation between the different GNSSs are assimilated into the ISBs. After removing the mean trend, the variations of the one day ISBs are shown in the right panel, in which the average STD is improved to 0.29 ns. Figure 4 shows the one day ISB variations for Galileo, in which the average STD is 0.25 ns.   Figure 5 shows the STD of each station's one day ISB series for BDS-2 and Galileo, respectively. The original BDS-2 ISB series included the changing trend, which obtained a relatively large STD for each station. The ISB series, including the datum bias of the satellite clock offset estimation, cannot show the stability of the inter-system bias caused by the receiver-dependent hardware delay between   Figure 3 shows the variation of ISBs with the means removed for BDS-2. A significant gradient trend is presented in the left panel for all ISB series. This changing trend enlarged the average STD of all ISB series to 1.32 ns, which indicates that the clock datum differences in the satellite clock correction estimation between the different GNSSs are assimilated into the ISBs. After removing the mean trend, the variations of the one day ISBs are shown in the right panel, in which the average STD is improved to 0.29 ns. Figure 4 shows the one day ISB variations for Galileo, in which the average STD is 0.25 ns.   Figure 5 shows the STD of each station's one day ISB series for BDS-2 and Galileo, respectively. The original BDS-2 ISB series included the changing trend, which obtained a relatively large STD for each station. The ISB series, including the datum bias of the satellite clock offset estimation, cannot show the stability of the inter-system bias caused by the receiver-dependent hardware delay between  Figure 5 shows the STD of each station's one day ISB series for BDS-2 and Galileo, respectively. The original BDS-2 ISB series included the changing trend, which obtained a relatively large STD for each station. The ISB series, including the datum bias of the satellite clock offset estimation, cannot show the stability of the inter-system bias caused by the receiver-dependent hardware delay between the different GNSSs. After a straightforward removal of the mean trend of the ISB series, the STD of each station's ISB is significantly reduced compared to the original STD, except for station "TID1." Comparing with Galileo, it is observed that they have similar stability, as indicated by the STD.
Remote Sens. 2020, 12, x FOR PEER REVIEW 9 of 20 the different GNSSs. After a straightforward removal of the mean trend of the ISB series, the STD of each station's ISB is significantly reduced compared to the original STD, except for station "TID1." Comparing with Galileo, it is observed that they have similar stability, as indicated by the STD.

FCB Results
After detailed analysis of the stability of ISB, we used the uncombined PPP method to estimate the FCB products with GPS, BDS-2, and Galileo. For this, 316 stations from the MGEX network were used to estimate the FCB products, of which, all stations support the tracking of GPS constellations. As shown in Figure 6, about 138 stations were tracking the BDS-2 satellite, whereas Galileo was tracked by 190 stations. The float ambiguities for GPS, BDS-2, and Galileo in the FCB estimation were obtained from the uncombined multi-GNSS PPP model processing with coordinate constraints.

GPS FCB Results
For PPP AR, the performance depends on the accuracy of the FCB estimation, which was evaluated with the a posteriori residuals of the ambiguities that participated in the FCB estimation. The residuals were expected to be zero, which would indicate that there are highly consistent FCBs between float ambiguities. The residual distributions in the cycles were calculated for all of the epochs of the FCB estimation over one day time span. Figure 7 presents the a posteriori distribution of the residuals of the integer linear-combined ambiguities for the GPS FCBs. For the two (NL and WL) combinations, the same magnitude of the root mean square (RMS) of the a posteriori residuals, which are less than 0.07 cycles, indicates highly consistent FCB estimation. The RMS of the a posteriori

FCB Results
After detailed analysis of the stability of ISB, we used the uncombined PPP method to estimate the FCB products with GPS, BDS-2, and Galileo. For this, 316 stations from the MGEX network were used to estimate the FCB products, of which, all stations support the tracking of GPS constellations. As shown in Figure 6, about 138 stations were tracking the BDS-2 satellite, whereas Galileo was tracked by 190 stations. The float ambiguities for GPS, BDS-2, and Galileo in the FCB estimation were obtained from the uncombined multi-GNSS PPP model processing with coordinate constraints.
Remote Sens. 2020, 12, x FOR PEER REVIEW 9 of 20 the different GNSSs. After a straightforward removal of the mean trend of the ISB series, the STD of each station's ISB is significantly reduced compared to the original STD, except for station "TID1." Comparing with Galileo, it is observed that they have similar stability, as indicated by the STD.

FCB Results
After detailed analysis of the stability of ISB, we used the uncombined PPP method to estimate the FCB products with GPS, BDS-2, and Galileo. For this, 316 stations from the MGEX network were used to estimate the FCB products, of which, all stations support the tracking of GPS constellations. As shown in Figure 6, about 138 stations were tracking the BDS-2 satellite, whereas Galileo was tracked by 190 stations. The float ambiguities for GPS, BDS-2, and Galileo in the FCB estimation were obtained from the uncombined multi-GNSS PPP model processing with coordinate constraints.

GPS FCB Results
For PPP AR, the performance depends on the accuracy of the FCB estimation, which was evaluated with the a posteriori residuals of the ambiguities that participated in the FCB estimation. The residuals were expected to be zero, which would indicate that there are highly consistent FCBs between float ambiguities. The residual distributions in the cycles were calculated for all of the epochs of the FCB estimation over one day time span. Figure 7 presents the a posteriori distribution of the residuals of the integer linear-combined ambiguities for the GPS FCBs. For the two (NL and WL) combinations, the same magnitude of the root mean square (RMS) of the a posteriori residuals, which are less than 0.07 cycles, indicates highly consistent FCB estimation. The RMS of the a posteriori

GPS FCB Results
For PPP AR, the performance depends on the accuracy of the FCB estimation, which was evaluated with the a posteriori residuals of the ambiguities that participated in the FCB estimation. The residuals were expected to be zero, which would indicate that there are highly consistent FCBs between float ambiguities. The residual distributions in the cycles were calculated for all of the epochs of the FCB estimation over one day time span. Figure 7 presents the a posteriori distribution of the residuals of the integer linear-combined ambiguities for the GPS FCBs. For the two (NL and WL) combinations, the same magnitude of the root mean square (RMS) of the a posteriori residuals, which are less than 0.07 cycles, indicates highly consistent FCB estimation. The RMS of the a posteriori residuals is 0.067 cycles for the NL combination and 0.069 cycles for the WL combination. For the distribution of the residuals, 98.6% of the NL residuals and 98.9% of the WL residuals are within [-0. 25   The WL FCB series is more stable than that of the NL combination. To further analyze the stability of the FCB series, Figure 9 shows the STDs of the one day GPS FCB series for each satellite NL and WL combination. Significantly, the STDs of the NL combination for most of the satellites are larger than that of the WL combination. Additionally, the STD of the G18 NL combination is about 0.115 cycles. The G18 signals have been broadcasted by the old BLOCK IIA (SVN34) satellite since March 20, 2018, which can account for the large variation. Excluding G18, the average STD is 0.030 cycles for NL combination FCB and 0.015 cycles for WL combination FCB.   The WL FCB series is more stable than that of the NL combination. To further analyze the stability of the FCB series, Figure 9 shows the STDs of the one day GPS FCB series for each satellite NL and WL combination. Significantly, the STDs of the NL combination for most of the satellites are larger than that of the WL combination. Additionally, the STD of the G18 NL combination is about 0.115 cycles. The G18 signals have been broadcasted by the old BLOCK IIA (SVN34) satellite since 20 March 2018, which can account for the large variation. Excluding G18, the average STD is 0.030 cycles for NL combination FCB and 0.015 cycles for WL combination FCB.   The WL FCB series is more stable than that of the NL combination. To further analyze the stability of the FCB series, Figure 9 shows the STDs of the one day GPS FCB series for each satellite NL and WL combination. Significantly, the STDs of the NL combination for most of the satellites are larger than that of the WL combination. Additionally, the STD of the G18 NL combination is about 0.115 cycles. The G18 signals have been broadcasted by the old BLOCK IIA (SVN34) satellite since March 20, 2018, which can account for the large variation. Excluding G18, the average STD is 0.030 cycles for NL combination FCB and 0.015 cycles for WL combination FCB.    Figure 10 shows the distribution of the residuals of the BDS-2 NL and WL combinations for the IGSO and MEO satellites. The RMS of the residuals for the NL combination is 0.063 cycles, and is 0.091 cycles for that of the WL combination. For the distribution of the residuals, it is 98.9% and 97.3% for the NL and WL combination residuals, which are located in [−0.25, 0.25] (cycles). The consistency of the NL combination is better than that of the WL combination, which indicates that the stability of the BDS-2 WL combination is worse than that of GPS. The BDS-2 FCB time series are presented in Figure 11. Considering the low precision of geostationary (GEO) satellite orbits, only the IGSO and MEO satellites were included in the FCB estimation. Similarly, as for GPS, the WL FCB series are more stable than the NL FCB series. For three MEO satellites, the fluctuations of FCB are larger than that of the IGSO. In a good IGSO coverage area, IGSO satellites are almost continuously tracked. Figure 12 shows the STD of each satellite FCB series for the NL and WL combinations. For the NL combination, the STD of the MEO satellites is significantly higher than that of the IGSO satellites. Consequently, the average STD of the NL FCB series for the grouped BDS-2 MEO satellites is 0.073 cycles, while it is 0.0156 cycles for the grouped BDS-2 IGSO satellites. The average STD of the WL FCB for all BDS-2 satellites is 0.013 cycles.  Figure 10 shows the distribution of the residuals of the BDS-2 NL and WL combinations for the IGSO and MEO satellites. The RMS of the residuals for the NL combination is 0.063 cycles, and is 0.091 cycles for that of the WL combination. For the distribution of the residuals, it is 98.9% and 97.3% for the NL and WL combination residuals, which are located in [−0.25, 0.25] (cycles). The consistency of the NL combination is better than that of the WL combination, which indicates that the stability of the BDS-2 WL combination is worse than that of GPS.  Figure 10 shows the distribution of the residuals of the BDS-2 NL and WL combinations for the IGSO and MEO satellites. The RMS of the residuals for the NL combination is 0.063 cycles, and is 0.091 cycles for that of the WL combination. For the distribution of the residuals, it is 98.9% and 97.3% for the NL and WL combination residuals, which are located in [−0.25, 0.25] (cycles). The consistency of the NL combination is better than that of the WL combination, which indicates that the stability of the BDS-2 WL combination is worse than that of GPS. The BDS-2 FCB time series are presented in Figure 11. Considering the low precision of geostationary (GEO) satellite orbits, only the IGSO and MEO satellites were included in the FCB estimation. Similarly, as for GPS, the WL FCB series are more stable than the NL FCB series. For three MEO satellites, the fluctuations of FCB are larger than that of the IGSO. In a good IGSO coverage area, IGSO satellites are almost continuously tracked. Figure 12 shows the STD of each satellite FCB series for the NL and WL combinations. For the NL combination, the STD of the MEO satellites is significantly higher than that of the IGSO satellites. Consequently, the average STD of the NL FCB series for the grouped BDS-2 MEO satellites is 0.073 cycles, while it is 0.0156 cycles for the grouped BDS-2 IGSO satellites. The average STD of the WL FCB for all BDS-2 satellites is 0.013 cycles. The BDS-2 FCB time series are presented in Figure 11. Considering the low precision of geostationary (GEO) satellite orbits, only the IGSO and MEO satellites were included in the FCB estimation. Similarly, as for GPS, the WL FCB series are more stable than the NL FCB series. For three MEO satellites, the fluctuations of FCB are larger than that of the IGSO. In a good IGSO coverage area, IGSO satellites are almost continuously tracked. Figure 12 shows the STD of each satellite FCB series for the NL and WL combinations. For the NL combination, the STD of the MEO satellites is significantly higher than that of the IGSO satellites. Consequently, the average STD of the NL FCB series for the grouped BDS-2 MEO satellites is 0.073 cycles, while it is 0.0156 cycles for the grouped BDS-2 IGSO satellites. The average STD of the WL FCB for all BDS-2 satellites is 0.013 cycles.

Galileo FCB Results
The results of the Galileo FCB estimation are illustrated in Figures 13-15, which present the distribution of the residuals, the FCBs time series, and the corresponding STD values. Compared to the GPS and BDS-2 results, the distribution of the residuals for Galileo is significantly close to zero, where the RMS is 0.058 and 0.045 for the NL and WL combinations, respectively. It is 98.9% and 99.8% for the NL and WL combination residuals, which are located in [−0.25, 0.25] (cycles). The stability of the FCB series for the NL and WL combinations, as shown in Figure 14, is also better than the other two GNSSs. Figure 15 presents the STDs of each satellite FCB time series, which indicates that the average STD of the NL combination is 0.0184 cycles and that of the WL combination is 0.006 cycles. In reference [59], it is confirmed that the WL FCBs of Galileo are quite stable over a long time and can be treated as constant over several months or even years.

Galileo FCB Results
The results of the Galileo FCB estimation are illustrated in Figures 13-15, which present the distribution of the residuals, the FCBs time series, and the corresponding STD values. Compared to the GPS and BDS-2 results, the distribution of the residuals for Galileo is significantly close to zero, where the RMS is 0.058 and 0.045 for the NL and WL combinations, respectively. It is 98.9% and 99.8% for the NL and WL combination residuals, which are located in [−0.25, 0.25] (cycles). The stability of the FCB series for the NL and WL combinations, as shown in Figure 14, is also better than the other two GNSSs. Figure 15 presents the STDs of each satellite FCB time series, which indicates that the average STD of the NL combination is 0.0184 cycles and that of the WL combination is 0.006 cycles. In reference [59], it is confirmed that the WL FCBs of Galileo are quite stable over a long time and can be treated as constant over several months or even years.

Galileo FCB Results
The results of the Galileo FCB estimation are illustrated in Figures 13-15, which present the distribution of the residuals, the FCBs time series, and the corresponding STD values. Compared to the GPS and BDS-2 results, the distribution of the residuals for Galileo is significantly close to zero, where the RMS is 0.058 and 0.045 for the NL and WL combinations, respectively. It is 98.9% and 99.8% for the NL and WL combination residuals, which are located in [−0.25, 0.25] (cycles). The stability of the FCB series for the NL and WL combinations, as shown in Figure 14, is also better than the other two GNSSs. Figure 15 presents the STDs of each satellite FCB time series, which indicates that the average STD of the NL combination is 0.0184 cycles and that of the WL combination is 0.006 cycles. In reference [59], it is confirmed that the WL FCBs of Galileo are quite stable over a long time and can be treated as constant over several months or even years.   For all GNSSs, the FCB estimations for the NL and WL combinations are quite accurate, which can support the uncombined PPP AR model. Although the FCB series of the NL combination is not as stable as that of the WL combination, it is quite stable over a short time span for the NL combination, such as for 15-30 min. The relatively small fluctuations between the adjacent sessions indicate that FCB estimation of real-time processing is easily achieved. For the BDS-2 MEO satellites, their NL FCB series is not as accurate as that of the IGSO satellites, which will be improved with the   For all GNSSs, the FCB estimations for the NL and WL combinations are quite accurate, which can support the uncombined PPP AR model. Although the FCB series of the NL combination is not as stable as that of the WL combination, it is quite stable over a short time span for the NL combination, such as for 15-30 min. The relatively small fluctuations between the adjacent sessions indicate that FCB estimation of real-time processing is easily achieved. For the BDS-2 MEO satellites, their NL FCB series is not as accurate as that of the IGSO satellites, which will be improved with the   For all GNSSs, the FCB estimations for the NL and WL combinations are quite accurate, which can support the uncombined PPP AR model. Although the FCB series of the NL combination is not as stable as that of the WL combination, it is quite stable over a short time span for the NL combination, such as for 15-30 min. The relatively small fluctuations between the adjacent sessions indicate that FCB estimation of real-time processing is easily achieved. For the BDS-2 MEO satellites, their NL FCB series is not as accurate as that of the IGSO satellites, which will be improved with the For all GNSSs, the FCB estimations for the NL and WL combinations are quite accurate, which can support the uncombined PPP AR model. Although the FCB series of the NL combination is not as stable as that of the WL combination, it is quite stable over a short time span for the NL combination, such as for 15-30 min. The relatively small fluctuations between the adjacent sessions indicate that FCB estimation of real-time processing is easily achieved. For the BDS-2 MEO satellites, their NL FCB series is not as accurate as that of the IGSO satellites, which will be improved with the new-generation BDS; BDS-3 is already functioning with more operational satellites than the BDS-2 used here.

Positioning Results
To further validate the accuracy of the FCB estimation, as well as to evaluate the positioning performance with the float and fixed solutions in the uncombined multi-GNSS PPP model, the data were processed in a static mode, collected from a network, as shown in Figure 1, on 3 September 2019. The 3 h PPP AR sessions of the 24 h station data were processed with GPS, BDS-2, and Galileo observations. The positioning accuracy is firstly presented; then, the advantages of the ambiguity success fixing rate and the convergence time are presented for the uncombined multi-GNSS PPP model. The coordinates from the International GNSS Service (IGS) combined weekly solutions were used as the reference. Four combinations of GNSSs were used in the processing, namely, only-GPS, GPS/BDS-2, GPS/Galileo, and GPS/BDS-2/Galileo, which were marked as "G", "GC", "GE", and "GCE" in this study. For convenience, the ambiguity-fixed solutions were marked as "AR" and the ambiguity-float solutions as "float".

PPP AR Accuracy
The PPP ambiguity-fixed solutions were processed with four type of systems and were compared to the ambiguity-float solutions. Figure 16 shows the average RMS of the ambiguity-float and ambiguity-fixed solutions in the eastward, northward upward, and three-dimensional directions for four systems, respectively. In the eastward component, "GCE" performs the best and has higher accuracy, followed by "GC" and "GE", while "G" presents the worst performance. The same trend is also found in the other components, though there are some slight differences in the dual-system solutions. Compared to the ambiguity-float solutions, significant improvements are presented in the eastward and upward directions for the AR results. For the eastward component, the RMS of the ambiguity-float solutions is 1.02, 0.87, 0.81, and 0.78 cm for "G", "GC", "GE", and "GCE", while that of ambiguity-fixed solutions is 0.43, 0.43, 0.36, and 0.36 cm, with 58%, 50%, 56%, and 53% improvements, respectively. With PPP AR, the RMSs of the ambiguity-fixed solutions for the eastward component of the four systems are as accurate as of that in the northward direction. For the three-dimensional component, the RMSs of the four systems' ambiguity-float solutions are 1.92, 1.76, 1.68, and 1.62 cm for "G", "GC", "GE", and "GCE", and 1.34, 1.19, 1.21, and 1.14 cm for the ambiguity-fixed solutions, with 30%, 32%, 28%, and 30% improvements, respectively. This indicates that the increase in satellites from the multi-GNSS constellations significantly improves the positioning accuracy, and the ambiguity resolution further enhances this improvement.
Remote Sens. 2020, 12, x FOR PEER REVIEW 14 of 20 new-generation BDS; BDS-3 is already functioning with more operational satellites than the BDS-2 used here.

Positioning Results
To further validate the accuracy of the FCB estimation, as well as to evaluate the positioning performance with the float and fixed solutions in the uncombined multi-GNSS PPP model, the data were processed in a static mode, collected from a network, as shown in Figure 1, on 3 September 2019. The 3 h PPP AR sessions of the 24 h station data were processed with GPS, BDS-2, and Galileo observations. The positioning accuracy is firstly presented; then, the advantages of the ambiguity success fixing rate and the convergence time are presented for the uncombined multi-GNSS PPP model. The coordinates from the International GNSS Service (IGS) combined weekly solutions were used as the reference. Four combinations of GNSSs were used in the processing, namely, only-GPS, GPS/BDS-2, GPS/Galileo, and GPS/BDS-2/Galileo, which were marked as "G", "GC", "GE", and "GCE" in this study. For convenience, the ambiguity-fixed solutions were marked as "AR" and the ambiguity-float solutions as "float".

PPP AR Accuracy
The PPP ambiguity-fixed solutions were processed with four type of systems and were compared to the ambiguity-float solutions. Figure 16 shows the average RMS of the ambiguity-float and ambiguity-fixed solutions in the eastward, northward upward, and three-dimensional directions for four systems, respectively. In the eastward component, "GCE" performs the best and has higher accuracy, followed by "GC" and "GE", while "G" presents the worst performance. The same trend is also found in the other components, though there are some slight differences in the dual-system solutions. Compared to the ambiguity-float solutions, significant improvements are presented in the eastward and upward directions for the AR results. For the eastward component, the RMS of the ambiguity-float solutions is 1.02, 0.87, 0.81, and 0.78 cm for "G", "GC", "GE", and "GCE", while that of ambiguity-fixed solutions is 0.43, 0.43, 0.36, and 0.36 cm, with 58%, 50%, 56%, and 53% improvements, respectively. With PPP AR, the RMSs of the ambiguity-fixed solutions for the eastward component of the four systems are as accurate as of that in the northward direction. For the three-dimensional component, the RMSs of the four systems' ambiguity-float solutions are 1.92, 1.76, 1.68, and 1.62 cm for "G", "GC", "GE", and "GCE", and 1.34, 1.19, 1.21, and 1.14 cm for the ambiguityfixed solutions, with 30%, 32%, 28%, and 30% improvements, respectively. This indicates that the increase in satellites from the multi-GNSS constellations significantly improves the positioning accuracy, and the ambiguity resolution further enhances this improvement.

Convergence Time
In general, multi-GNSS processing can accelerate the convergence of positioning accuracy. Figure 17 shows the positioning bias variations for the four systems with the float and AR solutions in different directions. Indeed, multi-GNSS positioning improves the accuracy, compared to the only-GPS solution or the dual-system solutions at the same epoch. Meanwhile, the positioning convergence is further improved by the resolution of the ambiguity. As shown in Table 2, an accuracy level of 20 cm is achieved for all of the solutions after 10 min of observation. However, about a 10 cm level accuracy is obtained for all ambiguity-fixed solutions. With 30 min of the observations, the accuracy of the ambiguity-float solutions is better than 10 cm, while this is about 2~3 cm for the ambiguity-fixed solutions. For the ambiguity-float and the ambiguity-fixed solutions, "GCE" presents the best performance, followed by "GE" and "GC", while "G" has the worst performance. Here, the threshold value of the three-dimensional RMS was 10 cm for calculating the convergence time. The convergence times are listed in Table 3, which shows 16 min for the "GCE" ambiguity-float solution to achieve a 10 cm accuracy level, and 7.5 min for the AR solutions of "GCE" with a 53% improvement. In the ambiguity-float solutions, the maximum convergence time is 22 min for "G"; however, it is only 10.5 min for the ambiguity-fixed solutions.

Convergence Time
In general, multi-GNSS processing can accelerate the convergence of positioning accuracy. Figure 17 shows the positioning bias variations for the four systems with the float and AR solutions in different directions. Indeed, multi-GNSS positioning improves the accuracy, compared to the only-GPS solution or the dual-system solutions at the same epoch. Meanwhile, the positioning convergence is further improved by the resolution of the ambiguity. As shown in Table 2, an accuracy level of 20 cm is achieved for all of the solutions after 10 min of observation. However, about a 10 cm level accuracy is obtained for all ambiguity-fixed solutions. With 30 min of the observations, the accuracy of the ambiguity-float solutions is better than 10 cm, while this is about 2~3 cm for the ambiguity-fixed solutions. For the ambiguity-float and the ambiguity-fixed solutions, "GCE" presents the best performance, followed by "GE" and "GC", while "G" has the worst performance. Here, the threshold value of the three-dimensional RMS was 10 cm for calculating the convergence time. The convergence times are listed in Table 3, which shows 16 min for the "GCE" ambiguity-float solution to achieve a 10 cm accuracy level, and 7.5 min for the AR solutions of "GCE" with a 53% improvement. In the ambiguity-float solutions, the maximum convergence time is 22 min for "G"; however, it is only 10.5 min for the ambiguity-fixed solutions. Figure 17. Convergence of the positioning bias for "G", "GC", "GE". and "GCE" systems in the ambiguity-float and ambiguity-fixed solutions with 1.5 hour data.  Figure 17. Convergence of the positioning bias for "G", "GC", "GE". and "GCE" systems in the ambiguity-float and ambiguity-fixed solutions with 1.5 h data.

Ambiguity Success Fixing Rate
We define the ambiguity success fixing rate as the ratio of the number of fixed solutions over all sessions to that of all solutions within the same time length. Figure 18 shows the series of the ambiguity success fixing rate for the four types of solutions. The notable differences are shown in the initial stage of positioning in Table 4. With 10 min of data, the "GCE" achieves a 99.0% ambiguity success fixing rate. In contrast, it is approximately 94.3% for the "GC" solution and 86.5% for the "GE" solution, while for the "G" solution, it is only 75.5%. After 1 h of data processing, the ambiguity success fixing rate is better than 99.5% for the dual-GNSS or triple-GNSS systems, while it is only 97.9% for the only-GPS system.  We define the ambiguity success fixing rate as the ratio of the number of fixed solutions over all sessions to that of all solutions within the same time length. Figure 18 shows the series of the ambiguity success fixing rate for the four types of solutions. The notable differences are shown in the initial stage of positioning in Table 4. With 10 min of data, the "GCE" achieves a 99.0% ambiguity success fixing rate. In contrast, it is approximately 94.3% for the "GC" solution and 86.5% for the "GE" solution, while for the "G" solution, it is only 75.5%. After 1 h of data processing, the ambiguity success fixing rate is better than 99.5% for the dual-GNSS or triple-GNSS systems, while it is only 97.9% for the only-GPS system.

Discussion
When multi-frequency signals and multi-GNSS observations are present, the uncombined PPP model is expected to flexibly and efficiently compute the PPP solutions and to estimate the FCB products. In this study, the ISB series of BDS-2 and Galileo relative to GPS were estimated to show the variations. Considering the clock datum time-varying differences being coupled with the satellite clock offset estimation, unlike for Galileo, the ISB series of BDS-2 exhibited gradient trends. In this work, we simply estimated the ISB parameters by the white noise model in the PPP processing and FCB estimation to assimilate the uncertain bias. The STD of the NL FCB series for the BDS-2 MEO satellites performed worst due to less available satellites and the uneven distribution of stations. With continuous construction of the BDS-3 and Galileo constellations, further comprehensive work on the uncombined PPP method with multi-GNSS observations should be undertaken. For the ambiguity-float or -fixed solutions, the uncombined PPP model with GPS/BDS-2/Galileo observations significantly achieved higher positioning accuracy and shorter convergence time than that of the only-GPS, GPS/Galileo, and GPS/BDS-2 solutions-noting that the PPP-AR solutions with GPS/BDS-2/Galileo observations performed best among all. In this model, the convergence time was 7.5 min and final positioning accuracy was 1.14 cm.

Conclusions
We comprehensively studied the performance of uncombined multi-GNSS PPP processing, as well as FCB estimation with GPS/BDS-2/Galileo observations. The performance of the positioning accuracy, the convergence time, and the ambiguity success fixing rate was evaluated. Firstly, the differences between the BDS-2 and Galileo ISB series were presented to exhibit the fluctuations of ISBs. Although ISB is mainly affected by the receiver-dependent bias, the clock datum of the satellite clock offset estimation also accounts for the variations of the one day ISB series. Hence, we selected the white noise model for ISB estimation, considering the estimation of precise satellite clock offset products among different GNSS constellations. With the GFZ precise satellite clock products, the averaged STD of the BDS-2 ISB series was 1.32 ns, while the value associated with Galileo was 0.29 ns.
Furthermore, the FCBs of GPS, BDS-2, and Galileo were estimated and evaluated with the ambiguity a posteriori residual distributions and the STDs. The RMSs of the residuals for the NL and WL combinations among all systems were better than 0.1 cycles, indicating that relatively high-precision FCB is more likely to achieve PPP AR. For the RMSs of the residuals for the FCB estimation of the NL combination, the values were 0.067, 0.063, and 0.058 cycles for GPS, BDS-2, and Galileo, respectively, while for the WL combination, they were 0.069, 0.091, and 0.045 cycles.
Finally, rapid convergence and high positioning accuracy were achieved in the implementation of PPP AR using the multi-GNSS observations. For the ambiguity-fixed solutions, the positioning accuracy of all type of GNSS combinations was improved by over 30% compared to that of the float solutions. The RMS of the positioning accuracy for "GCE" PPP AR was 1.14 cm and it exhibited a 15% improvement compared to the result for "G". Correspondingly, the convergence time was improved to 7.5 min for "GCE", resulting in a 29% improvement. A higher ambiguity success fixing rate was also achieved for "GCE" within only 10 min of data, and the percentages were computed as 75.5%, 94.3%, 86.5%, and 99.0% for "G", "GC", "GE", and "GCE", respectively.
With the rapid development of BDS-3, its orbits and positioning performance have been demonstrated with initial evaluations [60][61][62][63][64]. FCB estimation and PPP AR performance using BDS-3 satellites will be studied in future work.