GLONASS Aided GPS Ambiguity Fixed Precise Point Positioning

Abstract

The Precise Point Positioning (PPP) concept enables centimetre-level positioning accuracy by employing one Global Navigation Satellite System (GNSS) receiver. The main advantage of PPP over conventional Real Time Kinematic (cRTK) methods is that a local reference network infrastructure is not required. Only a global reference network with approximately 50 stations is needed because reference GNSS data is required for generating precise error correction products for PPP. However, the current implementation of PPP is not suitable for some applications due to the long time period (i.e. convergence time of up to 60 minutes) required to obtain an accurate position solution. This paper presents a new method to reduce the time required for initial integer ambiguity resolution and to improve position accuracy. It is based on combining GPS and GLONASS measurements to calculate the float ambiguity positioning solution initially, followed by the resolution of GPS integer ambiguities.

The results show that using the GPS/GLONASS float solution can, on average, reduce the time to initial GPS ambiguity resolution by approximately 5% compared to using the GPS float solution alone. In addition, average vertical and horizontal positioning errors at the initial ambiguity resolution epoch can be reduced by approximately 17% and 4%, respectively.

1. INTRODUCTION

The state-of-the-art method to obtain centimetre-level positioning accuracy is conventional Real Time Kinematic (cRTK), which is based on cancelling or reducing errors that are correlated between two Global Navigation Satellite System (GNSS) receivers by differencing measurements across the receivers and satellites. This method provides centimetre-level positioning accuracy with short convergence times. However, its main drawback is that it requires the use of a local reference network infrastructure or a single reference station, which may not always be available. In the case of cRTK, the baseline between receivers must typically be shorter than 50 km when employing dual-frequency GNSS receivers. However, in the case of high ionospheric activity, ionospheric de-correlation can be significant even with shorter baselines, leading to poor positioning accuracy. The cost of building dense reference networks can be high and potentially prohibitive in remote areas (Kaplan and Hegarty, 2006).

With Precise Point Positioning (PPP), a local reference infrastructure is not required (Zumberge et al., 1997). This leads to cost savings when high accuracy positioning is needed, for example, in remote areas. Instead of correcting errors by differencing measurements across receivers and satellites, they are corrected by employing correction products generated from data captured by global reference networks and, where appropriate, modelling the errors. The major error sources for PPP are the ionosphere, troposphere, satellite clock, and satellite orbit and receiver clock. In addition, site-displacement, antenna phase centre offset, antenna phase centre variation and phase wind-up errors should be corrected to obtain high positioning accuracy (Abdel-salam, 2005).

There are already many possible applications of PPP including land surveying, precision farming or crystal deformation monitoring. In general, using PPP is practical if there is no local reference network infrastructure available and the long PPP solution convergence time is not a problem for the application (Bisnath and Gao, 2007).

The traditional PPP model employs the ionosphere-free GPS L1/L2 observable to eliminate the first-order ionospheric error, which can cause metre-level range errors (Heroux et al., 1993). The tropospheric errors are corrected by employing, for example, the Saastamoinen model and mapping function (Saastamoinen, 1973). The International GNSS Service (IGS) precise satellite orbit and clock products (IGS, 2009) are used to correct satellite orbit and clock errors. The site-displacement effects are caused by the periodic movements of stations compared to the International Terrestrial Reference Frame (ITRF) (Kouba, 2009a). They can be corrected using site-displacement models to estimate correction terms and adding these to position estimates. The antenna offsets and variations can be corrected using corrections in the Antenna Exchange Format (Kouba, 2009a). The satellite antenna phase wind-up correction can be calculated based on the phase centre coordinates of the receiver and satellite antennae (Wu et al., 1992).

Carrier-phase ambiguities are estimated as float numbers. According to Bisnath and Gao (2007), it may take 30 minutes, or more, depending on the processing scenario, to obtain static 3D position solution accuracy better than 10 cm.

A key issue with traditional PPP models is that they do not fix carrier-phase ambiguities to integers because of Fractional Cycle Biases (FCB) or Un-Calibrated Phase Delay (UPD) errors in GNSS carrier-phase signals. This issue has recently been addressed by Geng et al. (2010) through the FCB estimation method and Collins (2008) and Laurichesse et al. (2011) through the Integer Recoverable Clock (IRC) method. To date these two methods are the main fixed ambiguity PPP approaches for the GPS-only scenarios. The principle of the FCB estimation method is to use a global reference network to estimate narrow-lane FCB corrections. The narrow-lane, which has a wavelength of 10·7 cm, is a combination of GNSS L1 and L2 measurements (Collins, 2008). The calculation of it is explained later in this paper. These FCB corrections are applied at the user-level when narrow-lane ambiguity resolution is carried out. In contrast, in the IRC method, FCB corrections are not needed for narrow-lane ambiguity resolution because a carrier-phase reference network solution is employed when estimating satellite clock corrections.

The main weakness of the current GPS fixed ambiguity PPP methods above is the long time period (30 to 60 minutes) required to obtain an initial ambiguity fixed position solution (Geng et al., 2010; Laurichesse, 2011). Other issues include: lack of an appropriate ambiguity validation method and of integrity monitoring, and insufficient performance for centimetre-level position accuracy. These additional issues are addressed by the GPS-only PPP model presented in Jokinen et al. (2012). However, the length of time required to obtain the initial ambiguity fix is still an open issue (Jokinen et al., 2012).

It is shown in Jokinen et al. (2011) that employing both GPS and GLONASS can reduce the PPP float solution convergence time compared to GPS alone in some cases. Reussner and Wanninger (2011) showed that GLONASS ambiguity resolution in the case of PPP could be used, in a geometry-dependent way, by using external ionosphere corrections. According to their work, solving GLONASS ambiguities is not possible in a geometry-free way as in the case of fixed ambiguity GPS PPP, because receiver specific code biases make it impossible to use the Melbourne-Wubbena combination to estimate wide-lane float ambiguities (Melbourne, 1985, Wubbena, 1985). Furthermore, Reussner and Wanninger (2011) showed that calibrating the code biases is difficult because of receiver specific variations. In addition, calibration of GLONASS inter-channel phase biases is required (Reussner and Wanninger, 2011). For these reasons, solving GLONASS carrier-phase ambiguities in the case of PPP is a difficult task. On the other hand, employing both GLONASS and GPS measurements when estimating a float position solution is straightforward and should provide performance benefits.

This paper builds on the work on the combined usage of GPS and GLONASS measurements for PPP to propose a new method for employing GLONASS to facilitate the resolution of GPS integer ambiguities. Section 2 describes the GPS/GLONASS PPP method developed, focussing on error models and corrections, as well as GPS ambiguity resolution and validation. Section 3 presents the test results and Section 4 concludes the paper.

2. PPP METHOD

2.1. Error models and corrections

The French Space Agency (CNES) satellite orbit and clock corrections are employed both for GPS and GLONASS satellites (CNES, 2012). There are no FCB errors in the GPS satellite clock corrections because they are generated by employing a carrier-phase network solution (Laurichesse, 2011). In contrast, the GLONASS satellite clock corrections are generated in the usual way by employing a code-phase network solution (CNES, 2012). The CNES satellite clock and orbit corrections are generated in real time. Therefore, the results presented in this paper, although post-processed, give a good indication of real-time (epoch by epoch) positioning performance.

The standard deviation of the CNES GPS satellite clock corrections is 0·12 ns, which is equivalent to 3·6 cm, and the 3D accuracy of CNES GPS orbit corrections is better than 5 cm on average (Laurichesse, 2011). The satellite orbit error in the radial direction is significantly smaller than the 3D orbit error (Laurichesse et al., 2008). According to Laurichesse (2011), the accuracy of CNES GPS orbit and clock corrections is sufficient to enable GPS narrow-lane ambiguity resolution. For the CNES GLONASS satellite orbit and clock corrections, accuracy estimates are not available in the public domain. Therefore, it is assumed that the accuracy of GLONASS orbit and clock corrections is worse than of GPS orbit and clock corrections.

The antenna phase centre offsets and variations can cause decimetre-level errors and therefore must be corrected in PPP (Abdel-salam, 2005). Both the receiver and satellite antenna phase centre offsets and variations are therefore, corrected by employing corrections in the ANTEX 1.4 format provided by IGS (IGS, 2009).

According to Kaplan and Hegarty (2006), the magnitude of the total tropospheric delay can typically vary between 2·4 m and 25 m, depending on the satellite elevation. The tropospheric delay is corrected using the UNB3m troposphere model (Leandro et al., 2006) and the Global Mapping Function (GMF) (Boehm et al., 2006). The dry component is corrected directly by employing the UNB3m model. The wet component is mapped to the range domain by employing GMF, and estimated as a state in an Extended Kalman Filter (EKF). Tropospheric gradients to the north and east directions are also estimated as states. The Chen mapping function as shown in (1), is used to map tropospheric gradient to the range domain (Chen and Herring, 1997).

(1)$$m(\varepsilon )_{azi} = \displaystyle{1 \over {\sin \varepsilon \,\cos \varepsilon + 0{\cdot} 0032}}$$

where ε is the elevation angle of the satellite.

The total tropospheric delay correction is calculated as:

(2)$$d_{trop} = m(\varepsilon )_h d_h + m(\varepsilon )_w d_w + m(\varepsilon )_{azi} (G_N \cos \phi + G_E \sin \phi) $$

Where:

ϕ

is the azimuth angle of the satellite

m(ε)_h

is the tropospheric dry delay mapping function

d _h

is the tropospheric dry delay

m(ε)_w

is the tropospheric wet delay mapping function

d _w

is the tropospheric wet delay

m(ε)_azi

is the tropospheric gradient mapping function

G _N

is the tropospheric gradient value to the north direction and

G _E

is the tropospheric gradient value to the east direction.

C1/P1 and C2/P2 Differential Code Biases (DCB) are biases between C1 and P1 or C2 and P2 code-phase measurements. When C1 or C2 measurements are used, the biases must be corrected because CNES precise orbit and clock corrections refer to P1 and P2 measurements. The correction of the biases is done using products from the Centre for Orbit Determination in Europe (CODE) (Dach et al., 2007). DCB corrections are provided as monthly files because they are fairly stable. For example, the approximate day by day reproducibility of C1/P1 DCB is 1·5 cm (RMS). In the case of PPP, correcting DCB is particularly important when employing the Melbourne-Wubbena combination (see section 2.2.1) to estimate wide-lane ambiguities, because the impact of DCBs can be up to 1·2 m or 1·4 wide-lane cycles. (Dach et al., 2007).

The satellite phase wind-up correction is calculated based on the receiver antenna and satellite antenna phase centre coordinates (Wu et al., 1992).The wind-up effect can cause errors smaller than one cycle. Satellites which are in the eclipsing phase are excluded because it is not possible to model their yaw attitude accurately (Kouba, 2009b).

The first order ionospheric error is eliminated by employing the ionosphere-free measurement combination, with higher order terms ignored due to their relatively small contributions to the range error. The higher order ionospheric error terms can cause an error of a few centimetres for satellites at low elevation angles, if the ionosphere is highly active (Dach et al., 2007).

To obtain centimetre-level positioning accuracy, site displacement effects need to be corrected (Dach et al., 2007). In particular, solid earth tides can cause significant position errors, typically 30 cm in the radial direction and 5 cm in the horizontal direction (Kouba, 2009a). Solid earth tides, ocean loading and pole tides are corrected by employing models presented in Kouba (2009a). The atmospheric pressure loading effect can occasionally cause centimetre-level vertical errors (Urquhart, 2009), which are considered insignificant and therefore ignored in this paper.

Multipath and Non-Line-Of-Sight (NLOS) errors as in Kaplan and Hegarty (2006) are dealt with using elevation-based measurement weighting and detection using integrity monitoring.

Cycle slip detection is carried out before using measurements to the position estimation, and is based on detecting changes in the ionospheric Total Electron Content (TEC) Rate (TECR) and the Melbourne–Wubbena Wide Lane (MWWL) ambiguity between the current and previous measurement epochs (Liu, 2011). When there is no cycle slip, both TECR and MWWL differences between the epochs are small.

If a cycle slip occurs, a cycle slip correction method, mainly based on work presented in Banville and Langley (2009), is applied. The magnitude on the L1 and L2 measurements is estimated based on the change of geometry-dependent wide-lane ambiguity and TECR between the current and previous epochs. Cycle slips are corrected based on these estimated values. The geometry-dependent wide-lane is used for cycle slip correction instead of MWWL because the geometry dependent wide-lane has smaller multipath and noise. Cycle slip corrections are validated by checking TECR and geometry-dependent wide-lane difference between the current and previous epochs. The cycle slip correction is rejected if the validation test is not passed.

2.2. GPS ambiguity resolution and validation

The Between-Satellite-Difference (BSD) operation is used to remove receiver clock error and receiver side FCB from the measurements by differencing the GPS code-phase measurements, as shown in equation (3), and carrier-phase measurements, as shown in equation (4). The highest elevation satellite is chosen as a base-satellite and the measurements of the base-satellite are subtracted from the measurements of other satellites. The BSD was previously used in the FCB estimation by Geng et al. (2010). In this paper, BSD measurements are used together with the IRC fixed-ambiguity PPP method. Un-differenced GLONASS measurements are used.

The terms and notation used in the formulation of the relevant subsequent mathematical expressions are:

F

is the frequency band index of the received GPS or GLONASS signals. The signal can be either on the L1 (index=1) or L2 (index=2) frequency band

f

is the frequency of the GNSS signal in Hertz

f ₁

is GPS or GLONASS L1 frequency and

f ₂

is GPS or GLONASS L2 frequency

C _F

is a GPS or GLONASS C-code observation on the frequency band F in metres

P _F

is a GPS or GLONASS P-code observation on the frequency band F in metres. In the case of receivers providing C-code observations, P-code observations can be obtained by employing DCB corrections

L _F

is a GPS or GLONASS carrier-phase observation on the frequency band F in metres

j

is the index of a base satellite

i

is the index of a satellite (i≠j)

(3)$$P_F^{ij} = \; P_F^i - \; P_F^j $$

(4)$$L_F^{ij} = \; L_F^i - \; L_F^j $$

The ionosphere-free combination is calculated by using (5) for code-phase measurements (P ₃) and (6) for carrier-phase measurements (L ₃) in metres. Ionosphere-free ambiguity terms are estimated as Extended Kalman Filter (EKF) states. The principle of Kalman filtering as originally presented in Kalman (1960) is to minimise Gaussian noise in the system. However, it is not possible to solve for ionosphere-free ambiguities directly because these ambiguities have a non-integer nature (Geng et al., 2010). This can be seen in (7) where the ionosphere-free ambiguity term (b _c) is calculated based on the L1 ambiguity (N ₁) and L2 ambiguity (N ₂) terms. The value of b _c is not an integer, even when N ₁ and N ₂ are fixed to integers. Therefore, GPS ambiguity resolution is carried out in two stages: the wide-lane ambiguity is fixed in the first stage and the narrow-lane ambiguity is fixed in the second stage. GLONASS ambiguities are kept as float.

(5)$$P_3 = \; \displaystyle{1 \over ({\,f_1^2 - f_2^2) (\,f_1^2 P_1 - f_2^2 P_2 )}}$$

(6)$$L_3 = \; \displaystyle{1 \over ({f_1^2 - f_2^2) (f_1^2 L_1 - f_2^2 L_2 )}}$$

(7)$$b_c = \; \displaystyle{1 \over ({\,f_1^2 - f_2^2) (\,f_1^2 N_1 - f_1 f_2 N_2 )}}$$

2.2.1. Wide-lane ambiguity resolution and validation

The Melbourne-Wubbena combination in (8) is employed to estimate GPS wide-lane ambiguity float values (Melbourne, 1985; Wubbena, 1985). The CNES wide-lane FCB corrections are applied to the estimates to recover the integer nature of the wide-lane ambiguities (Laurichesse, 2011). Thereafter, the estimates are fed to the EKF to filter out measurement noise. Using the EKF instead of averaging Melbourne-Wubbena wide-lane ambiguities can be beneficial because wide-lane ambiguities are correlated with each other. The correlation is due to the BSD operation on the measurements. It is important to apply the same wide-lane corrections on the client side as used in the generation of the IRC satellite clock corrections, otherwise wide-lane ambiguities could be fixed to the wrong integers.

(8)$$L_6^{i,j} = \displaystyle{{\,f_1 L_1^{i,j} - f_2 L_2^{i,j}} \over {\,f_1 - f_2 \;}} - \displaystyle{{\,f_1 P_1^{i,j\;} + f_2 P_2^{i,j}} \over {\,f_1 + f_2}} $$

Wide-lane ambiguity resolution is achieved by rounding the float wide-lane ambiguity estimate to the nearest integer value. The wide-lane integer ambiguity candidates are validated by using the probability-based test presented in Dong and Bock (1989). With this test, an ambiguity validation decision is based on the distance from the nearest integer and the standard deviation of the float wide-lane estimate. It is possible to use this simple method to fix wide-lane ambiguities because the wavelength of the wide-lane ambiguities is long (86 cm) compared to the typical error budget.

2.2.2. Narrow-lane ambiguity resolution with the minimum constellation method

The narrow-lane float ambiguity estimate (b_nl^i,j) shown in (9) is calculated based on the float ionosphere-free ambiguity estimate (b_c(i, j)), from (7), and the fixed wide-lane ambiguity (N_wl^i,j) Geng et al., 2010).

(9)$$b_{nl}^{i,j} = \displaystyle{{\,f_1 + f_2} \over {\,f_1}} b_c^{i,j} {\rm \;} - {\rm \;} \displaystyle{{\,f_2} \over {\,f_1 - f_2}} N_{wl}^{i,j} $$

The narrow-lane integer candidate vector is obtained by using the Least-squares AMBiguity De-correlation Adjustment (LAMBDA) method (Jonge and Tiberius, 1996). With this method, the integer ambiguity least squares search is carried out by using a transformation that de-correlates the ambiguities. At least four narrow-lane float ambiguities are required before the LAMBDA method can be used. The Minimum Constellation Method (MCM) developed in Schuster et al. (2012) and adopted in Jokinen et al. (2012) is used when there are more than four float ambiguities available.

The principle of the MCM is to attempt ambiguity resolution with all possible sub-set combinations of float ambiguities if the full-set cannot be used. For example, if there are five narrow-lane float ambiguities available, ambiguity resolution can be attempted for one five-satellite group and five four-satellite groups. The benefit of employing the MCM is the achievement of fast initial ambiguity resolution by exploiting measurement sub-sets that have the required quality and geometry conditions. (Jokinen et al., 2012).

To reduce the amount of computational power required, narrow-lane ambiguity groups consisting of low elevation satellites (below 15°) or narrow-lane ambiguity groups with high Position Dilution Of Precision (PDOP) can be excluded when there are more than nine narrow-lane ambiguities available. Nine is empirically selected as the critical value for the number of combinations to test because the required computation time increases significantly when the number of ambiguity combinations increases. Nine is chosen as the threshold because the data processing would otherwise take too long to be practical. For example, there are 126 four ambiguity groups to test when there are nine float narrow-lane ambiguities available. The number increases to 210 when there are ten float narrow-lane ambiguities available.

2.2.3. Narrow-lane ambiguity validation

Narrow-lane ambiguity validation is carried out using the ratio test in (10), where R₂ is the Sum of Squared Errors (SSE) of ambiguity residuals of the second most probable ambiguity candidate vector. R₁ is the SSE of the most probable ambiguity candidate vector. These vectors are calculated as shown in (11),

(10)$$\displaystyle{{{\bf R}_2} \over {{\bf R}_1}} \gt threshold$$

(11)$${\bi R_i} = (\hat a - \check{a} _i )^T G_{\hat a}^{ - 1} (\hat a - \check{a}_i )$$

where: â is the float ambiguity vector

$ {\check a}_i $ is an integer ambiguity candidate vector and

G_â is the variance matrix of the float ambiguity vector

The fixed ambiguity candidate vector is accepted if ${\textstyle{{{\bf R}_2} \over {{\bf R}_1}}} $ is larger than the test threshold. The ratio test does not test the absolute correctness of the fixed ambiguity candidate vector. Instead it provides information on the closeness of the integer ambiguity candidate vector to the float ambiguity vector. Therefore, the integer ambiguity candidate vectors accepted could be incorrect if the float ambiguity vector is incorrect, for example, as a result of multipath or other measurement errors. (Teunissen and Verhagen, 2008)

An empirically chosen constant ratio test threshold value such as 1·5 or 3·0 is used in many RTK and PPP software packages. However, there is no credible theoretical or practical justification for the use of a constant threshold. For higher reliability ambiguity validation, the ratio test threshold can be calculated based on the required confidence level and the degrees of freedom. (Feng et al., 2012)

It is assumed that the residuals of integer ambiguity candidate vectors obey the non-central chi-square distribution shown in (12) where:

• n is the number of degrees of freedom (i.e. the number of float narrow-lane ambiguities in the candidate vector)

• δ _i is the non-central parameter and

• i is the frequency index (i.e. either 1 or 2).

Both the numerator and denominator of the ratio test statistics obey the non-central chi-square distribution shown in (13). Therefore, the ratio test statistics obey the double non-central F-distribution shown in (14), where (δ ₁,δ ₂) are the parameters of the F-distribution, if it is assumed that R₁ and R₂ are independent. (Feng et al., 2012)

(12)$$\displaystyle{{R_i} \over {n\sigma ^2}} = \; \chi (n,\delta _i )$$

(13)$$\displaystyle{{{\bf R}_2} \over {{\bf R}_1}} = \; \displaystyle{{{\bf R}_2 /(n\sigma ^2 )} \over {{\bf R}_1 /(n\sigma ^2 )}}$$

(14)$$\displaystyle{{{\bf R}_2} \over {{\bf R}_1}} \sim \; F(n,n,\delta _2, \delta _1 )$$

The ratio test threshold value can be obtained by solving numerically the integral shown in (15) (Feng et al., 2012). The calculation is based on the number of degrees of freedom (n) and the required confidence level (P _c). In the case of a practical software implementation, there is no need to calculate this integer as shown in (15) in real-time. A table of ratio test thresholds for different combinations of confidence levels and degrees of freedom can be included in the software.

(15)$$P_c = \; \mathop \int \limits_0^{\displaystyle{{{\bf R}_2} \over {{\bf R}_1}}} F(x|n,n;\delta _2, \delta _1 )dx$$

2.2.4. Making wrong ambiguity resolution less likely

As presented in the previous section, the ratio test provides information about the closeness of an integer ambiguity candidate vector to the float ambiguity vector. However, it does not validate the correctness of the float ambiguity vector.

Figure 1 shows an example of wrong ambiguity resolution, which occurs when the float position solution is still converging. In this case, the ratio test alone is used with a confidence level threshold of 99·99% to validate the ambiguities. The values of fixed and float ambiguities at the initial ambiguity resolution epoch are shown in Table 1. It can be seen that the float ambiguity values are close to integers. Thus, the ratio test accepts the fixed ambiguity candidate vector with the test statistic value of 24·93, when the threshold is 24·06. However, the float ambiguity values are far from the correct values, because the float position solution is still converging. The correct float ambiguity values correspond to when the magnitude of the position error is zero. However, it not possible to know the exact correct float ambiguity values when using real data. Thus, it is assumed that float ambiguities are far from the correct values when the position error is large (3D error more than 10·7 cm). In summary, this case shows that wrong ambiguity resolution may occur when using the ratio test alone and the float ambiguities are far from the correct values.

Figure 1. 3D position error in the case of fixing ambiguities wrongly.

Table 1. Ambiguity values in a wrong ambiguity resolution case.

To reduce the likelihood of a wrong ambiguity resolution, additional validation is required to complement the ratio test.

When employing MCM, there may be more than one ambiguity candidate vector accepted by the ratio test. Therefore, the best ambiguity candidate vector has to be determined, based on the following criteria established in this paper: the largest number of fixed ambiguities in the vector and the highest ratio test statistic. Thereafter, it is required that the best ambiguity candidate vector is the same over a given number of epochs consecutively, and that the ambiguities belonging to the vector can be fixed to the same integers during the chosen time-period. In this paper, 30 epochs are used as an empirically chosen threshold in the decision to accept or not the ambiguity resolution. The theory behind this test is that the float ambiguity values are changing when the float solution is converging. Thus, if for 30 epochs (150 s with 5 s data-rate) the best ambiguity candidate vector is the same and ambiguities in the vector are fixed to the same values, then it is less likely that float ambiguities have not converged to the correct values and are close to the nearest integer. The value of 30 epochs is selected based on testing ambiguity resolution with different values such as 10, 20, 30 and 40 and selecting the value which gives a good compromise between the reliability and rapidity of ambiguity resolution.

With the LAMBDA method, the number of float narrow-lane ambiguities in a vector must be at least four. However, in this paper at least five ambiguities are required. This is because it is less likely that five ambiguities are close to wrong integers than four ambiguities when the float position solution is far from the correct position. The probability that all float ambiguities are close to wrong integers decreases as the number of the ambiguities tested increases.

A five degree elevation mask is used in this paper. However, ambiguity resolution is only attempted for an ambiguity term when the elevation angle for the satellite is at least ten degrees. The reason for this selection is that it is beneficial to use low elevation satellites in the calculation to improve tropospheric error estimation. However, ambiguity resolution for low elevation satellites, particularly for satellites below 10 degrees, is risky because of the increased likelihood of multi-path errors.

When employing the confidence level-based ratio test, the ratio test statistic must be larger than the threshold for the integer ambiguity candidate vector to be accepted. A confidence level threshold value of 99·99% is used during the float solution convergence period and 99·9% is used otherwise. The solution is defined to be in the convergence period when the longest carrier-phase lock time is smaller than 2000 seconds. The reason for using the higher ratio test confidence threshold during the float solution convergence period is that the position error is typically larger during this period. Thus, it is more likely that float ambiguities are close to wrong integers.

The choice of ambiguity validation parameters is a trade-off between the rapidity and reliability of ambiguity resolution. For example, it is possible to obtain ambiguity resolution in a short period of time (e.g. less than 500 seconds) under favourable conditions if using less strict ambiguity validation parameters than used in this paper. However, the probability of wrong ambiguity resolution is then also higher.

2.2.5. Calculating the ambiguity fixed position solution

After the integer ambiguity candidate is accepted by the ratio test, fixed ionosphere-free ambiguity (b _c^i,j) values are calculated based on the fixed wide-lane and narrow-lane ambiguities as shown in (16) (Geng et al., 2010). It is not possible to fix ionosphere-free ambiguities to integer values but the ionosphere-free ambiguities can be fixed to float values. In this paper, an initial fixed ambiguity position solution is defined as a fix of at least four ionosphere-free ambiguities.

(16)$$b_c^{i,j} \; = \; \displaystyle{{\,f_1} \over {\,f_1 + f_2}} (N_{nl}^{i,j} ) + \displaystyle{{\,f_1 f_2} \over {\,f_1^2 - f_2^2}} N_{wl}^{i,j} $$

2.2.6. The benefit of adding GLONASS

In this paper, both GPS and GLONASS measurements are used when estimating the float position solution. Using GLONASS measurements in addition to GPS measurements can improve float position solution convergence in many cases as presented in Jokinen et al. (2011). The reasons for the improvement are better satellite geometry and the larger number of measurements available, which can help, for example, tropospheric error estimation. The narrow-lane ambiguity resolution is largely affected by the accuracy of the float position solution. Therefore, ambiguity resolution performance can be improved by employing GLONASS measurements in those cases where the use of GLONASS improves the quality of the float position solution.

2.3. The Extended Kalman Filter (EKF)

The principle of the Kalman filter was originally presented in Kalman (1960). Kalman filters are commonly applied in navigation to filter out measurement noise, which is typically assumed to be Gaussian distributed. In this paper, an EKF similar to that in Feng et al. (2009) is employed. The estimated EKF states are latitude, longitude, altitude, GLONASS receiver clock error, tropospheric wet delay, tropospheric gradients to the east and north directions, GPS and GLONASS ionosphere-free ambiguities and GPS wide-lane ambiguities.

The EKF is updated separately with BSD GPS and un-differenced GLONASS measurements. The standard deviation of GPS carrier-phase measurements is obtained based on (17), where a is 0·003 m, b is 0·003 m and el is the elevation angle of the satellite. The formula and the values of the GPS measurement standard deviations are similar to those used in the RTKLIB software (Takasu, 2012). It is assumed that standard deviations are 100 times larger for GPS code-phase measurements than for GPS carrier-phase measurements. It is also assumed that the quality of GLONASS orbit and clock corrections is worse than the GPS clock and orbit corrections. Thus, the magnitude of the GLONASS measurement standard deviation is set as larger than for GPS measurements. The standard deviation of GLONASS carrier-phase measurements is set at three times the GPS carrier-phase measurements. For GLONASS code-phase measurements, the standard deviation is set at 20 times the GPS code-phase measurements. In addition, there are no separate code-phase and carrier-phase satellite clock corrections for GLONASS and there can be inter-channel code-biases in GLONASS code-phase measurements. This is another reason why GLONASS code-phase measurements are significantly underweighted against GPS code-phase measurements.

(17)$$std = a + \displaystyle{b \over {\sin (el)}}$$

3. TEST RESULTS

The effect of using both GPS and GLONASS compared to using GPS alone is analysed. The first aim is to test the effect of using a combined float position solution when attempting GPS ambiguity resolution against a float position solution obtained by GPS only.

All tests presented in this paper are carried out with the iNsight project (http://www.insight-gnss.org/) POINT software. The POINT software is programmed in C/C++ and is flexible in its configuration to test different positioning scenarios.

3.1. Test data

Static GNSS data recorded by IGS is used to test ambiguity resolution algorithms. Originally recorded at 1 Hz the data were resampled to 5 seconds to match the rate of the satellite clock correction products. A five degree satellite elevation mask is used for these tests. Data from 12 different stations recorded on 10 different days (30–31 December 2011, 2–6 January 2012 and 10–12 February 2012) are used. The codes of the IGS stations used are: BRST, HERT, KIR0, LAMA, MARS, MATE, SASS, SOFI, TITZ, UNB3, WTZR and ZIM2.

3.2. Ambiguity resolution

The average times required to obtain an initial ambiguity resolution with GPS-only and the combination of GPS and GLONASS are computed and compared. Initial ambiguity resolution is defined as the fixing of at least four integer narrow-lane carrier-phase ambiguities. A position error at the initial ambiguity resolution epoch is obtained by comparing coordinates estimated by the POINT software to the known International Terrestrial Reference Frame (ITRF) 2008 coordinates of the stations. For the ZIM2 and TITZ stations, the known coordinates are provided by the GNSS Data Centre (GDC) as in GDC (2012). For the remainder of the stations, the known coordinates are provided by ITRF (2012). The estimated accuracy of the station coordinates provided by GDC and ITRF is at the millimetre level, enabling the investigation of the accuracy of the position solutions provided by POINT.

Figure 2 shows the average time required to obtain ambiguity resolution on different days. Results suggest that combining GPS with GLONASS reduces, on average, the time required to obtain the initial ambiguity resolution compared to the GPS-alone scenario on all ten days. Depending on the day, the magnitude of reduction is between 15 and 230 seconds. Figure 3 shows the average time needed for initial ambiguity resolution at different stations. It can be seen that the average time is lower at ten different stations in the GPS/GLONASS case. The largest reduction of 589 s is obtained when processing data from the SASS station. The station specific differences in the average time may depend on the receiver and antennae used, as well as the local tropospheric conditions. It can be seen from Figure 4 that using GLONASS with GPS results in an initial ambiguity resolution time 5% shorter when averaged over all stations and days tested. In general, it can be seen that using both GPS and GLONASS gives a small reduction in the time required to fix ambiguities. Although the reduction is small it is beneficial given that adding GLONASS to the processing can be done at little additional cost, because most modern professional GNSS receivers already support GLONASS. The standard deviations of the time required to obtain an initial ambiguity resolution are shown by error bars. Ambiguities can be fixed for all tested days and stations when RINEX data are available, using both the GPS-only and GLONASS with GPS processing models.

Figure 2. Average time for initial ambiguity resolution for different days over all stations (standard deviation shown by bars).

Figure 3. Average time for initial ambiguity resolution at different stations over all tested time periods (standard deviation shown by error bars).

Figure 4. Average time for initial ambiguity resolution for all tests (standard deviation shown by error bars).

The average horizontal position errors on different days are shown in Figure 5 and for different stations in Figure 6. In terms of the average horizontal position error, using GPS with GLONASS gives a smaller or equal average error for the data on eight of the days analysed (between 0 and 0·3 cm) and for nine stations (between 0 and 0·2 cm). The largest increase of 0·2 cm to the average horizontal position error caused by using both GPS and GLONASS occurs on 31 June 2011. The horizontal position error is larger in the MARS, SASS, WTZR and ZIM2 stations when employing both GPS and GLONASS. There may be problems with GLONASS orbit and clock corrections for some of the satellites. Thus, the horizontal error is larger in some stations. In general, the difference in the average horizontal error in the different processing cases is small (between 0 and 0·3 cm).

Figure 5. Average horizontal position errors – different days – all stations (standard deviation is shown by error bars).

Figure 6. Average horizontal position errors – different stations – all time-periods (standard deviation shown by error bars).

The average vertical position errors on different days are shown in Figure 7. It can be seen that using GLONASS with GPS reduces the average horizontal error on all the days tested, varying between 0·1 cm and 1·4 cm. The vertical position error at different stations is shown in Figure 8. It can be seen that using GPS with GLONASS gives smaller or equal vertical errors in all stations except WTZR. In the WTZR station, using GPS with GLONASS causes a 0·4 cm increase in the average vertical error compared to using GPS alone. There are many possible reasons which may cause the increase: the worse quality of GLONASS orbit and clock correction products, accuracy of site-displacement corrections, partially wrong ambiguity resolution, incorrect tropospheric modelling and incorrect input parameters for the tropospheric model. In the other stations, using GLONASS with GPS reduces the average vertical error between 0 and 1·5 cm.

Figure 7. Average vertical position errors on different days over all stations (standard deviation shown by error bars).

Figure 8. Average vertical position errors at different stations over all tested time-periods (standard deviation shown by error bars).

Figure 9 shows the average 3D, horizontal and vertical position errors based on all data tested. It can be seen that using GLONASS with GPS reduces the average 3D error by 11%, horizontal error by 4% and vertical error by 17%.

Figure 9. Average 3D, horizontal and vertical errors based on all data tested (standard deviation is shown by error bars).

In general, it can be seen that using GLONASS with GPS gives benefits for fixed-ambiguity PPP processing by reducing the time required for the initial ambiguity resolution, and by increasing accuracy. The average reduction in time to initial ambiguity resolution is 5%, in the horizontal error (at initial ambiguity resolution epoch) 4% and in the vertical error (at initial ambiguity resolution epoch) 17%, based on the all data tested. It is likely that smaller vertical errors are obtained due to the increase in the number of satellites, which helps tropospheric error estimation, when using both GPS and GLONASS. On the other hand, the horizontal errors are not reduced as much as the vertical errors because the increased number of satellites due to adding GLONASS does not provide a significant benefit to the horizontal error estimation. In addition, the quality of the GLONASS satellite orbit predictions may not be as good as the quality of the GPS satellite orbit predictions, which can cause a negative impact to the horizontal error estimation.

Using both GPS and GLONASS for processing is particularly useful in terms of reducing the vertical position error. The reduction of the time needed to obtain the initial ambiguity resolution and horizontal position error is less significant. However, the cost of adding GLONASS to the processing is small, because most modern professional GNSS receivers can already provide both GPS and GLONASS measurements. Thus, the method of aiding GPS PPP ambiguity resolution with GLONASS is suitable for common use.

4. CONCLUSIONS

A novel PPP method which uses GLONASS to aid GPS ambiguity resolution is presented in this paper. Between-Satellite-Difference (BSD) GPS measurements and un-differenced GLONASS measurements are used. The resolution of GPS carrier-phase ambiguities is attempted and GLONASS carrier-phase ambiguities are kept as float numbers. The minimum constellation method is used to reduce the time required to obtain an initial ambiguity resolution. Ambiguity validation is carried out by employing the ratio test with a confidence level based variable threshold and additional validation methods, presented in this paper.

PPP ambiguity resolution is tested for two cases: GPS only and the combination of GPS with GLONASS. The results show that using both GPS and GLONASS can help reduce the time needed to achieve initial GPS ambiguity resolution by approximately 5% and reduce vertical and horizontal position errors at the initial ambiguity resolution epoch by approximately 17% and 4%, respectively.

Future work will investigate PPP with new GNSS signals. These may help reduce the time needed for initial ambiguity resolution. In addition, future work can include the developing and testing of integrity monitoring algorithms.