Ionosonde total electron content evaluation using International Global Navigation Satellite System Service data

In this work, a period of 2 years (2016–2017) of ionospheric total electron content (ITEC) from ionosondes operating in Brazil is compared to the International GNSS (Global Navigation Satellite System) Service (IGS) vertical total electron content (vTEC) data. Sounding instruments from the National Institute for Space Research (INPE) provided the ionograms used, which were filtered based on confidence score (CS) and C-Level flag evaluation. Differences between vTEC from IGS maps and ionosonde TEC were accumulated in terms of root mean squared error (RMSE). As expected, we noticed that the ITEC values provided by ionosondes are systematically underestimated, which is attributed to a limitation in the electron density modeling for the ionogram topside that considers a fixed scale height, which makes density values decay too rapidly above ∼800 km, while IGS takes in account electron density from GNSS stations up to the satellite network orbits. The topside density profiles covering the plasmasphere were re-modeled using two different approaches: an optimization of the adapted α-Chapman exponential decay that includes a transition function between the F2 layer and plasmasphere and a corrected version of the NeQuick topside formulation. The electron density integration height was extended to 20 000 km to compute TEC. Chapman parameters for the F2 layer were extracted from each ionogram, and the plasmaspheric scale height was set to 10 000 km. A criterion to optimize the proportionality coefficient used to calculate the plasmaspheric basis density was introduced in this work. The NeQuick variable scale height was calculated using empirical parameters determined with data from Swarm satellites. The mean RMSE for the whole period using adapted α-Chapman optimization reached a minimum of 5.32 TECU, that is, 23 % lower than initial ITEC errors, while for the NeQuick topside formulation the error was reduced by 27 %.

1 Introduction

The understanding of ionospheric behavior provides important information about the space weather. In addition, the electron content affects group and phase delays of radio waves passing through the ionosphere and impacts, among others, the Global Navigation Satellite System (GNSS). Different instruments are capable of evaluating electron density in the ionosphere, and validations among different sources of data can lead to interesting conclusions. While ionosonde instruments provide “ground truth” measures for the bottom side of the ionospheric profile and estimate the topside using an exponential decay function, ground GNSS stations receiving radio signals from orbiting satellites can provide large-scale details of the entire ionosphere structure and even plasmasphere (Huang and Reinisch, 2001; Reinisch and Huang, 2001; Jakowski, 2005; Reinisch and Galkin, 2011; Jin and Jin, 2011). The analysis proposed in this work is based on comparisons between total electron content (TEC) estimated using density profiles derived from ionograms and vertical total electron content (vTEC) from the International GNSS Service (IGS). While the IGS has its own intrinsic quality control through a ranking system (Hernández-Pajares et al., 2009), ionosonde data are evaluated by its auto-scaling system, and quality scores are assigned to each ionogram. The study was performed in the Brazilian region for a 2-year period (2016–2017), where ionosonde data from National Institute for Space Research (INPE) were available. Indeed, the plasmaspheric electron density has been considered using two different models: an adapted α-Chapman function (Jakowski, 2005) with a simple optimization and a corrected version of the NeQuick topside formulation (Pezzopane and Pignalberi, 2019).

1.1 IGS vTEC maps

IGS vTEC maps are considered to be a reliable ionospheric information product, which was achieved from integrating scientific community efforts (Hernández-Pajares et al., 2009). Such maps are generated by a combination of data from different research institutions within a method that is based on ranking different vTEC maps to compose the final product (Hernández-Pajares et al., 2009). The process begins with raw data from the GNSS ground network being acquired and sent to Ionospheric Associate Analysis Centers (IAACs) so the vTEC maps can be generated using the IONosphere map EXchange (IONEX) format (Schaer et al., 1998). To achieve a high level of quality, these vTEC maps are evaluated, and its ability to reproduce corresponding slant TEC (STEC) maps is checked. Next, a combination process takes place using a weighted mean of the available vTEC maps. The final step before making the maps available for access on the IGS server is a validation process. It compares the vTEC maps to an independent source: dual-frequency altimeters onboard the TOPEX, JASON, and ENVISAT satellites.

1.2 Ionosonde data

An ionosonde measures the returning echoes of pulse signals at a fixed location to estimate ionospheric characteristics, and the ionogram trace can be processed to result in a vertical electron density profile. The bottom side profile starts with measures at ∼90 km up to the peak of the F2 layer (foF2), around 350 km. The ionosphere topside profile, instead, is modeled using an exponential decay function. The integration of electron density in height produces an estimate for the TEC value.

Ionograms can be interpreted either manually by experts or automatically using software. The auto-scaling ionosonde data availability, rather than manual scaling, contributes to meeting practical applications (Jiang et al., 2015). Different systems were created and concentrated efforts have been applied to improve auto-scaling (Reinisch and Xueqin, 1983; Scotto and Pezzopane, 2002, 2007; Reinisch et al., 2005). Also, a standard archiving output (SAO) format was created by the initiative of the Ionospheric Informatics Working Group (IIWG) to store and disseminate auto-scaled data. Initially, SAO format considered only ionograms scaled by an automatic real-time ionogram scaler with true height (ARTIST); however, it evolved to hold scaled data from other sounder systems (Galkin, 2006).

Table 1 Correspondent CS values for C-Level flag representation, adapted from Galkin et al. (2013).

Figure 1 Locations of available ionosonde data during 2016–2017.

Ionogram quality

Several attempts had been made to verify ionogram auto-scaling system quality (Reinisch et al., 2005; Enell et al., 2016; Pezzopane et al., 2017). Early comparisons between manual and automatic scaled ionospheric parameters revealed limitations in ARTIST system performance due to the absence of quality metrics (Gilbert and Smith, 1988). Recent versions of the ARTIST system improved quality-proof methods, enhancing their results (Bamford et al., 2008; Galkin and Reinisch, 2008) and facing problems related to auto-scaling (Pezzopane and Scotto, 2007; Stankov et al., 2012).

Ionosonde data scaled by ARTIST 5 have a quality metric called the confidence score (CS). Such a metric is based on quality criteria supported by concepts of ionogram interpretation and algorithms that specify the uncertainty and confidence of scaled results (Galkin et al., 2013). The CS metric includes quality checking solutions introduced by confidence calculation schemes developed since the late 1980s: the auto-scaling confidence level (ACL) quality flag, the two-digit confidence level (C-Level), and the QualScan quality control (McNamara, 2006; Galkin et al., 2013). The estimation of the confidence score occurs during ionogram processing. Ionogram interpretation criteria consider not only analysis of extracted trace shapes, but also ionospheric conditions to compute per-point-error reduction. The CS starts with a value of 100. If an interpretation criterion is found, its per-point-error value is subtracted from the CS. To be considered acceptable for further use, auto-scaling records need to reach a CS above a predefined threshold value, which is generally 40 (Galkin et al., 2013).

The SAO format version 4 does not have the CS on its specifications, but has the C-Level representation. The two digits' range goes from 11 (highest confidence) to 55 (lowest confidence). The CS produced by ARTIST 5 can be converted to C-Level representation using Table 1 (Galkin et al., 2013).

Figure 2 Distribution of C-Level flags for ionograms during 2016–2017.

Figure 3 The same as Fig. 2 but for each ionosonde station separately.

Table 2 Ionosonde locations and correspondent closest IGS grid data location.

Figure 4 Hourly distribution of C-Level flags.

2 Methodology

Ionosonde data were obtained from the INPE database using files in SAO format (version 4) and scaled by the ARTIST auto-scaling system (version 5). Data from up to five instruments (see Fig. 1) were available at the same time for the period considered (2016–2017). Although ionosondes can generate ionograms in less than 10 min intervals, a 1 h interval between soundings was considered in this work, except for the comparisons with IGS data. In that case, a 2 h interval is used to match IGS data availability.

C-Level values were extracted from all ionograms, and despite auto-scaled ionosonde data with a CS above 40 able to be considered acceptable (Galkin et al., 2013), we chose to use only those achieving a CS above 60, corresponding to C-Levels 11 and 22. Figure 2 shows the total number of C-Level flag occurrences for the available data. Figure 3 shows the same distribution, however, considering each ionosonde station separately. It can be seen from Figs. 2 and 3 that the majority of C-Level flag occurrences are in classification levels 11 or 22. Also, more than half of C-Level flags achieve a CS above 80 (C-Level 11).

Figure 5 vTEC and ITEC daily mean variation for the period under analysis considering the mean value for each C-Level flag classification.

Considering the daily variation in ionosphere electron density, it would be interesting to analyze also the ionosondes' data quality variation within daytime hours taking into account all available data. In Fig. 4 the increase in the occurrence of C-Level 11 after sunset can be seen, while C-Level 22 decreases. The aggregation of C-Level flags 11 and 22 ionograms (used in this work) provides ionosonde data with over 1000 samples even for the period with low occurrences – see the red curve in Fig. 4.

TEC values from IGS maps and ionosondes were compared at the same date or time using the closest geographic correspondence as shown in Table 2, considering the IGS data grid (5^∘ in longitude per 2.5^∘ in latitude, every 2 h). The analysis is mainly based on the accumulation of TEC differences by applying the root mean squared error (RMSE) as defined in Eq. (1) (Chai and Draxler, 2014):

$$\begin{array}{}\text{(1)}& \mathrm{RMSE}=\sqrt{{\displaystyle \frac{\mathrm{1}}{n}}\sum _{i=\mathrm{1}}^n{e}_i^{\mathrm{2}}},\end{array}$$

where errors e_i are the differences (with $i=\mathrm{1},\mathrm{2},\mathrm{\dots },n$) between TEC values from the IGS and ionosonde, and n is the total number of values considered.

Figure 6 The same as Fig. 5 but considering only C-Level flags 11 and 22 and for each ionosonde station separately.

3 Experiments and results

Ionospheric total electron content (ITEC) daily variability for each ionosonde C-Level flag is shown in Fig. 5, and for each ionosonde station separately, considering only C-Level flags 11 and 22, it is shown in Fig. 6. Since the daily mean ITEC values from flags 11 and 22 follow the IGS vTEC data variation, they are coherent. On the other hand, it can be seen that the vTEC values are consistently higher than ITEC for the whole period and for every ionosonde. In Fig. 5, a noisy and incoherent TEC variation for flags greater than 22 in both vTEC and ITEC is noticeable. Obviously, since the results for flags greater than 22 have low confidence, they may have errors, but this reason can not be used to explain the noise in vTEC values. The data representative of higher flags is low; i.e., there is a reduced number of points and heterogeneous distributions during daytime and nighttime. Such unbalanced distribution can produce daily mean TEC representing only daytime or nighttime that during consecutive days leads to noisy curves.

Figure 7 Daily mean, considering all ionosondes: (a) ITEC and IGS vTEC differences in terms of RMSE; (b) F2 layer critical frequency; (c) maximum altitude of density profile; and (d) plasma frequency at maximum altitude.

In Fig. 6 we can observe that some ionosondes presented a lack of data for a few days or even entire months. The seasonal variation in ITEC was similar for all stations. During the autumn and winter seasons in the Southern Hemisphere we can notice a decrease in ITEC values for both years evaluated.

All the panels in Fig. 7 present daily mean values, considering the density profiles of all ionosondes. In panel (a) the RMSE when comparing ITEC and IGS vTEC is shown. The seasonal variability in TEC differences seems highly correlated with the ionization distribution along the analyzed period. Figure 7b shows the peak of plasma frequency (foF2), and we can observe the periods of high foF2 values corresponding to high RMSE. The maximum altitude used for electron density integration in ionosondes (Fig. 7c) does not change significantly, rarely reaching 900 km. Figure 7d shows the plasma frequency at the maximum altitude of the density profile, which indicates the level from where it is necessary to extend the ionosphere structure evaluation to higher altitudes. Considering the fixed scale height used in digisonde topside profile modeling, such a contribution has not been included in the ITEC calculation, since values of electron density decay too rapidly above ∼800 km, and the simple extension of maximum integration altitude is insufficient for proper comparisons to IGS vTEC. This is the main reason why the ITEC values from ionosondes are underestimated when compared to IGS values. It is well known that vTEC values from the IGS data represent the integrated electron density along the signal path between the receiver and the satellite altitude (∼ 20 000 km). Thus, this analysis is in agreement with what is shown in Figs. 5 and 6.

Figure 8 Variation of total RMSE with plasmaspheric basis density proportionality coefficient.

Figure 9 Daily mean, considering all ionosondes, of ITEC (blue), IGS vTEC (red), and density profile integration up to 20 000 km with topside reconstruction using adapted α-Chapman (orange) and NeQuick formulation (green).

Different analytical functions have been used to model the topside ionospheric density profile (e.g., exponential, Epstein, Chapman) (Nsumei et al., 2012; Pignalberi et al., 2018a; Reinisch et al., 2007). These functions and their variations may adopt fixed or variable scale height. In this work, the ionogram topside profiles were re-modeled using two different approaches that consider fixed and variable scale heights: an adapted α-Chapman exponential decay (Jakowski, 2005) that includes an optimized transition function between the F2 layer and plasmasphere and the NeQuick topside formulation with modeled scale height as a function of a corrected version of the empirical parameter H₀ (Pezzopane and Pignalberi, 2019).

3.1 Adapted α-Chapman

The adapted α-Chapman introduced by Jakowski (2005) defines the topside profile N_T as

$$\begin{array}{}\text{(2)}& \begin{array}{rl}{N}_{\mathrm{T}}\left(h\right)& =\mathit{\text{Nm}}\mathrm{F}\mathrm{2}\cdot \exp \left({\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}(\mathrm{1}-z-e^{-z})\right)+N_{P\mathrm{0}}\cdot \exp \left({\displaystyle \frac{-h}{H_p}}\right),\\ & \mathrm{where}z={\displaystyle \frac{h-\mathit{\text{hm}}\mathrm{F}\mathrm{2}}{H_T}}.\end{array}\end{array}$$

The ionograms provided F2 scale height H_T, the electron peak density NmF2 that can be derived from measured critical frequency foF2 using $\mathit{\text{Nm}}\mathrm{F}\mathrm{2}=(\mathrm{1}/\mathrm{80.6})\cdot (\mathit{\text{fo}}\mathrm{F}\mathrm{2}{)}^{\mathrm{2}}$, and peak height hmF2. According to Jakowski (2005), the plasmaspheric scale height H_p can be defined as 10 000 km and the plasmaspheric basis density N_P0 is assumed to be proportional to NmF2, i.e., $N_{P\mathrm{0}}=K\cdot \mathit{\text{Nm}}\mathrm{F}\mathrm{2}$. Using the topside reconstruction of the density profile shown above, the maximum integration height used to estimate ionosonde TEC values was defined as 20 000 km, corresponding to an approximation for the satellite orbit.

In this work, different values for the proportionality coefficient K are examined, and the optimal factor that minimizes the global RMSE is used. Figure 8 shows the ionosonde TEC differences to IGS vTEC in terms of the mean RMSE in the whole period, considering all ionosondes. When K is set to zero, Eq. (2) is reduced to regular α-Chapman decay, and the plasmaspheric slowly decaying exponential term is ignored. As K increases, the underestimated ionosonde TEC values move closer to IGS, hence reducing RMSE. However, after an optimal K, in this experiment equal to 1∕175, the plasmaspheric contribution is exceeded, increasing again RMSE.

Figure 10 Daily mean RMSE variation, considering all ionosondes, for ITEC (blue circle), adapted α-Chapman (orange triangle), and NeQuick topside formulation (green cross).

3.2 NeQuick topside formulation

The NeQuick topside analytical formulation (Nava et al., 2008; Pezzopane and Pignalberi, 2019) is based on a semi-Epstein layer describing the topside electron density profile N_T as

$$\begin{array}{}\text{(3)}& \begin{array}{rl}{N}_{\mathrm{T}}\left(h\right)& =\mathrm{4}\cdot \mathit{\text{Nm}}\mathrm{F}\mathrm{2}\cdot {\displaystyle \frac{\exp \left(z\right)}{(\mathrm{1}+\exp (z){)}^{\mathrm{2}}}},\\ & \mathrm{where}z={\displaystyle \frac{h-\mathit{\text{hm}}\mathrm{F}\mathrm{2}}{H_T}},\end{array}\end{array}$$

In this approach, the modeled scale height H_T is dependent on height h, hmF2, and also on the empirical parameter H₀:

$$\begin{array}{}\text{(4)}& H_T\left(h\right)=H_{\mathrm{0}}\cdot \left[\mathrm{1}+{\displaystyle \frac{\mathrm{100}\cdot \mathrm{0.125}\cdot (h-\mathit{\text{hm}}\mathrm{F}\mathrm{2})}{\mathrm{100}\cdot H_{\mathrm{0}}+\mathrm{0.125}\cdot (h-\mathit{\text{hm}}\mathrm{F}\mathrm{2})}}\right].\end{array}$$

H₀ can be calculated using NmF2, foF2, the propagation factor (M(3000)F2), hmF2, and the smoothed sunspot number (R₁₂) as presented by Nava et al. (2008). Pezzopane and Pignalberi (2019) proposed a new formulation for H₀ based on electron density measurements made by the Swarm satellite constellation. The formulation is

$$\begin{array}{}\text{(5)}& H_{\mathrm{0}}=\left\{\begin{array}{ll}& H_{\mathrm{0},\mathrm{AC}}+(H_{\mathrm{0},\mathrm{B}}-H_{\mathrm{0},\mathrm{AC}})\cdot \frac{h-\mathit{\text{hm}}\mathrm{F}\mathrm{2}}{\mathrm{600}},\\ & \mathrm{for}\mathit{\text{hm}}\mathrm{F}\mathrm{2}\le h<\mathit{\text{hm}}\mathrm{F}\mathrm{2}+\mathrm{600},\\ & H_{\mathrm{0},\mathrm{B}},\\ & \mathrm{for}h\ge \mathit{\text{hm}}\mathrm{F}\mathrm{2}+\mathrm{600},\end{array}\right.\end{array}$$

where two 2-D grids provide the values of H_0,AC and H_0,B as a function of foF2 and hmF2. Specifically, the grids have been calculated as the median values obtained by using the NeQuick topside formulation, IRI UP modeled values (Pignalberi et al., 2018b, c), Swarm A and C (for H_0,AC), or Swarm B (for H_0,B) electron density measurements (Pezzopane and Pignalberi, 2019). In our experiments, hmF2 and foF2 were obtained from the post-processed ionograms and H_0,AC and H_0,B grids were generously provided by Michael Pezzopane and Alessio Pignalberi from the Istituto Nazionale di Geofisica e Vulcanologia, Italy.

3.3 Comparative evaluation

Figure 9 presents a comparison, considering all ionosondes, among the daily mean ITEC (in blue), IGS vTEC (in red), and density profile integration up to 20 000 km with topside reconstruction using an adapted α-Chapman function (in orange) and a NeQuick formulation (in green). The results are for the years 2016 (top panel) and 2017 (bottom panel). All calculated values follow the reference, i.e., the IGS vTEC variations showing a semi-annual dependency with a minimum in June solstice as well as the day-to-day variability. Indeed, the TEC values obtained with the optimization of an adapted α-Chapman are very similar to the ones from the NeQuick topside procedure, and both are much closer to IGS vTEC than ITEC values.

To identify the best methodology to estimate the plasmaspheric TEC, the daily mean RMSE variation shown in Fig. 10 can be assessed. The differences to IGS vTEC using adapted α-Chapman (orange triangle) and NeQuick (green cross) approaches can be compared to the ITEC errors (blue circle) as shown in Fig. 7a. In general, the RMSE values calculated with adapted α-Chapman and NeQuick are similar. However, the lowest RMSE values along the 2 years 2016 and 2017 belong to the NeQuick criterion, as shown in Fig. 10. In addition, the mean RMSE for the whole period, using adapted α-Chapman reconstruction with the proportionality coefficient optimization, has a minimum value of 5.32 TECU, while using NeQuick topside reconstruction the error is 5.05 TECU.

4 Conclusions

This paper presented a 2-year period validation of ionosonde data, using IGS vTEC as a reference. Ionogram electron density profiles were first selected based on the confidence score and then integrated in height. As expected, ITEC values were systematically underestimated, which is consistent with an ionospheric topside modeling limitation that uses a fixed scale height, which almost neglects plasma above ∼800 km, while IGS data consider electron densities from the GNSS stations up to the satellite. This claim was supported by the examination of Figs. 5, 6, and 7. The ionogram topside profiles were re-modeled using two different approaches: optimization of an adapted α-Chapman exponential decay and the NeQuick topside formulation, based on a semi-Epstein layer with modeled scale height as a function of a corrected version of the H₀ empirical parameter. The electron density integration height was extended to an approximation of satellite orbits. Hence, as expected, for both topside reconstructions the plasmaspheric ionization contribution brought ionosonde TEC values closer to IGS observations. In our experiments, the improvement was significant for determining TEC using ionosonde data, as shown in Figs. 9 and 10. Although both procedures for calculating plasmaspheric TEC yield similar results, the NeQuick criterion shows lower RMSE values, as we can clearly see in Fig. 10.