Abstract
A ground-based augmentation system (GBAS) is a critical component in civil aviation that augments the Global Positioning System (GPS) in providing precision approach and landing capabilities with guaranteed accuracy and integrity. The GBAS ground facility broadcasts a parameter known as σvig to the aircraft, which is used to compute vertical protection levels for evaluating navigation integrity. σvig represents the standard deviation of the vertical ionospheric gradients, which bounds the spatial gradients under nominal conditions. Although the time-step method has been widely utilized to estimate ionospheric spatial gradients, this strategy suffers from temporal effects. In this paper, an improved time-step method is developed for separating temporal gradients from spatial gradients using observation data collected from the Hong Kong Satellite Positioning Reference Station Network. We investigated two parameters: σtg, which bounds the standard deviation of temporal gradients, and σvig. The results show that a constant value of 5.5 mm/km can serve as an upper bound for all σtg values. However, the results of σvig vary seasonally, with maximum and minimum values occurring at the equinoxes and summer, respectively. To reflect this seasonality, quadratic polynomial expressions, given as functions of the day of the year, were derived to provide an upper bound for all σvig values.
1 INTRODUCTION
The Global Positioning System (GPS) has been widely applied to various fields, including agriculture, construction, mining, and transportation. In civil aviation, GPS is combined with ground-based augmentation systems (GBASs) to provide enhanced levels of service that support aircraft precision approach and landing operations. A GBAS broadcasts differential corrections and integrity-related information to aircraft via a very-high-frequency transmitter to ensure accuracy and aviation integrity during aircraft approach. Although ionospheric errors are largely removed by differential corrections, a residual error remains because of spatial and temporal decorrelations between the GBAS reference station and the aircraft. As the larger of the two correlations (Lee et al., 2007), the spatial decorrelation, which is expressed in terms of ionospheric spatial gradients with units of mm/km, may degrade GBAS performance. Ionospheric spatial gradients vary significantly between nominal and anomalous ionospheric conditions. For example, in the conterminous United States (CONUS), the nominal ionospheric spatial gradient value is approximately 4 mm/km (Lee et al., 2007); in contrast, during anomalous conditions, ionospheric spatial gradients larger than 400 mm/km have been observed in the CONUS (Pullen et al., 2009).
To protect users against anomalous ionospheric conditions, ionospheric monitors and the geometry screening method are implemented in GBASs (ICAO, 2018; Lee et al., 2011; Pullen et al., 2009). Under nominal ionospheric conditions (RTCA DO-245A, 2008), spatial gradients are considered in the computation of protection levels, which bound the true position error. The vertical protection level (VPL) for a nominal or fault-free case is computed by the following (RTCA DO-253C, 2008):
1
where Kf fmd is a multiplier determined by the probability of fault-free missed detection (on the order of 10−7), which covers the entire integrity risk due to nominal or fault-free (known as H0 in GBAS) conditions (RTCA DO-245A, 2008). svert, i is the vertical component of the weighted least-squares projection matrix of satellite i, and N is the number of satellites used in the position solution. σi represents the pseudorange standard deviation of satellite i under the fault-free hypothesis, which is computed as follows (RTCA DO-253C, 2008):
2
where 
and 
indicate the GBAS ground subsystem and airborne error contributions, respectively. σtropo,i is the uncertainty associated with the tropospheric error. σiono,i is the residual ionospheric delay uncertainty (RTCA DO-253C, 2008):
3
where Fpp is the vertical-to-slant obliquity factor. xair is the distance between the aircraft location and the GBAS reference point. τ is the time constant used in the smoothing filter, and vair is the horizontal approach velocity of the aircraft. σvig is the nominal standard deviation of vertical ionospheric gradients, which must bound the ionospheric spatial gradients under nominal conditions. Equations (1)–(3) show that σvig is crucial in calculating the VPL.
For mid-latitude regions, the value of σvig has been determined to be several mm/km (Lee et al., 2007; Mayer et al., 2009). However, in low-latitude regions, the ionosphere is more active. Therefore, the value of 4 mm/km is not sufficient to bound nominal ionospheric conditions in low-latitude regions. For example, in Japan, σvig was estimated to be 6.02 mm/km based on observation data from a dense global navigation satellite system (GNSS) network called GEONET (Yoshihara et al., 2010). In India, an overbounding σvig of 14.4 mm/km is necessary under active ionospheric conditions (Reddy & Achanta, 2016). An overall bounding σvig value of 14 mm/km is needed to bound nominal ionospheric behavior during the daytime in Brazil (Chang et al., 2021). Because the ionospheric layer presents different behaviors depending on location, it is crucial to examine the ionospheric nominal conditions in Hong Kong, which is located at low latitudes between 22.13° and 22.58° north (the geomagnetic latitude of Hong Kong is between 12.65° north and 13.49° north, whereas the geomagnetic latitude for India is between 1° south and 24° north), to determine an appropriate value of σvig for the deployment of GBASs in this area.
Because ionospheric delays can be estimated by dual-frequency observations, three methods have been widely applied to estimate ionospheric spatial gradients: the dual-frequency station-pair method (Lee et al., 2007), the satellite-pair method (Supriadi et al., 2022), and the time-step method (Chang et al., 2021). The dual-frequency station-pair method utilizes the differential ionospheric delay between a pair of stations with the same satellite at the same epoch to calculate the spatial gradients. This method is intuitive because its configuration resembles the GBAS ground station and aircraft configuration. However, the remaining inter-frequency bias (IFB) in receivers, if not estimated accurately, affects the spatial gradient estimation results. A single-frequency station-pair method in which the IFB is mitigated based on single-frequency observations was proposed to estimate the ionospheric spatial gradient (Budtho et al., 2018; Saito et al., 2012). However, both the single-frequency and dual-frequency station-pair methods require a dense station network with short distances (typically less than 40 km) for measuring ionospheric spatial gradients. The satellite-pair method takes the differential ionospheric delay between a pair of satellites with the same receiver at the same epoch. This method cancels out the IFB of the receiver, and the IFB of the satellite can be removed by ionosphere exchange (IONEX) data collected from the International GNSS Service (IGS) (Supriadi et al., 2022). The separation distance between two ionospheric pierce points (IPPs) is generally larger than 100 km, and residual errors can be largely reduced after being divided by the large IPP separation distance (Supriadi et al., 2022). However, because the IPP separation distance in the satellite-pair method is much larger than the GBAS-applicable distance (10–40 km), this method is not suitable for representing spatial gradients in GBAS applications.
The dual-frequency time-step method estimates spatial gradients by comparing the ionospheric delay of a single satellite and a single receiver between different epochs. This method is advantageous because both the receiver and satellite IFBs are eliminated, and the IPP distance can be regulated by adjusting the time interval. Use of this method has also been proposed to detect threatening spatial gradients in real time (Marini-Pereira et al., 2022, 2023). In addition, a single-frequency time-step method was proposed with desired ionospheric gradients estimated by Kalman filtering and ambiguity fixed by the least-squares ambiguity decorrelation adjustment (LAMBDA) method (Budtho et al., 2020). However, ionospheric gradients estimated by either the single-frequency or dual-frequency time-step method are computed over a time interval instead of at a single epoch, thereby including both temporal and spatial gradients (Chang et al., 2021). The inclusion of temporal gradients might result in an overestimation or underestimation of the desired spatial gradients. A “geometric similarity” method has been proposed for evaluating temporal gradients using data collected from dense GNSS observation networks, e.g., GEONET, which comprises more than 1300 reference receivers and an average receiver separation of approximately 20 km (Chang et al., 2021). A large, dense network such as GEONET (more than 1300 stations, with coverage over all of Japan’s islands) is required to meet the IPP selection criterion, but is difficult to apply in Hong Kong, where the coverage radius of GNSS reference stations is less than 25 km. Therefore, in this study, we propose an improved time-step method in which ionospheric temporal gradients can be separated from spatial gradients based on data collected from the Hong Kong Satellite Positioning Reference Station Network (SatRef).
In this paper, the nominal ionospheric conditions in Hong Kong are examined, and the bounding values of ionospheric spatial and temporal gradients are determined by using the improved time-step method. In Section 2, an overview of the existing methods for estimating ionospheric spatial gradients is provided, and corresponding overbounding σvig results are compared. In Section 3, the improved time-step method, in which the temporal gradients are separated from the spatial gradients, is described. In Section 4, the overbounding results of σtg and σvig for each day from 2012 to 2020 are presented, and quadratic polynomial expressions are developed to serve as an upper bound for all σvig values with seasonality. Finally, in Section 5, conclusions and final remarks are provided.
2 IONOSPHERIC GRADIENT ESTIMATION METHOD
The dual-frequency station-pair, satellite-pair, and time-step methods are based on the ionospheric delay estimated by dual-frequency measurements. The slant ionospheric delay can be estimated by dual-frequency code or carrier-phase observations expressed as follows (Datta-Barua et al., 2010):
4
5
where 
, with fL1 and fL2 denoting the L1 and L2 frequencies, respectively. c is the speed of light. ϕ and ρ are the carrier-phase and code observations with contained noise indicated by εϕ and ερ. IFBk and IFBi indicate the IFB biases for receiver k and satellite i. NL1 and NL2 are integer ambiguities in the frequencies of the carrier-phase observations at L1 and L2, respectively. Because the carrier-phase observations have much less multipath and thermal noise, Iϕ is less noisy than Iρ. Iϕ with the unknown ambiguities NL1 and NL2 is leveled by Iρ (Kim et al., 2015):
6
where Îρ(t) is the leveled ionospheric delay estimation, t is the epoch index, and M is the smoothing weight. M = t when t < Ns, and M = Ns when t ≥ Ns. Ns is the maximum number of samples within the time at which the filter reaches a steady state (Marini-Pereira et al., 2022). The slant ionospheric delay varies with satellite elevation. In particular, the GPS signals of low-elevation satellites experience a larger delay than those of high-elevation satellites. To compute the equivalent vertical delay at the IPP, the slant ionospheric delay is divided by the obliquity factor Ob, assuming that the ionosphere is a thin shell. Ob is defined as follows (Oliveira et al., 2020):
7
where el is the elevation angle between the line of sight (LOS) and the local horizon. Re is the radius of the Earth. hl is the height of the thin-shell model, which is assumed to be 350 km.
GNSS data are collected from all SatRef stations. SatRef was established in 2001 with 6 stations (Kumar et al., 2016) and now has 18 stations, as shown in Figure 1. The approximate locations of these 18 GNSS reference stations are labeled as red points on the GeoInfo map obtained from https://www.map.gov.hk. The physical separation between any two stations ranges between 4.8 km and 50 km. Considering the purpose of this study, which is to investigate and evaluate ionospheric gradients under nominal conditions, the data sets during periods with a planetary K (Kp) of less than 3 and an L1-band amplitude scintillation index S4 of less than 0.2 were selected (Chang et al., 2021). The geomagnetic 3-h Kp index is an important measure representing the energy input from the solar wind to earth, expressed as an integer in the range of 0–9, with 1 being calm and 5 or more indicating a geomagnetic storm (Bartels et al., 1939). The Kp index is widely used to characterize ionospheric disturbance on a global scale (Matzka et al., 2021). In addition, the S4 index is a metric for indicating the intensity of amplitude scintillation (Van Dierendonck et al., 1993). Amplitude scintillation events frequently occur in the equatorial and low-latitude regions, mainly due to equatorial ionization anomalies (EIAs) and ionospheric bubbles formed after sunset (de Oliveira Nascimento Brassarote et al., 2018). The S4 index can be utilized to identify local ionospheric distortions for low latitudes.
Distribution of SatRef stations in 2023
For the sake of simplicity, the ionospheric delay in the following text refers specifically to the vertical ionospheric delay. The ionospheric delay estimate is divided by the IPP separation distance to compute the ionospheric spatial gradients. Based on different configurations of stations and satellites, three methods, i.e., the dual-frequency station-pair method, satellite-pair method, and time-step method, are designed for spatial gradient estimation. In the following subsections, a brief overview of these methods is provided, and the corresponding ionospheric spatial gradient results are presented.
2.1 Dual-Frequency Station-Pair Method
The dual-frequency station-pair method has been used to estimate ionospheric spatial gradients under nominal and anomalous ionospheric conditions (Kim et al., 2015; Lee et al., 2007; Saito et al., 2017). The satellite and station configuration of the dual-frequency station-pair method is illustrated in Figure 2.
Satellite and receiver configuration of the dual-frequency station-pair method
As shown in Figure 2, two IPPs (IPP1 and IPP2) indicate the intersections of the LOS of two stations (R1 and R2) observing the same satellite S1 and the thin-shell ionosphere. The ionospheric delay at IPP1 is compared with the delay at IPP2 at the same epoch t. The ionospheric spatial gradient at t is then computed by dividing the delay difference by the IPP separation distance:
8
where 
and 
are the ionospheric delays at IPP2 and IPP1, respectively. It should be noted that the estimated gradients are values projected to the baseline between a pair of IPPs. Thus, the gradient is generally underestimated because the gradient is a two-dimensional (2D) vector. Estimated gradients should accurately represent the true ionospheric anomaly in all directions, considering that the direction between a GBAS ground station and aircraft can take any orientation. Relying on a single pair of stations may lead to an underestimated 2D gradient magnitude. Therefore, utilizing multiple pairs of stations is essential to address this issue.
The dual-frequency station-pair method involves calculating the ionospheric delay difference between two stations viewing the same satellite, which cancels out the satellite IFB. The IFBs of the two receivers can be estimated and removed with a simpler and faster method by searching for the best bias candidate that provides the minimum derivation of the total electron content (TEC) (Ma & Maruyama, 2003). However, the error in the estimated receiver IFB may degrade the accuracy of the spatial gradient estimation results. The separation distance of the IPP is determined by the physical separation of the stations. Therefore, to obtain ionospheric spatial gradients at GBAS-applicable distances (10–40 km), a dense GNSS network with appropriate separation distances is needed. The sparse distribution of the GNSS reference stations in some areas limits the application of the dual-frequency station-pair method for estimating ionospheric spatial gradients over short distances.
Based on the dual-frequency station-pair method, the ionospheric spatial gradient results from the day of year (DOY) 103 (April 13) in 2014 are examined as an example, as shown in Figure 3. The differential ionospheric delay as a function of the IPP separation distance is shown in the left panel. The range of the IPP separation distance is 2–47 km, which is similar to the range of physical separation between SatRef stations. The vertical ionospheric delay difference is divided into bins, where the color scale of each pixel indicates the number of samples in each bin. With the ionospheric delay difference divided by the corresponding IPP separation distance, the ionospheric spatial gradient results can be obtained. The collected ionospheric spatial gradient results are normalized by subtracting the mean (μvig) from the results and dividing the results by the standard deviation (σvig). The folded cumulative distribution function (CDF) of the normalized ionospheric spatial gradient results is indicated by a black dashed line (computed by 1 minus the CDF for x-axis values larger than 0) in the right panel of Figure 3. A nominal Gaussian distribution with a mean of zero and standard deviation of σvig, indicated by the red line, is not able to overbound the non-Gaussian tail of the actual distribution of the normalized gradient results (Lee et al., 2007). Therefore, the nominal Gaussian must be inflated, and the required inflation factor (f) is obtained when the inflated nominal Gaussian distribution covers all of the actual distribution at a probability of 10−4 and above (Yoshihara et al., 2010). Finally, the overbounding σvig value is computed as |μvig| + f × σvig. In this example, μvig, σvig, and f are −0.05 mm/km, 14.41 mm/km, and 2.83, respectively. The obtained overbounding σvig value computed by |μvig| + f × σvig is 32.33 mm/km.
Dual-frequency station-pair results for April 13, 2014
The obtained σvig of 32.33 mm/km is much larger than the results obtained in previous studies, such as those by Chang et al. (2021), Supriadi et al. (2022), and Yoshihara et al. (2010). Thus, this value cannot represent real ionospheric nominal conditions. To explore the cause of this discrepancy, data from another day (DOY 257, i.e., September 14) in 2014 were used to obtain the dual-frequency station-pair σvig results shown in Table 1 for comparison. The number of samples on DOY 257 is larger because more stations were in operation. In addition, the results are compared for different IPP separation distance ranges. Table 1 shows that the magnitudes of σvig and f decrease sharply as the IPP separation distance increases. This trend indicates that the large σvig value is caused by the remaining uncalibrated receiver IFB and carrier-phase leveling error, which dominate the ionospheric gradient estimation results. Consequently, a large σvig with a large f is needed to bound the non-Gaussian tail of the actual distribution. Increasing the IPP separation distance can reduce the overbounding σvig value because the error effect is mitigated when the delay difference is divided by longer IPP distances. However, the number of samples decreases as the required IPP separation distance increases, and the number of samples might be insufficient to overbound σvig at a probability of 10−4. For example, as shown in the fourth row of Table 1, the sample number decreases to below 105 when the IPP separation distance is set to 30–47 km. Therefore, the dual-frequency station-pair method is not suitable for obtaining the desired overbounding σvig value in Hong Kong, as most station separations in SatRef are less than 30 km.
Overbounding σvig Results Estimated by the Dual-Frequency Station-Pair Method
2.2 Satellite-Pair Method
The satellite-pair method is another approach for estimating ionospheric spatial gradients (Supriadi et al., 2022). The satellite and station configuration of the satellite-pair method is illustrated in Figure 4.
Satellite and receiver configuration of the satellite-pair method
As shown in Figure 4, IPP1 and IPP2 represent the locations at which the LOS of a pair of satellites (S1 and S2) observed by the same receiver R1 intersect with the thin-shell ionosphere. The ionospheric delay difference is computed by comparing the ionospheric delays at IPP1 and IPP2 at the same epoch. The delay difference is then divided by the IPP separation distance to compute the ionospheric spatial gradient, as indicated by Equation (8). The satellite-pair method cancels out the receiver IFB by taking the ionospheric delay difference between a pair of satellites observed by the same receiver. The remaining satellite IFB can be removed by using the IFB estimation results calculated by the GNSS network of the IGS. Although the satellite IFB corrected by IONEX data may contain residual errors, it is estimated by many receivers around the world, thereby having a higher accuracy than the receiver IFB estimate based on local receivers (Supriadi et al., 2022).
For the satellite-pair method, the ionospheric spatial gradient results from DOY 103 (April 13) in 2014 are used as an example, as shown in Figure 5. The differential ionospheric delay as a function of the IPP separation distance and normalized ionospheric spatial gradients are shown on the left and right panels of Figure 5, respectively. The IPP separation distance ranges from 76 km to 1065 km, which is much larger than the IPP distance range of the dual-frequency station-pair method. The folded CDF of the normalized ionospheric spatial gradients of the satellite-pair method is thinner than that of the dual-frequency station-pair method, thereby requiring a smaller f to bound the non-Gaussian tails. In this example, the overbounding σvig value obtained by the satellite-pair method is 10.48 mm/km.
Satellite-pair results for April 13, 2014
The overbounding σvig results for another day (DOY 257, i.e., September 14) in 2014 are shown in Table 2. The obtained overbounding σvig value is smaller than that obtained via the dual-frequency station-pair method. This result is due to the elimination of the receiver and satellite IFBs, and the error effect is largely reduced by the long IPP separation distance. However, such a long IPP separation distance (larger than 100 km) is not suitable for representing ionospheric spatial gradients at the range of 10–40 km relevant to GBAS applications.
Overbounding σvig Results Estimated by the Satellite-Pair Method
2.3 Time-Step Method
The dual-frequency time-step method is an alternative method used for estimating ionospheric gradients on nominal days (Lee et al., 2007) and for validating anomalous ionospheric gradients when the density of ground stations is insufficient (Yoon, Kim et al., 2017; Yoon, Lee et al., 2017). The satellite and station configuration of the time-step method is illustrated in Figure 6.
Satellite and receiver configuration of the time-step method
As depicted in Figure 6, IPP1 and IPP2 indicate the locations at which the LOS of a single satellite S1 observed by a single receiver R1 at one epoch (t) and a later epoch (t + ΔT) intersect with the thin-shell ionosphere. The ionospheric delay difference is computed by comparing the ionospheric delays at IPP1 and IPP2. The delay difference is then divided by the IPP distance to compute the ionospheric gradient:
9
where 
and 
are the vertical ionospheric delays of IPP2 and IPP1, respectively. Similar to the dual-frequency station-pair method, the estimated gradient is constrained to the baseline direction and ignores the 2D nature of the gradient.
Because the gradient g is estimated over a time interval, the result contains both the ionospheric spatial gradient and the temporal gradient. In other words, the time-step method may overestimate or underestimate the desired results when applied to estimate ionospheric spatial gradients (Chang et al., 2021). Despite this limitation, the time-step method has several advantages. First, the receiver and satellite IFBs are eliminated when the delays of the same satellite and receiver are differentiated. Second, assuming that no cycle slip occurs during ΔT, the integer ambiguity can also be removed; therefore, Îρ can be replaced by Iϕ in Equation (9) to further reduce the noise effect in the gradient estimation. Third, the IPP separation distance is determined by the satellite elevation and ΔT. Therefore, the desired IPP separation distance can be obtained by adjusting ΔT.
As shown in Figure 7, the IPP separation distance varies from 3 km to 80 km as ΔT increases from 1 min to 10 min. For example, 90% of samples have an IPP separation distance of 15–30 km when ΔT is chosen as 5 min. However, one limitation of the time-step method is that the IPP separation direction is along the satellite track determined by the satellite movement. In Figure 8, a histogram of the IPP separation directions of all samples on April 13 is presented, where IPP direction azimuths larger than 180° are converted to be within 0°–180°. The results show that the direction of IPP separation primarily lies in the north-south direction, indicating that the estimated gradients are dominant in the north-south direction and result in a less accurate gradient estimation in the east-west direction. This result presents one of the drawbacks of the time-step method, which is unable to describe 2D ionospheric gradients (Budtho et al., 2020).
Box plot of the IPP separation distance for different ΔT values
Direction of the IPP separation on April 13 obtained using the time-step method
For the time-step method, the ionospheric spatial gradient results from DOY 103 (April 13) in 2014 are used as an example, as shown in Figure 9. The differential ionospheric delay as a function of the IPP separation distance and normalized ionospheric spatial gradients are shown in the left and right panels of Figure 9, respectively. ΔT is chosen as 1 min in this example, and the corresponding IPP separation distance ranges from 3 km to 9 km. In this example, the resultant overbounding σvig value obtained using the time-step method is 13.99 mm/km.
Time-step results for April 13, 2014
The overbounding σvig results of the time-step method for another day (DOY 257, i.e., September 14) in 2014 are shown in Table 3. Different ΔT values with corresponding IPP separation distance ranges are considered. A longer IPP separation distance has little impact on the overbounding σvig value because the residual error is negligible, as only carrier-phase measurements are used. It should be noted that small-scale or medium-scale fluctuations of the ionosphere can cause TEC variations and affect the ionospheric gradient. TEC perturbations (usually enhancements) may be caused by sporadic E layers or medium-scale traveling ionospheric disturbances (MSTIDs). The differences between sporadic E layers and MSTIDs lie in their spatial scales and durations. The spatial scale and duration of sporadic E layers are smaller and shorter, respectively, than those of MSTIDs in the F region (Saito et al., 2021; Maeda and Heki, 2014; Saito et al., 1998). In addition to sporadic E layers and MSTIDs, equatorial plasma bubbles are another common ionosphere disturbance source in low latitudes such as Hong Kong (Sousasantos et al., 2021). Therefore, the magnitudes of the spatial and temporal gradients can vary when these disturbances are present. The time-step method offers several advantages over the dual-frequency station-pair and satellite-pair methods. First, this method quantifies a single satellite and a single station as a pair and obtains targeted IPP separation distances by adjusting ΔT; for example, a ΔT value of 5 min corresponds to IPP separation distances ranging from 14.8 km to 40.8 km, which falls within the range of GBAS-applicable distances. Second, the time-step method eliminates both satellite and receiver IFBs, resulting in reduced noise and errors in ionospheric gradient estimates. The noise effect is further reduced because only carrier-phase observations are utilized. Therefore, the time-step method can provide a good estimate of ionospheric gradients with sufficient data within the IPP separation distances of interest. However, it is important to note that the gradient estimated by the time-step method is a combination of the spatial gradient and temporal effects, whereas the value of σvig reflects only the ionospheric spatial gradient. Therefore, it is helpful to separate the temporal effects from the spatial gradients to obtain a more accurate overbounding σvig value.
Overbounding σvig Results Estimated by the Time-Step Method
3 IMPROVED TIME-STEP METHOD
As mentioned previously, the time-step method involves comparing ionospheric delay differences between two epochs. This comparison can be divided into two parts: the temporal delay difference and the spatial delay difference (Chang et al., 2021):
10
Here, the first part, i.e., 
, represents the temporal delay difference and can be computed by the product of the temporal gradient and the time interval. The second part, 
, represents the spatial delay difference and can be computed by the product of the spatial gradient and the IPP separation distance. Therefore, the ionospheric delay difference can be rewritten as follows:
11
where gte,t and gsp represent the ionospheric temporal and spatial gradients, respectively. 
is the residual error due to noise and multipath. Dividing both sides by the IPP separation distance modifies Equation (11) as follows:
12
where the left side of the equation represents the time-step ionospheric gradients. The first term indicates the temporal effect caused by gte,t. For simplicity, this term is represented by gte. The units of gsp and gte are mm/km. The last term, i.e., 
divided by the IPP separation distance, is negligible compared with gsp and gte, as only carrier-phase measurements are used and the separation distance is chosen to be larger than 10 km.
Figure 10 shows time-step ionospheric gradients for the HKNP, HKSS, and HKWS SatRef stations between 3 p.m. and 8 p.m. local time (LT) on April 13, 2014. Because the separation distances between stations are sufficiently smaller than the station–satellite distances, the IPP separation distances for these three stations are the same as those shown in the middle panel. The time-step ionospheric delay differences and ionospheric gradients for these three stations show similar trends. However, although the fluctuations in the ionospheric gradient curve are the same for HKSS and HKWS, these fluctuations are not the same for HKNP and HKSS (or HKNP and HKWS). Because the spatial gradients are larger and dominate the time-step ionospheric gradients, the trends and fluctuations of the time-step ionospheric gradient curves indicate the spatial gradients and temporal gradients, respectively. The same trend is observed for HKSS, HKWS, and HKNP, but the fluctuations in these stations are not similar. This result indicates that all SatRef stations share the same ionospheric spatial gradient, whereas the ionospheric temporal gradients are the same only when two stations are close to each other. In this case, the HKSS and HKWS stations have a physical separation distance of 6.8 km, and the HKSS and HKNP stations have a physical separation distance of 43 km. The spatial gradient similarity between the SatRef stations is further illustrated in Figure 11, with HKNP-HKKS and HKNP-HKTK pairs (HKTK and HKKS have a physical separation distance of 22 km) shown as an example.
Time-step ionospheric gradient on April 13, 2014 for GPS pseudorandom noise (PRN) 5
Ionospheric spatial gradient obtained by the dual-frequency station-pair method on April 13, 2014 for GPS PRN 5
To further demonstrate that the fluctuations in the time-step ionospheric gradient curve are caused by the temporal gradient, the ionospheric spatial gradient computed by the dual-frequency station-pair method are compared in Figure 12. The station pair of HKNP-HKSS is chosen because the IPP separation distance of this station pair is larger than 40 km, thereby largely reducing the residual error and multipath effect. Because only carrier-phase measurements are used in the time-step ionospheric gradient computation and the IPP distance is 15–40 km during this period, the residual error and multipath in the time-step ionospheric gradients are negligible. Specifically, the largest noise effect in the time-step ionospheric gradients is smaller than 0.41 mm/km (computed by 
, considering a 2-mm noise level in carrier-phase measurements, independence between L1 and L2 measurements, independence between different epochs, and an IPP separation distance of 15 km). However, as shown in Figure 12, the fluctuations in the time-step ionospheric gradients and ionospheric gradient difference between adjacent epochs are much larger than those in the dual-frequency station-pair spatial gradients, with the gradient fluctuations exceeding 5 mm/km. Therefore, the fluctuations in the time-step ionospheric gradients are caused by temporal gradients instead of noise effects and must be separated from the ionospheric spatial gradients.
Comparison between the time-step ionospheric gradient of HKNP and the dual-frequency station-pair ionospheric spatial gradient of HKNP-HKSS on April 13, 2014 for GPS PRN5: (left) ionospheric gradient estimated by the time-step and station-pair methods; (right) ionospheric gradient change rates of the time-step and station-pair methods, calculated from the ionospheric gradient difference between adjacent epochs
The time-step method in which the ionospheric spatial and temporal gradients are separated is referred to as the improved time-step method. The separation process is based on two assumptions: (1) Overall, the temporal gradients are smaller than the spatial gradients. (2) The spatial gradients do not change as rapidly over time as the temporal gradients. The first assumption has been widely accepted in similar studies (Chang et al., 2021; Lee et al., 2007). The second assumption is in accordance with nominal ionospheric conditions, as shown in Figure 12. The spatial gradient estimated by smoothing the time-step gradient differs in magnitude from the spatial gradient estimated by the station-pair method. This difference is due to the different choice of IPP pairs between the time-step method and the station-pair method (Budtho et al., 2020). A simple way to separate the temporal gradient from the spatial gradients is to smooth the time-step ionospheric gradient. The smoothed time-step ionospheric gradient is the desired spatial gradient, and the temporal ionospheric gradient is then computed by subtracting the ionospheric spatial gradient from the time-step ionospheric gradient. Locally estimated scatterplot smoothing (LOESS) is used herein, with a span of 10% of the total number of data points to smooth the nonlinear time-step ionospheric gradient. The separation results for HKNP on April 13, 2014 are shown in Figure 13. After the ionospheric spatial gradient has been separated from the temporal gradient, the difference between the ionospheric spatial gradients is significantly reduced and is similar in magnitude to that of the dual-frequency station-pair method.
Separation of the spatial and temporal ionospheric gradients in the time-step ionospheric gradient of HKNP on April 13, 2014 for GPS PRN 5: (left) time-step ionospheric gradient and spatial and temporal gradients separated by the improved time-step method; (right) ionospheric gradient change rates among the station-pair, time-step, and improved time-step methods
The change rate is calculated as the ionospheric gradient difference between adjacent epochs divided by the time interval.
For the improved time-step method with a ΔT value of 5 min, the folded CDFs of the normalized ionospheric spatial and temporal gradients are shown in Figure 14 on a logarithmic scale. For the ionospheric spatial gradient, the obtained σvig, μvig, and f are 9.94 mm/km, −0.27 mm/km, and 1.37, respectively. The overbounding σvig value computed by |μvig| + f × σvig is 13.78 mm/km. For the ionospheric temporal gradient, the obtained standard deviation (σtg), mean (μtg), and f are 0.59 mm/km, 0.0014 mm/km, and 1.82, respectively. The overbounding σtg value computed by |μtg| + f × σtg is 1.07 mm/km.
Improved time-step results for April 13, 2014: (left) overbounding results for σvig; (right) overbounding results for σtg
The overbounding σvig results for DOYs 103 and 257 in 2014 are compared between the dual-frequency station-pair method, time-step method (ΔT = 5 min), and improved time-step method (ΔT = 5 min) in Table 4. On DOY 103, the IPP separation distance range is chosen as 20–40 km because not enough samples can be obtained when the range is set as 30–40 km. The magnitude of the overbounding σvig value of the improved time-step method is the smallest among these three methods. This result is attributed to the usage of carrier-phase measurements and the separation of the spatial and temporal gradients. Compared with the time-step method, the magnitude of the overbounding σvig value of the improved time-step method is closer to that of the dual-frequency station-pair method for DOY 257. Therefore, the improved time-step method outperforms the time-step method and is more suitable for estimating σvig when the dual-frequency station-pair method is not able to obtain a sufficient number of samples.
Comparison of Overbounding σvig Results Between the Dual-Frequency Station-Pair, Time-Step, and Improved Time-Step Methods
4 VARIABILITY OF IONOSPHERIC SPATIAL AND TEMPORAL GRADIENTS
Using the improved time-step method, both σvig and σtg can be obtained. σvig and σtg represent the standard deviation of the spatial and temporal variations of ionospheric delay, respectively. Both parameters provide information about ionospheric conditions under nominal conditions. In addition, the temporal gradient affects the positioning accuracy because the smoothing filter in the GBAS introduces an additional error caused by the time-varying effect. The additional error can be quantified by substituting the magnitude of σtg into the simulation or equation of the smoothing filter (Huang et al., 2008). To improve the GBAS performance by reducing the ionospheric impact on the smoothing filter, divergence-free (DFree) and ionosphere-free (IFree) techniques have been proposed based on dual-frequency measurements (McGraw & Young, 2005). DFree techniques can remove the time-varying delay effect of the smoothing filter, and IFree techniques can remove the first-order ionosphere effect (Konno, 2007). The DFree improvement can be quantified by comparing the DFree smoothing filter with the current smoothing filter while considering σtg.
Overbounding σvig and σtg results for 2014–2020 as a function of the DOY are shown in Figure 15 and Figure 16, where μvig and σvig are indicated by blue and red lines, respectively. Because solar activity plays a key role in ionospheric activity, the ionospheric spatial gradients and σvig tend to be larger around the peak of a solar cycle, i.e., 2012–2015, than in other years. In addition, because of the pronounced EIA effect during spring and autumn equinoxes (Chen et al., 2016), the overbounding σvig value tends to be greater during the March and September equinoxes than during other seasons. In contrast, the overbounding σtg results for 2014–2020 in Figure 16 do not show an apparent seasonal dependence.
Daily σvig overbound results for 2014–2020
Daily σtg overbound results for 2014–2020
Although a constant value such as 16 mm/km is sufficiently conservative to overbound all σvig values for the entire year, it fails to capture the seasonality of σvig. In addition, using an overly conservative bounding value for σvig negatively impacts the performance of GBASs (Silva & Monico, 2022). Therefore, to better reflect the seasonality of σvig and improve the GBAS performance, quadratic polynomial expressions that bound all σvig values as a function of the DOY were developed. Although a higher-order polynomial or trigonometric function could also be used, a quadratic polynomial function is chosen for simplicity.
The maximum overbounding σvig and σtg values for each day of 2012–2020 are plotted in Figure 17, as indicated by blue dots. Notably, the σvig results of the first peak in the spring equinox are larger than those of the second peak in the autumn equinox; thus, two independent quadratic polynomial functions were designed to account for this variability. Figure 17 illustrates the determination of the quadratic polynomial functions that serve as upper bounds for the σvig results. The experimental results of σvig obtained by the improved time-step method are fitted by quadratic polynomial curves (red dashed curves). Then, the overbounding curves are obtained by raising the fitted curves to cover all σvig results. A closed-form expression of the bounding quadratic polynomial for the σvig results can be given as a function of the DOY as follows:
13
where doy indicates the DOY. In contrast to σvig, no discernible pattern or trend is observed in the σtg results, as shown in Figure 17. Therefore, a conservative constant value of 5.5 mm/km is selected to bound all σtg results.
Closed-form expressions as functions of the DOY to bound all σvig and σtg values
5 CONCLUSIONS
The time-step method requires only a single satellite and station for ionospheric gradient estimation and allows for an adjustable IPP separation distance. However, the estimation results contain both spatial and temporal gradients, which increase σvig beyond its intent of covering the spatial gradient only. To address this issue, an improved time-step method was developed to separate spatial and temporal gradients based on SatRef data for Hong Kong. The spatial gradients are obtained by smoothing the time-step ionospheric gradient estimates, and the temporal gradients are then computed by subtracting the ionospheric spatial gradient from the time-step ionospheric gradient. Based on the improved time-step method, the nominal ionospheric spatial and temporal gradients for 2012–2020 were investigated. The overbounding values of σvig and σtg were independently determined for each DOY. The overbounding σvig values reveal a seasonal trend, with larger values occurring at the spring and autumn equinoxes. To account for this trend, two quadratic polynomial expressions that provide upper bounds for the overbounding σvig values were designed as functions of the DOY. These closed-form expressions are expected to achieve better GBAS performance. Because no discernible pattern or trend was observed in the overbounding σtg results, a conservative constant value of 5.5 mm/km was determined to bound all overbounding σtg values. This improved time-step method can also be used to investigate nominal ionospheric spatial and temporal gradients for other regions with similar reference station network configurations.
HOW TO CITE THIS ARTICLE
Li, W., & Jiang, Y. (2024). Improved time-step method for bounding nominal spatial and temporal ionospheric gradients for ground-based augmentation systems in Hong Kong. NAVIGATION, 71(2). https://doi.org/10.33012/navi.642
CONFLICT OF INTEREST
The authors declare no potential conflicts of interest.
ACKNOWLEDGMENTS
The work described in this paper was supported by the Research Grants Council of the Hong Kong Special Administrative Region, China (Project No. 25202520; 15214523), and the National Natural Science Foundation of China (Grant No. 42004029).
This is an open access article under the terms of the Creative Commons Attribution License, which permits use, distribution and reproduction in any medium, provided the original work is properly cited.
Jump to section
Related Articles
Cited By...
- No citing articles found.