Abstract
In ground-based augmentation system (GBAS) approach service type D (GAST-D), various ionospheric monitors are implemented in both aircraft and ground facilities to detect ionospheric anomalies. Additionally, the largest undetected differential range errors induced by ionospheric anomalies must be examined because these errors are used in geometry screening to identify potentially unsafe satellite geometries. Based on the ionospheric front threat model, a closed-form expression of the largest undetected ionospheric range error has been established for GBAS approach service type C (GAST-C), where only ground ionospheric monitoring is involved. This paper presents a closed-form expression for GAST-D, and both the ionospheric front model and plasma bubble threat model are taken into consideration. Based on exhaustive simulations among all possible ionospheric threat conditions, the expression is determined as a linear function of the relative speed and gradient magnitude of the ionospheric anomaly. Compared with the linear expression of ionospheric errors for GAST-C, the expression for GAST-D demonstrates that the use of additional ionospheric monitors and a smaller time constant for the code-carrier smoothing filter can effectively reduce the largest undetected ionospheric range error.
1 INTRODUCTION
Located at heights between approximately 50 km and 1000 km, the ionosphere is a shell of dispersive medium that contains electrons and electrically charged atoms (Misra & Enge, 2004). As global positioning system (GPS) signals travel through the ionosphere, they experience a delay, which is proportional to the density of total electron content. A ground-based augmentation system (GBAS) augments existing GPS signals by broadcasting differential corrections and integrity information to enable precision approaches for aircraft (Hoffmann & Walton, 2018). Under nominal conditions, the ionospheric delay is partly eliminated by differential corrections, and residual ionospheric errors caused by spatial decorrelation between the ground and aircraft are bounded by the integrity parameter σvig (Chang et al., 2021). However, in the case of an ionospheric anomaly, the spatial decorrelation might become too large to be bounded. For such cases, ionospheric monitors come into play. Still, undetected errors always exist, which may induce large vertical position errors, thus posing a risk to user safety.
One way to mitigate this threat is by screening out potentially unsafe satellite geometries, known as geometry screening (Pullen et al., 2009). A satellite geometry is unsafe if it has an unacceptably large position error induced by a hypothetical worst-case ionospheric anomaly, even if it has an acceptable vertical protection level (VPL). Unsafe satellite geometries are screened out by inflating σvig to ensure that the VPL exceeds the vertical alert limit and to thus ensure that all unsafe geometries are unusable (Lee et al., 2011). With the geometry screening procedure, the largest ionosphere-induced differential range error (Er) must be determined in order to compute the worst-case ionosphere-induced position error. Previously, Er has been obtained by exhaustive simulations under specific scenarios and parameter ranges. The results are then used to generate look-up tables for future reference (Lee et al., 2006). This process involves uncertainty and low efficiency, as the exhaustive simulations must be repeated if the scenario or parameter ranges change. Kim et al. (2021) proposed a closed-form expression as a bound for all possible scenarios and parameter ranges, provided in the form of a linear function of the gradient and relative speed of the ionospheric anomaly. The expression established for GBAS approach service type C (GAST-C) is based on the ionospheric front model used to describe mid-latitude ionospheric anomalies. In GAST-C, a code-carrier divergence (CCD) monitor is applied on the ground to detect ionospheric anomalies for category I precision approaches (Simili & Pervan, 2006). The geometry screening method is implemented at the ground station and must screen out all potentially unsafe satellite geometries that users might experience.
GBAS approach service type D (GAST-D), which has stricter performance requirements, was proposed to support more demanding CAT II/III precision approaches. This paper focuses on establishing a linear closed-form expression of the largest undetected ionosphere-induced differential range error for GAST-D. In GAST-D, the responsibility of detecting ionospheric anomalies is shared between the aircraft and ground. A CCD monitor and an ionosphere gradient monitor (IGM) are implemented on the ground, whereas dual-solution pseudorange ionospheric gradient monitoring (DSIGMA) is implemented on the aircraft (RTCA DO-253D, 2017). Moreover, geometry screening is transferred to the aircraft to make use of the known aircraft satellite geometries, which simplifies the screening procedure (Lee et al., 2011; ICAO, 2018). This transfer also improves availability because it avoids the conservatism implied by the ground, aiming to ensure that all possible subset geometries are safe for the aircraft (Marini-Pereira et al., 2021). The satellite geometry screening in GAST-D GBAS is completely different from the ground geometry screening previously described for GAST-C GBAS. The vertical projection factor for satellites is evaluated by comparing max(svert2) × EIG with max(Ev), where svert2 is the sum of the two largest absolute elements of the vertical projection matrix svert and EIG and max(Ev) are the maximum value of the ionospheric gradient error and user-specified values for the vertical position error limit, respectively (RTCA DO-253D, 2017). The satellite geometries are screened out when max(svert2) × EIG is larger than max(Ev). EIG is calculated on the aircraft by adding YEIG and the product of MEIG and the distance from the GBAS reference point to the runway threshold. The parameters YEIG and MEIG in message type 2 are broadcasted by the GBAS ground station. The first parameter, YEIG, is the maximum value of the ionospheric gradient error at the GBAS reference point, and the second parameter, MEIG, is the slope of the maximum ionospheric gradient error versus distance (ICAO, 2018). The values of YEIG and MEIG must be determined to ensure that all predetermined potentially undetected Er values remain below EIG.
Using a technique similar to that for GAST-C (Kim et al., 2021), this paper presents a method for generating the largest undetected Er for GAST-D in order to help determine EIG for each threshold and to determine YEIG and MEIG. A closed-form expression of the largest Er as a linear function of the ionospheric gradient and speed is given in this paper. Different from the expression for GAST-C, the derivation of the closed-form expression for GAST-D GBAS considers both the ionospheric front model and plasma bubble model; moreover, ionospheric monitors are distributed on both the ground and aircraft.
This paper is organized as follows. First, the geometric models of impact scenarios under an ionospheric front and plasma bubble are described. Then, an explicit expression of Er induced by an ionospheric anomaly is derived, followed by an analysis of ionospheric monitors to obtain undetected ionosphere-induced errors. Finally, the closed-form expression that bounds the largest undetected Er is derived based on exhaustive simulation results.
2 IONOSPHERIC IMPACT SCENARIOS
Researchers have studied the characteristics of ionospheric anomalies in mid-latitude and low-latitude regions, e.g., the continental United States (Lee et al., 2017), Europe (Robert et al., 2018), Brazil (Lee et al., 2015), and the Asia-Pacific region (Saito & Yoshihara, 2017). Two types of simplified threat models, i.e., the wedge model and trapezoid model, have been used to describe the ionospheric front and plasma bubble, respectively (Pullen et al., 2009; Saito et al., 2009; Saito et al., 2017). The wedge model used to describe an ionospheric front contains three parameters: the gradient, i.e., the linear change between the maximum ionospheric delay and nominal ionospheric delay, g; the front width, w; and the front propagation speed relative to the speed of the ionospheric pierce point of ground (IPPGF) Δvm. The trapezoid model for a plasma bubble is formed by a pair of wedges, which contains seven parameters: a pair of gradients, g1, g2; the width of each wedge, w1, w2; the width of plasma bubble depletion, wb; and the speed of the plasma bubble relative to the speed of IPPGF, Δvl. The total differential delay is assumed to be the same on both sides of the plasma bubble (Depth = g1w1 = g2w2).
The geometries of the aircraft and ground under the impact of a moving ionospheric front and plasma bubble are shown in Figures 1 and 2. To fully examine the impact of the ionospheric front and plasma bubble, an “ahead case” and a “behind case” are considered in the scenarios. In the ahead case, the aircraft and ionospheric anomaly are assumed to move in the same direction, and the ionospheric anomaly affects the aircraft prior to the ground. In contrast, in the behind case, the ionospheric anomaly and the aircraft are assumed to move in opposite directions, and the ground is affected by the anomaly prior to the aircraft.
Geometry model for the aircraft and ground under an ionospheric front (a) Ahead case with the aircraft and front moving in the same direction (b) Behind case with the aircraft and front moving in opposite directions
Geometry model for the aircraft and ground under a plasma bubble (a) Ahead case with the aircraft and bubble moving in the same direction (b) Behind case with the aircraft and bubble moving in opposite directions
During the approach, the aircraft is assumed to move at a constant velocity of vac (70 m/s) parallel to the runway toward the landing threshold point (LTP), which is assumed to be located at a distance of x = 5 km from the centroid of the ground. Because the distances from the aircraft and ground to their corresponding IPPs (IPPAC and IPPGF) are much smaller than their distances to satellites, the relative speed between IPPAC and IPPGF is assumed to be equal to vac. D denotes the initial distance between the ionospheric anomaly and the ground in the ahead case. Dair denotes the initial distance between the LTP and the aircraft in the behind case. The parameters shown in Figures 1 and 2 are assumed to remain constant during the aircraft approach.
The generation of undetected Er (Erundetected) for GAST-D is shown in Figure 3. The ionospheric delays for the aircraft (IAC) and ground (IGF) pass through carrier smoothing code (CSC) filters, and Er is calculated from the difference between two CSC filter outputs (ÎAC and ÎGF). When all ionospheric monitors fail to detect the ionospheric anomaly with the required probability, the Er at the LTP is recorded as Erundetected.
Diagram of Erundetected generation caused by an ionospheric anomaly
3 EXPLICIT EXPRESSION OF DIFFERENTIAL RANGE ERRORS
Er is computed by subtracting ÎGF from ÎAC:
1
IAC, IGF, and the CSC filter are analyzed first. The relevant variables and their definitions are listed in Table 1.
Nomenclature
3.1 Ionospheric Delay
IAC and IGF are determined by the position of IPPAC and IPPGF relative to the ionospheric front or plasma bubble. Because the aircraft and ionospheric anomaly move in the same direction in the ahead case, slow-moving and fast-moving scenarios are considered in order to ensure that the aircraft experiences the entire ionospheric delay during the approach. In the slow-moving scenario, the speed of the ionospheric anomaly relative to IPPGF is less than the speed of the aircraft, i.e., Δvm < vac, Δvl < vac, and the IPPAC at the beginning of the approach is assumed to occur at the rear edge of the ionospheric front model or plasma bubble model, as shown in Figure 4(a) and Figure 5(a). In contrast, the IPPAC at the beginning of the approach is assumed to be at the front edge of the ionospheric front model or plasma bubble model under the fast-moving scenario with Δvm > vac, Δvl > vac. The geometric models of fast-moving scenarios are shown in Figure 4(b) and Figure 5(b). In the behind case, IPPGF is assumed to occur at the rear edge of the ionospheric front model or plasma bubble model, as shown in Figure 6 and Figure 7. The explicit expressions of IAC and IGF for the ionospheric front model and plasma bubble model in the ahead and behind cases are derived as follows.
IAC and IGF under an ionospheric front in the ahead case (Kim et al., 2021) (a) slow-moving scenario (b) fast-moving scenario
IAC and IGF for the plasma bubble model in the ahead case (a) slow-moving scenario (b) fast-moving scenario (c) ionospheric delay for aircraft and GF
IAC and IGF under an ionospheric front in the behind case
IAC and IGF for the plasma bubble model in the behind case
The derivation of IAC and IGF under the ionospheric front model in the ahead case is the same as in GAST-C (Kim et al.,2021). The temporal changes in IAC and IGF are shown in Figure 4. gAC and gGF are the temporal gradients for the aircraft and ground, respectively; TAC and TGF are the time periods for the aircraft IPP and ground IPP to pass through the entire slope of the ionospheric front, respectively; tLTP is the time period for the aircraft to arrive at the LTP; tGF is the time period before the ionospheric front starts to affect the ground. Using the unit step function u(t), IAC for the slow-moving and fast-moving scenarios in the ahead case, 
and 
, can be expressed as follows (Kim et al., 2021):
2
where t indicates the current time, R(t) =tu(t) is a ramp function, the subscript m indicates the ionospheric front model, the superscripts sl and fa indicate the slow-moving and fast-moving scenarios, respectively, and the superscript ah indicates the ahead case. Moreover, IGF for the slow-moving and fast-moving scenarios, 
and 
, can be expressed as follows:
3
The relevant parameters are calculated as follows (Kim et al., 2021):
4
5
6
This paper also considers IAC and IGF under the ionospheric front model in the behind case, as shown in Figure 6. tAC is the time at which the ionospheric front starts to impact the aircraft. For the ionospheric front model in the behind case, IAC and IGF are denoted as 
and 
respectively:
7
where the superscript be indicates the behind case. The related parameters are calculated as follows:
8
9
The threat model of a plasma bubble is also studied in this paper, where the temporal change in IAC and IGF follows the same trend under slow-moving and fast-moving scenarios in the ahead case, as shown in Figure 5. TAC1, TAC2, and TAC3 are the time periods for the aircraft IPP to pass through the nearest slope, depletion, and farthest slope, respectively; TGF1, TGF2, and TGF3 are the time periods for the ground IPP to pass through the nearest slope, depletion, and farthest slope, respectively. gAC1 and gAC2 are the temporal gradients for the aircraft when its IPP is under the nearest and farthest slopes, respectively. gGF1 and gGF2 are the temporal gradients for the ground when its IPP is under the nearest and farthest slopes. For the slow-moving and fast-moving scenarios in the plasma bubble model, IAC is denoted as 
and 
, respectively. The expressions are defined with u(t):
10
where the subscript l indicates the plasma bubble model.
For the slow-moving and fast-moving scenarios in the plasma bubble mode, IGF is denoted as 
and 
, respectively. The expressions are defined with u(t):
11
The relevant ground parameters are calculated as follows:
12
13
The relevant aircraft parameters are calculated as follows:
14
15
16
IAC and IGF for the plasma bubble model in the behind case are shown in Figure 7. IAC in the behind case, 
, can be calculated from 
in Equation (10) by substituting t with t – tl,AC, because the plasma bubble starts to affect the aircraft after tl,AC in the behind case.
IGF for the plasma bubble model in the behind case, 
, can be expressed as follows:
17
The relevant parameters are calculated as follows:
18
19
20
3.2 CSC Filter
The same CSC filter is implemented on the ground and on the aircraft to reduce the multipath and noise errors in pseudorange measurements (Hatch, 1982):
21
where 
is the smoothed pseudorange. ω = τ /T is the filter weight with τ as the smoothing time constant, where τ is 30 s in GAST-D to reduce the time-delay effect on ionospheric errors compared with 100 s in GAST-C (Konno,2007; RTCA DO-253D, 2017). T is the sample interval, and ρ and ϕ are the input code and carrier measurements.
If we have only I in the input measurements with ρ(t) = I(t), ϕ(t) = −I(t), the Laplace transform of the CSC filter H(s) can be approximated as follows (Kim et al., 2021):
22
where I(s) and Î(s) represent I(t) at the input and output of the CSC filter in the Laplace domain, respectively.
For the ionospheric front model in the ahead case, explicit expressions of Î at the aircraft and ground under the slow-moving scenario are denoted as 
and 
, respectively. Explicit expressions of Î at the aircraft and ground under the fast-moving scenario are denoted as 
and 
, respectively:
23
where 
and 
computed by taking the Laplace transform of 
and 
, respectively. L−1 indicates the inverse Laplace transform. Ω(t) is used to simplify the expression of the CCD monitor output, expressed as follows:
24
In the same manner, explicit expressions of Î for the aircraft and ground in the behind case, 
and 
, are as follows:
25
Similarly, explicit expressions of Î for aircraft in the plasma bubble model with the ahead case under the slow-moving and fast-moving scenarios are denoted as 
and 
, respectively:
26
where 
and 
are obtained by taking the Laplace transform of 
and 
, respectively.
Explicit expressions of Î for the ground in the plasma bubble model with the ahead case under the slow-moving and fast-moving scenarios are denoted as 
and 
, respectively:
27
where 
is obtained by taking the Laplace transform of 
.
In the behind case, Î on the ground for the plasma bubble model, 
, can be expressed as follows:
28
where 
is obtained by taking the Laplace transform of 
. Î for the aircraft in the behind case under the plasma bubble model, 
, can be computed from 
by substituting t with t – tl,AC, because the plasma bubble starts to affect the aircraft after tl AC in the behind case.
4 IONOSPHERIC MONITORS
In GAST-D, a CCD monitor, DSIGMA monitor, and IGM are implemented for ionosphere monitoring. Compared with GAST-C, in which only a CCD monitor is installed on the ground, the overall performance in GAST-D is improved by monitoring ionospheric anomalies at both the ground and aircraft. Explicit expressions for these monitors are given as follows.
4.1 CCD Monitor
The CCD monitor implemented on the ground consists of a second-order cascaded filter to detect the ionospheric gradient (Simili & Pervan, 2006). In the Laplace domain, the output of the CCD monitor, Outccd, can be expressed as follows (Kim et al., 2021):
29
where τccd is the filter time constant of the CCD monitor, which is set as 25 s (Pullen et al., 2017). I(s) is the Laplace form of IGF, with Im,GF(s) denoted for the ionospheric front model and Il,GF(s) for the plasma bubble model.
By substituting 
or 
for I(s) in Equation (29), Outccd for the ionospheric front model in the ahead case, 
, can be explicitly expressed as follows:
30
where Ψ(t) is used to simplify the expression:
31
Outccd for the behind case with the ionospheric front model, 
, can be computed as follows:
32
Outccd for the plasma bubble model in the ahead case, 
, is obtained by substituting 
or 
for I(s) in Equation (29):
33
Similarly, Outccd for the plasma bubble model in the behind case, 
, can be computed by substituting 
into Equation (29):
34
The minimum detectable divergence rate of the CCD (MDDRccd) is computed based on the probability of a false alarm (PFA), probability of missed detection (PMD), and standard deviation of test statistics (Pullen et al., 2017):
35
where σccd is the standard deviation of test statistics, Kffd is the K-factor for a given PFA, and Kmd is the K-factor for a given PMD. The values are set as follows: σccd is bounded as 6.9 mm/s, Kffd is 5.91, and Kmd is 6.0 (Pullen et al., 2017). Because the simulation only involves deterministic parameters, the derived CCD monitor outputs are theoretical and free of noise. Meanwhile, the MDDRccd value is compared with the CCD monitor outputs.
The 
and 
results for different Δvm and Δvl values under the ionospheric front model and plasma bubble model in the ahead case are shown in Figure 8. The simulation results for the ionospheric front model were obtained for g = 500 mm/km and w = 25 km. For the plasma bubble model, g1 = 500 mm/km, w1 =25 km, g2 = 250 mm/km, w2 = 50 km, and wb = 25 km. Detection occurs when 
or 
exceeds MDDRccd, as indicated by the black ellipse.
Simulation of the CCD monitor output (a) Ionospheric front model (b) Plasma bubble model
As shown in Figure 8, the CCD monitor cannot detect an ionospheric front when Δvm is less than 85 m/s with g = 500 mm/km. Moreover, the CCD monitor cannot detect a plasma bubble when Δvl is less than 85 m/s with g1 = 500 mm/km and g2 =250 mm/km. For the plasma bubble model, two peaks with opposite signs are observed due to the presence of double slopes in the plasma bubble model. The absolute value of the peak is determined by g1 and g2. Moreover, 
and 
are larger with earlier detection when Δvm and Δvl are larger. However, the time period during which the CCD monitor is affected by the ionospheric front or plasma bubble tends to be shorter when Δvm or Δvl is larger.
4.2 DSIGMA
On an aircraft, DSIGMA monitors the difference between two CSC filter outputs with different time constants:
36
where the superscripts τ2 and τ1 represent two time constants, i.e., 100 s and 30 s, respectively. For the ionospheric front model, Outdsigma for the slow-moving and fast-moving scenarios in the ahead case, 
and 
, can be explicitly expressed by substituting 
and 
in Equation (33):
37
where Θ(t) is expressed as follows:
38
By substituting 
in Equation (33), Outdsigma for the ionospheric front model in the behind case, 
, can be computed as follows:
39
Similarly, Outdsigma in the ahead case for the slow-moving and fast-moving scenarios in the plasma bubble model are denoted as 
and 
, respectively. The expressions are obtained by substituting 
and 
into Equation (33):
40
Outdsigma for the plasma bubble model in the behind case, 
, can also be expressed by 
with t substituted by t – tl,AC.
The minimum detectable error (MDE) of DSIGMA (MDEdsigma) is computed from PFA, PMD, and the standard deviation of test statistics (Pullen et al., 2017):
41
where σdsigma is the standard deviation of test statistics. The values of Kffd, Kmd, and σdsigma are set as follows: σdsigma is bounded as 0.714 m, Kffd is 5.61, and Kmd is 6.0 (Pullen et al., 2017).
The output of DSIGMA for different Δvm or Δvl values in the ahead case is shown in Figure 9. The simulation results for the ionospheric front model were obtained for g = 500mm/km and w = 25 km. The simulation results for the plasma bubble model were obtained for g1, = 500 mm/km, g2 = 250 mm/km, w1 = 25 km, w2 =50 km, and wb = 25 km. Detection occurs when the output of DSIGMA exceeds MDEdsigma, as indicated by the black ellipses.
Simulation results for the DSIGMA output (a) Ionospheric front model (b) Plasma bubble model
As shown in Figure 9, the DSIGMA monitor cannot detect an ionospheric front when Δvm is between 42 m/s and 98 m/s with g = 500 mm/km. Moreover, the DSIGMA monitor cannot detect a plasma bubble when Δv1 is between 42 m/s and 98 m/s with g1 = 500 mm/km and g2=250 mm/km. Two peaks with opposite signs are observed because of the presence of double slopes in the plasma bubble model. The absolute value of the peak is determined by g1 and g2. Moreover, a larger speed difference between Δvm and vac results in a larger 
, whereas the period during which the DSIGMA monitor is affected by the ionospheric front tends to be shorter. Similarly, a larger difference between Δvl and vac results in a larger 
, whereas the period during which the DSIGMA monitor is affected by the plasma bubble tends to be shorter.
4.3 IGM
An IGM on the ground directly monitors the ionospheric gradient among multiple ground stations by using double-difference carrier phase measurements, where the baselines are on the order of several hundred meters (Khanafseh et al., 2012).
For the ionospheric front model, the IGM outputs in both the ahead and behind cases, 
and 
, are the gradient magnitudes when the ground is affected by an anomaly:
42
For the plasma bubble model, the IGM outputs in the ahead and behind cases are denoted as 
and 
, respectively:
43
The minimum detectable gradient (MDG) of the IGM (MDGigm) for the IGM proposed by Khanafseh et al. (2012) is computed from the PFA, PMD, and standard deviation of test statistics (Pullen et al., 2017):
44
where σigm,co and σigm,in are continuity and integrity standard deviations of the IGM, respectively. The values of σigm,co, σigm,in, Kffd, and Kmd are set as follows: σigm,co and σigm,in are bounded as 16.7 mm/km and 26.3 mm/km, respectively, Kffd is 5.91, and Kmd is 6.0 (Pullen et al., 2017).
The 
results for different Δvm and g values in the ionospheric front model and 
results for different Δvl, g1 and g2 values in the plasma bubble model are shown in Figure 10. The simulation results for the ionospheric front model were obtained for w = 25 km and D = 10 km. The simulation results for the plasma bubble model were obtained for w1 = 25 km, w2= 50 km, wb = 25 km, and D = 10 km. Detection occurs when 
or 
exceeds MDGigm, as indicated by the black ellipse.
Simulation results of the IGM output (a) Ionospheric front model (b) Plasma bubble model
As shown in Figure 10, 
is larger when g is larger, and 
is larger when g1 is larger. For a larger Δvm or Δv1, detection occurs earlier, but the period during which the IGM is affected by the ionospheric front or plasma bubble tends to be shorter.
5 CLOSED-FORM ANALYSIS OF UNDETECTED DIFFERENTIAL RANGE ERRORS
In this section, a closed-form expression for the largest ionosphere-induced Erundetected 
is derived as a function of threat parameters of ionospheric speed and gradient magnitude, i.e., Δvm and g for the ionospheric front model and Δvl, g1, and g2 for the plasma bubble model. Therefore, an exhaustive simulation must be conducted to search for 
among other parameters, e.g., w and D.
The parameter ranges for the ionospheric front model and plasma bubble model (shown in Tables 2 and 3) are defined similarly to those for GAST-C (Kim et al., 2017; Lee et al., 2017; Saito et al., 2017). The upper limit of D is generated as follows:
45
Range of Plasma Bubble Threat Model Parameters for Simulation (Lee et al., 2017; Saito et al., 2017)
Range of Ionospheric Front Threat Model Parameters for Simulation (Kim et al., 2017)
Normally, 200 km is sufficient to generate 
Dmax for the ionospheric front model is denoted as Dm,max which is needed only when Δvm is close to vAC:
46
Dmax for the plasma bubble model is denoted as Dl,max, which is needed when Δvl is close to vAC:
47
For each simulation round with a particular set of threat parameters, Er is calculated only when the aircraft arrives at the LTP to reduce the computation load. The maximum outputs of the CCD monitor 
, DSIGMA monitor 
, and IGM 
during the simulation are used because there are multiple opportunities to detect the ionospheric anomaly (Pullen et al., 2017). Erundetected is recorded when 
is smaller than MDDRccd, 
is smaller than MDEdsigma, and 
is smaller than MDGigm.
For a given w, Δvm, and g for the ionospheric front model in the ahead case, 
can be obtained by searching for the maximum Erundetected across all possible values of D:
48
For the ionospheric front model in the behind case, 
can be obtained by substituting Dair for D in Equation (42). For the plasma bubble model, 
can be obtained in the same way. When 
, and 
are computed, the recovery procedure should be considered. As an example, Figure 11 illustrates the recovery procedure for the CCD monitor.
Recovery procedure for the CCD monitor
As shown in Figure 11, recovery occurs when Outccd remains below the re-admittance level of the CCD monitor (Treccd) for a certain period (Bdelay). The CCD monitor and 
are reset at the time of recovery. The parameters for each ionospheric monitor related to the recovery procedure are listed in Table 4.
Recovery Procedure Parameters
Taking the ahead case of the ionospheric front model as an example, the determination of 
across all possible D values is shown in Figure 12. The results were obtained for g =500 mm/km and w = 25 km. Δvm is set as 40 m/s and 120 m/s for the slow-moving and fast-moving scenarios, respectively. 
is found by searching for the maximum Er within the range of D for which no detection occurs before tLTP for all ionospheric monitors.
for the ionospheric front model, given g = 500 mm/km and w = 25 km (a) Slow-moving scenario for the ahead case, Δvm = 40m/s (b) Fast-moving scenario for the ahead case, Δvm = 120 m/s
can be obtained by exhaustive offline simulations, as shown above. Based on the simulation results, a linear closed-form expression is developed to bound the results as an analytical way to represent the results. Although a polynomial expression can bound 
more tightly, a simplified, conservative linear expression is used.
for the ionospheric front model is shown in Figure 13. The results are divided into two parts: severe gradients with gradient magnitudes larger than or equal to 250 mm/km and moderate gradients with gradient magnitudes smaller than 250 mm/km. The results for different gradient values are shown by different colors.
Exhaustive simulation results of the largest undetected differential range errors induced by an ionospheric front (a) Ionospheric front model with severe gradients (b) Ionospheric front model with moderate gradients
In Figure 13, only the ahead case is presented for the severe gradient condition, while both the ahead case and the behind case are presented for the moderate gradient condition. This difference is due to the fact that the geometric model in the behind case assumes that IPPGF occurs at the rear edge of the ionospheric front model, which would be detected by the IGM when the simulation begins. When Δvm is small, the CCD monitor is unable to detect an ionospheric front. However, because of the presence of the DSIGMA monitor and the IGM, the ionospheric front can be detected before the aircraft arrives at the LTP with no Erundetected generated. For a larger Δvm, 
reaches a maximum because all ionospheric monitors fail to detect the ionospheric front before the aircraft arrives at the LTP. As Δvm becomes larger, 
starts to decrease because the ionospheric front is observed by the ionospheric monitors before 
reaches a maximum. When Δvm is sufficiently large, 
converges to a value proportional to the gradient magnitude. For the moderate gradient condition, both the behind case and ahead case must be considered. In addition, 
for the behind case is larger than that of the ahead case because the ionospheric front affects the ground earlier than the aircraft and cannot be detected by the ground IGM. 
continues increasing as Δvm increases until the ionospheric front can be detected by the DSIGMA monitor.
The linear expression of 
for the ionospheric front model is divided into two parts: severe gradients with g ≤ 250 mm/km and moderate gradients with g < 250 mm/km, as shown in Figure 14. The blue dashed line represents exhaustive simulation results with g = 500 mm/km and g = 200 mm/km for severe gradients and moderate gradients, respectively. The black solid line indicates the linear expression as a function of Δvm and g.
Example of the linear bound as a function of Δvm and g of the ionospheric front model (a) Ionospheric front model with a severe gradient (b) Ionospheric front model with a moderate gradient
As shown in Figure 14, three transition points, namely A, B, and C, are defined to bound 
under a severe gradient. Point A is used to indicate Δvm when the ionospheric front is fully mitigated without 
generated. In contrast, point B indicates the maximum 
. The upper bound of 
between point A and point B is conservatively set as the maximum 
. The x-coordinate of point C denotes the Δvm value when 
is the smallest, and the y-coordinate of point C indicates the converged 
value when Δvm is close to 500 m/s. To obtain a linear expression of 
, it is necessary to define the expressions of points A, B, and C as a function of Δvm and g. 
, and 
are used to denote the Δvm values at points A, B, and C, respectively, where the superscript se indicates the severe gradient condition. The 
values for points B and C are represented by 
and 
, respectively.
and 
as a function of g (mm/km) are shown in Figure 15, indicated by blue and black lines, respectively. The red dots represent the simulated results of Δvm when 
reaches the maximum value. For a given g, the speed difference between vAC and Δvm for the red dots can be approximately calculated from MDEdsigma, τ1 and τ2.
49
Determination of Δvm for points A and B
Then, 
and 
can be derived analytically as the upper and lower limits of the red dots:
50
is obtained when the physical separation of the ground IPP and aircraft IPP is equal to zero:
51
and 
are determined by a fitting method, as shown in Figure 16. The red stars in the top and bottom figures show the simulated results of the maximum and converged 
. The black line represents the fitted result expressed as a function of g:
52
Determination of 
for points B and C
With points A, B, and C determined, the linear closed-form expression for 
in the ionospheric front model with a severe gradient can be expressed as follows:
53
where Erde is a decreasing linear function:
54
The established linear expression has the advantage that the structure remains the same if we vary the aircraft speed as input. In other words, the parameters of the expression can be adjusted based on the input aircraft speed with no need to define a new expression. To demonstrate this feature, we used a vAC of 82.8 m/s as an example; this value is the aircraft speed at the LTP, as defined in a previously reported speed profile (ICAO, 2018). The simulation results and linear expression are shown in Figure 17.
Simulation results and linear expression for g = 500 mm/km and vAC = 82.8 m/s in the ionospheric front model for the ahead case
As shown in Figure 17, the linear expression can still bound the simulation results for the largest undetected ionosphere-induced differential range error after vAC changes.
In the same manner, two transition points, A and B, are defined to establish the linear expression to bound 
for the ionospheric front model under a moderate gradient. The x-coordinate of point A, i.e., 
, indicates the Δvm value beyond which the DSIGMA monitor is able to detect the ionospheric front. The upper bound of 
is conservatively set as the maximum 
when Δvm is smaller than or equal to 
. The y-coordinate of point B, i.e., 
, denotes the converged 
when Δvm is equal to 500 m/s.
The linear closed-form expression for 
in the ionospheric front model with a moderate gradient can be expressed as follows:
55
where the superscript mo indicates the moderate gradient condition and 
is 500 m/s.
and 
are determined as a function of g by a fitting method using the simulation results:
55
can also be expressed as a function of g with the fitting method:
56
The simulation results and closed-form linear expression as a function of Δvl, g1, and g2 for the plasma bubble model are shown in Figure 18. The established linear expression needs to bound 
induced by the two gradients of the plasma bubble model. Similar to the linear expression established for the ionospheric front model, the results and linear expression are divided into two parts. One component corresponds to severe gradients, when at least one of the two gradients is larger than or equal to MDGigm (250 mm/km). The second component corresponds to moderate gradients, when the magnitudes of both gradients are smaller than 250 mm/km.
Example of a linear expression as a function of Δvl, g1 and g2 for the plasma bubble model (a) Plasma bubble model with severe gradients (b) Plasma bubble model with moderate gradients
For severe gradients, a large 
is induced in the behind case, where the ground is affected by the plasma bubble earlier than the aircraft when Δvl is small. As Δvl becomes larger, the largest 
occurs in the ahead case when one of the two gradients is smaller than MDGigm. When Δvl becomes sufficiently large, 
converges to a value proportional to the magnitude of the larger of the two gradients. For moderate gradients, 
increases to its maximum and starts to decrease when the plasma bubble can be detected by the DSIGMA monitor.
Points A and B are defined as two transition points for the linear expression of 
:
57
where the subscript l indicates the plasma bubble model. 
is equal to 
defined for the ionospheric front model under severe gradients. 
and 
are determined by the fitting method:
58
where gmax = max{g1, g2} is the larger of the two gradients for the plasma bubble model.
is also determined by the fitting method:
59
Similarly, the explicit expression for 
under moderate gradients can be expressed as follows:
60
where 
and 
are determined by the fitting method using simulation results:
61
Here, 
is 500 m/s, and 
is determined as a function of gmax by the fitting method:
62
5.1 Comparison Between GAST-C and GAST-D
The linear expressions of 
for GAST-C and GAST-D are compared in Figure 19. The blue solid line and red dashed line represent the linear expressions for GAST-C and GAST-D, respectively, for the ionospheric front model with g = 500 mm/km. The green dashed line shows the linear expression of GAST-D for the plasma bubble model with g1 = 500 mm/km and g2=500 mm/km. The linear closed-form expression of GAST-C for the ionospheric front model was obtained from Kim et al. (2021).
Comparison of linear expressions between GAST-C and GAST-D
As shown in Figure 19, due to the small divergence rate induced by the near-stationary ionospheric front, it is difficult for the CCD monitor to detect these fronts in GAST-C (Luo et al., 2005). Therefore, 
reaches a maximum under slow-moving ionospheric fronts. In contrast, GAST-D applies an additional DSIGMA monitor, and thus, the near-stationary ionospheric fronts are well detected without 
generated. However, it should be noted that for the plasma bubble model, 
can still reach several meters when Δvl is small. Besides the DSIGMA and CCD monitors, GAST-D applies an IGM on the ground to directly detect the ionospheric gradient. Moreover, GSAT-D utilizes a smaller time constant of 30 s in the CSC filter. The additional ionospheric monitors and reduced CSC filter time constant largely reduce the delay effect and result in a decrease of 38% for 
when Δvm and Δvl are sufficiently large.
The results in Pullen et al. (2017) show that 
can be bounded by 2.75 m under a mid-latitude ionospheric front, whereas our linear expression shows that 
can be as large as 3.34 m. The main reason for this inconsistency is that a speed profile defined in ICAO Annex 10 is used in Pullen et al. (2017), while a constant aircraft speed is adopted in this research. With a constant aircraft speed, the DSIGMA monitor cannot take advantage of aircraft speed changes. This also explains why relatively large range errors occur when the speed of an ionospheric anomaly is close to the aircraft speed. This increase in range error is disadvantageous because it may impact the service availability. However, the availability reduction is limited to the condition in which the anomaly speed is close the aircraft speed.
6 CONCLUSION
This paper has established linear closed-form expressions to bound the 
induced by ionospheric anomalies for GAST-D under threat models for an ionospheric front and plasma bubble. In GAST-D GBAS, geometry screening is transferred from the ground station to the aircraft to obtain higher availability. The GAST-D geometry screening in the aircraft system requires EIG, which is computed from a linear combination of the up-linked parameters YEIG and MEIG. Because EIG must bound all predetermined potentially 
values, the established expressions can be used to help determine YEIG and MEIG for each runway threshold.
The expressions are divided into two parts to bound 
under severe and moderate gradients. Under severe gradients, near-stationary ionospheric fronts can be fully mitigated by the DSIGMA monitor, whereas slow-moving plasma bubbles are likely to induce an Erundetected of several meters. Under moderate gradients, an expression is determined as a decreasing linear function when the DSIGMA monitor is able to detect the anomaly. The established expression can be applied to the development of GAST-D operation for mid-latitude and low-latitude regions. Previous results with GAST-C show vulnerability in defending against near-stationary ionospheric fronts. The additional DSIGMA monitor in GAST-D greatly improves the mitigation of such cases under severe gradients. In addition, the smaller time constant used in the CSC filter reduces the time-delay impact on ionospheric errors. The derived linear expressions depend on the assumed ionospheric threat models and given monitor thresholds, and a constant aircraft speed is used in the derivation. Further study will involve aircraft speed profiles, as defined in ICAO Annex 10.
HOW TO CITE THIS ARTICLE
Li, W., & Jiang, Y. (2023). Closed-form study of undetected range errors induced by ionospheric anomalies for GAST-D GBAS. NAVIGATION, 70(4). https://doi.org/10.33012/navi.603
CONFLICT OF INTEREST
The authors declare no potential conflicts of interest.
ACKNOWLEDGMENTS
We gratefully acknowledge support from Research Grants Council (RGC) (grant no. 25202520) and National Natural Science Foundation of China (NSFC) (grant no. 42004029), which funded this work.
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.