WAAS and the Ionosphere – A Historical Perspective: Mitigating Mesoscale Irregularities

  • NAVIGATION: Journal of the Institute of Navigation
  • April 2026,
  • 73
  • navi.757;
  • DOI: https://doi.org/10.33012/navi.757

Abstract

To enable the use of global navigation satellite systems (GNSSs) for aircraft navigation, satellite-based augmentation systems have been implemented worldwide to guarantee the accuracy and integrity of aircraft position estimates derived from observations of GNSS signals. For over two decades, the United States’ Wide Area Augmentation System (WAAS) has protected users of the Global Positioning System from threats to position accuracy posed by ionospheric disturbances over North America. A prior companion paper (Sparks et al., 2022) reviews how WAAS has protected users from the disruptive impact of moderate and extreme iono spheric storms. The present paper addresses the methodology adopted by WAAS to protect users from the influence of ionospheric disturbances that are more modest in magnitude, both those well-sampled and those poorly sampled. A subsequent companion paper traces in greater detail the evolution of the WAAS undersampled ionospheric irregularity threat model used to augment the integrity confidence bounds that quantify position accuracy.

Keywords

1 INTRODUCTION

The ionosphere is the largest source of aircraft positioning error for singlefrequency users of global navigation satellite systems (GNSSs). Ionospheric disturbances can dramatically amplify the delay that a radio signal experiences as it travels from an emitting navigation satellite to an airborne receiver. Over the past two decades, satellite-based augmentation systems (SBASs) have been implemented worldwide to ensure the accuracy and integrity of position estimates computed from GNSS measurements for aircraft navigation. Each SBAS broadcasts navigation messages enabling a user (1) to estimate the ionospheric delays that affect GNSS position estimates and (2) to bound reliably the associated estimation error. Ionospheric threats to GNSS positioning accuracy that an SBAS must counter vary in magnitude depending upon the degree of ionospheric disturbance. Over North America, the level of ionospheric perturbation is monitored by the Wide Area Augmentation System (WAAS), the satellite-based augmentation of the United States’ Global Positioning System (GPS). Sparks (2012) presented an introduction to the influence of space weather on the operation of WAAS.

A prior companion paper (Sparks et al., 2022) addressed how WAAS, from its inception to the present, has protected users from the potentially dangerous impact of ionospheric storms. The present paper addresses the methodology used by WAAS to mitigate threats arising from more modest disturbances, i.e., irregularities in electron density that escape local detection either due to their limited spatial extent or because they form in regions that are sampled poorly by network measurements. This methodology relies upon an undersampled ionospheric irregularity threat model to augment the integrity confidence bounds that cover dependably the error in estimates of the ionospheric delay affecting GPS signals detected by a WAAS user. This paper provides an overview describing how the ionospheric threat model is constructed. As in the prior companion paper, throughout this paper, the terms local and irregularity refer to the mesoscale, i.e., [∼100-1000 km]. A subsequent companion paper (Sparks et al., 2026) traces the evolution of the threat model from its original incarnation at WAAS’s commissioning on July 10, 2003, to the threat model fielded in 2023. These three papers may be regarded as supplementing the broad review by Walter et al. (2018) reporting the improvements in WAAS service and integrity since 2003.

WAAS must necessarily bound with a sufficiently high probability all delay estimation errors that arise, independent of the state of the ionosphere. In the absence of a means to distinguish between nominal and disturbed conditions, WAAS would be required to assume that disturbed conditions are always present. The ability to detect ionospheric disturbances allows WAAS to provide a high level of service during periods of nominal, quiet conditions. A robust detection scheme guarantees that a reduced level of service need be imposed only during periods when significant ionospheric disturbances are observed.

WAAS relies on geostationary satellites to disseminate correction messages that allow a WAAS user to bound with this high degree of certainty the position error due to ionospheric delay. WAAS is regarded as broadcasting hazardously misleading information (HMI) whenever the actual error in a user’s vertical (or horizontal) position estimate exceeds a certain mandatory limit. The most restrictive integrity specification required of WAAS is that, for aircraft navigation using precision approach with vertical guidance, the probability of broadcasting HMI must not exceed more than one occurrence in every 10,000,000 runway approaches (resulting in either landings or missed approaches). This corresponds to a probability of broadcasting HMI below 10−7.

To verify that the probability of broadcasting HMI does not exceed this limit, analyses of WAAS operations rely on a fault tree in which each node of the tree is associated with a distinct threat to system integrity, and the probability of that threat occurring is allocated an upper bound. Since ionospheric disturbances constitute only one of several possible sources of HMI, the probability of failure allocated to the monitor responsible for ensuring that a user’s position error is adequately bounded must be less than the 10−7 value allocated to the entire system. Accordingly, the fault tree allocates to this monitor a value of 2.25×10−8 per approach as an upper limit on the probability of broadcasting HMI owing to the missed detection of an ionospheric irregularity. The extreme storm detector, discussed below and described by Sparks et al. (2022), protects WAAS users from the adverse impact of major ionospheric storms. Thus, the probability restricting the broadcast of HMI applies only at times when the extreme storm detector has not been triggered.

Transmitted at uniform time intervals, WAAS messages consist, in part, of data that reflect the local state of the ionosphere, defined at regularly spaced points on a grid over North America. These data are comprised of a set of ionospheric grid delays (IGDs) and a corresponding set of grid ionospheric vertical errors (GIVEs). At each ionospheric grid point (IGP), the broadcast IGD represents an estimate of the ionospheric delay that a signal propagating vertically through that point would experience. The GIVE assigned to that IGP provides an extremely conservative integrity confidence bound on the corresponding delay estimation error. Vertical delay estimates and their confidence bounds are derived from fits of equivalent vertical delay, i.e., slant delay measurements recorded by the WAAS network of dual-frequency receivers, converted to estimates of vertical delay using the well-known thin-shell model of the ionosphere. As required by the Minimum Operational Performance Standards (MOPS) for Global Psositioning System/Satellite-Based Augmentation System Airborne Equipment (RTCA, 2016), WAAS sets the model shell height to 350 km, a representative height of the peak of the ionosphere’s F2 layer. The fit methodology that generates the IGDs has been summarized by Sparks et al. (2022), and a more detailed description has been presented by Sparks et al. (2011a).

A regional mask of the global IGP grid stipulated by the MOPS (RTCA, 2016) identifies the IGPs at which IGDs and GIVEs are broadcast. Figure 1 shows the set of IGPs that comprise the IGP mask adopted in Release 62-CY23. Color-coded markers display the minimum GIVE value broadcast at each IGP. Drawn in blue is the corresponding pierce point filter line (PPFL) used to eliminate measurements from fits of ionospheric delay at low latitude, where the ionosphere tends to be disturbed even under nominal ionospheric conditions. Figure 3 of Sparks et al. (2022) shows how the set of IGPs that comprise the IGP mask has evolved over the life of WAAS. Since the commissioning of WAAS in 2003, this mask has been modified six times, sometimes in response to an expansion of the number of wide-area reference stations (WRSs) in the WAAS receiver network. (Note: in Figures 3(c) and 3(d) of Sparks et al. (2022), the location of the PPFL was erroneously plotted 5° higher in latitude than its actual fielded position.) The locations of the WRS sites currently in operation are listed in Table A1.

FIGURE 1

WAAS IGP mask for Release 62-CY23 (May 19, 2023)

The color of the square at each IGP indicates the minimum value (in meters) of the GIVE broadcast at that IGP. Plotted in blue is the PPFL used to remove measurements from fits of ionospheric delay at low latitude. Black squares identify IGPs near Hawaii where, owing to the isolation of the WRS in Honolulu from other WAAS WRSs, GIVEs are always set to “not monitored.”

To determine the ionospheric delay correction (and its integrity confidence bound) associated with a given GPS signal detected by a user’s single-frequency receiver, the WAAS user must first determine the ionospheric pierce point (IPP) at which the signal raypath intersects the ionospheric grid at its specified shell height of 350 km. An estimate of the vertical delay at this point is then obtained by performing bilinear interpolation of the IGDs at the nearest IGPs surrounding that point. In a similar fashion, an integrity confidence bound associated with the vertical delay estimate, designated the user ionospheric vertical error (UIVE), is calculated by interpolating the GIVEs at the same IGPs. To generate an estimate of the slant delay along the raypath from the satellite to the user, the interpolated vertical delay estimate is multiplied by the appropriate thin-shell obliquity factor, as determined by the elevation angle of the raypath at the receiver. To obtain an integrity confidence bound for this slant delay estimate, the UIVE is multiplied by the same obliquity factor. With a slant delay estimate and confidence bound for each detected GPS signal, the user can reduce the corrupting influence of ionospheric delay from his or her position estimate and derive a confidence bound on the error in that estimate. This confidence bound is used, in turn, to evaluate the user’s horizontal protection level (HPL) and vertical protection level (VPL). WAAS becomes unavailable for a user when the HPL or VPL exceeds, respectively, the horizontal alert limit (HAL) or the vertical alert limit (VAL) associated with a given navigation mode and level of aviation service (RTCA, 2016; Sparks et al., 2011b).

The level of aviation service determines the minimum height above the runway down to which WAAS can provide vertical guidance. If the required visual reference to continue a precision approach has not been established by this decision height, a missed approach must be initiated. The decision height assigned to each level of aviation service and the associated VAL are distinct. For localizer performance with vertical guidance (LPV) service, the decision height and VAL are 250 feet and 50 m, respectively; for LPV200 service, these values are 200 feet and 35 m, respectively.

Under quiet conditions, a planar model of the ionosphere has been found to describe accurately the spatial variation of vertical ionospheric delay on meso-scales (Walter et al., 2001). The WAAS GIVE monitor is designed to address a single ionospheric threat, namely, the threat that local mesoscale behavior in the vicinity of a given IGP does not conform well to the planar model assumed when estimating the vertical delay at that IGP. The GIVE monitor evaluates IGDs and GIVEs to be broadcast and ensures that the integrity of the GIVEs is sufficient to guarantee that the UIVE, derived by the user at his or her IPP, bounds vertical delay errors with the mandated probability. For WAAS purposes, the standard χ2 goodness-of-fit parameter associated with an IGD estimate serves as a reliable indicator of the local level of ionospheric disturbance near the IGP (Walter et al., 2001). Each GIVE is calculated, in part, from the formal error of the estimate of the vertical delay at the IGP, inflated to account for the statistical uncertainty of this error. Prior to being broadcast, the GIVEs determined by the GIVE monitor may be increased or set to “not monitored” by the range domain monitor (RDM) or the user position monitor (UPM), which, like the GIVE monitor, are part of the WAAS integrity data monitoring capability.

The GIVE monitor allocates to a set of ionospheric disturbance detectors part of its responsibility for providing protection from deviations of estimated vertical delay from planarity (see the work by Sparks et al. (2022)). A local irregularity detector (LID) at each IGP is triggered when an irregularity metric proportional to the χ2 goodness-of-fit parameter for a fit of equivalent vertical delay centered on the IGP exceeds a specified threshold, causing the GIVE monitor to set the GIVE at that IGP to a conservative value, GIVEmax, that safely and reliably bounds the maximum WAAS ionospheric vertical delay estimation error ever observed. When the GIVE at an IGP is GIVEmax, WAAS ceases to support vertical guidance at that location. To indicate the state of the ionosphere as a function of time on larger spatial scale (i.e., over North America), WAAS defines the ionospheric perturbation metric (IPM) to be the maximum value of the irregularity metric over the entire WAAS IGP grid. This metric serves as the basis for both the WAAS extreme storm detector (ESD) and the WAAS moderate storm detector (MSD). The ESD distinguishes extreme ionospheric storms from less intense disturbances and determines the length of the duration following extreme storm onset over which GIVEs are set to GIVEmax at all WAAS IGPs. The tripping of the MSD forces the GIVE monitor to broadcast GIVEs that are more conservative than would be broadcast under nominal conditions. Figure 2 of Sparks et al. (2022) presents a process flow diagram that shows how the ESD and the MSD fit into the WAAS scheme for monitoring ionospheric disturbances (the operation of the LIDs is effectively subsumed in the box labeled Compute IPM).

By reducing the degree of correlation between measurements of ionospheric delay recorded by WAAS WRSs and those observed by a user, the presence of well-sampled mesoscale irregularities generally causes the magnitude of the user’s formal delay estimation error to increase. Whether this formal error accurately reflects the true error, however, depends upon how well the electron density fluctuations in the ionosphere are sampled. This, in turn, depends upon the spatial distribution of the measurement raypaths that connect GPS satellites to WAAS receivers. Undersampling can occur in regions where the measurement coverage is either sparse or highly non-uniform.

To shield the user from threats due to mesoscale density gradients that escape detection by the GIVE monitor, the inflated formal error in the GIVE must be augmented by an additional amount. The GIVE monitor retrieves the amount of this GIVE enhancement from a table of adjustments provided by a branch of the WAAS undersampled ionospheric irregularity threat model (Sparks et al., 2001; Altshuler et al., 2001; Sparks et al., 2026). Currently, the threat model is comprised of four branches: (1) the quiet-time branch consists of a table of adjustments to be consulted under nominal, quiet ionospheric conditions, (2) the disturbed-time branch consists of a table to be accessed under disturbed conditions when the MSD has been triggered (Sparks et al., 2022), (3) the unknown extreme storm state branch consists of a table of relatively large adjustments to be used when the state of the ionosphere is unknown (e.g., after a system reset), and (4) the “realistic” branch consists of a table (similar to that of the disturbed-time branch) designed to model the difference between estimates of vertical delay based upon fits that include all measurements in an input data set and estimates of the same vertical delay based upon fits that exclude measurements whose IPPs lie in a region of generally enhanced ionospheric structure at low latitude (see the discussion of Sparks et al. (2026) concerning the PPFL shown in Figure 1).

An updated ionospheric threat model has been fielded with each system release that has expanded the number of network receivers or that has modified the algorithm used to construct the threat model. Each of these updated threat models has incorporated any threats found to have arisen from ionospheric disturbances occurring subsequent to the fielding of the prior threat model. It is, of course, crucial that an adjustment to the GIVE be large enough to bound safely all neighboring undersampled threats; however, if the threat model is too conservative, the resulting GIVEs will be large, and system availability will suffer. WAAS requires a method for assembling the threat model that will safely handle threats posed by the ionosphere without being overly pessimistic. In practice, the adjustment provided by the ionospheric threat model is generally the largest term defining the GIVE, and consequently, the ionosphere becomes the determining factor limiting WAAS performance (see Equation (2) below).

Each branch of the undersampled ionospheric irregularity threat model is constructed offline from sets of archived measurements of GPS signals that fail to sample (or poorly sample) disturbed regions of the ionosphere, where each set is recorded by WAAS receivers in a single epoch. An estimate of the vertical delay at the IPP of each test signal in a given epoch is derived from a fit of the equivalent vertical delay of signals with neighboring IPPs in that epoch. The fit is centered on the IGP nearest the test signal’s IPP. The difference between the estimated vertical delay and the equivalent vertical delay of the test signal then defines a fit residual. Fit residuals calculated in this fashion are tabulated to determine the maximum values they have assumed as a function of metrics characterizing the spatial distribution of measurement IPPs about the corresponding fit centers.

To extend the range of data sets examined, a technique designated data deprivation is used to simulate alternative samplings of the given set of ionospheric disturbances actually observed (as described by Sparks et al. (2026)). Two types of data deprivation masks have been defined: those that exclude measurements that lie in the interior of the coverage region and those that exclude measurements whose IPPs lie at the edge of coverage. The former masks permit evaluation of the impact of irregularities that are sufficiently compact to escape detection, whereas the latter masks enable an assessment of the influence of irregularities that avoid detection owing to their location at the fringes of the sampling region.

Figure 2 presents a flow chart illustrating how each of these branches is tabulated offline. The unknown extreme storm state branch, the disturbed-time branch, and the quiet-time branch are restricted to tabulation of delay fit residuals associated with IPPs near each IGP where the fit centered on that IGP has failed to trigger the LID (Sparks et al., 2022). The largest of these fit residuals found to occur in historical data sets determines the table of adjustments to apply to the GIVE when the ionospheric state is unknown. The quiet-time and disturbed-time branches exclude from tabulation fit residuals that arise when the ESD has been triggered. In addition, the quiet-time branch excludes from tabulation all fit residuals that occur when the MSD has been triggered. As discussed by Sparks et al. (2022), the system-wide ESD assumes responsibility for protecting the user from compact disturbances that have been observed to originate during extreme storms and have been found to persist many hours into nighttime. WAAS assumes that such highly localized events occur only following the onset of extreme storms. If responsibility for protecting the user from these events were assumed by the undersampled ionospheric irregularity threat model, the fit residuals that arise from these disturbances would require increasing the magnitude of the GIVE augmentation provided by the threat model, thereby severely reducing WAAS availability under nominal ionospheric conditions (see Section 8).

The primary objective of this paper is to describe how, over the first 20 years of WAAS operation, the undersampled ionospheric irregularity threat model and its contribution to defining the GIVE at each IGP have protected WAAS users from the detrimental influence of mesoscale ionospheric disturbances that fail to trip any of WAAS’s LIDs. This paper first reviews how the GIVE has inflated the variance of the formal error in the vertical delay estimate at each IGP to address threats associated with well-sampled irregularities. It then examines how the GIVE has relied on the WAAS undersampled ionospheric irregularity threat model to provide protection from threats posed by poorly sampled irregularities. (A detailed account of the considerations and decisions that have driven the evolution of the ionospheric threat model is deferred to a subsequent companion publication ((Sparks et al., 2026).) As noted in our first companion paper (Sparks et al., 2022)), prior publications concerning the technology discussed here have generally been limited to conference proceedings. Once again, our secondary objective is to provide a comprehensive bibliography for accessing these publications.

FIGURE 2

Process flow diagram showing how the branches of the WAAS undersampled ionospheric irregularity threat model are tabulated

The process consists of the following steps: (1) measurements are recorded at WAAS WRSs, (2) measurements are converted to equivalent vertical delay, (3) estimates of vertical delay are computed at fits centered on IGPs, (4) the “realistic” branch is tabulated by keeping all threats that arise when no data deprivation (see Section 4) is invoked, (5) the unknown extreme storm state branch is tabulated using only the LID to exclude data, (6) the disturbed-time branch is tabulated using both the LID and the ESD to exclude data, and (7) the quiet-time branch is tabulated using the LID, the ESD, and MSD to exclude data.

Section 2 of this paper describes how the GIVE at each IGP is defined. Section 3 examines how this GIVE protects the user from the influence of ionospheric irregularities that are well-sampled. Section 4 presents an overview of how the GIVE provides protection from the impact of irregularities that are poorly sampled. Section 5 investigates how the values of various fit parameters, averaged over time under nominal conditions, vary as a function of the position of the IGP at which an IGD is calculated. Section 6 explores threat characteristics that govern the structure of the current ionospheric threat model. Section 7 analyzes the geographic distribution of undersampled threats that have appeared over the lifetime of WAAS operations. Section 8 offers insight as to how the burden of protecting the user from these threats is distributed among the ESD, MSD, LIDs, and the undersampled ionospheric irregularity threat model. Section 9 provides a summary. Appendix A contains tables that identify (1) the current WAAS receiver sites, (2) terms defining the GIVE as they have been modified over time, and (3) the critical points that determine the structure of the quiet-time branch of the Release 62-CY23 threat model.

2 SPECIFICATION OF THE GRID IONOSPHERIC VERTICAL ERROR

The integrity confidence bound on the error in a user’s estimated vertical ionospheric delay must address the statistical uncertainty arising from both the measurement error and the failure of the sampled regions of the ionosphere to match exactly the assumed model for ionospheric delay. In addition, this bound must protect the user from the presence of ionospheric disturbances that are undersampled. Furthermore, it must account for the error that accrues from the spatial interpolation of broadcast IGDs. Finally, it must reckon with a possible increase in the estimation error occurring between the observation time of the measurements upon which the estimation of the vertical delay is based and the time at which the user evaluates the confidence bound to constrain his or her position error.

At the v-th IGP, the GIVE is specified in terms of a standard normal (Gaussian) distribution designed, in the offline evaluation of the ionospheric threat model, to overbound the tails of the actual distribution of the fit residual error for estimates of vertical delay near the IGP (Sparks et al., 2011b). By convention (RTCA, 2016), the GIVE at an IGP is defined to be the product of the standard deviation of this overbounding distribution and a constant, KGIVE = 3.29, that defines an interval corresponding to a confidence level of 99.9%. The overbounding distribution is specified in terms of a distribution whose standard deviation is designated (σGIVE,v. The distribution is constructed by requiring that it overbound the tails of the actual error distribution to a distance of at least KHMI σGIVE,v from the origin. As mandated by the MOPS (RTCA, 2016), the constant KHMI is set to 5.33 ensuring that the probability of broadcasting HMI does not exceed 10−7, i.e., assuming a Gaussian distribution, KHMI determines a confidence interval with a confidence level of (1–10−7). If the ionosphere were the sole source of user range error, the broadcast GIVE would be simply defined as KGIVE σGIVE,v. However, since there are other sources of user range error, a tighter restriction on the amount of range error due to the ionosphere must be imposed. This is achieved by broadening the distribution that overbounds the fit residual errors, i.e., by multiplying σGIVE,v by the ratio KHMI_GIVE/KHMI, where KHMI_GIVE is a constant that defines a tighter confidence interval. As noted in Section 1, the WAAS fault tree allocates to the GIVE monitor a value of 2.25×10−8 per approach as an upper limit on the probability of broadcasting HMI due to the missed detection of an ionospheric irregularity. This corresponds to a KHMI_GIVE value of 5.592. Thus, the GIVE at the v-th IGP becomes: GIVEv[KGIVEKHMIGIVEKHMIσGIVE,v]indexed1 where the brackets indicate that the GIVE values broadcast are quantized. In a user receiver ionospheric correction message, the GIVE is quantized to a discrete value by rounding its computed magnitude upward to the nearest quantization level. The MOPS (RTCA, 2016) specifies these quantization levels in meters (at L1) to be 0.3, 0.6, 0.9, 1.2, 1.5, 1.8, 2.1, 2.4, 2.7, 3.0, 3.6, 4.5, 6.0, 15.0, and 45.0. One additional quantization level may be used to indicate that the GIVE is “not monitored.”

The maximum quantization level of 45 m at the GPS frequency L1 (1575.42 MHz), designated GIVEmax, is sufficiently large to bound all residual errors at L1 under all ionospheric conditions. A computed GIVE that exceeds GIVEmax is reset to GIVEmax. Furthermore, WAAS broadcasts GIVEmax at any IGP where the LID has tripped. In practice, a floor value, stipulated at each IGP to increase system integrity, limits access to the lower quantization levels. The color-coded markers in Figure 1 display the minimum GIVE value broadcast at each IGP.

In each WAAS release after its Initial Operating Capability (IOC), the standard deviation σGIVE,v of the overbounding distribution upon which the GIVE at an IGP is based has been determined by a variance of the following form: σGIVE,v2=σGIVE,well-sampled,v2+σGIVE,undersampled,v22 where σGIVE,well-sampled,v2 establishes an error bound that mitigates the adverse impact of estimation error arising from well-sampled irregularities and σGIVE,undersampled,v2 augments this error bound to protect the user from estimation error generated by undersampled ionospheric threats. σGIVE,well-sampled,v2 is primarily based on the formal error variance estimated at each IGP, whereas σGIVE,undersampled,v2 is determined by a lookup table compiled from vertical delay fit residuals that have been observed during historical ionospheric disturbances, as discussed in Section 7 of the subsequent paper in this series (Sparks et al., 2026). Table A2 provides a summary of how the error variances used to evaluate σGIVE,v2 have evolved over the life of WAAS.

In addition to broadcasting GIVEs at IGPs, each proportional to σGIVE,v, WAAS broadcasts ancillary parameters in Message Type 10 that allow the user to inflate further each σGIVE,v prior to performing the interpolations that establish the UIVE for each GPS signal received by the user. The objective of this additional inflation is to compensate for possible growth in the level of ionospheric disturbance during the period over which the broadcast GIVEs are to remain in force. The MOPS (RTCA, 2016) provides two options (designated σionogrid,ν2) for degrading the σGIVE,v at an IGP: σionogrid,v2=(σ˜GIVE,v2+εiono)2,option1,orσ˜GIVE,v2+εiono2,option2,3 where: εionoCiono_ramp(tt0)+Ciono_step[tt0tmessage]floor.4

Here, Ciono_ramp determines the linear rate of continuous change in εiono, t is the current WAAS network time, t0 is the time at which WAAS begins transmitting an ionospheric correction message at the geostationary satellite, Ciono_step determines the size of a step-function jump in εiono, tmessage is the minimum update interval for ionospheric correction messages, and the brackets [x]floor indicate taking the greatest integer less than x. For each upgrade of the ionospheric threat model, the option selected in Equation (3) and the operational system parameters Ciono_step, Ciono_ramp, and tMessage have been determined such that εiono always produces an inflated σGIVE,v that overbounds the temporal threats in question over the life of the correction message. The temporal duration over which a correction message may be operative is twice tmessage (to cover the possibility of a missed message) plus the time required by the system to calculate IGDs and GIVEs for the maximum number of IGPs and to include these data in the broadcast message. For its initial operations, WAAS adopted option 2 and continued to use option 2 until WFO Release 3A in 2011 when it switched to option 1.

3 PROTECTING THE USER FROM WELL-SAMPLED IONOSPHERIC IRREGULARITIES

As discussed by Sparks et al. (2022), WAAS defines an ionospheric irregularity as any ionospheric behavior that cannot be accurately described by the assumed planar model (i.e., gradients in electron density that are well-represented by this model are not considered to be “irregularities”). The formal error in the vertical delay estimated at a given point provides a basis for formulating a confidence interval that reliably bounds the delay estimation error arising from well-sampled irregularities located near that point. Here, we designate σ2 to represent this formal error at a given observation point located in the neighborhood of an arbitrary fit center, and σv2 identifies the formal error of this point when the fit center is the v-th IGP. The term σIGP,v2 is the counterpart of σv2 after undergoing inflation. To define σGIVE,well–sampled,v2, it is first necessary to consider σ2 as the sum of two terms representing the variances of process noise due to the ionosphere and measurement noise originating in the observing instruments (see the discussion leading to Equation (20) of Sparks et al. (2022)): σ2(w)=σprocess2(w)+σmeasurement2(w),5 where: σprocess2(w)wTCw2wTc+c0,6 and: σmeasurement2(w)wTMw.7

Here, w is a vector of weights applied to a fit of N measurements centered on an IGP and determined by minimizing the estimation variance subject to the constraint that the estimator be unbiased, C is an N x N matrix describing the covariance of the detrended ionospheric delay between measurement locations, c is an N-vector whose elements specify the covariance between the detrended delay at the measurement locations and the scalar field describing small irregularities that are superimposed on the planar trend at a position ∆x near the IGP, c0 is the variance of the scalar field, assumed to be a constant independent of ∆x, and M is an N x N matrix describing the covariance of the measurement noise between measurement locations. All quantities in Equations (6) and (7) have been defined in the prior companion paper (Sparks et al., 2022).

Under nominal conditions, the irregularity detector at each IGP can be safely trusted to determine when random statistical fluctuations in the ionosphere exceed a threshold that renders the user’s computed integrity confidence bounds unreliable (Sparks et al., 2022). When ionospheric conditions become slightly worse than nominal, however, the probability that the detector will fail to react properly to a local, well-sampled ionospheric threat becomes finite. To compensate for this statistical uncertainty, the WAAS GIVE monitor inflates at each IGP the formal error upon which the GIVE is based. The magnitude of the inflation depends on the degree of ionospheric turbulence: the more disturbed the ionosphere, the greater the inflation. Magnifying the formal error in this fashion ensures that the irregularity detector at the IGP is triggered whenever warranted by a significant increase in the level of ionospheric noise neighboring the IGP, thereby rendering the broadcast values of the IGD and the GIVE at the IGP unusable.

The remainder of this section describes how the algorithm that serves to protect users from the potentially inimical influence of well-sampled ionospheric irregularities has evolved since the commissioning of WAAS in 2003.

3.1 Initial Operating Capability

Studies conducted prior to the commissioning of WAAS in July of 2003 (Hansen et al., 2000a; Hansen et al., 2000b) demonstrated that, under nominal conditions, the difference in delay imposed by the ionosphere on two radio signals propagating vertically is, to zeroth order, a linear function of the horizontal distance separating the signals. This behavior is consistent with the planar model of the vertical delay assumed by WAAS in each fit. It was initially anticipated that the first-order decorrelation function used to relate ionospheric measurements made at separate locations would be a function of (1) the measurement elevation angles and (2) the distance between the thin-shell pierce points. These studies concluded, however, that the first-order deviations of vertical delay from planarity are uncorrelated, implying that (σdecorrtotal)2, the delay covariance associated with widely separated IPPs, can be treated as a constant matching (σdecorrnominal)2, the delay covariance associated with nearly coincident IPPs. Setting (σdecorrtotal)2 equal to (σdecorrnominal)2 simplifies various terms that appear in the estimation of ionospheric delay as presented by Sparks et al. (2022): the covariance matrix C in Equation (6) becomes diagonal, the vector c vanishes, c0 reduces to (σdecorrnominal)2, and the weighting vector w at an IGP simplifies as follows: w=WG(GTWG)1s,8 where s = sIGP specifies the weighting at the IGP: sIGP[100]T,9

W is the weighting matrix: W[C+M]110

and the rows of the G matrix identify the locations of the IPPs of the measurements in the fit relative to the IGP position: G[1Δx1Te^eastΔx1Te^north1Δx2Te^eastΔx2Te^north1ΔxNTe^eastΔxNTe^north]11 with e^east and e^north representing east and north unit vectors, respectively, in an Earth-centered, Earth-fixed local coordinate system. When Equation (8) and C, c, and c0 as defined by Sparks et al. (2022) are substituted into Equation (6), the variance of the process noise evaluated at the v-th IGP then reduces to the following: σprocess,v2=sT(GTWG)1s+(σdecorrnominal)2.12

The first term on the right-hand side reflects the uncertainty in estimating the orientation of the ionospheric plane, and the second term describes the inherent uncertainty in vertical delay about that plane. Based on the studies cited above, IOC WAAS assigned a value of 35 cm to σdecorrnominal. It was later concluded that the first-order decorrelation function describing measurement errors results from multipath rather than from ionospheric perturbations; thus, it is more accurate to affirm that the ionospheric decorrelation at mid-latitude is bounded by 35 cm about a plane (Walter et al., 2001).

WAAS determines the state of the irregularity detector at an IGP according to an irregularity metric derived from the χ2 goodness-of-fit parameter generated by the fit of vertical delay at the IGP. The goodness-of-fit parameter is evaluated as follows (see Equations (25) and (32) of Sparks et al. (2022)): χ2IIPPT[WWG(GTWG)1GW]IIPP13 where: I¯IPP[I¯1I¯2I¯N]T14 is the vector of slant delays converted to vertical at each of the fit IPPs. In IOC, the irregularity metric was defined as follows: χirreg2χ2χ2(1PfaN¯)15 where χ2(PN¯) denotes the inverse of the chi-square probability function for N¯ degrees of freedom evaluated at P, i.e., the χ2 value that bounds the fraction P of the distribution below it (Abramowitz & Stegun, 1972), and Pfa is a permissible rate of false alarms set to 10−3 in IOC. The irregularity detector trips when the following condition holds: χirreg2>Tirreg,trip16 where Tirreg,trip is the value of χirreg2 at which the assumed planar model ceases to be an accurate representation of local ionospheric behavior. The irregularity detector trip threshold Tirreg,trip was specified in IOC to be 1.

To make certain in IOC that the computed GIVEs always satisfied a probability of broadcasting HMI below 10−7, WAAS adhered to the highly conservative assumption that the ionosphere is always in a state just below the level guaranteed to trigger each LID. An inflation factor Rirreg was designed to amplify the magnitude of the formal error at each IGP to correspond to this worst-case condition. To establish a numerical value for Rirrg, WAAS adopted a tractable, if pes-simistic, model in which all formal error variances were uniformly increased by Rirreg2 (Walter et al., 2001). The value of Rirreg was then determined according to the statistics of a chi-squared distribution. Given a threshold for triggering each irregularity detector based upon the permissible rate of false alarms Pfa, Rirreg was assigned the value required to achieve an allowable probability of missed detection Pmd. Since the Rirreg2 defined here is a constant that depends only upon Pfa, Pmd, and the number of degrees of freedom, but not upon the state of the ionosphere, it has been designated Static Rirreg2 or Rirreg–static2: Rirreg-static2χ2((1Pfa)N¯)χ2(PmdN¯)17 where N¯ here represents the degrees of freedom, N – 3, when N measurements are fit to retrieve the three model parameters representing the estimated vertical delay at the IGP and its gradients in the eastern and northern directions. In IOC, Pmd was assigned a value of 2.25×10−4.

To evaluate both the differential ionospheric correction and the integrity confidence bound associated with the measurement of a given GPS signal, the WAAS user must interpolate the IGDs and GIVEs broadcast at the corners of the grid cell in which the measurement IPP is located. The interpolation is performed according to the bilinear scheme defined in the MOPS (RTCA, 2016). Since the user interpolates IGDs and GIVEs at (generally) four IGPs to form an estimated vertical delay and a UIVE at the user IPP, confidence intervals must bound vertical delay errors in all regions between grid points. (Note: when the user IPP lies at the edge of the IGP grid, its estimated vertical delay and UIVE may be determined in WAAS from three-point interpolation, rather than four-point interpolation, provided the IPP lies within the right triangle formed by the three IGPs nearest the IPP.) Figure 3 displays the dependence of the weight assigned to the IGD or the GIVE at a given IGP according to its distance from the IPP in four-point interpolation. For example, when the IGP at the center of Figure 3 serves as the lower left-hand IGP of the interpolation cell in which a given IPP lies, then the interpolation weight assigned to this IGP is specified by the color of the point at the IPP location.

To establish in IOC a confidence bound on the statistical spatial error for estimates of vertical delay located in the vicinity of an IGP, Walter et al. (2001) relied on the fact that the formal error due to the assumed planar model exhibits a global minimum near the weighted centroid of the measurement IPPs in the fit (where each IPP weight is determined by the ionospheric and measurement noise) and that this error increases as one moves away from that point. The equation governing the increase in the variance that bounds the error at a distance ∆X from the v-th IGP is as follows:

σv2(δθ,δϕ)qvT(δθ,δϕ)(GTWG)1qv(δθ,δϕ)18

where the angles δθ and δϕ represent the change in latitude and longitude, respectively, from the IGP position, W and G are the matrices defined by Equations (10) and (11), respectively, and qv identifies the location of interest relative to the IGP: qv(δθ,δϕ)[1 ΔxT(δθ,δϕ)e^east ΔxT(δθ,δϕ)e^north].19

To bound the worst interpolated estimation error that a user may experience in the four cells surrounding each IGP, Walter et al. (2001) pointed out that it is sufficient to consider the deviations from the assumed planar model of vertical delay that occur in the threat domain of each IGP. At each IGP, this domain consists of the set of points situated closer to the given IGP than to any of the surrounding IGPs, i.e., the rectangular region bordered by the midpoints in latitude and longitude between the IGP and its nearest neighboring IGPs. In Figure 3, the magenta rectangle identifies this domain for the IGP at the center of the plot. The maximum in the bounding variance in this domain will occur at the farthest distance from the weighted centroid of the measurement IPPs, i.e., at one of the domain corners, as illustrated by Walter et al. (2001). Consequently, taking the maximum of the standard deviations of the formal error at the four corners of the magenta box in Figure 3 defines a confidence interval that, in IOC, bounded the statistical spatial error for an IPP near the central IGP, i.e.,σv2(σvmax corner)2, as follows: σvmax cornermax[σvcorner(δθN,δϕE)σvcorner(δθN,δϕW)σvcorner(δθS,δϕE)σvcorner(δθS,δϕW)],20 and δθN,δθS,δϕE, and δϕW are the angular distances from the IGP to the midpoints of the segments connecting the IGP to the nearest neighboring IGPs to the north, south, east, and west, respectively, i.e., the corners of the domain enclosed by the magenta rectangle in Figure 3.

For IOC, the inflated error variance adopted to protect the user from estimation error due to well-sampled irregularities was as follows:

σGIVE,well-sampled,v2Rirreg-static2[(σvmax corner)2+(σdecorrnominal)2].21

Rirregstatic2 was not allowed to drop below 1; if its calculated value was less than 1, it was reset to be 1. Multiplying (σvmax corner)2 at each IGP by this factor ensured that each irregularity detector would be triggered with the mandated probability of missed detection.

Rather than inserting Equation (21) directly into Equation (2) to evaluate σGIVE,v2 IOC defined σGIVE,v2 by combining terms from σGIVE,well-sampled,v2 and σGIVE,undersampled,v2 as follows:

σ˜GIVE,v2Rirreg-static2σ˜vmaxcorner2+maxRirreg-staticσdecorrnominal2,σundersampled,v2+σROT2,22

where σundersampled,v2 bounded those estimation errors occurring as a result of spatial undersampling (see Section 4) and σROT2 augmented εiono2 in option 2 of Equation (3) with Ciono_step = 0 in Equation (4). (Since the data used to fit Ciono ramp in Equation (4) did not conform to a function that vanished at t = t0 and increased quadratically in t, IOC was forced to introduce a non-zero σROT2 to ensure that σionogrid,v2 adequately addressed the threat posed by the growth of fit residuals between evaluations of planar fits; see Section 7.1.2 of Sparks et al. (2026)) This equation implicitly assumes that undersampled ionospheric conditions in the vicinity of the IGP are either nominal and covered by (Rirreg-staticσdecorrnominal)2 or disturbed and covered by σundersampled,v2. As a consequence, σundersampled,v2 here effectively serves as a lower bound for the second term on the right-hand side of Equation (22). This approach to combining protection from well-sampled threats with protection from undersampled threats was subsequently modified in LPV Release 6/7.

3.2 LPV Release 6/7 and Subsequent Updates

The approach to inflating the formal estimation error taken in IOC was extremely conservative. As discussed above, the applied inflation factor Rirreg-static2 was independent of the actual state of the ionosphere. Having the GIVE monitor treat the ionosphere as if its level of disturbance was, at all times, nearly sufficient to trip an LID neglected the fact that the ionosphere rarely experiences disturbances of this magnitude. While this approach helped to guarantee the integrity of the user’s computed confidence bounds, it resulted in large GIVE values that unnecessarily restricted system availability under nominal ionospheric conditions.

FIGURE 3

Dependence of the central IGP‘s interpolation weight on the distance to a user IPP lying within one of the four grid cells depicted

The magenta rectangle encloses the threat domain of the central IGP, namely, the set of points that lie closer to this IGP than to any of the neighboring IGPs. An estimate of the vertical delay at the IPP located in the upper right-hand grid cell can be determined by interpolating the vertical delay at the four corners of this cell. The interpolation weight for the IGP at the lower left-hand corner of this cell (i.e., the central IGP of the figure) is determined by the color at the location of the IPP.

To address this shortcoming, the first major upgrade of the GIVE monitor, occurring in September of 2007, implemented changes to the calculation of the formal error variance that reduced the magnitude of its contribution to the broadcast GIVE without impairing system integrity. First, the definition of Rirreg was modified in LPV Release 6/7 to reflect the current state of the ionosphere. The resulting inflation factor has been termed dynamic Rirreg (or Rirreg–dynamic) to distinguish it from the static Rirreg of IOC. As originally proposed by Altshuler et al. (2002), the χ2((1Pfa)N¯) in the numerator of Equation (17) was replaced with a functional dependence on the χ2 computed from the fit (see Equation (13)), indicating the measured level of disturbance. Furthermore, the denominator of Equation (17) was reformulated in terms of the probability of broadcasting HMI rather than the probability of the irregularity detector failing to be trip under disturbed conditions (Sparks et al., 2011b). Specifically, the χ2(PmdN¯) in the denominator was replaced by χ2N,lowerbound, defined such that if X is a standard normal random variable and Y is a χ2 random variable with N¯ degrees of freedom, then the following inequality holds:

P(|X|(KHMIGIVEYχN¯,lower bound))PHMI,iono,23

where PHMI,iono represents the bound on the probability of broadcasting HMI allocated to ionospheric errors and KHMI_GIVE is the scalar chosen to adjust the protection provided by the GIVE monitor against the broadcast of HMI (as described in the discussion of Equation (1) above), defining a confidence interval corresponding to a confidence level of (1–2.25×10−8). The equation actually used to determine χN¯,lower bound2 was as follows (Blanch, 2003):

PHMI,iono=2πϕ=0π/2(KHMI_GIVE2χN¯,lower bound21sin2(ϕ)+1)N¯/2dϕ24

A parameter study that examined how various choices of PHMI,iono affected system performance determined that PHMI,iono should be set to 10−10, a suballocation of the 2.25×10−8 probability of failure allocated to the GIVE monitor (Blanch, 2003).

As discussed by Sparks et al. (2022), IOC WAAS did not address the possibility that measurement noise could hinder the ability of the χ2 of a fit to reflect accurately disturbed conditions in the ionosphere. To compensate for this contingency, the irregularity detector at each IGP in LPV Release 6/7 incorporated into the irregularity metric χirreg2 an inflation factor designated Rnoise (see Equation (35) of Sparks et al. (2022) for a definition):

χirreg2Rnoiseχ2χ2(1PfaN¯)25

The same reasoning led LPV Release 6/7 to include Rnoise in the definition of Rirreg _dynamic:

Rirreg-dynamic2Rnoiseχ2χN¯,lower bound226

Fortunately, the inclusion of this Rnoise factor, greater than or equal to 1, did not negate the considerable reduction in GIVE magnitudes achieved by moving from static to dynamic Rirreg. As in IOC, Rirreg was not permitted to fall below unity.

Prior to LPV Release 6/7, it was determined that the evaluation of the user’s inflated formal estimation error variance did not require taking the maximum value at the corners of a threat domain centered on an IGP (see Equation (21)), i.e., that this approach was overly conservative since the interpolation scheme already accounted for the spatial variation of the error variance over the grid cell and that including the maximum value at the corners effectively double-counted this error variance (Blanch, 2003). Historically, the development of this algorithm had preceded the development of the undersampled ionospheric irregularity threat model and had been intended to mitigate threats at the edges of the threat domain of each IGP. Only later was it realized that such a threat model would be needed to mitigate these threats. It was then demonstrated that if the interpolation scheme used to evaluate the vertical delay at the user’s IPP were applied in a like manner to the GIVEs broadcast at the IGPs bounding this IPP, the resulting UIVE would satisfy the WAAS integrity requirements. This provided an incentive to simplify Equation (22). Finally, it was realized that the measurement noise is bounded independent of the state of the ionosphere (Blanch et al., 2003b). This implies that Rirreg–dynamic should be used to scale only those terms in the estimation variance that arise due to the ionospheric process noise and not those due to the measurement noise. Consequently, in LPV Release 6/7 (and in all later releases), the inflated formal estimation variance at the v-th IGP has taken the following formml:

σIGP,v2=Rirregdynamic2σprocess,v2+σmeasurement,v227

Implementing this change reduced the excessive conservatism that the formal estimation error had contributed to the broadcast GIVE.

Substituting Equation (12) into Equation (27) defines the inflated formal error used in LPV Release 6/7 and LPV Release 8/9. WFO Release 3, however, implemented a vertical delay estimation model based on Kriging (Blanch, 2002; Blanch, 2003; Blanch et al., 2003a; Sparks et al., 2010; Sparks et al., 2011a), thereby introducing correlations between vertical delay estimates (in σprocess2).(σdecorrtotal)2 was no longer assumed to be equal to (σdecorrnominal)2 and the process noise was as follows:

σprocess,v2wTCw2wTc+(σdecorrtotal)228

in accord with Equation (6).

In IOC, antenna biases were thought to be sufficiently small so as not to require a term accounting for them in the GIVE. Prior to LPV Release 6/7, this assumption was re-examined, and it was decided henceforth to include in σGIVE,well-sampled,v2 a term bounding the maximum error associated with bias in pseudorange measurements at a WRS due to group delay variation in an antenna element. For the κ-th IPP in the fit, the curves that bound the “average” antenna bias for the GPS frequencies L1 and L2 (1227.6 MHz) were specified to be functions of the elevation angle ακ of the raypath connecting the satellite to the receiver:

μL1,κ=cL1aL1+bL1expακ/dL1μL2,κ=cL2aL2+bL2expακ/dL229,

where aL1, bL1, cL1, dL1, aL2, bL2, cL2, and dL2 are constants. For each IPP, these equations were combined to form the maximum average antenna bias using a methodology that mimics the conversion of L1 and L2 pseudorange measurements to vertical delay estimates (as described by Sparks et al. (2022)):

μκ=(fL22fL12fL22)(μL1,κ+μL2,κ)F(ακ,hi).30

The total maximum antenna bias error term was then calculated as the sum of the products of the absolute value of each weighting vector coefficient Wv,κ and the individual maximum average antenna bias for the corresponding fit IPP:

μtot,v=κ=1N|wv,κ|μκ31

where each Wv,κ is an element of the vector w from Equation (8) with s set equal to qvT(Δθ,Δϕ). Here, the arguments ∆θ and ∆ϕ refer to the location of the IPP relative to the v-th IGP.

The total maximum antenna bias error term was incorporated into the variance mitigating the impact of estimation error arising from well-sampled irregularities as follows:

σGIVE,well-sampled,v2[σIGP,v+μtot,vKHMI_GIVE]2.32

The contribution of the total maximum antenna bias error term to σGIVE2 has not been modified in subsequent upgrades of the GIVE monitor. The specification of σGIVE,well-sampled,v2 did not change in LPV Release 8/9.2. The implementation of Kriging in WFO Release 3, however, altered the value of σIGP,v for the fit based on a given set of measurements, but the form of Equation (32) remained the same. Equation (32) continues to govern the contribution to the overbounding distribution σGIVE,v2 that serves to protect WAAS users from threats posed by well-sampled irregularities under current operations.

4 PROTECTING THE USER FROM POORLY SAMPLED IONOSPHERIC IRREGULARITIES

In addition to protecting the user from the influence of well-sampled ionospheric irregularities that fail to trip the nearest LID, the GIVE monitor must also shield the user from threats posed by irregularities that are undersampled. Ionospheric irregularities may be undersampled because (1) they are very localized or (2) they are situated in regions where the distribution of measurement IPPs is highly non-uniform. The onus of protecting the WAAS user from such threats is effectively shared between the ESD, the MSD, the LIDs, and the undersampled ionospheric irregularity threat model.

When WAAS performs a fit of equivalent vertical delay centered on a given IGP (as defined by Sparks et al. (2022)), the three fit parameters retrieved may be used to form an estimate of the vertical delay at any IPP located within the threat domain of the IGP. The estimation error that results from undersampling may be regarded as safely bounded by the inflated formal estimation error variance whenever the following inequality is satisfied:

|I¯κIκ|2>Kundersampled2σκ2,33

where I¯κ is the slant delay of the κ-th signal converted to vertical at the measurement IPP using the elevation-dependent thin-shell obliquity factor, Iκ is the corresponding estimate of vertical delay at the same IPP location, σκ2 is the inflated formal error variance of the delay estimate at the IPP, and Kundersampled2 determines an upper bound on the square of the residual in terms of the inflated formal error variance. I¯κ serves as the truth value for the vertical delay estimate at the IPP in question. In practice, Kundersampled is set to a value of 5.33, corresponding, as previously noted, to a bounding probability of (1–10−7) for a Gaussian distribution.

When the inequality in Equation (33) is violated and neither the ESD nor the LID at the IGP of the fit has been triggered, there may exist an undersampled ionospheric disturbance that constitutes a potential danger to the integrity of a user’s computed error confidence bound on his or her position estimate. Unlike triggering either the ESD or an LID, triggering the MSD is not guaranteed to protect the user from a trigger-inducing threat; this remains the responsibility of the disturbed-time branch of the ionospheric threat model. When the ESD and the LID have not been tripped, the magnitude of the threat will depend upon the spatial distribution of the slant delay measurements in the vicinity of the disturbance.

Under nominal ionospheric conditions, we can expect the probability that the inequality in Equation (33) is satisfied to meet the mandated requirement for any signal measurement. In the presence of ionospheric disturbance, however, this inequality may be more readily violated. To ensure that the broadcast GIVE bounds the residual error |I¯κIκ| with the required probability when Equation (33) does not hold, the square of the standard deviation of the overbounding distribution σGIVE,v is augmented by σGIVE,undersampled,v2 as indicated by Equation (2). This augmentation necessitates the development of a model for ionospheric threats. If WAAS were to set σGIVE,undersampled,v2 to a constant, this constant would need to be sufficiently large to bound the worst errors ever arising from the presence of undersampled irregularities, and this, in turn, would severely restrict WAAS availability. However, since the magnitude of such threats depends upon the spatial distribution of the measurements used to estimate Iκ, availability can be enhanced by making σGIVE,undersampled,v2 a function of one or more metrics that characterize the spatial distribution of the IPPs included in the fit used.

The value of σGIVE,undersampled,2 that the GIVE monitor assigns to the broadcast GIVE at a given IGP is determined by the undersampled ionospheric irregularity threat model. As discussed in Section 1, the threat model is currently comprised of four branches; each branch provides a table of GIVE adjustments to be utilized under various circumstances. The quiet-time branch provides GIVE adjustments under nominal, quiet ionospheric conditions. When the MSD is triggered, the GIVE monitor immediately switches to retrieving the GIVE adjustments from the disturbed-time branch; it returns to retrieving the adjustments from the quiet-time branch once the MSD transitions back to a quiet-time state. The warm-up (or “unknown extreme storm state”) branch of the threat model provides the adjustments during the warm-up period after the WAAS correction and verification (C&V) subsystem has been reset; consequently, the GIVE monitor does not know whether an ionospheric storm has occurred that would have tripped the ESD. At the end of this warm-up period, either the quiet-time branch or the disturbed-time branch provides the GIVE adjustments, depending on the state of the MSD. The “realistic” branch of the threat model is used to define GIVE values for use in the downstream WAAS UPM and is not used in the GIVE monitor.

Each branch of the threat model is tabulated offline in an analysis of a representative set of historical WAAS measurements recorded under disturbed ionospheric conditions. As dictated by the inequality in Equation (33), the tabulation of ionospheric threats used to define σGIVE,undersampled,v2 in Equation (2) must include the threat detected by the κ-th measurement δκ2>0, where δκ2 is defined as follows:

δκ2|I¯κIκ|2Kundersampled2σκ2.34

Under nominal, quiet-time ionospheric conditions, we can expect δκ2, so defined, to be negative. In the presence of ionospheric disturbance, however, it may become positive. For each threat model branch, the dependence of σGIVE,undersampled,v2 on the spatial distribution of fit IPPs is determined by tabulating values of a slight modification of δκ2, a variance difference designated σ¯undersampled,κ2:

σ¯undersampled,κ2Kinflate2|I¯κIκ|2Kundersampled2σκ2,35 where Kinflate is an additional inflation factor greater than or equal to unity. Each table of GIVE adjustments is constructed from values of σ¯undersampled,κ2 that satisfy σ¯undersampled,κ2>0. Setting Kinflate greater than 1 expands the set of threats included in the tabulation of the threat model branch. In IOC, Kinflate was set to 1.1. Subsequently, however, this assignment was deemed to generate overly conservative values of σ¯undersampled,κ2. From LPV Release 6/7 onward, Kinflate has been set to 1.

WAAS has used two metrics to characterize the distribution of measurement IPPs incorporated into each fit of vertical delay: (1) the Euclidean (straight-line) distance Rfit from the fit center IGP at xIGP to the fit domain boundary and (2) the relative centroid metric (RCM), i.e., the ratio of Rcentroid to Rfit whereRcentroid|xcentroidxIGP| is the distance from the fit center IGP at xIGP to the weighted centroid of the fit IPPs at xcentroid. When the positions of the NIPP fit IPPs are expressed in Earth-centered, Earth-fixed rectangular coordinates, the i-th component of xcentroid is calculated in the operational system as follows:

xcentroid,ixIPP,iTWdiag11TWdiag1,36

where:

xIPP,i[xi,1xi,2xi,NIPP]37

and 1 is a column vector of NIPP elements, each of which is 1. Here, xi,j is the magnitude of the i-th component of the position vector for the j-th IPP, and Wdiag consists of the diagonal elements of the weighting matrix defined by Equation (10). The rationale for including the weighting matrix in the calculation of xcentroid is that measurements with large measurement noise should have less impact on the metric. It should be noted, however, that, in the construction of the ionospheric threat model, Wdiag is set to the NIPP x NIPP identity matrix, i.e., all IPPs are weighted equally, since errors in the ionospheric truth data used to construct the threat model are assumed to be negligibly small.

Rfit may be regarded as a proxy for the mean IPP density, which necessarily decreases as the fit domain expands to encircle a targeted number Ntarget of fit IPPs. When fit IPPs are uniformly distributed throughout the fit domain, the RCM is nearly zero. This term approaches unity when the fit IPPs congregate near a single point at the edge of the fit domain. Thus, the RCM may be taken as a measure of the degree of uniformity in the distribution of IPPs throughout the fit domain: a larger RCM value corresponds to a greater likelihood that an irregularity will be undersampled. The measurements included in the fit are selected according to the same search algorithm used to estimate the vertical delay at an IGP, as described by Sparks et al. (2022) and summarized in Section 5.1.

The raw data for each branch of the ionospheric threat model are generated offline by determining the maximum values of σ¯undersampled,κ2, observed over the measurements κ in sets of test data, as a function of these two IPP distribution metrics:

σundersampledraw(Rfit,RCM)(maxκ,tfit,residual(σ¯undersampled,κ2))1/238

where the maximization is conducted over positive values of σ¯undersampled,κ2 within the time interval tfit,residual following each fit epoch. At any given IGP, σ¯undersampled,κ2 is evaluated only at IPPs that lie within the threat domain of the IGP and only when the ESD and the irregularity detectors have not been triggered. Furthermore, when the MSD has tripped, σ¯undersampled,κ2 is excluded as well from the tabulation of the quiet-time branch of the ionospheric threat model. In recent releases, threats have also been excluded from the tabulation of the ionospheric threat model using UIVE floor culling, i.e., excluding from the threat model each threat that is covered by the floor value assigned to the UIVE at the threat IPP, as determined by interpolating the GIVE floors at the corners of the cell in which the IPP is located. Section 7.6 of Sparks et al. (2026) describes in greater detail the methodology used to invoke UIVE floor culling.

Since the probability that a threat will be undersampled should increase as the distribution of fit IPPs becomes more sparse (larger Rfit) or increasingly non-uniform (larger RCM), the tabulated raw data are subjected to a two-dimensional overbound to guarantee that the contribution of the threat model to each GIVE will be a monotonically increasing function of each of the two fit IPP distribution metrics:

σundersampledoverbound(σundersampledraw(Rfit,RCM)).39

Figure 4 shows the quiet-time and disturbed-time branches for the Release 62-CY23 ionospheric threat model. Figure 4(a) displays σundersampledraw(Rfit,RCM) (Rfit, RCM) for the quiet-time branch. No data lie below a fit radius of 800 km, as this is the minimum fit radius used in the fits. Figure 4(b) shows the overbound of the quiet-time raw data. The undersampled irregularity term in Equation (2) can now be defined in terms of σundersampled as specified by Equation (39):

σGIVE,undersampled,v2=σundersampled,v2,40

where the subscript v on the right-hand side indicates that the value of σundersampled has been taken from a threat model table and applied to the v-th IGP. For a fit at this IGP, the table from which the value is taken depends upon whether the MSD has been triggered. Figures 4(c) and 4(d) show the raw data and overbound for the disturbed-time branch of the Release 62-CY23 ionospheric threat model.

Overbounds such as those shown in Figures 4(b) and 4(d) can be specified in terms of their critical points, defined as follows:

P(Rfit,RCM)PcritifP(Rfit,RCM)>P(RfitΔfit,RCM) andP(Rfit,RCM)>P(Rfit,RCMΔrcm)41

In other words, a critical point at the lower left-hand corner of a tabulation bin has a value of σundersampled that is strictly larger than that of the bin directly below it as well as that of the bin directly to the left of it. A listing of the critical points and their associated σundersampled uniquely determines an entire threat model table (see Table A3).

Threat models have been plotted using a color bar scale covering the interval [0, 2] to confine attention to the most relevant values of σundersampled. Values of σundersampled exceeding 2 m give rise to GIVEs of 15 m or greater, at which point a user’s VPL generally exceeds the VAL, making the LPV service unavailable.

FIGURE 4

The Release 62-CY23 spatial-temporal ionospheric threat model: (a) quiet-time raw data, (b) quiet-time overbound, (c) disturbed-time raw data, (d) disturbed-time overbound

In addition to its disturbed-time and quiet-time branches, the Release 62-CY23 ionospheric threat model is comprised of two other branches depicted in Figure 5. Figure 5(a) shows the threat model branch to be consulted when the state of the ionosphere, as monitored by the ESD, is unknown. This situation arises, for example, after a reset of the WAAS C&V subsystem. After such a reset, the ESD is initially assumed to be in a triggered state. The warm-up period tESD,warmup following a reset has always been set to 8 h to be comparable to the ESD recovery confirm interval, as discussed by Sparks et al. (2022). Once the duration of C&V subsystem operation extends beyond the warm-up period, the disturbed-time branch of the threat model is used if the MSD has been tripped; otherwise, the quiet-time branch is used. The computation of the warm-up branch is identical to that of the disturbed-time branch except that the ESD is turned off, i.e., all threats are retained in the tabulation that would otherwise be excluded due to the triggering of the ESD. Figure 5(a) testifies to how conservative the fielded undersampled ionospheric irregularity threat model would need to be to protect WAAS users from all ionospheric events in the absence of an ESD.

The final branch of the Release 62-CY23 threat model, displayed in Figure 5(b), is designated the “realistic” branch. The purpose of the “realistic” σundersampled2 is to model the difference between vertical delay estimation that uses all IPPs located in southern Mexico and the vertical delay estimation that does not use the Mexican IPPs that lie south of the PPFL. This difference affects the performance of the WAAS UPM when users process data broadcast at IGPs south of the PPFL. At the v-th IGP σundersampled–realistic,v is computed in the same fashion as the σundersampled,v of the disturbed-time branch except that (1) the threshold for tripping each LID is set to a large value, so that threats are excluded from tabulation only when the ESD has been triggered, and (2) the tabulation is conducted without data deprivation (only observed, not simulated, threats are tabulated). Other standard features of the disturbed-time branch tabulation, such as filtering IPPs with the PPFL and using the UIVE floor to cull threats, are retained. Values of σundersampled–realistic taken from the “realistic” branch of the threat model are used to compute modified GIVEs only at IGPs south of the PPFL (for an IGP north of the PPFL, σundersampled,v is set to Zero). These modified GIVEs are required for operation of the WAAS UPM (see the work by Walter and Blanch (2017) for a description of the UPM fielded in 2017).

FIGURE 5

The Release 62-CY23 spatial-temporal ionospheric threat model: (a) warm-up (“unknown extreme storm state”) branch, (b) “realistic” branch

5 IGD ESTIMATION PARAMETERS AVERAGED OVER TIME UNDER QUIET-TIME CONDITIONS

When the GIVE monitor calculates the IGD and GIVE to be broadcast for a given IGP, the magnitudes of various estimation parameters depend critically upon the location of the IGP relative to the edge of WAAS coverage. This section examines how these parameters vary spatially under nominal conditions. For each estimation parameter, a plot shows the average value assumed at each IGP over the nine-day interval of June 3–11, 2022, a period of low geomagnetic activity (the geomagnetic Kp index did not rise above 3). There are three parts to this section, each highlighting different aspects of the GIVE monitor’s IGD evaluation algorithm. The first concerns geometric parameters (fit radius Rfit, relative centroid metric RCM, and the number of fit points NIPP) that are determined by the IPP search algorithm, the second examines the actual computation of the IGD and its associated goodness-of-fit parameter χirreg2 as specified by Equation (25), and the third addresses the broadcast GIVE. These plots have been generated using prototype simulations based upon measurements recorded at WAAS stations during the quiet period in question. Disturbed ionospheric conditions result in larger vertical delays, increased magnitudes of the irregularity metric χirreg2, and higher GIVE values.

5.1 Geometric Fit Parameters

Details of the algorithm used by the WAAS GIVE monitor to search for the measurement IPPs to be included in an IGD fit at an IGP have been described by Sparks et al. (2022); a brief summary is provided here. Initially, the algorithm examines all IPPs within a domain centered on the IGP with a fit radius Rfit = Rmin (800 km). If NIPP, the number of IPPs in this fit domain, is less than the targeted number Ntarget = 30, the algorithm increases Rfit until the fit domain encompasses Ntarget IPPs; then, it evaluates the IGD and subsequently the GIVE. If Rfit reaches Rmax (2100 km) without the fit domain having encircled Ntarget IPPs, the IGD and GIVE are evaluated using the NIPP points that have been located, provided that NIPPNmin (where Nmin = 10). If the GIVE monitor fails to locate at least Nmin IPPs within the distance Rmax from the fit center, however, the GIVE to be broadcast is set to “not monitored.”

The values of the algorithmic parameters Rmin, Rmax, Ntarget, and Nmin assigned above were first proposed prior to the commissioning of WAAS in 2003. Early references to these values, cited in conference proceedings (see, for example, the work by Sparks et al. (2002)), are based on an unpublished sensitivity analysis, conducted by the Raytheon Company in 2000, which examined how varying these values affects system availability. These values were used to generate a baseline for system availability on December 12, 2000, consisting of the percentage of the conterminous United States (CONUS) achieving 95% or 100% availability. The values of the geometric fit parameters were then adjusted one parameter at a time, and the results of the study were specified in terms of the percentage increase or decrease of the area of CONUS that was found to maintain these levels of availability. Three different values for each parameter were used to assess performance. The point of this study was not to determine optimal values for these parameters but to demonstrate the sensitivity of the system to these parameters. The study concluded that adjusting Rmax and/or Nmin offered the greatest prospect for enhanced performance but that the level of performance improvement did not justify the substantial effort that would be required to rework the algorithm validation and the HMI analysis of the system. Consequently, the values of the geometric fit parameters have been left at their baseline values throughout the life of WAAS.

Figure 6 shows the disparity in size between a fit domain at the edge of coverage and a fit domain in the interior. At the edge of coverage, Rfit must rise significantly to encircle at least Ntarget IPPs. The fit radius for the IGP marked by the magenta X at [35°N, 130°W] is close to Rmax, whereas that of the IGP marked by the red X at [45°N, 95°W] approaches Rmin. Note: measurements recorded at the WAAS station in Honolulu are excluded here since they are never used to compute IGDs and GIVEs.

As shown in Figure 7(a), the average NIPP is greater than or equal to 30 for nearly all of the IGPs in the grid. Only those IGPs at the outermost edge of coverage give rise to fits with lower values of NIPP, and even here, the domains centered on these IGPs encompass a sufficient number of IPPs to estimate IGDs (NIPPNmin). The various updates of the WAAS IGP set designed to enhance performance and integrity have retained the ability to conduct fits at nearly all IGPs in all epochs.

Figure 7(b) reveals the increase in Rfit needed to capture 30 IPPs in the fit domain as the fit center IGP moves from the interior of the service volume to its edge. Note that over land, a fit radius Rfit = Rmin is sufficient to encircle NIPP = Ntarget IPPs, but as the fit center IGP moves beyond the coastline, a larger fit radius is required to encompass Ntarget points.

FIGURE 6

Example of two IGD fit domains, one centered on an IGP in the interior of the WAAS service volume (red) and a second at the edge of coverage (magenta)

The location of the IGP at the center of each fit domain is indicated by an X. Brown squares indicate WRS locations. Small magenta circles identify the locations of GPS satellites projected onto the thin shell where the measurement IPPs are located.

In addition to NIPP and Rfit, there is a third pertinent geometric parameter, as noted in Section 4, that plays a key role in the construction of the WAAS ionospheric threat model, namely, the RCM. The RCM is defined to lie between 0.0 and 1.0, and the magnitude of the RCM is larger at the edge of coverage. Figure 7(c) shows that the RCM magnitude is generally less than 0.4 for IGPs over land and over the bulk of the WAAS network. Canadian IGPs, however, each tend to have a slightly elevated average RCM, as IPPs are more sparse in the Canadian region. The RCM value increases greatly for IGPs over an ocean. Note that some IGPs have average RCM values of roughly 0.8. Even IGPs that have poor access to IPPs, however, can benefit performance.

5.2 Delay Estimation Parameters

Once the set of fit IPPs has been identified, the vertical ionospheric delay at an IGP is computed via linear regression based on Kriging, as discussed both in Section 3 and in greater detail in the work by Sparks et al. (2022). The values of these estimated delays can vary greatly depending on the time of day and local ionospheric conditions. The values shown in Figure 8(a) represent daily averages of vertical delay under nominal conditions. As expected, the ionospheric delay in Figure 8(a) is much larger near the geomagnetic equator than in the northern regions. It is interesting to note that delay values are slightly larger in the southwest corner of the WAAS grid than in the southeast corner, where the WAAS IGPs are farther from the magnetic equator.

FIGURE 7

Geometric fit parameters at IGPs, averaged over June 3-11, 2022: (a) number of IPPs in the fit, (b) fit radius Rfit in kilometers, (c) relative centroid metric RCM

When the GIVE monitor calculates the vertical delay at a given IGP, the estimation algorithm also returns the χ2 goodness-of-fit statistic (see Equation (13)). Before the χirreg2 irregularity metric can be evaluated at the IGP, however, the inflation factor Rnoise must be computed (see Equation (25)). This term grows as the magnitude of the measurement noise increases. As discussed by Sparks et al. (2022), the χ2 statistic is suppressed when measurement noise becomes sufficiently large, which can result in a failure of χ2 to reflect accurately the magnitude of the error in the estimated vertical delay. As a consequence, an unsafe operational condition can arise. The value of Rnoise computed by the GIVE monitor serves to mitigate such a condition. The time-averaged values that Rnoise assumes under nominal conditions are shown in Figure 8(b). Note that these values can exceed 2 throughout the WAAS IGP grid and rise above 3 at both low and high latitude.

The metric χirreg2 quantifies the quality of a particular fit. Since the underlying physical model used to estimate vertical delay is planar, χirreg2 can be expected to increase in regions where the behavior of the ionosphere is typically non-planar. As shown in Figure 8(c), such is the case at low latitudes, where the largest values of the averaged χirreg2 reflect the influence of the daily traversal of the non-planar equatorial anomaly at low geomagnetic latitude. Note again that Figure 8(c) has been generated under nominal conditions. It is interesting that, under these conditions, χirreg2 attains its lowest values at mid-latitude. (The low values of χirreg2 near the edge of coverage are due to the smaller number of points in the fits at the IGPs located there.)

FIGURE 8

Vertical delay fit parameters at IGPs under nominal conditions, averaged over June 3–11, 2022: (a) the IGD in meters at L1, (b) Rnoise, (c) χirreg2.

The irregularity detector trip threshold Tirreg,trip has been tuned to ensure that the irregularity detector at each IGP trips rarely, even in the southern regions. (Note that even the largest values of χirreg2 in Figure 8(c) are well below the irregularity detector trip threshold of 3.) The threshold value implemented in the construction of the ionospheric threat model is the same as that implemented in WAAS operations. As discussed in Section 4, values of σ¯undersampled,κ2 are tabulated in the threat model only when the LID and the ESD have not tripped. Thus, the magnitude of Tirreg,trip is a key parameter in determining the contribution of σGIVE,undersampled,v2 to σGIVE,v2. It is important to realize, however, that Tirreg,trip is not an integrity parameter, i.e., adjusting this threshold does not alter which threats are covered. When Tirreg,trip is increased to a higher value, some threats that would have caused an irregularity detector to trip (thereby setting the GIVE to its maximum value of 45 m) must now be incorporated into the threat model. A higher value of Tirreg,trip tends to increase the magnitude of σundersampled,v retrieved from the threat model but to decrease the frequency of irregularity detector trips. The optimal value of Tirreg,trip is determined by its impact on WAAS availability.

5.3 GIVE Parameters

Equations (2), (32), and (40) determine the variance σGIVE,v2 of the overbounding distribution upon which the broadcast GIVE at an IGP is based. To evaluate this term for a given IGP, a second inflation factor, Rirreg–dynamic, must first be calculated according to Equation (26). The values that Rirreg–dynamic assumes under nominal conditions, averaged over time, are shown in Figure 9(a). When the ionosphere is locally well-sampled near a fit center IGP, the value of χ2 in Equation (26) is generally small enough to counteract any inflation provided by Rnoise, reducing the overall GIVE inflation factor to unity. Thus, the impact of Rirreg–dynamic is generally negligible unless the measurement noise is unusually large or the ionosphere is disturbed, rendering it significantly greater than unity.

Figure 9(b) shows the time-averaged standard deviation of the delay formal error σIGP under nominal ionospheric conditions, inflated by Rirreg-dynamic2 according to Equation (27). Figure 9(c) shows the corresponding time-averaged values of σundersampled used to augment σGIVE,well-sampled,v2 to form the GIVE variance according to Equation (2). As discussed above, σundersampled is designed to augment the inflated formal error sufficiently so that the probability a threat generates, near an IGP, an ionospheric delay residual exceeding KHMI_GIVEσGIVE,ν/KHMI is less than 2.25X10−8. Note the stark contrast between Figures 9(b) and 9(c). Not only is σIGP smaller than σundersampled over much of the IGP grid, the former is largely uniform whereas the latter grows substantially as IGPs advance toward the edge of WAAS coverage.

Figure 9(d) shows the time-averaged values of σGIVE, the standard deviation of the GIVE error bound under nominal ionospheric conditions as determined by Equation (2). The broadcast GIVE is evaluated by multiplying σGIVE by constants according to Equation (1) and quantizing as described in Section 2. Figure 9(e) shows the time average of the quantized GIVEs actually broadcast to users. The color bar scale in Figure 9(e) is limited to the range [3, 6]. The lower bound of this range is determined by the fact that 3 m is the minimum floor value conservatively assigned to a GIVE at any IGP. The upper bound defines a critical limit affecting precision approach capabilities. Unpublished studies of historical data conducted at the outset of the WAAS program showed that a GIVE of 6 m is sufficient to support LPV and LPV200 in a local region. However, once the unquantized GIVE at a given IGP exceeds 6 m, it is quantized to 15 m (RTCA, 2016), which is too large to support a precision approach to an airport runway located in the ground region beneath that IGP. Thus, Figure 9(e) shows the broadcast GIVE behavior in the dynamic range that supports a precision approach, i.e., 3–6 m. Note that the information about ranging sources transmitted by monitored IGPs supports the overall performance of WAAS even when the IGPs are located over water some distance from the coastline.

6 ANALYSIS OF THE THREATS THAT DETERMINE THE RELEASE 62-CY23 THREAT MODEL

This section examines the circumstances that have influenced the form and content of the Release 62-CY23 threat model. Recall from above that a threat is defined to exist when a vertical delay fit residual |I¯κIκ| violates the inequality in Equation (33). To distinguish the impact of actual threats from that of threats simulated using data deprivation as described in Section 1, we tabulate (1) two-dimensional histograms revealing the distributions of counts for both sets of threats as functions of fit radius Rfit and fit relative centroid metric RCM and (2) the maximum vertical delay fit residuals associated with these threats, also binned as functions of Rfit and RCM. Subsequently, we investigate the specific threats that give rise to critical points in the threat model to elucidate the origins of the critical threats that establish its structure.

FIGURE 9

GIVE parameters at IGPs, averaged over June 3-11, 2022: (a) Rirreg, (b) σIGP in meters, (c) σundersampled in meters, (d) σGIVE in meters, (e) quantized GIVE in meters

For the storm data set used to construct the Release 62-CY23 threat model, Figure 10(a) shows the distribution of ionospheric threats as a function of Rfit and RCM when the tabulation includes both the threats simulated through data deprivation and the threats corresponding to actual delay measurement geometries. This tabulation ignores all threats that were excluded from the construction of the threat model by the ESD, the MSD, and the irregularity detectors (but not those that were excluded by UIVE floor culling). Figure 10(b) shows the results of restricting the tabulation of the threat distribution to only those threats stemming from actual delay measurements. Note that the range of the color bar scale in Figure 10(a) is much larger than that in Figure 10(b). Clearly, the structure of the threat model is largely determined by the simulated measurement geometries.

FIGURE 10

Tabulation of ionospheric threats generated during storms included in the Release 62-CY23 threat model (i.e., delay fit residuals |I¯I| that violate the inequality in Equations (33)), binned as a function of fit radius Rfit and RCM: (a) threat counts with data deprivation, (b) threat counts without data deprivation, (c) maximum vertical delay fit residuals with data deprivation, (d) maximum vertical delay fit residuals without data deprivation

Figures 10(c) and 10(d) show the maximum vertical delay fit residual as a function of Rfit and RCM, respectively, for the same sets of threats included in Figures 10(a) and 10(b). Comparing Figures 10(c) and 10(d) confirms that the worst threats occur for the simulated geometries.

For Release 62-CY23, the addition of four solar-cycle-24 storms resulted in one new critical point in both the threat model’s quiet-time branch and its disturbed-time branch; in both instances, the origin of the point was the storm data from March 6, 2016. The critical points of the quiet-time branch of the Release 62-CY23 threat model are shown in Table A3, listed in ascending order according to the magnitude of σundersampled.

Figures 11(a) and 11(b) display the values of σundersampled at the critical points of, respectively, the quiet-time and disturbed-time branches of the Release 62-CY23 threat model. Notably, the critical points for the two threat model branches match for RCM values less than 0.5; threats excluded from the quiet-time branch by the MSD are found to occur only for distributions of IPPs that are highly non-uniform in the vicinity of the IGP at the fit center. Figures 11(c) and 11(d) relate the critical points of these branches to the corresponding dates of the storms that generated them.

FIGURE 11

Critical points of the quiet-time and disturbed-time branches of the Release 62-CY23 WAAS undersampled ionospheric irregularity threat model

Plots (a), (c), (e), and (g) display data associated with critical points of the quiet-time branch; the data displayed are (a) the value of σundersampled for each critical threat, (c) the date of the observational data generating each critical threat, (e) the WAAS receiver site observing each critical threat, and (g) the data deprivation mask associated with each critical threat. The counterparts of these plots for critical points of the disturbed-time threat model branch are shown in plots (b), (d), (f), and (h), respectively.

It is interesting to note that the majority of critical points are not produced by extreme storms, i.e., those of unusually large magnitude. Only the disturbances of July 16, 2000, and November 20, 2003, are associated with storms whose magnitudes rank among the top five in the data set (Sparks et al., 2022), and these two storms generate critical points only at high values of fit radius Rfit and the fit RCM, where again the distributions of IPPs about the fit center are highly non-uniform. Figures 11(e) and 11(f) match the critical points of these threat model branches to the sites of the receivers that recorded the measurements sampling the threats. The receivers in question tend to be located at northern latitudes within the service volume. Figures 11(g) and 11(h) identify the data deprivation mask used to generate each critical point in these branches. For single-station deprivation, the first number in the mask label is 10, and the second number identifies the station removed (see Table A1); all other mask labels refer to directional station deprivation, where the first number identifies the direction (also indicated by the abbreviation in parentheses corresponding to a cardinal direction) and the second number indicates the number of stations removed.

Figure 12 identifies the locations of the IPPs corresponding to the signals that detect the ionospheric threats associated with the critical points of each threat model branch. Perhaps the surprise here is that there are not more critical points arising from disturbances in the southern part of the service volume, where the ionosphere tends to be more disturbed. Note that Figure 12(a) displays critical point IPPs over Florida. These IPPs all arise from measurements recorded by the WRS in Jacksonville (see Table A3). However, when threats excluded from the tabulation of the quiet-time branch of the threat model by the tripping of the MSD are included in the tabulation of the disturbed-time branch, these IPPs disappear. New critical point IPPs that appear in Figure 12(b) at higher latitude include some corresponding to critical points with values of σundersampled larger than those associated with the Jacksonville IPPs but also with smaller values of the fit radius Rfit and/or the RCM. Applying the overbound to the raw data of the disturbed-time branch now removes the threats associated with the Jacksonville IPPs from the list of measurements that generate critical points.

FIGURE 12

Locations of the IPPs of the raypaths detecting the ionospheric threats that generate each critical point for (a) the quiet-time branch and (b) the disturbed-time branch of the Release 62-CY23 WAAS undersampled ionospheric irregularity threat model

7 SPATIAL DISTRIBUTION OF UNDERSAMPLED IONSPHERIC IRREGULARITY THREATS

This section probes the worst undersampled ionospheric threats to have plagued WAAS over the first 20 years of its operation (Sparks & Altshuler, 2021). The set of worst threats is, by definition, comprised of each fit residual that produces the largest value of σ¯undersampled,κ2 greater than zero (see Equation (35)) in the threat domain of the fit. As discussed in Section 4, responsibility for protecting the WAAS user from such threats is apportioned among the ESD, the MSD, the LIDs, and the undersampled ionospheric irregularity threat model.

Figure 13 shows the spatial dependence of threat IPPs when various categories of IPPs are systematically removed. Results for solar cycle 23 are plotted separately from those of solar cycle 24, as the level of solar activity and, consequently, the magnitudes of geomagnetic storms were much larger in solar cycle 23 than in solar cycle 24. (Additionally, the WAAS receiver network during solar cycle 24 included 13 stations not in operation during the worst storms of solar cycle 23.) In each plot, we include lines of constant geomagnetic latitude (calculated according to a single iteration of the scheme of Bowring (1976)). Mesoscale ionospheric disturbances tend to become more pronounced (especially during daylight hours) as one approaches the geomagnetic equator and auroral latitudes.

Figures 13(a) and 13(b) show the spatial distribution of the IPPs associated with the worst threats to have occurred during each solar cycle – both actual threats (where the search for fit IPPs has examined all IPPs in the threat epoch) and simulated threats (where data deprivation has removed IPPs from the search). Note that the magnitudes of the worst threats in solar cycle 23 are considerably larger than those of solar cycle 24. Figures 13(c) and 13(d) remove from Figures 13(a) and 13(b), respectively, the threats that trigger the LID at the nearest IGP. Note that Figure 13(c) strongly resembles Figure 13(a), and that Figures 13(d) and 13(b) also resemble each other. This trend suggests that the formal error calculated when the ionosphere becomes disturbed generally bounds the accompanying increase in the magnitudes of fit residuals or, equivalently, that the tripping of an LID is a relatively rare event.

Figures 13(e) and 13(f) reveal the impact of the ESD. Figure 13(f) is identical to Figure 13(d), since the ESD did not trip during solar cycle 24. Figure 13(e) shows, however, that, had the ESD been operational during solar cycle 23, this detector would have provided protection from many of the worst threats to occur. Figures 13(g) and 13(h) illustrate the impact of the MSD. Since the ionospheric disturbances in solar cycle 24 were rarely large enough to trip this detector, the distribution of threats in Figure 13(h) again resembles those of Figures 13(d) and 13(f). However, the impact that the MSD would have had on solar-cycle-23 threats is discernible. Note that the geographic distributions for solar cycles 23 and 24 in Figures 13(g) and 13(h), respectively, resemble each other much more closely.

Figures 13(i) and 13(j) demonstrate the influence of supplementing threat exclusion based on the ESD, the MSD, and the irregularity detectors with threat exclusion based on UIVE floor culling. With this additional threat exclusion, the tendency of the geographic distributions for solar cycles 23 and 24 to resemble each other becomes even more pronounced.

FIGURE 13

Geographic distribution of IPPs associated with the worst ionospheric threats examined in the construction of the WAAS ionospheric threat model from storm data occurring over the period 2000–2019: worst threats in (a) solar cycle 23 and (b) solar cycle 24; worst threats excluding those mitigated by LIDs in (c) solar cycle 23 and (d) solar cycle 24; worst threats excluding those mitigated by LIDs and the ESD in (e) solar cycle 23 and (f) solar cycle 24; worst threats excluding those mitigated by LIDs, the ESD, and the MSD in (g) solar cycle 23 and (h) solar cycle 24; worst threats excluding those mitigated by LIDs, the ESD, the MSD, and UIVE floor culling in (i) solar cycle 23 and (j) solar cycle 24; the brown curves superimposed on each map indicate lines of constant geomagnetic latitude

Figures 13(i) and 13(j) demonstrate that the irregularity detectors, the ESD, the MSD, and UIVE floor culling tend to mitigate the impact of the more extreme solar events affecting the ionosphere, so that the threat model becomes responsible for providing protection only from ionospheric disturbances whose likelihood of occurrence is relatively independent of the overall magnitude of solar activity. Perhaps surprising is the nearly complete absence, at mid-geomagnetic latitudes in North America, of threats that require mitigation by the ionospheric threat model. These results provide reassurance that the WAAS methodology for mitigating ionospheric threats is robust.

8 IMPACT OF THE IONOSPHERIC THREAT MODEL ON AVAILABILITY

As described in Section 1, one may employ WAAS for aircraft navigation only when the locally computed HPL and VPL are bound by the corresponding HAL or VAL, respectively, for a specified level of service. WAAS availability for a given day at a specified user location is defined to be the fraction of the day that both protection levels remained below the corresponding alert limits. The dominant contribution affecting these protection levels generally comes from the GIVEs broadcast at the IGPs that define the vertices of the cell in which the user’s IPP is located. Since σGIVE,undersampled,v2 is usually much larger than σGIVE,well-sampled,v2 in Equation (2), the values assigned to σGIVE,undersampled,v2 by the ionospheric threat model play a critical role in limiting WAAS availability. We can measure the effects of making changes in a threat model by observing their impact on availability. This section analyzes a sequence of hypothetical and actual threat models to explore how the burden of protecting the user from ionospheric threats has been shared between the ESD, the MSD, the LIDs, and the undersampled ionospheric irregularity threat model.

Figure 14 displays plots of availability for a succession of ionospheric threat models, where responsibility for threat protection is increasingly transferred from the threat model to various detectors. The level of service is set to LPV200, which assigns the HAL and the VAL to 40 m and 35 m, respectively. Figure 14(a) presents a hypothetical threat model constructed to protect against all threats, i.e., the ESD and the MSD are turned off and the irregularity detector threat threshold is set to a sufficiently large value such that no irregularity detector ever trips. Figure 13(b) shows the corresponding availability computed for observations recorded on July 1, 2023, a day characterized by nominal ionospheric conditions. Only in the central region of the United States does availability exceed 99.9%, a testament to how well the formal error can bound the threat fit residuals when the threat IPPs are far from the edge of coverage.

FIGURE 14

Undersampled ionospheric irregularity threat model with various threats excluded and examples of the corresponding coverage for July 1, 2023: (a) no threats excluded, (b) coverage corresponding to (a), (c) excluding threats with LIDs, (d) coverage corresponding to (c), (e) excluding threats with LIDs and the ESD, (f) coverage corresponding to (e), (g) excluding threats with LIDs, the ESD, and the MSD, (h) coverage corresponding to (g), (i) excluding threats with LIDs, the ESD, the MSD, and UIVE floor culling, (j) coverage corresponding to (i)

Figure 14(c) shows how the threat model changes when the irregularity detector threshold is reduced to its operational value of 3.0, and Figure 14(d) displays a corresponding appreciable improvement in availability. Excluding threats from the threat model when the ESD has tripped (Figure 14(e)) results in a dramatic improvement in availability (Figure 14(f)). Figure 14(g) shows the quiet-time branch of the threat model when the MSD becomes operational. The improvement in availability apparent in Figure 14(h) is modest but significant at the edge of the coverage region, especially along the coast of California, where, at some IGPs, the difference in availability has been found to reach 75% (Sparks &Altshuler, 2014). Finally, Figure 14(i) removes from the threat model the threats covered by UIVE floor culling (and the threats whose IPPs lie in a cell where two or more of the IGPs at the cell’s vertices are not monitored). This situation again leads to a significant improvement in availability at the edge of coverage.

9 SUMMARY

WAAS broadcasts messages over North America that permit a user (1) to calibrate ionospheric delays that affect aircraft position estimates derived from GPS signals and (2) to bound reliably the error associated with these position estimates. A prior companion paper (Sparks et al., 2022) described how WAAS, during the first two decades of its operation, has protected the user from estimation errors that arise as a consequence of ionospheric storms. The present paper has reviewed the methodology adopted by WAAS to mitigate ionospheric threats posed by mesoscale disturbances whose magnitudes are smaller than those generated by storms.

WAAS broadcasts an estimate of the ionospheric vertical delay at each IGP and a GIVE that bounds the error in this estimate. When defining the GIVE at a specified IGP, the GIVE monitor inflates the formal error of the vertical delay estimate at the IGP to mitigate threats posed by well-sampled irregularities. The GIVE monitor then augments this bound using a term retrieved from the WAAS undersampled ionospheric irregularity threat model to address the threat posed by each local ionospheric irregularity that is poorly sampled, due either to the irregularity being highly localized or to the distribution of sampling measurement IPPs in its neighborhood being highly non-uniform. Following a missed message, the user applies further inflation to each GIVE to compensate for possible rapid growth of the irregularity over time.

A subsequent companion paper (Sparks et al., 2026) examines the considerations and decisions that have guided the evolution of the WAAS undersampled ionospheric irregularity threat model from its origin in IOC to the present.

HOW TO CITE THIS ARTICLE:

Sparks, L., Altshuler, E., Blanch, J., Walter, T., McCord, E., & Griffin Sanchez, R. (2026). WAAS and the Ionosphere – A Historical Perspective: Mitigating Mesoscale Irregularities. NAVIGATION, 73. https://doi.org/10.33012/navi.757

ACKNOWLEDGMENTS

The authors wish to dedicate this paper to the memory of Dr. Nitin Pandya, an esteemed colleague who is no longer with us. He contributed in many instrumental ways to the work described here.

The research of Lawrence Sparks was performed at the Jet Propulsion Laboratory / California Institute of Technology under contract to the National Aeronautics and Space Administration and the Federal Aviation Administration.

The research of Eric Altshuler was performed at the Sequoia Research Corporation under contract to Zeta Associates Incorporated and the Federal Aviation Administration.

APPENDIX A | HISTORY OF THE WAAS UNDERSAMPLED IONOSPHERIC IRREGULARITY THREAT MODEL

View this table:
TABLE A1

WAAS Receiver Sites

This table identifies the location of each WAAS station and its corresponding index.

Table A2 summarizes how error variances used to evaluate the broadcast GIVEs have evolved since the commissioning of WAAS on July 10, 2003. For each variable in the first column, the second column identifies the WAAS releases for which the corresponding expressions in the third column have served as the variable's definition. The final column identifies the equation in the text that determines this definition. The variables listed in the first column are as follows:

  • σGIVE,v2: the variance of a standard normal (Gaussian) distribution that overbounds the tails of the actual distribution of the fit residual error for estimates of vertical delay near the v-th IGP

  • σionogrid,v2: the variance used to degrade σGIVE,v2 to compensate for possible growth in the level of ionospheric disturbance during the period over which broadcast GIVEs remain in force

  • σGIVE,well-sampled,v2: the contribution to σGIVE,v2 that mitigates estimation error arising from well-sampled irregularities

  • σIGP,V2: the formal error variance σ2 for a vertical delay estimate evaluated at the v-th IGP, inflated to ensure that the probability of broadcasting HMI always falls below a mandatory limit

  • σprocess,v2: the component of the formal error variance for a vertical delay estimate arising from process noise

  • σmeasurement,v2: the component of the formal error variance for a vertical delay estimate arising from measurement noise

  • σvmax corner: the maximum of the standard deviations of the formal error for vertical delay estimates at the four corners of the threat domain of the v-th IGP

  • σv2: the formal error variance for a vertical delay estimate at an IPP location in the neighborhood of the v-th IGP

  • (σdecorrtotal)2: the vertical delay covariance associated with widely separated (uncorrelated) IPPs

  • (σdecorrnominal)2: the vertical delay covariance associated with nearly coincident IPPs

  • σROT2: augmentation of εiono2 ensuring that σionogrid,v2 (see option 2 in Equation (3)) overbounds fit residuals between evaluations of planar fits

  • σGIVE,undersampled,v2: the contribution to σGIVE,v2 that protects the user from estimation error generated by undersampled ionospheric threats

  • σundersampled: the overbound of σundersampledraw, a monotonically increasing function of fit radius and relative centroid metric (appending v to the subscript indicates a value, taken from a tabulation of σundersampled and applied to the v-th IGP)

  • σundersampledraw: the maximum value of σ¯undersampled,κ2 tabulated as a function of the fit radius and relative centroid metric

  • σ¯undersampled,κ2: the amount by which σκ2 fails to bound the square of the vertical delay fit residual at the κ-th IPP, inflated by Kinflate and normalized by Kundersampled2

View this table:
TABLE A2 Evolution of the Estimation Error Variances Used to Define the Broadcast GIVE at Each IGP
View this table:
TABLE A3 Critical Points for the Quiet-Time Branch of the Release 62-CY23 Threat Model, in Ascending Order According to the Magnitude of σundersampled,v

In Table A3, the columns present the date of the threat, the Coordinated Universal Time (UTC), σundersampled,v (in meters at L1), the latitude and longitude (in degrees) of the IGP at which the critical point fit was evaluated, the radius Rfit of the fit domain (in kilometers), the fit IPP relative centroid metric RCM, the fit irregularity metric χirreg2, the site of the observing receiver, the data deprivation mask ID, the deprivation mask domain ID, and the number of IPPs removed "maliciously" (see Sparks et al. (2026)). Integers identifying data deprivation masks are as follows:

0 -> no station deprivation: retain IPPs for all receiver sites

10 -> single-station deprivation: exclude IPPs for a single receiver site

11 -> directional station deprivation: exclude IPPs for southernmost receivers

12 -> directional station deprivation: exclude IPPs for northernmost receivers

13 -> directional station deprivation: exclude IPPs for westernmost receivers

14 -> directional station deprivation: exclude IPPs for easternmost receivers

15-> directional station deprivation: exclude IPPs for southwesternmost receivers

16 -> directional station deprivation: exclude IPPs for northeasternmost receivers

17 -> directional station deprivation: exclude IPPs for southeasternmost receivers

18 -> directional station deprivation: exclude IPPs for northwesternmost receivers

When the Data Mask is 10, the integer representing the Mask Domain identifies the receiver site as given in Table A1; otherwise, the Mask Domain integer identifies the number of receiver sites whose data are excluded in directional station deprivation.

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.

REFERENCES

  1. Abramowitz, M., & Stegun, I. E. (Eds.). (1972). Handbook of mathematical functions with formulas, graphs, and mathematical tables. National Bureau of Standards Applied Mathematics Series No. 55. https://doi.org/10.1115/1.3625776
  2. Altshuler, E., Cormier, D., & Go, H. (2002, September). Improvements to the WAAS ionospheric algorithms. In Proceedings of the 15th International Technical Meeting of the Satellite Division of the Institute of Navigation (ION GPS 2002), pp. 22562261. https://www.ion.org/publications/abstract.cfm?articleID=2245
  3. Altshuler, E., Fries, R., & Sparks, L. (2001, September). The WAAS ionospheric spatial threat model. In Proceedings of the 14th International Technical Meeting of the Satellite Division of the Institute of Navigation (ION GPS 2001), pp. 24632467. https://www.ion.org/publications/abstract.cfm?articleID=1921
  4. Blanch, J. (2002, September). An ionosphere estimation algorithm for WAAS based on Kriging. In Proceedings of the 15th International Technical Meeting of the Satellite Division of the Institute of Navigation (ION GPS 2002), pp. 816823. https://www.ion.org/publications/abstract.cfm?articleID=2085
  5. Blanch, J. (2003). Using Kriging to bound satellite ranging errors due to the ionosphere [Doctoral dissertation, Stanford University]. https://web.stanford.edu/group/scpnt/gpslab/pubs/theses/JuanBlanchThesis03.pdf
  6. Blanch, J., Walter, T., & Enge, P. (2003a, January). Adapting Kriging to the WAAS MOPS ionospheric grid. In Proceedings of the 2003 National Technical Meeting of the Institute of Navigation, pp. 848853. https://www.ion.org/publications/abstract.cfm?articleID=3831
  7. Blanch, J., Walter, T., & Enge, P. (2003b, January). Measurement noise versus process noise in ionosphere estimation for WAAS. In Proceedings of the 2003 National Technical Meeting of the Institute of Navigation, pp. 854860. https://www.ion.org/publications/abstract.cfm?articleID=3832
  8. Bowring, B. (1976). Transformation from spatial to geographical coordinates. Survey Review, 23(181), 323327. https://doi.org/10.1179/003962676791280626
  9. Hansen, A., Blanch, J., Walter, T., & Enge, P. (2000a, September). Ionospheric correlation analysis for WAAS: Quiet and stormy. In Proceedings of the 13th International Technical Meeting of the Satellite Division of the Institute of Navigation (ION GPS 2000), pp. 634642. https://www.ion.org/publications/abstract.cfm?articleID=1458
  10. Hansen, A., Peterson, E., Walter, T., & Enge, P. (2000b, January). Correlation structure of ionospheric estimation and correction for WAAS. In Proceedings of the 2000 National Technical Meeting of the Institute of Navigation, pp. 454463. https://www.ion.org/publications/abstract.cfm?articleID=55
  11. RTCA Special Committee 159. (2016). Minimum operational performance standards (MOPS) for Global Positioning System/satellite-based augmentation system airborne equipment (Technical Report No. RTCA/DO-229). RTCA.
  12. Sparks, L. (2012, February). Addressing the effects of space weather on airline navigation in Guidance and Control 2012, Advances in the Astronautical Sciences, Vol. 144, edited by Michael L. Osborne. In Proceedings of the 35th Annual Rocky Mountain Section Guidance and Control Conference of the American Astronautical Society.
  13. Sparks, L., & Altshuler, E. (2014, September). Improving WAAS availability along the coast of California. In Proceedings of the 27th International Technical Meeting of the Satellite Division of the Institute of Navigation (ION GNSS+ 2014), pp. 32993311. https://www.ion.org/publications/abstract.cfm?articleID=12452
  14. Sparks, L., & Altshuler, E. (2021, September). The spatial distribution of ionospheric threats to WAAS integrity, 2000 – 2019: A systematic analysis. In Proceedings of the 34th International Technical Meeting of the Satellite Division of the Institute of Navigation (ION GNSS+ 2021), pp. 39323944. https://doi.org/10.33012/2021.18062
  15. Sparks, L., Altshuler, E., Blanch, J., Walter, T., McCord, E., & Sanchez, R. G. (2026). WAAS and the ionosphere – a historical perspective: Threat model evolution. NAVIGATION, 73. https://doi.org/10.33012/navi.758
  16. Sparks, L., Altshuler, E., Pandya, N., Blanch, J., & Walter, T. (2022). WAAS and the ionosphere – a historical perspective: Monitoring storms. NAVIGATION, 69(1). https://doi.org/10.33012/navi.503
  17. Sparks, L., Blanch, J., & Pandya, N. (2010, September). Kriging as a means of improving WAAS availability. In Proceedings of the 23rd International Technical Meeting of the Satellite Division of the Institute of Navigation (ION GNSS 2010), pp. 20132020. https://www.ion.org/publications/abstract.cfm?articleID=9314
  18. Sparks, L., Blanch, J., & Pandya, N. (2011a). Estimating ionospheric delay using Kriging: 1. Methodology. Radio Science, 46(6). https://doi.org/10.1029/2011RS004667
  19. Sparks, L., Blanch, J., & Pandya, N. (2011b). Estimating ionospheric delay using Kriging: 2. Impact on satellite-based augmentation system availability. Radio Science, 46(6). https://doi.org/10.1029/2011RS004781
  20. Sparks, L., Komjathy, A., & Mannucci, A. J. (2002, May). Sudden ionospheric delay decorrelation and its impact on WAAS. In Proceedings of the 10th International Ionospheric Effects Symposium (IES2002), pp. 115122.
  21. Sparks, L., Pi, X., Mannucci, A. J., Walter, T., Blanch, J., Hansen, A., Enge, P., Altshuler, E., & Fries, R. (2001, June). The WAAS ionospheric threat model. In Proceedings of the International Beacon Satellite Symposium 2001.
  22. Walter, T., & Blanch, J. (2017). Improved user position monitor for WAAS. NAVIGATION, 64(1), 165175. https://doi.org/10.1002/navi.180
  23. Walter, T., Hansen, A., Blanch, J., Enge, P., Mannucci, T., Pi, X., Sparks, L., Iijima, B., El-Arini, B., Lejeune, R., Hagen, M., Altshuler, E., Fries, R., & Chu, A. (2001). Robust detection of ionospheric irregularities. NAVIGATION, 48(2), 89100. https://www.ion.org/publications/abstract.cfm?articleID=102299
  24. Walter, T., Shallberg, K., Altshuler, E., Wanner, W., Harris, C., & Stimmler, R. (2018). WAAS at 15. NAVIGATION, 65(4), 581600. https://doi.org/10.1002/navi.252
Loading
Loading
Loading
Loading