Signatures of Equatorial Plasma Bubbles and Ionospheric Scintillations from Magnetometer and GNSS Observations in the Indian Longitudes during the Space Weather Events of Early September 2017

: Scintillation due to ionospheric plasma irregularities remains a challenging task for the space science community as it can severely threaten the dynamic systems relying on space-based navigation services. In the present paper, we probe the ionospheric current and plasma irregularity characteristics from a latitudinal arrangement of magnetometers and Global Navigation Satellite System (GNSS) stations from the equator to the far low latitude location over the Indian longitudes, during the severe space weather events of 6–10 September 2017 that are associated with the strongest and consecutive solar flares in the 24th solar cycle. The night-time influence of partial ring current signatures in ASYH and the daytime influence of the disturbances in the ionospheric E region electric currents (D iono ) are highlighted during the event. The total electron content (TEC) from the latitudinal GNSS observables indicate a perturbed equatorial ionization anomaly (EIA) condition on 7 September, due to a sequence of M-class solar flares and associated prompt penetration electric fields (PPEFs), whereas the suppressed EIA on 8 September with an inverted equatorial electrojet (EEJ) suggests the driving disturbance dynamo electric current (Ddyn) corresponding to disturbance dynamo electric fields (DDEFs) penetration in the E region and additional contributions from the plausible storm-time compositional changes (O/N2) in the F-region. The concurrent analysis of the D iono and EEJ strengths help in identifying the pre-reversal effect (PRE) condition to seed the development of equatorial plasma bubbles (EPBs) during the local evening sector days’ levels and the timestamps of ROTI peak magnitudes substantiate the severity, latitudinal time lag in the peak of irregularity, and poleward expansion of EPBs and associated scintillations. The key findings from this study strengthen the understanding of evolution and the drifting characteristics of plasma irregularities over the Indian low latitudes.


Introduction
The study of spatiotemporal propagation of ionospheric irregularities during space weather events remains an interesting topic among the ionospheric community. Rapid fluctuations in ionospheric irregularities can introduce severe ionospheric scintillations to disrupt the satellite-based communication links, degrading the performance or even can result in the loss of lock in the global navigation satellite system (GNSS) receiver measurements. The scintillations are known to be the random fluctuations of the signal amplitude and phase caused by the rapid small-to-medium scale plasma irregularities in the ionosphere, whose occurrence characteristics may change during the disturbed space weather condition. Although most of the aspects of the large-scale delay errors corresponding to the total electron content (TEC) along the signal ray path through the ionosphere have been resolved in the past, the scintillation effects remain a challenging task for the ionospheric modeling community as they can cause a serious threat to dynamic systems, depending on the space-based navigation services. There have been several studies on diverse space weather events in the past that strengthen the understanding of variability and modeling to a greater extent. However, some of the obscured features associated with the plasma dynamics, particularly over the equatorial and low latitude sector yet remain unsolved. The day-to-day dependence of the lives, technology infrastructure, navigation, and positioning are now concerning the space-based navigation services, which in turn depends on accurate modeling and the mitigation of effects due to plasma irregularities. The quality of the received signals at the GNSS receiver may be severely damaged by the scintillations, affecting the accuracy and reliability of the position, velocity, and time (PVT) solution. On the other hand, the TEC and its derivatives, as well as the scintillation indices extracted from the GNSS observations, serve as the signals of opportunity for understanding and mitigating the effects in the real-time applications.
In general, scintillation is believed to be strongly influenced by local time, season, solar activity, geomagnetic conditions, and wave propagation from the underlying atmosphere. The coronal mass ejections and solar wind interaction with the earth's magnetosphere can cause severe ionospheric irregularities by introducing scintillations in the signals through plasma depletions associated with the equatorial plasma bubbles (EPBs) [1]. Geomagnetic storms may occur during the southward polarity of the interplanetary magnetic field (Bz) that manipulates the regular equatorial electric fields to trigger pre-reversal enhancements (PREs) to seed the generation and development of plasma bubbles [2][3][4]. Additionally, the substorms formed at the polar latitudes can alter the equatorial electric field during the recovery of geomagnetic storms resulting in an increase or decrease in the scintillation activity [5][6][7][8]. Generally, the ionospheric electrodynamics during the geomagnetic storms are affected by two sources: (1) the short-lived penetrating electric fields from high to low latitudes corresponding to the southward turning of the interplanetary magnetic field Bz to drive eastward (westward) the polarized disturbances at the day and evening sectors (nightside) [9][10][11][12]; and (2) the disturbance dynamo electric fields (DDEFs) resulting from the changes in neutral winds that develop a few hours after the onset of the storm and usually last for several hours via thermospheric wind dynamo action, often dominating in the recovery phase [13,14]. It has been observed that during an extended main phase of the storm, the prompt penetration phase may coincide with the local dusktime PRE at the equatorial ionosphere to augment the post-sunset ionospheric irregularities, even if this occurs with the coexistence of DDEFs [15]. The DDEFs play an important role in the occurrence of daytime irregularities during the recovery phase of a storm [16] and references therein. It has been established that the scintillations are more frequent during the solar maximum years in the equatorial and auroral regions, whereas their occurrences are sparse in the mid-latitudes [17,18]. In the equatorial regions, the scintillations are more prominent in the post-sunset and midnight hours due to the occurrence of plasma bubbles [19][20][21]. However, there are also shreds of evidence for daytime scintillations, mostly associated with the sporadic-E (Es) structures, whose effects are relatively mild and less frequent, compared to the night-time scintillations [22,23]. The scintillations are quantified in terms of two basic measurable quantities, namely the amplitude scintillation index (S4 index) and phase scintillation index (σΦ). In GNSS operations, amplitude scintillations refer to the rapid fluctuations in the signal intensity (or carrier-to-noise ratio) measured by the receiver, whereas phase scintillation refers to the rapid fluctuation in the carrier-phase measurements, determined by the standard deviation of the detrended carrier phase over a period [24,25]. While severe scintillations can cause a loss of lock in the GNSS receivers, often making it impossible to calculate the position for a period, a less severe scintillation condition may degrade the positioning and navigation accuracy by adding uncertainty to the signal. In the case of the amplitude scintillation, the signal strength is degraded resulting in the possible loss of lock to attempt the reacquisition of the signal, whereas the phase scintillation introduces a cycle slip or even the loss of signal lock by holding considerable amount of time for the reacquisition of the signal [17,18]. Hence, the dynamic systems relying on the GNSS carrier phase tracking measurements are very much vulnerable to the scintillation conditions [26].
Apart from the amplitude and phase scintillations, the intensity of ionospheric scintillations or TEC fluctuations can be quantitatively described by another index, called the rate of TEC index (ROTI) [27,28]. The amplitude and phase scintillation indices are usually sampled at a high frequency with a specially configured cost-effective GNSS receiver exclusively for space weather monitoring, whose availability is sparse across any region. However, the relatively abundant ordinary non-scintillation geodetic receivers sampled at a lower frequency could provide equivalent ROTI estimates, which are equally useful for understanding the ionospheric irregularities in the absence of ionospheric monitoring receivers [29][30][31]. ROTI refers to the standard deviation of the rate of change of TEC (ROT), whose values per minute can be obtained from the slant TEC estimations following the dispersive nature of refractive ionospheres in the dual-frequency phase and pseudorange observables. Over the last two decades, the extensive analysis of scintillation indices and ROTI variations demonstrated a close relationship among the parameters [25,32,33]. Hence, it has been used by several ionospheric groups to analyze the TEC fluctuations from an adequate number of GNSS stations across the globe. Moreover, the next generation international GNSS service (IGS) ROTI maps product can serve as a valuable tool for monitoring global ionospheric irregularities and reanalyzing the impact of plasma irregularities on the GNSS positioning in the past [34].
There are several studies presenting scintillation occurrence characteristics with the help of ROTI in the equatorial, low-latitude, and sub-auroral regions, confirming that ROTI can be taken as a parametric index to identify the spatiotemporal variation of scintillations [33,[35][36][37][38]. The ROTI observations from 21 reference stations in China show its consistent relation with scintillation characteristics [39]. The longitudinal variation in the percentage occurrence of ROTI is similar to the L-band variation [40]. Similarly, the generated ROTI and S4 maps from 26 GPS stations show an excellent agreement between the parameters, even during disturbed geomagnetic conditions [36]. The seasonal correlation studies between the ROTI and S4 at a low latitude station (Hyderabad) demonstrate a maximum correlation of 0.95, even during the disturbed days [1]. From the occurrence probability studies, it is observed that in most cases the log-normal distribution of S4 fits well in the overall range of ROTI [38]. The time series analysis of ROTI demonstrates that a value greater than 0.5 can be considered as a threshold for identifying the presence of ionospheric irregularities [41]. The scintillation characteristics were presented in the past by using the extracted ROTI from the GPS aided geo augmented system (GAGAN) and other data over the Indian low latitude region [42]. Furthermore, the observations at a low latitude location in Visakhapatnam in India reported strong fluctuations in ROTI during the St. Patrick's day (17 March 2015) geomagnetic storm [43].
The space weather event of 6-10 September 2017 was a complex event associated with multiple X-and M-class solar flares and severe geomagnetic storms, when the GNSS signals witnessed pronounced post-sunset scintillations due to strong EPB irregularities. The occurrences of plasma density irregularities, their drifting characteristics, and the mechanisms responsible for this occurring over different longitude sectors during the storm of September 2017 have been widely investigated in the earlier literature [3,8,[44][45][46][47][48][49][50][51]. From the investigation of the topside ionospheric conditions during this particular storm period, with uplooking GNSS TEC from the observations retrieved from TerraSAR-X, GRACE, Swarm, and MetOp-A, and the in situ electron density (Ne) from the Swarm satellite, Jimoh et al. [50] reported night-time ROTI enhancements across a wide latitudinal range during the main phase of the storm. With sufficient proof from the observations of the GPS network, SuperDARN radar, HF radar observations, and the in situ density structures in the Swarm and DMSP signals in North America, Nishimura et al. [48] reported the evolution of different density irregularity structures across a wide range of wavelengths in the middle-latitudes, during the main and recovery phases in this particular event. Concurrent observations of the network of GNSS receivers, ionosondes, and very high frequency (VHF) radars on 8 September 2017, reveal intense poleward drifting EPB irregularities in the American and Asian sectors attributing to the significantly enhanced PRE vertical drift, because of the in-phase prompt penetration under-shielding the eastward electric field enhancement during the local post-sunset hours [47,52]. Similarly, multi-instrumental analysis, including observations from ionosondes, GNSS ROTI maps, and the Swarm satellite in the South East Asia longitude of the sectors in China, indicated the existence of a pre-reversal vertical drift and night-time spread F (8 and 9 September 2017) accompanying spread Es (9 September 2017) at low and middle latitudes during the recovery stage of the storm [49]. However, ground-and satellite-based observations over Brazilian equatorial-low latitudes showed a substantial increase in the scintillations on 7 September 2017; the scintillation intensity was evidently inhibited on 8 September 2017 [51]. Concerning the Indian sector, with the help of EEJ strengths from ground-based magnetometers and L-band TEC and scintillation observations from GNSS signals during 8 September 2017, Alfonsi et al. [3] presented the evidence of freshly generated and inherited EPBs and their northward propagating structures.
Despite a good number of global and regional investigations on the event, we realize that there is still scope to improve the understanding of the driving mechanisms of and possible effects on a regional aspect, in particular, for the Indian longitudes that reported strong ionospheric disturbances/plasma irregularities during diverse geomagnetic storm events. Hence, in this paper, we probe the plasma irregularity characteristics that signify the strong scintillation effects on transionospheric signals during the space weather events of 6-10 September 2017. This is conducted by complementing the role of prompt penetration and ionospheric disturbance dynamo electric fields through an analysis of the ionospheric currents from local magnetometer recordings and real-time prompt penetration electric field model (PPEFM) runs, along with TEC and ROTI variations from the GNSS observables at latitudinally aligned locations from the equator to a higher low-latitude in the Indian longitude sector. From our observations in the literature and our own analysis, we realize that the earlier report on the Indian longitude is confined only to the event that took place on 8 September and does not highlight the EPB characteristics around the same time on the nearest day (10 September), which evidenced the visible plasma density irregularities at the equator and nearby latitudes [3]. Moreover, to the best of our knowledge, there are hardly any reports of this event with conclusive evidence of the signatures of ionospheric currents from a network of magnetometer observations along the Indian longitude sector. The probable influence of the enhanced asymmetric ring current (ASYH) signatures in the global magnetic data at the local post-midnight sector that present an indication of the prevailing ionospheric disturbance dynamo and changes in the thermospheric neutral compositions over the Indian longitude sector, has not been underlined in the earlier report on EPBs over the same region. We emphasize the spatio-temporal characteristics of ROTI to substantiate the understanding of the scintillations due to plasma bubble occurrences and their drifting behavior during the storm. The analysis from the results reveals strong plasma irregularities at the equator and their poleward expansion with a temporal lag that is discussed in the subsequent sections of this study.

Materials and Methods
The present work is undertaken for the space weather events of 6-10 September 2017 that are associated with the strongest and consecutive solar flares of the 24th solar cycle. The datasets used in this study are considered from four ground-based magnetometer observatories and four ground-based dual-frequency GPS stations, latitudinally aligned from the magnetic equator to the higher low-latitude beyond the northern anomaly crest in the Indian longitude sector, as shown in Figure 1 and Table 1. The ionospheric piercing point (IPP) locations for the selected GPS satellites (PRN-14 and PRN-32) during the storm day (8 September 2017), whose discussions are presented in the latter part of the study, are indicated by green (+) and blue (+) curves, respectively. Also, the interplanetary electric and magnetic field components (Ey and Bz), geomagnetic index (Kp), the components of the symmetric-and asymmetric-H indices (SYMH and ASYH), and auroral electrojet index (AE) are used in this study to support the analysis. The details about the datasets as well as the methodology involved in the data processing and analysis are described as follows. .

Global Geomagnetic Indices and Interplanetary Parameters
In order to understand the progress of the geomagnetic storm and support the analysis in the local observation data, the three-hourly Kp index during the storm period was taken from the GFZ German Research Centre for Geosciences (ftp://ftp.gfz-potsdam.de, accessed on 9 August 2021). The selective geomagnetic quiet days (day of the month: 3, 9, 10, 19, 21, 25, and 26) were chosen from the list of international quiet days, based on the relative indication of the day characteristics in the month. The 1-min averaged Ey, Bz, AE, SYMH, and ASYH data were obtained from the Omniweb interface (https://omniweb.gsfc.nasa.gov/, accessed on 9 August 2021).

Local Magnetometer Data Processing
The magnetometer recordings from four observatories along the Indian longitude were used to analyze the ionospheric E region response to the geomagnetic storm. The 1min datasets from three magnetometer observatories (Hyderabad (HYB), Alibag (ABG), and Jaipur (JAI)) were obtained from the International Real-Time Magnetic Observatory Network (INTERMAGNET; http://intermagnet.org, accessed on 14 October 2021), whereas the 1-min data for Tirunelveli (TIR) was obtained on request from the Indian Institute of Geomagnetism (http://iigm.res.in/, accessed on 6 October 2021). The geographic and geomagnetic coordinates of the magnetometer observatories along with their respective magnetic dip can be referred to in Table 1 and Figure 1. The variations of the horizontal component (H) of the magnetic fields at these observatories were considered to evaluate the disturbances in the ionosphere produced by the electric fields and currents at the E region heights by following two different approaches. The proxy equatorial electrojet (EEJ) strength to the EEJ current in the ionospheric E-layer was calculated by subtracting the local midnight baseline values from the magnetometer H variation at equatorial TIR and the off-equatorial ABG observatory, and then finding the difference between them, as shown in Equation (1). This approach was previously proposed by Nair et al. [53], and thereafter used by many researchers as a proxy to understand the characteristics of transmitted fields and the state of the equatorial ionosphere, especially during the daytime [2,54,55].
where ΔHTIR and ΔHABG are the variations of the magnetic field at TIR and ABG, respectively.
Further, in order to analyze the ionospheric E region response to the space weather event at different latitudes, we used the magnetic field data from four observatories, namely Tirunelveli (TIR), Hyderabad (HYB), Alibag (ABG), and Jaipur (JAI) in India, latitudinally aligned from the equator to the northern equatorial anomaly crest region, as shown in Figure 1. It is known from Biot and Savart's law that the integral of ionospheric electric currents and magnetospheric electric currents can be realized from the perturbation in the ground magnetic field [56,57]. The magnetic field perturbation corresponding to the ionospheric electric currents is denoted as Diono that is essentially restricted to the E region. However, at low and middle latitudes, Diono is effective only during the daytime when there are conductivities due to solar radiation. As per Le Huy and Amory-Mazaudier [57], during the storm, the H of geomagnetic field can be written as a superposition of components (see Equation (2)); Here, H0 is the baseline magnetic field due to Earth's core and SR is the regular variation of Earth's magnetic field on a given day, and Sq is the average of SR for selective magnetic quiet days. DM and Diono represent the variations in the magnetic fields associated with the magnetospheric and ionospheric currents, respectively.
ΔH can be obtained by subtracting H0 from the H component of the magnetic field. DM is calculated from the SYMH index normalized with the dip latitude (Φ), by assuming that the greatest contribution to the daytime H comes from the symmetric ring current, as presented in Equation (4).
Sq is estimated by taking the average value of quiet ΔHi from the 1-min resolution values of the selective quietest days in September 2017, as presented in Equation (5). The selective quiet days out of 10 international quiet days are chosen by following the quiet days' definition of the GFZ German Research Centre for Geosciences (ftp://ftp.gfz-potsdam.de/pub/home/obs/kp-ap/quietdst/, accessed on 9 August 2021).
At low latitudes, Diono predominantly constitutes the combined effect of the ionospheric disturbance polar currents (DP2) and the ionospheric disturbed dynamo currents (Ddyn) [59]. In terms of the magnetic disturbance, DP2 is the disturbance associated with the penetration of the magnetospheric convection (PPEF; prompt penetration electric field) [9,10], and Ddyn is the magnetic disturbance associated with the disturbed ionospheric dynamo (DDEF; disturbance dynamo electric field) [13,57].

GPS TEC Data Processing
The GPS observation data are used in this study for investigating the ionospheric total electron content (TEC) and the rate of change of TEC index (ROTI) fluctuations during the space weather events of 6-10 September 2017, with the GPS observables at four latitudinally aligned GNSS stations ranging from an equatorial location to a higher lowlatitude in the Indian longitude sector. All these stations are under the international GNSS service (IGS) network, of which one is situated in the vicinity of the magnetic equator in Colombo, Sri Lanka (SGOC), and the other three stations, Bangalore (IISC), Hyderabad (HYDE), and Lucknow (LCK4), are located in India. The geographic and geomagnetic coordinates of the stations, along with the respective magnetic dips, are presented in Table  1. The geographical location of the stations along with the approximate position of the magnetic equator and northern equatorial anomaly crest (EIA crest) contour line is shown in Figure 1. The dual-frequency GPS observation data of these stations were collected from the Crustal Dynamics Data Information System (CDDIS) (https://cddis.nasa.gov/archive/gnss/, accessed on 6 August 2021). The Almanac data required for extracting the TEC and ROTI from the GPS observables were accessed from the United States Coast Guard Navigation Center (https://www.navcen.uscg.gov/, accessed on 30 July 2021). To analyze the TEC and ROTI deviations during the disturbed days, we considered the average values of the parameters for the selected geomagnetic quiet days of the month.
Initially, the slant total electron content (STEC) was extracted from the dual-frequency GPS phase and pseudorange observables following the frequency-dependent dispersive characteristics of the ionosphere, as shown in Equation (6). The STEC refers to the total number of electrons along the ray path from the satellite to the receiver, measured with a unit of TECU where 1TECU = 10 16 e/m 2 .
where f1 and f2 are the dual GPS frequency signals and P1 and P2 are the respective pseudorange observations. The vertical equivalent (VTEC) of STEC at the ionospheric piercing point (IPP) was then determined using the single-layer ionospheric model mapping function, as shown in Equation (7).
where M(E) is the mapping function whose mathematical expression is presented in Equation (8).
Here, h is the ionospheric piercing point (IPP), considered to be at an altitude of 350 kilometres. Re is the radius of the earth (Re = 6378 km) and E is the elevation angle, chosen as 30° to eliminate the possible interferences due to the multipath and tropospheric effects while determining the VTEC data.
The ROTI corresponding to each satellite signal received is calculated as the standard deviation of the rate of change of TEC (ROT), as presented in Equations (9) and (10). The ROT is obtained by computing the relative TEC changes from one epoch to another.
where i indicates the visible satellite and n denotes the time of epoch. Hence, t − t is the time interval between the subsequent epochs. The standard deviation of the ROT is then calculated, as shown in Equation (10).
where ‹ › indicates the averaging over every nonoverlapping 5 min interval.
To analyze the severity of irregularity on the storm day, we calculated the diurnal variation of the maximum absolute error (MAE) and root mean square error (RMSE) of ROTI on the storm day (September 8), with respect to the average quiet days' reference levels, as presented in Equations (11) and (12), respectively.
where the ROTIS denotes the timestamp data of the storm day and ROTIQ refers to the corresponding temporal averages of the selected quiet days of the month. N is the number of data points considered in the study.

Results and Discussion
The highlighting feature of the space weather event of 6-10 September 2017 is that it was associated with multiple X-and M-class solar flare eruptions, including the strongest solar flare (X9.3) of solar cycle-24. Figure 2 shows the variations of the Kp, ASYH, SYMH, Bz, Ey, and AE parameters during the period, arranged from the bottom to the top panel. The arrow-headed, vertical, red, dashed lines indicate the onset of SSCs (SSC1 at 23:43 UT on 6 September 2017 and SSC2 at 23:00 UT on 7 September 2017) due to the shocks of the coronal mass ejections (CMEs). However, SSC-1 is not followed by a geomagnetic storm; rather, the episode started after SSC-2 with a two-stage main phase (stage-I and stage-II), manifesting two prominent minima in SYMH (−146 nT at 01:08 UT and −115 nT at 13:56 UT on 8 September 2017), followed by the recovery phase. It can be observed from Figure  2  Interestingly, there are also enhancements of the ASYH component coinciding with the downward SYMH, describing the partial ring current influences on the nightside. The enhanced ASYH with the intensified AE index confirms the possible thermal expansion of the atmosphere to signify the compositional changes in the atmosphere through the molecular mass transport from the auroral regions to the middle and low latitudes, and the transmission of DDEF due to the auroral activity and joule heating [63]. Moreover, main phase stage-II overlaps with the recovery period of main phase stage-I during the local afternoon hours, augmenting the complex dynamics during the episode. The recovery phase follows for almost three days with gradual improvement in SYMH to attain normalcy around 14:00 UT on 11 September 2017. A detailed description of the solar flare occurrence, the arrival of SSCs, and the various phases during this particular event, with the variations of solar, interplanetary, and geomagnetic indices, were reported in the earlier literature [2,64,65]. To ascertain the local time dynamics of the ionospheric currents during the storm event, we considered the local magnetometer recordings. In Figure 3, we present the 1-min variation of the horizontal magnetic field component H (in black), superimposed with the regular Sq variation (in blue) and Diono component (in red) that corresponds to the disturbances due to the ionospheric currents registered at all the four magnetometer observatories (top 4-panels) during the period from 6-10 September 2017. The bottom panel in the Figure depicts the variation of the EEJ index during the storm period and the corresponding quiet days' variation, determined by the ΔHTIR-ΔHABG method, as presented in Equation (1). All the stations exhibit the global pattern of a momentary sharp simultaneous increase in H at the SSC, followed by a two-step downfall in the magnitude corresponding to the ring current effects, and the recovery phase evidently lasted for several hours. However, there are discrepancies in the H variation at the four latitudinal stations due to the small local time deviation effects on the electric currents. As the onset of the storm corresponds to the local dawn sector (23: 00 UT on 7 September 2017), it is not Diono but the DP2 equivalent current (PPEF transmission) in line with the eastward flipping of Ey that is simultaneously observable at all the latitudes at the beginning of the storm, and generally occurs over a period of 1 to 3 h. On the other hand, the Ddyn is a perturbation in the local time, which is characterized by a diurnal period with an anti-Sq signature. By revisiting Figure 2, the strong ASYH component at the beginning of the storm confirms the signature of the Ddyn equivalent current (DDEF transmission). Additionally, the storm-time composition changes in the neutral atmosphere at higher altitudes in the F region, associated with the transportation of the enhanced thermospheric mean molecular mass from high latitudes towards the middle and low latitudes along the meridional and diurnal winds, could not be overlooked [2,63]. In general, the thermospheric neutral compositions are mainly related to the chemistry in the F region, though there is also the electric field and neutral wind, which influence the transport in the F region. Furthermore, this can be described through the ratio of the atomic oxygen (O) to molecular species (N2 and O2). During the storm, the upward vertical wind strengthens due to the excess joule heating over the auroral latitudes, raising the molecular rich air to higher altitudes. The enhanced thermospheric mean molecular mass is transported towards the middle and low latitudes through the meridional and diurnal winds to enhance (O/N2) at low latitudes, which results in the enhanced ionospheric density (through the photoionization of atomic oxygen, e.g., [2,60,63,66,67]). However, the molecules that arrive later in time, reduce the electron density by enhancing the recombination rate with ambient N2+ and O2+ (e.g., [42,68,69]). It should be noted that DDEF needs several hours to reach the low latitudes, which could be a result of SSC on the previous day (23:43 UT on 6 September 2017). The influence of ASYH corresponding to the partial ring currents was already discussed in the earlier literature [62,70]. Further, the reversal of the EEJ index variation to below the quiet-time reference level (about −54.34 nT at 00:45 UT) during the main phase stage-I, reinforces the signature of Ddyn, which should not be misrepresented as a counter electrojet (CEJ), as the latter is a phenomenon observed on the magnetically quiet days. The magnetospheric ring current and enhanced high-latitude ionospheric currents are responsible for the reversal of Sq and EEJ during the main phase of the geomagnetic storms. The reversal, due to the Ddyn, is strongly associated with the current systems driven by thermospheric storm winds related to the joule heating of the enhanced high-latitude currents [13,57]. It is important to note that, unlike September 8, a clear signature of Ddyn (anti-Sq) is remarkable at TIR on 9-10 September, which proves that the ionospheric dynamo is certainly active on 8 September, but there is also a very strong DP2 disturbance that weakens towards the recovery phase.
It can be seen from the Figure that during the main phase stage-I, the Diono moves downwards in the negative direction to attain the minimum values at the local morning sector with almost equal behavior, but different in magnitude and local time of occurrence as seen from Table 2 [13,57]. The occurrence of the second reverse EEJ with its peat around 09:40 UT, disagreeing with the Ey fluctuation, supports the prevailing Ddyn characteristics related to DDEF, whereas the abrupt flipping in the EEJ index from a negative to positive direction attributes to the PPEF. The dominance of PPEF over DDEF, resulted in an increase in an eastward electric field to enhance the vertical drift during the local evening sector, which could probably trigger a pre-reversal enhancement (PRE)-like scenario that provided the seed perturbation for the development of the equatorial plasma bubble (EPB) during the geomagnetic storm [3,15,59,71]. The impact of the disturbance electric fields on the evening-time PRE vertical drift and ESF/EPB irregularities was examined earlier and discussed through the Digisonde-based F layer heights and vertical drift measurements in Brazil [72]. The fundamental driver for PRE generation could be the evening eastward wind in the equatorial region through the F-region dynamo process [73][74][75][76][77]. Therefore, the PRE can occur on any night with the evening wind in the equator still blowing eastward, and the normal PRE generated by the F region dynamo process can be strengthened by the PPEF [77]. Moreover, the multi-instrumental and multi-station studies by Tulasi Ram et al. [15], during the different geomagnetic storms, suggests that the dominating eastward electric fields, due to a prompt penetration in the presence of strong ambient westward electric fields at the equator, can result in spread F and scintillations in both the VHF and L-band observations, even during moderate geomagnetic storms. The enhanced dusk-time zonal electric field and electrodynamic uplift of the equatorial F layer at the sunset terminator can facilitate the EPB growth, despite the existence of background disturbance dynamo electric fields. Tulasi Ram et al. [71] presented the underlying mechanisms of the effects of PPEF on the equatorial evening PRE to trigger the EPB occurrence, in terms of the conductivity gradients in the low latitude E region at the sunset terminator, during the enhanced eastward Pedersen current, driven by convection electric fields manifested by the prevailing southward Bz condition.
To date, we do not have access to the Canadian Advanced Digital Ionosonde (CADI) operated by the Indian Institute of Geomagnetism (IIG) at an equatorial location in Tirunelveli (8.73°N, 77.70°E), India to provide the evening vertical plasma drift as a proxy for the PRE strength. The other equatorial Digisonde 4D located at Trivandrum (8.54°N, 76.87°E), India failed to record the observations after 9:15 UT. However, in the personal communication with S. Sripathi, a Professor at IIG, it was confirmed that there was a strong PRE in the observed ionograms at the equatorial ionosonde (CADI), at Tirunelveli on 8 September, against a normal PRE on the nearest international quiet day (10 September). Alternatively, to observe the effect of the penetration electric fields on the dusk-time PRE on the storm day, we plotted the estimated equatorial electric fields (5 min temporal resolution) at 78°E longitude by running the real-time prompt penetration electric field model (PPEFM; https://geomag.colorado.edu/real-time-model-of-the-ionospheric-electric-fields, accessed on 9 August 2021) [78]. Using the real-time solar wind data from the ACE satellite, the PPEFM model uses a transfer function to determine the equatorial ionospheric eastward electric field (EEF) variation with the help of a climatological model to account for the quiet day variations of the EEF. The model takes the user input time and longitude to output the following: (a) the background quiet-time electric field; (b) the prompt penetration electric fields; and (c) the total electric field (the sum of quiet-time and prompt penetration electric fields), as shown in  It should be noted that the EEF records in the model are mostly unavailable for the period from 9-10 September, to postulate the role of PPEF for any alteration in regular PRE, but the trend of EEF from 5:00 UT to 9:45 UT on 10 September indicates the nonexistence of PPEF, leading to a regular duskside PRE. Moreover, the westward electric fields, due to DDEF, was certainly active on the storm day whose existence is affirmed by the inverted EEJ on 8 September and the anti-Sq signature of Diono at the equatorial magnetometer location (TIR) on 9-10 September. Using the extensive incoherent scatter radar observations from the Jicamarca Radio Observatory, Navarro et al. [82] showed that the disturbed F-region drifts during the equinox are strongly downwards (i.e., opposite to PRE) near dusk. The radar shows that, due to the disturbance dynamo effect, the nighttime drift can be largely upwards, but at dusk it is always downwards. Therefore, it is not possible to obtain the strong PRE due to the ionospheric disturbance dynamo on 8 September. Moreover, the physical mechanisms presented by Tulasi Ram et al. [71] suggests that, despite the existence of a strong ambient westward electric field due to DDEF, the eastward electric fields due to PPEF dominate at the duskside to generate large PREs for triggering the development of EPBs at low latitudes. The plausible existence of EPBs, due to the pronounced PRE on the storm day, produce scintillations in the GNSS signals, and the poleward drifting characteristics of the scintillations are analyzed from the latitudinal variation of TEC and ROTI in the subsequent sections in this study. In Figure  It can be noticed from the diurnal TEC variations that on 7 September, there is a sharp enhancement in TEC at LCK4, a location beyond the anomaly crest in the low latitudes. However, the low latitude stations, such as HYDE and IISC, witness a slight enhancement, whereas the equatorial station (SGOC) did not show any visible response in the TEC magnitude. Despite any significant negative fluctuations in SYMH, the substantial increase in TEC at the higher low-latitude location is due to the combined effect of multiple M-class flares and the eastward penetrating interplanetary electric field (Ey), driven by frequent fluctuations in the Bz component on the day [42,83]. This supports the perturbances around the equatorial ionization anomaly (EIA) crest latitudes, following an intensification in the fountain effect under about a two-fold increased EEJ strength on 7 September 2017 (123 nT at 05:04 UT), compared to the pre-flare condition (58 nT), as confirmed in Figure 3. On 8 September, however, an increase in TEC is prominently visible from 00:00 UT to 03:00 UT (local morning) at the far low latitude station LCK4 and a sparingly noticeable magnitude at both HYDE and IISC, though the equatorial station (SGOC) did not evident any responses. This is attributed to the effects of the storm-time PPEF and thermospheric neutral compositional changes in (O/N2) in the atmosphere. From the maps by Global Ultraviolet Imager (GUVI) on board the Thermosphere Ionosphere Mesosphere Energetics and Dynamics (TIMED) satellite, the significant increase in (O/N2) on 8 September, compared to its preceding day, was reported by Bagiya et al. [2] and Chakraborty et al. [83] over the Indian low latitude. The enhanced TEC at latitudes higher than the equatorial and near latitudes during this period suggests that the storm-induced thermospheric neutral compositional changes (O/N2) could not reach the equatorial and near latitudes by this time.
Further, the temporary decrease in TEC from 05:00 UT to 06:00 UT could be associated with the DDEF and late-arriving thermospheric neutral wind to reduce the electron density with an enhanced recombination rate of the ambient ions in the atmosphere [2]. The second enhancement diurnal peak in TEC around 08:00 UT, at the higher low-latitude location LCK4 and with delays of a few hours at the locations towards the equator, can be attributed to the increased EEJ current through the M8.1 solar flare around 07:40 UT, when SYMH was transitioning from stage-I to stage-II during the episode of the two-stage main phase on the day. Although the diurnal peak of TEC increased to some extent towards the equator (SGOC and IISC), the locations at far low latitudes (HYDE and LCK4) showed a decreased level in TEC, compared to the non-flare days. This confirms that the EIA is suppressed, resulting in an increased ionization density over the equator, compared to the higher low-latitudes. This indicates the presence of a westward electric field due to strong ionospheric disturbance dynamo corresponding to Ddyn over the equator that has, in fact, reversed the daytime magnetic field to show an inverted EEJ signature at the equatorial magnetometer location. The TEC variations during the post-storm periods (9 and 10 September), however, show a reduced magnitude, compared to the storm days. The comparative analysis of TEC at all stations shows that the maximum difference in the TEC magnitude between the storm and quiet period is on 8 September for the near-equatorial stations (SGOC and IISC), whereas at the higher low-latitude stations (HYDE and LCK4), the difference is more evident on 7 September.
The interesting phenomenon observed from the ROTI variation plots in Figure 5, is the significant increased level of scintillations at all stations around 12:00 UT to 18:00 UT on 8 September. The scintillation amplitude was the highest at equatorial SGOC, followed by the IISC and HYDE stations. The severity and number of occurrences beyond the anomaly crest location (LCK4) seem to be larger than those at the HYDE station. The occurrence time of the scintillations at all the stations, and the variation of other geomagnetic and interplanetary indices, as discussed earlier, confirms the effects of a storm. It is worth mentioning that, despite high EEJ currents and solar flare emissions on the previous day (7 September), no scintillation occurrences were marked, as the geomagnetic storm was triggered towards the end of the day. On 8 September, the seeds of EPBs in the post-sunset period were prompted by a sudden flipping of a reverse EEJ current to an enhanced positive EEJ current, peaking at around 13:00 UT, referring to the local sunset terminator over the Indian longitude sector. This was indeed stimulated by the enhancement in the Ey by the southward turning of Bz and the negative departure of the SYMH, imitating a PRE scenario during the local evening hours on the day [3,4,42]. In fact, the interplay between the PPEF (increased) corresponding to the symmetric ring current and the DDEF (decreased) corresponding to the partial ring current, resulted in a dominant eastward electric field to satisfy the PRE condition [70]. From the quiet-day averaged diurnal ROTI variations, it can be observed that the magnitude did not exceed 0.3 TECU/min in most of the cases, with a minor deviation at LCK4, restricted to under 0.5 TECU/min. However, on the storm day (8 September), the magnitude of ROTI was far above the quiet-time level with the highest magnitude at SGOC crossing a value of 2 TECU/min, whereas the other locations manifested the highest ROTI level below 1.5 TECU/min, during the post-sunset period. Another scintillation condition, but with relatively lesser spread over time, was also noticed on 10 September, nearly at the same central time as the 8 September scintillation event, which was relatively more prominent towards the equatorial region (SGOC and IISC), with faded occurrences at HYDE and non-existent of the same at LCK4 (the farthest station from the equator). However, the SYMH did not show any perceptible changes, though the diurnal peak TEC was roughly aligned with the scintillation occurrence time on the day. The supporting Ey, Bz, and modeled PPEF data were not available on the day, except for an hour of recordings from around 16:00 UT to 17:00 UT that indicated the minimal fluctuations in the amplitudes of the parameters. According to S. Sripathi, a Professor at IIG (on personal communication with), a stronger PRE was evident at their operated equatorial ionosonde (CADI) on 8 September (storm day), compared to the normal PRE on the nearest international quiet day (10 September). Hence, it can be speculated that the stronger PRE on 8 September caused the EPB to extend more poleward than that on 10 September [84,85]. Nevertheless, we concentrated on the plasma bubble event on 8 September in this study and analyzed the severity of the scintillations at different latitudinal locations during the whole episode.
In Figure 6, we classify the ROTI as no TEC irregularity (ROTI < 0.25), and weak (0.25 ≤ ROTI < 0.5), moderate (0.5 ≤ ROTI < 1), and strong (ROTI ≥ 1) TEC irregularities, following the classification scheme provided by Liu et al. [86]. As the ROTI < 0.25 does not really correspond to the scintillations, the values below 0.25 TECU/min are masked out from the analysis. It can be observed from Figure 6 that among the four locations, strong TEC fluctuations occurred at the equatorial location SGOC (2.4 TECU/min), followed by IISC (1.6 TECU/min), HYDE (1.5 TECU/min), and with exceptional increased ROTI values at LCK4 (2.4 TECU/min). The deviation at the latter station is due to the perturbations around the adjoining areas of the EIA crest latitude, associated with the storm-time disturbance effects. Usually, ROTI is relatively more sensitive to the large-scale fluctuations in the ionosphere with the scale size varying above a few kilometers, whereas the S4 index is sensitive to the small-to-medium scale irregularities with an order of a few hundred meters [3,39,87]. Basu et al. [88] suggested that, as the large-and small-scale irregularities often coexist in the equatorial irregularity structures during the post-sunset hours, the ROTI measurements can be used to predict the presence of the scintillations. In brief, ROTI, S4, and phase scintillation are all indices that characterize the state of the ionospheric perturbations. ROTI is most closely related to the phase scintillation as it expresses the sharpness of the GNSS phase fluctuations caused by the occurrence of ionospheric irregularities, and is characterized by the sharp TEC gradients as, for example, demonstrated by Jacobsen and Dahnn [89], Cherniak et al. [32], and Carrano et al. [90]. Mathematically, ROTI is explicitly described by using the differences/fluctuations in the phase measurements of two GNSS frequencies [89]. Carrano et al. [90] demonstrated that ROTI is a scaled version of the structure function of the phase fluctuations imparted to the wave by the irregularities. Empirically, Jacobsen and Dahnn [89] and Aquino et al. [91] obtained the same general pattern of similar risk statistics attained by respectively using ROTI and the phase scintillation measurements at Hammerfest. The analysis of the correlations between ROTI and S4/σΦ using GPS data in low latitude Hong Kong, demonstrates the consistent relationships among the indices for the geomagnetically disturbed days, emphasizing the feasibility of ROTI observations to characterize the ionospheric scintillations [39]. Moreover, the advantage of ROTI over the standard S4 is that it is calculated based on recordings from the regular dual-frequency GNSS receivers, which are relatively more abundant than the specialized scintillation monitoring (S4 and σΦ) GNSS receivers across any region [92]. Moreover, in the African longitude sector, Okoh et al. [20] showed that there was ~67% agreement in the occurrence of plasma bubbles (detected by an airglow imager) and significant ROTI fluctuations (captured by GNSS receivers within the field of view of the airglow imager). The structural similarities and correlations between the ROTI and scintillation indices (S4 and σΦ) prove that the ROTI could be considered as a proxy for S4 to serve as more than just an indicator of the irregularities in the ionosphere [30,90]. As the specialized scintillation monitoring receivers were not available, we tried to probe the event characteristics solely with the ROTI observations. We observed that, during the storm day (8 September), the moderate TEC fluctuation occurrences were comparatively minimal at HYDE, with an intermediate number of occurrences at LCK4 and a relatively higher percentage of occurrences at the IISC and SGOC stations. The increased number of occurrences at LCK4 than HYDE was probably because the LCK4 station is located near the outer edge of the EIA crest that often suffers more fluctuations than the temperate locations. The percentage occurrence rate of irregularities at all locations is presented in Table 3. It is realized that the occurrence rate of the strong irregularities is the highest at SGOC, followed by IISC, LCK4, and HYDE, whereas that of the moderate irregularities is at a maximum at LCK4, followed by SGOC, HYDE, and IISC. The occurrence rate of the weak irregularities is however at a maximum at HYDE, followed by IISC, LCK4, and SGOC.   Figure 7 depicts the diurnal variation of the maximum absolute error (MAE) and root mean square error (RMSE) of ROTI on the storm day (8 September) with respect to the average quiet days. The peak of MAE, as well as RMSE, occurred concurrently at around 16:00 UT at SGOC, 17:00 UT at IISC, 17:00 UT at HYDE, and 18:00 UT at LCK4, demonstrating a time delay in the progression of the irregularities from the equator towards the higher low-latitude locations. It is evident from Figure 7 that the time lag between the maximum value of MAE and RMSE between the equatorial (SGOC) and far low latitude station (LCK4) is about two hours. We also noticed that the diurnal MAE and RMSE values are at a maximum at SGOC and minimum at LCK4, elucidating the severity of the irregularity at the equator that is gradually weakening towards the far low latitude locations. To further understand the spatiotemporal movement of the plasma irregularities, we probed the timestamps of the peaks in the ROTI magnitudes extracted from the respective observables of two selective PRNs (PRN 14 and PRN32) at all the locations, as shown in Figure 8. These PRNs were chosen owing to their visibility during the storm period at all the stations considered in this study. We arranged the stations from bottom to top in the Figure, based on their latitudes to realize the latitudinal drifting characteristics. It is clear from the Figure that there are two prominent peaks in ROTI, which are replicated at all stations but with varying magnitudes and increasing time lags with respect to the equator. There is also a third peak in ROTI with a relatively lesser magnitude visible only at SGOC, as observed in PRN-32 and to some extent at the IISC station. We interpreted these as irregularities that were initially developed at the equator and transferred towards the higher latitudes. The vertical dashed lines in the ROTI plot at all stations indicate the timestamps of the corresponding peaks in magnitudes. It can be observed from the ROTI plots from PRN-14 at all stations that there are two peaks in the irregularities developed around 14:50 UT and 16:08 UT at SGOC, whose corresponding appearance at a slightly higher latitude station (IISC) is felt with a delay of around 20-30 min, and approximately equivalent time lags are observed between the ROTI peaks of the successive latitudinal locations. Similar trends in the time-varying peaks in ROTI are also noticeable in the PRN-32 plots. The systematic latitudinal time delay observed in the onset of the scintillations corresponds to the velocity of the plasma bubble rise at the magnetic equator that strongly depends on the F-region dynamo electric field in the east-west direction. As the fieldaligned plasma bubbles and the associated irregularities rise at the equator, the low-latitude extremities of the bubbles shift poleward, proportionately defining the latitudinal limit of scintillation occurrences in accordance with the upper height limit of the irregularities [93][94][95]. Our observations corroborate the apparent poleward expansion of the plasma irregularities in the form of bubbles, along with the magnetic field tubes, as argued in the previous literature, through the observations obtained from different instruments, such as ground-based GNSS, incoherent scatter radar, all-sky imager, and radar interferometer, in addition to the space-based defense meteorological satellite program (DMSP) and Swarm satellite observations in the American, Asian, and East Asian longitudes [3,21,44,52,96]. The poleward drift motion of the plasma bubbles usually decreases from local evening to dawn with the velocity being the smallest around the equinox, intermediate around summer, and the largest in winter [95]. Pertaining to the zonal irregularity drifts, at this moment, it is difficult to characterize the propagation characteristics in eastwest direction due to the lack of a dense network of ionospheric observations across the longitude sector. However, we wish to indicate that earlier global and regional studies with networks of GNSS TEC/scintillation monitoring receivers, very high frequency radars, and ionosondes on the same event, demonstrated the prevailing westward drifting of the poleward portion of irregularities in the Asian and Brazilian low latitudes on the storm day (8 September 2017), as opposed to its preceding days (eastward drifting) [44,47,51], including the observations at the North-East sector of the Indian region [3].

Summary and Conclusions
In this study, the ionospheric TEC and rate of change of TEC index (ROTI) variations during the space weather events of 6-10 September 2017 were investigated with a latitudinal arrangement of four magnetometer observatories and four GNSS stations in the Indian longitude sector, spreading from the equator to the higher low-latitudes. Previous reports discuss the response of global and regional ionospheric TEC, electron density, and other ionospheric profile parameters to the combined effects of the solar flares and geomagnetic storms during this event, through instrument observations as well as model simulations. We focused on the investigation of ionospheric currents (Diono) corresponding to DP2 variation (PPEF) and Ddyn variation (DDEF) from the ground magnetometer observations, and the analysis of TEC and ROTI to probe the storm-time ionospheric TEC irregularities and scintillations associated with the occurrences of equatorial plasma bubbles during the disturbed period. The important findings from the present study are as follows: 1. The concurrence of ASYH enhancement with the SYMH/local magnetometer H component depressions indicates joule heating at the auroral zone, resulting in the probable DDEF transmission and molecular exchange in conjunction with the PPEF transmission related to magnetospheric convection, making it a complex event in the Indian local time sector.
2. The enhanced ASYH signature influenced the local post-midnight to dawn sector, while the large decrease in Diono influenced the daytime ionospheric current during the storm.
3. The sharp enhancement in the diurnal TEC variations at the higher low-latitude location (LCK4), almost no visible TEC response at the equator (SGOC), and slight enhancements at the intermediate stations (IISC and HYDE) on 7 September are associated with the disturbed equatorial ionization anomaly (EIA), due to multiple M-class flares and prompt penetration electric fields (PPEFs).
4. The significant decrement in the diurnal TEC at the higher low-latitudes and enhancements at the equatorial and nearby sectors on 8 September, confirms the delayed DDEF penetration and reduced EEJ current to suppressed EIA that resulted in the increased ionization over the equator. Additionally, contributions from the storm-time compositional changes (O/N2) in the F-region are also important to characterize the suppressed EIA at low latitudes.
5. On 8 September, the cumulative effects of the southward turning of Bz, the negative departure in SYMH, and the flipped EEJ current conceived a pre-reversal enhancement (PRE)-like scenario. This indeed resulted in a more dominant eastward electric field during the combined effects of PPEF and DDEF during the local evening sector, which was complemented by the penetrating electric field calculations through the real-time PPEFM model. Thus, the PRE seeded the development of the equatorial plasma bubbles (EPBs) in the post-sunset period, which was captured in the ROTI variations at all the stations in our study.
6. The relatively stronger PRE on 8 September caused the EPB to extend more poleward than the movement observed on 10 September, the nearest geomagnetically quiet day.
7. The higher magnitude of ROTI at the equatorial location (SGOC), reaching a level of 2 TECU/min, compared to the other low latitude region, confirms the severity of the scintillations at the equator. This was substantiated from the analysis of the % occurrence rate of the strong, moderate, and weak TEC fluctuations in the ROTI data at the locations.
8. Moreover, the largest maximum absolute error (MAE) and root mean square error (RMSE) of ROTI at the equator (SGOC) and its temporal shifts towards higher latitudes suggest the latitudinal movement of irregularities on the day. 9. The analysis of ROTI variations from two selected GPS PRNs (PRN-14 and PRN-31) suggests the severity of plasma irregularities at the equator and its temporal poleward expansion with a lag between consecutive stations, corroborating the drifting of EPBs towards farther latitudes.
In brief, we stress the analysis of the disturbance time ionospheric electric fields/currents from the ground magnetometer observations, along with the TEC and ROTI variations from the standard geodetic GNSS receiver observations, to characterize the evolution and development of the plasma bubbles and their severity on the scintillations to meet the sparse density of the specialized scintillations monitoring the GNSS receivers across the Indian longitude sector. The research is in line with the regional efforts towards exploiting scintillation evolution and mitigations over the equatorial and low latitude Indian region.

Conflicts of Interest:
The authors declare no conflicts of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript; or in the decision to publish the results.

Abbreviations
The following abbreviations are frequently used in this manuscript: