Assessment of ionospheric corrections for PPP-RTK using regional ionosphere modelling

This paper presents an analysis of the ionospheric corrections required to get a significant improvement in PPP-RTK performance. The main aim was to determine the improvement in the position precision and time-to-first-fix in the PPP-RTK user side using ionospheric corrections computed from a network. The study consists of two main steps. The first one includes an empirical investigation of the ionosphere model precision necessary to greatly improve the PPP-RTK performance in a simulated environment in terms of precision and convergence time. In the second one, an optimal ionosphere representation was developed to provide precise ionospheric corrections by parameterizing the ionospheric slant delays after the PPP-RTK network processing in terms of ionosphere model coefficients and differential code biases using real GNSS measurements. Experimental results demonstrate that the proposed methodology can be used for reliable regional ionosphere modeling and satellite code bias estimation, due to the consistency of the satellite code bias estimates with those provided from the International GNSS Service Analysis Centres, the high stability of the estimated receiver and satellite code biases and the low least-squares residuals of the network-based ionosphere modeling solution. Finally, it has been shown that the precision of ionospheric corrections at zenith needs to be better than 5 cm to enable faster PPP-RTK solutions.


Introduction
The integer ambiguity resolution (IAR) enabled precise point positioning (PPP) method, the so-called PPP-RTK (Wubbena et al 2005), is a state-of-the-art global navigation satellite systems (GNSS) technique that allows to determine high-acc uracy positions with short convergence time. The main idea behind PPP-RTK is to extend the PPP technique (Zumberge et al 1997) by providing single-receiver users, apart from precise orbits and clocks, with additional corrections (satellite phase fact that the carrier-phase ambiguities need time to converge. These ambiguities are not estimable as integers, because they are lumped with the receiver and satellite phase biases. In relative positioning techniques, such as with real-time-kinematic (RTK), these biases are eliminated with double-differenced measurements and, as such, the double-differenced ambiguities can be fixed to their integers.
Several methods have been formulated in the past to recover the integerness of the user ambiguities (Ge et al 2005, Collins 2008, Mervart et al 2008, Laurichesse et al 2009, Geng et al 2012, therefore enabling the PPP-IAR method realization. They usually employ between-satellite single differencing on the raw observations or use the ionosphere-free linear combination of the latter in order to eliminate receiver-related parameters or ionospheric slant delays, respectively. The differences among these implementations lie in the choice of parameterization, in the corrections applied and, on several cases, in the estimation method. Although at first glance one would say that different corrections are provided to the user due to the different S-basis choice, it has been shown that their information content is the same and can achieve the same goal, namely enabling the construction of a system of observation equations at the user component in mixed-integer form (Teunissen and Khodabandeh 2015).
As a result, these methods are able to resolve the phase ambiguities in the single-receiver observations and lead to RTK-like (mm-cm level) positioning accuracy. The elimination of the ionospheric error, however, is unfavorable since such corrections are required for the transition to PPP-RTK mode which can achieve significant shortening in the convergence times of the PPP-IAR positioning results (Odijk et al 2016a). An undifferenced and uncombined PPP-RTK model for mulation, as used in Odijk et al (2016b) and Henkel et al (2018), shows an obvious advantage over differences and linear combinations, as it contains all GNSS estimable parameters and, therefore has the benefit of providing (biased) ionospheric slant delays of high-precision that can be used for measuring the total electron content (TEC) of the Earth's ionosphere.
Over the years, there has been an extensive research on measuring the Earth's TEC using GNSS data (Schaer et al 1995, Manucci et al 1998, Gao and Liu 2002, Ciraolo et al 2007, Brunini and Azpilicueta 2009, Li et al 2015. The ionospheric observables, usually derived from the widely used geometry-free (GF) code or phase measurements and the carrier-to-code levelling (CCL) method, do not represent the unbiased slant TEC, due to the presence of the unknown carrier phase ambiguities or code hardware delays. Although one should not base one's precision analysis of TEC on that of the ionospheric observable (Khodabandeh and Teunissen 2017), the CCL method has been proven inaccurate to levelling errors (code noise and multipath effects) which sometimes exceed a couple of TECU (Banville 2014) (where 1 Total Electron Content Unit = 10 16 electrons m -2 and corresponds to 16 cm at the L1 frequency).
In the last few years, the GNSS ionosphere research community started turning its attention to alternative approaches for retrieval of TEC measurements (TECM). A network-based geodetic processing was employed from the UPC (Technical University of Catalonia) to retrieve the ambiguity term (Rovira-Garcia et al 2016), which is lumped in the GF phase measurements. Then, one is able to obtain the undifferenced ambiguity-fixed carrier-phase ionospheric observables, which are affected only by the code hardware delays.
Further, several recent studies have used the PPP approach with raw observations to retrieve ionospheric observables (Tu et al 2013, Zhang 2016, Liu et al 2017, Liu et al 2018. Although PPP relies on precise orbit/clock products and includes a more complicated data processing than that of the CCL method, the PPP-derived ionospheric slant delays are not affected by levelling errors (which might have undesired effects on ionosphere modeling) (Ciraolo et al 2007) and are more precise (Zhang 2016). Due to its capability to resolve the integer ambiguities, PPP-IAR is expected to provide much more precise slant ionospheric observables, as shown in the current work, since IAR is the key to fast and high-precision GNSS parameter estimation (Teunissen and Montenbruck 2017).
As already stated, the variety of ionospheric observables is vast, the interpretability of which is important to take into account in TEC determination (Khodabandeh and Teunissen 2017). Regardless of the combination in use, it is easily understandable that there is lack of information content in the undifferenced GNSS data to obtain unbiased ionospheric delays. Therefore, in order to retrieve the unbiased TEC and the lumped biases, the rank-deficiency of the GNSS observations needs to be identified and removed using the S-system theory (Baarda 1973, Teunissen 1985. A brief introduction to the singularity-system theory is given in this paper. Then, in order to separate the TEC from the hardware delays, a mathematical representation function is necessary to describe the ionosphere in the spatial and temporal domain, assuming it as a single-layer model. It has already been shown that the spherical harmonic (SH) functions are suitable for VTEC (Vertical TEC) modeling on a global scale (Schaer 1999). In the regional scale which we are interested for, loworder SH functions (Zhang et al 2015), adjusted SH functions (Liu et al 2018), bi-quadratic basis functions (Brunini and Azpilicueta 2009) or the combination of SH functions with generalized trigonometric series functions (Li et al 2015) are usually used.
If one is aiming at improving the PPP-IAR user performance, in terms of precision and most importantly convergence time, one has to study how good the ionospheric corrections need to be to enable the realization of PPP-RTK. Therefore, in this work, we present the method to assess the precision of the ionospheric corrections required to improve the PPP-IAR performance at the user side by means of the necessary time to fix ambiguities to their integers. Our approach consists of a design computation scheme, where a user is simulated to process GPS-only dual-frequency undifferenced and uncombined code and carrier phase measurements using the PPP-RTK technique. By presenting this methodology, we propose an ionosphere modeling strategy to improve the TEC precision at a regional scale and obtained accurate satellite code biases that are useful for positioning, navigation and timing (PNT) applications. This is the main contribution of this work.
The structure of this work is as follows. Section 2 reviews in brief S-system theory for a general model formulation. In section 3, the methodology of the design computation is analyzed to obtain the required precision of ionospheric corrections for PPP-RTK. We close this section by describing in detail the real data and ionosphere representation used for retrieving the unbiased TEC and satellite code biases at a regional scale. Section 4 presents the results of the design computation and our ionosphere modeling approach based on PPP-IAR derived TECMs. We conclude in section 5.

Brief review of S -system theory
Let us start with a linear model: where the observation and parameter vectors of dimensions k and l are denoted by y and x, respectively. Here E(·) and D(·) denote the expectation and dispersion operators, respectively. The design matrix A ∈ R k×l is rank-deficient with rank(A) = q l, while the measurement variance-covariance matrix (VCM) Q yy is assumed positive definite. A rank-deficient design matrix implies that not all the unknown parameters can be unbiasedly determined, given the information content in y due to linear dependence of some of the columns of A. This rank deficiency is of size dim where N (·) denotes the null space and R(·) denotes the range space; these two spaces are complementary.
Due to this rank deficiency, the parameter vector can be decomposed into its estimable x S and non-estimable part x V , using the l × q and l × (l − q) basis matrices S and V, respectively: where α denotes the q-vector containing the estimable parameter functions, while β denotes the (l − q)-vector containing the non-estimable parameter functions. Although there is not a unique S, the choice of S determines which estimable parameters are solved for and what their interpretation is. By inserting equation (2) into the rank-deficient linear observation model equation (1), one obtains the full-rank model: where Ã denotes the k × q full-rank design matrix of rank q. From equation (3), one can easily observe that the parameter changes leave the observations invariant. It should be highlighted here that the estimable parameters x S and x S , based on the basis matrices S and S respectively, cannot be directly compared due to the different choice of singularitybasis. If one wants to compare the two solutions, one out of the two needs to be transferred to the other's S-basis using an S-transformation matrix (Baarda 1973, Teunissen 1985.

Methodology
This section starts with the functional model of the PPP-RTK user and ends with the models and algorithms used in the regional VTEC modeling.

PPP-RTK user design computation
Faster convergence times are expected if the PPP-RTK user corrects a priori for the ionospheric delays, which are computed and modeled at the network side. The reason is that their presence significantly affects the resolution of integer ambiguities and, therefore, the solution's convergence time.
Prior to the ionosphere modeling step, the first thing that one should do is to investigate how precise the ionospheric corrections need to be in order to enable faster integer ambiguity resolution, and therefore enable a significant reduction in the convergence time of the PPP-RTK solutions.
For this reason, a design computation is performed in this study, simulating a GPS-only dual-frequency PPP-RTK user environment in order to assess the effect of the ionospheric corrections precision on the time-to-first-fix (TTFF), instead of the convergence time since no real data are used in this case. TTFF here refers to the required amount of time needed to achieve successful integer ambiguity resolution based on a pre-defined success rate. The success rate is an important measure, since it indicates the probability that the ambiguities have been fixed to the correct integers. Once the ambiguities are resolved and TTFF is obtained, the estimable GNSS parameters will converge faster due to the stronger functional model.
The ionospheric corrections can be estimated and modeled within a PPP-RTK network component. The basis of the PPP-RTK network system in this study consists of the set of un-differenced and uncombined carrier phase and pseudorange observation equations. For a receiver-satellite combination r − s at frequency j, they are defined as (Teunissen and Kleusberg 1998): (4) where φ s r,j and p s r,j denote the phase and code measurements, ρ s r the receiver-satellite range, τ r the tropospheric zenith delay, m s r the tropospheric mapping function, dt r and dt s the receiver and satellite clock offsets, ι s r the (first-order) slant ionospheric delays on the first frequency, µ j the frequencydependent ionospheric coefficient, d r,j and d s ,j the receiver and satellite code biases, δ r,j and δ s ,j the receiver and satellite phase biases, a s r,j the integer phase ambiguity, and λ j the wavelength at frequency j.
The above variables have a receiver index r = 1, . . . , n, with n the number of receivers, a frequency index j = 1, . . . , f , with f the number of frequencies ( f = 2 in this paper), and a satellite index s = 1, . . . , m, with m the number of satellites. All variables are time-dependent and are expressed in meters, except for the phase biases and ambiguities, which are expressed in cycles while the latter remain constant over time unless a cycle slip occurs. The ionospheric coefficient is defined as the squared ratio of frequencies: µ j = ( f 1 /f j ) 2 . In a dual-frequency GPS-only case, µ 1 = 1 and µ 2 = (77/60) 2 .
The above un-differenced observation equations cannot be used directly to estimate all the unknown parameters, since the design matrix is rank-deficient. To solve for the rank-deficient system of observation equations in the PPP-RTK network side, the S-system theory (Teunissen 1985) is applied, according to which several parameters are mapped to others in order to allow for a full-rank system of observation equations. In this paper, a common clocks (pivot) receiver (CC-R) S-basis (Odijk et al 2016a) is used to overcome the rank deficiencies, and the estimable parameters (denoted using the ( .) symbol) at a single receiver r are presented in table 1.
Assuming the receiver and satellite positions are known and precise enough, the observed-minus-computed observation equations become as follows: The network-derived satellite clock offsets and satellite phase biases comprise the key for the single-receiver PPP-IAR users to enable integer ambiguity resolution. Although these estimable parameters are biased, they can still do the job for the PPP-IAR user if the latter employs the same functional model with the same parameter mapping that was used in the network component. In that case, the interpretation of the estimable user parameters is the same as this in the network component.
Linearizing the observation equations with respect to the unknown user position and applying the precise satellite orbits and network-derived corrections for the satellite clock offsets and phase biases, the user's dual-frequency code and carrier phase measurements (with user index u) are as following: where φs u,j and p s u,j denote the observed-minus-computed phase and code measurements; g s T u denotes the unit vectors pointing from the satellites to the receiver.
It can be deduced that the user's receiver phase biases and integer carrier phase ambiguities are separable now, leading to a full-rank system of observation equations, due to the fact that the integer ambiguities vanish for the pivot satellite while the receiver phase biases do not, unlike equation (4) where the phase biases and ambiguities are lumped into one frequencydependent ambiguity term. This separation is a direct consequence of the S-basis used in this study to overcome the rank deficiencies between the phase biases and ambiguities, as can be seen from table 1. As a result, the user is able to perform IAR when the estimated ambiguities are precise enough and meet a pre-defined success rate threshold.
3.1.1. Ionosphere-float model. The PPP-IAR user model consisting of the observation equations equations (7) and (8), in which the (biased) ionospheric slant delays are estimated as unknown parameters, is the so-called ionosphere-float model (Teunissen 1997).
The undifferenced and uncombined code and carrier phase measurements are described by the following stochastic model: Estimable dual-frequency network parameters, including their interpretation and conditions using the CC-R S-basis (the symbol p denotes the pivot satellite/receiver if it is used as superscript/subscript).

Estimable parameter Interpretation Conditions
where Qφφ and Qpp denote the dual-frequency measurement VCMs for the observed-minus-computed phase and code measurements, respectively. It is usually assumed that no correlation exists between frequencies and between code and phase measurements. In reality, however, there exists correlation (a) between code and phase measurements since both of them have been corrected for the satellite orbits and the network-derived satellite clock offsets, and (b) between frequencies for the phase measurements because of the applied satellite phase biases. Only the first type of correlation was taken into account in our study. The satellite phase biases and satellite clock offsets transmitted to the PPP user are the key for IAR-enabled precise point positioning. To achieve that, the integer ambiguities need to be fixed correctly. However, it is known that the ionosphere-float PPP-IAR model is rather weak in terms of integer ambiguity resolution, since the estimable parameters for the unknown ionospheric delays affect the solution's convergence time. Therefore, a great shortening in the convergence time is expected in case ionospheric corrections are available to PPP-IAR users (Odijk et al 2016a).
3.1.2. Ionosphere-fixed model. In case ionospheric corrections are provided to PPP-IAR users, by either spatial interpolation or function-based modeling, faster integer ambiguity resolution than by using the ionosphere-float model is expected, since unknown parameters for the ionosphere do not need to be estimated. This is the ionosphere-fixed model, in which such precise ionospheric corrections are provided to the PPP-IAR users that can be assumed to be deterministic. As a result, a combined parameter of the GF receiver and satellite code biases becomes estimable: where d u,GF and d s ,GF are scaled versions of the satellite and receiver differential code biases (DCB). It is, therefore, intentional to estimate the satellite DCBs (SDCBs) at the network side, in order to provide them to and allow the user to solve for less parameters making the used observational model stronger. The receiver DCBs (RDCBs) and SDCBs can be separated by selecting a proper S-basis. Therefore, the provision of ionospheric corrections and SDCBs to the user allows the estimability of a scaled version of RDCB.
3.1.3. Ionosphere-weighted model. The aforementioned ionosphere-fixed model changes to an ionosphere-weighted model, firstly introduced by Bock et al (1986), in case the provided ionospheric corrections are assumed to be stochastic parameters, rather than deterministic. The functional model remains the same as in the ionosphere-fixed model, but with a different VCM. In particular, the uncertainty σ ι of the provided ionospheric corrections is taken into account as following: where ⊗ denotes the matrix Kronecker product (Henderson et al 1983), and µ = (µ 1 , µ 2 ) T is the 2-vector containing the wavelength ratios.
Since the same ionospheric corrections are applied to both code and phase measurements, the VCM is no longer a blockdiagonal matrix. Instead of applying the stochastic ionospheric corrections directly to the phase and code measurements, the user can use them as pseudo-observations and weight them based on their standard deviations.
Therefore, the ionosphere-weighted model is the general model, from which the aforementioned other two models can be produced. If unknown parameters for the ionospheric delays are estimated or prior ionospheric corrections do not contribute to the solution, the model is transformed into the ionosphere-float one (σ ι = ∞). If, on the other hand, the ionospheric corrections are precise enough to be assumed deterministic, the ionosphere-float model is transformed into the ionosphere-fixed one (σ ι = 0).

Ionosphere modeling
GNSS-based measurements have proven to be capable of remotely sensing the Earth's dynamic ionosphere. As already stated, there are several methods to extract slant ionospheric slant delays with varying interpretation and precision. Out of all of them, we selected the PPP-IAR technique in order to process undifferenced and uncombined code and phase measurements and obtain the estimable parameters as shown in table 1.
The estimable ionospheric slant delays (expressed in meters) are biased by the receiver and satellite DCBs, using the chosen S-basis: In order to estimate both the VTEC and the satellite and receiver DCBs simultaneously, one needs to make use of the thin-layer ionosphere model. According to the latter, the ionosphere is assumed to be a spherical shell at a height of 450 km above the Earth's surface. The slant total electron content ι s r (STEC) is mapped to its vertical counterpart v s r at the points where the satellite-to-receiver signal paths intersect the ionospheric shell, the so-called ionospheric pierce points (IPPs), using the following mapping function M s r (i) (Brunini and Azpilicueta 2009): where R is the mean Earth's radius, h is the height of the ionospheric shell (450 km in our case), and Z s r (i) is the zenith angle of satellite s observed from receiver r at epoch i.
Then, considering that ι s r (i) = M s r (i) · v s r (i), the VTEC at an IPP can be mathematically modeled by a wide variety of representation functions in the time and space domain. In this study, the generalized trigonometric series function (Li et al 2015) (sum of a polynomial function and a finite Fourier series) was used to model VTEC on a regional scale: where φ IPP and φ REC denote the geomagnetic latitude of the IPPs and the receivers, respectively; Λ IPP denote the solar longitude of the IPPs; A, B and K are the maximum orders of expansion; E ab , C k and S k are the model coefficients to be estimated as functions of time.
Thus, if the model coefficients are stored in a vector w and their corresponding scaling factors at epoch i are stored in the design matrix A s r , the VTEC in units of meters is described as v s In addition to the VTEC, which is expressed as a function of ionosphere model coefficients, the biased ionospheric observables contain the RDCBs and the SDCBs, which also need to be estimated. The system of observation equations shows a rank-deficiency, though, since the receiver and satellite code biases cannot be unbiasedly estimated. For this reason, a proper S-basis was employed according to the S-system theory in order to form a full-rank system of observation equations.
From the observation equation equation (12) at a single receiver r , it can be deduced that the ionospheric observables determined by a single receiver r are not enough to determine both the ionosphere model coefficients and the receiver and satellite hardware delays. As such, the ionospheric observables from n receivers of a regional network shall be used. This also strengthens the observation model, since the satellites are tracked by multiple instead of only one receivers, which allows for reliable estimation of satellite DCBs. Therefore, given the mapping function equation (13), the functional model reads where A(i) is the partial design matrix containing the scaling factors of the ionosphere model coefficients w(i) for all receivers r = 1, . . . , n, and M(i) is a diagonal matrix of order m · n having the corresponding mapping functions at epoch i as its entries. The estimable receiver and satellite DCBs are denoted as d r and ds , respectively. The n-vector having 1's as its entries is denoted as e n , and the unit matrix of order n is denoted as I n .

Results and analysis
We begin this section by describing the experimental setup both for the PPP-RTK user design computation and the ionosphere modeling using real GNSS data. These are then followed by numerical results, from which the major findings are described in detail.

Experimental setup
In the first part of our study, we simulated a dual-frequency GPS-only PPP-RTK user component in various receiver locations around the world, and employed the common clocks (pivot) receiver (CC-R) S-basis, shown in section 3.1, to overcome the rank deficiencies. In this part, a formal analysis was performed and, therefore, no real data was used. The undifferenced code and carrier phase observations, sampled every 30 s during DOY 046/2014, were assigned with zenith-referenced a priori standard deviations (σ 0 p and σ 0 φ ) of 30 cm and 3 mm, respectively. An elevation-dependent scheme (el denotes the elevation angle) was used with the variances of the code and phase observations being calculated as (Dach et al 2015) An elevation mask of 10 • was used in this study to avoid measurements acquired from satellites close to the horizon, while the precision of the precise orbit information was assumed to be equal to 2.5 cm, instead of considering it as a deterministic quantity. The ionosphere-float and ionospherefixed GNSS models were initially used to get the extreme cases for the obtained TTFFs, since (a) in the first case, unknown parameters for the ionospheric delays are estimated by the user making the model weak in terms of IAR, while (b) in the second case it is assumed that deterministic ionospheric corrections are provided to the user enabling fast IAR. The ionosphere-weighted model was then employed using a varying precision for the ionospheric corrections in order to find the optimal stochastic ionospheric corrections that can enable a shortening in the TTFF and, consequently, the convergence time.
The GNSS parameter estimation is performed in a Kalman filter. In this regard, the process noise of the parameters linked in time are listed in table 2. The parameters not listed are estimated as unlinked parameters in time. The estimable parameters share the same process noise values as the unbiased ones. The first step of the mixed-integer GNSS model solution results into the so-called float solution, if one ignores the integer property of the carrier phase ambiguities: where a is the 2(m − 1) ambiguity vector, and b is the vector containing the rest of the estimated parameters (estimable receiver clock offset, tropospheric zenith delay, estimable ionospheric slant delays, estimable receiver code and phase bias). The parameters of interest for the TTFF evaluation are the estimated float double-difference (DD) ambiguities, since their successful fixing depends on their precision, contained in Qââ. This is the input of the second step of the mixed-integer GNSS model solution, which focuses on the integer constraint a ∈ Z 2(m−1) , i.e. the mapping of the float ambiguities â into their corresponding integer ones ǎ with an integer mapping I : R 2(m−1) → Z 2(m−1) such as ǎ = I(â).
In order to evaluate whether the integer ambiguities can be estimated successfully, the integer least squares (Teunissen 1993) lower bound was used, i.e. the bootsrapped success rate P s,B . Since the bootstrapped estimator performs almost optimally after decorrelating the ambiguities using the Z -transformation of the LAMBDA method (Teunissen 1995), the success rate P s,B is evaluated for the decorrelated ambiguities ẑ = Zâ: Based on the decorrelated ambiguity VCM, the success rate can be evaluated (Verhagen 2005): with Φ denoting the cumulative normal distribution, and σ i|I the standard deviation of the ith least-squares ambiguity obtained through a conditioning on the previous I = 1, . . . , i − 1 ambiguities. Given a user-defined minimum threshold for the ambiguity success rate, which we set to 99.5% in our study, successful integer ambiguity resolution occurs when the estimated bootstrapped success rate is larger than this threshold.
Apart from the simulations, we also used real GNSS data to validate the performance of the PPP-RTK technique with raw code and phase measurements for regional ionospheric VTEC modeling. For this reason, CORS geodetic-grade receivers in North Carolina (US) of the NGS (National Geodetic Survey) network were selected in order to form a regional network (see figure 1), for the scope of this study. The dual-frequency GPS dataset (L1C, L2C, C1C, C2W) was sampled every 30 s by 45 geodetic-grade receivers of the same receiver type (TRIMBLE NETR5) over the DOY 046, 2014. The data were processed in our PPP-RTK engine, according to the parameterization given in table 1. The code and phase measurements were weighted according to their elevation, with an undifferenced zenith-referenced standard deviation of 30 cm and 3 mm, respectively. We also used a cut-off elevation angle of 5° to discard noisy measurements. A Kalman filter is used for the GNSS parameter estimation, using the precise orbit and clock products distributed by the European Space Operation Center (ESOC) of the European Space Agency (ESA). Unknown parameters for the ionospheric delays are estimated for every receiver-satellite link.
As stated earlier, the PPP-IAR-derived biased ionospheric slant delays serve as input for modeling the VTEC in the selected regional network. In order to allow for a rigorous and reliable determination of the ionospheric model coefficients and the RDCBs and SDCBs, a pre-processing of the TECM observables was applied. In particular, the data of the first and last 50 epochs (25 min) of each observable arc were excluded  in order to avoid estimates computed during the convergence period and the satellite setting interval. The unknown parameters are the ionosphere model coefficients, receiver and satellite DCBs. At this stage, a Kalman filter was used to determine these parameters epoch by epoch, using a cut-off elevation angle of 12 • to discard noisy measurements.

Results of the design computation
In this section, the impact of the zenith-referenced ionospheric corrections precision on ambiguities, and as such, on the achieved position precision and TTFF is investigated at the PPP-RTK user component. The TTFF is defined as the number of epochs required to obtain reliable integer ambiguity fixing, based on the pre-defined probability of correct integer ambiguity fixing which we set to 99.5%.
The achieved formal precision of the horizontal and vertical position components is illustrated in figures 2 and 3, using GPS dual-frequency measurements and making use of the ionosphere-float, -fixed and -weighted models. It can be easily seen that the TTFF is almost 30 min in case of the ionosphere-float model, while the ionosphere-fixed model achieves instantaneous IAR with the formal precision of the horizontal and vertical components reaching the 2 and 5 cm, respectively. As expected, the achieved formal precision in the vertical position component is worse than in the horizontal component for all examined cases.
Moreover, one can observe the effect of the ionospheric corrections precision on the TTFF of the solutions. In particular, in case the precision of the corrections is 16 cm (almost 1 TECU), then the improvement in TTFF is negligible. It can be seen that although an ionospheric precision ranging from 8 to 16 cm shortens the TTFF of the solutions, this improvement is in general small. Ionospheric corrections with precision of 4 cm lead to a TTFF equal to 17 min for the PPP-RTK model, almost half of the time required for the ionospherefloat PPP-IAR model to achieve successful integer ambiguity resolution.
After an extensive data analysis with various ionospheric error precisions and for 1440 initialization times during the day, i.e. for every minute of the day, in order to take into account the effect of satellite geometry, we present the final results in figure 4. One can deduce that a significant improvement in the TTFF for the PPP-RTK user is observed in case the ionospheric corrections have a precision better than 5 cm.

Results of the ionosphere VTEC modeling
First, we focus on figure 5, illustrating the ambiguity-fixed estimable ionospheric slant delays from a single CORS receiver in North Carolina for all the observed satellites throughout a day. Ambiguity-fixed estimate is a parameter that has been estimated after successful integer ambiguity resolution. The magnitude of the ionospheric delays is relative since they are biased by the receiver and satellite DCBs. However, one can easily observe the typical signature of ionosphere due to the higher variation and larger magnitudes of the ionospheric observables at daytime than at night (North Carolina: UTC = LT-4 h).
In figure 6, the formal precision of the ambiguity-float and ambiguity-fixed ionospheric delay estimates are depicted. The convergence process of the estimates at each arc beginning is obvious, in which formal precisions are within 10 TECU. Formal precisions converge to the 0.20 TECU level in 30 min at minimum and in 2 h at maximum for the ambiguity-float case, depending on the observational session duration, while in 1 min at minimum and in 10 min at maximum for the ambiguity-fixed case. After the convergence time, the ambiguityfloat estimates have a formal precision still larger than 0.10 TECU, whereas the ambiguity-fixed ones can reach the 0.06 TECU precision level, clearly showing the dramatic improvement in precision after IAR. Moreover, it can be seen that formal precision becomes worse at the end of the arcs, which is due to the satellite setting towards the horizon.
Apart from the formal precisions, the precision of the ionospheric slant delay estimates can be validated with the between-receiver (BR) differences of short or zero baselines (Ciraolo et al 2007, Liu et al 2017. The BR differences of STECs eliminate the lumped satellite DCBs and most of the ionospheric errors, with the BR RDCBs and the ionospheric delay residuals being the remaining parameters.
The BR STEC differences of two CORS stations in North Carolina, with inter-station distances of 45 km and 500 km, are shown in figures 7 and 8, respectively, on DOY 064, 2014. All four receivers share the same receiver type, meaning that the BR RDCBs are eliminated. In the 45 km baseline, it can be seen that the BR differences show a fast convergence, with most of the differences not exceeding the 0.01 TECU level. Although the same level is achieved in the 500 km baseline, one can easily observe the much longer convergence time, probably due to small-scale changes in the ionosphere within this distance. Through this analysis, we conclude that PPP-IAR with raw observations can extract high-precision TEC measurements and avoid the levelling errors that are present in the CCL method.
The performance of the proposed regional ionosphere VTEC modeling algorithm was first evaluated based on the least-squares residuals of the ambiguity-fixed ionospheric observables over the selected day, as shown in figure 9. It can be observed that most of the measurement residuals do not exceed the 1.00 TECU level, while 90% of them range within 0.50 TECU. However, a few measurement residuals exceeding the 2.00 TECU level can be observed. These residuals correspond to the measurements acquired by the newly tracked and lost satellites observed from the CORS network receivers and, therefore, their corresponding estimates need some time to converge. The root mean square (RMS) of the residuals for all the receiver-satellite links throughout the selected day is equal to 0.48 TECU, which indicates that the selected representation function can fit well the ionosphere on the selected day. Then, an assessment of the modeled ionospheric corrections followed. For this reason, the self-consistency test (Orus et al 2005) was used, which analyzes the slant ionospheric delay variations along a continuous arc (satellite pass) over each station. The epoch in which the satellite is at its highest elevation was assigned as the reference epoch (Hernandez-Pajares et al 2017). This is an internal consistency test, providing a quality measure for the STEC computed by the used ionosphere model. The self-consistency metric is defined by the daily root mean square of the STEC variation ∆ι: where ι 0 is the PPP-IAR-derived ionospheric observable, ι m is the ionospheric slant delay derived from the estimated model coefficients, el max denotes the highest satellite elevation, i the epoch. The receiver and satellite code biases are assumed to be constant over time and are, therefore, cancelled in the differencing over a continuous arc. Figure 10 illustrates the self-consistency RMS measure for all receiver-satellite pairs, where the receiver investigated for the self-consistency test was excluded from the modeling step to avoid over-optimistic results. One can easily observe that most of the RMS values do not exceed the 1.50 TECU level, while the overall RMS equals to 1.10 TECU. However, it seems that there exist a few outliers, since RMS values greater than 2.00 TECU are observed for a few receiver-satellite pairs. The estimated ionosphere model was externally validated using the IGS global ionosphere map (GIM) over the selected region for the selected day of year. In particular, the global Center for Orbit Determination in Europe (CODE) GIMderived gridded VTEC values were compared to the modeled VTEC in the regional area, resulting to a mean offset of 0.91 TECU and an RMS equal to 4.50 TECU. Therefore, one concludes to the fact that there is a bias between the modeled and IGS-derived TEC models, which is due to a variety of factors. First of all, the CODE-derived global TEC map was assumed to be our ground-truth, although it is known that its accuracy ranges from 2 to 8 TECU. Moreover, a further VTEC interpolation from the grid points to points of interest result in general to larger errors. In addition, the ionosphere representation function used in our study may not be able to model the medium-scale variations of the ionosphere in our regional network, causing the detected bias.
Another performance indicator of our ionosphere modeling methodology is the behavior of the satellite and receiver DCBs. The estimates for the ambiguity-fixed GPS SDCBs on DOY 046 are shown in figure 11. A stable behavior is easily observed for the code biases of almost all GPS satellites after the convergence process. The convergence time ranges from a few to several hours and, therefore, the SDCBs that are observed only for a short amount of time do not converge to a constant value. Actually, this is a disadvantage of the regional ionosphere modeling, whereas within a global network (for global VTEC modeling) the satellites are observed without gaps due to the global distribution of the stations. Moreover, it seems that the DCBs of several GPS satellites do not remain stable over the time interval during which the satellites are observed. This is attributed to the fact that those satellites are not tracked sufficiently well from all the CORS receivers of the regional network we used. Figure 12 illustrates the formal precision of the GPS SDCBs. It is observed that most of the satellite code biases can reach a precision between 3 and 20 ps, or 1 and 6 mm respectively, depending on the observational session duration. It is deduced, therefore, that the longer a satellite is observed from the network, the more precise its DCB estimate becomes.
In order to further validate the performance of our proposed methodology for regional VTEC modeling and satellite DCB estimation, our estimable SDCBs were compared to those provided by the IGS. Within the Multi-GNSS Experiment (MGEX) , the German Aerospace Center (DLR) provides satellite DCB products for multi-GNSS signals, including the GPS C1C-C2W which we are interested for. Due to the fact that our estimable SDCBs are estimated based on a different S-basis than that of the IGS (zero-mean condition of satellites), their direct comparison is not possible. In order to allow for their comparison, they have to be transferred to a common S-basis. Given that, our estimable SDCBs were transferred to the S-basis of the IGS using an S-transformation (Baarda 1973). Figure 13 shows the errors of our estimable SDCBs (averaged after convergence), with respect to the DLR-derived SDCBs, based on a single-day dataset. In contrast to figure 12 which shows the standard deviation of the satellite DCB estimates (derived from their variance-covariance matrix), the errors shown in figure 13 serve as a measure of accuracy for the satellite DCB estimates. The satellite DCB errors show a zero mean, which was expected since both SDCBs are referred to the zero-mean condition. In addition, it can be deduced that the estimable DCBs of most of the GPS satellites do not deviate more than 2.0 ns from the published products, while 70% of them show errors lower than 1.5 ns. Overall, the RMS of the errors is equal to 1.3 ns.
In addition to the SDCBs, the temporal behavior of the RDCBs is another performance indicator of the proposed algorithm. Figure 14 illustrates the estimates for the ambiguity-fixed GPS RDCBs of the CORS network receivers on DOY 046. A stable behavior can be observed for the code biases of all the used CORS receivers, which is more visible when their average is removed ( figure 15). Their standard deviation ranges from 0.02 to 0.17 ns within a day, and this observed receiver DCB stability is another indicator that DCB  estimation is feasible in regional networks and can achieve high-precision results.

Conclusions
The main idea behind PPP-RTK is to extend the PPP technique by providing single-receiver users, apart from precise orbits and clocks, with external information (satellite phase biases, ionospheric and tropospheric corrections) so as to enable fast integer ambiguity resolution and, therefore, short convergence time. Although the undifferenced and uncombined GNSS observation model shows great flexibility for a potential multi-GNSS integration for strengthening the model with dynamic constraints on all the parameters, the unknown spatially correlated ionospheric errors still affect the GNSS observables and the convergence time, since the ionospherefloat PPP-IAR model is rather weak in terms of integer ambiguity resolution.
Faster convergence times are expected if ionospheric corrections are provided to the PPP-IAR users. Although the GNSS community has conducted a thorough research on measuring Earth's TEC, they extensively used the geometry-free code or phase measurements and the carrier-to-code leveling method, which are prone to levelling errors. As a result, the alternative approach of undifferenced and uncombined PPP was recently introduced as a means of extracting more accurate TEC measurements, although still biased by hardware delays. Due to its capability to resolve the integer ambiguities, PPP-IAR is the key to obtain high-precision TEC observables which are still biased by hardware delays, but unaffected by code noise and multipath. In this study, an analysis of the ionospheric corrections required to get a significant improvement in PPP-RTK performance was investigated. The main aim was to determine the improvement in the positioning precision and TTFF in the PPP-RTK user side using ionospheric corrections modeled from a network. The performed design computations clearly showed that faster PPP-RTK solutions are expected in case ionospheric corrections of 5 cm (~0.31 TECU) precision are available to the users, since the carrier-phase ambiguities can be fixed to their integer values faster.
Then, we proposed a methodology to model the PPP-IAR derived (biased) ionospheric delays on a regional scale, by parameterizing the ionospheric slant delays in terms of ionosphere model coefficients and DCBs using real GNSS measurements. PPP-IAR processing can provide high-precision ionospheric slant delays to be used for measuring the Earth's TEC. It was deduced that the proposed methodology can be used for reliable regional ionosphere modeling (RMS equal to 1.10 and 4.50 TECU for internal and external validation, respectively) and estimation of satellite DCBs (RMS of errors equal to 1.30 ns with respect to DLR products), which are required for the ionosphere-weighted PPP-RTK model. Although the zenith-referenced precision of our modeled VTEC reached the 5 cm level, a further investigation is needed to evaluate our modeled ionospheric corrections at the PPP-RTK user side in terms of convergence time reduction. Our method can be used for both real-time and post-processing, since in our study the measurements were processed epoch by epoch with a Kalman filter. The accuracy of the proposed methodology is expected to improve when a two-layer model is used for better modeling the ionospheric structure, and alternative ionosphere regional representation functions are employed.