Ionospheric tomographic common clock model of undifferenced uncombined GNSS measurements

Abstract

In this manuscript, we introduce the Ionospheric Tomographic Common Clock (ITCC) model of undifferenced uncombined GNSS measurements. It is intended for improving the Wide Area precise positioning in a consistent and simple way in the multi-GNSS context, and without the need of external precise real-time products. This is the case, in particular, of the satellite clocks, which are estimated at the Wide Area GNSS network Central Processing Facility (CPF) referred to the reference receiver one; and the precise realtime ionospheric corrections, simultaneously computed under a voxel-based tomographic model with satellite clocks and other geodetic unknowns, from the uncombined and undifferenced pseudoranges and carrier phase measurements at the CPF from the Wide Area GNSS network area. The model, without fixing the carrier phase ambiguities for the time being (just constraining them by the simultaneous solution of both ionospheric and geometric components of the uncombined GNSS model), has been successfully applied and assessed against previous precise positioning techniques. This has been done by emulating real-time conditions for Wide Area GPS users during 2018 in Poland.

1 Introduction

In this paper we introduce the Ionospheric Tomographic Common Clock (ITCC) model of undifferenced uncombined GNSS measurements. It allows, in particular, the precise positioning service in Wide Area networks based on undifferenced and uncombined modelling, which is specially suitable for multi-GNSS and multifrequency signals, due to its simplicity and bounded computational burden. Furthermore, it does not need to be continuously provided with external real-time products (like precise satellite clocks), e.g. in Wide Area networks. The model has been implemented on an updated version of the tomographic model TOMION for precise ionospheric sounding, tropospheric sounding and GNSS Navigation, Hernández-Pajares et al (2000, 2001); Colombo et al (2002); Hernández-Pajares et al (2003); Graffigna et al (2019b, 2019a).

TOMION v1 was initially designed for solving ionospheric tomographic models based on dual-frequency GPS measurements and for generating Global Ionospheric Maps, Hernandez-Pajares et al (1997); Hernández-Pajares et al (1999). Afterwards, this ionospheric model was further developed and synergistically connected with a precise GNSS non-frequency-dependent terms model, which was able to estimate position, clock errors and zenith tropospheric delay: TOMION v2. It implemented not only Precise Point Positioning (PPP) but also innovative techniques like the Wide Area Real-Time kinematics (WARTK), at both Central Processing Facility (CPF) and roving user levels Hernández-Pajares et al (2000, 2003). Recently, TOMION v2 has been updated by incorporating: (1) the most up-to-date tropospheric modelling, which provides precise estimation of gradients that reflect the dynamics of the troposphere even during hurricanes Graffigna et al (2019b); and (2) the further realistic models in case of very stable receiver clocks and the earth’s crust elastic response of hydrological load variations Graffigna et al (2019a).

TOMION v1 and v2 computations are based on combination of observables with different frequencies. The combination of observables removes common errors and unknown parameters. However, proper dynamic models of the parameters may provide methods that do not require combination, nor differencing of observables, thus yielding observables and thereby estimables with lower uncertainty. As an alternative to the combined-and-differenced approach, new GNSS-based positioning and navigation models based on undifferenced and uncombined (UU) GNSS data have been developed by the geodetic community, e.g. Zhang et al (2011); Teunissen and Khodabandeh (2015); Odijk et al (2016); Zhang et al (2019). The main advantage of UU GNSS models (hereinafter UUM) is that their implementation is simpler with respect to combined and differenced (CD) methods. For example, the main structure of the UUM equations does not depend on the GNSS constellation, thus easing their implementation in multi-GNSS scenarios. Moreover, by not combining GNSS data in time and/or in frequency, the UU estimables are less noisy when stochastic methods are used to solve the GNSS equations. Furthermore, combined methods, such as standard PPP, remove almost all the ionospheric delay, which turns out to be a bottleneck for the convergence time of the ambiguities, up to the point that it may take an hour to achieve cm-level accuracy. However, the provision of external ionospheric corrections increases the convergence speed by improving the ambiguity success-rate, and, consequently, it reduces the time for the user to achieve a given threshold in positioning accuracy, Olivares-Pulido et al (2019).

There are as many UUMs as linear combinations of estimables. For example, one of the first introduced UUM is the Common Clocks (CC) model, which provides estimated receiver and satellite clocks parameters that are common for all phase and code observables in the model, De Jonge (1998); Odijk et al (2016, 2017). High-precision positioning models, such as, for example, PPP-RTK (which is integer ambiguity resolution enabled PPP) have been developed based on the UU approach, e.g. see Zhang et al (2011); Teunissen and Khodabandeh (2015); Zhang et al (2019).

To the extent of our knowledge, all UUMs present in the current literature estimate the ionospheric delay as either proxies of the slant total electron content (STEC) estimates, or as vertical total electron content VTEC estimates, Zhang et al (2012); Odijk et al (2016, 2017).

It is well known that the STEC introduces a rank defect in the system of positioning equations when it is chosen as the ionospheric estimable, Odijk et al (2016). We will show later on that such rank defect can be fixed by computing the STEC from the VTEC. Models that compute the STEC by mapping the VTEC are based on the hypothesis that most of the ionospheric electron density is on a thin-layer surface at around $450\,km$ height, Klobuchar (1987). However, the thin-layer modelling is detrimental for the ionospheric modelling, especially for large geodetic networks and in low-latitude regions.

Another aspect to consider when either STEC or VTEC estimates are computed is that the dimension of the matrices involved in the computation, e.g. covariance, Jacobian (a.k.a. design matrix, model matrix), scale geometrically with the number of receivers and transmitters. Therefore, the computational burden of the model increases geometrically as well.

An alternative to both aforementioned models is the tomographic approach. It is well known that ionospheric tomography increases the accuracy and resolution of the ionospheric modelling with respect to the thin-layer model, Hernández-Pajares et al (1999). Additionally, later on, we will prove that the tomographic model also fixes the rank, and it makes the design matrix scale linearly with the data set, thus lightening the computational loading.

Note that, in the tomographic approach, the ionospheric delay is computed as the numerical integration of the ionospheric electron density at some points (pixels in 2D, voxels in 3D) along the line of sight from the transmitter to the receiver. Therefore, although the integrals, i.e. the ionospheric delays, are different for each satellite-transmitter pair, the integrand, the ionospheric electron density model, is common to all users and receivers. In other words, all users and receivers see the same ionospheric model. This is the tomographic method to estimate the ionospheric delay, which consists in estimating the ionospheric electron density at each voxel.

The tomographic model provides a dynamic model for the ionospheric electron density, which allows for the estimation of the ionospheric delay in the state-space, instead of the observational space. The tomographic method is a further extension of geodetic methods based on space-state variables, such as, for example, PPP. Furthermore, the STEC is, as opposed to the ionospheric electron density, a parameter relative to the receiver and the transmitter, hence it may be considered as an observational variable, rather than a state-space one.

This paper is organized as follows: The second section introduces the new ITCC model, implemented in TOMION v3. Then, the third section presents and discusses some results with field data from a geodetic network in Poland. Finally, section four presents the conclusions and future work.

2 ITCC model

The original system of equations for UUMs is rank deficient. Consequently, not all the original parameters can be estimated. Some strategies can overcome this problem by combining and/or differencing observables (PPP and RTK, respectively). The model presented in this paper is based on the linear combination of estimables, not observables, thereby fixing the rank deficiency. Namely, the model consists in computing the satellite clocks, the receivers clocks and their hardware biases with respect to a common reference receiver (i.e. the pivot) in the network, alongside with a tomographic sounding of the ionosphere. The system of equations for the network (i.e. the CPF) consists in modelling the code $P(t)_{i,j}^k$ and carrier-phase $L(t)_{i,j}^k$ observables as follows:

$$\begin{aligned} \begin{array}{ccl} P(t)_{i,j}^k &{} = &{} \rho (t)_j^k+c\cdot \left[ dt(t)_j-dt(t)^k\right] +M_T(t)_j^k\cdot Z_T(t)_j\\ &{}&{}+\frac{\beta }{f_i^2}S(t)_j^k+\frac{\beta }{f_i^2}\left[ D(t)_j-D(t)^k\right] \\ L(t)_{i,j}^k &{} = &{} \rho (t)_j^k+c\cdot \left[ dt(t)_j-dt(t)^k\right] +M_T(t)_j^k\cdot Z_T(t)_j\\ &{}&{}-\frac{\beta }{f_i^2}S(t)_j^k+B_{i,j}^{k,m}+\frac{c}{f_i}\phi (t)_j^k \\ \end{array} \end{aligned}$$

(1)

where the indexes i, j and k refer to the ith frequency, the jth receiver and the kth satellite, respectively. The geometric distance from the jth receiver to the kth satellite, $\rho (t)_j^k$, can be linearized in the usual way in terms of an approximate position of the receiver. $dt_j$ and $dt^k$ are the receiver and satellite clock errors, respectively. The tropospheric delay is computed as the product of the tropospheric mapping function $M_T(t)_j^k$ and the vertical tropospheric delay $Z_T(t)_j$ for both the hydrostatic and nonhydrostatic components. And the simultaneous estimation of the corresponding tropospheric gradients can be activated, specially suitable in extreme weather scenarios, see for instance Graffigna et al (2019b). The ionospheric delay consists in the STEC, $S(t)_j^k$, rescaled with the parameter $\beta /f_i^2$, where $\beta \simeq 40.3$ (S.I.). Then, the receiver and transmitter Delay Code Biases (DCB) are $D_j$ and $D^k$, respectively; the carrier phase bias (hereinafter phase ambiguity) B for the continuous (i.e. without cycle-slips) phase arc “m” including the unknown integer number of cycles and the non-integer parts of the receiver and transmitter at the given frequency; and, finally, the wind-up angle $\phi$.

Since we have analysed $P_1,P_2$-equivalent GPS measurements in the first ITCC application, presented below, we have taken as code bias “datum” in equation 1 the typical one: the ionospheric-free combination of them equals to zero, i.e. we estimate the receiver and satellite delay code biases of the geometry-free combination ($D_j$ and $D^k$, respectively). In future ITCC applications, to more than two signals and from different GNSS, this condition can be replaced by taking as zero both receiver and transmitter code biases for the first frequency processed ($D_{1,j}\equiv 0$, $D_1^k\equiv 0$) and estimating directly (except for the reference receiver) the k transmitter and j receiver biases for the remaining signals ($D_{i,j}$, $D_i^k$).

To better understand the ITCC model it should be taken into account that during the short “life-time” of a given carrier phase ambiguity unknown (i.e. during the phase-lock time, up to some hours for the predominant MEO GNSS transmitters) the ambiguity is assumed to be time-constant, similarly to the standard PPP method. And the idea of “constraining the ambiguities” means the same that it meant in the previous combined approach (WARTK) but implemented in a most straightforward way in ITCC:

• 1. In the previous WARTK combined approach, the geometric (ionospheric-free) and ionospheric (geometry-free) equations were coupled, i.e. mutually constrained, by the Melbourne–Wubbena combinations, once it is properly smoothed. This produces as benefit the acceleration of the convergence of the ionospheric model computed from the permanent network of stations Hernández-Pajares et al (2002) and the acceleration of the convergence of the geometric model computed by the roving user (e.g. in positioning in Hernández-Pajares et al (2000) and in zenith tropospheric delay determination in Hernández-Pajares et al (2001)).

• 2. In the uncombined ITCC approach presented in this manuscript, the simultaneous set of equations contain both type of unknowns (ionospheric and geometric ones) directly combined in a natural way. As a consequence the connection between both geometric and ionospheric models played in WARTK by the Melbourne–Wubbena combination, is not needed anymore in ITCC. And the mutual constrain between both “submodels” (ionospheric or geometric) is implicit and depends on the corresponding uncertainty: usually the geometric sub-model is benefitted from the ionospheric (received) one by the roving user, and the reversal for the permanent network of GNSS receivers.

2.1 Identification and removal of the rank deficiency

As it stands now, Eq. (1) is rank-deficient. In other words, not all the parameters can be projected onto the space of observables by the Jacobian matrix. This means that when such parameters are shifted by an arbitrary constant, the observables remain invariant. For example, by adding the same constant a in the receiver and satellite clocks, $dt_j+a$ and $dt^k+a$, respectively, the code and carrier-phase do not change. Similarly, any arbitrary constant on the receiver and satellite DCBs does not modify the observables. There are many ways to remove these two types of rank deficiency. Several examples are given in Odijk et al (2016) and Teunissen and Khodabandeh (2015). In this paper, we have chosen the so-called Common-Clock Model (CCM), whereby the clock and DCB of one receiver are taken as reference. Then, the rank defect is fixed by referring the clock errors and pseudorange hardware biases to the corresponding values at the reference receiver:

$$\begin{aligned} \begin{array}{ccc} dt(t)_0 &{} \equiv &{} 0\\ D(t)_0 &{} \equiv &{} 0\,. \end{array} \end{aligned}$$

(2)

There is still another rank defect: that between the ionospheric delay STEC, clocks and DCBs. The size of this rank defect is $m+n-1$, where n and m are the number of receivers and satellites, respectively. If not removed, this rank defect would yield estimated STEC parameters biased by their satellite and receiver DCBs, see for example Section 5.1 in Odijk et al (2016). In combination with Eq. (2), then the rank defect of the ionospheric delay would be with either the clocks or the DCB parameters.

In this case, an arbitrary constant a on the STEC can be nullified by shifting, for example, the satellite clocks by the term $(\beta /f_i^2)\cdot a$, hence the rank defect.

In general, the ionospheric delay I is obtained by rescaling an ionospheric parameter ${\mathcal {I}}$ as follows:

$$\begin{aligned} I \equiv {\mathcal {J}}\cdot {\mathcal {I}}\,, \end{aligned}$$

(3)

where ${\mathcal {J}}$ is the Jacobian coefficient that maps ${\mathcal {I}}$ onto the observational space. Should the ionospheric parameter ${\mathcal {I}}$ be shifted by a constant a, then, in order to cancel it, the clock parameter should be shifted by ${\mathcal {J}}\cdot \displaystyle \frac{\beta }{f_i^2}\cdot a$, thus introducing the following shift term r in Eq. (1):

$$\begin{aligned} r\equiv \left( {\mathcal {J}}-1\right) \cdot \displaystyle \frac{\beta }{f_i^2}\cdot a\,. \end{aligned}$$

(4)

Note that when $I=STEC$, then ${\mathcal {I}}=STEC$ and ${\mathcal {J}}=1$, therefore, $r=0$ as expected. Consequently, the modelled observationals P and L remain invariant under such shift a on both parameters STEC and $dt^s$. This is equivalent to say that both parameters STEC and $dt^s$ are correlated.

However, when $I\equiv M_I\left( {\mathcal {E}}\right) \cdot VTEC$, which is the case for the thin-shell ionospheric model, then the shift term is the following:

$$\begin{aligned} r\equiv \left( M_I\left( {\mathcal {E}}\right) -1\right) \cdot \displaystyle \frac{\beta }{f_i^2}\cdot a\,, \end{aligned}$$

(5)

where ${\mathcal {J}}\equiv M_I\left( {\mathcal {E}}\right)$ is the ionospheric mapping function, and ${\mathcal {E}}$ the elevation angle of the satellite above the local horizon. In such case, $r\ne 0$, and the modelled observables P and L would be affected as well. Consequently, the ionospheric thin-shell model fixes the rank, as expected as well, Odijk et al (2016). In terms of correlation, we say that they are decorrelated due to the change in the relative geometry between receiver and satellite.

Although the thin-shell model fixes the rank, its implementation requires short inter-receiver distances of about $70\,km$, thus hampering its application for larger networks (e.g. wide area networks, or continental networks of low receiver density in Australia, or either worldwide networks such as the IGS network of permanent receivers).

The mismodelling that the thin-layer hypothesis introduces in positioning is too large for models that aim at cm-level accuracy, Banville and Langley (2013). Furthermore, this mismodelling increases as the user approaches to the geomagnetic equator Hernández-Pajares et al (1999).

The tomographic approach has proven to be a more accurate ionospheric model than the thin-shell one, Juan et al (1997). Instead of estimating STEC for each satellite-receiver pair, the ionospheric electron density is estimated inside small regions, i.e. voxels. The STEC is the integration of the ionospheric electron density $N_e$ along the line of sight between the satellite and the receiver:

$$\begin{aligned} S(t)_j^k = \sum _{p=1}^{n(t)_k^j} N_e(t)_{m(p)}\cdot \delta l(t)_{m(p),j}^k\,, \end{aligned}$$

(6)

where $n(t)_k^j$ is the index number of the illuminated voxels along the line of sight between the receiver j and the satellite k at time t, $\delta l(t)_{m(p),j}^k$ is the section of the LOS that goes through the $p^{th}$ illuminated voxel, which index is m(p).

Now, following Eq. (3), in the tomographic model, the ionospheric coefficients of the Jacobian are as follows

$$\begin{aligned} {\mathcal {J}} = \displaystyle \left( \frac{\beta }{f_i^2}\right) \cdot \delta l(t)_{m(p),j}^k\,, \end{aligned}$$

(7)

and the set of ionospheric parameters ${\mathcal {I}}$ consists of the ionospheric electron density values inside each voxel:

$$\begin{aligned} {\mathcal {I}} = N_e(t)_{m(p)}\,. \end{aligned}$$

(8)

Consequently, the shift term r is the following expression:

$$\begin{aligned} r = \left( \delta l(t)_{m(p),j}^k-1\right) \cdot \displaystyle \frac{\beta }{f_i^2}\cdot a \,. \end{aligned}$$

(9)

Due to the time-dependence of $\delta l(t)_{m(p),j}$, the shift term r cannot be permanently cancelled by the constant a. Consequently, the modelled observables will not remain invariant under such shift. Thus, the tomographic modelling also fixes the rank defect between the ionospheric and hardware delays.

Similarly to the thin-layer model, the evolution of the relative geometry between satellites and receivers is what allows for the decorrelation between the ionospheric delay, clocks parameters, and DCBs.

2.2 The ITCC model

The combination of Eqs. (1) and (6) yields the Ionospheric Tomographic CC (ITCC) model as follows:

$$\begin{aligned} P(t)_{i,j}^k= & {} \rho (t)_j^k+c\cdot \left[ dt(t)_j-dt(t)^k\right] +M_T(t)_j^k\cdot Z_T(t)_j\nonumber \\&+\frac{\beta }{f_i^2}\left[ \sum _{p=1}^{n(t)_k^j} N_e(t)_{m(p)}\cdot \delta l(t)_{m(p),j}^k\right] +\frac{\beta }{f_i^2}\left[ D(t)_j-D(t)^k\right] \nonumber \\ L(t)_{i,j}^k= & {} \rho (t)_j^k+c\cdot \left[ dt(t)_j-dt(t)^k\right] +M_T(t)_j^k\cdot Z_T(t)_j \nonumber \\&-\frac{\beta }{f_i^2}\left[ \sum _{p=1}^{n(t)_k^j} N_e(t)_{m(p)}\cdot \delta l(t)_{m(p),j}^k\right] +B_{i,j}^{k,m}+\frac{c}{f_i}\phi (t)_j^k \,, \nonumber \\ \end{aligned}$$

(10)

where the ionospheric electron densities within each voxel, $N_e(t)_{m(p)}$, are the new ionospheric estimables, thus replacing the STEC in the parameters set.

With the tomographic approach, the ionospheric parameters $N_e(t)_{m(p)}$ can be geometrically decorrelated each other accordingly to the relative geometry of the satellite-receiver pair and its local contribution for each ionospheric voxel via the Jacobian coefficient $\delta l(t)_{m(p),j}^k$.

Tomographic models have proved to perform well for inter-receiver distances ranging from $100\,km$ up to $3000\,km$, Hernández-Pajares et al (2000).

The Jacobian matrix of CC models that estimates either the STEC (SCC) or the VTEC (VCC) scales geometrically with the number of receivers and satellites, n and m, respectively (see Eq. (35) in Odijk et al (2016)). However, for the ITCC model, once the number of voxels p is set up by the spatial resolution, the addition of new receivers and/or constellations would increase the number of clocks and hardware delays to be estimated, but not the number of ionospheric parameters for the ITCC model. Therefore, once the tomographic model achieves the desired spatial resolution for the ionospheric voxels, the number of parameters to be estimated scales as $n+m$. Consequently, for large GNSS networks, the computational burden is lower for ITCC than for SCC and VCC models.

All together considered, the ITCC model fixes the rank defect for the ionospheric parameters via geometric decorrelation, it reduces the ionospheric gradient mismodelling present in the thin-shell models, reduces the computational loading of the Jacobian matrix with respect to SCC and VCC models, and allows for the estimability of DCBs. Table 1 summarizes the advantages of the ITCC model over SCC and VCC models.

However, if the size of the network is not large enough, about $70\,km$ interstation distance, then the ionospheric coefficients in the Jacobian would be too similar, up to the point that it would introduce a rank-defect. This problem also affects SCC and VCC models, Odijk et al (2016).

Table 1 Comparison of models in terms of their algebraic features and applicability in large networks (✓and ✗indicates whether the given condition is fulfilled or not, respectively)

2.3 The GNSS equations for the roving users

Given that the satellites’ clocks and hardware delays, and ionospheric parameters are computed by the network, the system of equations to be solved for roving users will be as follows:

$$\begin{aligned} \begin{array}{ccl} {\tilde{P}}(t)_{i,rov}^k &{} = &{} -{\hat{\rho }}(t)_{rov,0}^k\cdot \Delta \mathbf {X}+c\cdot dt(t)_{rov}+M_T(t)_{rov}^k\\ &{}&{}\cdot Z_T(t)_{rov}+\frac{\beta }{f_i^2}D(t)_{rov} \\ {\tilde{L}}(t)_{i,rov}^k &{} = &{} -{\hat{\rho }}(t)_{rov,0}^k\cdot \Delta \mathbf {X}+c\cdot dt(t)_{rov}+M_T(t)_{rov}^k\\ &{}&{}\cdot Z_T(t)_{rov}+B_{i,rov}^{k,m}+\frac{c}{f_i}\phi (t)_{rov}^k\,, \\ \end{array} \end{aligned}$$

(11)

where ${\hat{\rho }}(t)_{rov,0}^k$ is the unity vector pointing the GNSS transmitter from the approximate position $\mathbf {X}_0$ of the rover and $\Delta \mathbf {X}=\mathbf {X}-\mathbf {X}_0$ is the position vector correction to be estimated. And ${\tilde{P}}(t)_{i,rov}^k$ and ${\tilde{L}}(t)_{i,rov}^k$ are the corrected code and carrier-phase (AKA pre-fit residuals), respectively. Namely,

$$\begin{aligned} \begin{array}{ccl} {\tilde{P}}(t)_{i,rov}^k &{} \equiv &{} P(t)_{i,rov}^k - \rho (t)_{rov,0}^k + c\cdot dt(t)^k\\ &{}&{} - \frac{\beta }{f_i^2}S(t)_{rov}^k + \frac{\beta }{f_i^2}D(t)^k\\ {\tilde{L}}(t)_{i,rov}^k &{} \equiv &{} L(t)_{i,rov}^k - \rho (t)_{rov,0}^k + c\cdot dt(t)^k + \frac{\beta }{f_i^2}S(t)_{rov}^k\,, \\ \end{array} \end{aligned}$$

(12)

where $\rho (t)_{rov,0}^k$ is the distance from the GNSS transmitter to the initial position guess $\mathbf {X}_0$, $S(t)_{rov}^k$ is computed using Eq. (6) with the ionospheric electron density parameters computed by the network, $Ne(t)_p$.

Fig. 1

Map representing the permanent receivers in Poland, containing the X-axis and Y-axis the longitude and latitude, respectively, in degrees. Data from days 249–255 in year 2018 (green squares) has been analysed in this work, as well as the reference receiver (red square) and permanent GPS ones treated as rover users (magenta circles)

Fig. 2

Cumulative distribution of distances between the permanent receivers in Poland being considered during days 249–255 of year 2018

Fig. 3

Cumulative distribution of distances between the receiver SZAM, treated as roving user in this study, regarding to the permanent receivers in Poland, considered during days 249–255 of year 2018

Fig. 4

Error in the East coordinate for the further away roving GNSS user SZAM during the first day (first column) and last day (second column) of the analysed period, 249–255 of year 2018, after forcing user cold start every 7200 s. It is shown for WARTK without constraining nor fixing the carrier phase ambiguities, e.g. somehow equivalent to PPP (first row), ITCC without ambiguity fixing (second row) and WARTK with double-differenced ambiguities fixed (third row)

Additionally, the roving users will use the ultrarapid precise predicted orbits, available hours in advance from the International GNSS Service (IGS, Johnston et al (2017)). Therefore, there is no need for providing external real-time products to the network, making the service of the wide area network autonomous for long periods of time. All the geodetic and atmospheric parameters are simultaneously provided by the network, thus reducing the precise positioning convergence time for roving GNSS users with respect to users that compute their own ionospheric delays, Olivares-Pulido et al (2019). In this regard, the following section presents an example for a Wide Area geodetic network in Poland.

Fig. 5

Error histograms in the vertical coordinate, for the further away roving GNSS user SZAM during the overall 7-days period, 249–255 of year 2018 (first column) and just 300 seconds after every cold start (second column). It is shown for WARTK without constraining nor fixing the carrier phase ambiguities, e.g. somehow equivalent to PPP (WAPPP, first row), ITCC without ambiguity fixing just constraining (second row) and WARTK with double-differenced ambiguities fixed (third row)

Fig. 6

Cumulative distribution function of the convergence time for roving user SZAM to reach horizontal error of 10 cm (first column), 20 cm (second column) and 30 cm (third column) with WAPPP (first row), ITCC (second row) and WARTK with double-differenced ambiguity fixing (third row)

3 A first example of ITCC: Wide Area precise GNSS positioning in Poland

In this section, we present a summary of the first ITCC application, performed over a geodetic network in Poland, in terms of Wide Area positioning error and convergence time, during days 149-155 of year 2018, emulating real-time conditions. For a better assessment, the same data set was previously analysed by WARTK method, its CD counterpart model.

The network consists in the permanent receivers represented under green squares in Fig. 1, including the reference receiver (j=0) WLOC close to the centre of the network, and treating the permanent receiver SZAM as ITCC real-time roving user. It can be seen in Fig. 2 that the 95% of the distances are greater than about 130 km, with less than 3% of baselines shorter than 100km. And in Fig. 3 it is shown the receiver SZAM, treated as real-time kinematic user, is placed at typical Wide Area user distances from the permanent sites (around 120 km to the nearest one).

The overall receiver state of roving user, SZAM, is fully reset every 2 hours (7200 s), emulating a new cold start, up to a total of 84 convergence events in the analysed period comprising 7 days. This high number of independent full user resets allows to properly characterize both the positioning error and the convergence time, comparing: (1) A similar undifferenced but combined and not using ionospheric corrections technique (hereinafter Wide Area PPP, WAPPP); (2) ITCC, and (3) WARTK, a combined and undifferenced geodetic model, along with the tomographic ionospheric model, used only in the ambiguity fixation events through the Melbourne–Wübbena combination Hernández-Pajares et al (2000). WARTK should provide the best performance, in particular than this first implementation of ITCC equivalent to constraining only, and not fixing, the carrier phase ambiguities and code instrumental delays.

Focusing on the users on land (e.g. precise farming), a first glance of the behaviour of ITCC on horizontal positioning error can be seen in Fig. 4, for the East (E) coordinate of the roving user SZAM in Poland at +120 km from the permanent receivers, the most sensitive component to carrier phase ambiguity fixing or constraining in general (see for example Ge et al (2008) and Geng et al (2019)). It can be seen as the performance of the plain ITCC (without fixing ambiguities for the time being in this formulation) is at short term (less than 3000 sec after the cold start) in between the standard precise positioning (WAPPP) and the WARTK after fixing double-differenced carrier phase ambiguities. The higher noise of ITCC at the end of the convergence period, compared not only with WARTK under ambiguity fixing, but also with WAPPP, can be understood for the higher number of unknowns in the CPF model. As it has been discussed above this noise can be significantly reduced in future applications, by processing other GNSS systems beyond GPS, because under ITCC the number of ionospheric unknowns (electron density values of voxels, predominant in the set of unknowns) should not suffer significant modifications.

Table 2 Approximate percentile 95% of the convergence time in seconds of receiver SZAM treated as roving user to horizontal positioning error of 10 cm (first column), 20 cm (second column) and 30 cm (second column), performing a cold start every two hours during days 149-155 of year 2018

Fig. 7

Comparison of Zenith Tropospheric Delay computed with the ITCC technique implemented in TOMION (red) with the global uncombined processing fixing the carrier phase ambiguities by using FUSING software, in blue (day 149 of year 2018)

In Fig. 5 it can be seen the distribution of vertical positioning errors for the same roving-like receiver SZAM but for the whole studied period of 7 days, treated in the same way as roving user, compared for the same three techniques in the same row order, for all the epochs (left-hand column), and just from 300 second and after each reset performed every 2 hours (right-hand column). It can be seen in the error RMS values that the coordinate most affected by the geometric dilution of precision on ground receivers improves from 0.8 m of vertical error RMS under traditional Wide Area Precise Positioning (WAPPP) to 0.5 m of vertical error RMS when the plain ITCC is applied, due to the very quick reduction in larger errors in ITCC. Nevertheless when we consider only the position after 5 minutes of each cold start, then the WAPPP improvement (0.1 m of error RMS) is higher than for ITCC (0.2 m of error RMS), consistent with a certain higher noise after convergence already seen in previous Fig. 4 for E coordinate. The origin, likely related with the high number of unknowns associated with the direct usage of tomographic voxel model for the STEC term S, and its potential mitigation by using the four GNSSs, can be confirmed in next studies.

However, in terms of convergence time of the horizontal coordinate (hereinafter CTH), ITCC performs clearly better than WAPPP. For instance in the 95% of cases the CTH is smaller than 1350s for ITCC vs less than 4000s for WAPPP to reach 10 cm; and 250s vs 700s to reach 30cm, respectively (see Fig. 6). As expected, WARTK with double-differenced ambiguity fixing provides a significantly better performance to reach 10 cm (95% of cases with CTH smaller than 750s) and somehow better to reach 30 cm (95% of cases with CTH smaller than 125s), compared with ITCC (250s). These and other details are summarized in Table 2.

We have performed an additional assessment, focused in this case on the permanent receivers, and on the zenith tropospheric delay (ZTD) estimation, a parameter which can be considered a quality factor as it is described below.

It can be seen in Fig. 7 the ZTD estimated for some representative permanent GNSS receivers by using: (1) TOMION implementation of the plain ITCC, i.e. floating the ambiguities in an undifferenced uncombined approach at Wide Area scale estimating the ionospheric delay with an ionospheric tomographic model; and (2) the ZTD computed with the independent software FUSING, Gu et al (2018), by implementing the global uncombined undifferenced ambiguity-fixing processing by estimating directly the STEC as random walk, initialized with the values provided by the CODE GIM. It can be seen an agreement, typically at 1cm-level and with larger deviations at some periods, but noisier in ITCC consistently with the positioning results already discussed. This is the case, taking into account as well that Fig. 7 corresponds to the first day of the time series, i.e. including the convergence of the Central Processing Facility filter during the first half to one hour approximately.

4 Conclusions

The Ionospheric Tomographic Common Clock (ITCC) model of undifferenced uncombined GNSS measurements has been introduced and tested. The model has been implemented for a Wide Area network and roving users.

It has been proved that the tomographic approach fixes the rank deficiency of the ionospheric estimates, thus allowing for the computation of ionospheric electron density, clocks and DCBs. Although ionospheric thin-shell models also fix the rank deficiency, they do so at the expense of increasing the ionospheric gradient mismodelling. Additionally, the number of parameters in the ITCC model increases linearly with the number of satellites and receivers, as opposed to SCC and VCC models, which scale geometrically.

The ITCC model, without fixing the carrier phase ambiguities for the time being (just constraining them), is assessed and compared with previous precise positioning techniques during days 249-255 of year 2018 in Poland. On the one hand the proposed approach is 3 times faster than similar approach with no ionospheric corrections, and on the other it requires twice the time than WARTK with ambiguity fixing, to achieve 10 cm of horizontal error. And in despite the final accuracy is worse than in both previous approaches, coinciding with the larger number of unknowns, which includes as well the electron densities of voxels, it can be significantly mitigated in next step. Indeed, this can be achieved by including the additional GNSS measurements, not just the GPS ones, taking profit that in the present formulation the number of unknowns will remain the same. A final boost of both, convergence time and final accuracy, is expected in a further step by incorporating the carrier phase ambiguity fixing.

The tomographic approach allows the network model to simultaneously compute geodetic and ionospheric parameters. Namely, the network can estimate and then provide satellites clocks and hardware delays, and ionospheric corrections to the rover model, thus further increasing its convergence speed for precise positioning. Both models for the rover and the network can be implemented for a Wide Area GNSS network without needing external real-time products.

Finally, it has been shown that, without fixing the ambiguities, this model is faster than a model based on differenced and combined observables. Furthermore, it is expected to further increase the positioning accuracy by fixing the ambiguities in ITCC in future works, characterizing as well the ITCC performance dependence on different regions (e.g. vs geomagnetic latitude) and periods (e.g. vs solar cycle phase and seasons).