Space Weather Services for Civil Aviation—Challenges and Solutions

Abstract

This paper presents a review on the PECASUS service, which provides advisories on enhanced space weather activity for civil aviation. The advisories are tailored according to the Standards and Recommended Practices of the International Civil Aviation Organization (ICAO). Advisories are disseminated in three impact areas: radiation levels at flight altitudes, GNSS-based navigation and positioning, and HF communication. The review, which is based on the experiences of the authors from two years of running pilot ICAO services, describes empirical models behind PECASUS products and lists ground- and space-based sensors, providing inputs for the models and 24/7 manual monitoring activities. As a concrete example of PECASUS performance, its products for a post-storm ionospheric F2-layer depression event are analyzed in more detail. As PECASUS models are particularly tailored to describe F2-layer thinning, they reproduce observations more accurately than the International Reference Ionosphere model (IRI(STORM)), but, on the other hand, it is recognized that the service performance is much affected by the coverage of its input data. Therefore, more efforts will be directed toward systematic measuring of the availability, timeliness and quality of the data provision in the next steps of the service development.

1. Introduction

The Sun controls the conditions in near-Earth space and in the upper parts of the atmosphere, just as it controls tropospheric weather and, thus, our everyday life. The concept of space weather is used to describe solar driven processes in Earth’s magnetosphere–ionosphere–thermosphere system, which can occasionally challenge man-made technologies both in space and on ground [1]. Today, aviation is in many ways dependent on technology that is prone to space weather storms, causing strong spatio-temporal variations in the geospace environment. The most obvious examples of such technologies are global navigation satellite systems (GNSS) and over-the-horizon (OTH) communication by radio waves in the high-frequency range (HF; 3–30 MHz). Examples of space weather events that have had undisputed impacts on aviation procedures are the storms on 23 January 2012 and on 3–13 September 2017. During the former storm, some cross polar flights were rerouted due to uncertainties in HF communication [2], while during the latter storm, a significant drop in the availability of the European GNSS satellite-based augmentation system took place [3]. Space weather storms can also disturb activities at airports. An exceptionally strong solar burst at radio wavelengths around 1 GHz disturbed secondary air traffic control radars at several European airports on 4 November 2015 [4].

The International Civil Aviation Organization (ICAO) develops principles and techniques for enhanced safety in aviation. Specifications of ICAO-endorsed systems for risk mitigation are given in the organization’s Standards and Recommended Practices (SARPs) that are regularly updated with a three-year cycle. Space weather centers and advisories were introduced in the SARPs for Meteorological Service in 2018. The advisories provide information on space weather phenomena that can affect HF radio communication, communication via satellite, GNSS-based navigation and surveillance systems, or pose radiation risk to aircraft occupants [5]. A space weather center with its back-up arrangements is expected to monitor space weather by 24/7 operations with similar reliability, availability and maintainability requirements, as used in ICAO-endorsed meteorological services.

The impact areas of ICAO interest are related to a multitude of phenomena in the complicated system of solar–terrestrial interactions. Solar flares in X-rays and extreme ultraviolet (EUV) affect directly the radio wave propagation conditions in the ionosphere [6]. Consequently, HF communication problems can appear as changes in the optimal communication frequency or as signal fading. GNSS-based services may suffer from scintillation or unexpected group delays in the signal. Solar X-ray flares and EUV eruptions are often accompanied by bursts of energetic particles [7], which can strengthen ionospheric disturbances and, in extreme cases, also enhance radiation levels at flight altitudes. Solar radiation and particle bursts require prompt reactions from space weather centers as they reach the Earth’s atmosphere with delay times of minutes to hours after their initiation. Besides direct solar bursts, the centers need to have the capability of estimating the time evolution and effectiveness of coronal mass ejections (CMEs, ref. [8]) and other structures in the solar wind (e.g., high speed streams—HSS [9,10]). These structures have typical propagation times of 2–5 days [11] and, depending on their magnetic topology, they may have the potential to trigger geomagnetic storms [12]. During, and occasionally also some days after a geomagnetic storm, the magnetosphere–ionosphere–thermosphere system can generate similar risks for aviation as direct solar radiation and particle bursts do. CME/HSS–driven storms challenge space weather centers with their great variability in spatial extent and duration.

The Manual on Space Weather Information in Support of Air Navigation published by ICAO [13] specifies a set of parameters to be monitored by space weather centers, and provides also thresholds of advisory issuances for two intensity levels: moderate and severe. The parameters and related thresholds that are listed in the first edition of the Manual are summarized in Table 1. Potential space weather impacts on GNSS-based services should be continuously observed with indices that characterize scintillation in the GNSS signal phase and amplitude and with the vertical total electron content (TEC) of the ionosphere. Spatio-temporal variations in TEC can cause uncertainties in the GNSS positioning performance [14]. The statutory procedures of aircrew radiation risk monitoring should be supported with estimates of effective dose rates at a multitude of different flight levels. Space weather–related problems in HF-communication can be triggered by four different types of events: absorption at auroral latitudes, polar cap absorption (PCA), solar X-ray flares and post-storm depression. The three first events are associated with electron density enhancements in the ionospheric D-layer, causing HF signal fadings. In the fourth case, the capability of the ionospheric F-layer to reflect radio signals decreases, due to significant drops in its electron density. The parameters to be followed are the level of global geomagnetic activity (for auroral absorption) as characterized by the Kp-index, attenuation in cosmic noise in the HF regime (for PCA), solar flux in the range of 0.1–0.8 nm (for X-ray flares), and maximum usable frequency (MUF, for post-storm depression).

Information on enhanced space weather activity causing potential problems in the above described impact areas is disseminated by advisories with a standardized structure. An example of such an advisory in its basic form is given in Figure 1. Besides some standard information (time stamp, service provider, advisory numbering, etc.), the advisory describes the space weather impact by its type, intensity and spatial extent. At the time of this writing, the extent is characterized by using a geographic grid with a resolution of 30° in latitude and 15° in longitude (hereafter called the ICAO standard grid). The advisory template includes also rows for forecasts with lead times of 6, 12, 18 and 24 h.

Three global space weather centers were found to be compliant with the ICAO specifications and were asked to start official operations in November 2019. PECASUS is one of these three global centers, and is operated as a joint effort of research institutes from ten ICAO member countries. In this article, we provide a comprehensive review on the instrumentation and data analysis methodologies that are used in PECASUS in order to ensure timely issuance of validated ICAO advisories. The review is based on the authors’ experiences gathered during years 2019–2021, when the PECASUS system was run in its pilot phase for ICAO. The core instruments and the data analysis methods behind PECASUS services are described in Section 2. Section 3 and Section 4 provide short overviews of the PECASUS operational and self-evaluation activities. Section 5 showcases PECASUS performance during a post-storm depression situation, where the ICAO advisory thresholds were temporarily exceeded. Finally, Section 6 summarizes the lessons learned by the PECASUS team during its first years of operation and presents some prospects for future improvements in ICAO space weather services.

2. Observations and Methods

In the following subsections, the PECASUS approaches to monitor space weather activity in the three impact areas of ICAO interest—GNSS, radiation, and HF communication—are described. For each impact area, we first give some general information about space plasma processes that are behind the harmful impacts, and about some standard methods used in the space weather community to observe those processes. Subsequently, we describe the algorithms that PECASUS uses in its near-real-time services.

2.1. GNSS

In the impact area of GNSS, ICAO advisories are disseminated in two subcategories: in the cases of enhanced scintillation (Section 2.1.1) and of enhanced TEC (Section 2.1.2).

2.1.1. Scintillation

Small-scale irregularities in the ionospheric electron density cause scintillation in GNSS signals. During space weather storms, enhanced signal phase or amplitude scintillation can complicate the signal reception, e.g., by loss-of-lock events [15]. At high latitudes, phase scintillation is more typical than amplitude scintillation, and the disturbances are often related to auroras and so-called polar cap patches with localized gradients in electron density [16]. At equatorial latitudes, plasma bubbles building at local sunset cause amplitude scintillation [17]. Scintillation is monitored with specific GNSS receivers that have higher sampling rates than the standard receivers. As the coverage of scintillation, receivers are limited particularly in global scales; actual scintillation measurements are often supported with a proxy called rate of TEC (ROT) derived from networks of standard GNSS receivers with better coverage. The ROT index (ROTI) is retrieved per each receiver–satellite link, according to the formula of [18]: $ROTI=\sqrt{<RO{T}^2>-<ROT{>}^2}$ in which $ROT=dTEC/dt$ and indicates the temporal average calculated, e.g., over 5 min. As time derivatives are used, the challenging task of robust TEC calibration for absolute values can be avoided.

In PECASUS services, scintillation is monitored with high-cadence GNSS receivers (up to 50 Hz) whose carrier signal at L1 frequency (1575.42 MHz) is used to derive standard deviations in one-minute time windows. The amplitude scintillation index $S_4$ is normalized with the signal amplitude, and the phase scintillation index ${\sigma }_{\varphi }$ is derived from the de-trended signal and given in radians. The geographic locations of the PECASUS scintillation receiver stations are shown in Figure 2.

2.1.2. Total Electron Content

Besides scintillation, the ionosphere affects GNSS signal propagation by its bulk electron content. The signal gets delayed, which causes some error in its range used as input in the positioning procedures. In the first approximation, the range error (d) is directly proportional to the integrated electron density, along the signal path in the ionosphere, i.e., $d=K{f}^{-2}\int N_e\left(s\right)ds=K{f}^{-2}STEC$, where K is a constant 40.3 m${}^3$s${}^{-2}$, f is the signal frequency, $N_e\left(s\right)$ is the electron density, and s defines the signal path. $STEC$, i.e., slant total electron content, can be converted to the vertical total electron content ($VTEC$) with the following formula, where the ionosphere is described with a single layer approximation [14]:

(1)$VTEC=STEC{(1-{\left(\frac{R_{e}cosϵ}{R_e+h_{sp}}\right)}^2)}^{1/2},$

where $R_e$ is the Earth’s radius, $ϵ$ is the signal elevation angle, and $h_{sp}$ is the altitude of the single-layer ionosphere (typically 400 km). A common standard is to express $VTEC$ values in TEC units (TECUs), with 1 TECU corresponding to 10${}^{16}$m${}^{-2}$ (i.e., 16.2 cm in the range error for GPS L1 signal). While ionospheric delay can be reduced from services using dual-frequency GNSS receivers, applications relying on single frequency information need additional support from ionospheric models for adequate positioning performance. Such applications are still used in aviation and, therefore, estimates on global $N_e$ are critical to enable safe traffic management and landing procedures. Statistical ionospheric models of volumetric $N_e$, such as NeQuick [19] or $VTEC$, such as Klobuchar [20], provide convenient ways to provide such information as being embedded in the GNSS signals. It is obvious, however, that during space weather storms, range errors cannot be estimated reliably with statistical models. For improved error estimates, statistical models can be supported by TEC information from ground-based networks of dual-frequency GNSS receivers.

For redundancy reasons, PECASUS hosts two parallel systems for global TEC monitoring: one developed under the German Ionosphere Monitoring and Prediction Center (IMPC; https://impc.dlr.de/, accessed on 27 August 2021) and the other developed under the Italian electronic space weather upper atmosphere system (eSWua; https://doi.org/10.13127/eswua/tec, accessed on 16 September 2021). Both systems use as inputs freely accessible near-real-time (NRT) GNSS receiver data from global networks, such as EUREF (https://www.epncb.oma.be/, accessed on 27 August 2021) and the International GNSS Service (IGS, https://www.igs.org/, accessed on 27 August 2021). The resulting TEC maps have spatial resolutions of 2.5° in latitude and 5.0° in longitude and are updated once every 5 min by IMPC and once every 15 min by eSWua. Besides global maps, both systems also generate higher resolution TEC maps for the European region. These maps have a time resolution of 5 min, and the spatial resolution is 2° in latitude and longitude in the IMPC version and 0.5° in the eSWua version. For fine-scale investigations, eSWua is able to produce TEC maps over the Mediterranean territory with a spatial resolution of 0.1° and with a time resolution of 10 min. The exceptionally high resolution is achieved by the combination of IGS and RING (Rete Integrata Nazionale GPS; http://ring.gm.ingv.it/, accessed on 27 August 2021) GNSS receiver networks. Details about STEC calibration (removing hardware based biases and ambiguity in the carrier phase observable) for the IMPC TEC maps can be found in [14] and for the eSWua maps in [21,22,23].

In the IMPC approach, calibrated $STEC$ values are converted to VTECs with Equation (1) and are combined with the empirical Neustrelitz TEC model (NTCM) [24] describing average TEC conditions. NTCM represents the spatio-temporal variations of TEC based on a polynomial consisting of linear terms corresponding to the diurnal and semi-diurnal variation, the annual and semi-annual variation, dependence on geomagnetic latitude, solar zenith angle and dependence on solar activity coefficients, multiplied by factors defined by linear least squares with high-precision data from the IGS long-term archives. In NRT operations, a streamlined version of the model is used that ignores seasonal and solar cycle dependencies, but allows rapid model adjustment according to the current day’s GNSS data and thus, generation of TEC maps in near real time. For this purpose, the NTCM model coefficients and the corresponding TEC calibration of all input data streams are performed on the basis of a weighted least squares solution every 24 h, using the complete data set of the previous day, recursively adjusted to the TEC observations of the last 5 min.

Since ground-based GNSS data are often unevenly distributed, the inclusion of the NTCM in the TEC reconstruction helps to overcome such data gaps, which occur mainly over the oceans. Therefore, the assimilation uses observed VTEC values at their ionospheric piercing points to further fit the TEC values of the NTCM model with two-dimensional Gaussian weighting factors, depending on the distances between the piercing points and the model grid points. The final results represent measured TEC values near the measurement points, while slightly modified model values are provided at larger distances from the measurements.

In the eSWua approach, $VTEC$ values at their piercing points are combined with the NeQuick model [19]. As the first step, data ingestion includes the calculation of differences between the observed $VTEC$s and those by NeQuick at the corresponding piercing points. Subsequently, the differences are interpolated by kriging in order to create a global map of differences that is matching with the observations. Finally, the TEC nowcast is achieved by summing the map of differences with the TEC map by NeQuick. The second step of this process is under development, as standard kriging methods assume an isotropic distribution of data points, which obviously is not true for GNSS receiver networks. This deficiency can be partly overcome by handling differences instead of absolute TEC values, but as a further improvement the algorithm will, in its next release, manage the anisotropy by creating a variogram surface (for more details, see [25]).

2.2. Radiation

Earth’s atmosphere and magnetic field protect us from the detrimental impacts of galactic cosmic rays (GCR) and solar proton events (SPEs). The shielding, however, is partial. According to the Lorentz force law, energetic charged particles propagating toward the Earth get deviated from their original trajectory by the geomagnetic field component parallel to the Earth’s surface. As this component gets smaller at higher latitudes, also the shielding by the geomagnetic field becomes weaker there. Therefore, assessment of effective dose rates for cross polar flights is of particular interest to airliners. The shielding effect by the atmosphere occurs due to collisions that slow down penetrating particles, but on the other hand, they also generate showers of secondary particles that can enhance radiation levels at flight altitudes. Secondary showers are ubiquitous in Earth’s atmosphere, due to the always-present GCR. Modeling the atmospheric radiation level due to GCR is easier than estimating the impact by SPEs because variations in GCR are slow, and the spectrum of primary particles is well known. On the other hand, additional radiation by SPEs is characterized by its abrupt appearance, short duration and energy spectra, varying significantly from case to case. In the context of very strong SPEs, secondary showers are detected as ground level enhancements (GLEs), i.e., enhanced count rates in ground-based neutron monitors.

In PECASUS, effective doses ($\mu$Sievert/hour) at flight altitudes are estimated with the AVIDOS tool [26], which can handle both GCR and SPE impacts on radiation levels. The contribution of GCR on radiation levels at flight altitudes is described elsewhere [27]. For the impact of SPEs, the approach uses, as its primary data source, NRT count rates from the neutron monitor in Oulu (Finland, 65.05°N 25.37°E; [28]) with support by the ANeMoS GLE alerting service [29]. In the case of a GLE, AVIDOS first constructs two estimates for the primary proton energy spectrum—one soft and one hard version. In the method, the energy spectra of the primary protons are assumed to have a power law form in particle rigidity (particle momentum divided by its charge):

(2)$J_{SCR}(R,t)=A_0\left(t\right)R^{-\gamma \left(t\right)},$

where R is the primary proton rigidity expressed in GV, t is time, $A_0\left(t\right)$ is amplitude of the proton flux (in particles/(cm${}^2·$ s · sr · GV)), and $\gamma \left(t\right)$ is the spectral index (softness of the spectrum). $A_0\left(t\right)$ is derived from the neutron monitor data. The relationship between the primary spectra and neutron monitor count rates, $N(R,t,z)$, can be described by the following formula [30]:

(3)$N(R,t,z)={\int }_{R_{min}}^{\infty }S(R,z)J_{SCR}(R,t)dR,$

where $R_{min}$ is the minimum considered rigidity ($R_{min}=0.8$ GV for Oulu neutron monitor station), z is the atmospheric depth of the station, and $S(R,z)$ is the so-called yield function for the considered neutron monitor. For the yield function, AVIDOS uses the parametrization developed by [31]:

(4)$log\left(S(R,z)\right)=\sum _{m,n=0}^3{C}_{mn}{z}^m{\left(log\left(R\right)\right)}^n,$

where $C_{mn}$ are parametrization coefficients based on Monte Carlo simulations of neutron monitor responses (see [31] for more details and exact values of $C_{mn}$) and $log$ stands for decimal logarithm.

When the two estimates (hard and soft) of the proton spectrum are ready, the spectra are folded with results of Monte Carlo simulations of radiation transport through the atmosphere performed with the Geant4 application PLANETOCOSMICS ([32,33]), which results in a large set of secondary particle fluences. Finally, to obtain the radiation protection quantity effective dose, the particle fluences are folded with appropriate fluence-to-dose conversion coefficients and are summed up. In AVIDOS, two types of such conversion coefficients can be used: that based on the older recommendations of the International Commission of Radiological Protection (ICRP) from 1990 [34] or that based on the more recent ones from 2007 [35].

AVIDOS performance was validated in an international comparison campaign performed by EURADOS’ (the European Radiation Dosimetry Group, www.eurados.org, accessed on 27 August 2021) Working Group 11 on high energy radiation fields. In this exercise, nine codes estimating radiation exposure due to GLE were employed to calculate radiation doses for GLE 42 (September 1989) and GLE69 (January 2005) for three representative flight routes. The survey showed that estimated flight route doses may differ by up to a factor of ten [36]. Further research and verification by on-board aeroplane measurements are needed for more general level understanding on the model discrepancies. The NRT version of AVIDOS tackles this issue partly by the set-up of two estimates for the primary proton population spectra (hard and soft ones). Consequently, as output, the code yields maximum and minimum boundaries of the expected effective dose. By deducing the minimum–maximum range instead of a single value, the code gives implicit information about its performance (with the assumption that large ranges imply increased uncertainty in the estimates).

Complementary to ground observations, PECASUS also makes use of NRT space data to generate SPE forecasts, which are internally used to alert the operator for increased risks of radiation storms and possible impacts on aviation. Two predictor tools are adopted: HESPERIA UMASEP-500 ([37,38]) and COMESEP SEP Forecast ([39]). The HESPERIA UMASEP-500 tool produces real-time predictions of the intensity of the first hour of >500 MeV SPEs near-Earth on the basis of a time-delayed correlation analysis between the GOES soft X-ray flux and the GOES differential and integral proton fluxes in the energy channels above 500 MeV. When the correlation estimation surpasses a threshold, and the associated flare is greater than a specific X-ray peak flux, a >500 MeV SPE warning is issued before (up to tens of minutes) detection by neutron monitors. Particles of that high energy can easily reach the top of the atmosphere within minutes and create enhanced radiation levels at flight altitudes.

The COMESEP alert system provides (in addition to geomagnetic storm alerts) SPE alerts in the lower energy range, i.e., >10 and >60 MeV. The SPE alerts are triggered by the automated detection of solar activity (>M1 solar X-ray flares). Based on a statistical analysis of historical data correlating flare and CME observations with SPE characteristics, a risk level estimate is derived from NRT data based on the expected occurrence probability and severity. The alerts are continuously updated along the timely evolution of event(s) in the case that new information on the solar eruption is received. As explained in the next paragraph on HF communication, these solar energetic particles mainly spiral-in into the polar regions, causing polar cap absorption events. In the PECASUS service, the COMESEP system is used as an advanced warning for increased risk of radiation impacts confined to the transpolar regions.

2.3. HF Communication

In the impact area of HF communication, ICAO advisories are disseminated in four subcategories: in the cases of shortwave fadeout, polar cap absorption, auroral absorption (all described below in Section 2.3.1) and of post-storm depression (Section 2.3.2).

2.3.1. Shortwave Fadeout, Polar Cap and Auroral Absorption

Solar X-ray flares and energetic particle precipitation increase the electron density in the lower part of the ionosphere (D-layer). While the shortwave fadeout impact by X-ray flares lasts from some minutes to a few hours, an absorption event due to energetic particle precipitation can last for days. In both cases, the interaction between radio waves and D-layer plasma enhances when compared to normal conditions. As this interaction is affected strongly by the ionization of neutral atmospheric constituents and by electron collisions with atmospheric neutral particles, radio waves are absorbed instead of reflected, like in the case of F-layer interaction. D-layer absorption can be a problem particularly in the low end of the HF-spectrum because electron density in the D-layer is too low to affect higher frequencies, even during disturbed conditions. A solar X-ray flare causing D-layer absorption is known as a shortwave fadeout, a legacy term from radio wave applications using frequencies lower than HF. Airliners need to be aware of shortwave fadeouts because they challenge HF-communications in the sunlit ionosphere. At the time of this writing, PECASUS advisories on fadeouts are based on X-ray observations by the XRS instrument (0.1–0.8 nm; [40]) of the geostationary GOES-16 satellite—maintained by the U.S. National Oceanic and Atmospheric Administration (NOAA).

The D-RAP model (https://www.swpc.noaa.gov/content/global-d-region-absorption-prediction-documentation, accessed on 27 August 2021) is widely used in space weather services—including PECASUS—to estimate the spatial extent and frequency impact of D-layer absorption, due to solar X-ray flares and SPEs. As inputs, the model uses X-ray and proton fluxes by the GOES-16 satellite and the $Kp$ index (more details about $Kp$ are given later in this section). For X-ray flares, the model provides an estimate of the highest affected frequency in the dayside ionosphere and spatial extent of the impacted area, the latter of which is used as supporting information in ICAO services. However, the intensity of impact (MOD/SEV) is determined from direct X-ray flux measurements (c.f. Table 1). For SPE-driven polar cap absorption, D-RAP gives estimates of 30 MHz signal attenuation in decibels (dB) based partly on empirical relationships as described by [41]. The relationship between proton flux and day- and nightside absorption levels (in dB) can be expressed as follows:

(5)$A_d{=0.155(J(E>Ec\left(InvLat\right))}^{1/2},$

and

(6)$A_n{=0.020(J(E>Ec\left(InvLat\right))}^{1/2},$

where $A_d$ and $A_n$ are the dayside and nightside absorption levels (with sun elevation angles >10° and <−10°, respectively), J(E>Ec(InvLat)) is the proton flux (in (cm${}^2$· s · sr)${}^{-1}$) for energies larger than $E_c$. $E_c$ (in MeV) is the cut-off energy above which protons are able to overcome the geomagnetic shield and, thus, precipitate into the ionosphere at the given location, whose position is characterized with the McIlwain L-parameter [42]. In the twilight zone, absorption levels are estimated with linear bi-interpolation from $A_d$ and $A_n$ according to the sun elevation angle and InvLat.

D-RAP performance in the polar cap absorption assessment was evaluated with riometer data from six high-latitude stations at different longitudes for 13 space weather storms during the years 1998–2005 [43]. The survey reveals that the model performs reasonably well for high-latitude stations and can provide useful qualitative information also during very disturbed periods. However, for more accurate quantitative estimates, more work is still needed, as D-RAP seems to under-estimate absorption events, particularly in the European sector.

The D-RAP approach cannot capture energetic electron precipitation of magnetospheric origin, which can cause HF-signal absorption at auroral latitudes. Energy spectra and spatial distribution of electrons is linked with magnetospheric dynamics in such a complicated way that they cannot be characterized using the single-point measurement approach (at geostationary satellites) used for proton fluxes. This complication can be partly overcome with the knowledge that the energetic electron precipitation being at high levels is particularly common during geomagnetic storms. Therefore, ICAO has chosen to use the $Kp$-index as the indicator of enhanced auroral absorption in its advisories. $Kp$-index is composed from measurements of 13 well-established geomagnetic observatories located at sub-auroral latitudes. It describes geomagnetic activity as the deviation from the regular daily variation with a three-hour resolution [44] and with a scale from 0 to 9 [45]. Statistical analysis of $Kp$-archives since 1933 reveals that the occurrence of high enough $Kp$-values for ICAO advisories on auroral absorption varies, according to the solar cycle being on average seven per year during years around the solar maxima and once per year around the solar minima. Typical durations of geomagnetic storms with peak $Kp\ge 8$ are in the range of 1–3 days.

2.3.2. Post-Storm Depression

Geomagnetic storms are associated with post-storm depression of electron density in the ionospheric F2 layer. Depressions are usually more intense, longer-living and widespread during the main and recovery phases of intense geomagnetic storms. According to a common recommendation between the ICAO space weather centers, post-storm depression advisories may be issued up to 96 h after $Kp\ge 6$ has been reached. With decreased electron density, the ionosphere’s capability to reflect radio waves lowers, in particular, at the higher end of the HF-regime. Post-storm depression is caused by changes in the thermosphere–ionosphere interaction as a consequence of enhanced Joule heating and particle precipitation at auroral latitudes [46]. Heating at high latitudes causes up-welling, where neutral molecules are lifted to typical F2 layer altitudes (300–400 km) and cause the loss of electrons due enhanced recombination. Heating generates also horizontal pressure gradients and advection of molecule-rich plasma from auroral latitudes to lower latitudes. The spatial distribution and intensity of a post-storm depression depends strongly on season and local time. During summer time, the background neutral wind pattern has an equatorial component, which supports equatorward propagation of molecule-rich plasma fronts. Therefore, post-storm depression events penetrate to lower latitudes and are stronger there during summer than during winter months. Furthermore, observers in early local time morning hours at mid-latitudes experience more severe depression effects than those in the afternoon sector. The recovery of ionosphere–thermosphere conditions back to normal may take 1–2 days after a strong post-storm depression event.

Ionosonde Measurements

Post-storm depressions can be detected with ionosondes, which are vertically sounding pulsed radars operated in a range of HF frequencies (typically 1–20 MHz). According to this measurement concept, the ionosphere is probed with gradually increasing frequencies. Signals transmitted from the ground are reflected back from the ionosphere at the altitude where the ionospheric plasma frequency ($f_p,f_p\approx 9\sqrt{N_e}$) is greater than or equal to the signal frequency. The highest frequency that is reflected back gives a good estimate for the maximum $N_e$ in the ionosphere. Approximation of the reflection altitude can be obtained from the signal travel times, with the assumption of the signal propagation speed being that of light. The maximal electron density is typically in the F2 layer. URSI (International Union of Radio Science) specifies the ionospheric parameter codes. $foF2$ is used for the F2 layer plasma frequency peak, and $h^{\prime }F2$ is used for the virtual altitude of $foF2$. The latter is the height obtained directly from the measured radio wave travel times; it can be converted to true height $hpF2$ or $hmF2$ only by using a theoretical fitting curve [47].

The maximum usable frequency listed in the ICAO advisory threshold parameters (c.f. Table 1) can be related with URSI parameters $MUF3000F2$ and $M3000F2$, the former of which is the F2 layer maximum usable frequency for a 3000 km path and the latter of which is $MUF3000F2/foF2$. Obtaining $foF2$ from ionosonde measurements is straightforward, as it can be easily scaled directly from the ionogram. Estimating $MUF3000F2$ or $M3000F2$ requires either fitting the ionogram traces with a semi-empirical curve [47] or knowledge of the true peak height $hmF2$. The latter, in turn, needs to be calculated from the measured $h’F2$, using some theoretical model curve to take into account the delay of the signal due to ionization below the reflection altitude.

The PECASUS Approach

In PECASUS, both $foF2$ and $MUF3000F2$ are used to monitor post-storm depression. As pointed out above, $foF2$ is available as a standard parameter for all ionosonde types, which makes it a robust source for global monitoring. The network of ionosondes that is used in the PECASUS approach at the time of this writing is presented in Figure 3. $foF2$ depression, $\delta foF2$, is estimated for each available ionosonde as follows:

(7)$\delta foF2=-\frac{foF2-<foF2>}{<foF2>}$

where $<foF2>$ is the 30-days’ median value of $foF2$. Positive values of $\delta foF2$ correspond to depression of $foF2$. $foF2$ values from each available ionosonde are checked and compared with the 30-day median. The outputs of the $foF2$ algorithm-$\delta foF2$ values and excesses of the ICAO MOD/SEV thresholds are given in the resolution of the ICAO standard grid (c.f. Section 1) and with a time resolution of 15 min. If, for a specific ionosonde, there are fewer than 20 values for the past 30 days, the 30-day median is estimated, using the NeQuick model [19]. The $\delta foF2$ value from a particular ionosonde is naturally considered to be representative for its own cell in the ICAO standard grid, but also for the neighboring cells if the central angle between their boundaries and the ionosonde location is less than 10°.

As stormy conditions can, in some cases, challenge reliable ionosonde measurements, the strategy here presented is to apply some conservative rules in the $\delta foF2$ algorithm in order to avoid false warnings. If several ionosondes contribute to $\delta foF2$ for a given cell and all those values are below the ICAO advisory thresholds, then the median of the individual values is used to represent the cell conditions. If only one of the contributing ionosondes is suggesting excess of a ICAO threshold, that measurement is ignored, and the cell obtains the median value from the remaining ionosonde measurements. If more than one of the contributing ionosondes are giving higher $\delta foF2$ than the thresholds, then the median is determined from those ionosondes, and the other ionosondes are ignored. A similar rule is applied also when creating the maps of ICAO threshold excesses. A cell is denoted with yellow or red shading only if two or more ionosondes contributing to its content are showing strong enough depressions for MOD or SEV advisory issuances. The color is determined according to the stronger activity level. Examples of $\delta foF2$ maps created with these principles are shown in Figure 4.

The method used to produce $MUF3000F2$ global maps is composed of two models for $foF2$: global dynamic model (GDMF2) [48] and EUROMAP [49]. $MUF3000F2$ maps are created from $foF2$ maps by using propagation factor $M3000F2$ estimates from the empirical international reference model (IRI; [50]). The space weather impact on MUF is characterized by the ratio between $MUF3000F2$ and background conditions, which come from the empirical model by [51]. Maps are created with a one hour time resolution, based on NRT data and indices. As many of the indices are available also as forecasts, the methodologies described below are suitable for forecasts with some moderate modifications.

GDMF2 is a global empirical model particularly for mid-latitude post-storm depression conditions, which relates $\Delta foF2$ with a neutral composition and temperature taken from any thermospheric model, for instance, NRLMSISE-00 [52]. The impact of solar and geomagnetic activity is taken into account with the solar radio flux at the wavelength of 10.7 cm (F10.7) and with the $Kp$-index converted to $ap$ (for conversion; see, for example, ref. [48]), respectively. $ap$ is a linear variant of $Kp$[53] with values in the range of 0–400. In order to describe the inertia in thermospheric response to enhanced Joule heating, GDMF2 uses instead of $ap$ its time-weighted accumulation as defined by [54]:

(8)$ap(t,\tau )=(1-\tau )\sum _{n=0}^{}{\tau }^{n}ap(t-3n),$

where t is time and $\tau$ is a parameter characterizing the thermospheric response time. Typical values for $\tau$ are around 0.7–0.8, which correspond to response times of 8–13 h [48].

The EUROMAP approach is applied to well-established ionosonde stations that can provide NRT data and have long-term archives of polished products available for the analysis. EUROMAP includes regression models to describe strong negative disturbances ($ap(\tau ,t)>30$) and training models over the previous 28 days. Both these components are specifically tailored for each ionosonde station used in the model. Constructing a storm regression model typically requires a data archive from the station of interest of more than 20 years in order to ensure good enough statistics on storm conditions. The F10.7 is used to characterize the impact of seasonal and 11-year cycle variations of solar activity on $foF2$, while T-index is used to construct the $foF2$ background for each ionospheric station (for more details on T-index derivation, see [55]). The training model is based on the knowledge that $\Delta foF2$ correlates well with its previous value during quiet conditions. Therefore, $\Delta foF2\left(t\right)$ for the ongoing hour can be estimated with the following formula:

(9)$\Delta foF2\left(t\right)=C_0+C_1\Delta foF2(t-1)+C_{2}ap(t,\tau ),$

where $\Delta foF2(t-1)$ is the $\Delta foF2\left(t\right)$ of previous hour and coefficients $C_0$, $C_1$, and $C_2$ are defined with the least squares fit to the data of the previous 28 days (one solar rotation period).

3. PECASUS Operational Activity

The first and central task for PECASUS is to operate a framework capable of monitoring all space weather activity that is relevant for aviation, and to issue and disseminate advisories as the ICAO requirements specify. The PECASUS framework, to which all PECASUS consortium members contribute, consists of a number of central hubs and expert groups.

The main Advisory Dissemination Hub, responsible for the provision of advisories to the aviation community, is hosted by the Finnish Meteorological Institute, Finland (FMI), which is also the consortium lead. The Advisory Consolidation and Production Hub is operated by the Solar Terrestrial Center of Excellence, Belgium (STCE). The STCE operators are responsible for the monitoring of aviation-relevant space weather through a specifically developed web-application tool, which supports them to produce the space weather advisories when the pre-defined conditions are met (see Table 1). The content and format of the draft advisories are checked in the Advisory Dissemination Hub before they are pushed to the Aeronautical Fixed Telecommunication Network. The Met Office, United Kingdom (UKMO) has the role of the Resilience Hub to support the other two hubs in standby mode with the capability of taking over either role on short notice. This ensures the continuity of service in the case of technical problems with the advisory production or with the dissemination processes. The three hubs are manned with over a dozen operators, working in shifts, each equipped with a dedicated training in space weather forecasting and PECASUS operations. The provision of data, observations and models are supported by several expert groups, members of the PECASUS consortium: Seibersdorf Laboratories, Austria (SL), German Aerospace Center, Germany (DLR), Space Research Centre of Polish Academy of Sciences, Poland (SRCPAS), Istituto Nazionale di Geofisica e Vulcanologia, Italy (INGV), Frederick University, Cyprus (FU) and South African National Space Agency (SANSA), UKMO, and STCE. Self-evaluation is conducted by the Royal Netherlands Meteorological Institute, the Netherlands (KNMI).

The official ICAO space weather monitoring that was initiated with three global centers was organized with a bi-weekly rotational principle, where a center serves either as an On Duty Center (ODC), a Primary Backup Center (PBC) or a Secondary Backup Center (SBC). The official handovers between the centers is set on Tuesdays, every two weeks, at 08:00 UTC. While PECASUS serves as the official ODC only for two weeks out of every six weeks, the PECASUS operational services are run continuously, with the only difference being that, during the four weeks of the back-up role, the PECASUS advisories are circulated only among the consortium. At the time of this writing, PECASUS issued only a small number of advisories, as the operational services started close to the solar minimum of activity. Nevertheless, several short-lived events had already been observed and a sequence of advisories had been issued for GNSS scintillation, HF-COM PCA, HF-COM post-storm depression, and HF-COM short-wave fadeout. An example post-storm depression event is discussed in Section 5.

Data and models directly related to the ICAO requirements for issuing advisories are continuously monitored. Such critical inputs (described in Section 2) are provided by the PECASUS expert groups and are centralized at STCE and FMI. For resilience of the continuous operations, redundancy is implemented in the data streams when possible: PECASUS is using and monitoring multiple vertical TEC global maps built independently by DLR and INGV, two D-RAP models run independently by UKMO and the NOAA Space Weather Prediction Center, as well as various MUF and $foF2$ maps provided by SRCPAS and INGV. PECASUS also monitors the following K-indices: $Kp$-index from NOAA and Potsdam, K-index measured in the local station in Dourbes (Belgium), as well as NRT $Kp$-index from NOAA. As a proxy for GOES X-ray measurements, X-ray fluxes at wavelengths of 6–20 nm and 17–80 nm from the LYRA instrument on-board PROBA2 ([56]) can be used. This solution was tested successfully on 1 June 2021, when data interruption in the GOES X-ray flux was observed around 15:00 UTC and 20:00 UTC. However, it is good to recognize that the PROBA2 back-up arrangement is available only for exceptional cases. LYRA data are usually downloaded with several hours of delay and measurements may contain gaps due to missing down-links or calibration campaigns.

The incoming real-time space weather data are monitored partially automatically. Special routines are implemented to alert the operator on duty when one or more of the critical data or model outputs exceeds specific thresholds. Once the operator is alerted, she/he analyzes the data manually to decide whether an advisory needs to be issued or not. The driving strategy in the automatic routine is to minimize the number of missed events, thereby accepting the potential triggering of false alerts that the operator on duty needs to disregard. Table 2 provides the set of criteria used to alert the operator on duty, who, later on, performs further manual checks to analyze the event.

Indirectly related data and models are also intensively used during the operations activities. For example, increases in the GOES-16 proton flux are monitored as possible precursors of HF-COM polar cap absorption and used for surveillance of radiation events in addition to the ground based neutron monitors (see [27,57,58]). Additionally, the solar wind speed, density, and Interplanetary Magnetic Field (IMF) measurements of CMEs and HSS are closely monitored since they are the most relevant parameters in terms of their impact at the Earth, or geoeffectiveness ([59,60]). STCE forecasters use several methods for obtaining the CME parameters. One of them is CACTus ([61]), developed at $STCE$, which detects coronal mass ejections in images sequences from coronagraphs (such as SOHO/LASCO and STEREO/COR2). Furthermore, it allows a real time estimation of the initial speed and angular width of CMEs.

4. Self Evaluation, User Feedback Assessment, and Training

PECASUS self evaluation includes the verification and validation of products, using independent data sources, analyses of anomalous situations, and annual post-processing of PECASUS data archives and error logs. The work contains also coordinated reporting procedures for off-nominal situations. These are evaluated as part of the ISO quality management system, and the total performance is reported to the Meteorological Operations Group (MOG) that is part of the Meteorological Panel of ICAO. Based on the lessons learned, performance will be a subject of continuous attention.

During the forthcoming years, there will be the important but difficult task of post-processing the issued advisories. This includes the validation of nowcasts and verification of forecasts. This will provide insight for the correctness and timeliness of advisories, and, in some cases, may also reveal cases where an advisory has incorrectly not been issued. The challenges for validation are caused by the fact that many space weather data (observations and models) are not yet available in a strictly standardized, quality-controlled way, like meteorological data. An important task for self-evaluation is to find potential systematic discrepancies between actual observations and their estimates based on the methods described in Section 2.

As the aviation world is multifaceted, organizing user feedback is not trivial. We hope, however, to organize user engagement events or present space weather at aviation meetings. This will be more effective in the forthcoming years when the Sun’s new cycle progresses and the number of advisories issued are increased.

Finally, a space weather service needs trained staff, especially operators. PECASUS duty officers receive dedicated training provided by the Space Weather Education Center (SWEC, https://www.stce.be/SWEC, accessed on 27 August 2021) of the STCE. PECASUS duty officers learn to analyze and interpret the PECASUS data and to navigate through the PECASUS dashboard. As a standard, the STCE organizes ’Space Weather Introductory Courses-SWIC’ several times a year. The SWIC includes live lectures and exercises, written course material with references and an evaluation with certificate issued by the STCE. New PECASUS duty officers follow this course to become acquainted with space weather and the STCE space weather bulletins. Tailored courses for aviation end-users are organized on request by the STCE. An online User’s Guide to Space Weather Advisories for Pilots is also available (https://www.stce.be/PECASUS_guide4pilots, accessed on 27 August 2021 or http://pecasus.eu/?page_id=593, accessed on 27 August 2021).

5. PECASUS Advisories on Post-Storm Depression during 28–30 September 2020

The first space weather event with official ICAO advisories took place from 28 September to 1 October 2020. Activity was triggered by high speed solar wind ($\approx 600$ km/s) emanating from a solar coronal hole. High-wind speed combined with a southward component of IMF enhanced the energy input from the solar wind to the magnetosphere, which launched a geomagnetic storm with provisional $Kp=6-$ during the three last UT-hours of 27 September. High $Kp$-values were observed also earlier, around noon of 24 September. The first signs of post-storm depression were observed around 05:30 UT on 28 September by the Australia–Canada–France–Japan consortium (ACFJ), which served as the on-duty center at that time. During the event, ACFJ sent a total of 13 advisories about subsequent depressions in different regions of the globe over a period of 3.5 days. The analysis below is concentrated on the period of 28–30 September UT, when depression was observed in the European sector.

As Figure 5 shows, enhanced geomagnetic activity (with $Kp$-values 4–5) continued still throughout the day on 28 September (UT hours 0–24 in Figure 5), but activity ceased over the next two days (UT hours 24–72 in Figure 5). Over the next three days (28–30 September), moderate post storm depression took place during daytime at all stations. In Nicosia, which is the lowest latitude ionosonde station among those shown, the effect was smallest, with a clear depression only between 5:00 and 11:00 UT on 28 September (UT hours 5–11 in Figure 5). During night times, however, clear depression was observed at all the other stations. The day–night difference in depression intensity is due to the different thermospheric circulation during daytime and nighttime hours. The poleward component in the daytime neutral wind restricts depression to higher latitudes during daytime hours. A comparison between EUROMAP and IRI(STORM) models and observations at Pruhonice and Rome is shown in Figure 6. It is encouraging to notice that during this moderate negative MUF depression event, EUROMAP follows better the observations than the IRI(STORM) model [62].

Figure 7 shows the spatial extent of night time depression on 28 September as characterized by the ratio of $MUF3000F2$ to its background, suggesting ICAO threshold excesses in the Northern mid-latitude band at longitudes ca. from 30° W to 105° E. This result can be compared with the advisory recommendations by the $foF2$ algorithm, which are shown for hours between 19:00 and 24:00 UT in Figure 8. The $foF2$ service recommends advisories on MOD level for the European region, but besides mid-latitudes, high latitudes are impacted. On the other hand, the longitudinal extent of the impacted area in the $foF2$ map is narrower than that in the $MUF3000F2$ map. These differences between the two products can be explained by their different approaches: $foF2$ maps are mostly controlled by the availability of NRT ionosonde data, while $MUF3000F2$ maps reflect information also from long-term statistics. Both aspects are considered to be valuable assets in PECASUS operations.

6. Conclusions and Future Prospects

The present work describes the PECASUS system of data analysis methods and models, which provides advisories on enhanced space weather, compliant with the Standards and Recommended Practices of the International Civil Aviation Organization (ICAO). The present review is based on the authors’ experiences gathered during the years 2019–2021, when the PECASUS system was run in its pilot phase for ICAO. In the ICAO framework, advisories are disseminated in three impact areas: radiation levels at flight altitudes, GNSS based navigation and positioning, and HF communication. Space weather activity levels and thresholds for advisory issuances are characterized by a set of geophysical parameters, which are widely used among the space weather community both in research and services. The PECASUS machinery with its core models, their data inputs and augmenting models are summarized in Table 3.

PECASUS space weather services rely on software solutions and NRT observations that were originally developed for research purposes, and often the development work originated outside the PECASUS consortium. Both these factors naturally pose a risk for the reliability, availability and maintainability of the PECASUS system, which, according to the ICAO requirements, must be above $95\%$. This challenge is partly tackled by duplicating some of the critical products—such as GNSS TEC, scintillation indices and post-storm depression estimates—with more than one independent methodology. The critical elements are typically assets of extensive use by their developers and by other well-established national and international space weather services, such as that in the NOAA Space Weather Prediction Center (SWPC) and in the Space Safety program of the European Space Agency. For example, GOES-16 X-ray flux and GOES-16 Proton flux are crucial external data taken from SWPC that are directly monitored by the PECASUS operators and are also used as input for models, such as the D-RAP model, and AVIDOS maps. Wide use in operational services motivates the continued maintenance and development of these critical products, although in some cases, sustained funding mechanisms supporting them are missing. ICAO’s joining to the stake holders of space weather services will hopefully improve this situation in future.

An important part of PECASUS development is to combine the actual service provision with careful quality management, which includes a clear division of responsibilities, systematic self-evaluation with user feedback collection, and education activities. The responsibility to maintain automated data and product provision routines is shared in PECASUS among the original developers of the algorithms. The 24/7 semi-automatic and manual operations are conducted in the PECASUS Advisory Consolidation Hub, Dissemination Hub and Resilience Hub. Product verification and validation and service performance are part of the self-evaluation activity of PECASUS, as well as the development and setting up of an education program.

The next steps in the PECASUS development will include improvements in its spatial coverage. Particularly in scintillation and post-storm depression, the service performance is constrained by the coverage of input data used in the analysis. In scintillation measurements, information by the station network of 50 Hz receivers (c.f. Figure 2) can be supported with the ROTI proxy (c.f. Section 2.1.1). PECASUS already has routines to create global ROTI maps with time and space resolutions of 15 min and 2.5° in latitude and 5° in longitude. More work is still needed to gain a better understanding of the relationship between ROTI values and ICAO scintillation thresholds, which most likely will vary, according to latitude and level of activity. As improvements in the global TEC mapping algorithms, the PECASUS research teams are investigating possibilities to use direction-based variograms for kriging interpolation and novel methods for differential code bias estimation [63].

For the case of post-storm depressions, no proxies exist to replace actual ionosonde measurements. A way forward would be to add more stations with historical archives to the process that is currently used to create $MUF3000F2$ maps for the European region [64]. The algorithm combines NRT measurements and historical information from different locations by the ordinary kriging method for improved spatial coverage and better resilience against occasional gaps in observations.

The template of ICAO advisories include also entries for space weather forecasts. The forecast lead times (6, 12, 18 and 24 h, c.f. Figure 1), which come as heritage from meteorological services, are, in most cases, too challenging for space weather services with the current physical knowledge of the phenomena that are of interest for ICAO. However, the common understanding among ICAO, airliners and service providers is that, even with nowcasting information alone, the space weather advisories are a welcome step forward in enhancing safety in civil aviation against space weather related risks. It is clear, however, that demand for forecasts will become stronger once the service users become more familiar with the space weather discipline.

PECASUS partners’ work for forecasts of post-storm depressions is a natural starting point to tackle the challenges on space weather forecasts. EUROMAP is already now a working $foF2$ forecast method for the European region based on historical and real-time $foF2$ observations. Severe and moderate MUF depressions related to high geomagnetic activity with forecasts of 1–24 h lead times can be achieved with higher accuracy than those by IRI(STORM) [49], while no advantage over IRI(STORM) is obtained for forecasts with lead times greater than 3 h during low geomagnetic activity. The next step will be to achieve forecasts of moderate and severe MUF depressions during low geomagnetic activity for all lead times with better accuracy, and to use this model globally for all ionosonde stations providing NRT data. NRT $MUF3000F2$ maps over Europe are also generated as a PECASUS service [64] and planned to be improved in several ways—by means of data assimilation for a wider network of ionosondes, by better definition of the background ionosphere, and by an interpolation procedure aimed to reduce undesired local and border effects.

PECASUS will continue also the work for improved TEC forecasts. The eSWua TEC services (https://doi.org/10.13127/eswua/tec, accessed on 27 August 2021) include a novel empirical model for global TEC forecasts with a 24 h lead time [65]. The approach applies a nonlinear autoregressive neural network and is able to reproduce ionospheric conditions especially during disturbed periods. This product by [66] is implemented as a real-time service in the Ionosphere Prediction Service (IPS) prototype (http://ips.gsc-europa.eu, accessed on 27 August 2021). Machine learning and data mining covering observations from several solar cycles may be also more generally an attractive way to make progress in space weather forecasts. Long-term archives enable also probabilistic forecasts, whose interpretation is familiar for airliners from meteorological services.

As a step towards physics based modeling, the European Heliospheric Forecasting Information Asset (EUHFORIA, ref. [67]) was implemented at the PECASUS Advisory Consolidation and Production Hub with the goal to support space weather forecasting. EUHFORIA is a 3D Magneto–Hydro–Dynamic solar and heliospheric model that simulates solar wind (including HSS) and the evolution of CMEs from 0.1 AU to 2 AU. The code was specifically designed to estimate arrival times of CMEs and HSSs, and their impact in space weather operational activity. For the estimation of the potential impact of a CME, EUHFORIA needs input data of CME initiation in order to launch a simulation. Such input is available from the CACTus tool [61] which is already in routine use at the Advisory Consolidation and Production Hub. After its validation in an operational environment, EUHFORIA can be used by PECASUS as a supporting asset when considering inputs for the forecast entries of ICAO advisories.

Finally, some effort was directed toward measuring the availability, timeliness and quality of the data provision used in PECASUS and the issued advisories. However, this is an ongoing process that needs to be developed in the future in order for PECASUS to be able to evaluate and quantify its performance and the validity of the issued advisories. For example, validation campaigns will be organized, and statistical studies involving the issued advisories (which will become more numerous as the solar activity increases) will be carried out. While the PECASUS system is still in its infancy, it has already demonstrated that a self-organized community of research institutes is capable of providing a 24/7 operational service for a world organization as ICAO. We expect that—as more space weather–affected industries formulate their requirements—similar concepts and ways of working can be expanded to an increasing number of services.

Author Contributions

Conceptualization, A.-M.H., D.B. and M.G.; methodology, software, investigation, and formal analysis, J.A., P.B., J.B., C.C., E.D.D., J.d.P., M.D., E.D., K.H., H.H., E.H., M.K., T.L., M.L., Y.M., L.P., E.P., L.R., D.S., L.S., I.S., L.T., M.T., K.v.D., P.V., T.V., V.W. and A.Z.; validation, E.D., K.v.D., B.v.d.O. and K.Ö.; resources, funding acquisition and supervision, A.-M.H., D.B., M.G., J.B., I.S., V.R.; data curation, J.A., J.d.P., K.H., E.H., T.L., Y.M., L.T. and P.V.; writing—original draft preparation, K.K.; writing—review and editing, D.B., C.C., E.D.D., J.d.P., E.D., K.H., H.H., E.H., M.K., M.L., Y.M., L.P., L.R., V.R., D.S., L.S., I.S., B.V., P.V., T.V. and A.Z.; visualization, K.K., L.P., D.S. and A.Z.; project administration, A.-M.H., J.A., M.G., K.K. and K.Ö. All authors have read and agreed to the published version of the manuscript.

Funding

The ICAO member countries of the PECASUS consortium have agreed to maintain the above described space weather services at their own cost during the period from November 2019 to November 2022.

Data Availability Statement

PECASUS services during the project pilot phase (2019–2022) have used the following publicly available datasets:

• https://www.swpc.noaa.gov/, accessed on 27 August 2021): NRT GOES, Kp, solar wind

• https://www.gfz-potsdam.de/en/kp-index/, accessed on 27 August 2021): The official Kp index

• https://cosmicrays.oulu.fi/, accessed on 27 August 2021): Oulu NM data

• https://www.sws.bom.gov.au/HF_Systems, accessed on 27 August 2021): Ionospheric T-index

• https://www.igs.org/, accessed on 27 August 2021): NRT GNSS data

• https://www.epncb.oma.be/, accessed on 27 August 2021): NRT GNSS data

• https://impc.dlr.de/, accessed on 27 August 2021): TEC and scintillation indices

• https://doi.org/10.13127/eswua/tec, accessed on 27 August 2021): TEC

• https://doi.org/10.13127/eswua/gnss, accessed on 27 August 2021): Scintillation indices

• https://doi.org/10.13127/eswua/hf, accessed on 27 August 2021): Ionosonde data

• https://www.ngdc.noaa.gov/stp/spaceweather.html, accessed on 27 August 2021): Ionosonde data

Acknowledgments

The above listed instances are acknowledged for their data services. ICAO Meteorology Panel with its Meteorological Information and Service Development and the other ICAO space weather centers—NOAA/SWPC, Australia–Canada–France–Japan consortium, and China–Russia consortium—are acknowledged for fruitful collaboration in planning and implementation of the ICAO services. EUHFORIA was created as a joint effort between KU Leuven and the University of Helsinki and is being developed further by the the project EUHFORIA 2.0, a European Union’s Horizon 2020 research and innovation programs under grant agreement No 870405. K.K. thanks Atte Kapanen and Theresa Hoppe for their valuable assistance in preparing the figures for this paper and Petri Koskimaa for the statics on Kp-index.

Conflicts of Interest

The authors declare no conflict of interest.

Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Figures and Tables

Figure 1. Fictitious example of a space weather advisory. The latitudinal extent of the impacted area is specified with three-letter acronyms, which for the Northern hemispheric latitude bands are HNH, MNH, and EQN for 60–90°, 30–60°, and 0–30°, respectively. For the Southern hemispheric latitude bands, the corresponding acronyms are HSH, MSH, and EQS.

Figure 2. Locations of the GNSS scintillation receiver stations that are used in the PECASUS services. Stations maintained by the Ionosphere Monitoring and Prediction Center (https://impc.dlr.de/, accessed on 9 September 2021) and by the electronic space weather upper atmosphere (https://doi.org/10.13127/eswua/gnss, accessed on 9 September 2021) are marked with the green circles and blue squares, respectively.

Figure 3. Network of ionosonde stations that are used to support the PECASUS δfoF2 service during the project pilot phase. Stations maintained by the PECASUS consortium are marked with triangles. Stations, whose NRT data are fetched from the U.S. NOAA National Centers for Environmental Information, are marked with blue and gray dots for 15 min and 60 min time resolutions, respectively. Stations whose observations are used in the analysis of Section 5 are marked with orange triangles and cyan dots.

Figure 4. Examples of outputs from the PECASUS δfoF2 algorithm (28 September 2020, 20:00 UTC): (top) median δfoF2 values from NRT ionosonde data in the ICAO standard grid and (bottom) recommendations for post-storm depression advisories with the information on their intensity level (MOD/SEV).

Figure 5. Kp-index values (1st panel), and observed foF2 (black line) together with background foF2 (green line) of ionosonde stations at Juliusruh (2nd panel), Moscow (3rd panel), Pruhonice (4th panel), Rome (5th panel), and Nicosia (6th panel) during 28–30 September 2020. For Rome, Nicosia, Juliusruh, Pruhonice and Moscow locations, see the two orange triangles and three cyan dots in Figure 3.

Figure 6. MUF variations during 28–30 September by EUROMAP and IRI (STORM) models (red and blue curves, respectively) compared with the background model [51] (green curve) and with the ionosonde data (black curve) from Pruhonice (top panel) and Rome (bottom panel).

Figure 7. Ratio of MUF3000F2 to the background model at 20:00 UT on 28 September 2020. Regions with values below 0.7 in the map are those for potential ICAO advisory issuances.

Figure 8. Recommendations by the PECASUS δfoF2 algorithm for post-storm depression advisories for the evening of 28 September 2020.

Table 1

Parameters and thresholds that are used to issue space weather advisories according to ICAO regulations.

Impact Area	Parameter (Unit)	Moderate	Severe
GNSS	Amplitude scintillation S4 (dimensionless)	0.5	0.8
	Phase scintillation ${\sigma }_{\varphi }$ (radians)	0.4	0.7
	Vertical TEC (TEC Unit)	125	175
Radiation	Effective dose ($\mu$Sievert/hour)	30	80
HF	Auroral absorption (Kp)	8	9
	PCA (dB from 30 MHz riometer data)	2	5
	Solar X-ray (W/m${}^2$) (0.1–0.8 nm)	${10}^{-4}$	${10}^{-3}$
	MUF (%)	30	50

Table 2

Criteria for alerting the PECASUS operator on duty. Based on these criteria, the operator investigates and possibly sends an advisory. The moderated thresholds are given for each data type in Table 1, respectively.

Monitored Effects	Alert for the Operator
$S_4$ and ${\sigma }_{\varphi }$ Scintillation	Over the last 15 min, 50% of the scintillation
	observations recorded at one station have reached
	the moderated threshold.
Vertical TEC maps	A single data point anywhere on the map reached
	the moderated threshold.
Effective Dose FL ≤ 460	A single data point in the AVIDOS model at FL ≤ 460
	reached the moderated threshold.
Effective Dose FL > 460	A single data point in the AVIDOS model reached the
	sever threshold.
Auroral Absorption (Kp)	One of the monitored K and Kp indices reached the
	moderated threshold.
Polar Cap Absorption	A single data point anywhere on the maps reached
	the moderated threshold.
Shortwave Fadeout	GOES-16 X-ray flux reached the moderated threshold.
Post-storm depression	One or several areas in $\delta foF2$ maps remain
	above the moderated threshold for 45 min and
	one of the K and $Kp$ index was greater or equal to
	6 over the past 5 days.

Table 3

A summary of PECASUS methods and data inputs used to create ICAO space weather advisories.

Impact Area	PECASUS Solution	Augmenting Models	Inputs
GNSS	$S_4$ and ${\sigma }_{\varphi }$ from 16 stations		50 Hz GNSS data
	TEC maps	NeQuick [19];	IGS and EUREF
	Neustrelitz TEC [24]
Radiation	Avidos	ANeMos	Neutron Monitors
		COMESEP [39]	GOES protons
		HESPERIA	GOES X-rays
		UMASEP-500 [37]; [38]
HF	EUROMAP [49]	IRI [50]	ionosondes
	GDMF2 [48]	MUF background [51]	F10.7
	foF2 warnings	NeQuick [19]	Kp
			T-index
	D-RAP [41]		GOES protons
			GOES X-rays
Overall	24/7		solar wind
Space Weather	Monitoring		density,
	CACTus [61]		velocity,
			IMF
			Kp
			NRT Kp (NOAA)
			Local K
			GOES X-rays
			PROBA2 X-rays