On the Electron Temperature in the Topside Ionosphere as Seen by Swarm Satellites, Incoherent Scatter Radars, and the International Reference Ionosphere Model

Abstract

The global statistical median behavior of the electron temperature (T_e) in the topside ionosphere was investigated through in-situ data collected by Langmuir Probes on-board the European Space Agency Swarm satellites constellation from the beginning of 2014 to the end of 2020. This is the first time that such an analysis, based on such a large time window, has been carried out globally, encompassing more than half a solar cycle, from the activity peak of 2014 to the minimum of 2020. The results show that Swarm data can help in understanding the main features of T_e in the topside ionosphere in a way never achieved before. T_e data measured by Swarm satellites were also compared to data modeled by the empirical climatological International Reference Ionosphere (IRI) model and data measured by Jicamarca (12.0°S, 76.8°W), Arecibo (18.2°N, 66.4°W), and Millstone Hill (42.6°N, 71.5°W) Incoherent Scatter Radars (ISRs). Moreover, the correction of Swarm T_e data recently proposed by Lomidze was applied and evaluated. These analyses were performed for two main reasons: (1) to understand how the IRI model deviates from the measurements; and (2) to test the reliability of the Swarm dataset as a new possible dataset to be included in the underlying empirical dataset layer of the IRI model. The results show that the application of the Lomidze correction improved the agreement with ISR data above all at mid latitudes and during daytime, and it was effective in reducing the mismatch between Swarm and IRI T_e values. This suggests that future developments of the IRI T_e model should include the Swarm dataset with the Lomidze correction. However, the existence of a quasi-linear relation between measured and modeled T_e values was well verified only below about 2200 K, while for higher values it was completely lost. This is an important result that IRI T_e model developers should properly consider when using the Swarm dataset.

1. Introduction

The ionosphere is a plasma medium consisting of ions and electrons with also a neutral component that, for most of the cases, are not in thermal equilibrium. In fact, since the early satellites and rockets missions in the 1950s and 1960s [1], a bulk of evidence has been accumulated demonstrating the non-thermal equilibrium between ions, electrons, and neutrals in the ionosphere. Plasma temperature (ions plus electrons) usually exceeds the neutral one by a factor strongly dependent on altitude, time, and location; moreover, large differences among ion and electron temperatures are also observed. As a consequence, when talking about the temperature of the ionosphere, we need to differentiate between ion (T_i), electron (T_e), and neutral (T_n) temperature.

The non-thermal equilibrium of the ionosphere has a number of physical and chemical implications for the ionospheric dynamics and composition. For example, the temperature of charged particles deeply influences the distribution of electron density by driving the recombination rate in the bottomside ionosphere, and affects the diffusive equilibrium state in the topside ionosphere [2,3]. The main source of energy heating the ionosphere is the Sun: either directly, through the EUV and X-ray radiations illuminating the Earth, or indirectly, at high latitudes, through the precipitation of charged particles induced by the interactions between the solar wind and the magnetosphere-ionosphere system and the consequent Joule heating [1,4,5,6]. In the ionization processes occurring in the ionosphere, a certain amount of energy is released as an excess of kinetic energy of photoelectrons. Such photoelectrons dissipate their energy through collisions with charged and neutral particles in a cascade process that first raises the temperature of the electrons (because the transfer of energy via elastic collisions is more efficient between particles of similar mass). In a second stage, the established temperature difference drives the transfer of energy between electrons and the surrounding medium. Above the F2 region, most of the energy transfer is from electrons to ions because ion-electron Coulomb collisions are more efficient than neutral-electron ones; then, the temperature of the ions also increases above the neutrals one. As a further step, the ions lose their excess energy through collisions with neutrals and other ions, thereby raising the temperature of neutrals. As a result of these complex interactions, the entire photoelectron energy dissipates in the neutrals and it is finally transferred down in the lower layers of the atmosphere by conduction and radiation.

Among the three mentioned temperatures, T_e is the most variable because of the very low heat capacity of electrons [4]; moreover, electrons are the most affected by the solar radiation energy input and by the interaction with charged particles entering the Earth’s ionosphere following the geomagnetic field lines at auroral latitudes. The close coupling with the geomagnetic field lines and the plasmasphere produces very distinct spatial patterns in the T_e distribution [6]. As a consequence, T_e exhibits very distinctive diurnal, seasonal, spatial, and solar activity trends. This large variability makes it difficult to obtain a precise and reliable description of T_e through in-situ or remote sensing measurements. In fact, in-situ measurements taken by Langmuir Probes (LPs) and Retarding Potential Analyzers (RPAs) [1,7] on-board Low-Earth-Orbit (LEO) satellites provide a global spatial description for different local times, but they are highly dependent on the satellite orbital configuration and altitude and are limited to the satellite mission lifetime. Conversely, ground-based remote sensing techniques such as Incoherent Scatter Radars (ISRs) [7,8] can better describe the altitudinal, diurnal, and seasonal variability but for a fixed geographic location and for the very few existing ISR facilities. Moreover, the solar activity variability can be described when long time series become available. Furthermore, the cost of this measurement technique strongly limits the number and distribution of ISRs worldwide. As a consequence, for modeling purposes, both in-situ satellites and remote sensing ISRs observations are very important and complement each other.

An outstanding example is represented by the International Reference Ionosphere (IRI) model [9], which is the empirical climatological model of the ionosphere used as reference by the ionospheric community, and recently recognized as the official ISO (International Organization for Standardization) standard for the ionosphere [10]. Through the exploitation of different measurement datasets of ionospheric parameters collected by ionosondes, ISRs, LEO satellites, Global Positioning System (GPS) radio occultation, and rockets, IRI describes the hourly monthly median behavior of the electron density, electron and ion temperature, and ion composition in the ionospheric altitude range, on a global basis, for different levels of solar and magnetic activity. Due to its empirical nature, IRI evolves and improves when new datasets become available and are used to derive the empirical coefficients of its analytical formulations. In this framework, the ability to model the changes of the ionospheric ground state as a function of the different climatological conditions of both the interplanetary medium and the Earth’s seasons is also fundamental in relation to space weather science.

Despite a long series of measurements acquired by both satellites and ISRs, the major current gap in modeling and studying ionospheric T_e consists in the lack of an accurate description of its dependence on the most important geographic and geophysical parameters. There is still insufficient knowledge of the global distribution and relevant description of T_e longitudinal variations as well as of its seasonal and hemispheric asymmetries, despite some attempts [11,12,13]. A difference between long-term and short-term T_e variations with F10.7 has also been noted [14], which further complicates the modeling of T_e. Empirical models still struggle to describe small-scale structures such as the morning overshoot [15,16,17] and its dependence on solar and magnetic activity or the mid-latitude T_e enhancement [18,19]. The dependence on the magnetic activity and the choice of an appropriate magnetic activity proxy for T_e modeling is still an open problem; for example, Brace and Theis [20] found no clear dependence of T_e on the ap index, with exception of perhaps a slight systematic increase of T_e with ap at polar latitudes. Apart from modeling, there are still gaps in the fundamental theory. For instance, Pavlov et al. [21] introduced an anomalous conductivity to obtain an agreement between theory and EXOS-D data in the plasmasphere. Even very advanced physics-based models still encounter problems with reproducing observational data (e.g., the post-sunset T_e enhancement) and must be driven by measured data [22]. Another tricky point concerns the existing problems associated with the frequent T_e differences among satellite data from different missions, and between in-situ data and ground-based ISR measurements [23,24], which require in many cases a cross calibration and/or a correction of these data. Because of this, studying and global modeling of T_e has become a very challenging problem.

From this point of view, T_e observations collected by LPs on-board the European Space Agency (ESA) Swarm satellites [25] are very promising for a future inclusion in the IRI T_e model [26]. In fact, the three Swarm satellites launched at the end of 2013 and still operating provide T_e data in the topside ionosphere region on a global basis, for different local times, at a 2-Hz rate. As a consequence, the present dataset of Swarm LP T_e measurements is that large to cover different latitudes, local times, seasons, and also solar activity levels. Contrary to the electron density, T_e data have not been extensively exploited so far; only a few and rather narrowly focused studies have been published, such as [17,24]. In this context, we aimed to investigate the main climatological trends of T_e in the topside ionosphere by taking advantage of seven full years (from the beginning of 2014 to the end of 2020) of Swarm T_e observations by comparing the corresponding spatial, diurnal, and seasonal trends with those modeled by IRI. Recently, to improve the accuracy of Swarm T_e data, Lomidze et al. [24] proposed a correction to observations, based on the comparison with three ISRs located at different latitudes. In this study, we applied the Lomidze correction to Swarm data to verify its performance when compared against IRI-modeled T_e data and T_e data observed by Jicamarca (12.0°S, 76.8°W), Arecibo (18.2°N, 66.4°W), and Millstone Hill (42.6°N, 71.5°W) ISRs over the last decades.

The paper is organized as follows. In Section 2, the Swarm and ISRs T_e datasets, along with the IRI T_e model and the Lomidze correction, are discussed. In Section 3, an overall comparison between Swarm-measured (original and after applying the Lomidze correction) and IRI-modeled T_e values is first given based on the entire dataset composed by the three Swarm satellites; in addition, Swarm B and IRI T_e data are climatologically compared with the corresponding values recorded by the Jicamarca, Arecibo, and Millstone Hill ISRs. In the same section, we also describe and discuss the geographic, diurnal, and seasonal trends shown by both Swarm B and IRI data. Section 4 is devoted to the final discussion of our results and is supported by the conditioned probability density function analysis of Swarm-measured and IRI-modeled data. The conclusions and the future developments are the subject of Section 5.

2. Materials and Methods

2.1. ESA Swarm In-Situ Electron Temperature Data and Application of the Lomidze Correction

The ESA Swarm mission is a constellation of three LEO satellites launched at the end of 2013, and still operating, with the aim of studying the geomagnetic field, the electric currents in the magnetosphere and ionosphere, and the impact of the solar wind on the dynamics of the upper atmosphere [25]. The three Swarm satellites are labeled A, B, and C, and share the same design and on-board instrumentation. All three satellites were put in a circular near-polar orbit. Swarm A and C have the same orbit configuration (inclination of 87.35°, initial altitude of about 460 km, east-west separation of about 1–1.5° in longitude), while Swarm B has a different one (inclination of 87.75°, initial altitude of about 510 km). As a consequence, Swarm A and C fly in tandem, while Swarm B moves away from the couple A and C by covering different local times. Each satellite takes about 130–140 days to cover all the local times.

Among the instruments carried by Swarm satellites, we focused on LPs data [27]. LPs provide measurements of in-situ electron density (N_e) and T_e along the satellite orbit with an original sampling frequency of 2 Hz in the harmonic mode [24,28]. Moreover, 1-Hz data are also provided by interpolation of the 2-Hz dataset; it is this dataset containing time series of in-situ N_e and T_e data at 1-Hz rate that is used in this work. Specifically, LP data collected by Swarm A, B, and C from 1 January 2014 to 31 December 2020 (seven full years) are taken into account. Swarm’s data are freely downloadable at ftp://swarm-diss.eo.esa.int (accessed on 24 September 2021).

Swarm’s LP data are provided with three flags [28]: Flags_LP contains information on the source of LP measurements, and Flags_Ne and Flags_Te characterize N_e and T_e measurements, respectively. Detailed information about Swarm data quality is available at https://earth.esa.int/web/guest/swarm/data-access/quality-of-swarm-l1b-l2cat2-products (accessed on 24 September 2021). In this study, only the most reliable data were considered; specifically, those recorded by LPs in High-Gain mode. LP data were then selected with the following choice of the flags: Flags_LP = 1, Flags_Ne ≤ 29, and Flags_Te = 10 or 20.

Lomidze et al. [24] developed a correction of both N_e and T_e Swarm LP data based on the analysis of Swarm measurements matching coincident values retrieved by ISRs, COSMIC/FORMOSAT-3 GPS radio occultation, and ionosondes. For what concerns the Swarm’s T_e data correction, Lomidze et al. [24] separately considered the measurements collected when the LPs are in High-Gain mode (HG in the following) from those collected with the LPs in Low-Gain mode. Since only HG data were considered in this study, we applied only the HG Lomidze correction to T_e data. As a consequence, the Equations (4)–(6) in Lomidze et al. [24] were applied on the Swarm A, B, and C T_e datasets, namely:

$$T_{\mathrm{e}}{}_{,\mathrm{Lom}}=1.2815\cdot T_{\mathrm{e}}{}_{,\mathrm{LP}}-1167+7.293\cdot N_{\mathrm{e},}{}_{\mathrm{LP}}/{10}^4\mathrm{Swarm}\mathrm{A},\mathrm{HG}$$

(1)

$$T_{\mathrm{e}}{}_{,\mathrm{Lom}}=1.2248\cdot T_{\mathrm{e}}{}_{,\mathrm{LP}}-1047+8.548\cdot N_{\mathrm{e}}{}_{,\mathrm{LP}}/{10}^4{\mathrm{Swarm}\mathrm{B},\mathrm{HG}}_{}$$

(2)

$$T_{\mathrm{e},}{}_{\mathrm{Lom}}=1.1334\cdot T_{\mathrm{e},}{}_{\mathrm{LP}}-762+4.088\cdot N_{\mathrm{e},}{}_{\mathrm{LP}}/{10}^4{\mathrm{Swarm}\mathrm{C},\mathrm{HG}}_{}$$

(3)

In Equations (1)–(3) T_e,Lom is the Lomidze-corrected T_e, T_e,LP is the original Swarm-measured value, N_e,LP is the corresponding N_e value measured by Swarm LP (i.e., unadjusted). In Equations (1)–(3), and throughout the paper, T_e values are in K and N_e values are in cm⁻³. Equation (1) is valid for the dataset of HG Swarm A LP T_e dataset, Equations (2) and (3) for the corresponding Swarm B and C datasets, respectively.

As outlined by Equations (1)–(3), the Lomidze correction is based on a regression analysis of Swarm-measured T_e,LP with data retrieved by the three ISRs located at Jicamarca, Millstone Hill, and Arecibo from December 2013 to June 2016 (for more information on the dataset selection on which Equations (1)–(3) are based see [24]). The shortness of the dataset did not allow a characterization of possible diurnal, seasonal, geographic, and solar activity dependences of the corrections introduced in Equations (1)–(3). As a consequence, Equations (1)–(3) were applied on the entire Swarm datasets although some residual diurnal, seasonal, geographic, and solar activity dependences were somewhat expected. The application of these ISRs-based corrections allowed greatly improving the accuracy of Swarm T_e data when compared to ISRs ones (stepping from median differences between 300 and 400 K to near zero) without affecting the corresponding precision and correlation [24]. As outlined by Lomidze et al. [24], the inclusion of the term N_e,LP in Equations (1)–(3) produced a worsening in the precision of the correction (less than 60 K) but with the advantage of removing some certain non-physical patterns from Swarm T_e data.

In this work, both Swarm T_e,LP original data and Lomidze-corrected T_e,Lom data were employed. In Section 3, a statistical comparison between T_e data measured by Swarm, calibrated after applying the Lomidze correction, and the corresponding ones modeled by IRI is presented. In the same section, a statistical comparison with T_e data recorded by Jicamarca, Arecibo, and Millstone Hill ISRs is also performed.

2.2. Electron Temperature Observations by Incoherent Scatter Radars

ISRs exploit the Thomson backscatter from ionospheric electrons to retrieve the value of several ionospheric parameters including T_e [8]. Specifically, T_e profiles are recorded by ISRs in the ionospheric altitude range from about 100 to 1000 km.

T_e observations collected by Jicamarca, Arecibo, and Millstone Hill ISRs were used in this work for validation purposes. Jicamarca ISR is located right above the magnetic equator at a Quasi-Dipole (QD, [29]) latitude of 0.2°N, Arecibo (QD latitude 27.0° N) is a low-latitude station, while Millstone Hill (QD latitude 51.8°N) is a mid-latitude/sub-auroral station. Jicamarca data are from 1996 to 2020, Arecibo ones from 1974 to 2015, and Millstone Hill ones from 1976 to 2020. ISRs data are not continuous in time but the length of the datasets here used allow for a quite uniform representation of the different diurnal and seasonal trends.

The ISRs T_e data (including corresponding T_e error) were downloaded from the http://cedar.openmadrigal.org website (accessed on 24 September 2021). The altitudes of Swarm satellites fall into standard remote sensing altitude ranges of all three ISRs, which provide reliable data for decades. A similar dataset referring to T_i, for Millstone Hill was recently employed in Pignalberi et al. [30].

To perform a comparison with Swarm B-measured and IRI-modeled values, we selected T_e values recorded by ISRs in the range 510 ± 20 km. Only the most reliable ISRs’ data were selected, i.e., those with a percentage error lower than 10%. This filtering mostly affected data recorded during the night and early-morning, which are those characterized by the largest dispersion.

2.3. Electron Temperature Description by the IRI Model

In the last version of IRI (IRI-2016, [9]), two different models for T_e are implemented; they are conventionally named BIL-1995, and TBT-2012 (http://irimodel.org/, accessed on 24 September 2021). The BIL-1995 is the oldest option and is based on the works by Bilitza [31,32,33] who, through the application of the Booker approach [34] to the Brace and Theis [20] and Spenner and Pluge [35] models, describes the T_e vertical profile for the ionospheric altitude range between 60 and 2000 km. Since the IRI-2001 version [36], IRI has included an additional option based on Intercosmos data [37], that since the IRI-2012 version [38] was updated with the most recent TBT-2012 model [26] which incorporates all available satellite T_e data. The solar activity dependence, also described in Truhlik et al. [26], is mainly based on results from Truhlik et al. [39]. In the current IRI version, this additional option is called TBT-2012+SA and represents the option related to the TBT-2012 model with the solar activity term switched on. Since the TBT-2012+SA is the current default option for the T_e modeling in IRI, it has been applied in this work to perform a comparison with Swarm and ISRs measured data and will be briefly described here.

In the TBT-2012 model [26], T_e values at five fixed anchor points at 350, 550, 850, 1400, and 2000 km of altitude are modeled through a spherical harmonic expansion in a system of associated Legendre polynomials (up to the 8th order) in terms of the magnetic local time and a latitudinal coordinate based on a combination of invariant and dip latitude (invdip, [40]). Through the use of T_e measurements from LPs on-board several satellite missions (for a detailed description of the dataset see Bilitza et al. [23] and Truhlik et al. [26]), the coefficients of the spherical harmonic expansions were calculated for summer and winter solstices, and for the combined equinoxes (independently of the hemisphere); then, interhemispheric differences in magnetic coordinates were not considered in the model. The solar activity dependence was then included by selecting three solar activity ranges and describing the corresponding variability as a function of the PF_10.7 solar index [26].

The complete vertical T_e profile was obtained by applying the Booker profile function formalism [33,34]. First, a piece-wise function was built making a linear interpolation of modeled T_e values between anchor points, which produced a vertical profile divided in regions with a constant T_e gradient. Second, Epstein-step functions were applied to realize the transition from regions with different gradients and obtain the skeleton profile, which produced a continuous analytical representation of the T_e vertical profile. Finally, the skeleton profile was integrated to obtain the Booker profile function of T_e, which output a continuous and smooth description of the entire modeled T_e vertical profile. For a detailed mathematical description of the Booker profile function formalism, refer to the Appendix section of Pignalberi et al. [30].

In IRI, T_e is constrained to be equal to T_n described by the NRLMSISE-00 model [41] for altitudes lower than 120 km. An example of IRI (TBT-2012+SA default option) modeled vertical T_e profile is given in Figure 1 (green curve) where the anchor points that IRI uses to describe the entire profile are highlighted. An anchor point at a non-fixed altitude around 270 km (T_e,m in Figure 1) was introduced on the basis of Jicamarca and Arecibo ISRs data to describe a local daytime maximum in T_e characterizing mostly the low latitudes [33].

In this work, the IRI model (TBT-2012+SA option) was run on the entire Swarm satellites dataset as if it were collocated and synchronized with LPs on-board Swarm satellites. As a consequence, 1-Hz time series of IRI-modeled T_e values were produced for the time and locations probed by Swarm satellites. In this way, a one-to-one comparison between Swarm-measured and IRI-modeled T_e values was guaranteed.

3. Results

3.1. Overall Statistical Comparison between Measured and Modeled Electron Temperature Values

First of all, we aimed to investigate the effect of the application of the Lomidze correction on Swarm T_e measurements through a comparison with IRI-modeled values for the period 2014–2020 for all the three Swarm satellites.

An overall comparison between Swarm-measured (original and Lomidze-corrected) and IRI-modeled T_e values is given in Figure 2 and Figure 3 as joint histograms and distributions of residuals between measured and modeled data. From the residuals between measured and modeled values, the following statistical metrics were calculated: mean and standard deviation of the residuals, Root Mean Square Error (RMSE, Equation (4)), and Relative Root Mean Square Error (RRMSE, Equation (5)).

$$\mathrm{RMSE}[\mathrm{K}]=\sqrt{\frac{{{\displaystyle \sum }}_{i=1}^N{\left(T_{\mathrm{e},}{}_{\mathrm{measured},i}-T_{\mathrm{e},}{}_{\mathrm{modeled},i}\right)}^2}{N}}$$

(4)

$$\mathrm{RRMSE}[\%]=\sqrt{\frac{{{\displaystyle \sum }}_{i=1}^N{\left(\frac{T_{\mathrm{e}}{}_{,\mathrm{measured},i}-T_{\mathrm{e}}{}_{,\mathrm{modeled},i}}{T_{\mathrm{e}}{}_{,\mathrm{modeled},i}}\cdot 100\right)}^2}{N}}$$

(5)

In Equations (4) and (5), T_e,measured is the T_e measured by Swarm (i.e., T_e,LP) or corrected according to Lomidze et al. [24] (i.e., T_e,Lom), while T_e,modeled is the corresponding value modeled by IRI. Temperatures are in K. N is the number of values on which the statistics are calculated.

Figure 2 represents the joint histograms between T_e values measured by Swarm and modeled by IRI and gives information on the distribution of measured and modeled data as a function of the corresponding magnitude. Most of the values cluster in the reddish region, i.e., where more than 10⁴ values are counted inside each bin (the scale is logarithmic). It is interesting to note that IRI-modeled data never go beyond about 4300 K, while Swarm data (both original and Lomidze-corrected ones) show a long tail well beyond 5000 K (not shown in the figure), even though the application of the Lomidze correction reduces their magnitude. To characterize the accordance between modeled and measured values and the global effect of the application of the Lomidze correction, corresponding residuals were studied.

Figure 3 represents the statistical distributions of the residuals between values measured by Swarm and modeled by IRI. By looking at the statistical metrics reported in the upper left corner of each panel, it emerges how the application of the Lomidze correction to Swarm data improved the accuracy when compared to IRI-modeled data. In fact, the mean of the residuals decreased from about 560–650 K to about 150–250 K, with a 400 K average improvement. Conversely, the application of the Lomidze correction did not affect the precision of Swarm data, as highlighted by the standard deviation values that remained in the range 380–470 K. Nevertheless, the application of the Lomidze correction generally improved the agreement with IRI data, as demonstrated by the RMSE and RRMSE values that decreased, respectively, from about 690–780 K to about 440–530 K, and from about 27–29% to about 20–23%.

To sum up, the results shown in Figure 3 demonstrate how the application of the Lomidze correction to Swarm T_e data can improve the agreement with IRI-modeled values, without impacting the corresponding precision. However, we are aware that the comparison with a model, IRI, cannot guarantee the reliability and efficacy of the Lomidze correction; the only comparison that can truly tell whether this correction improves real T_e observations is that with other independent measurements. This is why we performed a comparison against T_e data observed by ISRs.

3.2. Statistical Comparison against ISRs Data

We implemented a statistical comparison among Swarm (original and corrected values), IRI, and ISRs T_e data for the Jicamarca, Arecibo, and Millstone Hill locations to verify the ability of Swarm and IRI to describe the diurnal and seasonal trends of T_e shown by ISRs data, and their accuracy.

Unlike Lomidze et al. [24], we did not look for correspondences (in both time and space) between ISRs and Swarm data because of the limited number of correspondences achievable (of the order of some hundreds from 2014 to 2020, for each satellite). Instead, we applied a statistical approach by comparing the diurnal and seasonal trends exhibited by both T_e values collected by ISRs and by Swarm B for the ISRs locations. In this way, we could investigate how much Swarm, and also the IRI model, can describe the diurnal and seasonal T_e trends at equatorial, low, and mid latitudes by making a comparison with Jicamarca, Arecibo, and Millstone Hill data, respectively. This comparison and those of the following section were limited to the Swarm B satellite because it is the nearest to the IRI vertical profile anchor point at 550 km of altitude. Moreover, as follows from Figure 2 and Figure 3, the statistical analysis for all three satellites was very similar, with slightly larger differences for Swarm B.

Specifically, T_e values collected at Jicamarca, Arecibo, and Millstone Hill in the range (510 ± 20) km were binned as a function of the magnetic local time (MLT, in QD coordinates) and season, and corresponding boxplots were generated. The binning in 15-minutes wide (in MLT) bins allowed describing the diurnal trend of T_e at ISRs location. The seasonal trend was studied by binning data in four periods (expressed in day of the year, doy) around solstices and equinoxes and covering the entire year:

• March Equinox:                  35 ≤ doy ≤ 125;
• June Solstice:                       126 ≤ doy ≤ 217;
• September Equinox:           218 ≤ doy ≤ 309;
• December Solstice:             doy ≤ 34 OR doy ≥ 310.

No solar activity data selection was performed because the quantity of both ISRs and Swarm data does not allow further split, and because the investigation of T_e dependence on solar activity is outside of the scope of this study. As a consequence, the results in Figure 4, Figure 5 and Figure 6 are representative of average solar activity conditions of the years covered by ISRs datasets (see Section 2.2), and of the years 2014–2020 for Swarm and IRI.

In summary, every ISR dataset was split in 384 (96 diurnal x 4 seasonal) bins and the following statistical metrics were calculated for the boxplot representation:

• 5th percentile, i.e., the lower whisker;
• 25th percentile, i.e., the first quartile;
• 50th percentile, i.e., the second quartile, representative of the median;
• 75th percentile i.e., the third quartile;
• 95th percentile, i.e., the upper whisker.

For this kind of comparison, the median is the more appropriate statistical metric because it cuts off the outliers of the distribution. The difference between 75th and 25th percentiles defines the inter-quartile range (IQR), which provides a measure of the data dispersion around the median value, while the whiskers highlight the tails of the distribution.

Figure 4, Figure 5 and Figure 6 show the boxplots for the Jicamarca (Figure 4), Arecibo (Figure 5), and Millstone Hill (Figure 6) sites. For the comparison with ISRs results, both Swarm B-measured (original and corrected with Lomidze) and IRI-modeled T_e data were selected and binned as done for ISRs data. Specifically, from the Swarm B dataset, we selected those measurements falling inside a region centered at ISRs’ geographic coordinates and extending ±5° in latitude and ±20° in longitude. We tested regions of different extension in both latitude and longitude around the ISRs’ location. The choice we made is the one balancing the amount of data needed for the statistical comparison and the preservation of the characteristic trends at the ISRs locations. Then, the selected Swarm B (original and corrected with Lomidze) and IRI data were binned as ISRs data and the corresponding median values are superposed to ISRs’ boxplots in Figure 4, Figure 5 and Figure 6.

Figure 4 shows how equatorial values at Jicamarca are characterized by the presence of a remarkable morning peak between 5 and 8 MLT, with median values from ISR exceeding 3000 K; this phenomenon is often referred as morning overshoot [15,16,17,42]. Overall, daytime values stand at about 1000–1300 K, while nighttime ones at about 700–1000 K, in a quite steady fashion. Differently from the dawn solar terminator hours, the ones related to the dusk solar terminator do not exhibit a peak but only a small ramp with T_e smoothly stepping from the daytime to the nighttime behavior. Of course, the MLTs of both the morning overshoot and the dusk ramp slightly change according to the season. The MLTs characterized by the morning overshoot exhibit the highest dispersion of data as described by the IQRs and whiskers values. The seasonal variability is really weak at Jicamarca due to the illumination conditions that do not change significantly over the year. Both Swarm and IRI can catch the main diurnal features of T_e for each season but with remarkable differences in the absolute values. IRI represents reliably nighttime values, while it slightly overestimates daytime values and underestimates values in the morning sector when the overshoot occurs. Differently, Swarm (both original and corrected with Lomidze) performs a general overestimation independently of the hour. Moreover, Swarm values do not exhibit differences between day and night, with the exception of the morning peak. The application of the Lomidze correction reduces the original overestimation of Swarm measurements, but not enough to make them match with ISR data. When comparing with ISRs data, we have to take into account that the IRI TBT-2012+SA option is based on data recorded from satellites [26], that are often limited in both space and time. In fact, the differences between IRI and Jicamarca ISR values might be related to the limited spatial and temporal resolution of the model (20 degrees in latitude and 80 minutes in MLT), resulting in the inability of the model to accurately describe the equatorial daytime dip in T_e and thus resulting in an overestimation of daytime T_e at the magnetic equator by IRI. At MLTs of the narrow and high morning peak, this limited resolution is most likely at the base of both the underestimation and some mismatches of the position of the morning peak maximum made by the IRI model. Another reason that could explain the observed differences is the fact that in-situ probes can produce higher T_e values in a non-equilibrium plasma [43] or if some sort of probe contamination occurs [44]. A possible combination of both effects thus could result in both an overestimation of about 30% during daytime and an underestimation of about 20% of the peak maximum made by the IRI model (as evidenced by Figure 4). During nighttime, when plasma is in thermal equilibrium and the limited spatial and temporal resolution of the model is not so important, the agreement between IRI and ISR is very good.

Figure 5 shows the results for Arecibo. Compared to Jicamarca, at Arecibo the seasonal variability is more marked above all for MLT hours from noon to dusk. Like Jicamarca, Arecibo values also show a remarkable morning increase with temperatures abruptly increasing from about 1000 to 2500–3000 K. However, Arecibo daytime values do not decrease as for Jicamarca but remain at around 2000 K (except for the December solstice) with significant differences in the sector 14–19 MLT, depending on the season. In fact, March and September equinoxes show a very similar trend where the first hump corresponding to the morning increase is followed by a small second hump before sunset, while for the June solstice (summer), the second hump does not appear and values steadily decrease from daytime to nighttime. The December solstice (winter) is instead characterized by quite constant values (at about 2500 K like the morning peak one) during daytime and a very distinct decrease at sunset. IRI can reliably represent the nighttime T_e at Arecibo, while some departures are evident during daytime, above all in the morning for the June solstice (underestimation) and in the central hours of the day at equinoxes (overestimation). Compared to Jicamarca, at Arecibo, the accordance between Swarm and ISR data is better, especially in the daytime sector, while nighttime values are still affected by a remarkable overestimation (even after applying the Lomidze correction). A very good accordance is found for the winter season (December solstice) during daytime. Swarm and IRI can describe both the diurnal and seasonal trends of T_e at Arecibo.

Figure 6 summarizes the results for Millstone Hill, where T_e manifests a very distinct diurnal behavior characterized by a sharp division between nighttime and daytime values. Nighttime values are of the order of 1100–1500 K, while the daytime values range between 2500 and 3000 K. Differently from Arecibo, in this case, daytime values are always quite constant with a slight tendency to maximize at noon during summer (June solstice). Winter (December solstice) nighttime values are remarkably higher than summer ones due to conjugate heating [23,45]. Except for June solstice, IRI can reliably describe the diurnal behavior of T_e at Millstone Hill. It is the same for Swarm after applying the Lomidze correction; in this case the agreement with ISR’s data is remarkable for both the representation of the diurnal and seasonal trends and the accuracy of the absolute values; only a slight overestimation is visible during nighttime particularly at December solstice.

It has to be considered that the Lomidze correction, aside from being based on the same three ISRs here studied, is heavily biased toward Millstone Hill ISR data that were more represented in the Lomidze et al. [24] dataset than Arecibo data and above all Jicamarca. As a consequence, since Lomidze et al. [24] cumulated the data of the three ISRs to derive a single calibration curve, it should be expected that the correction is more effective at Millstone Hill than at Arecibo and Jicamarca. This agrees very well with the results shown in Figure 4, Figure 5 and Figure 6. Moreover, the Lomidze correction does not discriminate between daytime and nighttime or between different seasons. Our results show that the Lomidze correction has the best performance for daytime hours and the worst for nighttime hours, with no distinct seasonal difference. This suggests that a latitudinal and day/night dependence should be included in the Lomidze correction when enough data become available in the future.

As we can see from Figure 4, Figure 5 and Figure 6, the daily T_e profile varies with the ISR geographic location. For example, at Jicamarca, which is at the magnetic equator, T_e shows a single sharp maximum around dawn that can be associated with the well-known morning overshoot. This feature is characteristic of low latitudes, and in fact is present also in the daily profiles acquired at Arecibo. Here, and at Millstone Hill (middle and sub-auroral latitudes), the increase of T_e in the dayside is due to the increased sunlight and a second hump at dusk marks the decrease of T_e in passing from day to night. At Millstone Hill, an increased T_e during nighttime is recorded at the December solstice, mainly due to the decrease of plasma density in the F layer and the following reduction of the collisional cooling effect, and also to the conjugate heating. To better highlight and describe all of these different spatial features exhibited by the topside T_e, in the next section, the main spatial trends of T_e will also be investigated, in both geographic and magnetic coordinates.

3.3. Statistical Trends of Swarm-Measured Electron Temperature Values and Comparison with IRI-Modeled Ones

The results of Section 3.1 and 3.2 demonstrate that the application of the Lomidze correction improves the agreement between Swarm data and those modeled by IRI and measured by different ISRs. This is why, in the following analyses, we focus only on Swarm B data corrected according to Lomidze et al. [24].

Swarm B data from 2014 to 2020, and corresponding IRI-modeled values, were binned as a function of the geographic coordinates with bins 2.5°-wide in latitude and 5°-wide in longitude. For each bin, the median value was considered, as representative of the main statistical trend, after discarding the outliers through the Median Absolute Deviation statistical procedure (see Section 3.3 by Pignalberi et al. [46] for a detailed description). The binning was first performed on the entire dataset and then by selecting MLT (in QD coordinates) sectors representative of the dawn (5–7 MLT), daytime (12–15 MLT), and nighttime (0–3 MLT) conditions. This was done because at a fixed altitude, the T_e diurnal variations are known to be the most important along with the latitudinal ones [20,26,39].

Figure 7 shows the corresponding results considering all the MLTs and the different MLT sectors, for both Swarm-measured and IRI-modeled data, and the corresponding percentage of the normalized residuals between them. The figure highlights how the T_e spatial distribution follows the configuration of the Earth’s magnetic field lines. In fact, the characteristic shape of the magnetic equator and of the auroral regions stand out, and this is somewhat expected because the plasma in the topside ionosphere is generally constrained by the geomagnetic field [47]. This behavior has also been recently confirmed by Giannattasio et al. [48,49] that, by using Swarm data, identified patterns in the parallel electrical conductivity (which is crucially affected by T_e) related to features of the geomagnetic field such as regions R1 and R2 [50] or the magnetic cusp [51]. They also identified variations in the parallel electrical conductivity due to MLT, QD latitude, season, solar, and geomagnetic activity. The MLTs around dawn highlight the fact that the morning peak is a phenomenon affecting only the low latitudes, which follows the magnetic equator shape and reinforces in the longitudinal sector between 80°W and 80°E (including South America and Africa) where the departure between the magnetic equator and the geographic one is the largest. At the morning peak, Swarm B values exceed 3000 K in many longitudinal sectors, well above those recorded at mid latitudes and even higher than high latitude ones. In contrast, during daytime and nighttime, the high latitudes experience the highest values; in both cases, the increase of T_e from mid to high latitudes is limited to a quite narrow latitudinal range where the spatial gradient is concentrated: at mid latitudes during daytime and at the sub-auroral and auroral boundaries during nighttime. When compared to IRI-modeled values, the largest differences emerge at magnetic equator latitudes at dawn, and at low and mid latitudes during daytime and nighttime, respectively. Overall, Swarm measurements exceed IRI values for most cases and this is in line with the results of Figure 2 and Figure 3. The excess is particularly clear at the morning peak and at low and mid latitudes at nighttime. If we recall the analysis shown in the previous section, based on the comparison with ISRs data (Figure 4, Figure 5 and Figure 6), we can say that Swarm overestimates the nighttime values; at the morning peak, Jicamarca ISR results (Figure 4) point out that the difference between Swarm and IRI is mainly due to a joint action of both the Swarm overestimation and the IRI underestimation.

As highlighted by Figure 7, the Earth’s magnetic field has a strong influence on the diurnal trend of Swarm and IRI values. Therefore, we binned T_e data also as a function of QD latitude and MLT with bins 2.5°-wide and 15 minutes-wide, respectively. The binning was first performed on the entire dataset and then by selecting different days of the year representative of the March equinox, June solstice, September equinox, and December solstice conditions, as done in Section 3.2. for ISRs data.

Figure 8 shows that the overall diurnal trend of T_e is characterized by different well-defined features at different latitudes. At low latitudes, the morning peak is the main feature and small day–night differences are visible, which confirms the results from Jicamarca ISR (Figure 4). At mid latitudes (30–60° QD latitude), no distinct morning peak is present and the main feature is the remarkable difference of absolute values between daytime and nighttime, which confirms the results from Millstone Hill ISR (Figure 6). At high latitudes, the diurnal trend is less pronounced with high values throughout the day that intensify during daytime at the auroral boundaries. The agreement between Swarm and IRI is very good for what concerns the overall diurnal pattern. The main differences characterize the morning peak hours and the nighttime hours for low and mid latitudes. It should be noted that IRI cannot describe the very well-located morning peak maximum just over the magnetic equator, and the T_e increase at the auroral boundaries during daytime. For both cases the reason has to be found in the limited spatial resolution of the model. As expected, no remarkable differences characterize the equinoxes, while the comparison between the solstices highlights seasonal differences and a different description of these by Swarm and IRI. In fact, while IRI shows a clear seasonal dependence for daytime hours, with higher T_e values in the winter hemisphere than in the summer one, Swarm data show a less marked seasonal difference.

The diurnal behavior of T_e shown in Figure 4, Figure 5, Figure 6, Figure 7 and Figure 8 is also consistent with previous results in the literature, where the observed trend was recognized to be mainly due to solar UV radiation (in the dayside) and the precipitation of particles (both in the dayside and nightside) caused by magnetosphere-ionosphere coupling [52]. Apart from the obvious contribution due to the excess of energy of photoelectrons in the dayside especially at mid latitudes, electrons, protons and alpha particles of solar origin may be injected into the ionosphere during open magnetosphere conditions, i.e., when reconnection at the magnetopause occurs, at QD latitudes of about 80° and around 12:00 MLT [53,54]. The precipitation of such particles is expected to enhance T_e in the cusp regions [51] and at lower latitudes in the pre-noon sector [55]. T_e also increases on the nightside due to particle precipitation. In more detail, plasma acceleration occurs in response to reconnection phenomena in the far-away geomagnetic tail regions, generating a three-dimensional current system modeled by the current wedge paradigm [56] that connects the tail to the nightside ionosphere in correspondence of the R2 region [50]. Interestingly, the peak of such enhancement is observed in correspondence with the boundary between R2 and the main ionospheric trough in the pre-dawn sector due to the joint action of particle precipitation and the decrease of N_e [49,57].

The diurnal variation of T_e also shows a clear dependence on the local season, which emerges from Figure 8 when considering the six-months shift of local seasons in Northern and Southern hemispheres in correspondence with solstices and equinoxes. Specifically, the increase of T_e is remarkable above all in the nightside winter, from about 20 MLT to pre-dawn at QD latitudes of about 60° and especially in the Southern hemisphere, and in the dayside summer. Such a behavior can be easily interpreted in light of the previous literature. In fact, previous studies concluded that field-aligned currents density changes from winter to summer by about a factor of 2 in the dayside due to the increase of photoelectron flux and the following increase of T_e [58,59,60,61,62]. In this case, the increase of T_e is driven by both the increase of EUV photoionization and, at very high latitudes (~80°), by precipitations in the cusp region around noon. In the nightside, particle precipitation is the dominant process at auroral and sub-auroral latitudes and becomes more and more important when proceeding from summer to winter. In particular, during winter the enhancement of T_e reaches the maximum expansion from dusk to about 08:00 MLT, especially in the Southern hemisphere. The explanation of this behavior is the simultaneous collapse of both electron and ion density in the ionospheric F layer, and the following collapse of the contribution of collisional cooling to T_e. Therefore, in the absence of mechanisms that reduce the energy of precipitating particles T_e is free to increase due, e.g., to energy exchanges with the nightside magnetosphere. Not surprisingly, there is an anticorrelation between N_e and T_e, and T_e is significantly enhanced in regions of depleted N_e, such as between region R2 and the main trough in the nightside winter; while, in contrast, it decreases in regions of enhanced N_e due to energy loss to the ions [19,48,49,57,63]. Some studies in the past tried to establish a link between N_e and T_e, since N_e is a major factor in electron energy loss due to the ionospheric plasma quasi-neutrality (density of ions is nearly equal to N_e). Generally, an anti-correlation was found. Brace and Theis [20] established an interrelationship between N_e and T_e using Atmosphere Explorer C data. Kakinami et al. [64] found a U shape in the T_e–N_e dependence using Hinotori satellite data. Their results suggest that an additional heat source has to be considered, even though the mechanism is unknown. Su et al. [65] dealt with the investigation of T_e–N_e correlation using data from the DEMETER satellite, by highlighting the variability of this correlation depending on the latitude, local time, and season. Hu et al. [66] used N_e measured from COSMIC as a driver for T_e in their neural network T_e model. Generally, they found a better agreement between their model and ISR data than between ISR data and IRI values. Thus, it is evident that the T_e–N_e relationship could help in the investigation of T_e behavior in the ionosphere but it is a complex problem itself. This is why, in this study, we did not investigate the T_e–N_e correlation using the Swarm data and we focus only on T_e.

Compared to the spatial and diurnal trends, the interpretation of the seasonal trends exhibited by Swarm data is more difficult and the differences with IRI are significant. Overall, the main seasonal differences concern the summer and winter seasons for specific MLT sectors. Nighttime T_e values at high latitudes have been found to be deeply influenced by the solar illumination conditions [67]; in fact, high values are found in the summer hemisphere lit by the Sun. In contrast, low latitudes are not affected by seasonal changes at nighttime. T_e at mid latitudes during local winter can be affected by strong heating due to conjugate photoelectrons if the opposite hemisphere is already sunlit.

4. Discussion

The results of Section 3.1. highlighted how the application of the Lomidze correction to Swarm-measured data reduces the mismatch with IRI-modeled data but, as pointed out by Section 3.3. results and the comparison with ISRs data (Section 3.2), several differences between measured and modeled data are still present for what concerns both T_e magnitude and the corresponding spatial, diurnal, and seasonal variations.

To better highlight the relation between measured and modeled data, the conditioned probability density function between them has been calculated. Specifically, we calculated the probability density of measured T_e values, both original and Lomidze-corrected, conditioned by IRI-modeled ones in 50 K-wide range around a given value T₀ (ranging from 25 to 4975 K in steps of 50 K):

$$p\left(T_{\mathrm{e}}{}_{,\mathrm{measured}}|T_{\mathrm{e},}{}_{\mathrm{modeled}}\right)|T_{\mathrm{e},}{}_{\mathrm{modeled}}\in \left[T_0-25;T_0+25\right)\mathrm{K}$$

(6)

In Equation (6), T_e,measured is the T_e measured by Swarm (i.e., T_e,LP) or corrected according to Lomidze et al. [24] (i.e., T_e,Lom), while T_e,modeled is the corresponding one modeled by IRI and p is the conditioned probability density. Equation (6) was calculated for each T₀ equal to the T_e,modeled values used to define the ranges for the conditioned probability density calculation. Conditioned probability density values are represented in Figure 9 as colored bins, while bins with p < 10⁻⁶ are in white. For each [T₀ − 25; T₀ + 25) K range, the median is calculated and represented as a black circle, and corresponding error bars are the median absolute deviation values [46].

Median values of Figure 9 were used to perform a linear regression analysis and to obtain the corresponding slope, intercept, and Pearson correlation coefficient, namely:

$$R =\frac{\mathrm{cov}\left(T_{\mathrm{e}}{}_{,\mathrm{measured}},T_{\mathrm{e},}{}_{\mathrm{modeled}}\right)}{{\sigma }_{T_{\mathrm{e},}{}_{\mathrm{measured}}}{\sigma }_{T_{\mathrm{e},}{}_{\mathrm{modeled}}}}\in \left[-1,1\right]$$

(7)

where cov is the covariance matrix and σ the variance.

Linear regression best fit lines are represented in Figure 9 in magenta, with corresponding slope, intercept, and R values reported in each panel. The analysis was restricted to T_e,modeled ≤ 2200 K because the linear relation between measured and modeled data fully breaks above this threshold, as clearly highlighted by the median values of Figure 9.

The results of Figure 9 confirm and complement those shown by Figure 2 and Figure 3. Specifically, the slopes of the linear best fit, which range between 0.80 and 0.89 for the original Swarm data, after applying the Lomidze correction range between 0.97 and 1.00, while the intercept values decrease from 920–990 K to 140–320 K. This confirms the accuracy improvement provided by the Lomidze correction. Differently, the correlation coefficient (R) was slightly affected by the Lomidze correction because this correction is substantially linear, apart from the term dependent on N_e,LP (see Equations (1)–(3)); anyway, the correlation is very high in both cases (R > 0.94). However, we have to take into account that these results are valid only in the range T_e,modeled ≤ 2200 K, where a linear dependence with a slope near 1 between measured and modeled data is evident. A non-linear dependence is clearly visible for T_e,modeled > 2200 K and in principle could be indifferently due to the Swarm LPs instrumentation, possible shortcomings of the Lomidze correction, and the IRI representation of the highest T_e values.

Recently, Pignalberi et al. [68] studied the occurrence of very high Swarm-measured T_e values (T_e,measured ≥ 6000 K), and corresponding spatial, diurnal, and seasonal variations. As expected, very high T_e values were recorded above all at high latitudes but with well-defined trends. Specifically, the occurrence of very high T_e,measured values maximizes in winter and is minimum in summer, with intermediate values at both equinoxes. Moreover, the occurrence is asymmetric between hemispheres, being highest in the Southern hemisphere, where it can reach 40% of the total observations in winter. These extremely high T_e values are not represented at all by IRI, which is limited to values below about 4300 K (see Figure 2 and Figure 9). As a consequence, the departure from the linear dependence between T_e,measured and T_e,modeled for the highest values cannot be related to the application of the Lomidze correction. The fact that the occurrence of T_e,measured ≥ 6000 K values is quite remarkable at high latitudes and that it usually happens at small spatial and temporal scales [68] drives us to the conclusion that the median climatological representation given by the IRI model cannot properly represent such very fast and steep T_e variations. This is mainly due to the limited spatial (20 degrees in invdip latitude) and temporal (80 minutes in MLT) resolution of its current spherical harmonic expansion approach [26]. Nevertheless, we point out that the Lomidze correction is based on a dataset of T_e measurements from Arecibo, Jicamarca, and Millstone Hill ISRs, which are all low- and mid-latitude stations. This means that no high-latitude measurement was considered and, most importantly, only values in the range T_e,measured < 3500 K (in HG mode) were used to derive the calibration curves (see Figure 7 in [24]). With regard to this, it would be of great benefit to include high-latitude ISR stations (e.g., Resolute Bay (74.7°N, 94.9°W), Poker Flat (65.1°N, 147.5°W,), Irkutsk (52.9°N, 103.3°E), EISCAT Tromsø (69.3°N, 19.1°E), EISCAT Kiruna (67.5°N, 20.3°E), EISCAT Sodankylä (67.2°N, 26.4°E), EISCAT Longyearbyen (78.1°N, 16.0°E)) in the Lomidze dataset to also characterize the high T_e values usually recorded at auroral latitudes. Finally, the possible overestimation of T_e values measured by Swarm LPs due to instrumental/environmental effects cannot be completely ignored. This is why, a more in-depth analysis, based on different datasets and measurement techniques, is needed to clarify the origin and reliability of the very high T_e values recorded by Swarm satellites.

5. Conclusions

In this study, the main statistical median features of the topside ionosphere T_e were investigated on a global basis through data collected by LPs on-board Swarm satellites from the beginning of 2014 to the end of 2020. Swarm T_e measured data were compared with corresponding values modeled by IRI, and with data measured by Jicamarca, Arecibo, and Millstone Hill ISRs, to evaluate the reliability of both Swarm data and IRI model in the description of the spatial, diurnal, and seasonal trends of T_e in the topside ionosphere.

We also investigated the effectiveness of the correction to Swarm T_e values proposed by Lomidze et al. [24] through comparison with IRI-modeled data and ISRs-measured ones. Specifically, we found that the application of the Lomidze correction:

• improves the agreement between Swarm data and corresponding ones modeled by IRI when the entire dataset is considered, for every Swarm satellite. This is attested by the average 400 K improvement in the mean residuals between Swarm and IRI after the Lomidze correction application;
• does not alter either the dispersion of Swarm data around the mean and the correlation between Swarm and IRI, due to the linear character of the correction;
• reduces the Swarm data RMSE from about 690–780 K to about 440–530 K, and RRMSE from about 27–29% to about 20–23%, when compared to IRI data;
• generally improves the agreement between Swarm and ISR data. The improvement is particularly evident at Millstone Hill while it is lower at Arecibo and even lower at Jicamarca. Moreover, the correction is more effective during daytime than nighttime;
• does not alter the linear relation trend between measured and modeled T_e values in the range below 2200 K, but it improves the corresponding slope (closer to 1) and intercept (closer to 0) values.

The results obtained with ISRs data suggest that the calibration of Swarm T_e data made by Lomidze et al. [24] is heavily biased towards mid latitudes and daytime hours due to the ISRs dataset they used, and to the fact that they cumulated the data of the three ISRs to derive a single calibration curve. As a consequence, the accuracy of Swarm T_e data would benefit from the inclusion of the latitudinal and diurnal dependence in the calibration procedure. Swarm T_e data are very valuable to study the relative changes of T_e, but a more detailed correction is needed to obtain more reliable absolute values of T_e especially at low and equatorial latitudes. Moreover, to better describe the high T_e values, the inclusion of data recorded by high-latitude ISRs in the Lomidze dataset is advisable.

The orbital configuration of Swarm satellites and the extension of the dataset are very tailored for the study of the global features of T_e with a spatial and time resolution hardly accessible in the past. Because of this, we were able to study the median spatial, diurnal, and seasonal trends of T_e in the topside ionosphere in a very accurate fashion. Moreover, this allowed implementing an important statistical comparison between measured data and corresponding ones modeled by IRI. This comparison is beneficial both to validate the IRI model and to unveil the most important departures between measured and modeled data for a future inclusion of the Swarm T_e dataset in the IRI model. The statistical comparison between Swarm and IRI evidenced that:

• the largest differences emerge at magnetic equator latitudes at dawn, and at low and mid latitudes during daytime and nighttime, respectively;
• IRI needs to be improved in the description of the morning peak at low latitudes. This can be achieved by increasing the order of spherical harmonics underlying the IRI Te description;
• IRI needs to be improved in summer daytime for both hemispheres, and in the description of the high Te values characteristic of daytime auroral oval latitudes;
• IRI data never go beyond about 4300 K, while Swarm data show values well beyond 5000 K, and in general Swarm values are higher than IRI ones.

The main statistical features of T_e in the topside ionosphere have been pointed out and discussed in the paper and linked to the physical phenomena affecting the T_e behavior at these altitudes, and to the existing literature on the topic. One of the most relevant results of this study is the characterization of the high T_e values recorded by Swarm LPs, which are customary above all at high latitudes. These high T_e values are not represented by the climatological description given by the IRI model. However, IRI has been developed to represent the median behavior of T_e at climatological scales, and extreme (off-median) values or very fast and steep T_e variations are out of the scope of the model. Nevertheless, a careful characterization of the high T_e values recorded by Swarm LPs deserves particular attention for both modeling and scientific purposes.

In this paper, we did not deal with T_e variations with the solar and magnetic activity. Since Swarm satellites are still in orbit and the mission will be supported for the next years, the continuous growth of the dataset will allow also studying the effects of the solar and magnetic activity variations on T_e in the years to come. Moreover, we focused on the large-scale spatial and temporal variations of T_e in order to make a fair comparison with the IRI model that describes the climatological behavior of the ionosphere. However, the very-high resolution of Swarm LPs data also allows the study of small-scale T_e features [68,69] which complement the results of this study and can be of help in the comprehension of the T_e behavior in the topside ionosphere.