Sensing the ionosphere with the Spire radio occultation constellation

Radio occultation (RO) provides a cost-effective component of the overall sensor mix required to characterise the ionosphere over wide areas and in areas where it is not possible to deploy ground sensors. The paper describes the RO constellation that has been developed and deployed by Spire Global. This constellation and its associated ground infrastructure are now producing data that can be used to characterise the bulk ionosphere, lower ionosphere perturbations, and ionospheric scintillation.


Introduction
The ionosphere is a dynamic environment that varies over scale lengths ranging from metres to thousands of kilometres, and on scale times of fractions of seconds to decades. Furthermore, many radio systems operate through or via the ionosphere whose performance can be affected by changing ionospheric conditions. Ionospheric impacts are generally limited to systems operating at L band (1-2 GHz) and below but cover a wide range of applications, including high frequency (HF, communications and radar, very high frequency (VHF, 30-300 MHz) satellite communications, and ultra-high frequency (UHF, 300-3000 MHz) space track radar and space-based radar (Angling et al., 2012).
Although many ground-based techniques have been developed to monitor the ionosphere, space-based methods are being increasingly used since they can measure the ionosphere over regions where ground-based sensors cannot easily be located, such as the oceans. This paper is concerned with a particular space-based method known as Global Navigation Satellite System (GNSS) radio occultation (RO) (Jakowski et al., 2009). This is a form of atmospheric limb sounding that uses transmissions from GNSS satellites in medium Earth orbit (MEO) and corresponding receivers on low Earth Orbit (LEO) satellites. GNSS RO was first demonstrated for remote sensing the Earth's ionosphere (and the neutral atmosphere) by the GPS/MET instrument (e.g., Hajj et al., 1994;Kursinski et al., 1996). Further individual missions such as CHAMP (Jakowski et al., 2002) and SAC-C (Hajj et al., 2004) provided opportunities to improve data processing methods. However, the ability to provide global ionospheric monitoring has come about through the use of constellations of RO satellites such as the Constellation Observing System for Meteorology, Ionosphere and Climate (COSMIC) (Hajj et al., 2000) and COSMIC2 (Hsu et al., 2018). Furthermore, commercial operators are now flying constellations of RO satellites. One such operator is Spire Global (https://spire.com), and this paper describes how Spire collects, processes, and uses ionospheric RO data.
The following parts of this section provide brief descriptions of how GNSS can measure the ionosphere (Sect. 1.1) and the RO method (Sect. 1.2). Then Section 2 describes the Spire measurement system (constellation, satellite platform, GNSS receiver, etc.), and Section 3 describes three example uses of the data global ionospheric data assimilation, high-resolution E region perturbation detection, and ionospheric scintillation detection. Finally, Section 4 describes future developments.

GNSS measurements of the ionosphere
Dual-frequency GNSS signals can be used to determine the total electron content (TEC) between a GNSS transmitter (Tx) and a GNSS receiver (Rx) (Jin et al., 2014). The TEC is the line integral of electron density on a path (l) between the Tx and the Rx. It is related to the frequency (f in Hz), excess phase (S in metres), phase refractive index (n) and electron density (N e in electrons/m 3 ) by Gorbunov & Kornblueh (2001): TEC is often expressed in TEC units (TECu), where 1 TECu equals 10 16 electrons/m 2 . A dual-frequency (f 1 , f 2 ) GNSS receiver can record pseudorange (P 1 , P 2 , the distance measured by the GNSS code between Tx and Rx, assuming propagation at the velocity of light in a vacuum) and phase (L 1 , L 2 ) for both frequencies.
The path separation between the two signals is usually small, and the assumption is often made that they both travel along the same straight line. The TEC can be estimated by forming the geometry-free combination (Teunissen & Montenbruck, 2017). This is advantageous as the combination cancels errors in the satellite position and clock.
The TEC estimated from the pseudorange is noisy due to multipath and to limitations of the receiver front end bandwidth. It is, however, an absolute measurement except for a bias: where DCB is the differential code bias. In contrast, the TEC estimate from the phase measurements has less noise but includes an unknown phase ambiguity (b): Since L1 and L2 signals can suffer from cycle slips, these need to be detected and corrected (Blewitt, 1990;Savastano et al., 2017). One approach to obtain calibrated TEC is to level the phase TEC to the pseudorange TEC, either by simply modifying the mean value or by weighting each point in the TEC arc according to the impact of multipath as measured by the divergence between code and phase (Pedatella, 2011). The final step is to estimate the Tx and Rx DCBs. For the GNSS transmitters and a network of ground receivers, these can be solved for as part of a thin shell ionosphere mapping process, and values are published by the International GNSS Service (IGS; Beutler et al., 1999) and other institutions such as the Chinese Academy of Science (CAS) (Wang et al., 2020). For space-based receivers, other methods are required for DCB estimation (e.g. Yue et al., 2011;Zhong et al., 2016).

Radio occultation
Radio occultation (RO) methods are a form of atmospheric limb sounding that was originally developed at Stanford University and the Jet Propulsion Laboratory (JPL) to sense planetary atmospheres (Fjeldbo et al., 1971). GNSS-RO can be applied to the Earth's atmosphere by monitoring transmissions from GNSS with a receiver on an LEO satellite. As the LEO satellite moves in its orbit, the GNSS satellites rise above or set below the horizon. Thus, the GNSS signals traverse the atmosphere at a range of heights. Currently, there are over a hundred active GNSS satellites transmitting radio-navigation signals, thereby offering a great opportunity to collect RO atmospheric soundings from LEO.
Comprehensive descriptions of GNSS-RO techniques can be found in Hardy et al. (1994), Kursinski et al. (1997), and Hajj et al., (2002). In summary, when electromagnetic radiation passes through the atmosphere, it is refracted. The magnitude of the refraction (bending angle a in Fig. 1) depends on the gradient of the atmosphere refractivity normal to the path, which depends on the composition of the atmosphere (see Eq. (6)). The bending angle cannot be measured directly; rather, the bending can be estimated using the Doppler shift of the signal, given precise knowledge of the satellites' positions and velocities.
An Abel Transform is often then used to invert the GNSS-RO measurements. Given an assumption of spherical symmetry, the bending angle (a) of the ray between the GNSS satellite and a LEO is (e.g., Hajj et al., 2002): where a is the impact parameter (= nrsin(/) = const); n is the refractive index; x is the refractional radius (= rn(r)). Then the refractive index can be expressed as: Formally, this requires the extension of the bending angle profile to infinite values of impact parameter. In practice, the bending angles above the LEO height are extrapolated from lower heights or derived from a model. The atmospheric refractivity (N) comprises terms dependent on the neutral atmosphere and the ionospheric electron density (n e ): where a 1 = 77.6 (K/mbar), a 2 = 3.73 Â 10 5 (K 2 /mbar), P (mbar) is the total pressure, P w (mbar) is the water vapour partial pressure, T (K) is the temperature, n e (m À3 ) is the electron density, f (Hz) is the operating frequency, k = 40.3 Â 10 6 (m 3 s À2 ) and the final terms represent higher-order terms related to propagation in a magnetised medium which are often neglected (Hoque & Jakowski, 2008). At ionospheric altitudes, the terms in pressure and temperature can also be dropped so that the refractivity is only related to the electron density and the frequency. Consequently, the refractive index profile derived from the Abel Transform can be converted to electron density. An alternative approach is to estimate the TEC using the method described in Section 1.1. Then, by assuming spherical symmetry, it is possible to relate the TEC to the electron density profile (e.g., Schreiner et al., 1999): where r 0 is the straight-line impact parameter (i.e., it is the shortest distance from the assumed straight-line ray to the centre of the Earth). This expression can be inverted thus: Both these approaches provide vertical electron density profiles with the high vertical resolution (on the order of the width of the first Fresnel zone) but with poor horizontal resolution since the geometry of the measurement results in the measured bending angles or TECs containing information from a large region of the ionosphere. Furthermore, the presence of horizontal structures will limit the performance of methods that assume spherical symmetry. Many methods have been developed to attempt to overcome this (e.g. Hernández-Pajares et al., 2000;Wu et al., 2009;Mannucci et al., 2011;Pedatella et al., 2015). A further alternative is to dispense with the Abel Transform altogether and assimilate the slant TEC (or potentially excess phase) directly into an ionospheric model (e.g. Angling, 2008). This will be discussed in Section 3.1.
2 The Spire ionospheric measurement system 2.1 Overview Figure 2 illustrates Spire's GNSS-RO observing system architecture. Spire owns the entire stack of the observing system, including the design of the GNSS payload hardware and software, the ground station network, mission operations, and ground-side processing (Cappaert, 2018;Bryce & Cappaert, 2019). This gives Spire the ability to modify one or more of the components of the observing system to optimise for ionospheric measurement quality, quantity, and latency.

The LEMUR satellite platform
Spire Low Earth Multi-Use Receiver (LEMUR) secondgeneration satellites (LEMUR2) use the 3U CubeSat form factor (Johnstone, 2020). An overview of the platform capability is given in Table 1. From flight model 91 (FM91) onwards, LEMUR2 v3.4 has a third stowable solar panel to increase the power generating surface (Fig. 3). Due to the increase in the available power that this bus provides, the GNSS-RO payload can be operated at high duty cycles; 98% duty cycle has been tested, though in practice lower duty cycles are usually used.
A new major revision of the LEMUR2 bus (v 4.0) has been designed to incorporate internal architectural changes and additional features such as an X-band downlink and S-band uplink capability to improve data quantity and latency. The LEMUR2 v4.X will become the standard Spire satellite platform in the future.

The STRATOS GNSS payload
The Spire GNSS-RO payload comprises the STRATOS GNSS receiver and antennas. The STRATOS payload has three antennas (Fig. 4). The precise orbit determination (POD) antenna is a single patch, right hand circularly polarised (RHCP)   antenna aligned in normal operation to point to the zenith. It can collect signals from all GNSS constellations, but only GPS signals are currently processed to form precise orbit and ionospheric observables. The RO antennas mounted on the side of the LEMUR2 each consist of three RHCP patches to create a higher gain antenna for observing occulting and low-elevation GNSS signals. These antennas are mounted on opposite faces of the LEMUR2 and nominally are oriented in the ram (forward RO, FRO) and anti-ram (backward RO, BRO) velocity directions to observe rising and setting occultation events. These antennas have a radiation pattern covering local satellite elevations from approximately +10°to À40°and with an azimuth of ±40°around the velocity vector direction. These antennas can be used for observing ionospheric radio occultations from +10°elevation above the local horizontal down to +60 km above the Earth's surface (and neutral atmosphere RO below that, Fig. 5).
The GNSS signals from all three antennas are separately collected, down-converted to complex baseband, and then digitised (Fig. 6). The digital signal processing (DSP) backend picks up those digitised streams, tracks the GNSS signals on them, and creates raw observation data files, in which it records the observed two frequency propagation delays. The raw data is split into the following streams: 1 Hz sampled, closed-loop tracked, phase and pseudorange observations. These are used to obtain precise orbits and to retrieve TEC. Currently, these are of GPS signals only; 50 Hz sampled, open-loop tracked (Ao et al., 2009), phase observations of occulting GNSS signals. 50 Hz sampled, closed-loop tracked, phase observations from a reference satellite at high elevation. This is used in post-processing for removing receiver clock errors.
A summary of the current STRATOS measurement setup is given in Table 2.

The Spire constellation and coverage
Spire's constellation consists of 3U LEMUR2 CubeSats. As of January 2021, more than 110 LEMUR satellites have been launched, and more than 40 can make GNSS-RO measurements. The satellites are spread over various LEO orbits, between 400 and 600 km in altitude and between equatorial inclinations and sun synchronous orbits (SSO) (Fig. 7).   There are currently over 100 active GNSS satellites transmitting radio-navigation signals. Over a day, occultations made using these transmitters are evenly geographically distributed (Fig. 8). However, due to the current configuration of the Spire constellation, the occultations are not yet distributed evenly in local time.

Total electron content data
The dual-frequency phase and pseudorange observations can be combined to estimate the TEC. These TEC measurements contain transmitter and receiver differential code biases that must be removed to obtain absolute TEC measurements. The transmitter (GNSS) DCBs are obtained daily from an external source (currently CAS) and subtracted from each TEC arc. For measurements made through the POD antenna, the minimum value of each TEC arc after the application of the GNSS DCB is stored. Then the lower quartile of the minimum TEC values, occurring at elevation angles greater than 40°and on the nightside of the orbit (1800 to 0600 local time), over the past seven days is used as an estimate for the receiver DCB (Zhong et al., 2016) and removed from the measurement. The TEC measurements made via the BRO and FRO antennas are not DCB corrected since it is impossible to apply the same type of DCB correction when high elevation (i.e., low TEC) arcs are not available.
The daily number of calibrated TEC measurements collected via the POD antenna has increased greatly from January 2018 to date (Fig. 9). TEC production is currently approximately 9 Â 10 6 observations/day (at 1 Hz sampling) across the constellation; this includes both high elevation measurements of the topside ionosphere and plasmasphere, as well as some low elevation radio occultation measurements of the bottom side ionosphere (Fig. 10). Approximately 4000 relative TEC RO profiles are also collected via the BRO/FRO antenna per day. In the future, these will be leveled with the POD data to provide a combined product.
Individual TEC arcs show the expected increase in TEC as the ray path transits the F-Region peak. Figure 11 shows an example of a calibrated TEC arc collected between FM101 and the GPS G30 satellite (labelled antPOD in the figure). It also shows a shorter, uncalibrated arc from the BRO antenna. This has been roughly aligned with the POD arc. The BRO antenna provides much better gain in an RO geometry, and this is reflected in the better signal-to-noise ratio (SNR) (Fig. 12).

Data latency
Low latency data is important for effective specification of ionospheric conditions, given the short timescales it can vary. Currently, the median time between an RO observation being made and the data is available within the Spire data centre is approximately 90 min across the constellation (e.g., Fig. 13). Note that different individual satellites have different latencies largely due to the relationship between their orbit and the location of the ground stations. With its large constellation of satellites and network of ground stations, Spire performs an ongoing optimization process to solve for the most efficient utilization of ground station resources while maximizing data volumes and minimizing data latencies. Optimization includes: Last-in-first-out (LIFO) data ordering during downlinking that ensures that, within the available contact times, the latest data are always downlinked before older data; Contact optimization that favours scheduling contacts with satellites that have larger quantities of lower-latency data in preference to those with older data.
In the future, latency will be further reduced by a number of different means, including intersatellite communication links that will allow more efficient data flow from the satellites to the ground. M.J. Angling et al.: J. Space Weather Space Clim. 2021, 11, 56 3 Example uses of Spire data Ionospheric data collected by the Spire constellation can be used in multiple ways. In this section we describe three typical uses: Global data assimilation of total electron content. Detection of low altitude ionospheric perturbations. Detection of ionospheric scintillation.

Global data assimilation of total electron content
Data assimilation (DA) is a mathematical method that combines observations with a numerical model taking into account the error statistics of both information sources (Twomey, 1977;Rodgers, 2000). For geophysical applications, observations are usually more accurate than the numerical model forecasts. However, observations are often sparse in space over short time windows. On the other hand, numerical models provide a solution at all model grid points but often suffer due to unknown model bias, lack of spatial resolution, or limited physics approximations. The combination of the two information sources can produce an optimal estimate of the state of the ionosphere provided suitable background model, and observation error statistics are used.
DA systems can produce 3D images of the ionosphere using data provided by a range of instruments such as ground-based GNSS receivers (e.g., operated by the International GNSS Service, IGS; Beutler et al., 1999), RO receivers (Komjathy et al., 2010;Limberger, 2015), and ground-based ionosondes (Angling & Jackson-Booth, 2011;Galkin et al., 2012). One such model is STEAM (the Spire TEC Environment Assimilation Model). In STEAM, a Local Ensemble Transform Kalman Filter (LETKF, Elvidge and Angling, 2019;Hunt et al., 2007) is used to assimilate ground and space-based observations. The LETKF is a resilient data assimilation method, providing useful corrections to a background state even in the presence of forecast model and/or observation operator non-linearities. It uses an ensemble of ionospheric models to represent the background  M.J. Angling et al.: J. Space Weather Space Clim. 2021, 11, 56 error covariance matrix, and it is necessary that the constructed ensemble is consistent with the model dynamics. In STEAM, this is achieved by generating the ensemble using a version of NeQuick (Nava et al., 2008) modified to be run with sets of perturbed model parameters that have the desired covariances in latitude and longitude. Since the STEAM ensemble is based on an empirical model there is no need to forecast each ensemble member to the next time step; instead, only the mean of the ensemble electron density and statistics of the perturbed F2 layer peak parameters forecasted from the current to the next timestep. These are then reflated to create a new ensemble of electron densities. STEAM is currently operated on a regular 5°Â5°regular latitude-longitude grid and with 50 height levels. Figure 14 shows an example vTEC map derived from a 3D STEAM electron density analysis. It also shows the locations of the input data and the increments made to the NeQuick background model. Note that at the current low solar flux conditions, NeQuick has a significant positive bias in the daytime ionosphere.   M.J. Angling et al.: J. Space Weather Space Clim. 2021, 11, 56

E Region perturbations
Perturbations (i.e., sporadic E and travelling ionospheric disturbances, TIDs) in the lower ionosphere can be measured using a "bottom-up" TEC estimation approach (Wu, 2018). This method uses high-rate (e.g., 50 Hz sampling) single frequency excess phase observations from a radio occultation receiver to derive profiles of relative slant TEC (sTEC) in the E region. The 50 Hz sampling of the relative sTEC profile corresponds to a vertical resolution of around 100 m, and this, therefore, allows the detection of fine vertical structures in the E region ionosphere. Figure 15 shows two observations made almost simultaneously (within 30 s) by two different Spire satellites (FM075 and FM080). Each receiver was making observations of transmitters from two different GNSS constellations (GPS and Glonass). Therefore, it is likely that the enhanced TEC observed at~106 km is due to a sporadic E layer in the ionosphere rather than an artefact of the collection or processing system. This sTEC measurement method and a means of automatically categorising the perturbations will be fully described in a future paper.

Ionospheric scintillation
Scintillation is characterised by rapid variations of a signal's amplitude, phase, and direction of arrival, and modification of its coherence time and bandwidth (Cannon et al., 2006) due to propagation through small-scale electron density irregularities in the ionosphere (Fejer & Kelly, 1980). The ionospheric irregularities vary in different regions around the globe regarding shape, size, plasma drift velocity, and orientation with respect to the geomagnetic field. They can be small-scale structures (from centimetres to meter scale) or large-scale ambient structures embedded in the ionosphere (Aarons & Guidice, 1966).
The occurrence of scintillation is frequency, location, and solar cycle-dependent (Aarons, 1982). Historically scintillation measurements have been made using stellar sources (Sofieva et al., 2013, and references therein), multifrequency radio beacon transmitters, such as the Defense Nuclear Agency's (DNA) Wideband satellite (Fremouw et al., 1978) and the USAF's Hilat and Polar Bear satellites (Basu et al., 1986). Following the advent of GNSS, these satellites have been used for scintillation measurements and monitoring (Kintner et al., 2007). As well as ground-based GNSS receivers, LEO receivers can also be used to measure scintillation; e.g., the CORRIS receiver on C/NOFS has been used to study equatorial scintillation (Huang et al., 2012a(Huang et al., , 2012b, and the COSMIC constellation  M.J. Angling et al.: J. Space Weather Space Clim. 2021, 11, 56 has been used to determine scintillation climatology (Dymond, 2012).
The Spire RO constellation can also be used to observe ionospheric scintillation. Figure 16 shows an example of a scintillation event detected on 8 March 2020 on the L1 frequency. The signal-to-noise time series is divided into two segments corresponding to scintillation (segment shaded in light blue) and non-scintillation (segment shaded in light orange) conditions. Even in the scintillation segment, the fade depths are relatively small (i.e., less than 3 dB), which is expected at this point in the solar cycle.
The power spectrum density (PSD) of the signal intensity and phase can be used to provide information about the irregularity characteristics causing ionospheric scintillation. Under weak scintillation conditions, the slope of the power-law spectrum (p) helps determine the wavelength dependence of the scintillation and provides useful information about the scale sizes of the ionospheric irregularities (Crane, 1976). Furthermore, the change in slope indicates the Fresnel frequency, and this provides information on the predominant scale size of ionospheric irregularities that cause amplitude scintillation. Figure 17 shows the intensity power spectrum for the two highlighted segments of Figure 16. The PSD for the section corresponding to scintillation has a slope of À2.77 and shows a change in the gradient at a Fresnel frequency of 45 mHz

Future developments
The Spire data collection system will continue to be developed across areas such as constellation design, receiver design, and new types of measurements. In the future, with a more planned constellation and greater duty cycle operation, greater measurement density can be achieved. Figure 18 shows the locations of simulated RO events for a sample hour over Europe given a constellation of 48 satellites distributed across eight Sun-synchronous orbital planes, running at 100% duty cycle, and using signals from GPS, Galileo, Glonass, and QZSS. Given this setup, coverage is broadly uniform across the Earth.   M.J. Angling et al.: J. Space Weather Space Clim. 2021, 11, 56 The STRATOS receiver is also undergoing enhancements to increase the number of GNSS channels that it can simultaneously track and provide self-calibration for DCBs. The first versions of the new receiver are now in orbit, and more will be deployed through 2021.
The STRATOS receiver can also be used to make GNSS reflectometry measurements (GNSS-R) (Jales et al., 2020;Nguyen et al., 2020). These measurements open new opportunities to make additional measurements of both scintillation (Molina & Camps, 2020) and TEC (Wang & Morton, 2021).

Conclusion
Radio occultation provides a cost-effective means of observing the ionosphere. It can form an important part of the overall sensor mix required to characterise the ionosphere over wide areas, especially in areas where it is impossible to deploy ground sensors. Spire has developed and deployed an RO constellation, which is now producing data that can be used to characterise the bulk ionosphere, lower ionosphere perturbations, and ionospheric scintillation. It is expected that Spire will continue to grow its capacity to make these measurements over the coming years.
Acknowledgements. The COSMIC-2 data used in Figure 8 was obtained via the COSMIC Data Analysis and Archive Center (CDAAC), https://doi.org/10.5065/t353-c093. The editor thanks two anonymous reviewers for their assistance in evaluating this paper.

Data availability
The data files used to produce Figures 11,12,15,16,and 17 are available on request from the corresponding author. Additional data is available for research purposes via the European Space Agency's (ESA) Earth Online portal, which offers European research institutes access to earth observation data from satellite missions operated directly by ESA and from Third Party Missions (TPM). Furthermore, researchers funded by the US Federal Government can access Spire data through the NASA Commercial Smallsat Data Acquisition Program (CSDAP), which evaluates and procures data from commercial sources that support NASA's Earth science research and application goals.