Abstract
Single-station PPP-RTK is a special case of PPP-RTK in that corrections are computed, instead of a network, by only one single GNSS receiver. The present contribution aims to develop a framework to generate multi-epoch, single-station corrections, thereby providing PPP-RTK users the capability to time-predict corrections that are subject to time delay or latency. By presenting analytical expressions of the user ambiguity variance matrix, we address how the ambiguity resolution performance is driven by the correction latency and therefore by the uncertainty involved in the time-prediction of single-station PPP-RTK corrections. Supported by numerical results, our analytical study shows that the number of satellites and number of frequencies work in tandem to enable one to increase the correction latency, yet ensuring successful single-receiver ambiguity resolution.
Similar content being viewed by others
1 Introduction
Integer ambiguity resolution enabled precise point positioning (PPP-RTK) extends the PPP concept (Heroux and Kouba 1995; Zumberge et al. 1997) by providing single-receiver users, next to the network-derived orbits and clocks, also information about the satellite biases and (optionally) atmospheric delays (Wubbena et al. 2005; Teunissen et al. 2010). This information, when properly provided, enables recovery of the integerness of the user ambiguities, thus enabling single-receiver ambiguity resolution. In this contribution we discuss a special case of PPP-RTK, single-station PPP-RTK, in which the network consists of only one single reference receiver. Particular attention will be given to the role of the correction latency, i.e. the delay in time after the PPP-RTK corrections are generated and the time they are applied to the user data.
The concept of single-station PPP-RTK was introduced by Khodabandeh and Teunissen (2015), in that a proper set of estimable parameters of a single GNSS station were shown to act as if they are the satellite clocks and phase/code biases, respectively. Similar to its network-based variants presented in, e.g. (Ge et al. 2008; Laurichesse et al. 2009; Collins et al. 2010; Bertiger et al. 2010; Odijk et al. 2015), single-station PPP-RTK has the potential of benefiting from an efficient state-space packaging of the corrections to reduce their transmission rate, i.e. the frequency with which the corrections are to be provided to the user. A reduction in the transmission rate comes at the cost of delayed corrections. Consequently, the user is required to time-predict the corrections so as to bridge the gap between the corrections’ generation time and the user positioning time. It is therefore the goal of the present contribution to develop a framework to generate and time-predict multi-epoch, single-station PPP-RTK corrections, thereby addressing the impact the correction latency has on the user ambiguity resolution performance.
The importance of such single-station framework is threefold, which is well appreciated in the light of the recent proliferation of low-cost GNSS sensors (Banville et al. 2019). Firstly, as the PPP-RTK corrections are network-adjusted versions of their single-station counterparts (Khodabandeh and Teunissen 2015), a careful study of their multi-epoch, single-station formulation provides insight into the quality of such corrections. Secondly, as the rapid growth in the number of network receivers demands an efficient way of data processing in terms of computational capacity, one may choose to decentralize the underlying data processing and distribute the task of the computing centre across individual network receivers (Khodabandeh and Teunissen 2019). This, in turn, paves the way for various cooperative and peer-to-peer applications, see e.g. (Deambrogio et al. 2010). Thirdly, the single-station PPP-RTK concept offers the capabilities of network-based correction providers to individual single-station receivers, thereby realizing low-cost positioning infrastructure.
This contribution is organized as follows. In Sect. 2, we discuss the fine distinction between single-station PPP-RTK and single-baseline RTK. How their respective corrections are related to one another is addressed. The advantage of using single-station PPP-RTK comes to the fore when one exploits the distinct temporal properties of the GNSS individual parameters to send them in their estimable form, each with its own transmission rate, to the user. In Sect. 3, we show how this can be realized through an application of \({\mathcal {S}}\)-system theory (Baarda 1973; Teunissen 1985). For a class of polynomial dynamic models, single-receiver estimable parameters of both the PPP-RTK provider and user are identified. These estimable parameters are shown to inherit the temporal properties of their original counterparts. In Sect. 4, we present closed-form expressions for the user ambiguity variance matrix and analyse the impact the latency has on the single-receiver ambiguity resolution performance. Our analytical study benefits from the single-epoch building blocks of the ambiguity variance matrix. In Sect. 5, we quantify the extent to which the precision of the user ambiguity solutions differs from that of their minimum-variance versions, when the user ignores the corrections’ uncertainty by treating them as nonrandom, i.e. deterministic. Finally, a summary with conclusions is provided in Sect. 6.
2 RTK versus single-station PPP-RTK
The ultimate goal of both network-RTK and PPP-RTK is to deliver integer ambiguity-resolved solutions for single-receiver users. The fine distinction of these two methods, and their single-station versions, can be best explained through the effect of their products serving as corrections on the single-receiver user observations. We hereby start with the user observation equations and highlight why single-receiver user ambiguity resolution cannot stand on its own, requiring corrections from an external provider.
2.1 Single-receiver observation equations
Let \(\varDelta \phi _{u,j}^s\) and \(\varDelta p_{u,j}^s\), respectively, denote the observed-minus-computed phase and code observables of satellite s (\(s=1,\ldots ,m\)) on frequency j (\(j=1,\ldots ,f\)) that are collected by the user receiver u. If we make use of the multi-frequency vector notation \(\varDelta \phi _u^s\!=\![\varDelta \phi _{u,1}^s\!,\!\ldots \!,\!\varDelta \phi _{u,f}^s]^T\) and \(\varDelta p_u^s\!=\![\varDelta p_{u,1}^s,\ldots ,\varDelta p_{u,f}^s]^T\), the user observation equations can be expressed as (Khodabandeh and Teunissen 2015)
with \(\mathsf {E}(\cdot )\) being the expectation operator. The receiver-satellite range increment \(\varDelta \rho _{u}^{s}\) is parameterized in terms of the \(\nu \)-vector \(\varDelta x_u\) containing the user’s position increments and/or the zenith tropospheric delay (ZTD). Parameter \(\nu \) can take the values \(\nu =3\) (position-only model), \(\nu =1\) (ZTD-only model) or \(\nu =4\) (position-plus-ZTD model). Thus, \(\varDelta \rho _{u}^{s}={g_u^{s}}^T\varDelta x_u\), with the \(\nu \)-vector \(g_u^{s}\) containing the receiver-satellite direction vector and/or the tropospheric mapping function. Here and in the following, the precise orbital corrections are assumed included in the observed-minus-computed observables. The diagonal matrix \(\varLambda =\mathrm{diag}(\lambda _{1}, \ldots , \lambda _{f})\), formed by the wavelengths \(\lambda _j\), links the phase observables to the integer ambiguities \(z_u^s\!=\![z_{u,1}^s,\ldots ,z_{u,f}^s]^T\) as well as the receiver and satellite phase biases \(\delta _u\!=\![\delta _{u,1},\ldots ,\delta _{u,f}]^T\) and \(\delta ^s\!=\![\delta _{,1}^s,\ldots ,\delta _{,f}^s]^T\), respectively. The common receiver and satellite clocks are denoted by \(dt_u\) and \(dt^{s}\), while the receiver and satellite code biases are denoted by \(d_u\!=\![d_{u,1},\ldots ,d_{u,f}]^T\) and \(d^s\!=\![d_{,1}^s,\ldots ,d_{,f}^s]^T\). The term \(\iota _u^{s}\) denotes the first-order slant ionospheric delay experienced on the first frequency. The corresponding coefficients are defined as the ratios \(\mu _j=(\lambda _j^2/\lambda _1^2)\) forming the f-vector \(\mu =[\mu _1,\ldots ,\mu _f]^T\). The f-vector of ones is denoted by \(e=[1, \ldots ,1]^{T}\). Apart from \(z_{u}^{s}\), \(\delta _{u}\) and \(\delta ^{s}\) that are expressed in cycles (and the dimensionless direction vector \(g_u^s\)), the rest of the quantities are all in units of length.
The user observation equations (1) do not contain enough information to solve for an integer ambiguity-resolved user position. The reason being that (i) the range increments \(\varDelta \rho _{u}^{s}\) are biased by the satellite clocks \(dt^s\) prohibiting unbiased estimation of the user position \(\varDelta x_u\) and that (ii) the integer ambiguities \(z_u^s\) are biased by the satellite phase biases \(\delta ^s\) prohibiting integer ambiguity resolution. Note that the presence of the receiver clock \(dt_u\) and phase biases \(\delta _u\) is not of concern, as they can be eliminated by forming between-satellite differences of the observables. The user thus takes recourse to an external provider so as to acquire information on the satellite clocks \(dt^s\) and phase biases \(\delta ^s\).
For both network-RTK and PPP-RTK, the role of such provider is taken by a network of receivers (Odijk et al. 2015). Special cases are single-baseline RTK and single-station PPP-RTK, i.e. when the network consists of only one single reference receiver. The observation equations of such reference receiver, say receiver r, follow by replacing the user index u with r in (1) and setting \(\varDelta \rho _{r}^{s}=0\), assuming that the reference receiver position (ZTD) is a priori known (corrected), i.e. \(\varDelta x_r=0\). Thus, the reference receiver observation equations
form the basis for the provision of information on the satellite clocks \(dt^s\) and phase biases \(\delta ^s\).
2.2 State-space representation of the corrections
The goal is now to find a class of corrections with which the user can correct his GNSS observables and obtain ambiguity-resolved positioning solutions. As the reference receiver observables, in (2), contain satellite clock and phase bias information, it seems reasonable to take the vector \([\varDelta \phi _{r}^{s^T},\varDelta p_{r}^{s^T}]^T\) as a member of such class. Subtracting the stated vector from the user observables in (1) gives
in which use is made of the between-receiver differencing notation \((\cdot )_{ru}=(\cdot )_u-(\cdot )_r\). As the corrections cancel the satellite-specific terms \(dt^s\), \(\delta ^s\) and \(d^s\), an integer ambiguity-resolved user position solution would follow by solving, e.g. the between-satellite differences of the user-corrected observation equations (3). The correction vector \([\varDelta \phi _{r}^{s^T},\varDelta p_{r}^{s^T}]^T\) is, however, not the only vector that can bring the user observation equations (1) into the same structure as that of (3). In fact, if \(\varDelta \phi _{r}^{s}\) and \(\varDelta p_{r}^{s}\), on the left-hand side of (3), are replaced by their arbitrary weighted least-squares adjusted versions, the right-hand side remains unchanged. This can be understood as follows. Let \(\varDelta \hat{\phi }_{r}^{s}\) and \(\varDelta \hat{p}_{r}^{s}\), indicated by the \(\hat{.}\)-symbol, denote least-squares adjusted phase and code observables of the reference receiver r, respectively. The claim follows then from the least-squares unbiasedness properties \(\mathsf {E}(\varDelta \hat{\phi }_{r}^{s})=\mathsf {E}(\varDelta \phi _{r}^{s})\) and \(\mathsf {E}(\varDelta \hat{p}_{r}^{s})=\mathsf {E}(\varDelta p_{r}^{s})\).
The difference between the adjusted observables \(\varDelta \hat{\phi }_{r}^{s}\) and \(\varDelta \hat{p}_{r}^{s}\) and their original counterparts \(\varDelta \phi _{r}^{s}\) and \(\varDelta p_{r}^{s}\) is due to the presence of single-station model’s redundancy (Teunissen 2000a). The redundancy is defined here as the number of model’s observables minus the number of its estimable parameters (which will be treated in Sect. 4.1). In the absence of redundancy, however, the adjusted observables reduce to the original ones, showing that the adjusted observables encompass the original observables as special cases. A class of single-station corrections is therefore characterized as follows
Since the clock and phase-bias information, upon (4), is conveyed to the user in terms of (adjusted) observables, the above expression of the corrections is referred to as ‘observation-space representation’ (OSR). Such representation is commonly employed in RTK (Wubbena et al. 2005).
As with any real-time GNSS applications, RTK users need to have instant access to the corrections (4), a situation that often demands high data transmission. The data transmission can, however, be optimized through the minimization of the corrections’ transmission rate, i.e. the frequency with which the corrections are to be provided. This entails careful modelling of the corrections’ temporal behaviour so as to enable the user to time-predict the corrections, thereby bridging the gap between the corrections’ generation time and the positioning time (Wang et al. 2017). To get insight into the temporal characteristic of the corrections (4), let us first consider how they are formed by the individual parameters given in (2). From (2) and the unbiasedness properties of the adjusted observables follows that
Accordingly, the temporal behaviour of the corrections (4) is not merely governed by the time-stable parameters (i.e. the ambiguities and biases), but by the highly time-varying parameters (i.e. the clocks and ionospheric delays). As a consequence, the transmission rate of the corrections (4) is required to be rather high, often higher than 0.1–1 Hz for sub-cm positioning (Wubbena et al. 2005; Wang et al. 2017).
Provision of state-space corrections over time. The bias and clock/ionospheric solutions are, respectively, provided to the user, \(k_o\) and k epochs after the filter’s initialization (\(k\ge k_o\)). Collecting GNSS measurements at epoch \(l\ge k\), the user is required to time-predict the corrections (cf. Sect. 3.4)
The distinct temporal properties of the parameters, in (5), suggests sending individual parameter solutions, each with its own transmission rate, to the user. The user can then form the corrections (4) as combinations of such individual solutions (see Fig. 1). The information content of the observation equations (2) and (5) is, however, not sufficient to unbiasedly determine the individual parameters (Odijk et al. 2015). Only estimable combinations of them can be solved for. By constraining a minimum set of parameters as an \({\mathcal {S}}\)-basis, a full-rank version of the observation equations is formed, delivering such estimable parameters. Different choices of \({\mathcal {S}}\)-basis lead to different sets of estimable parameters. The adjusted observables (4) do, however, remain invariant for the choice of \({\mathcal {S}}\)-basis, as they are estimable by definition (Teunissen 1985). Let the estimable parameters be distinguished from their original counterparts by the \(\tilde{.}\)-symbol. In terms of the estimable parameter solutions, the corrections (4) can be expressed as (compare with 5)
Of the constituent solutions, forming the corrections (6), the clock and phase bias solutions \(d\hat{{\tilde{t}}}^{s}\) and \((\hat{\tilde{\delta }}^{s}-\hat{\tilde{z}}_{r}^{s})\) are needed for computing an ambiguity-resolved user position (Khodabandeh and Teunissen 2015). With the provision of the code bias solutions \(\hat{\tilde{d}}^{s}\) on the other hand, the user is able to make full use of the ‘multi-frequency’ code data. Also when the distance between the user receiver u and the reference receiver r is sufficiently short so that \(\iota _r^s\approx \iota _u^s\), the provision of the ionospheric solutions \(\hat{\tilde{\iota }}_r^{s}\) would strengthen the user model (3), i.e. by improving the user ambiguity resolution performance through constraining the unknowns \(\iota _{ru}^s=\iota _{u}^s-\iota _{r}^s\) to zero (Odijk 2002; Wielgosz et al. 2008). As long as the user positioning performance is concerned, the provision of the receiver-specific terms \(d\hat{{\tilde{t}}}_r\), \(\hat{\tilde{\delta }}_r\) and \(\hat{\tilde{d}}_r\) is, however, not required. This suggests that the ‘packaging’ of the corrections (6) can be further optimized in the sense of avoiding the transmission of unnecessary pieces of information to the user. Excluding the receiver-specific terms from (6) gives another class of single-station corrections
in which the phase and code corrections \(\hat{c}_{\phi }^{s}\) and \(\hat{c}_{p}^{s}\) are formed by the following constituent solutions
Since the corrections are provided to the user in terms of parameter (state) solutions, the above expression of the corrections is referred to as ‘state-space representation’ (SSR). Such representation is employed in PPP-RTK (Wubbena et al. 2005). SSR offers the flexibility to provide corrections of the time-stable and time-varying parameters separately, allowing the user to time-predict the individual corrections (8) on the basis of their distinct temporal behaviours. This generally leads to lower data transmission rates than those in OSR.
We conclude this section by comparing the single-station PPP-RTK user’s corrected observation equations with their RTK counterparts (3). After applying the PPP-RTK corrections (7), the user observation equations (2) can be shown to read (Khodabandeh and Teunissen 2016)
Compare (9) with (3). Both are identical in structure, with a difference, that the user unknown parameters in (3) are replaced by their PPP-RTK estimable versions indicated by the \(\tilde{.}\)-symbol in (9). In the next section we discuss how the estimability of the PPP-RTK user parameters is dictated by the choice of \({\mathcal {S}}\)-basis made by the provider, i.e. by the reference receiver r. A brief summary comparing single-station PPP-RTK with single-baseline RTK is given in Table 1.
3 Single-receiver estimable parameters
In the previous section, we discussed that single-station PPP-RTK and RTK differ from one another in their ways to provide corrections to the user (Table 1). Comparison of the between-satellite differences of the corrections (6) and (7) does, however, reveal the following identity
where the between-satellite differencing notation is defined as \((\cdot )^{1s}=(\cdot )^{s}-(\cdot )^{1}\). According to (10), the user positioning performances aided by the corrections (4) and (8) would be the same, if both the RTK and PPP-RTK providers employ a same estimation method for generating the corrections. From a data transmission perspective, however, PPP-RTK has the advantage over RTK in that estimable versions of the individual parameters are solved and disseminated to the user with the aim of reducing the underlying transmission rate. Characterizing such an advantage is not to degrade the RTK performance, but rather to highlight the flexibility offered by SSR. In this section we show how this can be realized through an application of \({\mathcal {S}}\)-system theory to a multi-epoch formulation of the reference receiver observation equations (2).
3.1 Multi-epoch formulation
Suppose that the reference receiver r collects a time series of the observations \(y_{r}^s=[\varDelta \phi _{r}^{s^T},\varDelta p_{r}^{s^T}]^T\), in (2), over k successive epochs \(i=1,\ldots ,k\). As the time-varying corrections (8) need to be recursively generated over time, the provider may process the stated time-series using minimum mean squared error filtering techniques such as the Kalman filter, see e.g. (Sanso 1980; Herring et al. 1990; Teunissen and Khodabandeh 2013). The standard assumptions on which the Kalman filter is based are placed on two types of models: the filter’s measurement model and dynamic model. We first consider the measurement model. Defining the notations
and \(e_o=[e^T,e^T]^T\), the filter’s measurement model is structured by the system of observation equations (2) as follows
for \(i=1,\ldots ,k\). The satellite and receiver biases form the vectors
while the satellite clock and ionospheric delays form the vectors
In the above equations and in the following, we use the epoch argument i to emphasize the time dependency of the quantities. In the case of the ambiguity parameter \(z_r^s\), however, the argument i is omitted. This is because the ambiguities are constant in time if no slips occur.
The measurement model (12) links the parameters to the GNSS observables. It does not describe the variability of the parameters in time though. It is the filter’s dynamic model that captures the temporal characteristics of the time-varying parameters. We consider the following class of dynamic models
for \(i=2,\ldots ,k\). The term \(\partial ^{^{q}}\!b_{r,i}^s\) denotes the qth-order time-derivative of \(b_{r,i}^s\), with \(\partial ^{^0}\!b_{r,i}^s=b_{r,i}^s\). \(\tau \) is the measurement sampling period, i.e. the reciprocal of the sampling rate. The highest degree of the polynomial equations (the highest power of \(\tau \)) is given by the integer h. The zero-mean process noises \(n_{a_r,i}\), \(n_{a^s,i}\) and \(n_{\partial ^{^{h}}\!b^s\!,i}\) are assumed to take zero samples (Teunissen 2001).
The first two equations, in (15), describe constant-state or random-walk processes (Odijk et al. 2015). They state that the biases \(a_{r,i}\) and \(a_{,i}^s\) (\(i=1,\ldots ,k\)) are time-constant with a degree of uncertainty. The uncertainty in the time-stability of \(a_{r,i}\) and \(a_{,i}^s\) is modelled through the choice of variances of their process noises \(n_{a_r,i}\) and \(n_{a^s,i}\). The smaller values the variances take, the more time-stable the biases are assumed. In the extreme case where the biases are assumed to behave constant in time, the stated variances should be set to zero.
The choice of random-walk processes for modelling the GNSS biases’ temporal behaviour is justified since such biases are reported to behave rather stable over time, see e.g. (Komjathy et al. 2005; Zhang et al. 2018). Random-walk processes may, however, not fully capture the clock and ionospheric parameters’ time-variation as they can rapidly change in time (Wang et al. 2017). As an extension of the random-walk process, we therefore use a polynomial dynamic model (the third set of the equations in 15) to describe the temporal behaviour of the satellite clock and ionospheric parameters \(b_{r,i}^s\) (\(i=1,\ldots ,k\)). As for the receiver clocks \(dt_r(i)\), no dynamic model is given in (15), meaning that the receiver clocks are assumed unlinked in time.
To appreciate the choice for the polynomial dynamic model, consider the last equation of (15) by setting \(q\mapsto h\), i.e.
It states that the hth-order time-derivative of \(b_{r,i}^s\) is time-constant with a degree of uncertainty. The polynomial model therefore reduces to a random-walk (constant-state) model for \(h=0\). When \(h\ge 1\), the remaining h equations of (15), for \(q=0,\ldots ,h-1\), follow by repeated integration of the constant derivative \(\partial ^{^h}\!b_{r,i-1}^s\) with respect to the sampling period \(\tau \) (Teunissen 2001). In the following we present full-rank systems of equations underlying the Kalman filter of the single receiver r for the two cases \(h=0\) and \(h=1\). Thus, the time-variations of the satellite clock and ionospheric parameters are modelled by ‘constant-state’ (\(h=0\)) and ‘constant-velocity’ (\(h=1\)) processes.
3.2 Constant-state setup (\(h=0\))
As shown by de Jonge (1998) and Odijk et al. (2015), the system of undifferenced GNSS observation equations is rank-deficient. This implies that the filter’s measurement and dynamic models (12) and (15) cannot directly serve as input for Kalman filtering. To remove the underlying rank-deficiencies, a set of the parameters (\({\mathcal {S}}\)-basis) need to get lumped into the remaining parameters (estimable parameters), recovering the system of equations to one of full rank. Thus, it is the full-rank versions of (12) and (15) that are applicable for Kalman filtering. For the constant-state case \(h=0\), such full-rank models can be expressed as
for \(i=2,\ldots ,k\). The constant-state setup (17) is divided into three sets of equations. The first set \({\mathrm{I}}\) is required to initialize the filter by estimating the parameters at the first epoch \(i=1\). Not all the parameters can, however, be solved for. The receiver clock \(dt_r(1)\) is lumped into the satellite clocks \(dt^s(i)\), whereas the receiver biases \(a_{r,1}\) are lumped into the satellite biases \(a_{,i}^s\). That is why their \(\tilde{.}\)-versions are constrained to zero to constitute the \({\mathcal {S}}\)-basis. Also not all the satellite biases \(a_{,1}^s\) can be solved for. The ionosphere-free (IF) and geometry-free (GF) combinations of the codes biases \(d_{,j}^s(1)\) on the first two frequencies (\(j=1,2\)) are fully absorbed by the satellite clocks \(dt^s(i)\) and the ionospheric delays \(\iota _r^s(i)\), respectively (Khodabandeh and Teunissen 2015). The parameters \(\tilde{d}_{,j}^s(1)\) (\(j=1,2\)), forming the vector \(F^T\!\tilde{a}_{,1}^s\), are therefore constrained to zero as well. That the ambiguities \(z_r^s\) are absent in (17) is due to the fact that they get lumped into the satellite phase biases \(\delta _{,j}^s(i)\). Thus, \(\tilde{z}_r^s=0\). As a result, what the filter delivers as solutions do not represent unbiased estimation of the original parameters, but that of their estimable counterparts. The interpretation of the single-station estimable parameters, together with the definition of the IF and GF combinations, is provided in Table 2.
Table 2 also presents the chosen \({\mathcal {S}}\)-basis formed by the parameters at the first epoch only. This is due to the presence of the second set \({\mathrm{I}\!\mathrm{I}}\), i.e. the dynamic models which provide information on the time-difference of the parameters. If such dynamic models are absent, the stated parameters at the remaining epochs (\(i=2,\ldots ,k\)) would need to be chosen as \({\mathcal {S}}\)-basis as well. From the second epoch onwards, all the parameters in the third set \({\mathrm{I}\!\mathrm{I}\!\mathrm{I}}\) can be recursively estimated in their estimable form. It should be remarked that the estimable parameters, in (17), inherit the temporal properties of the original parameters. This is because they are biased by time-constant offsets only, i.e. the \({\mathcal {S}}\)-basis parameters at the first epoch \(i=1\) (cf. Table 2).
If the variance of the process noises \(n_{b^s\!,i}\) is close to zero, the satellite clock and ionospheric parameters \(\tilde{b}_{r,i}^s\) are then considered ‘time-stable’ under the constant-state setup (17). This is an implicit assumption one makes in baseline positioning to eliminate the satellite clock and ionosphere delays (for short baselines) by double-differencing. In other words, the difference between the signal reception times, experienced by the reference and rover receivers, can exceed 10 milliseconds within which one ‘common’ satellite clock parameter is assumed for the two receivers (de Jonge 1998, pp. 27). Although such assumption seems to be realistic for short time intervals (several milliseconds), it does not necessarily hold for longer intervals (over several seconds). In such cases, the constant-state setup (17) needs to be extended so as to accommodate the rapid time-variation of the clock and ionospheric delays.
3.3 Constant-velocity setup (\(h=1\))
The constant-velocity setup follows by extending (17) so as to include the velocity parameters \(\partial dt^{s}(i)\) and \(\partial \iota _r^{s}(i)\) as additional unknowns. The corresponding full-rank version can be expressed as
for \(i=3,\ldots ,k\). Similar to (17), the constant-velocity setup (18) is divided into three sets of equations. Instead of one epoch, however, the filter’s initialization requires the observation equations of the first two epochs (set \({\mathrm{I}}\)). This is due to the inclusion of the velocity parameters \(\partial b_{r,i}^s\) (\(i=1,\ldots ,k\)) as unknowns present in the dynamic model (i.e. in set \({\mathrm{I}\!\mathrm{I}}\)). These unknowns are, however, not unbiasedly estimable. The extra rank-deficiency, caused by the presence of \(\partial dt^{s}(i)\), can be removed by extending the earlier \({\mathcal {S}}\)-basis to the receiver clock at the second epoch, i.e. \(dt_r(2)\). As a consequence, the clocks \(dt^{s}(i)\) and \(dt_r(i)\) as well as the rates \(\partial dt^{s}(i)\) are biased by the receiver velocity term
Table 3 presents changes in the earlier single-station estimable parameters by switching from the constant-state (17) to constant-velocity setup (18). As shown, the ionospheric rates \(\partial \tilde{\iota }_r^{s}(i)\) (\(i=1,\ldots ,k\)) can be unbiasedly determined under the single-station model (18). From the third epoch onwards, all the parameters in the third set \({\mathrm{I}\!\mathrm{I}\!\mathrm{I}}\) are recursively estimated in their estimable form.
With the class of dynamic models (15) provided, one can formulate further full-rank single-station setups, similar to (17) and (18), to account for the higher-order time-derivatives of the satellite clock and ionospheric parameters \(b_{r,i}^s\) like, e.g. the acceleration parameters \(\partial ^{^{2}}\!b_{r,i}^s\). This involves an extension of the earlier \({\mathcal {S}}\)-basis to the receiver clocks at further epochs, i.e. epoch \(i=3\) and beyond. For a detailed discussion on \({\mathcal {S}}\)-system theory applied to GNSS observation equations, see (Odijk et al. 2015).
3.4 User estimable parameters
Given the single-station setup (17) or (18), the PPP-RTK provider, i.e. the reference receiver r, is in a position to recursively generate the corrections (8) over time. As stated previously, the distinct temporal properties of the original parameters \(a_{,i}^s\) and \(b_{r,i}^s\) carry over to their estimable versions \(\tilde{a}_{,i}^s\) and \(\tilde{b}_{r,i}^s\). Thus, the transmission rate of the individual corrections (8) can be tuned on the basis of the time-variations of the original parameters. While the time-stability of the bias parameters \(a_{,i}^s\) can hold within several minutes up to few hours, solutions of the time-varying parameters \(b_{r,i}^s\) often need to be provided with the transmission periods not higher than 10–30 seconds (Wang et al. 2017).
Let us now revisit the scenario sketched in Fig. 1 in which the bias solutions \(\hat{\tilde{a}}_{,k_o}^s\) are provided to the user, \(k_o\) epochs after the filter’s initialization of the PPP-RTK provider. Due to a low transmission rate taken for \(\hat{\tilde{a}}_{,k_o}^s\), the updated bias corrections are not provided until \(k_o\) epochs afterwards. During the interval \([k_o,~2k_o]\), the solutions of the time-varying parameters (i.e. \(\hat{\tilde{b}}_{r,k}^s\)) need to be provided, say, k epochs after the filter’s initialization (\(k_o\le k< 2k_o\)). Collecting GNSS measurements at epoch \(l\ge k\), the user is therefore required to time-predict the PPP-RTK corrections using the individual corrections
Thanks to the dynamic models employed by the single-station provider (i.e. set \({\mathrm{I}\!\mathrm{I}}\) in 18), the time-prediction of the corrections at epoch \(i=l\) follows as
Here and in the following we refer to the time difference \((l-k)\tau \) as latency. Thus upon employing the constant-velocity setup (18) by the provider, the user is able to take care of the latency by time-predicting the parameters \(\hat{\tilde{a}}_{,l}^s\) and \(\tilde{b}_{r,l}^s\). Provision of the OSR corrections’ variation rate, for bridging the latency, can also be made, see e.g. (Vollath et al. 2012). This advantage comes at the expense that the velocity terms \(\partial \hat{\tilde{b}}_{r,k}^s\) need to be provided as extra corrections, i.e. the number of time-varying corrections becomes doubled. In view of the reduction of the corrections’ transmission rates, doubling the number of time-varying corrections will only be plausible if indeed it reduces the transmission rate of \(\hat{\tilde{b}}_{r,k}^s\)/\(\partial \hat{\tilde{b}}_{r,k}^s\) by a factor of more than two.
Making use of the notations \(\hat{c}^s=[\hat{c}_{\phi }^{s^T},\hat{c}_{p}^{s^T}]^T\) and (11), together with \(\hat{\tilde{z}}_{r}^{s}=0\), the corrections in combined form (7) read then
By applying the predicted corrections (21) to the user observations \(y_{u}^s=[\varDelta \phi _{u}^{s^T},\varDelta p_{u}^{s^T}]^T\), the user’s corrected observation equations (9) take the following form
for \(s=2,\ldots ,m\), where \(\varDelta \tilde{\rho }_{u}^{s}={g_u^{s}}^T\varDelta \tilde{x}_u\) and \(\mu _o=[-\mu ^T,\mu ^T]^T\) (cf. 11). The vector of user estimable biases is given as \(\tilde{a}_{u}=[\tilde{\delta }_{u}^T,\tilde{d}_{u}^T]^T\).
If the user collects measurements from at least \((\nu +1)\) satellites (i.e. \(m\ge \nu +1\)), the system of observation equations (22) is full-rank (Khodabandeh and Teunissen 2015). Due to the corrections’ dependency on the provider’s \({\mathcal {S}}\)-basis (cf. Table 1), the estimability of the user parameters will be driven by the \({\mathcal {S}}\)-basis. Their interpretation is presented in Table 4. As shown, the user position/ZTD increments \(\varDelta x_u\) are unbiasedly estimable. Importantly, the estimable ambiguities \(\tilde{z}_{ru}^{s}\) are straightforward double-differenced ambiguities and thus integer-valued. The ambiguities of the reference satellite \(s=1\), i.e. \(z_{ru}^1\), are lumped into the estimable receiver phase biases \(\tilde{\delta }_{u}\). That is why they are absent in the first equation of (22).
Now consider the case in which the distance between the user u and provider r is short enough to assume \(\iota _u^s=\iota _r^s\). Under this assumption, the estimable ionospheric delays \(\tilde{\iota }_{ru}^{s}\) reduce to (Table 4)
since \(\iota _{ru}^{s}(l)=0\). Substitution of (23) into (22) gives
for \(s=2,\ldots \!,m\). Thus, the m unknowns \(\tilde{\iota }_{ru}^{s}\) (\(s=1,\ldots ,m\)) are replaced by only one unknown \(d_{ru, \scriptscriptstyle G\!F}\). This increase in the user model’s redundancy raises the potential for instantaneous ambiguity resolution, where the integer ambiguities \(\tilde{z}_{ru}^{s}\) can be resolved by only a single epoch of user data (Odijk 2002). In the following we therefore focus our attention on the ionosphere-corrected model (24) to study the role played by the latency \((l-k)\tau \), and thus by the prediction error of (20), in the user single-epoch ambiguity resolution performance.
4 User ambiguity resolution performance
In this section, we present analytical expressions for the PPP-RTK user ambiguity variance matrix and analyse the impact the correction latency \((l-k)\tau \) has on the single-receiver ambiguity resolution performance. As here the focus is to be restricted to the double-differenced (DD) ambiguities \(\tilde{z}_{ru}^{s}\), we will be working with the between-satellite single-differenced version of (24). Using the notation \({\hat{y}}_{r}^s=[\varDelta \hat{\phi }_{r}^{s^T},\varDelta \hat{p}_{r}^{s^T}]^T\) and the identity (10), such observation equations follow by subtracting the first set from the second set of the equations (24) and substituting \({\hat{y}}_{r}^{1s}=\hat{c}^{1s}\),
With the above DD-like user-corrected observation equations, the user can obtain ambiguity-resolved solutions for the position \(\varDelta \tilde{x}_u\), providing that ‘successful’ integer ambiguity resolution is applied. The decision whether or not ambiguity resolution is deemed successful is determined by the probability of correct integer estimation, the so-called ambiguity ‘success-rate’ (Teunissen 1999b). The success-rate is driven by the ambiguity variance matrix, which on its turn is governed by the precision of both the user data \(y_{u}^{1s}\) and time-predicted provider data \({\hat{y}}_{r}^{1s}\) in (25). We therefore first discuss the redundancy of the constant-velocity model (18) and address how it contributes to the precision of \({\hat{y}}_{r}^{1s}\).
4.1 Redundancy of the constant-velocity setup
Since the precision of the least-squares adjusted observations \({\hat{y}}_{r}^{1s}\) depends on the condition equations of the underlying model (Teunissen 2000a), it is important to identify them for the case of the constant-velocity model (18). The number of condition equations is equal to the model’s redundancy. In Table 5 we outline how the redundancy is built up from the conditions imposed by the dynamic model in (18), i.e. the second set \({\mathrm{I}\!\mathrm{I}}\). To distinguish those conditions that contribute to the DD ambiguity solutions, the condition equations are divided into two parts: the satellite-differenced and satellite-averaged types (Khodabandeh and Teunissen 2018). Let us first consider the condition equations of the satellite-differenced type indicated by the superscript \((.)^{1s}\). These conditions are built up from the between-satellite single-differenced version of (18) in which the receiver-specific parameters \(d{\tilde{t}}_r(i)\) and \(\tilde{a}_{r,i}\) are absent. At the first epoch \(i=1\), only \(2f-2\) (or \(2[f-1]\)) out of 2f biases \(\tilde{a}_{,1}^{1s}\) are estimable per satellite pair 1s. Thanks to the constant-state dynamic model taken for the biases \(\tilde{a}_{,i}^{1s}\) (\(i=1,\ldots ,k\)), these \(2(f-1)\) estimable biases must behave constant over k epochs, thereby giving \(2(f-1)\) times \((k-1)\) condition equations per satellite pair 1s. Likewise, because of the constant-velocity dynamic model taken for the clock/ionospheric parameters \(\tilde{b}_{r,i}^{1s}\), the difference \(\tilde{b}_{r,i}^{1s}-\tilde{b}_{r,i-1}^{1s}\) (\(i=2,\ldots ,k\)) must behave constant over \(k-1\) epochs, giving 2 times \((k-2)\) condition equations per satellite pair 1s. Thus, the total number of satellite-differenced condition equations is equal to \(2[f(k-1)-1](m-1)\). It is these equations that contribute to the user DD ambiguity solutions, as their satellite-averaged counterparts are uncorrelated with the DD measurements (Khodabandeh and Teunissen 2018).
We now turn our attention to the satellite-averaged conditions indicated by the superscript \((.)^{\bar{s}}\) in Table 5. Due the presence of the estimable receiver clock \(d{\tilde{t}}_r(i)\) from the third epoch \(i=3\) onwards, the constant-state assumption on \(\tilde{a}_{,i}^{\bar{s}}\) only creates \(2(f-1)\) condition equations over the first two epochs \(i=1,2\). For the same reason, the difference \(\tilde{b}_{r,i}^{\bar{s}}-\tilde{b}_{r,i-1}^{\bar{s}}\) does not create condition equations. Instead, that the estimable receiver clocks \(d{\tilde{t}}_r(i)\) are assumed frequency-common among the 2f observables \(y_{r,i}^{\bar{s}}\), in (18), gives \((2f-1)\) condition equations per epoch i (\(i=3,\ldots ,k\)), increasing the model’s redundancy by \((2f-1)(k-2)\).
In order to evaluate the user ambiguity variance matrix, in the following the stochastic model of the between-satellite single-differenced measurements \(y_{r,i}=[y_{r,i}^{{12}^T},\ldots ,y_{r,i}^{{1m}^T}]^T\) is assumed given as
in which \(C=D_{m}^{T}C_{S}D_{m}\), with the cofactor matrix \(C_{S}\!=\!\mathrm{diag}(w_1^{-1}\!,\ldots ,\!w_{m}^{-1})\) whose positive elements \(w_{s}\) are the satellite elevation-dependent weights. The \(m\times (m-1)\) matrix \(D_m\) forms between-satellite differences (Teunissen 1997). The \(f\times f\) positive-definite matrices \(Q_\phi \) and \(Q_p\) are the cofactor matrices of the phase and code observable types, respectively. The dispersion operator is denoted as \(\mathsf {D}(.)\), while the symbol \(\otimes \) denotes the Kronecker product (Henderson et al. 1983).
In case of the dynamic model, the satellite biases \(\tilde{a}_{,i}^{1s}\) are assumed constant in time, and thus, the process noises \(n_{a^s,i}\), in (18), are set to be identically zero. If the zero-mean acceleration of the clock/ionospheric parameters \(\tilde{b}_{r,i}^{1s}\), i.e. \(\partial ^{^{2}}\!b_{r,i}^{1s}\), are normally distributed, the (co)variance matrices of the process noises \(n_{b,i}=[n_{b^{12},i}^T,\ldots ,n_{b^{1\!m},i}^T]^T\) and \(n_{\partial b,i}= [n_{\partial b^{12}\!,i}^T,\!\ldots \!,n_{\partial b^{1\!m}\!,i}^T]^T\), linking the parameters from epoch \(j\le i\) to epoch i, read then (Teunissen 2001)
where
The spectral density of the clock and ionospheric acceleration parameters \(\partial ^{^{2}}\!b_{r,i}^{s}\) are given by \(\sigma _{\partial ^{2}\!\rho }^2\) and \(\sigma _{\partial ^{2}\!\iota }^2\), respectively. In case of the satellite clocks, the spectral density \(\sigma _{\partial ^{2}\!\rho }^2\) may vary between 1 and 7 \([\mathrm{mm}^2/\mathrm{sec}^3]\), whereas the ionospheric spectral density \(\sigma _{\partial ^{2}\!\iota }^2\) may vary from fractions of \([\mathrm{mm}^2/\mathrm{sec}^3]\) to several \([\mathrm{mm}^2/\mathrm{sec}^3]\), depending on the location and time of the receiver data collection (Wang et al. 2018b).
4.2 Two single-epoch building blocks
The model’s redundancy of the single-station setup (18), outlined in Table 5, enables the PPP-RTK provider to compute the adjusted observations \({\hat{y}}_{r}^{1s}\). Would one discard this redundancy, then the combined correction \(\hat{c}^{1s} = {\hat{y}}_{r}^{1s}\) reduces to the single-epoch measurements \({\hat{y}}_{r}^{1s} = y_{r}^{1s}\). Substitution into (25) gives
which is the standard DD model. This indeed shows that single-station PPP-RTK encompasses single-epoch RTK as a special case. Let the minimum-variance least-squares solution for the vector \(z= [\tilde{z}_{ru}^{{12}^T},\ldots ,\tilde{z}_{ru}^{{1m}^T}]^T\), under the RTK model (29), be given by \({\hat{z}}^\mathrm{RTK}\). We will show that its variance matrix \(Q_{{\hat{z}}{\hat{z}}}^\mathrm{RTK}\) can serve as a building block of the user ambiguity variance matrix that is aided by the multi-epoch corrections \(\hat{c}^{1s}\).
To form the second building block, for argument-sake suppose that only the single-epoch bias solutions \(\hat{\tilde{a}}^{1s}\) are provided to the user. Due to the absence of the clock/ionospheric solutions \(\hat{\tilde{b}}_{r}^{1s}\), the corrections (21) are changed to the term \(\hat{c}^{1s}\!\mapsto \! (y_{r}^{1s}\!-\!B\hat{\tilde{b}}_{r}^{1s})\). Replacing \(y_{r}^{1s}\) by this term in (29) gives
with the biased range increment \(\varDelta \varrho _{u}^{1s}=\varDelta \tilde{\rho }_{u}^{1s}-d{\tilde{t}}^{1s}\). The model (30) is referred to as geometry-free (GF) since the information about the relative receiver-satellite geometry \(\varDelta \tilde{\rho }_{u}^{1s}=g_u^{1s}\varDelta \tilde{x}_u\) is hidden by the presence of the unknown satellite clocks \(d{\tilde{t}}^{1s}\), thereby being disregarded (Teunissen 1997). It is not the position/ZTD increment \(\varDelta \tilde{x}_u\) (of size \(\nu \)) that are involved as unknowns in the model, but instead the biased range increment \(\varDelta \varrho _{u}^{1s}\) (of size \(m-1\)). Note also that the ionospheric delays \(\tilde{\iota }_r^{1s}\) are included as extra unknowns. Although the GF model (30) has more unknowns than the RTK model (29) (thus being a weaker model), it is full-rank and allows one to compute the corresponding minimum-variance least-squares ambiguity solution \({\hat{z}}^\mathrm{GF}\). The variance matrix \(Q_{{\hat{z}}{\hat{z}}}^\mathrm{GF}\) will be shown to serve as the second building block. The precision of the solutions \({\hat{z}}^\mathrm{GF}\) and \({\hat{z}}^\mathrm{RTK}\) is characterized below.
Lemma 1
(RTK and GF ambiguity variance matrices) Let the variance matrix of the user data \(y_{u}\) be equal to that of the reference receiver \(y_{r}\) as given in (26). Then, the variance matrices of the ambiguity solutions, under the RTK and GF models (29) and (30), are given as
with the projector \(\mathcal {P}_{G}=G(G^TC^{-1}G)^{-1}G^TC^{-1}\), where \(G=[g_{u}^{12},\ldots ,g_{u}^{1m}]^T\) and \(\mu _e=\mu -\frac{\sigma _{\hat{\rho }\hat{\iota }}}{\sigma _{\hat{\iota }}^2}e\). The definition of the variance-type scalars \(\sigma _{\hat{\rho }|\iota }^2\) and \(\sigma _{\hat{\iota }}^2\) is provided in Table 6.
Proof
The single-epoch ambiguity variance matrices (31) follow from the results of (Teunissen 1997). \(\square \)
The factor 2, in the expressions (31), highlights the dependency of the ambiguity precision on the variance of the between-receiver single-differenced data \(y_{u}^{s}\!-\!y_{r}^{s}\) that is twice its undifferenced version. The expressions have been decomposed such as to clearly show the contribution of the phase and code precision. Next to the phase cofactor \(Q_{\phi }\) which is common to both the RTK and GF variance matrices, the ambiguity precision is also governed by the rather poor precision of the code data through the code-driven variances \(\sigma _{\hat{\rho }|\iota }^2\) and \(\sigma _{\hat{\iota }}^2\) (cf. Table 6). While the RTK variance matrix is only affected by the range variance \(\sigma _{\hat{\rho }|\iota }^2\), its GF counterpart is affected by both the range and ionospheric variances \(\sigma _{\hat{\rho }|\iota }^2\) and \(\sigma _{\hat{\iota }}^2\). This is due to the inclusion of the extra unknowns \(\tilde{\iota }_r^{1s}\) present in (30). In contrast to its GF counterpart, however, the RTK variance matrix depends on the \((m-1)\times \nu \) receiver-satellite geometry matrix G through the projector \(\mathcal {P}_G\). When a minimum number of satellites (i.e. \(m=1+\nu \)) are tracked by the user, G becomes then a square matrix. For square and invertible G matrices, the projector reduces to an identity matrix \(\mathcal {P}_G=I_{m-1}\) (thus \(\mathcal {P}_G C=C\)). This is the weakest case for the single-epoch short-baseline RTK model (29) in terms of ambiguity resolution capability (Teunissen 1997). The strongest case, known as geometry-fixed, follows when \(\mathcal {P}_G=0\), i.e. when the position/ZTD increment \(\varDelta \tilde{x}\) is absent in the model (29). In the geometry-fixed case, the RTK variance matrix is fully governed by the very precise phase data.
To get an indication of the ambiguity resolution strength of the single-epoch models (29) and (30), one can use the Ambiguity Dilution Of Precision (ADOP). The ADOP is defined as (Teunissen 1997)
where |.| denotes the determinant. From the ADOP, one can infer an upper bound for the bootstrapped ambiguity success-rate (Teunissen 2000b). The smaller the ADOP, the higher the upper bound of the success-rate becomes. For ADOP smaller than 0.14 cycles, the stated upper bound is always higher than 99%. The following lemma provides closed-form ADOP expressions for the variance matrices (31).
Lemma 2
(Single-epoch RTK and GF ADOPs) The ADOPs of the RTK and GF models (29) and (30) are, respectively, given as
with \(\bar{\sigma }_\phi =\sqrt{|Q_\phi |}^{\frac{1}{f}}\), \(\bar{\lambda }=\prod _{j=1}^f\lambda _j^{\frac{1}{f}}\), and
where \(\bar{w}=(1/m)\sum _{s=1}^m w_s\) and \(\bar{\bar{w}}=\prod _{s=1}^m w_s^{\frac{1}{m}}\) (cf. 26). The variance-type scalars \(\sigma _{\hat{\rho }|\iota }^2\), \(\sigma _{\check{\rho }|\iota }^2\), \(\sigma _{\hat{\iota }}^2\) and \(\sigma _{\check{\iota }}^2\) are defined in Table 6.
Proof
The ADOP expressions (33) follow from an application of the results of (Odijk and Teunissen 2008). The proof of the inequality (34) is in Appendix. \(\square \)
To get a better insight into the ADOP expressions (33), let us first make some approximation. The inequality (34) implies that the term \(\sqrt{2}\,w_o\) is bounded from above by 2 for equal elevation-dependent weights \(w^s=1\) (\(s=1,\ldots ,m\)), i.e. when \(\bar{\bar{w}}=\bar{w}=1\). The precision of current GNSS phase data is also about 1 per cent of their wavelength, i.e. \((\bar{\sigma }_\phi /\bar{\lambda })\approx 0.01\). This gives the approximation \(\sqrt{2}\,w_o(\bar{\sigma }_\phi /\bar{\lambda })\!\approx \! 0.02\), showing that the ADOPs (33) would have been indeed very small if the variance ratios \((\sigma _{\hat{\rho }|\iota }/\sigma _{\check{\rho }|\iota })\) and \((\sigma _{\hat{\iota }}/\sigma _{\check{\iota }})\) were absent. As shown in Table 6, the magnitude of these ratios is about 100, since the code data are almost two orders of magnitude less precise than their phase counterparts. We therefore arrive at the following rather crude approximation
Accordingly, the GF ADOP is roughly about \(0.02\times 100=2\) and \(0.02\times 100^{\frac{2}{3}}\approx 0.4\) cycles for the dual- (\(f=2\)) and triple-frequency (\(f=3\)) scenarios. These values are way larger than 0.14 cycles, showing that instantaneous ambiguity resolution is not feasible under the GF model (30). Hence, although the GF model is in mixed-integer form, it is unlikely to successfully map the ambiguity solutions \({\hat{z}}^\mathrm{GF}\) to their correct integers. In case of the RTK model (29), however, the ambiguity resolution performance benefits from the number of satellites m (see the dependency of the RTK ADOP on m). When the user receiver tracks a minimum number of satellites (i.e. \(m-1=\nu \)), the dual-frequency RTK ADOP is about \(0.02\times 100^{\frac{1}{2}}\approx 0.2\) cycles which is 10 times smaller than its GF counterpart. With a modest increase in the number of satellites, say \(m-1=2\nu \), the RTK ADOP reduces to about \(0.02\times 100^{\frac{1}{4}}\approx 0.06\) cycles. With the single-epoch RTK model (29), the user can therefore successfully resolve ambiguities instantaneously.
While the ADOP is a useful indicator of the overall model’s ambiguity resolution strength, it does not provide information on the precision of individual ambiguities. Such information can be provided by the conditional standard deviations of the LAMBDA-decorrelated ambiguities whose geometric average is equal to the ADOP (Teunissen et al. 1997). In Fig. 2 we present spectra of the stated conditional standard deviations for Galileo dual- and triple-frequency cases where the number of visible satellites is \(m=8\), with three positioning parameters, i.e. \(\nu =3\). For both dual- and triple-frequency cases, all the RTK spectra are flat and well below 0.1 cycles (green triangles). This is, however, not the case with the GF spectra. The GF spectra (black dots) have a ‘staircase’ pattern. When looking from left to right of the figure, the first 7 (i.e. \(m-1\)) spectra represent the rather large standard deviations of the GF ambiguities on one frequency (several cycles). There is a large discontinuity from the first to second 7 spectra for both dual- and triple-frequency cases (0.2–0.3 cycles) and ditto from the second to third 7 spectra for triple-frequency case (0.05 cycles). This is because the second and third 7 spectra represent the relatively small standard deviations of the frequency-differenced GF ambiguities (e.g. widelane and extra widelane ambiguities). Although the large value of the triple-frequency GF ADOP (0.31 cycles) shows that instantaneous full ambiguity resolution is not possible, it is still possible to reliably resolve a subset of the ambiguities (Jonkman et al. 2000), e.g. the last 7 ambiguities having spectra of similar magnitude as that of their RTK versions.
4.3 Precision of the minimum-variance ambiguities
With the single-epoch RTK model (29), the user can thus perform, in principle, instantaneous ambiguity resolution. One may then be inclined to conclude that the need for establishing multi-epoch filter setups like (18) is redundant. One should, however, be remarked that the spectacular ambiguity resolution performance of single-epoch RTK, presented in Fig. 2, relies on a stringent assumption, namely that the single-epoch corrections are instantaneously (i.e. with no latency) assumed given to the user. In practice, the corrections are delivered with nonzero latency. A more realistic description of the user ambiguity precision is therefore required when the aiding PPP-RTK corrections are subject to time delay.
Theorem 1
(Precision of the minimum-variance ambiguities) Given the scenario sketched in Fig. 1, suppose that the single-station bias and clock/ionospheric corrections of epochs \(k_o\) and \(k\!\ge \!k_o\) are given to the user at epoch \(l\!\ge \! k\), respectively. The variance matrix of the minimum-variance least-squares ambiguity solutions, under the PPP-RTK model (24), can then be expressed as
where
with the projector \(\mathcal {P}_G^\bot = I_{m-1}-\mathcal {P}_{G}\). The \(2\times 2\) cofactor matrix \(Q_{{\hat{b}}_l}\) depends on the latency \((l-k)\tau \) and is given by
The \(2\times 2\) matrices \(([k+1]/k)Q_{{\hat{b}}{\hat{b}}}\) (Table 6) and \(Q_{\partial {\hat{b}}_{l}}\) (cf. Appendix) denote the time-averaged and time-differenced components of \(Q_{{\hat{b}}_l}\), respectively. The variance-type scalar \(\sigma _{\hat{\iota }_{l}}^2\) is the second diagonal entry of \(Q_{{\hat{b}}_l}\).
Proof
see Appendix. \(\square \)
The expression (36) holds for the whole class of polynomial dynamic models (15) underlying the single-station model of the PPP-RTK provider. To facilitate its interpretation, the user ambiguity matrix (36) has been written as the sum of three terms. The first term is a scaled version of the single-epoch RTK variance matrix \(Q_{{\hat{z}}{\hat{z}}}^\mathrm{RTK}\). The factor \((k_o+1)/(2k_o)\) indicates how the user ambiguity precision is improved when the multi-epoch bias solutions \(\hat{\tilde{a}}_{,k_o}^s\) are provided, \(k_o\) epochs after the filter’s initialization at the provider side. The factor is equal to 1 at the initialization (\(k_o=1\)) and decreases the more the number of epochs \(k_o\). In the extreme case when \(k_o\rightarrow \infty \), the factor reaches its minimum value 0.5. This shows, for sufficiently large number of epochs, that one can at most tackle the uncertainty of the corrections, thereby leaving the precision of the user data to solely govern user ambiguity resolution (Khodabandeh and Teunissen 2015). The second term is due to the difference between the epochs \(k_o\) and k at which the bias and clock/ionospheric solutions \(\hat{\tilde{a}}_{,k_o}^s\) and \(\hat{\tilde{b}}_{r,k}^s\) are provided. This term is driven by the difference \((Q_{{\hat{z}}{\hat{z}}}^\mathrm{GF}-Q_{{\hat{z}}{\hat{z}}}^\mathrm{RTK})\). If both the solutions \(\hat{\tilde{a}}_{,k_o}^s\) and \(\hat{\tilde{b}}_{r,k}^s\) are provided at the same epoch \(k_o=k\), the second term would then vanish. The third (last) term Q(l) depends on the latency \((l-k)\tau \) through the \(2\times 2\) cofactor matrix \(Q_{\partial {\hat{b}}_{l}}\). For the zero-latency case \(l=k\), the magnitude of the entries of \(Q_{\partial {\hat{b}}_{l}}\) is never larger than the phase-driven matrix \(Q_{\check{b}\check{b}}\) (Table 6). Therefore, the approximation \(Q_{\partial {\hat{b}}_{l}}\approx 0\) leads to \(Q(l)\approx 0\) for the zero-latency case \(l=k\). However, the magnitude of \(Q_{\partial {\hat{b}}_{l}}\) (and thus Q(l)) gets larger, the longer the time delay \((l-k)\tau \) becomes. A general expression for the cofactor \(Q_{\partial {\hat{b}}_{l}}\) follows from the multi-epoch recursion of the provider’s filter (see 61 in Appendix). For the constant-velocity model (18), an easy-to-compute upper bound of \(Q_{\partial {\hat{b}}_{l}}\) is presented below.
Lemma 3
(Upper bound of \(Q_{\partial {\hat{b}}_l}\)) Let the \(2\times 2\) cofactor matrix S, in (28), be approximated as a scaled version of the phase-driven matrix \(Q_{\check{b}\check{b}}\) (Table 6) by \(\tau ^3S\approx \gamma \, Q_{\check{b}\check{b}}\), where \(\gamma \) is a nonnegative scalar. Under the single-station constant-velocity setup (18), the cofactor matrix \(Q_{\partial {\hat{b}}_{l}}\) in (38) is then bounded from above as follows
with \({\scriptstyle \varDelta }=l-k\) and \(\bar{\gamma }=\gamma /6\).
Proof
see Appendix. \(\square \)
Figure 3 shows the factor \(\gamma _{_l}\) (dotted lines) as well as its upper bound as functions of the latency \((l-k)\) for three different values \(\gamma =0.5,~1\) and 2 (on a logarithmic scale). As shown, for the zero-latency case \(l=k\), the factor \(\gamma _{_l}\) is always smaller than 1, confirming that the entries of \(Q_{{\hat{b}}_l}\) are never larger than those of the phase-driven cofactor matrix \(Q_{\check{b}\check{b}}\) (Table 6). The factor \(\gamma _{_l}\) is also shown to be monotonically increasing with latency. The factor is, however, not so much dependent on the epoch k (see the zoom-in of the figure). The larger the scalar \(\gamma \), the larger the factor \(\gamma _{_l}\) becomes. The magnitude of \(\gamma \) shall be chosen on the basis of the temporal properties of the satellite clock/ionospheric delays. As stated previously, their spectral densities \(\sigma _{\partial ^{2}\!\rho }^2\) and \(\sigma _{\partial ^{2}\!\iota }^2\), in (28), are of the order of \([\mathrm{mm}^2/\mathrm{sec}^3]\), whereas the GNSS phase variance is of the order of \(\mathrm{mm}^2\). For the approximation \(\tau ^3S\approx \gamma \, Q_{\check{b}\check{b}}\), in the following we therefore set \(\gamma =1\) to analyse the impact of the latency on the user ambiguity resolution performance. The latency impact on the user ADOP is characterized below.
Factor \(\gamma _{_l}\) (dotted lines) as well as its upper bound (solid lines), in (39), as functions of the latency \((l-k)\) [epoch] for three different values \(\gamma =0.5,~1\) and 2. The grey lines, shown in the zoom-in, correspond to different epochs \(k=10,\ldots ,200\)
Lemma 4
(ADOP of the PPP-RTK model) The ADOP of the PPP-RTK model (24) can be expressed as
in which the latency factor \(\mathcal {L}\) is given by
where
The variance-type scalar \(\omega _{\hat{\iota }}^2\) is the second diagonal entry of \(\varOmega _{{\hat{b}}}\). The definition of the cofactor matrices \(Q_{{\hat{b}}{\hat{b}}}\) and \(Q_{\check{b}\check{b}}\) as well as their entries \(\sigma _{\hat{\iota }}^2\) and \(\sigma _{\check{\iota }}^2\) is provided in Table 6.
Proof
see Appendix. \(\square \)
PPP-RTK user ADOP (40) (solid lines) as well as the corresponding latency factor \(\mathcal {L}\) (41) (dashed lines) as functions of the latency \((l-k)\tau \) [second] for different numbers of satellites m: (top) Galileo dual-frequency case E1/E5a; (bottom) triple-frequency case E1/E5a/E5b. The values follow by setting \(k_o=k=100\)
The PPP-RTK user ADOP (40) has been expressed as a scaled version of the single-epoch RTK ADOP. The first term (i.e. the multi-epoch factor \(\sqrt{(k_o+1)/2k_o}\)) is never larger than 1, while the second term (i.e. the latency factor \(\mathcal {L}\)) is always larger than 1. The role of the multi-epoch factor is to show how much improvement in the user ambiguity resolution performance is brought by the provider’s multi-epoch setup, whereas that of the latency factor is to indicate the degradation of the performance due to the time delay of the corrections. The multi-epoch factor is bounded from below by \(1/\sqrt{2}\), thereby decreasing the user ADOP at most by a factor of about 1.4. Hence, the multi-epoch filter setup is not so much to improve the user ambiguity resolution performance, but rather to enable the user to time-predict the corrections.
Figure 4 shows the user ADOP and its latency factor \(\mathcal {L}\) as functions of the latency \((l-k)\tau \) for Galileo dual- and triple-frequency cases. For both cases, the latency factor is shown to increase as both the latency and the number of satellites m increase. Contrary to the latency factor, the ADOP itself decreases the larger the number of satellites m. This is because the decrease in the single-epoch RTK ADOP (33) is larger than the increase in the latency factor (41) as a function of m.
Spectra of conditional standard deviations of LAMBDA-decorrelated ambiguities for the PPP-RTK model (24) for different correction latencies of 0, 3, 10, 15 and 20 seconds: (top) dual-frequency case E1/E5a; (bottom) triple-frequency case E1/E5a/E5b; (left) number of satellites \(m=8\); (right) number of satellites \(m=13\). Their corresponding single-epoch RTK counterparts (green triangles) are shown as reference
According to the results shown in Fig. 4, the correction latency should not be higher than 3 seconds to have ADOP smaller than 0.15 cycles when only 5 satellites are tracked by the dual-frequency user, i.e. when \(f=2\) and \(m=5\). For the case where more satellites are tracked (e.g. \(m=10\)), this maximum allowable latency can, however, be increased to around 10 seconds to guarantee similar ADOP values. By switching from the dual- to triple-frequency case, the stated latency can be even higher than 20 seconds to ensure such small ADOP values. Therefore, both the number of satellites m and number of frequencies f can work in tandem to enable one to increase the maximum allowable latency of the PPP-RTK corrections as far as the user-ADOP performance is concerned.
As stated previously, the ADOP does not convey information about the precision of individual ambiguities though. To obtain further insights into the PPP-RTK user ambiguity capabilities, we therefore also present spectra of the conditional standard deviations of the LAMBDA-decorrelated ambiguities in Fig. 5. The results have been shown for the latencies of 0, 3, 10, 15 and 20 seconds. Their single-epoch RTK counterparts (green triangles) are also shown as reference. For the zero-latency cases, all transformed spectra are shown to be flat and at a slightly lower level than their single-epoch RTK counterparts. As the correction latency gets higher, the discrepancy between the spectra gets larger, thus inheriting ‘staircase’ patterns more similar to that of their single-epoch GF versions (cf. Fig. 2). Also note that an increase in either of the number of satellites m and number of frequencies f can push the ambiguity spectra into a lower level. Importantly, one should be remarked that the cause of having small ADOP values for the triple-frequency case, in Fig. 4, must mainly be attributed to the small conditional standard deviations of the ‘frequency-differenced’ ambiguities. As shown in the bottom of Fig. 5, these conditional standard deviations remain small even when the latency increases to 20 seconds, just like their single-epoch GF versions (Fig. 2). This can be understood as follows. The frequency-differenced ambiguities (e.g. of widelane and extra widelane types) are functions of the combined GNSS measurements that have a less dependency on the geometric range, clock and ionospheric delays. Thus, the precision of such combined measurements (and linear functions thereof) is not affected by much, if the latency of the clock/ionospheric corrections increases. In such cases, small ADOP values do not indicate that full ambiguity resolution is feasible. Instead, one can go for the partial ambiguity resolution option and yet obtain positioning solutions of reasonable precision (Wang et al. 2018a).
Single-epoch user positioning solutions using multi-epoch single-station PPP-RTK corrections with latency of 10 seconds (Perth, Australia, Trimble receiver CUTA, 28 October 2019, with 1-Hz sampling rate): (left) E1/E5a Galileo-only, (middle) E1/E5a/E5b Galileo-only, (right) E1/E5a Galileo-plus-GPS. At the bottom are shown the ADOP (magenta lines) and number of tracked satellites (black lines) over time. The instantaneous, empirical success-rates, from left to right, are 99.3%, 99.9% and 100%, respectively
Example 1
To provide numerical insight into the user positioning performance, 1-Hz Galileo and GPS data-sets on frequencies L1-E1, L5-E5a and E5b are utilized to generate single-station PPP-RTK corrections under the multi-epoch model (18). The corrections are then given to a user (several meters apart) with latency of 10 seconds. The ionospheric delays of the reference and user receivers are therefore almost identical, \(\iota _u=\iota _r\). Using (20), the user time-predicts the corrections and applies them to about 5,400 epochs of data. The corresponding single-epoch positioning solutions are shown in Fig. 6. When looking from left to right, it can be observed how an increase in the number of frequencies and satellites reduces the number of incorrectly fixed solutions (red dots). For the dual-frequency Galileo-plus-GPS case, there is no incorrectly fixed solutions despite the larger ADOPs than those of the triple-frequency Galileo-only case. This can be understood from the flat and rather small dual-frequency ambiguity spectra (\(m=13\)) for the latency of 10 seconds in Fig. 5. It is also important to remark—unlike single-epoch RTK—that the user positioning solutions do not necessarily experience two orders of magnitude precision improvement for non-zero latency cases. This is because the user ambiguity-resolved carrier phase data are affected by the uncertainty of the time-predicted reference data \({\hat{y}}_r^{1s}\) (cf. 25). \(\boxtimes \)
From the user ambiguity spectra depicted in Fig. 5, we have learned that the spectra’s signature is similar to that of their single-epoch RTK versions for the zero-latency case, but will approach the staircase pattern of their single-epoch GF versions (Fig. 2) as the correction latency increases. This is what one should expect since the higher the latency, the larger the amount of uncertainty involved in the time-prediction of the clock/ionospheric corrections \(\hat{\tilde{b}}_{r,l}^s\) becomes (cf. 27). In the cases where the correction latency is extremely high, no information about \(\hat{\tilde{b}}_{r,l}^s\) can be retrieved by the user, just like the single-epoch GF model (30). One would thus expect the PPP-RTK user ambiguity performance to be a compromise between its single-epoch RTK and GF counterparts. This notion is made precise in the following corollary.
Corollary
(RTK and GF building blocks) Let the scalar \(\epsilon \) be the phase-to-code variance ratio. With \(\tau ^3S\approx \gamma \, Q_{\check{b}\check{b}}\) and the crude approximation \(Q_{\check{b}\check{b}}\approx \dfrac{\epsilon }{1+\epsilon }Q_{{\hat{b}}{\hat{b}}}\), the user ambiguity variance matrix (36) and ADOP (40) take the following forms
and
where
with the scalar \(\epsilon _{_l} = (\frac{\epsilon }{1+\epsilon })\gamma _{_l}\) (cf. 39).
Proof
see Appendix. \(\square \)
Thus, the PPP-RTK user ambiguity variance matrix (43) can be interpreted as a weighted average of its single-epoch RTK and GF counterparts that is scaled by \((k_o+1)/(2k_o)\). The weights \((1-\alpha )\) and \(\alpha \) tell us how close the PPP-RTK user ambiguity precision would be to that of the single-epoch RTK and GF ambiguities. For the zero-latency case (i.e. when \(\epsilon _{_l}\approx 0\) and \(k=k_o\)), the weight \(\alpha \) becomes zero, upon which the PPP-RTK ambiguity precision is never worse than its single-epoch RTK version. As the correction latency (and thus \(\epsilon _{_l}\)) increases, the weight \(\alpha \) gets larger. In the limit, when \(\epsilon _{_l}\rightarrow \infty \) (cf. Fig. 3), the weight \(\alpha \) tends to 1, upon which the stated ambiguity precision becomes of the order of its single-epoch GF version.
From the comparison of the ADOP expressions (40) and (44) also follows that
The above approximation shows the dependency of the latency factor (41) on various contributing factors, i.e. the number of satellites m, number of frequencies f, phase-to-code variance ratio \(\epsilon \) and weight \(\alpha \). The latency factor is shown to be monotonically increasing with \((\alpha /\epsilon )\) and m, but decreasing with f (Fig. 4).
5 Uncertainty of the corrections ignored
So far we based our analysis on an implicit assumption that the PPP-RTK user takes the uncertainty of the singe-station corrections \(\hat{c}^{s}\), in (21), into account. This means that the error variance matrix of the time-predicted corrections (20) is required to be made available to the user. By applying the error propagation law to the corrected data \((y_{u}^s\!-\!\hat{c}^{s})\) in (24), the user would then be in a position to employ a properly weighted least-squares adjustment, i.e. BLUE-estimation, to compute the minimum-variance ambiguity solutions (Teunissen 2000a).
In practice, however, the information about the corrections’ error variance matrix is often not provided to the user (Odijk et al. 2014). After all, the purpose of using SSR (8) is to limit the amount of information to be transmitted to the user. In the absence of such information, the corrections \(\hat{c}^{s}\) may be considered precise enough to be treated as nonrandom. For this case, the user only takes the variance matrix of his own data \(y_{u}^s\) for the weighting of the corrected data \((y_{u}^s\!-\hat{c}^{s})\) in (24). The consequence of this practice is that the weighted least squares parameter solutions, including those of the ambiguities, may loose its optimality property of ‘minimum-variance’. This has the following two implications for user integer ambiguity resolution:
-
Implication 1: Less precise float ambiguities. Such least-squares ambiguity solutions are not the most precise ones, thereby having variances larger than those described by the variance matrix (36). The user therefore cannot achieve the highest ambiguity success-rate.
-
Implication 2: Suboptimal integer estimation. The inverse of the corresponding reduced normal matrix, computed by such least squares adjustment, does not represent the variance matrix of the ambiguity solutions. Given such incorrect variance matrix, the corresponding integer least-squares estimation of the ambiguities is not optimal in the sense that it fails to deliver the highest success-rate (Teunissen 1999a). This again results in a reduction in the user ambiguity success-rate.
In the following theorem we address the extent to which the variance matrix of such suboptimal ambiguity solutions differs from its minimum-variance version (36).
Theorem 2
(Precision of the suboptimal ambiguities) Given the scenario sketched in Figure 1, suppose that the user ignores the uncertainty of the single-station corrections. Under the PPP-RTK model (24), the variance matrix of such least-squares ambiguity solutions, say \(Q_{{\hat{z}}{\hat{z}}}^\mathrm{PRK}\), is then related to its minimum-variance counterpart \(Q_{{\hat{z}}{\hat{z}}}\) as follows
where
The \(2\times 2\) cofactor matrices \(Q_{{\hat{b}}_l}\) and \(Q_{\partial {\hat{b}}_{l}}\) are given in (38), with scalars \(\sigma _{\hat{\iota }_{l}}^2\) and \(\sigma _{\partial \hat{\iota }_{l}}^2\) being their second diagonal entries, respectively. Moreover, the inverse of the corresponding reduced least-squares normal matrix is given as
Proof
see Appendix. \(\square \)
The expressions (47) and (48), together with (38), show that the difference between the precision of the suboptimal and minimum-variance user ambiguity solutions is driven by
which on its turn is governed by the \(2\times 2\) cofactor matrix \(Q_{\partial {\hat{b}}_{l}}\). The smaller the entries of \(Q_{\partial {\hat{b}}_{l}}\), the smaller this precision difference becomes. As previously illustrated in Fig. 3, the entries of \(Q_{\partial {\hat{b}}_{l}}\) are small (phase precision-level) for the zero-latency case and monotonically increase with latency. For low-latency cases, one can thus conclude that ignoring the PPP-RTK corrections’ uncertainty does not degrade the precision of the user ambiguities by much. With reference to Fig. 3, the variance matrix \(Q_{{\hat{z}}{\hat{z}}}^\mathrm{PRK}\) can, however, unboundedly deviate from its optimal version \(Q_{{\hat{z}}{\hat{z}}}\) when the correction latency is high. In the context of ambiguity resolution performance, the ADOP may be employed to quantify the terms ‘low’ and ‘high’ latencies.
Lemma 5
(Suboptimal-to-optimal ADOP-ratio) Let \(\mathrm{ADOP}^{PRK}\) denote the ADOP of the PPP-RTK model (24) where the uncertainty of the single-station corrections is ignored. This ADOP is related to its optimal version, in (40), as follows
where
with \(\omega _{\partial \hat{\iota }}^2\) being the second diagonal entry of \(\varOmega _{\partial {\hat{b}}}\).
Proof
see Appendix. \(\square \)
PPP-RTK user ADOP-ratio (51) as a function of the latency \((l-k)\tau \) [second] for different numbers of satellites m: Galileo dual-frequency case E1/E5a (solid lines); triple-frequency case E1/E5a/E5b (dashed lines). The values follow by setting \(k_o=k=100\)
Figure 7 shows the ADOP-ratio (51) as a function of the latency \((l-k)\tau \) [second] for different numbers of satellites m and for both the dual- and triple-frequency cases. As expected, this ratio increases with latency. Note also that the ratio does not depend on the number of satellites m by much, but the number of frequencies f. Importantly, it can be observed that the ADOP-ratio is close to 1 (less than 1.15) for low-latency cases (below 10 seconds). With reference to Implication 1, the user may therefore treat the PPP-RTK corrections as nonrandom and yet his ADOP-performance remains comparable to that of the minimum-variance case. This indeed would have held true if 1) the latency is low, and 2) the ambiguity variance matrix (47) is taken as input to the integer least-squares estimation. The latter concerns Implication 2. With the method of integer least-squares ambiguity resolution, the user achieves the highest possible ambiguity success-rate, providing that the correct variance matrix of the float ambiguities serves as input (Teunissen 1999a). Would the user treat the PPP-RTK corrections as nonrandom (deterministic), the weighted least-squares then incorrectly delivers the positive definite matrix \(M_{{\hat{z}}{\hat{z}}}^\mathrm{PRK}\), given in (49), as the ambiguity variance matrix. Comparison of (47) and (49) clearly shows that \(M_{{\hat{z}}{\hat{z}}}^\mathrm{PRK}\ne Q_{{\hat{z}}{\hat{z}}}^\mathrm{PRK}\). Thus, matrix \(M_{{\hat{z}}{\hat{z}}}^\mathrm{PRK}\) does not provide a realistic precision description of the user ambiguity solutions.
Examples of the 2D integer least-squares pull-in regions (hexagons) for a single set of 50,000 float ambiguity solutions (dots). The pull-in regions are formed by the correct (top) and an incorrect (bottom) ambiguity variance matrices. The number of to-be-correctly-fixed ambiguities, in the corresponding pull-in region (blue hexagons), is maximized when the correct variance matrix serves as input. The unique 99% confidence ellipse of the solutions (in black) is depicted in both figures
To gain a better understanding of the role played by the variance matrix in driving the user ambiguity success-rate, consider the two-dimensional (2D) examples given in Fig. 8. For a single set of 50,000 float ambiguity solutions, we have plotted integer least-squares pull-in regions (hexagons) that are formed by the correct and an incorrect ambiguity variance matrices. For the correct variance matrix case (top panel of the figure), the corresponding pull-in region (blue hexagon) captures a maximum number of the to-be-correctly-fixed ambiguities (green dots). This is because the shape and orientation of the pull-in region is governed by the variance matrix which describes the dispersion of the solutions. For the incorrect variance matrix case, however, (bottom panel of the figure), the pull-in region fails to capture the maximum number of the solutions, since the input variance matrix does not give a correct description of the solutions’ dispersion. As a consequence, the success-rate drops from 99.5% to 96.9%, leading to more to-be-wrongly-fixed ambiguities (red dots).
Single-epoch user east-north scatterplots using multi-epoch single-station PPP-RTK corrections with latencies of 5 (top) and 10 (bottom) seconds (cf. Fig. 6). From left to right (and top to bottom), the first three figures refer to the ‘incorrect variance matrix’ case where the user ignores the corrections’ uncertainty, thus using the wrong ambiguity variance matrix (49). The fourth figure (bottom-right) refers to the ‘correct variance matrix’ case where the user still ignores the corrections’ uncertainty, but uses the correct ambiguity variance matrix (47): E1/E5a Galileo-only (top-left) and E1/E5a Galileo-plus-GPS (the three remaining figures). The instantaneous, empirical success-rates are 99.5% (top-left), 100% (top-right), 98.6% (bottom-left) and 100% (bottom-right)
Example 2
(Example 1 continued) Consider again the 1-Hz Galileo and GPS data-sets as discussed in Example 1. We now make an extra assumption that the corrections’ variance matrix is absent and the user ignores the corrections’ uncertainty by treating them as nonrandom. The user ambiguity-resolved positioning performance is therefore suboptimal in the sense that it may fail to maximize the number of correctly fixed solutions. The corresponding 5,400 single-epoch positioning solutions are shown in Fig. 9. For the 5-second latency cases (top panel of the figure), it can be observed that the Galileo-GPS integration does indeed increase the success-rate from 99.5% to 100%. Thus despite the usage of an incorrect ambiguity variance matrix, the user success-rate remains high when the number of tracked satellites is sufficiently large. For the 10-second latency cases, however (bottom panel of the figure), the ambiguity success-rate of the same Galileo-plus-GPS scenario drops to 98.6%. As shown at the bottom-right of the figure, this reduction in the ambiguity success-rate can be avoided if the user makes use of the correct ambiguity variance matrix (47) as the input of his method of integer least squares ambiguity resolution.
Single-epoch user east-north scatterplots of station UWA0 (cf. Fig. 9). The instantaneous, empirical success-rates are 97.9% (top-left), 100% (top-right), 98.6% (bottom-left) and 99.9% (bottom-right)
To highlight the role played by the distance between the user and the provider, in Table 7 we present values of the empirical ambiguity success-rate for two different inter-station distances. Station CUTB is just a few meters away from the provider, whereas the distance between station UWA0 and the provider is 8km. For both stations, the zero-latency cases lead to success-rates of 100%. However, when the corrections’ uncertainty is ignored for the nonzero-latency cases, the reduction in the success-rate is more pronounced for station UWA0. As with Fig. 9, Figure 10 shows user east-north scatterplots, but now for station UWA0. As shown, the reduction in the success-rate can be avoided if the correct ambiguity variance matrix is taken. \(\square \)
The results, shown in Figs. 9 and 10, together with Table 7, are in agreement with our analytical results in (47) and (51). The conclusion reads therefore that the PPP-RTK user ambiguity resolution performance does not differ too much from its minimum-variance counterpart if 1) the correction latency would be low, and 2) the precision description of the float ambiguities are correctly given as input to the user’s method of ambiguity resolution. The latter can, for instance, be handled through evaluating the corrections’ variance matrix by the user himself (Odijk et al. 2014).
6 Summary and conclusions
In this contribution we developed a framework to generate and time-predict multi-epoch PPP-RTK corrections from the GNSS data of a single station. Analytical expressions of the user ambiguity variance matrix were presented to study the role played by the correction latency, and thus by the corrections’ prediction error, in the user single-epoch ambiguity resolution performance. Our main findings are summarized as follows:
-
SSR linked to OSR. We started off by discussing the fine distinction between single-station PPP-RTK and RTK. The combined PPP-RTK corrections (7), with state-space representation (SSR), were shown to differ from their RTK versions (4), with observation-space representation (OSR), only in their \({\mathcal {S}}\)-basis dependent receiver-specific terms (Table 1). Their respective user positioning performances would therefore be the same, if both the single-station PPP-RTK and RTK providers employ a same estimation method for generating the corrections.
-
Single-receiver estimable parameters. The distinct temporal properties of the GNSS individual parameters suggest sending parameter solutions, each with its own transmission rate, to the user. This is the advantage PPP-RTK has over RTK in that estimable versions of the individual parameters are solved and disseminated to the user with the aim of reducing the data transmission rate. To realize this, we applied \({\mathcal {S}}\)-system theory to a class of polynomial dynamic models (15) and identified single-receiver estimable parameters of both the provider and the user (Tables 2 to 4). It was shown that the estimable parameters inherit the temporal properties of their original counterparts.
-
Marginal ambiguity precision improvement by ‘multi-epoch’ corrections. The multi-epoch condition equations that are imposed by the single-station constant-velocity setup (18) were characterized (Table 5). It is these equations that contribute to the precision of the PPP-RTK corrections and therefore to that of the user ambiguity solutions. The provider’s multi-epoch filter setup is, however, not so much to improve the user ambiguity resolution performance, but rather to enable the user to time-predict the corrections.
-
Building blocks of the user ambiguity variance matrix. The single-station PPP-RTK user ambiguity performance was shown to be a compromise between its single-epoch RTK and geometry-free (GF) building blocks. The former follows when the corrections are provided to the user instantaneously (i.e. with zero-latency), whereas the latter follows when no information on the satellite clock/ionospheric delays is given. With the provision of \(k_o\)-epoch bias corrections, the user ambiguity matrix \(Q_{{\hat{z}}{\hat{z}}}\) was therefore shown to approximate
$$\begin{aligned} Q_{{\hat{z}}{\hat{z}}} \approx \left\{ \begin{array}{ll} \frac{k_o\!+\!1}{2k_o} Q_{{\hat{z}}{\hat{z}}}^\mathrm{RTK}, &{} ~\text {zero-latency}\\ \frac{k_o\!+\!1}{2k_o} Q_{{\hat{z}}{\hat{z}}}^\mathrm{GF}, &{} ~\text {high-latency}\end{array}\right. \end{aligned}$$(53) -
Role of the number of satellites and number of frequencies. Supported by the numerical results, our ADOP-analysis showed that both the number of satellites and number of frequencies can work in tandem to enable one to increase the correction latency of the PPP-RTK corrections. To have ADOP values smaller than 0.15 cycles, the latency should be below 5 seconds for dual-frequency data of 5 satellites. The stated latency can, however, be increased to over 10 seconds when dual- or triple-frequency data of more than 10 satellites are tracked.
-
Suboptimal vs. minimum-variance ambiguity solutions. In practice, the corrections’ variance matrix is often not provided to the user. We thus quantified the extent to which the precision of the user ambiguity solutions differs from that of their minimum-variance versions, when the user ignores the corrections’ uncertainty, thereby treating them as nonrandom. Such suboptimal user ambiguity matrix \(Q_{{\hat{z}}{\hat{z}}}^\mathrm{PRK}\) was shown to approximate
$$\begin{aligned} Q_{{\hat{z}}{{\hat{z}}}}^\mathrm{PRK} \approx \left\{ \begin{array}{ll} Q_{{{\hat{z}}}{{\hat{z}}}}, &{} ~\text {low-latency} \\ \infty , &{} ~\text {high-latency}\end{array}\right. \end{aligned}$$(54)Thus for low-latency cases, the precision of the user ambiguities is not degraded by much. For high-latency cases, the variance matrix \(Q_{{{\hat{z}}}{{\hat{z}}}}^\mathrm{PRK}\) can, however, unboundedly deviate from its optimal version \(Q_{{{\hat{z}}}{{\hat{z}}}}\). In order to do proper justice to the user ambiguity resolution strength, the user is required to achieve the highest possible ambiguity success-rate. This can be realized if the precision description of the float ambiguities is correctly given as input to the user’s method of ambiguity resolution.
Data Availability
The GNSS data used in the paper are freely accessible via the website http://saegnss2.curtin.edu.au/ldc/.
References
Baarda W (1973) S-transformations and Criterion Matrices. Tech. rep., Netherlands Geodetic Commission, Publ. on Geodesy, New Series, Vol. 5(1), Delft
Banville S, Lachapelle G, Ghoddousi-Fard R, Gratton P (2019) Automated processing of low-cost GNSS receiver data. In: Proceedings of Institute of Navigation GNSS+ 2019 conference, Miami
Bertiger W, Desai SD, Haines B, Harvey N, Moore AW, Owen S, Weiss JP (2010) Single receiver phase ambiguity resolution with GPS data. J Geod 84(5):327–337
Collins P, Bisnath S, Lahaye F, Heroux P (2010) Undifferenced GPS ambiguity resolution using the decoupled clock model and ambiguity datum fixing. Navigation 57(2):123–135
Deambrogio L, Palestini C, Bastia F, Gabelli G, Corazza GE, Samson J (2010) Impact of high-end receivers in a peer-to-peer cooperative localization system. In: Ubiquitous Positioning Indoor Navigation and Location Based Service, IEEE, https://doi.org/10.1109/UPINLBS.2010.5654054
Ge M, Gendt G, Rothacher M, Shi C, Liu J (2008) Resolution of GPS carrier-phase ambiguities in precise point positioning (PPP) with daily observations. J Geod 82(7):389–399
Henderson HV, Pukelsheim F, Searle SR (1983) On the history of the kronecker product. Linear Mult Algebra 14(2):113–120
Heroux P, Kouba J (1995) GPS precise point positioning with a difference. Paper presented at Geomatics 95. Canada, Ottawa, pp 13–15
Herring TA, Davis JL, Shapiro II (1990) Geodesy by radio interferometry: The application of Kalman Filtering to the analysis of very long baseline interferometry data. J Geophys Res 95(8):12,561–12,581
de Jonge PJ (1998) A processing strategy for the application of the GPS in networks. PhD thesis, Delft University of Technology, Publication on Geodesy, 46, Netherlands Geodetic Commission, Delft
Jonkman NF, Teunissen PJG, Joosten P, Odijk D (2000) GNSS long baseline ambiguity resolution: impact of a third navigation frequency. In: IAG symp 121, Geodesy Beyond 2000, Springer, pp 349–354
Khodabandeh A, Teunissen PJG (2015) An analytical study of PPP-RTK corrections: precision, correlation and user-impact. J Geod 89(11):1109–1132
Khodabandeh A, Teunissen PJG (2016) PPP-RTK and inter-system biases: the ISB look-up table as a means to support multi-system PPP-RTK. J Geod 90(9):837–851
Khodabandeh A, Teunissen PJG (2018) On the impact of GNSS ambiguity resolution: geometry, ionosphere, time and biases. J Geod 92(6):637–658
Khodabandeh A, Teunissen PJG (2019) Distributed least-squares estimation applied to GNSS networks. Measurem Sci Technol, 30(4): 044,005
Koch KR (1999) Parameter estimation and hypothesis testing in linear models. Springer, Berlin
Komjathy A, Sparks L, Wilson BD, Mannucci AJ (2005) Automated daily processing of more than 1000 ground-based GPS receivers for studying intense ionospheric storms. Radio Sci. https://doi.org/10.1029/2005RS003279
Laurichesse D, Mercier F, Berthias J, Broca P, Cerri L, CNES F, (2009) Integer ambiguity resolution on undifferenced GPS phase measurements and its application to PPP and satellite precise orbit determination. Navig 56(2):135–149
Odijk D (2002) Fast precise GPS positioning in the presence of ionospheric delays. Ph.D. thesis, Delft University of Technology, Publication on Geodesy, 52, Netherlands, Geodetic Commission, Delft
Odijk D, Teunissen PJG (2008) ADOP in closed form for a hierarchy of multi-frequency single-baseline GNSS models. J Geod 82(8):473–492
Odijk D, Teunissen PJG, Khodabandeh A (2014) Single-frequency PPP-RTK: theory and experimental results. IAG Symp 139:167–173
Odijk D, Zhang B, Khodabandeh A, Odolinski R, Teunissen PJG (2015) On the estimability of parameters in undifferenced, uncombined GNSS network and PPP-RTK user models by means of S-system theory. J Geod 90(1):15–44
Sanso F (1980) The minimum mean square estimation error principle in physical geodesy (stochastic and non-stochastic interpretation). Boll Geod Sci Affi 39(2):112–129
Teunissen PJG (1985) Generalized inverses, adjustment, the datum problem and S-transformations. In: Grafarend EW, Sanso F (eds) Optimization and design of geodetic networks. Springer, Berlin, pp 11–55
Teunissen PJG (1997) A canonical theory for short GPS baselines. Part I: The baseline precision, Part II: The ambiguity precision and correlation, Part III: The geometry of the ambiguity search space, Part IV: Precision versus reliability. J Geod 71(6): 320–336
Teunissen PJG (1999a) The effect on GNSS ambiguity resolution of using an approximate weight matrix. Artif Satell 34(2):59–64
Teunissen PJG (1999b) An optimality property of the integer least-squares estimator. J Geod 73(11):587–593
Teunissen PJG (2000a) Adjustment theory: an introduction. Delft University Press, series on Mathematical Geodesy and Positioning
Teunissen PJG (2000b) ADOP based upper bounds for the bootstrapped and the least squares ambiguity success. Artif Satell 34(4):171–179
Teunissen PJG (2001) Dynamic data processing. Delft University Press, Delft
Teunissen PJG, Khodabandeh A (2013) BLUE, BLUP and the Kalman filter: some new results. J Geod 87(5):1–13
Teunissen PJG, de Jonge PJ, Tiberius CCJM (1997) The least-squares ambiguity decorrelation adjustment: its performance on short GPS baselines and short observation spans. J Geod 71(10):589–602
Teunissen PJG, Odijk D, Zhang B (2010) PPP-RTK: results of CORS network-based ppp with integer ambiguity resolution. J Aeronaut Astronaut Aviat 42(4):223–229
Vollath U, Talbot NC, Glocker M, Chen X, Leandro R (2012) GNSS Signal Processing with Delta Phase. Patent Pub. No.: WO2012151006 A1, Trimble Navigation Limited, November 2012
Wang K, Khodabandeh A, Teunissen PJG (2017) A study on predicting network corrections in PPP-RTK processing. Adv Space Res 60(7):1463–1477
Wang K, Khodabandeh A, Teunissen PJG (2018a) Five-frequency Galileo long-baseline ambiguity resolution with multipath mitigation. GPS Solut 22(3):75
Wang K, Khodabandeh A, Teunissen PJG, Nadarajah N (2018b) Satellite-clock modeling in single-frequency PPP-RTK processing. J Surv Eng 144(2):04018,003, 10.1061/(ASCE)SU.1943-5428.0000252
Wielgosz P, Krankowski A, Sieradzki R, Grejner-Brzezinska DA (2008) Application of predictive regional ionosphere model to medium range RTK positioning. Acta Geophysica 56(4):1147–1161
Wubbena G, Schmitz M, Bagg A (2005) PPP-RTK: Precise Point Positioning using state-space representation in RTK networks. In: Proceedings of ION GNSS, pp 13–16
Zhang B, Liu T, Yuan Y (2018) GPS receiver phase biases estimable in PPP-RTK networks: dynamic characterization and impact analysis. J Geod 92(6):659–674
Zumberge JF, Heflin MB, Jefferson DC, Watkins MM, Webb FH (1997) Precise point positioning for the efficient and robust analysis of GPS data from large networks. J Geophys Res 102:5005–5017
Acknowledgements
This work benefited from author’s discussions with Dr Safoora Zaminpardaz from RMIT University. The GNSS data are provided by the Curtin GNSS Research Centre. This support is gratefully acknowledged.
Author information
Authors and Affiliations
Contributions
The author contributed to the results and writing of the paper.
Corresponding author
Appendix
Appendix
Proof of the inequality (34) in Lemma 2. Substitution of the identities \(\sum _{s=1}^m w_s=m\bar{w}\) and \(\prod _{s=1}^m w_s = \bar{\bar{w}}^m\) into the first expression of (34) gives
The above inequality holds if the term \((\frac{\bar{w}}{\bar{\bar{w}}})^{\frac{1}{2(m-1)}}\) is a monotonically decreasing function of \(m\ge 2\) when \(\bar{w}\ge \bar{\bar{w}}\). This is always the case because of the ‘arithmetic-geometric-mean inequality’. This, together with the inequality \(\sqrt{2}\ge m^{\frac{1}{2(m-1)}}\), gives (34). \(\square \)
Proof of Theorem (1)
To prove (36), we apply the least-squares conditional adjustment to the user GF ambiguity solutions. Such solutions follow from the GF model (30) by replacing the reference data \(y_{r}^{1s}\) with its time-averaged version \(y_{r,\bar{k}_o}^{1s}\) (Khodabandeh and Teunissen 2015),
Application of the error propagation law to (56) gives
since \(\mathsf {D}(\phi _{r}(\bar{k}_o))=(1/k_o)\mathsf {D}(\phi _{u}(l))\) and \(\mathsf {D}({\hat{b}}_{u,\bar{k}_o})=(1/k_o)\mathsf {D}({\hat{b}}_{u,l})\). When the ‘time-prediction’ of the satellite clock/ionospheric delays \({\hat{b}}_{u,l}\) are provided as \({\hat{b}}_{r,l|k}\), the ambiguity solutions (56) can be adjusted through the resultant misclosure vector
with \(u_1=[1,0]^T\), \(u_2=[0,1]^T\), and \(G^{\bot }\) being an orthogonal complement basis matrix of G. The first set of misclosure vector (58) states that the ionospheric delays of the user and provider must have the same mean, while the second set states that clock solutions must have the same mean. The expression (36) follows then from the least-squares conditional formula (Teunissen 2000a)
and the equality
with \(Q_{tt}\) and \(Q_{{{\hat{z}}}_{ru,\bar{k}_o}^\mathrm{GF}\!,t}\) being the variance matrix of t and covariance matrix between \({{\hat{z}}}_{ru,\bar{k}_o}^\mathrm{GF}\) and t, respectively. \(\square \)
Proof of Lemma (3)
Since the proof is lengthy, we only give a sketch. Let the vectors \(n_{\partial ^{^{q}}\!b,i}\) \((q=2,\ldots ,h)\) be defined in the same way as the process noise vectors \(n_{b,i}\) and \(n_{\partial b,i}\) in (27). Also let us define the higher-dimension vectors \(\vartheta _{r,i}=[\partial b_{r,i}^T,\ldots ,\partial ^{^{h}}\!b_{r,i}^T]^T\) and \(\eta _{i}=[n_{\partial b,i}^T,\ldots ,n_{\partial ^{^{h}}\!b,i}^T]^T\). Then, the last expression of the dynamic model (15) takes the following multi-epoch structure
where use is of the multi-epoch notation \((\cdot )_{[k]}\!=\![(\cdot )_{1}^T,\ldots ,(\cdot )_{k}^T]^T\). The \(k\times (k-1)\) matrix \(D_k\) forms between-epoch differences. The entries of the coefficient matrices \(F_1\) and \(F_2\) follow from (15) and are specified once the highest degree of the polynomial equations h is chosen. Solving the last set of equations \(i=k,k-1,\ldots , k-p\) of (61) (with \(0\le p\le k-1\)) gives an upper bound for the unknown parameter \(\partial b_{r,k}\). The upper bound gets sharper the more the epochs p. For the constant-velocity model \(h=1\), the upper bound (39) is obtained based on the last two epochs k and \(k-1\). \(\square \)
Proof of Lemma (4)
Application of the determinant factorization rule (Koch 1999) to (36) gives
The second term in (62) can be simplified as
since \(B^TQ_\phi ^{-1}\!B = Q_{\check{b}\check{b}}^{-1}-\!Q_{{\hat{b}}{\hat{b}}}^{-1}\). Likewise, the third term in (62) can be simplified as
since \(\mu _e^T[Q_\phi \!+\!\sigma _{\hat{\rho }|\iota }^2ee^T]^{-1}\mu _e = \sigma _{\check{\iota }}^{-2}\!-\!\sigma _{\hat{\iota }}^{-2}\). Substitution of (63) and (64) into (62), together with (32), gives (40). \(\square \)
Proof of Corollary
With \(\tau ^3S\approx \gamma \, Q_{\check{b}\check{b}}\) and approximation \(Q_{\check{b}\check{b}}\approx \dfrac{\epsilon }{1+\epsilon }Q_{{\hat{b}}{\hat{b}}}\), the matrix Q(l) in (37) takes the following form
The result (43) follows then by substituting (65) into (36). To prove (44), consider the approximation \(\varOmega _{{\hat{b}}}\approx \alpha \, Q_{{\hat{b}}{\hat{b}}}\) giving
Combining the above result with (62), along with (63) and (64), yields (44). \(\square \)
Proof of Theorem (2)
The proof goes along the same lines as the proof of Theorem (1). The suboptimal ambiguity solutions are obtained from an application of the weighted least-squares conditional adjustment (Teunissen 2000a) to the user GF ambiguity solutions (56). The variance matrix (47) follows then by an application of the error propagation law to the weighted adjusted ambiguity solutions. The incorrect variance matrix (49) is reported by the least-squares adjustment, when the variance matrix of \({\hat{y}}_{r}^{1s}\) in (25) is ignored (i.e. set to zero). The model becomes then identical to the RTK model (29), with a difference, that the data variance matrix is wrongly assumed to be reduced by half. \(\square \)
Proof of Lemma (5)
As with (62), an application of the determinant factorization rule to (47) gives
Dividing (67) by (62), together with (32), gives (51). \(\square \)
Rights and permissions
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.
About this article
Cite this article
Khodabandeh, A. Single-station PPP-RTK: correction latency and ambiguity resolution performance. J Geod 95, 42 (2021). https://doi.org/10.1007/s00190-021-01490-z
Received:
Accepted:
Published:
DOI: https://doi.org/10.1007/s00190-021-01490-z