Ionospheric TEC can be calculated using the different delays of two transmitted frequencies from multi-frequency GNSS, including GPS. Recent years have witnessed an increased interest in investigating the ionospheric density variations during geomagnetic disturbances due to the fast-growing number of ground-based GPS receivers and the resulting regional or even global scale 2-D GPS TEC maps. We use the world-wide GNSS Vertical TEC (VTEC) data distributed in the publicly available CEDAR Madrigal database. The processed and binned VTEC data have a spatial resolution of 1° × 1° in geographical coordinates and a temporal resolution of 5 min. Detailed information about the data processing procedure has been described in Rideout and Coster (2006).

The Madrigal TEC data were ingested into the recently developed TEC matrix completion tool called VISTA (Video Imputation with SoftImpute, Temporal smoothing and Auxiliary data) (Sun et al., 2021). The VISTA algorithm is based on the SoftImpute software and contains two extensions of temporal smoothing and auxiliary data. In this case, the auxiliary data is the spherical harmonic fitting of the Madrigal TEC data. The VISTA model is able to provide completed TEC map and preserves the meso-scale TEC features in the final output, which is important for tracking high-density structures in our study.

The SuperDARN is an international collaboration operating high-frequency (HF) coherent radars in both the Northern and Southern Hemispheres (Greenwald et al., 1995; Chisham et al., 2007). They enable us to image the large-scale ionospheric convection and thus electric fields with high temporal resolution (∼1–2 min). The SuperDARN radars operate 24 h a day, 365 days a year and have been generating a large amount of data.

Based on measurement of magnetic field perturbations from the Iridium Communications constellation of more than 70 satellites at low Earth orbits (∼780 km), Active Magnetosphere and Planetary Electrodynamics Response Experiment (AMPERE) can provide a global view of Birkeland FACs in the high-latitude ionosphere. The orbits of the satellites form 12 longitudinal planes equally spaced by 2 h, and the satellites along each track are separated by 9 min. The magnetic field perturbation data are calibrated and then used in spherical harmonic inversion (Anderson et al., 2000). The FACs are derived by taking the curl of the fitted perturbations and are provided in the Altitude Adjustment Corrected GeoMagnetic (AACGM) coordinates on a 1 h local time by 1° latitude grid. More detailed description of the data product can be found in Anderson et al. (2014) and references therein.

DMSP satellites are equipped with comprehensive instruments that can measure the precipitating particles from the solar wind and magnetosphere as well as the ionospheric thermal plasma characteristics. DMSP SSIES (Special Sensor for Ions, Electrons, and Scintillation) plasma instrument package is mainly used in this study and can measure the ionospheric thermal ion and electron temperatures, ion density and composition, as well as cross-track and vertical drifts. The DMSP horizontal velocities are shown in the inertial reference frame. In the auroral and polar region, vast majority of the ions at the DMSP altitude is O⁺ and the total ion density or plasma density is essentially the same as the O⁺ density. The product of the plasma density and vertical drift is used to study the vertical flux in the polar region (e.g., Coley et al., 2006). DMSP SSJ (Special Sensor J) precipitating particle instrument are also used to analyzing the generation mechanism of upflows. They can measure precipitating electrons and protons from 30 eV to 30 keV.

The PFISR radar is part of the NSF-supported advanced modular incoherent scatter radar (AMISR) facility that is used to conduct studies of the upper atmosphere and to observe space weather events. It can simultaneously measure the altitude profiles of important ionospheric parameters, such as electron density, electron/ion temperatures, and line-of-sight ion velocity. PFISR was in the four-beam low-duty cycle mode during this period. Two types of pulses were transmitted, a long pulse and an alternating code pulse. The former is appropriate for F-region measurements (>175 km) and the latter for the E region (<175 km). Electron density measurements from these two pulses were combined to produce the altitude profiles in this study. The convection flows are calculated using the method described in Heinselman and Nicolls (2008).

The Global Ultraviolet Imager (GUVI) onboard NASA TIMED can yield global maps of thermospheric composition, including maps of the ratio between the O and N₂ column densities, which can monitor the neutral composition changes during geomagnetic activities (Strickland et al., 1995; Zhang et al., 2004). During this storm, the TIMED satellite was in the early afternoon to postmidnight orbital plane, which is ideal to observe the thermosphere composition change impact on the SED.

Figure 1 shows the solar wind and IMF conditions for this geomagnetic storm together with the AE and Sym-H indices at the bottom. During this storm, the Sym-H minimum reached ∼−110 nT, and thus this storm should be classified as an intense geomagnetic storm. During these 2 days, the solar wind speed gradually increased from ∼360 km/s to ∼600 km/s but no shock formed, and the IMF had large fluctuations during the velocity increasing period, which are classical signatures of corotating interaction region (CIR). The three vertical lines highlight the three important times. The first vertical line indicates the sudden negative IMF By increase (t1: ∼11:25 UT on March 06, 2016) and the beginning of the enhanced geomagnetic activities signaled by the gradually elevated AE. The IMF southward turning occurred at ∼15 UT and led to the rapid Sym-H decrease and supersubstorm with AE reaching 1,500 nT. The second vertical line indicates the time when the Sym-H reached minimum (t2: ∼21:20 UT on March 06, 2016). The last vertical line highlights the decline of the enhanced solar wind driving (t3: 06 UT on March 07), signaled by the solar wind dynamic pressure decrease and the reduction in the IMF fluctuation magnitude.

The TEC, SuperDARN convection and AMPERE FACs were combined and plotted in Figure 2 at a cadence of every 3 h starting at 11 UT on March 06. These combined data are shown in the magnetic local time (MLT) and magnetic latitude coordinates and reveal the mutual evolution of the ionospheric TEC and the electrodynamics resulted from the magnetosphere-ionosphere coupling processes. A supplementary movie in the same format as Figure 2 is also provided at a time cadence of every 2 min. Figure 2A shows the conditions just before the enhanced solar wind driving t1, and the TEC map shows the quiet time typical diurnal pattern with higher TEC on the day side and lower TEC on the night side. The convection and FACs were consistent with the quiet-time patterns confined near the pole. After t1 (Figure 2B), the convection and FACs started to expand to lower latitudes, and the SED began to develop in the early afternoon sector. Following the IMF southward turning, solar wind energy continuously entered the geospace system and the convection and FACs kept expanding to low latitudes (Figure 2C). Near the peak of the storm time t2 (Figure 2D), the Region-2 FAC system reached below 60° mlat in the afternoon sector and close to 55° mlat in the premidnight sector. During this period, the SED further developed and multiple polar cap patches and the SED plume entered the polar cap and drifted anti-sunward. Figures 2E–H show the decay of solar wind driving signaled by the shrinking of convection and weakening of the FACs. The impact of the negative ionospheric storm also became clear first in the morning sector and then gradually expanded to later local times. In particular, the SED plume narrowed, the TEC values decreased, and the plume eventually completely disappeared.

In order to capture the FACs and SED dynamics in a single plot, the time series of the TEC and FACs at 14 MLT from 40° mlat to 90° mlat are shown in Figure 3. The three vertical lines indicate the same important times as shown in Figure 1. As one can see, the FACs at 14 MLT strengthened and expanded equatorward after t1 and retreated to higher latitudes after t2 during the recovery phase of the storm. After that, their strength and location were relatively steady until the end of March 07. Several SED plume and polar cap patches can be seen extending to the pole, including two wide plumes drifting anti-sunward at ∼17 UT, 20 UT.

The SED and plume disappeared during the ionospheric negative storm phase due to the large-scale thermosphere composition change. The daily TIMED GUVI O/N₂ ratio of March 05–07 is shown in Figures 4A–C in geographic latitude and longitude coordinates. The quiet time values on March 05 were subtracted from the storm day values, and the delta O/N₂ ratios are plotted in Figures 4D,E. In each panel, the trajectories of the TIMED satellite are shown and the UTs when the satellite crossed the 60° geographic latitude in the Northern Hemisphere are labeled at the top horizontal axis. The quiet time map reveals the classical hemispheric asymmetry of the O/N₂ ratio with higher values in the winter hemisphere, i.e., the Northern Hemisphere. During storm time, in the low and mid-latitude regions, the O/N₂ ratio increased during the storm, while in the auroral latitudes, the thermosphere O/N₂ ratio change became negative after ∼17 UT on March 06, first in the Northern Hemisphere and then extended to low latitude in both hemispheres on March 07. This equatorward expansion of reduced O/N₂ ratio was seen for all sectors except near 120° W in the Southern Hemisphere due to the dipole tilt effect and low magnetic latitudes comparing with other regions with similar geographic latitudes. The expansion of the negative O/N₂ ratio change to low latitude is a classical signature for thermosphere composition change during a storm and the major reason for the negative ionospheric storm in the auroral and mid-latitude regions (Prölss, 2008; Fuller-Rowell, 2011; Wang et al., 2021). The impact of the thermosphere composition change on the ionospheric TEC can be seen clearly in Figures 2E–H and the supplementary movie.

PFISR was right underneath the dayside SED plume during the early storm recovery phase (Figure 2E) and observed the altitude profiles of the SED plume during its decay phase. As shown Figures 2E,F and in the supplementary movie, the width of the SED plume rapidly reduced after 00 UT on March 07. This width reduction mainly happened near the poleward boundary of the plume and was associated with the horizontal transport of lower density plasma from later MLTs due to the enhanced convection driven by the short-lived but intense southward IMF Bz after 23:30 UT on March 06. Figure 5 shows the PFISR observations from 17 UT on March 06 to 06 UT on March 07 and the TEC/FAC map at 01 UT on March 07. From top to bottom, Figures 5B-G show the direction, magnitude and vector of convection flows, vertical flow due to a combination of the convection and field-aligned flow, and altitude profiles of electron densities from beam 4 and beam 2. In Figure 5A, the magenta circle indicates the location of PFISR at 01 UT, right near the poleward edge of the SED plume. In Figure 5G, the density profiles of the SED plume, the low-density plasma from later MLTs at the poleward shoulder of the SED plume (∼01–03 UT) and the mid-latitude trough (after ∼03 UT) revealed large density gradients near the poleward edge of the SED plume. The PFISR beam looking at higher latitudes and with lower elevation angle observed the poleward edge earlier at ∼00:30 UT. The peak SED electron density exceeded ∼10¹² m⁻³ at ∼400 km at ∼00 UT on March 07.

During this decay phase of the SED plume, the accompanied plasma convection flows were weakly northwestward (Figures 5A,C), while the vertical flows due to the combined E × B convection and field-aligned flow component turned negative near the low latitude boundary of PFISR and weakly positive near 67° mlat (Figure 5D). These observations are consistent with the SED decay phase observations shown in Zou et al. (2014). Downward plasma flows can push plasma to lower altitudes with denser neutrals and thus speed up the loss process, more specifically, the charge exchange reactions during the first step of the dissociative recombination. The large westward convection flows observed after 03 UT were subauroral polarization streams (SAPS) accompanied with proton precipitations based on the DMSP F17 observations (not shown).

Four DMSP satellites (F15, F16, F17, F18) operated during this period and provided the precipitating particles and thermal plasma measurements in the polar ionosphere. In particular, DMSP F15 was in the 14-02 MLT orbital plane, which is the most suitable for detecting the SED contribution to ion upflow fluxes among the available satellites. Figure 6 shows two periods of DMSP F15 measurements during the large negative IMF By period, i.e., 11:35-12:01 UT (A-E) and 13:17-13:41 UT (F-J). Two selected TEC/convection/FAC polar view plots at 11:46 UT and 13:26 UT are also shown at the top. These two periods were at the beginning of the enhanced solar wind driving, and the SED TEC in the dayside mid-latitude and subauroral regions just started to increase. In Figures 6B, G, the topside SED plasma, characterized by enhanced density and relatively low electron temperature, were measured by DMSP F15 at ∼11:47 UT and ∼13:27 UT. The accompanied anti-sunward horizontal convection flows in Figures 6D, I were larger than ∼2 km/s and sandwiched by the large negative By related upward and downward FACs near noon (Figures 6A, F), while the vertical flows were elevated to several hundred m/s and exceeded ∼1 km/s at 11:47 UT. The combination of the high density and large vertical flow produced the upward flux reaching ∼0.8–1 × 10¹⁴ m⁻²s⁻¹ during these periods. The second DMSP orbit shown in Figure 6 also observed relatively enhanced densities at ∼13:35 UT at 2.5 MLT and produced ∼3 × 10¹³ m⁻²s⁻¹ upward fluxes. Based on the movie, these enhanced densities were likely the remnant of the SED plasma transported to this location.

After the IMF Bz southward turning, clear Region-1 and Region-2 FACs and two-cell convection patterns developed. The TEC values within the SED continued to increase and then multiple patches and the SED plume were seen moving anti-sunward. Figure 7 shows the DMSP F15 measurements during the peak of the storm from 20:01 to 20:27 UT (A-E) and from 21:43 to 22:09 UT (F-J) on March 06. In the same format as Figure 6. In Figures 7A, F, the high-latitude convection and FACs had already expanded to beyond 60° mlat. The large SED density was observed near 65.7° mlat at ∼20:07 UT, and the plume was observed at ∼ 71.8° mlat (Figure 7B). The horizontal convection flows associated with SED were sunward, peaking at ∼2.5 km/s, and the vertical flows exceeded ∼1.2 km/s. The peak upward flux again reached ∼1 × 10¹⁴ m⁻²s⁻¹ near the poleward shoulder of the SED. Interestingly, the ion fluxes within the SED plume at ∼ 71.8° mlat were downward, similar to previous observations deep in the polar cap (Ren et al., 2019). One orbit later, DMSP F15 observed the fully grown SED plasma, and the topside density at 850 km reached ∼2 × 10¹¹ m⁻³, tripled the value at the beginning of the storm in Figure 6B. In addition, the largest upflow fluxes during this storm was observed, exceeding 3 × 10¹⁴ m⁻²s⁻¹ at ∼12 MLT.

Figure 8 shows the DMSP F15 measurements during the negative ionospheric storm phase from 02:56 to 03:14 UT (A-E) and from 04:34 to 04:56 UT (F-J) on March 07. In the same format as Figures 6, 7. Evident in the TEC polar view plot, the TEC values in the Northern Hemisphere gradually reduced first in the noon sector and later in the pre-midnight sector. The topside ionosphere density shown in Figures 8B,G reduced to about 30% of the peak topside density shown in Figure 7. Despite the large vertical flows, i.e., several hundred m/s, due to the reduced ionospheric density, the upflow fluxes dropped to ∼2-3 x 10¹³ m⁻²s⁻¹. When comparing the second DMSP orbit in Figure 8 with the first orbit in Figure 6, one can see that the convection flows (Figure 8I and Figure 6D) and the electron temperature (Figure 8H and Figure 6H) during these two time periods were nearly identical, while the upward fluxes differed significantly due to the status of the ionospheric storm and the change of the topside ionosphere condition.

Based on the plasma temperature associated with the ion upflow in the F region and the topside ionosphere, ion upflow events have been conventionally classified into two categories, type 1 and type 2 (Wahlund et al., 1992). The type 1 ion upflows are related to strong perpendicular electric fields, enhanced and anisotropic ion temperatures, and low electron densities below 300 km. The enhanced ion temperature results in pressure gradients that propel the ions to higher altitudes. The type 2 ion upflows are related to electron precipitation, electron temperature increase, and reduced or unaffected topside electron density. The largely increased electron temperature in the topside ionosphere leads to an increased ambipolar electric field and thus ions are pulled upward together with the expanding electrons. The type 2 ion upflows can sometimes be accompanied by enhanced ion temperature as well.

In order to analyze the generation mechanism of the intense ion upflow fluxes and the associated precipitating particles, selected DMSP F16 passes are shown in Figures 9–11. DMSP F16 was in the 15-03 MLT orbital plane. The format of these three figures is the same, and from top to bottom shows the TEC/FAC map (a), precipitating electron (b) and proton (c) energy fluxes in log scale, in situ density (d), plasma temperature €, horizontal and vertical velocities (f), and vertical fluxes (g). In Figure 9, DMSP F16 crossed the afternoon sector SED base before 19:37 UT and observed the increased in situ density of ∼10¹¹ m⁻³, with a small peak near the equatorward edge of the proton precipitation and the SAPS flow region (highlighted by a dashed line). The density peak was clearly associated with the peak upflow flux, ∼6 × 10¹³ m⁻²s⁻¹, and the SAPS. The westward flux transported by this SAPS channel was ∼2 × 10¹⁴ m⁻²s⁻¹, comparable to the previous reported values (Erickson et al., 2011).

Both electron and ion temperatures associated with this SAPS upflow were lower than the surrounding region. The fact that this peak upflow flux was collocated with SAPS suggests that this is likely a type 1 upflow. The type 1 ion upflows are usually associated with strong perpendicular electric field (e.g., Sellek et al., 1991; Heelis et al., 1993). Sellek et al. (1991) used an ionosphere and plasmasphere model to study the effect of frictional heating due to a 2 km/s westward drift. Their results showed that the O⁺ temperature profile increases from ∼1000 to ∼3200 K between 200 and 500 km, decreases rapidly from ∼3200 to ∼2000 K between 500 and 750 km, and no enhancement above 750 km. Similarly, Heelis et al. (1993) studied the effect of frictional heating of a 2 km/s horizontal drift on O⁺ temperature and upflow velocity with a focus on their transient dynamic evolution. Their simulation results showed a negative temperature gradient between 300 and 1,000 km about several minutes after reaching the peak velocity. Therefore, it is not unexpected that there was no clear ion temperature enhancement signature associated with this SAPS upflow at the DMSP altitude.

Figures 10, 11 show the DMSP passes crossing the SED plume during its peak and its decaying phase, respectively. In Figure 10, the peak upflow fluxes increased to ∼10¹⁴ m⁻²s⁻¹ at ∼23:05 UT, when the SED plume met with intense soft electron and proton precipitations. These soft particle precipitation signatures are classic cusp precipitations that directly coming from the equatorial magnetosheath. In Figure 11, the SED plume continued to contribute to the large ion upflow fluxes even during its decaying phase. DMSP F16 encountered the SED plume at ∼00:46–00:48 UT and observed the narrow plume contribution to the large upflow fluxes there. Like Figure 10, the precipitating soft electron and ion signatures indicate again this was near the dayside cusp region. Together with the signature of largely elevated electron temperature, the large ion upflowes shown in Figures 10, 11 should be classified as type 2 ion upflow seeded by the SED plume plasma. It is interesting to note that in Figures 10, 11, the ion vertical fluxes associated with the SED at lower latitudes were downward (before 23:02 UT in Figure 10 and at ∼00:46 UT in Figure 11), consistent with the PFISR observations described in the previous section and earlier studies (Zou et al., 2014; Ren et al., 2020).

In order to show the temporal evolution of integrated upward fluxes, we integrated the upward fluxes above 45° mlat in the Northern Hemisphere for all DMSP F15 and F16 satellites and created the time series of integrated upflow flux in Figure 12. DMSP 15 and DMSP 16 have similar orbits, i.e., covering afternoon and postmidnight sectors, and their observations are also comparable. The integrated fluxes gradually increased during the positive ionospheric storm, due to a combination of the elevated SED plasma density and the increased vertical flow speed. The hemispheric integrated flux reached its peak at ∼22-23 UT, shortly after the ring current attainded its peak and the Sym-H approached its minimum. After the ionospheric negative storm initiated due to the thermospheric composition change, the integrated flux decreased significantly, despite the fact that the vertical flow speed was comparable to that during the earlier positive storm condition. These observations clearly demonstrate the importance of the ionospheric storm phase and status in determining the total available ion upflow fluxes. DMSP F17 and F18 were close to each other in the dusk and dawn orbital plane, and the SED contributions to the upflow fluxes there are difficult to quantify. Because at these two local times, the SED plasmas usually have already traveled a long distance from the dayside convection throat to the nightside auroral zone and on the way back to the dayside via the return convection flows. High quality convection flow pattern would be required to trace the SED plasma parcel to ensure that the enhanced densities observed at the dusk and dawn auroral zone are indeed part of the SED remnant. Therefore, the DMSP F17 and F18 observations are not shown here and further studies regarding the contribution of SED plasmas to the ion upflow fluxes in the return flow region of the auroral zone should be performed.

A back-of-the-envelope estimation of the SED contribution to the integrated upflow fluxes in the dayside cusp region is performed and compared with the integrated upflow fluxes in the nightside auroral zone. In Ren et al. (2020) paper, the average ion upflow fluxes observed by PFISR in the nightside auroral zone is ∼ 2-3 × 10¹³ m⁻²s⁻¹. So, the SED-associated upflow fluxes in the cusp can be ∼5–10 times larger than the average upflow fluxes observed in the nightside auroral zone. If we assume the typical dayside cusp is located at 75° mlat and covers 2 h of MLT and 2° mlat, as suggested by the statistical picture shown in Newell, (2004), and the SED contributed upflow flux is ∼1-2 × 10¹⁴ m⁻²s⁻¹, the integrated upflow flux in the cusp with SED contribution could be ∼ 2-4 × 10²⁵ s⁻¹. Similarly, if assuming the nightside auroral zone is located at ∼65° mlat and covers 6 h of MLT and 5° mlat and the upflow flux is ∼2-3 × 10¹³ m⁻²s⁻¹, the total auroral zone upflow flux could be ∼4-6 × 10²⁵ s⁻¹. Based on the above estimations, the SED-associated upflow fluxes in the dayside cusp is comparable to the total upflow fluxes in the nightside auroral zone despite the much-limited spatial coverage of cusp.