Intersystem Bias in GPS, GLONASS, Galileo, BDS-3, and BDS-2 Integrated SPP: Characteristics and Performance Enhancement as a Priori Constraints

Global navigation satellite systems (GNSSs) have experienced rapid development in recent years. Currently, both GPS and GLONASS are in full operation. As for BDS, the official declaration of positioning, navigation, and timing (PNT) services by regional BDS (BDS-2) for the users over the Asia-Pacific regions and by global BDS (BDS-3) for global users was made on 27 December 2012 and on 27 December 2018, respectively. The deployment of Galileo constellation is also progressing rapidly and will be finished in the next few years. As of August 2021, there are 31 GPS satellites, 26 GLONASS satellites, 15 BDS-2 satellites, 34 BDS-3 satellites (including four experimental satellites), and 26 Galileo satellites in orbits. The multisystem combination has been an emerging trend, as it can significantly enhance the solutions of satellite-based positioning technology, such as real-time kinematic (RTK) positioning, precise point positioning (PPP) and single point positioning (SPP). This benefits from the increased number and the improved geometry of the visible satellites. Among the satellite-based position services, SPP belongs to the open position service and is widely used in satellite navigation and low-accuracy positioning due to its flexibility, low cost, and simplicity.

Despite the many benefits of multisystem combination, there are still many issues to be solved in the multisystem integrated data processing, especially in the aspect of compatibility among different GNSSs. A bias between observations from different GNSSs can be caused by the distinct time scales employed by them. In addition, the receiver-specific hardware delays are also distinct for different GNSSs. To achieve rigorous multi-GNSS position solutions, these intersystem biases (ISBs) must be properly handled. The introduction of additional parameters to estimate the ISBs is a widely used approach. However, this approach will reduce the redundancy of observations. The observations from at least two satellites for each satellite system can theoretically contribute to the solutions of multi-GNSS integrated positioning, which will limit its performance in the extremely constrained visibility environments. When only a satellite for a specific satellite system is available in the multi-GNSS combined processing, the observation residuals of this satellite will equal to zero and the ISB estimates will encounter a deviation from the true value [1].

Many researchers focused on the mitigation of ISB in precise relative positioning. The GPS/Galileo intersystem double-differenced model based on overlapping frequencies L1/E1 and L5/E5a was derived in Odijk and Teunissen [2]. The results indicated that the between-receiver ISB could be ignored for the baseline with two receivers of identical types, while it was not the case for that characterized by different receive types. When the differential ISB was a priori corrected to recover the integer nature of double-differenced ambiguities between different satellite systems, the success rate of short-baseline ambiguity resolution (AR) with different receiver types could be comparable to the case with the same receiver types. Paziewski and Wielgosz [3] also demonstrated that the obvious ISB could not be ignored for the baselines with receivers of different types in GPS/Galileo tightly combined precise relative positioning, and they also validated the benefits for phase AR from the introduction of the known ISBs. In addition, the results indicated that the single-epoch solutions could be achieved at a noise level of 1–2 mm and decimeters in terms of phase and code ISB, respectively. Tian et al. [4] conducted the ISB estimation for the GPS/BDS-2/Galileo tightly combined RTK. For the code- and phase-specific ISB, the influencing factors included the adopted overlapping frequencies, the firmware versions, and the receiver types. Similarly, the tight combination had the improved position accuracies and successful ambiguity-fixing rates in comparison with the loose combination (differencing the observations within a single satellite system). To identify the estimability of ISB between GPS and Galileo, and GPS and BDS-2, Mi et al. [5] employed the between-receiver single-differenced (SD) multi-GNSS observation equations to formulate the ionospheric-float, ionospheric-fixed, and ionospheric-weighted models. Based on this SD method, the ISB could be estimated for both overlapping and nonoverlapping frequencies in the RTK processing. Attributed to the reasonable ISB calibration, the multi-GNSS RTK accuracy of the intersystem differencing could be improved by 20–35%, as compared to that of the classical differencing.

Actually, more studies contributed to the estimation and analysis of ISB in PPP. Liu et al. [6] investigated the effects of receiver types, antenna types, and receiver software versions as well as the sources of precise satellite products (namely analysis centers) on the stability of ISB estimates from GPS/BDS-2 combined PPP solutions. Hong et al. [7] adopted the four-system combined PPP to investigate the ISB between GPS and one of GLONASS, Galileo, and BDS-2. The ISB estimates showed systematic deviations (which were nearly identical for different stations) when using the precise products from different analysis centers. Zhang et al. [8] analyzed the long-term characteristics of ISB among GNSSs based on uncombined PPP. The standard deviations (STDs) of the relatively stable epochwise ISB over a day for GPS with respect to GLONASS, Galileo, or BDS-2 were less than 0.6 ns, while the variations of the single-day ISB solutions (average value of the epochwise ISB over a day) could be up to 60 ns, which might be attributed to the satellite reference clock of different systems. Zhou et al. [9] comprehensively evaluated the influence of ISB stochastic modeling (constant, random walk process, and white noise process) on the performance of undifferenced and uncombined multisystem (GPS/GLONASS/Galileo/BDS-2/QZSS) integrated PPP with the precise products from different analysis centers. The results indicated that both the positioning accuracy and convergence time were comparable regardless of which ISB strategy was employed when using the precise products from CODE, CNES, and WHU, while the position solutions with constant ISBs showed worse performance among the three ISB strategies when using GFZ products. With the availability of BDS-3 data, several researchers focused on its compatibility with BDS-2. Cao et al. [10] analyzed the ISB estimates between BDS-3 and BDS-2 derived from BDS-3/BDS-2 integrated PPP solutions. No significant day-boundary discontinuities could be found for the ISB estimates, but a difference of several nanoseconds was noted among the stations with the receivers from the same brand. The phase ambiguity estimates, ionospheric delay estimates, and code residuals were obviously affected when ignoring the ISB between BDS-3 and BDS-2, while it was not the case for the positioning accuracies and phase residuals. Zhao et al. [11] also demonstrated that an additional ISB parameter should be introduced into the BDS-2/BDS-3 integrated PPP processing due to the different clock references for the precise satellite clock products and the inconsistent code hardware delays between them. Similar to the precise relative positioning, several researchers investigated the contribution of a priori ISB estimates to PPP solutions. Jiang et al. [12] focused on the short-term modeling and prediction of piecewise ISB estimates (every 30 min) from GPS/BDS-2 combined PPP. The modeling accuracy was about 0.7 ns, and the prediction accuracy was 0.57–1.21 ns with the 1-day period. When taking the derived ISBs as a priori constraints, the convergence time of PPP could be reduced by 2.4–19.6%. A tightly coupled model, in which the phase-specific ISB was estimated station-by-station and the resolvable ambiguity between different satellite systems was formed, was proposed by Geng et al. [13] for multi-GNSS PPP AR. The day-to-day scattering was within 0.05 cycles for 85% of phase-specific ISBs. Compared with the intra-system PPP AR, the initialization time of inter-system PPP AR was reduced by 10%.

Several studies were related to the ISB in SPP processing. Gioia and Borio [14] employed the Allan deviation to validate the stability of ISB in GPS/Galileo combined SPP solutions. In addition, the multi-GNSS navigation performance could degrade when mitigating the ISBs (instead of estimating them) with the time offset parameters provided in the broadcast ephemeris. Zeng et al. [15] focused on the ISB between GPS and BDS-2 derived from their integrated SPP processing. The results indicated that the ISB estimates were receiver-dependent and found to be stable over a day and repeatable over adjacent days based on statistical hypothesis testing. Moreover, the estimated ISBs of BDS-2 geostationary orbit (GEO), inclined geosynchronous orbit (IGSO), and medium earth orbit (MEO) satellites with respect to GPS were consistent. Torre and Caporali [16] performed an analysis of the ISB for multisystem combined SPP using the data sets at nine European stations from four days in 2013. The results indicated that the ISB estimates between GPS and Galileo and between GPS and QZSS were dependent on the employed broadcast or precise ephemeris, while it was not the case for those between GPS and GLONASS.

According to above description, most of the research concentrated on the mitigation and estimation of ISB in precise relative positioning and PPP, while only several studies focused on the ISB in SPP. In addition, the improved accuracy of broadcast ephemeris and the increased number of in-orbit satellites in recent years may leave the conclusions on the SPP-related ISB summarized in the existed studies obsolete. Most importantly, the existed SPP studies did not cover the ISB about BDS-3 satellites. Furthermore, the contribution of a priori ISB estimates to SPP solutions was not rigorously investigated. In this study, we carefully characterize the ISB estimates in GPS/GLONASS/Galileo/BDS-3/BDS-2 integrated SPP, and the performance enhancement of ISB estimates as a priori constraints for multisystem SPP position solutions with a prediction period from a day to a week under a cutoff elevation angle from 10° to 50° (to simply simulate the harsh observation environments) is rigorously evaluated. The paper starts with the multisystem integrated SPP model taking ISB as a priori constraints. Next, the results are exhibited and analyzed. Subsequently, the discussion is conducted. Finally, the main conclusions are summarized.

In this section, the multisystem integrated SPP model taking ISB as a priori constraints is described. Only the code observations on E1 frequency for Galileo, on B1 frequency for BDS, on G1 frequency for GLONASS, and on L1 frequency for GPS are adopted here, since most of the SPP users are still using the low-cost single-frequency receivers. The single-frequency code observations from a GNSS satellite can be modeled as: where $P$ is the measured pseudorange, $\rho$ is the geometric range between the phase centers of the satellite and receiver antennas, $c\mathrm{d}{t}_r$ and $c\mathrm{d}t$ are the physical clock errors of the receiver and the satellite, respectively, $I$ is the slant ionospheric delay, $T$ is the slant tropospheric delay, and $d_r$ and $d$ are the stable code hardware delays at the receiver and the satellite, respectively.

The alignment of time reference and coordinate reference is one of key issues when performing multi-GNSS integrated data processing. The coordinate transformations can be omitted for the SPP processing with a meter-level accuracy, as the difference of coordinate references among different satellite systems (i.e., GTRF for Galileo, CGCS2000 for BDS, PZ90.11 for GLONASS, and WGS-84 for GPS) is confined to several centimeters [16]. By contrast, the difference of the time scales employed by the four GNSSs (i.e., Galileo System Time for Galileo, BDS Time for BDS, GLONASS System Time for GLONASS, and GPS Time for GPS) must be properly considered in the multi-GNSS combined SPP, in view that the effect of the time reference is significant. To mitigate the effect of distinct time scales in the multi-GNSS combined SPP, one receiver clock offset parameter can be estimated for each satellite system. Alternatively, by taking the receiver clock offset parameter of a selected satellite system as the reference, for example, GPS, the additional ISB parameters can be introduced into the code observation equations for other satellite systems. In addition to the distinct time scales, the ISB parameters can also absorb the effect of different receiver-specific code hardware delays among different GNSSs. Due to the inconsistent code hardware delays between BDS-2 and BDS-3 satellites [10,11], an ISB parameter for BDS-2 and an ISB parameter for BDS-3 are estimated. After correcting the ionospheric delay and tropospheric delay with a priori model, fixing the satellite position with broadcast satellite orbits, and correcting the satellite clock offset with broadcast satellite clocks and time group delays (TGDs), the linearized observation model for four-system integrated SPP can be expressed as: where the superscripts $G$, $R$, $E$, $C2$, and $C3$ refer to GPS, GLONASS, Galileo, BDS-2, and BDS-3 satellites, respectively, $p$ is the observed-minus-computed (OMC) code observable, $\mu$ is the unit vector in line-of-sight direction, $X$ is the receiver coordinates in three dimensions (corrections referring to the approximate values), $c\mathrm{d}{t}_{r,est}$ is the receiver clock estimate grouped with the receiver-specific code hardware delays of GPS, and ${\mathrm{ISB}}^{R,G}$, ${\mathrm{ISB}}^{E,G}$, ${\mathrm{ISB}}^{C2,G}$, and ${\mathrm{ISB}}^{C3,G}$ denote the ISB parameter between GLONASS and GPS, between Galileo and GPS, between BDS-2 and GPS, and between BDS-3 and GPS, respectively. Actually, a priority list should be made for the selection of reference satellite system for the ISB estimation. When the GPS data are absent, the receiver clock offset parameter of an available satellite system is taken as the reference based on the priority list, and then the ISB parameters are reformulated accordingly for other available satellite systems.

In view that the ISB has a stable characteristic during a short period of time, for example, several days, the estimated ISBs in advance can be used to enhance the multi-GNSS integrated SPP solutions by prediction, especially in the much harsh visibility environments. Four pseudo-observation equations are introduced into the traditional observation model of four-system integrated SPP, that is: where ${\delta }_{\mathrm{ISB}}^{R,G}$, ${\delta }_{\mathrm{ISB}}^{E,G}$, ${\delta }_{\mathrm{ISB}}^{C2,G}$, and ${\delta }_{\mathrm{ISB}}^{C3,G}$ denote the a priori ISB observable between GLONASS and GPS, between Galileo and GPS, between BDS-2 and GPS, and between BDS-3 and GPS, respectively. In this study, the a priori ISB observable is taken as the average value of the epochwise ISB estimates over a day (with a prediction period of several days).

The estimated parameters include the three-dimensional (3D) receiver coordinates, the receiver clock offset, and the four ISBs, that is: where $U$ is the vector of estimates. As shown in Equation (1), the code observation equation is nonlinear. The estimates in Equation (4) are actually the corrections with respect to the approximate values through an iterative procedure.

Besides the rigorous observation model, a proper stochastic model is also vital for the optimal SPP solutions. Since the quality of observations is highly relevant to the satellite elevation angles, the elevation-dependent weighting scheme is adopted here. Assuming that there is no correlation among the code observations from different satellites (all the nondiagonal elements are equal to zeroes in the covariance matrix of code observations), the variance of code observations from a GNSS satellite can be computed as: where ${\sigma }_0$ is the STD of code observations at zenith, which is set to 0.45 m for GLONASS satellites (mitigating the negative effect of inter-frequency code hardware delay at the receiver) and 0.3 m for Galileo, GPS, BDS-2 and BDS-3 satellites, and $\sigma (ele)$ is the STD of code observations at an elevation angle $ele$. As to the variance of a priori ISB observable, it is taken as the average value of the epochwise variances of ISB estimates over a day (with a prediction period of several days). In this study, the estimator is the least squares adjustment.

The multisystem integration is expected to ensure the positioning performance in constrained visibility environments. However, the benefits from multisystem combination for positioning and navigation are limited under extremely harsh environments due to the distinct time scales and receiver-specific hardware delays among different GNSSs. For example, the position solutions are unsolvable if the visible satellites only include a GPS satellite, a GLONASS satellite, a Galileo satellite, a BDS-2 satellite, and a BDS-3 satellite, as the additional ISB parameters need to be simultaneously estimated. The characteristic analysis indicates that the daily ISBs have a good long-term stability; thus, the predicted high-accuracy ISBs (obtained in advance) can be used as a priori constraints to improve the benefits from multisystem combination for real-time positioning and navigation under extremely harsh environments. The four-system integrated SPP considering a priori ISB constraints with a prediction period of a day can provide a positioning accuracy better than 1.6 and 3.8 m and 2.9 and 7.2 m in the horizontal and vertical directions at an availability of approximately 100.0% and 90.0% under the cutoff elevations of 40° and 50°, respectively, which is of great significance for the satellite navigation and low-accuracy positioning users.

With the rapid development of GNSSs, especially for BDS (BDS-3/BDS-2) and Galileo, the multisystem combination has been an emerging trend in recent years. The SPP technology, which is the open position service of GNSSs and is widely used in satellite navigation and low-accuracy positioning, can also benefit from the multi-GNSS integration. The compatibility among different GNSSs is a key issue in the multisystem integrated data processing. To solve this issue, the additional ISB parameters should be introduced into the observation model. In this study, we carefully characterize the ISB estimates for GLONASS, Galileo, BDS-2, and BDS-3 with respect to GPS in the four-system integrated SPP, and the performance enhancement of predicted ISBs as a priori constraints for multisystem SPP position solutions with a prediction period from a day to a week under a cutoff elevation angle from 10° to 50° (to simply simulate the harsh observation environments by increasing the elevation masks) is rigorously evaluated. The data sets from 120 globally distributed MGEX stations (covering seven different receiver types from three manufacturers) spanning a month from DOY 41 to 70 of 2020 are employed.

The results indicate that the epochwise ISB estimates are relatively stable with a varying range of approximately 10 ns, and those of BDS-3 and BDS-2 show obvious differences by several nanoseconds. The daily ISBs do not present obvious fluctuations when the receiver firmware versions are changed, but the step changes can be observed for the daily ISBs when encountering the replacement of receiver types. The short-term stability of epochwise ISBs for GLONASS, Galileo, BDS-2, and BDS-3 with respect to GPS can be 2.335, 1.262, 1.741, and 1.532 ns, respectively, whereas the corresponding long-term stability for daily ISBs can be 1.258, 1.288, 2.713, and 2.566 ns, respectively. The effects of receiver types on the short-term and long-term stability of ISBs are not significant. The daily ISBs are roughly consistent for the stations equipped with the same type of receivers. The numerical values range from 20 to 45 ns for the daily ISBs between GLONASS and GPS, from 30 to 70 ns for the daily ISBs between BDS-2 and GPS, and from 25 to 65 ns for the daily ISBs between BDS-3 and GPS, respectively, while the numerical values of daily ISBs between Galileo and GPS range from −13 to −7 ns for the stations with the TRIMBLE ALLOY receiver and from −5 to 3 ns for the stations with the other six types of receivers, respectively.

The single-day prediction residuals of daily ISBs (within the same group) for the stations with various receiver types are comparable to each other. The single-day prediction accuracy of daily ISBs for GLONASS, Galileo, BDS-2, and BDS-3 with respect to GPS can be 1.055, 0.640, 1.242, and 0.849 ns, respectively. The prediction accuracy of daily ISBs degrades with the increasing time span of prediction and a prediction accuracy of 1.767, 1.954, 2.982, and 2.580 ns can still be achieved for the four groups of daily ISBs when the time span of prediction is increased to seven days, respectively. The correction rates of predicted daily ISBs between Galileo and GPS are only 78.0% even when the time span of prediction is set to a day and decrease to 32.8% with a time span of prediction of a week, whereas the corresponding correction rates for the other three groups of daily ISBs keep over 93.0% with the increasing time span of ISB prediction.

The accuracy improvements that benefit from the a priori ISB constraints are marginal (less than 10%) when the cutoff elevations are lower than or equal to 30°. The accuracy improvements coming with the a priori ISB constraints become significant when the elevation masks are further increased to 40° and 50°. The accuracy improvements after introducing the a priori ISB constraints can be 22%, 27%, 24%, and 24% and 25%, 26%, 23%, and 23% at an elevation mask of 40° and 50° with a time span of ISB prediction of a day in the east, north, up, and 3D directions, respectively. Under a cutoff elevation of 40° and 50°, the improvement rates decrease with the increasing time span of ISB prediction, and the corresponding accuracy improvements attributed to the a priori ISB constraints with a time span of ISB prediction of seven days can still be 8%, 20%, 12%, and 13% and 12%, 17%, 6%, and 7% in the four directions, respectively. The availability starts to decrease at an elevation mask of 40° and is reduced to only 64.0% at an elevation mask of 50° for the four-system SPP ignoring a priori ISB constraints (traditional four-system SPP). By contrast, the availability still keeps 100.0% when the cutoff elevation is increased to 40° and can be up to approximately 90.0% at an elevation mask of 50° for the four-system SPP considering a priori ISB constraints.

In addition, we also characterize the ISB estimates derived from the low-cost u-blox M8T receiver and two Xiaomi Mi8 smartphones and also investigate the performance enhancement of a priori ISB constraints for the multisystem SPP solutions. Although both the short-term stability of epochwise ISBs and the single-day prediction accuracy of daily ISBs with these devices are worse than those with the geodetic-type receivers, the position accuracy as well as the availability of the multisystem integrated SPP can be improved after introducing the a priori ISB constraints, especially for the SPP results with the smartphones. Usually, the receiver hardware delay is considered to change with the temperatures. However, the epochwise ISBs obtained with the u-blox M8T receiver and the Xiaomi Mi8 smartphone do not show obvious fluctuations in the first half an hour (with increased temperatures), which may be attributed to the worse ISB estimates with the low-cost SPP receivers.