Precision analysis of partial ambiguity resolution-enabled PPP using multi-GNSS and multi-frequency signals

A single-receiver integer ambiguity resolution-enabled precise point positioning (PPP-RTK) user experiences a long convergence time when the rather weak single-constellation dual-frequency ionosphere-ﬂoat model is used. Nowadays, the rapid development of Global Navigation Satellite Systems (GNSS) provides a multitude of available satellites and frequencies that can serve in improving the user’s model strength and, therefore, its ambiguity resolution and positioning capabilities. In this study, we provide insight into and analyze the global impact of a multi-GNSS (GPS, Galileo, BeiDou-3) multi-frequency integration on the expected ambiguity resolution and positioning performance of the ionosphere-ﬂoat uncombined PPP-RTK user model, and demonstrate whether it is the increased number of satellites or frequencies, or a combination thereof, that speeds up ambiguity-resolved positioning. Moreover, we explore the capabilities of both full (FAR) and partial (PAR) ambiguity resolution, considering the full ambiguity information content with the LAMBDA method, and investigate whether PAR is an eﬃcient solution to the multi-dimensional ambiguity case. The performance of our solutions is assessed in terms of the ambiguity success rate (ASR), the number of epochs (TTFA) to achieve both an ASR criterion and a horizontal positioning precision better than 10 cm, as well as the gain in precision improvement. Based on multi-system multi-frequency simulated data from nine globally distributed stations and a large number of kinematic solutions over a day, we found that the increase in number of frequencies enhances the ambiguity resolution performance, with PAR achieving a TTFA reduction of 70% when ﬁve instead of two Galileo frequencies are used, while the ambiguity-ﬂoat solution is only slightly improved. Further, our numerical results demonstrated that the increase in number of satellites leads to an improvement in both the positioning and ambiguity resolution performance, due to the improved geometry strength. It is shown that the GPS+Galileo+BeiDou solutions can achieve a TTFA of 6.5 and 4.5 min (at 90%) on a global scale when two and three frequencies are used, respectively, without any a priori information on the ionospheric delays. Finally, we analyzed the sensitivity of the PPP-RTK user’s performance to changes in the precision of the measurements. (cid:1) 2020 COSPAR. Published by Elsevier Ltd. This is an open access article under the CC BY license (http://creativecommons.org/licenses/by/4.0/).


Introduction
Integer ambiguity resolution (IAR) is the key to fast and precise Global Navigation Satellite System (GNSS) parameter estimation (Teunissen, 1995;Teunissen and Montenbruck, 2017). The purpose of resolving the carrier-phase ambiguities as integers is to gain a significant improvement upon the precision of the remaining model https://doi.org/10.1016/j.asr.2020.08.010 0273-1177/Ó 2020 COSPAR. Published by Elsevier Ltd. This is an open access article under the CC BY license (http://creativecommons.org/licenses/by/4.0/). parameters, with the position components being usually the main parameters of interest. This improvement takes place because one can take full advantage of the carrier-phase data that act as ultra-precise pseudoranges, once the integer ambiguities are resolved successfully.
The reliability of ambiguity resolution, usually expressed by the probability of correct integer estimation or the so-called ambiguity success rate (ASR), depends on the applied integer estimator and on the strength of the underlying observation model. Assuming that the selected method for integer estimation is the one that maximizes the success rate, it is then the model strength that purely drives the ASR. Given that there are no misspecifications in the applied functional and stochastic models, one can use the model-driven ASR as a diagnostic measure for the expected IAR performance. This measure is solely based on the precision of the least-squares estimated float ambiguities, which is captured by the variancecovariance (vc-) matrix of the ambiguities. Therefore, the stronger the model, the more precise the ambiguities and the higher the ASR will be. However, it is well known that the single-constellation dual-frequency ionosphere-float model is weak in the sense of its IAR capabilities due to the presence of ionospheric delays. It is expected that additional observations, e.g. due to inclusion of data over multiple frequencies and from multiple GNSS systems or a priori atmospheric information, will strengthen the model and lead to higher ASRs.
In this contribution, we analyze the performance of IAR in the context of multi-system multi-frequency integration in the IAR-enabled precise point positioning (PPP; Zumberge et al. (1997)) method, namely PPP-RTK (Wü bbena et al., 2005), and in particular the gain in position precision improvement one should expect.
The rapid development and modernization of multiple satellite systems along with a plethora of new signals and frequencies provide an improved satellite geometry, thus stronger positioning model, and higher redundancy than the until recently traditional GPS dual-frequency model. It is, thus, expected that both positioning and ambiguity resolution capabilities will be improved in light of these developments. On the other hand, such an integration also leads to an increased dimension of the vector of to-beresolved ambiguities, implying a probable decrease of the ASR, mainly because of rising satellites, and an IAR slowdown. In case of a strong model, either by accumulating data over many epochs or by integrating multi-system and/or multi-frequency data, it might not be necessary to resolve the complete vector of integer ambiguities, but instead a sufficiently large subset such that the gain in position precision is significant. We will, therefore, explore the capabilities of both full (FAR) and partial (PAR) ambiguity resolution, and investigate whether PAR is an efficient solution for such a high-dimensional problem.
Diverse studies have shown the benefit of using multiconstellation and/or multi-frequency data to obtain more precise and reliable solutions, claiming, also faster conver-gence in the absence of precise ionospheric corrections. Using simulated GPS and Galileo dual-frequency data, Tiberius et al. (2002) have shown that the combinedsystem model can achieve instantaneous IAR for short and medium baselines, unlike long baselines, with formal ASRs above 99.9% and 95.0%, respectively. They also found that the inclusion of a third frequency brings a only a slight improvement in the performance, which was later observed by Verhagen et al. (2007), Odijk et al. (2014a) and Xiao et al. (2019).
According to Zhang et al. (2003), the mean time-to-fixambiguities (TTFA) is about 70 and 35 s for 50 km baselines in GPS-only and Galileo-only triple-frequency solutions, respectively, using simulated 1 Hz data. Moreover, Ji et al. (2013) showed that it takes over 39 s to resolve the ambiguities for long baselines using simulated Galileo 1 Hz data on four frequencies with cascade ambiguity resolution. After analyzing triple-frequency GPS and Galileo 30 s data from a single baseline, Odijk et al. (2014a) showed that the mean TTFA is about 32 and 25 min for FAR and PAR, respectively, in a GPS-only dual-frequency mode, while a slight improvement was observed when using a third frequency. In the dual-system solution, the mean TTFA was about 7 min for FAR and 2 min for PAR. Based on real five-frequency Galileo data in Australia, Wang et al. (2018) concluded that a high number of frequencies is helpful to achieve high ASRs within a short time at the network side, which can be higher than 99.9% within the first 5 epochs of processing using all five frequencies. Li et al. (2018) analyzed the IAR performance of static PPP using dual-frequency GPS, BeiDou, GLONASS and Galileo data from a global network of stations. Based on solely wide-lane combinations, it was found that the four-system solution enables a TTFA of 10 min compared to the GPS-only solution which required 18 min on average to achieve successful IAR. Using real triple-frequency GPS, Galileo and BeiDou 1 s measurements from 17 stations around the world and wide-lane ambiguity combinations, Duong et al. (2020) found that the average ASR using the multi-system multi-frequency kinematic PPP model was about 15% higher than the triple-frequency GPS-only counterpart, with the former showing a TTFA of 199 s and the latter of 553 s. In addition, according to the experimental results for receivers mainly located in the Asia-Pacific area, Li et al. (2019) concluded that the average TTFA of a combined BeiDou and Galileo dual-frequency static PPP solution can be reduced by using the third frequency from 34 to 29 min.
Despite different combinations of systems and frequencies, underlying models, sampling rates and receiver locations, the above reported results are indicative of the great importance of a multi-system multi-frequency integration for IAR, in terms of higher ASRs, reduced TTFAs and shorter positioning convergence times. However, it has been a common practice to use the traditional TTFA as a sole indicator of the PAR capabilities, see e.g. Duong et al., 2020;Li et al., 2019), without considering whether there is a gain in precision improvement or a criterion on the position precision has been met, as done for instance by Odijk et al. (2014a) and Brack (2017). In addition to this, the effect of a high-dimensional ambiguity vector on the ASR has not been studied satisfactorily yet, and is addressed in our study by providing numerical insight into the PAR capabilities as a solution to the dimensional curse. Moreover, we believe that a thorough insight into the expected multi-dimensional ambiguity resolution performance has not been yet provided, since most of the recent studies are restricted to static positioning applications and/or to the characterization of the performance in local or regional areas, see (Odijk et al., 2014a;Li, 2018;Li et al., 2018;Deo and El-Mowafy, 2019;Liu et al., 2019;Xiao et al., 2019;Li et al., 2019). An analysis of the IAR performance on a global scale and under kinematic mode would provide potential multi-system multifrequency users with a baseline indicator regarding the performance level they should expect.
In this contribution, we study and assess the expected performance of multi-GNSS (GPS, Galileo, BeiDou-3) multi-frequency ambiguity resolution on the basis of the ionosphere-float uncombined PPP-RTK user functional and stochastic models, thus without the need to process or simulate real data. Based on realistic functional and stochastic assumptions, a formal analysis allows us to gain clear insight into the factors contributing to ambiguity resolution and conclude on whether the improved satellite geometry or the high number of frequencies, or a combination thereof, is mainly driving the ambiguity resolution and positioning performance. Compared to other studies, the kinematic PPP-RTK user performance is predicted based on globally distributed stations and on a large number of samples. Instead of forming a priori linear combinations of multi-frequency system-specific ambiguities, as commonly used in recent studies, we use the Least-Squares AMBiguity Decorrelation Adjustment (LAMBDA) method (Teunissen, 1995) to determine the best-resolvable ambiguity (subset) combinations without any loss of information. We close our analysis by a systematic comparison of the performance of different models in terms of the TTFA, the achieved position precision and the expected gain in position precision after IAR.
This contribution is organized as follows. In Section 2, we briefly review the theory of integer ambiguity resolution, the success rate and the partial ambiguity estimator chosen in this study. Then, Section 3 introduces the basis of our uncombined ionosphere-float PPP-RTK underlying model, its estimable parameters, the experimental setup and processing settings. This is then followed, in Section 4, by an analysis of the numerical results for various multisystem multi-frequency models using both FAR and PAR. Finally, the work is concluded in Section 5.

Integer ambiguity resolution
In this section, a brief review of the principles of integer ambiguity resolution is given, along with a discussion on the selected PAR estimator.

Mixed-integer GNSS model
Any carrier-phase based GNSS model can be cast into the following general linear(ized) system of observation equations: where EðÁÞ and DðÁÞ denote the expectation and dispersion operators, respectively; y the vector of uncombined observed-minus-computed (O-C) carrier-phase and code observations, a 2 Z n the vector of integer doubledifferenced carrier-phase ambiguities, b 2 R q the vector of real-valued parameters such as position components, clocks, atmospheric delays and hardware delays, A and B their respective partial design matrices, and Q yy the vcmatrix of the observations. This mixed-integer GNSS model is solved mainly in three steps. In the first step, the best linear unbiased estimation takes place, ignoring the integer property of the carrier-phase ambiguities, resulting in the so-called float solution of all parameters (with:-symbol): whereâ andb are the float estimators of the ambiguity and the real-valued parameters, respectively. Qââ; Qbb and Qâb ¼ Q T bâ are the corresponding (co-) variance matrices. The second step focuses on the integer constraint a 2 Z n . Its objective is to map the float ambiguity solutionâ into the integer solution a ¼ IðâÞ using an integer mapping I : R n # Z n from the n-dimensional space of reals to the n-dimensional space of integers. Various integer estimators can be selected for this step, with the most popular ones being the following: integer rounding (IR), integer bootstrapping (IB) and integer least-squares (ILS). ILS is the optimal method for integer estimation in the sense that it maximizes the probability of correct integer estimation, i.e. the success rate (Teunissen, 1999). The ILS estimator, a ILS ¼ arg min a2Z n kâ À ak 2 Qââ , is efficiently mechanized in the LAMBDA method, in which the real-valued ambiguities are transformed and decorrelated with the Z-matrix, such that the ambiguity search space turns from an elongated hyper-ellipsoid into spheroid-like, enabling a fast integer search in the transformed search space: whereẑ and Q^z^z denote the transformed float ambiguities and their vc-matrix, respectively. Once the integer outcomes a ¼ Z ÀT z are accepted, the ambiguity-float solution of the real-valued parameters,b, is corrected by virtue of their correlation with the ambiguities, obtaining the so-called ambiguity-fixed solution (with :symbol): b ¼b À QbâQ À1 aâ ðâ À aÞ ð 5Þ

Ambiguity success rate
The ambiguity-fixed solution will enjoy a precision that is in accordance with the high precision of carrier-phase data, due to the imposed integer ambiguity constraints. However, this is based on the assumption that the ILS integer solution corresponds to the correct solution. In any other case, a wrong integer solution can cause significant biases in the position solution that may exceed the error of the ambiguity-float solution.
The probability of correct integer estimation is driven by the chosen integer estimator and the precision of the float ambiguity solution, Qââ, which depends on the strength of the underlying model at hand. Therefore, to infer whether the integer outcomes can be reliably used, one requires a diagnostic measure in order to decide on their acceptance or rejection. In this study, we use the formal bootstrapping success rate as a measure for successful ambiguity resolution, which lower bounds the success rate of the optimal ILS estimator (Teunissen, 1999), since an exact easy-to-compute expression is available: with UðxÞ denoting the cumulative normal distribution function: where P ð z ILS ¼ zÞ and P ð z IB ¼ zÞ denote the ILS and IB ASRs of the decorrelated ambiguities z, respectively, n the number of the decorrelated ambiguities, and expðÁÞ is the natural exponential function. rẑ ijI represents the standard deviations of the i th decorrelated ambiguities conditioned on the previous I ¼ i þ 1; . . . ; n ambiguities, and are provided by the square roots of the entries of the diagonal matrix D after an L T DL-decomposition of the decorrelated ambiguity vc-matrix Q^z^z. The reason why we evaluate the ASR based on the decorrelated ambiguities is that the IB estimator depends on the ambiguity parametrization and, by applying the Z-transformation of LAMBDA, the IB ASR becomes a sharp lower bound to the ILS ASR (Verhagen et al., 2013). If the ASR is high and close to 1, the integer ambiguities can be assumed to be deterministic and, then, the vc-matrix of the ambiguity-fixed parameters from Eq. (5) turns into:

Partial ambiguity resolution
The capability to resolve the full set of ambiguities with a high ASR is not always feasible, as the observation model might not be strong enough. It is well known that the presence of ionospheric delay parameters in the ionospherefloat model makes the latter weak in terms of instantaneous (single-epoch) FAR. In such a case, one would need to accumulate data over several epochs to ensure that a high ASR is achieved. In addition, the rising of new satellites would require even more epochs to be accumulated, due to the introduction of estimated float ambiguities of low precision that bring down the ASR because of the multiplicative definition of the IB ASR.
The same situation occurs in a multi-system multifrequency integration, where the number of to-beresolved ambiguities can be very large. This can lead to a dimensional curse, because the event that each of the ambiguities is correct should have a probability close to 1 (Verhagen et al., 2012). As the dimension of the ambiguity vector increases, this probability tends to get smaller due to the multiplication of probabilities, which by definition are smaller than or equal to 1.
Despite the incapability to rapidly resolve the complete ambiguity vector in these cases, it may still be possible to resolve a subset of ambiguities, referred to hereafter as PAR. Several PAR methods have been proposed in literature, including e.g. the fixing of (extra) wide-lane ambiguities in multi-frequency models (Cao et al., 2007;Li, 2010), the fixing of ambiguities that are identical in the LAMBDA-based best and second-best solution (Dai et al., 2007;Lawrence, 2009), or the fixing of only the ambiguities that have been individually accepted based on fixedfailure-rate critical values after FAR (Brack and Gü nther, 2014).
In this contribution, we describe and use the modeldriven approach of Teunissen et al. (1999), which is easy to implement and defines the ambiguity subset to be resolved based on a minimum required success rate. The method commences with the decorrelation of ambiguities using the LAMBDA Z-transformation and then selects the largest possible subset, starting from the last decorrelated ambiguity, such that the user-defined success rate requirement is met: where P 0 is the minimum required success rate and z nÀkþ1 is the subset containing the last n À k þ 1 decorrelated ambi-guities. The reason why one starts from the n-th decorrelated ambiguity and continues until selecting the last n À k þ 1 entries is that the LAMBDA algorithm orders the entries of the decorrelated ambiguity vector in ascending order in terms of precision. This means that only the last n À k þ 1 decorrelated ambiguities,ẑ nÀkþ1 , will be fixed to their integers using an integer estimator, while the remaining subset,ẑ kÀ1 , will be kept as float. In case k ¼ 1, then the selected subset is identical to the full set and, thus, PAR would coincide with FAR. Once the n À k þ 1 ambiguities are fixed to their integers in the transformed ambiguity domain, the partially ambiguity-fixed solution of the real-valued parameters follows as: It can be seen from Eq. (10) that the improvement of the real-valued parameters after successful PAR is based on the decorrelated integer ambiguity subset. Unless k ¼ 1, the back-transformation of the full set of decorrelated ambiguities z PAR ¼ ½ẑ T kÀ1jK¼k;...;n ; z T nÀkþ1 T into the original set a PAR will not contain integer entries anymore, since these entries will be linear functions of all decorrelated ambiguities that contain both the integer-valued and the conditional real-valued ambiguities.
In general, the PAR-based solution will be less precise than the FAR-based counterpart, at least until the identified subset corresponds to the full set, but still more precise than the ambiguity-float solution. This does not necessarily imply that the PAR solution will be significantly more precise than the float one, since this is inextricably linked to the number of fixed ambiguities. Therefore, the availability of a PAR solution should not be confused with the availability of a high-precision position solution.

PPP-RTK processing strategy
This section provides a detailed description of the observational model and parameter estimability for the PPP-RTK concept, defines the performance measures and gives an overview of the data and experimental setup considered in our study.

Full-rank observation model
The observation equations for the single-system singleepoch uncombined carrier-phase (/ s H r;j ) and code (p s H r;j ) measurements between the receiver r and the satellite s H of system H on frequency j can be formulated as (Teunissen and Montenbruck, 2017): Eðp sH r;j Þ ¼q sH r þðdt r À dt sH Þþm sH r s r þ l j i sH r þðd H r;j À d sH ;j Þ where j ¼ 1; . . . ; f is the frequency index with f being the number of frequencies, and s H ¼ 1 H ; . . . ; m H is the system-specific satellite index with m H being the number of satellites per H-system, with H 2 fG; E; Cg. The letters G; E and C denote the GPS, Galileo and BeiDou systems, respectively. In the following, we make systematic use of the satellite index s H to discriminate between the satellites of different GNSSs and keep a generalized form of the equations. The term q s H r denotes the geometric receiversatellite range. The receiver clock and system-specific satellite clock parameters are represented by dt r and dt s H , respectively. s r and m s H r represent the zenith tropospheric delay parameter and the troposphere mapping function, respectively. The slant ionospheric delay for a receiver r and a system-specific satellite s H is denoted by i s H r and is linked to the observations through the ionospheric coefficient l j for frequency j. The system-specific receiver and satellite phase biases are denoted with d H r;j and d s H ;j , respectively, while d H r;j and d s H ;j denote those for code observations, respectively. The integer phase ambiguity is represented by a s H r;j and is linked to the phase data through the wavelength k j at frequency j. The phase biases and ambiguities are expressed in cycles, while the other parameters in units of range.
In this contribution, by means of uncombined measurements we refer to measurements that have not undergone any differencing and/or linear combinations. Using such a formulation provides the flexibility of having all parameters available for a possible further model strengthening and of easily extending the model to any number of frequencies (Odijk et al., 2016).
In our PPP-RTK network processing, the carrier-phase and code measurements from Eq. (11) are corrected for the receiver and satellite positions (lumped in the geometric range), assuming that they are both precisely known. This network system of equations is rank-defect as the information content is not sufficient to determine the absolute parameters, but only estimable functions of them. The rank-deficiencies of the network model can be identified and removed by defining a minimum set of S-basis parameters according to the S-system theory (Baarda, 1973;Teunissen, 1985). After reformulation, the O-C terms of the uncombined phase (D/ s H r;j ) and code (Dp s H r;j ) measurements for a full-rank single-system multi-frequency model read as: where r ¼ 1; . . . ; n with n being the number of network receivers, and p denotes the pivot receiver/satellite depending on whether it is used as subscript/superscript. The interpretation of the estimable parameters, denoted with theÁ-symbol, and the chosen S-basis are given in (f > 2). After bringing together and applying these corrections along with the orbital corrections at the observational level, the user-corrected single-system linearized observation equations (replacing the index r by the user index u) for the O-C phase (D/ s H u;j ) and code (Dp s H u;j ) data follow as: EðDp sH u;j þ dt sH þd sH ;j>2 Þ ¼g where g s H u is the receiver-satellite direction vector and Dx u is the position increment vector. The interpretation of the estimable parameters follows from the user version of those in Table 1, with r replaced by u. Note that the user estimable ambiguities,ã s H u;j ¼ a s H pu;j À a p H pu;j , are now in doubledifferenced form and thus integer.
Assuming that the network corrections are sufficiently precise, the stochastic model of the single-epoch singlesystem ionospheric-float user model is given as: where T denote the phase and code measurement vectors, respectively, per GNSS system. The terms r D/ u;j and r Dp u;j denote the zenith-referenced formal precision of the phase and code data, respectively, while the m Â m matrix W u ¼ diagðw 1 u ; . . . ; w m u Þ contains the weights for every satellite. The symbol denotes the Kronecker product.
The above model formulation can be easily extended when multiple GNSS systems are employed, as the rankdeficiency removal concept is applicable in the same manner for every H-system. In a multi-system integration, one has to be aware that the receiver code and phase biases are not experienced in the same way from system to system in common frequencies (Hegarty et al., 2004;Montenbruck et al., 2011;Odijk and Teunissen, 2013) 1 . This is the reason why the full-rank system model Eq. (12) results in system-specific estimable receiver code/phase biases and receiver clock offsets. One could take advantage of the overlapping frequencies among systems and parameterize the network model in terms of the inter-system biases (ISB) that can be estimated and then used for strengthening the multi-system PPP-RTK user model as shown in Khodabandeh and Teunissen (2016). In this contribution, we chose to parameterize the system of equations in system-specific parameters treating each H-system independently, with only the coordinate and troposphere parameters being common for all GNSSs. A consequence of this approach is that one pivot satellite should be taken per system, and not a common one across systems. The variance-covariance (vc) matrix of the single-GNSS measurements is extended by appending the measurement vc-matrices of the additional systems in a block-diagonal manner.

Performance measures
In our study, we use several indicators to analyze the expected positioning capabilities of multi-system multifrequency user models and show the impact of ambiguity resolution. Using the recursively estimated parameter solutions as a basis for our measures, we evaluate the formal ASR, the TTFA, the formal position precision and the precision gain after successful IAR. It is worth mentioning at this point that there is no need of simulating the code and phase data themselves or the integer ambiguities that are derived thereof, because our computations are based on a formal analysis making use of the model's design matrix and the measurement vc-matrix.
The formal ASR is based on Eq. (6) and is used to infer whether the ambiguities can be reliably fixed to their Table 1 Estimable parameters and chosen S-basis of the system-specific multi-frequency network model (the symbol p denotes the pivot satellite/receiver if it is used as superscript/subscript). correct integers. Its estimation is solely based on the vcmatrix of the estimated float ambiguities, Qââ, which is available at every epoch of our Kalman-filter-based processing. When FAR is attempted, the TTFA demonstrates the required time span to exceed a minimum ASR criterion, which was set equal to 99.5% in our computations. However, TTFA in PAR would be the required time to meet the ASR criterion but at the same time to achieve a horizontal positioning precision of better than 10 cm. For our ambiguity resolution computations, we used the LAMBDA software that has been recently released in Python (Psychas et al., 2019). The positioning precision along the north, east and up components was derived from the corresponding entries of Qbb which were available at every epoch, with the ambiguity-fixed counterparts from Q b b being available only after successful IAR. The contribution of IAR into the positioning domain was measured by means of the gain numbers (Teunissen, 1998): where f k¼1;2;3 are the gain vectors along the three directions. In our study, we analyze the precision gain only in the horizontal (radial) position.

Experimental setup
To analyze the ambiguity resolution and positioning performance of the ionosphere-float PPP-RTK user model based on multiple combinations of systems and frequencies, we selected 9 globally distributed International GNSS Service (IGS; Dow et al. (2009)) stations for our simulation, the locations of which are shown in Fig. 1. Due to differences in the satellite visibility from site to site, the selection of user stations across the globe allows us to get a clear insight into the expected performance, that is globally applicable and not location-restricted. The receiversatellite geometries for GPS, Galileo and BeiDou-3 were reconstructed based on the precisely known station coordinates and the satellite positions, derived from the IGS merged broadcast ephemeris files 1 , on a randomly chosen day of 2019. To conduct a global-scale numerical analysis, we used only the BeiDou-3 Medium Earth Orbit (MEO) satellites, as they are the ones that provide global service and can transmit triplefrequency signals.
The real-time parameter estimation in the multi-GNSS multi-frequency PPP-RTK user models was performed with the recursive minimum-mean-squared-error Kalman filter using a dynamic model. Moreover, in our analysis the user positioning is performed in kinematic mode, therefore treating the user's position components as unlinked parameters in time. Details concerning the user's data, filter and processing settings are given in Table 2. It is worth mentioning that we chose the same code and phase zenith-referenced standard deviation for all GNSS systems and frequencies. We also assumed that the network corrections, provided to the user for realizing single-receiver IAR, are deterministic quantities, while no correlation between frequencies or in time was assumed. Moreover, the ASR criterion we used to decide on whether or not the ambiguities have been reliably fixed to their correct integers was set to 99.5% for both FAR and PAR.
Finally, we emphasize that our numerical results are linked to the selected sampling rate of 30 s. An increase in the sampling rate will have a positive effect on the achieved performance. Although a low sampling rate is beneficial for ambiguity resolution due to the greater change in the receiver-satellite geometry from epoch to epoch, a higher sampling rate leads to a higher model strength within the same time span.

Results and analysis
In the following section, we will present the performance results achieved with the ionosphere-float Kalman-filtered PPP-RTK user model using multiple GNSSs and multiple frequencies on a multitude of user stations.

Single-station analysis
This subsection provides an initial analysis of the expected ambiguity resolution and positioning performance in a single user station, namely DLF1 (in the Netherlands), in order to get a first insight into the factors that contribute mostly to IAR and positioning. Fig. 2 shows the number of observed GPS, Galileo and BeiDou-3 satellites as tracked by DLF1 on the selected day. It can be easily observed that the number of satellites for a combined GPS+Galileo+BeiDou system almost triples compared to a single GNSS separately. The singleand multi-GNSS results are shown in Sections 4.1.1 and 4.1.2, respectively. Fig. 3 presents the GPS-only ambiguity resolution and positioning performance results, using both FAR and PAR, during the first 2 h (240 epochs) on DOY 274 of 2019. The results indicate that 56 epochs since the start are needed to exceed the ASR criterion of 99.5% in order to resolve the full set of ambiguities and achieve centimeter-level positioning. Despite the gain that the user experiences in the position precision compared to the ambiguity-float solution, the ASR shows a sharp decrease at the very next epoch that restricts the availability of a FAR-based solution. This is due to a newly tracked GPS satellite at the 57th epoch. We recall here that the formal bootstrapped ASR is defined as a multiplication of probabilities of correct integer estimation, with the one of the newly tracked satellite being poorly determined. This, practically, means that with every introduction of a new satellite into our filter solution and with the attempt to resolve the full set of ambiguities, the ASR will show large fluctuations and will not meet the ASR criterion consistently. This is where the prominent advantage of PAR lies in, since the PAR technique attempts to fix a subset of ambiguities, instead of the full vector, that meets the ASR criterion. In this way, the user has an automatic method of discarding ambiguities of poor precision that are usually linked to newly tracked satellites, thereby solving the effect of rapid fluctuations due to rising satellites. The ASR criterion for PAR is met in 4 epochs in this case. However, this does not necessarily mean that the partially-ambiguityfixed solution achieves centimeter level position precision, as it can be observed from the figure, since a large enough subset is needed for that. The results show that the PAR solution achieves a 10 cm horizontal position precision after 39 epochs (19.5 min) where 75% of the decorrelated ambiguities have been reliably fixed. The maximum precision gain of about 6 is observed a few epochs later, when the full set of ambiguities has been reliably fixed, as expected. At this epoch, the FAR-based and PAR-based results are identical, with an achieved accuracy of 2 cm. Based on this result, we conclude that a GPS-only dualfrequency user will not experience the ultimate cm-level precision in an instant with PAR, but a gradual improvement compared to the ambiguity-float solution, as has also been shown by Odijk et al. (2014a). Despite the continuous availability of IAR in the next epochs, it can be seen that the precision gain decreases exponentially towards the value 1, which marks the point where the ambiguity-float solution shares the same quality with the ambiguity-fixed counterpart. This is expected because our ambiguity-float PPP solution is based on an implicit accumulation of data  over epochs due to the employed Kalman filter, which allows the ambiguity-float solution to become more precise over time.

Single-GNSS performance
Further, our intention is to investigate whether it is the increased number of frequencies or satellites, or a combination thereof, that contributes mostly to the ambiguity resolution and positioning performance. We start our analysis with the introduction of the third GPS frequency, namely L5. Fig. 4 depicts the GPS-only triple-frequency results of DLF1 during the same time interval. Two main conclusions can be drawn in this case. First, it can be observed that there is only a slight improvement in the ambiguity-float position solution by incorporating the L5 measurements. The reason that the triple-frequency model provides almost the same ambiguity-float positioning performance with the dual-frequency model lies in the fact that the receiver-satellite geometry remains invariant when data in more frequencies are used, since the number of observed satellites remained constant.
The impact of the multi-frequency integration is visible, though, in the ambiguity-fixed solution. It can be seen that the TTFA for FAR is shortened compared to the dualfrequency case and is equal to 47 epochs (23 min), while PAR can instantaneously meet the ASR criterion. Despite the PAR availability, it can be seen that only a small subset of ambiguities was reliably fixed at the first epochs, with the maximum gain of 8.9 being acquired at the 38th epoch where 89% of the full decorrelated ambiguity vector has been fixed. Thus, after a time span of 38 epochs (18.5 min) the user enjoys a 2 cm position precision using PAR, which in case of FAR it would be delayed by 10 epochs (5 min). It can be concluded, therefore, that although the receiver-satellite geometry is not improved by using data in more frequencies, IAR is successfully achieved in a shorter time span that leads to highprecision positioning solutions.

Multi-GNSS performance
The combined GPS+Galileo dual-frequency ambiguity resolution and positioning performance results of station DLF1 are illustrated in Fig. 5. In this case, the initial TTFA for FAR is 77 epochs (38 min) and, as expected, longer than the GPS-only dual-frequency TTFA. This is mainly due to multiplicative nature of the bootstrapped ASR, because in a multi-system integration the dimension of the ambiguity vector increases and, therefore, the probability of correct integer estimation has the tendency to get smaller.
However, an improved positioning performance can be seen compared to the GPS-only dual-and triplefrequency solutions, due to mainly two reasons. First, the user enjoys an ambiguity-float position solution of higher precision that reaches the 10 cm level in 34 epochs (16.5 min), which is due to the improved satellite geometry after using both GPS and Galileo satellites in our processing. The second reason lies in the improved PAR performance. Despite the higher dimension of the ambiguity vector, it is observed that within only 14 epochs (6.5 min) the PAR strategy was able to identify 71% of all the decorrelated ambiguities that meet the pre-defined ASR criterion. At this time instance, the positioning precision experiences a significant gain of about 16, with the precision improving from about 24 to 1.5 cm. If one attempts FAR, a time span of 77 epochs (38 min) is needed.
To understand the underlying reason of the better PAR performance in the GPS+GAL solution compared to the GPS-only solution, we did an inspection of the conditional standard deviations (STDs) of the ambiguities at the 14th epoch before and after decorrelation, which are shown in Fig. 6. Looking at the conditional STDs before decorrelation, one can see that there are several discontinuities in the GPS+GAL ambiguity spectrum, with the number of conditional ambiguities of a few centi-cycles STD being larger than this of the GPS-only spectrum. Therefore, the presence of satellite redundancy and of many more small conditional STDs in the GPS+GAL solution allowed, after the Z-transformation, to pull down more large STDs than in the GPS-only case. And as it can be seen, a subset of 23 decorrelated ambiguities was identified in the combined solution that meets the ASR criterion and achieved a 1.5 cm precision, while the subset of only 5 ambiguities in the GPS-only solution did not improve the ambiguityfloat solution.
The GPS+Galileo triple-frequency results are shown in Fig. 7. It is shown that the incorporation of the third frequency brings only a marginal improvement in the ambiguity-float position precision, due to the invariant satellite geometry, that reaches 10 cm in 32 epochs (15.5 min) compared to 34 epochs (16.5 min) in the dualfrequency case. On the other hand, a significant shortening of the TTFA for FAR is observed, which reduces from 77 epochs (38 min) to 18 epochs (8.5 min) using the additional L5/E5a measurements. In analyzing the PAR-based results, it was found that 86% of decorrelated ambiguities are identified as the best-resolvable ambiguity subset within 11 epochs to meet the ASR criterion, compared to 14 epochs in dual-frequency case. By resolving this subset, the user's position precision increases from 28 to 1.5 cm, with a gain of about 19. Therefore, it is again demonstrated that the increase in number of frequencies slightly improves the ambiguity-float performance, but improves the IAR performance, which in turn improves the ambiguity-fixed positioning precision.
To further investigate the integration of more than 2 systems, we evaluated the ambiguity resolution and positioning performance of the GPS+GAL+BDS dualfrequency PPP-RTK user model at station DLF1, with the results being illustrated in Fig. 8. It can be observed that the ambiguity-float position precision reaches the 10 cm level in 29 epochs (14 min) which is an improvement of 15% compared to the GPS+GAL dual-frequency solution. The improvement is smaller compared to the addition of GAL satellites in the GPS-only solution, as the GPS +GAL model was already strong enough. The TTFA for FAR is 76 epochs (37.5 min) and is therefore almost equal to the one of the GPS+GAL dual-frequency model, while the PAR solution achieves a 1.3 cm horizontal position precision within the first 10 epochs (4.5 min) having 70% of the decorrelated ambiguities fixed and a gain of 21.5. Therefore, using the BDS measurements reduces the convergence time to 10 cm by 2 min.
The triple-GNSS triple-frequency solution is shown in Fig. 9. As it was previously concluded, it can be seen that the addition of the third frequency of BDS does not contribute in improving considerably the ambiguity-float position precision. Then, although the TTFA shows a decrease compared to the triple-GNSS dual-frequency solution, it does not improve compared to the dual-GNSS triplefrequency solution. It is seen, therefore, that the addition of an additional frequency has a great impact in reducing the FAR TTFA when the receiver-satellite geometry remains invariant. The partially-ambiguity-fixed solution surpasses the 10 cm level and reaches a position precision of 1.5 cm within the first 7 epochs (3 min), when 83% of ambiguities have been fixed and the corresponding gain is 23. Compared to all the previous solutions, the  triple-GNSS triple-frequency solution was able to achieve centimeter-level positioning in the shortest time span, equal to 3 min.

Global analysis
In this section, we present and analyze the multi-system multi-frequency user performance on a global scale. To conduct a rigorous analysis of the global IAR performance, one needs a sufficient number of sample solutions in order to infer realistic conclusions about the expected performance. To this end, we performed our formal analysis, based on kinematic PPP-RTK user positioning, for the 9 selected globally distributed IGS stations with a 2-h processing window being re-initialized every 1 min, in order to capture all possible receiver-satellite geometry changes, having in total 1320 sample solutions per model per station. The measures we used in this case to characterize the performance are the achieved positioning precision and expected gain after IAR. Thus, we sorted the horizontal position precision of all solutions for every epoch and identified the maximum precision that do not exceed 90% of all the sorted solutions, while the same procedure was followed for the gain numbers as well. Since it was shown that PAR performs better than FAR and can provide centimeter-level positioning results in a significantly shorter time span, we considered only PAR in this analysis.
The average number of observed satellites for GPS, Galileo and BeiDou-3 in the selected day is depicted in Fig. 10. The average number of tracked satellites in a single-system scenario is 9 for GPS, followed by 7 for Galileo and 5 for BeiDou. When satellites from multiple systems are combined, the number of visible satellites ranges from 18 to 23 on a global scale, leading to improved geometry for the users.
The convergence behavior and TTFA of the multisystem multi-frequency user positioning results, along with the precision gain after IAR, are illustrated in Fig. 11. The reported results are, therefore, indicative of the performance one should expect at 90% of the cases on a global scale.
From the single-GNSS results, it is clear that the GPS solutions can achieve better performance compared to the Galileo and BeiDou counterparts, using either two or three frequencies, which is due to the complete GPS constellation with 9 observed satellites on average, unlike Galileo and BeiDou with 7 and 5 observed satellites on average, respectively. Moreover, the GPS-only solution shows a TTFA of 25.5 min using PAR. Due to their incomplete constellations, the Galileo-only solution shows a TTFA of 83 min, while the BDS-only solution requires more than 2 h to resolve a sufficiently large enough subset to achieve high precision. Further, it can be seen that the addition of a third frequency in the single-constellation models does not improve the ambiguity-float position precision, as expected, which is due to the invariant satellite geometry. However, the ambiguity-fixed positioning precision enjoys a great increase due to the improved ambiguity resolution performance, which is more profound in the weaker Galileo-only and BeiDou-only solutions. In particular, the triple-frequency GPS-only and Galileo-only solutions show a TTFA of 23 and 51.5 min, respectively, while the BeiDou-only solution shows a considerable improvement compared to the dual-frequency case but still requires more than 120 min. This can be also seen in the bottom figures of panels (a) and (b), where the gain numbers of Galileo and BeiDou seem to increase earlier than in the dual-frequency case.
When at least two GNSS systems are combined, a tremendous improvement in the TTFA and positioning performance can be observed in both ambiguity-float and ambiguity-fixed solutions. Due to the increase of number of satellites in the solution, the ambiguity-float positioning precision shows a dramatic increase, with the best performance being reported from the combined GPS+Galileo +BeiDou solution. The worst performance, compared to the rest multi-GNSS solutions, is the one of the GPS +BDS solution due to the incomplete BeiDou-3 MEO constellation. In particular, the 10 cm precision level can be achieved within 20, 23 and 16.5 min (TTFA) for the GPS +Galileo, GPS+BeiDou and GPS+Galileo+BeiDou dual-frequency ambiguity-float solutions, while the addition of an extra frequency for all GNSSs does not improve the precision, as expected.
The PAR-based positioning solutions converge faster and are of higher precision for both dual and triplefrequency solutions. The TTFA for the dual-frequency GPS+Galileo, GPS+BeiDou and GPS+Galileo+BeiDou solutions are equal to 8.5, 11 and 6.5 min, respectively, with the triple-GNSS solution showing a horizontal precision gain above 10, while the two dual-GNSS solutions have a gain between 5 and 8. As expected, the incorporation of the third frequency per GNSS improves IAR and, therefore, the partially-ambiguity-fixed positioning. The TTFA for the triple-frequency GPS+Galileo, GPS+BeiDou and GPS+Galileo+BeiDou solutions are equal to 7.5, 9.5 and 4.5 min, respectively, with the improvement being on the order of about 13%, 16% and 44% compared to the dualfrequency counterparts. The TTFAs of the user's PARbased positioning results are summarized in Fig. 12.
At this point, it is important to emphasize that these results are sensitive to our measurement uncertainty assumptions, since the stochastic model affects both the ambiguity resolution and the positioning performance. In particular, the code and phase measurement STDs used were equal to 30 cm and 3 mm, respectively, with these choices being inappropriate in case a low-cost receiver is used. In addition, our stochatic model was also based on the simplifying assumption, as is commonly done in literature, that the PPP-RTK user-provided network corrections are sufficiently precise such that they can be treated as deterministic. This assumption might be valid when a network continuously generates such corrections that allows them to gain high precision over time, but this might not hold in sparse networks where data over short time spans are used and IAR is not successfully realized.
To investigate the impact of the stochastic model on the user's performance, we computed the user's positioning   11. PPP-RTK user horizontal positioning precision and gain (90th percentile) as function of time since the processing start, using PAR and the initially defined measurement vc-matrix Q yy . The processing window is re-initialized every 1 min within the selected day for all stations and system/ frequency combinations. G, E, C stand for GPS, Galileo and BeiDou-3, while the numbers within the parentheses denote the number of frequencies. The number of observed GPS, Galileo and BeiDou-3 satellites is 9, 7 and 5 on average. solutions based on varying measurement precision, with the PAR-based TTFAs being shown in Fig. 12. In general, one can observe that an increase of the observation STDs in the stochastic model deteriorates the user's performance in terms of TTFA, gain in position precision and achieved positioning precision over time, as expected. The reason for the increased TTFAs lies in the fact that it takes longer to achieve successful PAR with a positioning precision of better than 10 cm due to the increased measurement noise that directly affects the parameter estimation. In this way, one will experience less wrongly fixed solutions, as an overoptimistic stochastic model will lead to faster PAR realizations and, therefore, a larger number of wrongly fixed solutions. When the initial measurement vc-matrix is scaled by a factor of 4 or 8, implying the use of low-cost receiver data with a code observation STD of 60 or about 85 cm at zenith, respectively, the best-performing-single-GNSS GPS dualfrequency model achieves a TTFA of 56 or 93 min, respectively. Adding the third frequency signals reduces the TTFAs by 12% and 16%, respectively. As it can be seen, it is the integration of multi-GNSS signals that contributes the most to the reduction of TTFAs. The GPS+Galileo +BeiDou solutions achieve successful PAR with a positioning precision better than 10 cm within 17 and 27.5 min when two frequencies are used, respectively, with a 13--15% improvement by adding the third frequency signals. This sensitivity for varying observation STD seems to be less pronounced in the multi-GNSS multi-frequency models, especially for the triple-frequency GPS+Galileo +BeiDou model. Thus, as the user's model strength increases, by using multi-satellite and/or multi-frequency data, the TTFA sensitivity for varying code and phase precision becomes less.

Five-frequency Galileo PPP-RTK
The European satellite navigation system Galileo undergoes a rapid development in the last years, with the declaration of its initial services in December 2016 (GSA, 2018). The capability of Galileo to transmit signals on five Fig. 12. TTFA (90th percentile) for PAR-based kinematic PPP-RTK user solutions for varying measurement precision. The results are based on the sample solutions computed at the selected 9 globally distributed IGS stations with a 2-h processing window being re-initialized every 1 min for several system/frequency combinations. frequencies allows for faster ambiguity resolution compared to the rather weak dual-frequency model. It is, therefore, of great interest to evaluate the IAR and positioning performance of a multi-frequency Galileo user. To this end, we computed the PPP-RTK user's positioning solutions for the 9 globally distributed stations with every-1-min initializations of a 2-h processing window, having in total about 12000 sample solutions. We emphasize here that in our computations the Galileo observables of different frequencies were assumed to be uncorrelated, based on the fact that no significant correlation has been empirically found, even for the Galileo E5a/E5/E5b observables that are derived from the same wideband signal (Odijk et al., 2014b). The extracted 90th percentiles of the ambiguityfloat and ambiguity-fixed positioning precision over time, as well as of the gain in horizontal position precision, are shown in Fig. 13.
One can observe that the dual-frequency ambiguity-float position solution requires more than 2 h to reach the 10 cm level, while the increase in number of frequencies brings only a slight improvement in the performance, as it was observed in the previous sections. However, the increasing number of frequencies seems to bring a significant improvement in the PAR-based solutions. In particular, the dual-frequency ambiguity-fixed solution (90%) achieves a TTFA of more than 1 h, while the triple-and fourfrequency solutions bring the TTFA to 51 and 48 min, respectively. Although no significant improvement is observed when adding the fourth frequency, the fivefrequency user's model is shown to be greatly strengthened in terms of its PAR capabilities, having a TTFA of 25 min that is a reduction of almost 60 min compared to the dualfrequency counterpart. Therefore, it is concluded that a high number of frequencies is critical in achieving high ASRs, ambiguity-resolved positioning of high precision and large precision gains in a short time span.
In summary, and based on all of the above numerical results, it can be concluded that multiple GNSSs and multiple frequencies can bring a significant benefit to IAR, with the first being the main driving factor. Even though no a priori information was assumed for the slant ionospheric delays, the incorporation of multiple signals strengthens the model in such a way that IAR can be achieved in a short time span. This, of course, comes to a cost. The high-dimensional vc-matrix of the estimated float ambiguities seems to be a computational bottleneck. From our experiments, it was observed that the higher the dimension of the ambiguity vector, the higher the computational burden is in the decorrelation process and the construction of the Z-transformation matrix of LAMBDA. To speed up the decorrelation, one can make use of the Z-matrix from previous epochs to decorrelate Qââ in successive epochs assuming that no significant change exists in the satellite geometry, as proposed by Nardo et al. (2016).

Conclusions
In this contribution, we provided insight into and analyzed the expected ambiguity resolution and kinematic positioning performance of the ionosphere-float PPP-RTK user's model based on several multi-GNSS multifrequency models without any a priori ionospheric information. The user's performance was formally evaluated at several globally distributed stations and measured in terms of the formal ambiguity success rate, the required time to fix ambiguities, the achieved horizontal position precision and the expected gain after successful ambiguity resolution. Our first finding is that the increase in number of frequencies brings only a small improvement in the ambiguity-float position precision, due to the invariant satellite geometry, but improves the ambiguity resolution performance by providing a shorter TTFA for both FAR and PAR and, therefore, centimeter-level positioning can be achieved in a shorter time span. This was also demonstrated by analyzing the five-frequency Galileo user's performance, which shows a significant improvement of 70% compared to the dual-frequency model and a TTFA of 25 min. We also numerically demonstrated that PAR outperforms FAR for all solutions, in the sense that a Fig. 13. Galileo multi-frequency PPP-RTK user horizontal positioning precision and gain (90th percentile) as function of time since the processing start, with the processing window being re-initialized every 1 min within the selected day for all stations. The ambiguity-fixed results are obtained using PAR. sufficiently large subset of ambiguities can meet the ASR criterion and provide a position precision of better than 10 cm in a shorter time span and in a consistent manner. Although no considerable improvement in the position precision was found for single-constellation solutions using PAR, it was shown that significant gains can be expected when at least two GNSS systems are combined. Further, we showed that the increase in number of satellites used in the model improves both positioning and ambiguity resolution capabilities due to the strengthened geometry. Therefore, it is concluded that the satellite and frequency redundancy work in tandem to improve the user's performance, with the former being the main driving force for speeding up ambiguity-resolved positioning in the absence of ionospheric information. Despite the increase in the dimension of the ambiguity vector and the potential decrease of the FAR ASR, PAR was shown to be an efficient solution to the upcoming dimensional curse such that a large enough subset can be identified to both meet the ASR criterion and surpass the 10 cm precision level. Based on an extensive sample of solutions computed on a global scale, it was found that the dual-and triple-frequency GPS +Galileo+BeiDou PPP-RTK user solutions can reliably achieve a TTFA of 6.5 and 4.5 min (90th percentile), respectively, within which the user's position precision enjoys an improvement by more than one order of magnitude and gets better than 10 cm, based on the current status of constellations. Both triple-GNSS solutions provide significant gains and much shorter TTFAs compared to the single-GNSS solutions, which would be 23 min in case of the best-performing triple-frequency GPS solution. Finally, we provided numerical evidence on the sensitivity of the user's performance for varying code and phase precision, characterizing both high-grade geodetic and low-cost receivers, as well as considering the network corrections' precision. It was found that a longer time span is required to achieve PAR with significant gains due to the increased measurement noise, with the sensitivity being less pronounced for multi-GNSS multi-frequency models.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.