Improved precise positioning with BDS-3 quad-frequency signals

Abstract

The establishment of the BeiDou global navigation satellite system (BDS-3) has been completed, and the current constellation can independently provide positioning service globally. BDS-3 satellites provide quad-frequency signals, which can benefit the ambiguity resolution (AR) and high-precision positioning. This paper discusses the benefits of quad-frequency observations, including the precision gain of multi-frequency high-precision positioning and the sophisticated choice of extra-wide-lane (EWL) or wide-lane (WL) combinations for instantaneous EWL/WL AR. Additionally, the performance of EWL real-time kinematic (ERTK) positioning that only uses EWL/WL combinations is investigated. The results indicate that the horizontal positioning errors of ERTK positioning using ionosphere-free (IF) EWL observations are approximately 0.5 m for the baseline of 27 km and 1 m for the baseline of 300 km. Furthermore, the positioning errors are reduced to the centimetre level if the IF EWL observations are smoothed by narrow-lane observations for a short period.

Introduction

The BeiDou Navigation Satellite System (BDS) is one of the four global satellite navigation systems (GNSSs). The establishment of the space segment of BDS is divided into three steps (Yang et al. 2019). The first constellation of BDS was constructed in 2003, namely BeiDou Navigation Demonstration System. The second constellation, which is called BeiDou Regional Navigation Satellite System (BDS-2), was completed and put into service in July 2012. The BDS-2, which currently consists of five geostationary earth orbit (GEO) satellites, seven Inclined GeoSynchronous Orbit (IGSO) satellites, and three Medium Earth Orbit (MEO) satellites, can provide Positioning, Navigation, and Timing (PNT) services in the Asia Pacific Region. The third constellation, which is called BeiDou Global Navigation Satellite System (BDS-3), currently consists of three GEO satellites, three IGSO satellites, and 24 MEO satellites; one of the GEO satellites is in testing (Zhang et al. 2019b). BDS-3 can provide PNT service globally with positioning, timing, and velocity accuracies of 10 m, 20 ns, and 0.2 m/s, respectively (CSNO 2018). Additional BDS satellite information can be found at http://www.csno-tarc.cn/system/constellation.

BDS-2 satellites broadcast triple-frequency signals, i.e. B1I, B2I, and B3I, which are centred at 1561.098, 1207.140, and 1268.520 MHz, respectively. BDS-3 satellites broadcast five-frequency signals, i.e. B1I, B3I, B1C, B2a, and B2b; B1C, B2a, and B2b are centred at 1575.420, 1176.450, and 1207.140 MHz, respectively. However, among these five signals of BDS-3 satellites, only B1I, B1C, B2a, and B3I are officially available for public service as of June 21, 2020. This paper examines the benefits of these quad-frequency signals with regard to Ambiguity Resolution (AR) and positioning. For the simplification of the following analysis, the four public available signals are ordered by their frequencies, as shown in Table 1.

Table 1 Quad-frequency signals of BDS-3

The benefits of the triple-frequency signals of BDS-2 have been well-discussed in previous studies. Compared with dual-frequency signals, triple-frequency signals are more beneficial with regard to AR (Feng 2008), cycle slip detection (Zhang and Li 2016), high-precision positioning (Feng and Li 2010), ionosphere inversion (Spit 2011), etc. Regarding the triple-frequency AR, the three-carrier ambiguity resolution and cascading-integer-resolution methods are commonly used to improve the AR efficiency (Henkel and Günther 2012; Zhao et al. 2015; Zhang and He 2016). In addition to the improved AR efficiency, the positioning performance is improved for the Real-Time Kinematic (RTK) technique (Li et al. 2015; He et al. 2014) and the Precise Point Positioning (PPP) technique (Gu et al. 2015; Guo et al. 2016; Geng and Bock 2013) by applying triple-frequency signals. Besides BDS-2, the triple-frequency positioning of Galileo Navigation Satellite System (Galileo) (Liu et al. 2019), BDS-3 (Li et al. 2020b), and multi-GNSS (Li et al. 2019b, 2020a; Geng et al. 2020) is promising.

Beyond the triple-frequency, Zhang et al. (2020a) investigate the quad-frequency AR of BDS-3. Wang et al. (2018) and Li et al. (2020c) investigates the five-frequency AR of Galileo satellites for long-baseline RTK and PPP, respectively. Compared with triple-frequency signals, quad- and five-frequency signals can in principle improve the positioning (Hatch et al. 2000). The additional signals increase the number of redundant observations, which can enhance the strength of the positioning model, as well as the positioning accuracy and robustness (Weber and Karabatic 2009). Additionally, using more signals is beneficial for estimating or eliminating the frequency-dependent ionospheric delays, which can improve the positioning performance. Furthermore, the additional signals provide more choices of EWL and WL combinations, which may be beneficial for AR and positioning (Li et al. 2015). However, no researchers have comprehensively investigated the benefits of quad-frequency AR and positioning. Hence, in this study, we investigated the aforementioned benefits for quad-frequency signals of BDS-3 with regard to the float resolution, AR, and high-precision positioning with EWL/WL observations.

In the following sections, we first theoretically analyse the quad-frequency improvement for resolving the float solutions based on the between-satellite standalone PPP model. The analysis is also applicable to the RTK model. Then, useful EWL/WL combinations of quad-frequency signals are selected according to the total noise level for the cascading AR (Feng 2008; Feng and Li 2008). The EWL AR performance of the selected combinations is investigated according to the success rates, and the Narrow-Lane (NL) AR performance with EWL/WL ambiguities fixed is investigated according to its derived variance. Additionally, the ambiguity-fixed EWL/WL observations can be directly used for EWL RTK (ERTK) positioning (Li et al. 2017), which also benefits from the quad-frequency signals. Regarding ERTK positioning, the Ionosphere-Free (IF) positioning model is applied with three EWL/WL combinations. Then, the NL observations are used to smooth the EWL/WL observations with a Hatch filter to improve the positioning performance. Finally, two experiments are conducted to evaluate the performance of quad-frequency ERTK positioning.

Fast quad-frequency PPP convergence

In this study, we formulate the PPP model based on the between-satellite Single-Differenced (SD) observations. By applying the corrections of the precise satellite products, i.e. orbits, clock offsets, and differential code biases, the SD observation equations of the code and phase are given as follows:

$$\begin{aligned} \phi_{j} & = \varrho + T - \mu_{j} \iota - \lambda_{j} a_{j} + \varepsilon_{{\varPhi_{j} }} \\ p_{j} & = \varrho + T + \mu_{j} \iota + \varepsilon_{{P_{j} }} , \\ \end{aligned}$$

(1)

where the subscript $j$ represents the frequency number; $\varrho$ and $T$ represent the SD satellite-to-receiver range and tropospheric delay, respectively; and $\iota$ represents the SD ionospheric delay at the first frequency with $\mu_{j} = f_{1}^{2} /f_{j}^{2}$ ($f_{j}$ represents frequency j); $\lambda_{j}$ represents the wavelength of the signal; and $a_{j}$ represents the SD ambiguity in cycles (real value), including the satellite hardware delays of the code and phase. The inclusion of the satellite code hardware delay was explained by Li et al. (2019a). The variations of the differential code biases (Feng et al. 2017) are ignored in this study, as they are ignorable for BDS-3 satellites (Zhang et al. 2019a, b). $\varepsilon_{{\varPhi_{j} }}$ and $\varepsilon_{{P_{j} }}$ represent the observation noises of the phase and code, respectively. The other effects not explicitly shown in (2), such as the relativistic effect, the phase centre offset and variation, and the phase wind-up, have been precisely corrected in advance.

By collecting all the $f$-frequency single-epoch observations, the linearised static PPP observation model is expressed as follows:

$${\text{E}}\left[ {\begin{array}{*{20}c} {\varvec{p}_{k} } \\ {\phi_{k} } \\ \end{array} } \right] = \left[ {\begin{array}{*{20}l} {\varvec{e}_{f} \otimes \varvec{A}_{k} } \hfill & {\varvec{\mu}\otimes \varvec{I}_{s} } \hfill & 0 \hfill \\ {\varvec{e}_{f} \otimes \varvec{A}_{k} } \hfill & { -\varvec{\mu}\otimes \varvec{I}_{s} } \hfill & {\varvec{\varLambda}\otimes \varvec{I}_{s} } \hfill \\ \end{array} } \right]\left[ {\begin{array}{*{20}c}\varvec{\xi}\\ {\begin{array}{*{20}c} {\varvec{\iota}_{k} } \\ \varvec{a} \\ \end{array} } \\ \end{array} } \right] ,$$

(2)

where $k$ represents the epoch number, and $s$ represents the number of SD satellite pairs tracked in epoch $k$. $\phi_{k} = \left[ {\phi_{1,k} , \ldots , \phi_{f,k} } \right]^{\text{T}}$ and $\varvec{p}_{k} = \left[ {\varvec{p}_{1,k} , \ldots , \varvec{p}_{f,k} } \right]^{\text{T}}$ are the observation vectors of the phase and code, respectively, for all frequencies, where the observation vector of each frequency contains the SD observations of $s$ satellite pairs. $\varvec{A}_{k}$ is the design matrix for the parameters $\varvec{\xi}$ consisting of three coordinate unknowns $\varvec{x}$ and one residual zenith tropospheric delay $\tau$ corrected by the UNB3 model (Collins and Langley 1997), together with the Niell mapping function (Niell 1996). $\varvec{\varLambda}= {\text{diag}}\left( {\left[ {\lambda_{1} , \ldots , \lambda_{f} } \right]} \right)$ and $\varvec{a} = \left[ {\varvec{a}_{1} , \ldots , \varvec{a}_{f} } \right]^{\text{T}}$ denote the $f$-frequency SD ambiguities. $\varvec{e}_{f}$ is an f-column vector in which all the elements are 1. $\varvec{I}_{s}$ is an identity matrix with dimension $s$. $\varvec{\mu}= [\mu_{1} , \ldots ,\mu_{f} ]^{\text{T}}$ is the scalar vector for the SD ionosphere parameters ($\varvec{\iota}_{k}$).

The stochastic model of (2) can be expressed as follows:

$$\varvec{Q}_{\text{SD}} = {\text{diag}}\left( {\left[ {\sigma_{p}^{2} ,\sigma_{\phi }^{2} } \right]} \right) \otimes \varvec{I}_{f} \otimes \varvec{P}_{k}^{ - 1} ,$$

(3)

where $\sigma_{p}^{2}$ and $\sigma_{\phi }^{2}$ are the variance scalars of undifferenced observations of the code and phase in the zenith direction, respectively. $\varvec{P}_{k} = \varvec{Q}_{k}^{ - 1}$ is the SD weight matrix, and $\varvec{Q}_{k}$ is the corresponding cofactor matrix that captures the elevation-dependent dispersions of observations (Li et al. 2019a).

By applying the Least-Squares (LS) criterion, we obtain the normal equations. After reducing the epoch-wised ionospheric parameter $\varvec{\iota}_{k}$, the normal matrix of $\varvec{\xi}$ and a over $K$ epochs is given as follows:

$$\left[ {\begin{array}{*{20}c} {\left( {f\sigma_{p + \phi }^{ - 2} - \frac{{\sigma_{p - \phi }^{ - 4} \varvec{e}_{f}^{\text{T}} \varvec{\mu \mu }^{\text{T}} \varvec{e}_{f} }}{{\sigma_{p + \phi }^{ - 2}\varvec{\mu}^{\text{T}}\varvec{\mu}}}} \right) \otimes \mathop \sum \limits_{k = 1}^{K} \varvec{N}_{k} } & {\left( {\sigma_{\phi }^{ - 2} \varvec{e}_{f}^{\text{T}}\varvec{\varLambda}+ \frac{{\sigma_{\phi }^{ - 2} \sigma_{p - \phi }^{ - 2} \varvec{e}_{f}^{\text{T}} \varvec{\mu \mu }^{\text{T}} {\varvec{\Lambda}}}}{{\sigma_{p + \phi }^{ - 2}\varvec{\mu}^{\text{T}}\varvec{\mu}}}} \right) \otimes \mathop \sum \limits_{k = 1}^{K} \varvec{A}_{k}^{\text{T}} \varvec{P}_{k} } \\ {\left( {\sigma_{\phi }^{ - 2} \varvec{\varLambda e}_{f} + \frac{{\sigma_{\phi }^{ - 2} \sigma_{p - \phi }^{ - 2} {\varvec{\Lambda}}\varvec{\mu \mu }^{\text{T}} \varvec{e}_{f} }}{{\sigma_{p + \phi }^{ - 2}\varvec{\mu}^{\text{T}}\varvec{\mu}}}} \right) \otimes \mathop \sum \limits_{k = 1}^{K} \varvec{P}_{k} \varvec{A}_{k} } & {\left( {\sigma_{\phi }^{ - 2}\varvec{\varLambda}^{2} - \frac{{\sigma_{\phi }^{ - 4} {\varvec{\Lambda}}\varvec{\mu \mu }^{\text{T}} {\varvec{\Lambda}}}}{{\sigma_{p + \phi }^{ - 2}\varvec{\mu}^{\text{T}}\varvec{\mu}}}} \right) \otimes \mathop \sum \limits_{k = 1}^{K} \varvec{P}_{k} } \\ \end{array} } \right] ,$$

(4)

where $\varvec{N}_{k} = \varvec{A}_{k}^{\text{T}} \varvec{P}_{k} \varvec{A}_{k}$, $\sigma_{p + \phi }^{ - 2} = \sigma_{p}^{ - 2} + \sigma_{\phi }^{ - 2}$, and $\sigma_{p - \phi }^{ - 2} = \sigma_{p}^{ - 2} - \sigma_{\phi }^{ - 2}$. We further reduce the ambiguity parameters in the normal matrix as follows:

$$\varvec{N}_{\xi ,K} = \left( {f\sigma_{p + \phi }^{ - 2} - \frac{{\sigma_{p - \phi }^{ - 4} \varvec{e}_{f}^{\text{T}} \varvec{\mu \mu }^{\text{T}} \varvec{e}_{f} }}{{\sigma_{p + \phi }^{ - 2}\varvec{\mu}^{\text{T}}\varvec{\mu}}}} \right)\mathop \sum \limits_{k = 1}^{K} \varvec{N}_{k} - \left( {\frac{{\left( {f - 1} \right)\sigma_{\phi }^{ - 2} + \left( {f + 3} \right)\sigma_{p}^{ - 2} }}{{1 + \sigma_{p}^{ - 2} \sigma_{\phi }^{2} }}\frac{{\varvec{e}_{f}^{\text{T}} \varvec{\mu \mu }^{\text{T}} \varvec{e}_{f} }}{{\varvec{\mu}^{\text{T}}\varvec{\mu}}}} \right)\varvec{A}_{\Sigma }^{\text{T}} \varvec{P}_{\Sigma }^{ - 1} \varvec{A}_{\Sigma } ,$$

(5)

where $\varvec{A}_{\Sigma } = \mathop \sum \nolimits_{k = 1}^{K} \varvec{P}_{k} \varvec{A}_{k}$. For a short timespan, it is adequate to assume a time-invariable geometry for each satellite, i.e. $\varvec{A}_{k} = \varvec{A}$ and $\varvec{P}_{k} = \varvec{P}$. Then, the covariance matrix of $\varvec{\xi}$ for epoch $K$ is

$$\varvec{Q}_{\xi ,K} = \frac{1}{K}\frac{{\mu_{\text{sq}} }}{{f\sigma_{p + \phi }^{ - 2} \mu_{\text{sq}} - \left[ {\sigma_{p}^{ - 2} + f\sigma_{\phi }^{ - 2} } \right]\mu_{\Sigma }^{2} }}\left( {\varvec{A}^{\text{T}} \varvec{PA}} \right)^{ - 1} ,$$

(6)

where $\mu_{\text{sq}} =\varvec{\mu}^{\text{T}}\varvec{\mu}$ denotes the sum of squared ionospheric coefficients and “sq” denotes “squared”. $\varvec{ }\mu_{\Sigma } = \varvec{e}_{f}^{\text{T}}\varvec{\mu}$ denotes the sum of ionospheric coefficients. Regarding the analytical formula, the positioning precision gain from quad-frequency signals can be expressed with respect to the dual- and triple-frequency signals as follows:

$$\frac{{\varvec{Q}_{\xi ,K} \left( {f = 2} \right)}}{{\varvec{Q}_{\xi ,K} \left( {f = 4} \right)}} = \frac{{0.2439\sigma_{p}^{ - 2} - 11.0243\sigma_{\phi }^{ - 2} }}{{0.0002\sigma_{p}^{ - 2} - 1.9997\sigma_{\phi }^{ - 2} }}$$

(7)

$$\frac{{\varvec{Q}_{\xi ,K} \left( {f = 3} \right)}}{{\varvec{Q}_{\xi ,K} \left( {f = 4} \right)}} = \frac{{0.2439\sigma_{p}^{ - 2} - 11.0243\sigma_{\phi }^{ - 2} }}{{0.2341\sigma_{p}^{ - 2} - 5.2976\sigma_{\phi }^{ - 2} }} ,$$

(8)

where B1C and B1I signals are used for the dual-frequency case and $\varvec{\mu}^{\text{T}} = \left[ {1, 1.0184} \right]$; B1C, B1I, and B2a signals are used for the triple-frequency case and $\varvec{\mu}^{\text{T}} = \left[ {1, 1.0184, 1.7933} \right]$; and B1C, B1I, B3I, and B2a signals are used for the quad-frequency case and $\varvec{\mu}^{\text{T}} = \left[ {1, 1.0184, 1.5424, 1.7933} \right]$. For example, $\frac{{\varvec{Q}_{\xi ,K} \left( {f = 2} \right)}}{{\varvec{Q}_{\xi ,K} \left( {f = 4} \right)}} \approx 5.5$ and $\frac{{\varvec{Q}_{\xi ,K} \left( {f = 3} \right)}}{{\varvec{Q}_{\xi ,K} \left( {f = 4} \right)}} \approx 2.1$ for $\sigma_{p} = 100\sigma_{\phi }$. The results roughly give the intuitive improvements of how PPP convergence and precision can be gained from quad-frequency signals. Equations (7) and (8) are also applicable for the RTK positioning model, while Double-Differenced (DD) observations and stochastic model are used in the positioning models of (2) and (3) with the same design matrix and cofactor matrix. Hence, the same improvements can be obtained in the float resolution of RTK positioning by applying quad-frequency signals.

EWL/WL combinations and quad-frequency AR

The multi-frequency observations allow the formation of more useful combinations-particularly EWL/WL combinations-to enhance the AR efficiency. In this section, we first select the useful EWL/WL combinations according to the total noise level in the cycle. Then, we present two procedures for cascading quad-frequency AR.

Selection of EWL/WL combinations

Omitting the epoch subscript, the combined observation equations of quad-frequency signals are given as follows:

$$\phi_{{\left( {i,j,k,h} \right)}} = \varrho + T - \mu_{{\left( {i,j,k,h} \right)}} \iota - \lambda_{{\left( {i,j,k,h} \right)}} a_{{\left( {i,j,k,h} \right)}} + \varepsilon_{{\phi_{{\left( {i,j,k,h} \right)}} }} ,$$

(9)

where

$$\phi_{{\left( {i,j,k,h} \right)}} = \frac{{i \cdot f_{1} \cdot \phi_{1} + j \cdot f_{2} \cdot \phi_{2} + k \cdot f_{3} \cdot \phi_{3} + h \cdot f_{4} \cdot \phi_{4} }}{{i \cdot f_{1} + j \cdot f_{2} + k \cdot f_{3} + h \cdot f_{4} }}$$

(10)

$$\mu_{{\left( {i,j,k,h} \right)}} = \frac{{f_{1}^{2} \left( {i/f_{1} + j/f_{2} + k/f_{3} + h/f_{4} } \right)}}{{i \cdot f_{1} + j \cdot f_{2} + k \cdot f_{3} + h \cdot f_{4} }}$$

(11)

$$\lambda_{{\left( {i,j,k,h} \right)}} = \frac{c}{{i \cdot f_{1} + j \cdot f_{2} + k \cdot f_{3} + h \cdot f_{4} }} .$$

(12)

Here, $p_{{\left( {i,j,k,h} \right)}}$ has a similar expression to $\phi_{{\left( {i,j,k,h} \right)}}$. $c$ represents the velocity of light, and the coefficients of the combinations $i, j, k,$ and $h$ are all integers.

Many methods have been proposed to find the optimal combinations for triple-frequency signals for different purposes (Richert and El-Sheimy 2007; Feng 2008; Li et al. 2010a, 2015). Overall, the identified useful combinations are similar or even equivalent among all the methods (Li 2018; Li et al. 2015). In this study, we employ the total noise level relative to the combined wavelength in cycles to find the useful quad-frequency combinations for the AR based on either a Geometry-Based (GB) or Geometry-Free (GF) model.

The total noise level in cycles is defined to compensate the ionospheric and tropospheric biases and phase noises, as follows (Feng 2008):

$$\sigma_{\text{TC}} = \sqrt {\delta_{\tau }^{2} + \mu_{{\left( {i,j,k,h} \right)}}^{2} \delta_{\iota }^{2} + \sigma_{{\phi_{{\left( {i,j,k,h} \right)}} }}^{2} } /\lambda_{{\left( {i,j,k,h} \right)}} ,$$

(13)

where $\delta_{\iota }$ and $\delta_{\tau }$ represent the ionospheric and tropospheric biases, respectively. The subscript “TC” denotes the Total noise level in Cycles. $\sigma_{{\phi_{{\left( {i,j,k,h} \right)}} }}^{2} = \frac{{\left( {i \cdot f_{1} } \right)^{2} + \left( {j \cdot f_{2} } \right)^{2} + \left( {k \cdot f_{3} } \right)^{2} + \left( {h \cdot f_{4} } \right)^{2} }}{{\left( {i \cdot f_{1} + j \cdot f_{2} + k \cdot f_{3} + h \cdot f_{4} } \right)^{2} }}\sigma_{\phi }^{2}$. We empirically give the values of $\delta_{\iota }$, $\delta_{\tau }$, and $\sigma_{\phi }$ for calculating the total noise level. Equation (13) can be applied to either SD or DD observations, while the biases and noises correspond to SD and DD observations, respectively.

Table 2 presents several useful EWL/WL combinations and their total noise levels under different given error budgets. The phase observation noises are equally set as 5 mm for different frequencies. Because the sum of the coefficients always equals 0 for the EWL/WL combinations (Li et al. 2015; Zhang et al. 2020b), only three EWL/WL combinations are linearly independent. For each error budget, the three linearly independent combinations with the lowest total noise levels are in bold. The combinations $\phi_{{\left( {0,0,1, - 1} \right)}}$ and $\phi_{{\left( {0,1, - 3,2} \right)}}$ are preferable when the ionospheric bias is small; otherwise, $\phi_{{\left( { - 3,4, - 3,2} \right)}}$ and $\phi_{{\left( {2,0, - 7,5} \right)}}$ are more suitable, as they are less sensitive to the ionospheric bias. The combination $\phi_{{\left( {1, - 1,0,0} \right)}}$ is always the best owing to its longest wavelength. Overall, the total noise levels of the selected combinations presented in Table 2 are lower than those of the triple-frequency combinations (Li et al. 2010a), indicating a better EWL/WL AR performance.

Table 2 Total noise level $\sigma_{\text{TC}}$ for useful combinations under different given error budgets

Cascading AR in parameter domain

The cascading AR in the parameter domain starts with solving the uncombined float ambiguity solutions of quad-frequency signals using the GB model. The EWL/WL and NL ambiguities are obtained by transforming the uncombined float solutions into the specific combinations and then successively fixed to integers. Thus, the combinations are formed according to the resolved parameters instead of the observation. The uncombined DD observation equations of quad-frequency signals are as follows:

$${\text{E}}\left[ {\begin{array}{*{20}c} {\bar{\varvec{p}}} \\ {\bar{\phi }} \\ \end{array} } \right] = \left[ {\begin{array}{*{20}l} {\varvec{e}_{4} \otimes \varvec{A}} \hfill & {\varvec{\mu}\otimes \varvec{I}_{s} } \hfill & 0 \hfill \\ {\varvec{e}_{4} \otimes \varvec{A}} \hfill & {\varvec{\mu}\otimes \varvec{I}_{s} } \hfill & {\varvec{\varLambda}\otimes \varvec{I}_{s} } \hfill \\ \end{array} } \right]\left[ {\begin{array}{*{20}c} \varvec{x} \\ {\begin{array}{*{20}c} {\bar{\varvec{\iota }}} \\ \varvec{z} \\ \end{array} } \\ \end{array} } \right] ,$$

(14)

where $\bar{\varvec{p}}$ and $\bar{\phi }$ represent the quad-frequency DD observations of the code and phase, respectively. $\bar{\varvec{\iota }}$ represents the DD ionospheric parameter vector. $\varvec{z} = \left[ {\varvec{z}_{1}^{\text{T}} , \ldots , \varvec{z}_{4}^{\text{T}} } \right]^{\text{T}}$ represents the uncombined quad-frequency DD ambiguities. The other symbols have the same meanings as those in (2), except that $f = 4$, and the epoch number is omitted. Here, we neglect the DD residual zenith tropospheric parameter owing to its correlation with height component (Li et al. 2010b). The covariance matrix of the DD observations is given as follows:

$$\varvec{Q}_{\text{DD}} = 2{\text{diag}}\left( {\left[ {\sigma_{p}^{2} ,\sigma_{\phi }^{2} } \right]} \right) \otimes \varvec{I}_{4} \otimes \varvec{P}^{ - 1} .$$

(15)

After necessary derivations, the covariance matrix of float ambiguity solutions is given as

$$\varvec{Q}_{{\varvec{\hat{z}\hat{z}}}} = 2\left( {\varvec{\varLambda}^{ - 1} \left( {\sigma_{\phi }^{2} \varvec{I}_{4} + \sigma_{p}^{2} \frac{{\varvec{\mu \mu }^{\text{T}} }}{{\mu_{\text{sq}} }}} \right)\varvec{\varLambda}^{ - 1} } \right) \otimes \varvec{Q} + 2\sigma_{p}^{2}\varvec{\varLambda}^{ - 1} \frac{{\mu_{\text{sq}} \varvec{e}_{4} \varvec{e}_{4}^{\text{T}} + \mu_{\Sigma } \varvec{\mu e}_{4}^{\text{T}} + \mu_{\Sigma } \varvec{e}_{4}\varvec{\mu}^{\text{T}} + \frac{{\mu_{\Sigma }^{2} }}{{\mu_{\text{sq}} }}\varvec{\mu \mu }^{\text{T}} }}{{f\mu_{\text{sq}} - \mu_{\Sigma }^{2} }}\varvec{\varLambda}^{ - 1} \otimes \varvec{A}\left( {\varvec{A}^{\text{T}} \varvec{PA}} \right)^{ - 1} \varvec{A}^{\text{T}} .$$

(16)

In the cascading AR in the parameter domain, one can transform the float ambiguities $\hat{\varvec{z}}$ into three EWL/WL ambiguities and one uncombined ambiguity:

$$\left[ {\begin{array}{*{20}c} {\hat{\varvec{z}}_{E} } \\ {\hat{\varvec{z}}_{1} } \\ \end{array} } \right] = \left[ {\begin{array}{*{20}c} {\varvec{Z}_{E}^{\text{T}} } \\ {\varvec{c}_{1}^{\text{T}} } \\ \end{array} } \right]\hat{\varvec{z}},\,\left[ {\begin{array}{*{20}c} {\varvec{Q}_{{\hat{\varvec{z}}_{E} \hat{\varvec{z}}_{E} }} } & {\varvec{Q}_{{\hat{\varvec{z}}_{E} \hat{\varvec{z}}_{1} }} } \\ {\varvec{Q}_{{\hat{\varvec{z}}_{1} \hat{\varvec{z}}_{E} }} } & {\varvec{Q}_{{\hat{\varvec{z}}_{1} \hat{\varvec{z}}_{1} }} } \\ \end{array} } \right] ,$$

(17)

with $\varvec{Q}_{{\hat{\varvec{z}}_{E} \hat{\varvec{z}}_{E} }} = \varvec{Z}_{E}^{\text{T}} \varvec{Q}_{{\varvec{\hat{z}\hat{z}}}} \varvec{Z}_{E} , \varvec{ Q}_{{\hat{\varvec{a}}_{1} \hat{\varvec{a}}_{1} }} = \varvec{c}_{1}^{\text{T}} \varvec{Q}_{{\varvec{\hat{z}\hat{z}}}} \varvec{c}_{1} , \varvec{Q}_{{\hat{\varvec{a}}_{E} \hat{\varvec{a}}_{1} }} = \varvec{Z}_{E}^{\text{T}} \varvec{Q}_{{\varvec{\hat{z}\hat{z}}}} \varvec{c}_{1}$, and $\varvec{c}_{1}^{\text{T}} = \left[ {1,0,0,0} \right] \otimes \varvec{I}_{s}$. Here, the transformation matrix is defined by three EWL/WL combinations as $\varvec{Z}_{E}^{\text{T}} = \left[ {\begin{array}{*{20}c} {\begin{array}{*{20}c} 1 \quad & { - 1} \\ \end{array} } \quad & {\begin{array}{*{20}c} 0 \quad & 0 \\ \end{array} } \\ {\begin{array}{*{20}c} 0 \quad & 0 \\ \end{array} } \quad & {\begin{array}{*{20}c} 1 \quad & { - 1} \\ \end{array} } \\ {\begin{array}{*{20}c} 0 \quad & 1 \\ \end{array} } \quad & { - \begin{array}{*{20}c} 3 \quad & 2 \\ \end{array} } \\ \end{array} } \right] \otimes \varvec{I}_{s}$. One can then fix EWL/WL ambiguities according to the float solutions $\hat{\varvec{z}}_{E}$ and $\varvec{Q}_{{\hat{\varvec{z}}_{E} \hat{\varvec{z}}_{E} }}$. Usually, the rounding method can fix them with a success rate of nearly 100%. Once the EWL/WL ambiguities are fixed, which is denoted as $\varvec{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\smile}$}}{z} }_{E}$, the remained uncombined ambiguities are updated as follows:

$$\tilde{\varvec{z}}_{1} = \hat{\varvec{z}}_{1} - \varvec{Q}_{{\hat{\varvec{z}}_{E} \hat{\varvec{z}}_{1} }} \varvec{Q}_{{\hat{\varvec{z}}_{E} \hat{\varvec{z}}_{E} }}^{ - 1} \left( {\hat{\varvec{z}}_{E} - \varvec{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\smile}$}}{z} }_{E} } \right)$$

(18)

$$\varvec{Q}_{{\tilde{\varvec{z}}_{1} \tilde{\varvec{z}}_{1} }} = \varvec{Q}_{{\hat{\varvec{z}}_{1} \hat{\varvec{z}}_{1} }} - \varvec{Q}_{{\hat{\varvec{z}}_{1} \hat{\varvec{z}}_{E} }} \varvec{Q}_{{\hat{\varvec{z}}_{E} \hat{\varvec{z}}_{E} }}^{ - 1} \varvec{Q}_{{\hat{\varvec{z}}_{E} \hat{\varvec{z}}_{1} }} .$$

(19)

One can further try to fix these ambiguities according to $\tilde{\varvec{z}}_{1}$ and $\varvec{Q}_{{\tilde{\varvec{z}}_{1} \tilde{\varvec{z}}_{1} }}$ by using the LAMBDA method (Teunissen 1995).

Cascading AR in measurement domain

Owing to the extra-long wavelength, one can often fix EWL/WL ambiguities instantaneously by simply applying the GF model. This model is preferable for EWL AR compared to the GB model because all the geometric errors are completely eliminated, and the satellite positions do not need to be computed. The EWL/WL ambiguities are simply estimated by applying the linear combinations of phase and code observations, as follows:

$$\hat{z}_{{\left( {i,j,k,h} \right)}} = \frac{{p_{{\left( {a,b,c,d} \right)}} - \phi_{{\left( {i,j,k,h} \right)}} }}{{\lambda_{{\left( {i,j,k,h} \right)}} }},\,\sigma_{{\hat{z}_{{\left( {i,j,k,h} \right)}} }}^{2} = \frac{{\sigma_{{p_{{\left( {a,b,c,d} \right)}} }}^{2} + \sigma_{{\phi_{{\left( {i,j,k,h} \right)}} }}^{2} }}{{\lambda_{{\left( {i,j,k,h} \right)}}^{2} }}.$$

(20)

The float ambiguities may contain the ionospheric bias, which is given as

$$b_{{\hat{z}_{{\left( {i,j,k,h} \right)}} }} = \frac{{\mu_{{\left( {a,b,c,d} \right)}} + \mu_{{\left( {i,j,k,h} \right)}} }}{{\lambda_{{\left( {i,j,k,h} \right)}} }}\delta_{\iota } .$$

(21)

The total noise level in cycles of such a float ambiguity solution is defined to include the ionospheric biases, phase noises, and code noises, as follows:

$$\sigma_{\text{TCN}} = \sqrt {b_{{\hat{z}_{{\left( {i,j,k,h} \right)}} }}^{2} + \sigma_{{\phi_{{\left( {i,j,k,h} \right)}} }}^{2} + \sigma_{{p_{{\left( {a,b,c,d} \right)}} }}^{2} } /\lambda_{{\left( {i,j,k,h} \right)}} .$$

(22)

The subscript “TCN” denotes the Total noise level including Code Noises. The useful combination schemes can be selected by examining their total noises. For different given error budgets, Table 3 presents several useful combination schemes and their total noise levels. The precision of uncombined DD observation noises is equally set to 5 mm and 0.5 m for the phase and code of different frequencies, respectively. To intuitively show the capability of instantaneous EWL/WL AR, the bias-affected rounding success rates are presented in Table 3 (Teunissen 1997).

Table 3 Total noise level $\sigma_{{{\text{TC}}2}}$ and success rate under different given ionospheric biases of 10 cm, 20 cm and 100 cm, respectively

Because in the GF model, the geometric errors are completely eliminated by involving code combinations, the three selected EWL/WL schemes probably changed depending on the magnitudes of the ionospheric biases. Fortunately, the instantaneous AR success rate is 100% for almost all the selected schemes except for one scheme under the extreme ionospheric bias of 1 m. Compared with the GB model, the involvement of code combinations eliminates the geometric errors and thus makes the AR significantly simpler from the viewpoint of programming. For the GF model, when one EWL/WL ambiguity is fixed, it can be used to solve the next EWL/WL ambiguity (Li et al. 2010a). However, in such a case, incorrect fixing of the previous EWL/WL ambiguity will definitely affect the following fixing.

Li et al. (2015, 2017) proved that for an f-frequency system, only (f-1) EWL/WL combinations are independent, and all the other EWL/WL combinations can be recovered by the linear combinations of these (f−1) combinations. They further advised to first fix (f−1) EWL ambiguities according to the GF model and then the NL ambiguities. For the quad-frequency signals, three EWL/WL ambiguities are first resolved according to the GF model (20) and fixed to their integers. Then, all the EWL/WL integer ambiguities are available through linear combination and can be used to correct the corresponding EWL/WL observations. For the NL AR, owing to its relatively short wavelength, the IF GB model is typically applied because the IF operation can eliminate the ionospheric effect, and the GB model has better strength than the GF model. Alternatively, as suggested by Li et al. (2010a), one can form a Geometry- and Ionosphere-Free (GIF) model for triple-frequency signals to solve NL ambiguity. For the quad-frequency case, we can analogously form the GIF model as follows:

$$\hat{z}_{1} = \frac{{b_{1} \overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\smile}$}}{\phi }_{{{\text{E}}1}} + b_{2} \overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\smile}$}}{\phi }_{{{\text{E}}2}} + b_{3} \overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\smile}$}}{\phi }_{{{\text{E}}3}} - \overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\smile}$}}{\phi }_{1} }}{{\lambda_{1} }},\,\sigma_{{\hat{z}_{1} }}^{2} = \frac{{4\gamma^{2} \cdot \sigma_{\phi }^{2} }}{{\lambda_{1}^{2} }}.$$

(23)

The coefficients satisfy

$$\begin{aligned} & b_{1} + b_{2} + b_{3} = 1 \\ & b_{1} \mu_{{\left( {1, - 1,0,0} \right)}} + b_{2} \mu_{{\left( {0,0,1, - 1} \right)}} + b_{3} \mu_{{\left( {0,1, - 3,2} \right)}} = \mu_{1} \\ & \gamma^{2} = \left( {\frac{{b_{1} f_{1} }}{{f_{{\left( {1, - 1,0,0} \right)}} }} - 1} \right)^{2} + \left( {\frac{{b_{3} f_{2} }}{{f_{{\left( {0,1, - 3,2} \right)}} }} - \frac{{b_{1} f_{2} }}{{f_{{\left( {1, - 1,0,0} \right)}} }}} \right)^{2} + \left( {\frac{{b_{2} f_{3} }}{{f_{{\left( {0,0,1, - 1} \right)}} }} - \frac{{3b_{3} f_{3} }}{{f_{{\left( {0,1, - 3,2} \right)}} }}} \right)^{2} + \left( {\frac{{2b_{3} f_{4} }}{{f_{{\left( {0,1, - 3,2} \right)}} }} - \frac{{b_{2} f_{4} }}{{f_{{\left( {0,0,1, - 1} \right)}} }}} \right)^{2} = { \hbox{min} } .\\ \end{aligned}$$

(24)

By solving this minimisation problem, we obtain the coefficients $b_{1} = 0.2357$, $b_{2} = - 1.5049$, $b_{3} = 2.2692$, and $\sigma_{{\hat{z}_{1} }} = 1302\sigma_{\phi }$. If we take the precision of the undifferenced phase observation as 2 mm, $\sigma_{{\hat{z}_{1} }} = 2.60$ cycles, which is significantly better than the triple-frequency case of BDS-2 (Li et al. 2010a).

Long-distance ERTK with quad-frequency observations

Once the EWL ambiguities are fixed, the ambiguity-corrected EWL observation plays the role of a pseudo-range, but with a higher precision. The RTK positioning is immediately started with ambiguity-corrected EWL observations, which is referred to as ERTK (Li et al. 2017). The quad-frequency observations can also benefit the ERTK positioning by increasing the number of available EWL observations. In this section, we investigate the quad-frequency ERTK performance by directly using IF EWL/WL observations and their smoothed counterparts.

ERTK with EWL observations

For the long baseline, ionospheric delays cannot be completely cancelled in DD observations and must be considered in the positioning model. One can simply eliminate the residual ionospheric delays by using the IF observations. The IF ERTK observation equations are given as follows:

$$\left[ {\begin{array}{*{20}c} {\bar{\phi }_{{{\text{IF}},{\text{E}}13}} } \\ {\bar{\phi }_{{{\text{IF}},{\text{E}}23}} } \\ \end{array} } \right] = \left[ {\begin{array}{*{20}c} {\frac{{\mu_{{{\text{E}}3}} }}{{\mu_{{{\text{E}}3}} - \mu_{{{\text{E}}1}} }}\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\smile}$}}{\phi }_{{{\text{E}}1}} - \frac{{\mu_{{{\text{E}}1}} }}{{\mu_{{{\text{E}}3}} - \mu_{{{\text{E}}1}} }}\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\smile}$}}{\phi }_{{{\text{E}}3}} } \\ {\frac{{\mu_{{{\text{E}}3}} }}{{\mu_{{{\text{E}}3}} - \mu_{{{\text{E}}2}} }}\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\smile}$}}{\phi }_{{{\text{E}}2}} - \frac{{\mu_{{{\text{E}}2}} }}{{\mu_{{{\text{E}}3}} - \mu_{{{\text{E}}2}} }}\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\smile}$}}{\phi }_{{{\text{E}}3}} } \\ \end{array} } \right] = \left[ {\begin{array}{*{20}c} \varvec{A} \\ \varvec{A} \\ \end{array} } \right]\varvec{x}, 2\sigma_{\phi }^{2} \left( {\varvec{R}_{\text{IF}}^{\text{T}} \varvec{R}_{\text{IF}} } \right) \otimes \varvec{Q} ,$$

(25)

where the subscripts ${\text{E}}1$, ${\text{E}}2$, and ${\text{E}}3$ denote the EWL/WL combinations $\left( {1, - 1,0,0} \right)$, $\left( {0,0,1, - 1} \right)$, and $\left( {0,1, - 3,2} \right)$, respectively. The superscript ‘check’ is applied to denote the ambiguity-corrected EWL/WL observations. $\varvec{R}_{\text{IF}}^{\text{T}}$ is the transformation matrix from the uncombined quad-frequency DD observations to the IF EWL observations,

$$\varvec{R}_{\text{IF}}^{\text{T}} = \left[ {\begin{array}{*{20}c} {\frac{{\mu_{{{\text{E}}3}} }}{{\mu_{{{\text{E}}3}} - \mu_{{{\text{E}}1}} }}} & 0 & {\frac{{ - \mu_{{{\text{E}}1}} }}{{\mu_{{{\text{E}}3}} - \mu_{{{\text{E}}1}} }}} \\ 0 & {\frac{{\mu_{{{\text{E}}3}} }}{{\mu_{{{\text{E}}3}} - \mu_{{{\text{E}}2}} }}} & {\frac{{ - \mu_{{{\text{E}}2}} }}{{\mu_{{{\text{E}}3}} - \mu_{{{\text{E}}2}} }}} \\ \end{array} } \right]\left[ {\begin{array}{*{20}l} {f_{1} /f_{{{\text{E}}1}} } \hfill & { - f_{2} /f_{{{\text{E}}1}} } \hfill & 0 \hfill & 0 \hfill \\ 0 \hfill & 0 \hfill & {f_{3} /f_{{{\text{E}}2}} } \hfill & { - f_{4} /f_{{{\text{E}}2}} } \hfill \\ 0 \hfill & {f_{2} /f_{{{\text{E}}3}} } \hfill & { - 3f_{3} /f_{{{\text{E}}3}} } \hfill & {2f_{4} /f_{{{\text{E}}3}} } \hfill \\ \end{array} } \right] .$$

(26)

Suppose that the weight matrix $\varvec{P}$ is diagonal, for simplicity. The covariance matrix of the LS estimate in the IF positioning model is given as follows:

$$\varvec{Q}_{{\varvec{\hat{x}\hat{x}}}}^{\text{IF}} = \alpha \left( {\varvec{A}^{\text{T}} \varvec{PA}} \right)^{ - 1} ,\,\alpha = \frac{{2\sigma_{\phi }^{2} }}{{\varvec{e}_{2}^{\text{T}} \left( {\varvec{R}_{\text{IF}}^{\text{T}} \varvec{R}_{\text{IF}} } \right)^{ - 1} \varvec{e}_{2} }} \approx 10057.9\sigma_{\phi }^{2} .$$

(27)

For $\sigma_{\phi } = 2$ mm, $\sqrt \alpha = 0.2006$ m, which empirically means that the positioning precision is decimetre- to meter-level and better than that of the triple-frequency case even without code observations (Li. et al. 2017).

ERTK with smoothed EWL observations

As indicated by (27), the positioning precision is low although the ionospheric biases are eliminated in the IF ERTK model. The main reason for this is that the IF EWL/WL observations are too noisy, and the precisions of IF EWL/WL observations are as high as $\sigma_{{\bar{\phi }_{{{\text{IF}},{\text{E}}13}} }} = 467.8\sigma_{\phi }$ and $\sigma_{{\bar{\phi }_{{{\text{IF}},{\text{E}}13}} }} = 149.2\sigma_{\phi }$. Hence, to improve the IF ERTK precision, we can smooth the IF EWL/WL observations by using IF DD observations $\bar{\phi }_{{{\text{IF}},14}} = \frac{{f_{1}^{2} \bar{\phi }_{1} - f_{4}^{2} \bar{\phi }_{4} }}{{f_{1}^{2} - f_{4}^{2} }}$ with $\sigma_{{\bar{\phi }_{{{\text{IF}},14}} }} = 5.177\sigma_{\phi }$. For instance, if we smooth the IF EWL observations $\bar{\phi }_{{{\text{IF}},{\text{E}}13}}$ with IF observations $\bar{\phi }_{{{\text{IF}},14}}$, the smoothed observation after K epochs is given as follows (Hatch 1982):

$$\tilde{\phi }_{{{\text{IF}},{\text{E}}13}} \left( K \right) = \frac{1}{K}\mathop \sum \limits_{k = 1}^{K} \bar{\phi }_{{{\text{IF}},{\text{E}}13}} \left( K \right) + \frac{1}{K}\left( {\left( {K - 1} \right)\bar{\phi }_{{{\text{IF}},14}} \left( K \right) - \mathop \sum \limits_{k = 1}^{K - 1} \bar{\phi }_{{{\text{IF}},14}} \left( k \right)} \right) .$$

(28)

If we ignore the correlations between $\bar{\phi }_{{{\text{IF}},{\text{E}}13}}$ and $\bar{\phi }_{{{\text{IF}},14}}$ and between epochs, the precision of smoothed observation is given as follows:

$$\sigma_{{\tilde{\phi }_{{{\text{IF}},{\text{E}}13}} \left( K \right)}} = \sqrt {\frac{{\sigma_{{\bar{\phi }_{{{\text{IF}},{\text{E}}13}} }}^{2} + \left( {K - 1} \right)\sigma_{{\bar{\phi }_{{{\text{IF}},14}} }}^{2} }}{K}} \approx \sigma_{{\bar{\phi }_{{{\text{IF}},14}} }} .$$

(29)

Clearly, the precision of smoothed EWL/WL observations is significantly improved and is nearly equal to $\sigma_{{\bar{\phi }_{{{\text{IF}},14}} }}$ if a certain period of smoothing is applied. The smoothed IF EWL observations, i.e. $\tilde{\phi }_{{{\text{IF}},{\text{E}}12}}$ and $\tilde{\phi }_{{{\text{IF}},{\text{E}}23}}$, can then be used instead of $\bar{\phi }_{{{\text{IF}},{\text{E}}13}}$ and $\bar{\phi }_{{{\text{IF}},{\text{E}}23}}$ to form the positioning model (25). Considering that the precision $\sigma_{{\bar{\phi }_{{{\text{IF}},14}} }}$ is significantly lower than $\sigma_{{\bar{\phi }_{{{\text{IF}},{\text{E}}13}} }}$ and $\sigma_{{\bar{\phi }_{{{\text{IF}},{\text{E}}23}} }}$, the positioning precision can be significantly improved after a period of convergence.

Experiment and analysis of ERTK positioning

The quad-frequency data of BDS-3 were collected using Trimble Alloy receivers with a sampling interval of 1 s at four stations. Two baselines were formed, with lengths of 300 and 27 km. The long-baseline data were collected on Day of Year (Doy) 007, 2020, and the short-baseline data were collected on Doy 158, 2020. The coordinates of the four stations were precisely known and served as references. The proposed ERTK models were implemented in the TJRTK software, which was developed at Tongji University for multi-frequency multi-GNSS RTK processing and relevant engineering and scientific applications.

In the data processing, only quad-frequency data of the BDS-3 system were used, and the cut-off elevation was set as 10. The elevation-dependent stochastic model $\sigma = \frac{1.02}{\sin \theta + 0.02}\sigma_{90^\circ }$ was applied for the undifferenced measurements with a zenith precision of $\sigma_{90^\circ } = 2$ mm for the phase and 0.2 m for the code. These parameters were obtained according to the evaluation of the BDS stochastic characteristics based on the study of Li (2016). Although the data were post-processed, the processing was completely analogous to the real-time processing; i.e. the data loading and all computations were implemented epoch-by-epoch, and only the broadcasting ephemeris was used.

Figure 1 shows the number of satellites observed at the rover station of the short baseline and the corresponding Dilution of Precision (DOP) values. The number of tracked satellites ranged from 7 to 10 and was 9 on average. The DOP values varied within a range of 3, indicating that the BDS-3 system was able to constantly provide high-quality positioning service independently.

Fig. 1

Number of observed satellites and corresponding DOP values of the short baseline

Figure 2 shows the positioning errors along north (N), east (E), and up (U) directions of the short baseline with IF ERTK model. The positioning errors were unbiased but very noisy (several decimetres) due to the large noises of the IF EWL observations. The ERTK positioning errors obtained when the smoothed IF EWL observations were used are shown in Fig. 3. The positioning precision was significantly improved to the centimetre level after a short period of convergence. The Root-Mean-Square (RMS) values of the positioning errors are presented in Table 4. The positioning accuracy was comparable to that of normal RTK solutions with all the frequency ambiguities fixed.

Fig. 2

ERTK positioning errors with IF EWL observations for the short baseline

Fig. 3

ERTK positioning errors of smoothed IF EWL observations for the short baseline

Table 4 RMS values of the positioning errors for the short baseline

We also evaluated the positioning performance of ERTK for the 300-km baseline. The number of satellites was slightly smaller than that for the short baseline, as shown in Fig. 4, because the observation date was earlier, and some BDS-3 satellites were not in service at the time. The DOP values were slightly larger, but overall, they were comparable to those for the short baseline. Figures 5 and 6 show the ERTK positioning errors for the application of IF and smoothed IF EWL observations, respectively. Table 5 presents the RMS values of the positioning errors. The positioning errors for the IF model and the smoothed model were a few meters and a few centimetres, respectively. The positioning results for the long baseline were significantly worse than those for the short baseline because the unmodelled errors, such as the orbit error, residual tropospheric delay, and high-order ionospheric delay, increased with the baseline length. Additionally, the reduction in the number of observed satellites may have degraded the positioning performance. The results indicate that smoothing IF EWL observations can significantly improve the positioning accuracy from the decimetre level to the centimetre level, even for a long baseline.

Fig. 4

Number of observed satellites and corresponding DOP values for the long baseline

Fig. 5

Positioning errors of the IF ERTK model for the long baseline

Fig. 6

Positioning errors of the ionosphere-smoothed ERTK model for the long baseline

Table 5 RMS values of the positioning errors for the long baseline

Concluding remarks

All the BDS-3 satellites are currently in service, except for the last one, which is in orbit but still in testing. The BDS-3 can independently provide positioning service globally with quad-frequency signals. The quad-frequency signals benefit the high-precision positioning with regard to several aspects. They can theoretically shorten the convergence of the float solutions of PPP and RTK and improve the precision. Additionally, more EWL combinations are available owing to the quad-frequency signals, which improve not only the efficiency of AR but also the positioning performance of ERTK. It is also more efficient to estimate or eliminate the ionospheric delays with the quad-frequency signals, increasing the efficiency of AR and the positioning performance of ERTK.

ERTK with IF EWL observations can provide decimetre-level positioning service immediately after instantaneous EWL AR. If one uses the smoothed IF EWL observations, the ERTK solutions can be significantly improved and comparable to those of normal RTK. Moreover, the positioning performance is degraded with an increase in the baseline length owing to the increased unmodelled errors. To overcome this limitation, other GNSS signals can be incorporated to handle the unmodelled errors.