Characterization of the precision of PPP solutions as a function of latitude and session length

1 Introduction

Global Navigation Satellite Systems (GNSS) Precise Point Positioning (PPP; Zumberge et al., 1997) is a data processing technique that calculates absolute coordinates in the orbital reference frame with millimeter-level precision (Kouba et al., 2017) using pseudo-range and carrier phase measurements from a single, double-frequency GNSS receiver. Unlike the double differences (DD) positioning technique, PPP relies on precise satellite orbits and clock products to mitigate errors in GNSS signals (Kouba and Héroux, 2001).

In recent decades, PPP has gained popularity in surveying applications as an alternative to the relative positioning technique, especially in areas lacking geodetic infrastructure. This surge in popularity is attributed to the advantages over relative positioning because it reduces labor and logistics costs as it doesn’t require the use of local base stations, and the loss of accuracy associated with PPP is now relatively minor. Moreover, the system of equations necessary to solve for a PPP coordinate is simpler than that of a DD solution, and therefore, the computational power needed for PPP is much lower compared to DD.

Although PPP offers these advantages, there are some disadvantages related to the minimum duration of an observation session. Users typically need to collect data for at least 15 min to a few hours to resolve initial phase ambiguities and obtain a precise position solution (Kouba et al., 2017). Furthermore, although ultra-rapid and rapid satellite products exist, users requiring high-precision solutions need to wait for precise satellite orbital products in order to achieve millimeter-level precision. Hence, real-time PPP applications have limited precision compared to short- or medium-baseline (up to ∼100 km) DD solutions.

Several studies have been conducted on the understanding of PPP solutions’ scatter across various sub-daily durations relative to a 24-h session. Previous investigations focused on comparing PPP (using various software) with relative positioning, exploring ambiguity resolution effects on station parameters (position, tropospheric delay, etc.), and assessing the impact of incorporating observations from multiple GNSS constellations. For example, Ghoddousi-Fard and Dare (2006) presented an analysis of GNSS online processing services. Although their data set was rather limited, the study concluded that a 10-h session can achieve results comparable to those obtained using a 24-h observation session. Tsakiri (2008) indicated repeatable solutions at the 1–2 cm level for 24-h files and decreasing accuracy for shorter observation times. Soycan and Ata (2011) compared Global Positioning System (GPS)-only PPP solutions to Bernese software solutions achieving centimeter-level precision for durations exceeding 3 h. Gandolfi et al. (2017) conducted an analysis of PPP solutions focusing on short intervals and attained horizontal repeatability within 10 cm for datasets spanning half-an-hour durations. Barbarella et al. (2018) proposed an empirically based mathematical formula relating session duration to PPP accuracy. El-Shouny and Miky (2019) emphasized the time savings and reliable performance for sessions exceeding 2 h. Studies with small datasets have also been conducted by El-Mowafy (2011), Jamieson and Gillins (2018), Bulbul et al. (2021), and Bilgen et al. (2022).

Previous studies on the effects of session length have, in general, two main weaknesses. First, the datasets used in these analyses are small and may not be statistically representative of the actual outcomes for different durations. Also, studies with limited-size datasets cannot account for seasonal effects which can be observed if an entire year is analyzed. Second, all previous studies known to us focused on a particular latitudinal region, neglecting the effects on precision induced by GNSS orbit inclination (Swaszek et al. 2018). This study aims to address these omissions and explore the precision of geodetic-grade coordinates obtained using PPP. Our dataset comprises 365 days of 24-h GPS-only Receiver INdependent EXchange (RINEX) files for the year 2021 from 516 continuous GNSS stations across 13 networks. The study area was designed to include all possible latitude ranges from 90°N to 90°S. We achieved this wide latitude range by using data within the Americas, spanning longitudes 30°–130°W. The 24-h RINEX files were windowed into smaller, non-overlapping sessions of 1-, 2-, 3-, 4-, 6-, and 12-h which, together with the original 24-h observing sessions, were processed using the GPSPACE PPP software (https://github.com/CGS-GIS/GPSPACE). The results were then compared against each stations’ reference coordinate derived from a daily modeled position from extended linear trajectory models as described by Bevis and Brown (2014) and Bevis et al. (2019). Notably, this study boasts a vastly larger dataset compared to previous studies, consisting of over 11 million sessions. Our results show that the dispersion (or scatter) as a function of latitude and session length in the horizontal and vertical components can be modeled using a power-law function. We obtained parameters for each latitude bin for the north, east, and up components of our PPP analysis and, additionally, also provided fits for the total horizontal component. We present these results so that they can be used to estimate the expected scatter, as a proxy for observation uncertainty in a PPP session, based on its duration and latitude of the station.

2 Methods

2.1 Dataset

We gathered a one-year dataset consisting of 24-h RINEX files for the year 2021. At least one year of data was assembled to capture any seasonal effects on the PPP results. Although these effects are not directly analyzed in our study, including an entire year accounts for known quasi-periodic effects observed in GNSS time series. The 24-h RINEX files were acquired from the GNSS RINEX file archive within the Division of Geodetic Science at The Ohio State University (OSU), which continuously collects and stores permanent GNSS data from around the world for various studies.

Stations were selected within the longitudinal range of 30°–130°W, covering all latitudes from 90°N to 90°S (SI, Figure S1). This region roughly encompasses most of Canada, the lower 48 states of the United States, Mexico, Central America, and the entirety of South America. This selection provides the largest north-south land mass distribution with GNSS data.

In total, the selection comprises 516 continuous stations distributed across 13 networks (SI, Table S1 for the network distribution). Station elevations in this dataset range from −31 to 3,735 m. Figure 1 shows a histogram with the number of stations as a function of latitude.

Figure 1

Station histogram as a function of latitude used in the study.

To assess PPP solutions for short-duration observation periods, the daily RINEX files were divided into smaller time windows using the GFZRNX software package (Nischan, 2016). For every station and each day, this approach resulted in the creation of twenty-four 1-h files, twelve 2-h files, eight 3-h files, six 4-h files, two 12-h files, and one 24-h file.

3 Results

The effect of the GPS orbit configuration on dilution of precision (DOP) and, therefore, on the precision of a PPP solution has been already documented by other studies (e.g., Swaszek et al., 2018). To understand any latitude-dependent effects created by the GPS orbit configuration, we analyzed the dispersion dividing the data into latitude bins. Since we noticed, as expected, a symmetry between northern and southern latitudes, we grouped the stations based on their absolute value of latitude, with each bin separated by 10°. The size of the latitude bins was selected to maximize the resolution of the results and, at the same time, increase the number of stations per bin (so as to strengthen the statistics). We note, however, that wrapping the latitude about the geographic equator ignores some effects that are produced by ionospheric activity, which is approximately symmetric about the geomagnetic equator, which does not coincide with the geographic equator. Nevertheless, we analyzed the latitude-unwrapped results and found that our dataset did not include any significant signals that can be attributed to ionospheric activity. This might be, in part, since during the selected year of our analysis (2021) solar activity was not yet at its peak. Although the ionospheric effects are known to introduce noise in GPS solutions, isolating the ionospheric component of the scatter is outside the scope of our work.

We estimated the average standard deviation based on each station’s individual standard deviation in the east, north, and up components within each latitudinal group and as a function of duration bin, as shown in Figure 2. We also obtained an overall horizontal deviation which was computed as the norm of the east and north components. This provides a measure of the expected overall standard deviation in the horizontal component, which might be a more practical measure of uncertainty for certain users. The midpoint of each latitude range (e.g., 5°, 15°, 25°…, 85°) serves as the representative value for each group.

Figure 2

Average standard deviations of PPP solutions of continuous GPS stations grouped by duration and latitude range: (a) east component scatter for each latitude-session duration bin. Numbers in white represent the average scatter (in cm) for all stations within the bin; (b) same as (a) but for the north component; (c) same as (a) but for the horizontal component; and (d) same as (a) for the vertical component.

As expected, longer observation sessions show lower levels of scatter, which are graphically represented in Figure 2. A few interesting observations can be made from these plots, especially for short-duration sessions. In the east component (Figure 2a), the scatter decreases for increasing latitudes, while in the north (Figure 2b), the scatter peaks at mid-latitudes to ∼2.1 cm at 40°–50°. We attribute the observed behavior in the east scatter to the inability of a short PPP solution to correctly estimate the station’s clock, which is more notably manifested in the east component scatter due to the Sagnac effect (Ashby, 2002). As session duration increases, this effect is reduced but it should be noticed that the scatter is uniform across all latitudes only in the 24-h session. Interestingly, the east component shows a sudden improvement in scatter at the 60°–70° bin in Figure 2a. We attribute this sudden improvement to the additional satellites appearing over the horizon on the northern or southern end of the satellites’ sky-plots, which help constrain the station’s clock behavior for short sessions (similar to what is observed for the horizontal DOP in Figure 18 of Swaszek et al., 2018).

For the north component, we attribute the increased scatter at mid-latitudes to satellite orbit geometry which tends to create a “hole” in satellite sky plots in the north or south directions of the station (depending on which hemisphere the station is on). Nevertheless, this scatter represents about half of that in the east direction. The weak satellite coverage improves as latitude increases or decreases relative to a mid-latitude location (similar to what is observed in the east component around the 60°–70° bin). Moving towards the equator, the southern or northern satellite hole is significantly reduced. Moving towards the poles, the satellite hole tends to be closer to the zenith (or centered exactly at the zenith for stations at the geographic poles), allowing a station to observe satellites across all azimuth directions and thus improving the scatter in the north component. This increase in scatter is significantly reduced for sessions durations above 1 h. In the total horizontal component scatter (Figure 2c), it is clearly observed that the east component dominates over the north component scatter pattern.

For the vertical (Figure 2d), a surprising pattern emerges for 1 h duration sessions, where a larger scatter is observed at lower latitudes compared to high latitudes. This appears to contradict the orbit geometry (DOP) effect, previously discussed for the north component, which should worsen the vertical uncertainty at high latitudes, as described by Swaszek et al. (2018). Nevertheless, DOP estimates do not take into account that low-latitude stations are located in regions with higher water vapor content, while high-latitude regions are nearly always drier (Teng et al., 2013 – Figure 1). Water vapor is highly variable in space and time, and this heterogeneity makes the reliable estimation of delay parameters more difficult in humid areas. Figure 2d implies that the water vapor component has a stronger impact, in terms of scatter, than the geometric component.

To generalize the observed patterns across the east, north, and up components, we parameterized the average scatter shown in Figure 2. These parameterizations will allow us to directly calculate the expected scatter as a function of latitude and session duration. The next section discusses the proposed model for parameterizing the solution scatter for each latitude bin.

3.1 Parameterization

We introduce a model to describe the relationship between the precision of a PPP solution (σ) and the duration of the observation session (t), of the type,

(2) σ = a · t b ,

where parameters a and b are determined through a least squares optimization. To simplify the least squares procedure, we substitute α = e ^a to linearize equation (2) as

(3) ln ( σ ) = α + ln ( t ) · b .

This allowed us to employ the classical Gauss-Markov model to estimate the vector of unknown parameters, α and b. We note that while we can assign units to a ( cm / h b with b being dimensionless), this parameter does not represent a particular process other than the relationship between time and scatter. Hence, because these fits are purely kinematic, we don’t assign specific units to a and directly assume that the output, σ, is in cm.

We performed the adjustment for each latitude-duration pair bin, as shown in Table 1, where we report the value of a instead of α. The observed scatter and their fits are shown for the horizontal in Figure 3 and for the vertical in Figure 4 (the east and north component fits can be found in the SI).

Table 1

The estimated parameters of the proposed function

|Latmin|	|Latmax|	East		North		Horizontal		Vertical
		a	b	a	b	a	b	a	b
80	90	0.015	−0.598	0.017	−0.579	0.014	−0.601	0.035	−0.576
70	80	0.016	−0.638	0.018	−0.664	0.015	−0.688	0.032	−0.562
60	70	0.024	−0.631	0.016	−0.625	0.018	−0.656	0.032	−0.540
50	60	0.039	−0.745	0.016	−0.593	0.027	−0.746	0.038	−0.523
40	50	0.049	−0.924	0.018	−0.720	0.034	−0.919	0.050	−0.715
30	40	0.053	−0.913	0.017	−0.717	0.036	−0.937	0.059	−0.738
20	30	0.050	−0.857	0.014	−0.721	0.034	−0.882	0.068	−0.769
10	20	0.049	−0.793	0.014	−0.673	0.033	−0.801	0.064	−0.731
0	10	0.052	−0.809	0.014	−0.684	0.034	−0.817	0.059	−0.716

Figure 3

The observed scatter and power-law fits for the horizontal component. Each panel shows the central absolute latitude range with session duration in the x-axis (in h, log scale) and scatter in the y-axis (in cm, log scale). Blue line represents the data from Figure 2d and orange line the power-law fit.

Figure 4

The observed scatter and power-law fits for the vertical component. Each panel shows the central absolute latitude range with session duration in the x-axis (in h, log scale) and scatter in the y-axis (in cm, log scale). Blue line represents the data from Figure 2d and orange line the power-law fit.

For clarity, the data were linearly interpolated in intervals of 30 min for occupation duration and one degree for latitude using the estimated parameters for each latitude bin. Contours at specific precisions were then overlaid onto the interpolated data. The final contour plot for all components being analyzed is shown in Figure 5.

Figure 5

Expected precision of a PPP solution based on the duration and latitude of an observation as a continuous surface: (a) east, (b) north, (c) horizontal, and (d) vertical.

4 Discussion

In this work, we estimated PPP solution uncertainties as a function of session duration. By including data that spans over the entire global latitude range, we also provide estimates of uncertainty as a function of station latitude. Furthermore, we adjusted a series of power-law parameters to each latitude bin (as a function of duration) to obtain a continuous estimate of uncertainty based on latitude and duration. Unsurprisingly, longer observation sessions reduce positioning errors, as shown in many previous studies. Our work, however, presents results that were obtained with several orders of magnitude more data than previous works. Our uncertainty estimates also include data for an entire year, which makes them more robust against potential seasonal effects in the northern and southern hemispheres. Therefore, our presented method provides an overall more robust estimate of PPP uncertainty as a function of the two main variables that affect the quality of a solution: session duration and latitude.

We also discussed the observed effects in precision as a function of latitude for a short session duration. We noted that positioning errors in the vertical tend to decrease with distance from the equator, as opposed to the DOP value which increases with latitude. This is somewhat surprising given that favorable conditions for satellite geometry decrease in this direction. We argue that this behavior is a result of the latitude dependence of total column atmospheric water vapor (Teng et al., 2013), which is known to affect the quality of GNSS solutions due to its high temporal variability. Further analysis correlating satellite availability, dilution of precision, and uncertainty as a function of time of the year could provide deeper insights into this and other effects. Given that results are GPS-only, potential increases in precision could be made by incorporating other constellations and integer ambiguity resolution PPP techniques which were not considered in the present work.

Because we used GPS observations from geodetic reference networks, the sky view at nearly every station was very good to excellent. Nevertheless, surveyors are often required to position points with much poorer sky view due to the presence of buildings, trees, or major topography. So, our model relating attainable positioning accuracy to station latitude and the time span of observations will provide only a lower bound on positioning error for stations with relatively poor sky view.