Empirical Modeling of the Total Electron Content of the Ionosphere

Abstract

With the appearance of such satellite systems as GPS, GLONASS, Galileo, and others, the total electron content TEC measured by means of navigational satellites became a key parameter characterizing a state of the ionized space. In turn, functioning of navigational and telecommunication systems needs models of TEC for an estimation of accuracy of positioning, for the short-term and long-term prediction of this parameter. In this Chapter, empirical models of the total electron content are presented. The new result is their comparison. It is shown that the majority of them provide an adequate accuracy and reliability. As the basic application of TEC measurements, the problem of determination of maximum concentration NmF2 of the ionosphere with use of its equivalent slab thickness τ is considered. It is shown that existing models of τ are not global and do not provide sufficient accuracy in determining NmF2. An approach for new global model is offered.

1. Introduction

All processes on the Earth are related to the influence of the sun. Under the influence of solar radiation, the Earth is surrounded by an ionized shell, which is called the ionosphere. The role of the ionosphere in ensuring mankind activity cannot be overestimated: It softens the blow of the solar wind and provides wave propagation of various frequency ranges. The simplest example is the variety of communications systems that are affected by the ionosphere and are described in detail in [1]. Among them may be selected satellite communications, satellite navigation, including systems such as GPS, GLONASS, Galileo and others, space-based radars and imaging, terrestrial radar surveillance and tracing, and others. For the operation of navigation and communication systems, the most important parameter is the ionospheric total electron content TEC modeling capabilities and the use of which is the subject of this Chapter. TEC parameter is defined as the number of electrons in the atmospheric column of 1 m² and is measured in units of TECU, where TECU = 10¹⁶ electrons/m². Methods for measuring the TEC are described in detail in [2]. Due to the complexity and diversity of the ionospheric processes, different approaches to the modeling of ionospheric parameters were developed. Empirical (or statistical) models based on statistical analysis of the results of measurements in different parts of the globe for a long period of time are widely used. Empirical models describe some mean states of the ionosphere, so they cannot be used to describe, for example, ionospheric disturbances. However, such models are widely used because they are easy and convenient way to describe and predict the behavior of the ionospheric parameters. Considering the disturbed conditions is possible by adaptation of models to parameters of current diagnostics. The big need for such models leads to the development of various new options. In this Chapter, two methods of modeling the TEC will be considered: (1) the integration of theoretical or empirical N(h)-profile (Section 2) and (2) empirical models (Section 3). It will focus on assessing the proximity of new models to the experimental data. The presence of well-known advantages of monitoring TEC (a large number of receivers, continuous global monitoring, and data availability on the internet) has made TEC appealing to determine NmF2 (same foF2). To do this, we need to know the proportionality factor—the equivalent slab thickness of the ionosphere τ. Section 4 is devoted to simulation methods of τ.

2. Methods based on the integration of N(h)-profiles

Such methods are considered by the example of the most widely used model of the International Reference Ionosphere (IRI), which developed from the late 60s [3] under the auspices of Committee on Space Research (COSPAR) and International Union of Radio Research (URSI). Model IRI constantly modified, in particular, to improve the definition of the TEC, it has been modified three times in this century: in 2001, 2007, and 2012 [4–6], however, a satisfactory compliance with the experimental values failed, as illustrated by several examples. This paper uses a new version of the IRI-IRI-Plas [7], which includes elements not found in previous versions: (1) a new scale height of the topside ionosphere, (2) expansion of the IRI model to the plasmasphere, (3) adapting the model to measured value of the TEC. Section 2.1 includes a brief description of the model. Any new model should be tested on experimental data, so in Section 2.2, the results of testing this model according to the incoherent radar sounding, data of satellites CHAMP and DMSP, tomographic reconstructions are presented. In Section 2.3, the TEC values for new and previous versions of IRI are compared to experimental values and conditions in which modeling results are the best specified.

2.1. Description of IRI and IRI-Plas models

At present, the IRI model is the international standard for determining ionospheric parameters [8]. This is the statistical average model based on the huge amount of data of ground and satellite measurements. For the problems of wave propagation, its most important parameters are as follows: critical frequency foF2 of the F2 layer (or the maximum concentration NmF2, a linear relation with the square of the critical frequency), maximum height hmF2, propagation coefficient M3000F2 determining the maximum usable frequency MUF for the path length of 3000 km, altitude profile of the electron density N(h), the total electron content. Defining the parameters is made using coefficients CCIR and URSI, obtained by Fourier expansion according to the “1960s,” 1980s. Start parameters are the indices of solar activity. The input parameters are the date, latitude, and longitude of points on the globe. The adaptation of the model to the current diagnostic parameters (foF2, hmF2) and correction of disturbed conditions using the storm-factor SF [9] are provided. There are several basic versions of the model reflecting the most important stages of its modification: IRI79, IRI90, IRI95, IRI2001, IRI2007, IRI2012 [3–6]. The 2007 modification has two options [5]: IRI2007corr and IRI2007NeQ. The first option is a correction factor for the model IRI2001. The second option is a model of the topside ionosphere NeQuick [10]. At present, there is a new version IRI-Plas [7], which can be considered as a new modification of the model IRI, although in fact, it exists more than 12 years [11]. The main distinguishing features of this model are as follows: (1) the introduction of a new scale for the height of the topside ionosphere, (2) expansion of the IRI model to the plasmasphere, (3) ingestion of experimental values of TEC.

2.2. Testing the model IRI-Plas according to various experiments

Since one of the reasons for the discrepancies of measured and model TEC is the shape of the profile, this section presents the results of testing the model IRI-Plas according incoherent sounding radar ISR and satellites CHAMP and DMSP. Data of ISR is very seldom. We managed to gather them for the some stations on the globe. Results have been obtained for European stations StSantin, Tromso, Svaldbard, and for the American station Millstone Hill, Japanese Shigaraki, station Arecibo in Puerto Rico, from [12]. Figure 1 shows the results for the station StSantin. The first panel includes the N(h)-profile of the initial model, that is, profile, calculated by the model values of foF2 and hmF2. It is represented by symbol IRI (black circles). The symbol foF2 (squares) indicates N(h)-profile obtained by adapting the model to the experimental values only foF2. Triangle (symbol TEC) shows the profile obtained by adapting the model to the experimental values only TEC. The crosses show the profile for the model, adapted to the experimental values of the two parameters foF2(obs) and TEC(JPL). The hollow circles show the values measured by radar. One valuable source of information is the measurement of plasma frequency on satellites, flying at various altitudes. In the second panel, N(h)-profiles are compared with plasma frequency of satellite CHAMP (h ~ 400 km), in the third panel—with DMSP (h ~ 840 km).

Figure 1.

Comparison of model and experimental N(h)-profiles above the station StSantin.

The initial IRI model and its adaptation to only the TEC do not always provide a match with the profile of ISR. Coincidence is achieved only when adapting models to both parameters TEC and foF2. Similar results were obtained for the remaining stations. Reference [13] presents N(h)-profiles of Kharkov radar for conditions of low solar activity. The results for the two profiles of this series are presented in [14]. Increasing the statistics show that there may be differences, but in most cases this applies to the bottomside profile, which does not give a large contribution to TEC. Thus, despite the limited amount of data, we can conclude that the adapted profiles are quite close to the radar and satellite data at various points of the globe. The results for satellites CHAMP and DMSP are compared for the original IRI model and the model adapted to an experimental values foF2 together with TEC of one of the global maps (JPL, CODE, UPC, ESA). Square shows the plasma frequency. In cases where the flight time does not coincide with the time of TEC observation, this is indicated in parentheses. All the results show that the model and the experimental critical frequency can vary greatly, but the most important result is that through the point with the plasma frequency can pass multiple profiles, that is, measurement on separate low-flying satellites do not provide unambiguous profile. Unambiguity can be provided by use of data of simultaneous flights of two satellites [15].

2.3. Comparison of model and experimental values of the TEC

Methods for determination of the TEC have both similarities and differences. These differences lead to the differences of the TEC values for different methods. Reference [16] gives an example of the differences in the specific days on 25 and 28 April 2001 for the station Kiruna. Below in Figure 2, the TEC values are given for these days and other stations in various parts of the globe, as well as a comparison with the model values for the medians, because the models provide the medians. In the graphs representing the results for specific days, black circles show the values of the map JPL, squares—TEC of the map CODE, triangles correspond to the map UPC, crosses—the map ESA. In addition, asterisks show values for medians of the model IRI2001, circles and pluses present values of two options of model IRI2007 (corr and NeQuick), rhombs—values of the model IRI-Plas.

Figure 2.

Comparison of TEC according to the stations Juliusruh and Goosebay.

Significant differences may be seen from day to day, for example, of two days, the maximum value may be either for the map JPL (in most cases), and maps ESA or UPC. Quantitative assessment of conformity of experimental and model values can be illustrated with the help of absolute and relative standard deviation (SD) for the monthly median, considering the value of the map JPL as a reference. The results are given in Tables 1 and 2 for stations Juliusruh (54.6°N, 13.4°E), Moscow (55.5°N, 37.3°E), Manzhouli (49.4°N, 117.5°E), Goosebay (53.3°N, 60.4°W), Thule (77.5°N, 69.2°W), Ascension Island (7.9°S, 14.4°W), Grahamstown (33.3°S, 26.5°E), Port Stanley (51.7°S, 57.8°W). In the Table 1, the absolute standard deviation is given, in Table 2—the relative standard deviations.

JPL	CODE	UPC	ESA	IRI01	cor	NeQ	Plas
Julius	6	1.67	7.56	3	4.76	8.39	2.64
Moscow	4.69	3.02	7.17	2.16	3.66	6.64	2.60
Manzh	6.27	5.37	4.97	4.31	7.41	7.45	5.60
Goose	5.85	1.86	6.19	10.05	2.20	5.55	3.52
Thule	9.22	3.36	7.82	11.07	2.13	11.50	5.67
AscIs	3.81	8.02	8.67	10.57	12.82	12.61	21.04
Grah	6.11	3.50	8.35	4.52	4.45	4.19	5.26
PortS	3.41	6.50	7.02	10.09	6.81	7.60	9.34

Table 1.

Absolute RMS deviations of the different values of TEC from TEC (JPL), TECU.

JPL	CODE	UPC	ESA	IRI01	cor	NeQ	Plas
Julius	20	5.80	26.33	9	16.57	29.24	9.20
Moscow	15.68	10.09	23.97	7.23	12.23	22.21	8.68
Manzh	17.22	14.73	13.64	11.82	20.34	20.44	15.38
Goose	24.35	7.74	25.77	41.84	9.16	23.10	14.66
Thule	36.95	13.47	31.32	44.36	8.54	46.07	22.72
AscIs	1.93	7.78	9.03	10.09	13.36	13.08	17.93
Grah	19.00	10.88	25.96	14.05	13.83	13.04	16.36
PortS	12.98	24.71	26.70	38.37	25.88	28.89	35.49

Table 2.

The relative standard deviations from the values of TEC(JPL), %.

RMS differences for different maps when compared with the map JPL in a large range of latitudes and longitudes do not exceed 10 TECU, and the smallest differences were obtained between JPL and UPC. It makes 5–35%. Comparison of absolute deviations for different models shows that the best fit with the map JPL was provided by version “corr” of the IRI2007 model, for which the standard deviation does not exceed 10 TECU. The IRI-Plas model gives better results than IRI2001, except the equatorial station Ascension Island.

Thus, with a few exceptions model can provide values of TEC differences not exceeding the difference between the maps.

3. Methods of the empirical modeling

The empirical modeling of TEC, to which Section 3 is devoted, plays a huge role both for the prediction of TEC, and for testing models of type described in Section 2. For modeling TEC, basically, the method of orthogonal components [17, 18] is used; however, authors do not submit corresponding coefficients and functions. In Section 3.1, the simplest model of Klobuchara [19] is brief stated as it was unique for updating of delay of signals in an ionosphere many long years and till now is widely used for systems with single-frequency receivers though the authors using her have identified several weaknesses, for example [20]. Section 3.2 describes model [21] as an example of a model for a particular station, which should have a high degree of accuracy. The model is based on the values of biases given by the Laboratory JPL. This paper presents the results of an additional test showing that there are difficulties and for this type of models. Section 3.3 describes a new model NGM **(the Neustrelitz Global Model) [22], which in addition to the TEC model includes models of other parameters (NmF2, hmF2) [23, 24]. The authors of this model have conducted their own testing, but for definite conclusions about the effectiveness of the model, it is not enough, so the results of more extensive testing will be presented in Section 3.3. Section 3.4 describes the latest models of the TEC [25].

3.1. The model of Klobuchar

The model of Klobuchar was developed in the mid-seventies and includes one layer with infinitesimal thickness at height of 350 km. Slant TEC is calculated in a cross-point of a ray with this height. The model provides a delay estimation (in sec) for a day and night ionosphere along a vertical direction, using eight coefficients transmitted in the navigational message. The night correction is supposed to equal constant DC, fair on a global scale, in five nanoseconds (~1.5 m). The day delay is defined in the form of a cosine T^V_iono = DC + A cos[2π(t − Φ)/P] where A is amplitude, P is period, Ф is a phase depending on the geomagnetic latitude of under ionospheric point, T^V_iono is a vertical delay. Eight transmission coefficients of two polynomials of 3° include four coefficients for A and four coefficients for P. Controlling ground segment of GPS updates these coefficients according to the season and the level of solar activity. Phase Ф in the argument of the cosine is constant and equal to 14 h. If the argument [2π(t − Φ)/P] is greater than π/2, the cosine becomes negative, and T^V_iono includes only a constant DC. Delay along the line is calculated as T_iono = F * T^V_iono where F = 1 + 16(0.53 − El)³, El—the angle of elevation. Taylor expansion of the equation for T^V_iono gives an expression for the model of Klobuchar.

This model serves as a standard when comparing the effectiveness of the correction of the ionospheric delay.

3.2. Taiwan empirical model of TEC

The majority of empirical TEC models of new generation are statistical. In reference [21], some models were built for a single point (24°N, 120°E) using the biases of JPL laboratory from 1998 to 2007 for quiet geomagnetic conditions (Dst > −30 nT). Input parameters are local time (LT), day of the year (DOY), the index of solar activity (F10.7 or EUV). Since the choice of the best index from their huge number is not obvious, the authors [21] investigated the effect of this choice on the final result. Set of indexes included the average values of F10.7 and EUV for the period from 1 to 162 days. It most closely matches the model and experimental values of the daily TEC caused EUV, which provided standard deviation RMS = 9.2TECU compared with 15-day moving medians with their RMS = 10.4TECU and evaluation for IRI2007 version NeQuick RMS = 14.7TECU. Daily values of index EUV (0.1–50 nm), obtained by Solar Heliospheric Observatory SOHO, were taken on a site http://www.ngdc.noaa.gov/stp/SOLAR/ftpsolarradio.html. The functions have periods of variations in 6, 8, 12, and 24 h with a dominant period of 24 h. Synodic period, causing variations in solar index about 27 days, was clearly identified in the spectrum of the TEC variation, as well as semiyear variations of 183 days, year (332 days), and longer (609 days). TEC is the product of three functions of three parameters (EUV, DOY, and LT). The function describing the dependence on solar activity uses a cubic approximation. The factor of the seasonal dependence includes three harmonic multipliers, daily course includes four harmonics. DOY parameter is normalized by the number of days in a year. The coefficients αn are presented in [21]. It should be noted that these coefficients are given in truncated form in the article, and this can lead to errors. Examples of correspondence between model and experimental values are given in Figure 3 (calculations were performed using the full set of factors, kindly provided by one of the authors [21]). The results for August 2002 presented in [21] and our calculations coincide. This makes it possible to obtain the results for other months of 2002 and for the same months of low activity.

Figure 3.

Comparison of model and experimental TEC for the Taiwan model near the peak of solar activity.

It is perfectly visible seasonal variations of TEC at the given latitude and full compliance for autumn and winter months. In the spring and in the summer, the model underestimates values. RMS range is 4–14 TECU. The relative standard deviation amounts to 6–18%. For a minimum of solar activity, TEC values were 2–3 times less than at the maximum of solar activity. The model can both underestimate and overestimate the experimental values. The range of the absolute deviation was 1–10 TECU. If we compare these results with a 50% rating for Klobuchar model [19], we get improvement in 2–5 times. Traditionally, the comparison is made for the medians, because the model is median, and the definition of instantaneous values is not possible. But the model [21] provides instantaneous values. Figure 4 gives a comparison of the daily model and experimental values for August 2002.

Figure 4.

Comparison of daily model and experimental values of TEC for August 2002.

Good correspondence of dynamics of TEC variations that are confirmed by quantitative estimations of absolute deviations 6.4 TECU is visible. RMS of absolute deviations is 8.3 TECU, and relative deviations are 16.4%.

These results show high efficiency of the model and a way of its construction. It can be used for testing of other models.

3.3. Empirical model NGM

The NGM unlike the Taiwan model is global. Its structure can be described as follows. Model TEC (NGM) is given by product of five multipliers: TEC = Ф1 * Ф2 * Ф3 * Ф4 * Ф5 [22]. Each multiplier reflects dependence on the certain physical factor and is calculated with use from two to six coefficients. Coefficients are defined by a method of least squares superposition on experimental data for some years. Multiplier Ф1 describes dependence on local time LT, that is, on an zenit angle of the Sun, and includes daily, semidiurnal, 8-day variations. It is calculated with use of five coefficients. Multiplier Ф2 describes annual and semi-annual variations, using two factors. Multiplier Ф3 includes dependence of TEC on a geomagnetic latitude. The model includes equatorial anomaly in latitudinal course of TEC. Dependence on the solar activity is described by index F10.7. The model for NmF2 [23] includes 13 factors. The maxima of a daily course of TEC and NmF2 are fixed at LT = 14. The model for hmF2 [24] includes four factors. Data-ins are: doy—number of day in a year, D(21.3)—number of day on 21 March in a year (80 for not leap, 81—for leap), F10.7—monthly average value of index F10.7 for the concrete day, ϕ—a geographical latitude of a point, λ—a geographical longitude of a point, ϕm—a geomagnetic latitude of a point, sign σ = ϕ/|ϕ|, LT(array)—an array of local times. TEC in various latitudinal zones strongly differ on the properties; therefore, results are presented separately for each zone. Comparisons for a middle-latitude zone are illustrated on an example of European station Juliusruh. As all models are median, comparison is performed for monthly medians. Typical examples are given in Figure 5 for the conditions close to a maxima (2001) and minimum (2007) of solar activities. The first drawing shows absolute deviations |ΔTEC(med)| for 2001. In this case, comparison is carried out for two versions of the IRI model: IRI2001 and IRI-Plas to estimate, whether can improve model IRI-Plas results of the previous versions. The second drawing gives relative deviations σ(TEC(med)). Next drawings concern to 2007.

Figure 5.

Examples of comparison of results in the conditions of a maxima (2001) and a minimum (2007) of solar activities for middle-latitude station.

There are months when the NGM model provides the better results than both IRI models; however, in winter months, all models do not provide necessary correspondence with experimental data. The particular interest is represented by results for high and equatorial zones. In some papers, for example [26], the possibility of use of the IRI model in high latitudes was shown. If in middle latitudes, the results of comparison can be similar for several stations, in high latitudes due to a strong variability it is possible to expect differences; therefore, results in Figure 6 are given for several stations with various coordinates. It has appeared that results for high-latitude stations not strongly differ from results of middle-latitude station with some increase of deviations with a latitude.

Figure 6.

Comparison of daily courses of foF2 and TEC medians for high-latitude stations in the conditions of low (2007) and high (2001) solar activities.

Maximum deviations concern to the IRI2001 model, illustrating advantages of models NGM and IRI-Plas before this model. At comparison of results for models IRI-Plas and NGM, advantage has the IRI-Plas model. In the conditions of low solar activity for all stations, there are periods when deviations for the NGM model are less than for the IRI model. Absolute deviations are lower in maxima of solar activity, and relative deviations are higher. The big deviations are inherent in all models in winter months. For a low-latitude zone, results are illustrated on an example of the data of station Athens (Figure 7), for equatorial—Ascension Island (Figure 8).

Figure 7.

Comparison of annual dependences of TEC medians for various models in 2001 and 2006 for station Athens.

Figure 8.

Comparison of annual dependences of TEC medians for various models in 2001 and 2006 for station Ascension Island.

For low-latitude station Athens, the NGM model has not advantages before remaining models, but for the equatorial station Ascension Island, the big advantages are visible; however, it is not obvious that the same results will be for other equatorial stations. More detailed results are presented in [27]. Results for separate stations yet do not give an overall picture. It is interesting to reveal behavior of deviations depending on a latitude. Results are given in Figure 9. They concern to certain month and a longitudinal zone: European (April 2002 and July 2004) and American (April 2002 and November 2003). Cases were selected on the basis of the greatest number of stations.

Figure 9.

Examples of latitudinal dependences of medians for various conditions.

Graph shows ranges of latitudes in which this or that model has advantages; however, for other conditions results can be others. The best results in most cases concern to the IRI-Plas model. It is important that in most cases relative deviations do not exceed 20%. This is comprehensible result.

3.4. The Bulgarian global empirical model of TEC

Process of a model development goes continuously. This is an additional confirmation of an urgency of this process. The model [25, 28], on the one hand, is most physically justified, on the other hand, by estimations of authors of [25], their model is two times more exact, than the NGM model. In references [25, 28], it was developed not only the TEC model, but also the model of its error [28]. Difference from the NGM model is the taking into consideration not only the components caused by sunlight, but also regular wave structure of the tidal nature acting from the lower atmosphere. The model is constructed according to the map CODE for 1999–2011. Sliding medians are calculated by means of a 31-day window, and the median is assigned to central day of a window, that is, 16 numbers. Sliding medians are calculated independently for each point of the chosen grid. Daily data sets for each modified geomagnetic latitude, a geographical longitude, and time UT are obtained. One of the reasons of use of the modified geomagnetic latitude instead of geographical just also is the account of influence of the lower atmosphere and a thermosphere as this influence depends on a configuration of force lines of a magnetic field. The difference between geomagnetic and geographical frames generates an additional tidal response of the ionosphere. Spatial-temporary structure of TEC is represented in the form of [29]: TEC = Φ1 * Φ2 * Φ3. Function Φ1 is represented in the form of expansions in Taylor series, Ф2 and Ф3—in Fourier series. As parameter of solar activity, it is chosen not only index F10.7, but also its linear velocity KF. The seasonal factor includes 4 harmonics: the annual, semi-annual, 4 and 3 monthly. The daily variability includes three components: mean value TEC, a part describing solar components, and a part describing stationary planetary waves. The model includes 4374 constants which are defined by a method of least squares. The number of included components in Taylor’s and Fourier’s expansions is defined by a trial and error method with use of the following criterion: Components of higher order are rejected if their inclusion improves an error only in the third sign. In papers [25, 28], detailed investigation of deviations of model TEC values from observational ones by means of estimations of an average (regular) error (ME), a mean squared error (RMSE), standard deviation errors (STDE) was conducted. For all array of the used data, the following estimations are obtained: ME = 0.003TECU. For such value of ME, the other values are RMSE = STDE = 3.387TECU. These estimations are compared to estimations for the NGM model of TEC [22]: ME = −0.3TECU, RMSE = 7.5TECU. Thus, the Bulgarian model has a smaller error in two times. However, it is noticed that both models are climatological, that is, describe an average condition in quiet geomagnetic conditions, and the difference in number of coefficients (12 against 4374) is underlined. Authors [25] absolutely fairly do not consider a higher number of coefficients as a model shortage as these factors are calculated once; however, they are unavailable. Coefficients of the NGM model were published and can be used by any user. In turn, we can notice that in an error distribution of any model there are “tails” and it is important to define, which latitudinal zones and which conditions of solar activity they concern to. As any model cannot work equally well in all latitudinal zones and meet the possible requirements because of limitations of the approaches, the used data, distinction of physical processes, testing of models does not cease to be an actual problem.

In conclusion of this section, we will note reference [30] in which some methods were compared at an estimation of positioning accuracy. One of them is based on the TWIN model [31]. This model was used in [32] for correction of ionospheric delays in single-frequency receivers and has yielded results of positioning accuracy better than the Klobuchar model and standard global maps of TEC. Figures lay within 1–10 m. These figures and other results of the reference [32] show that basic distinctions between accuracies of positioning for these models are not present, but the TWIM model is constructed by data for low solar activity. The example for high activity is given in the paper [30] mentioned in [32]. In it, results of six methods were compared: (1) not corrected delays, (2) model [19], (3) IRI2001, (4) the prediction for 40 min by results of tomographic reconstructions, (5) a method of tomographic reconstructions MIDAS, (6) a two-frequency delay (it was used as a true delay). The basic emphasis was made on an estimation of a possibility to use a tomographic method for increase of positioning accuracy. As advantages, it is indicated a possibility of an obtaining of the data in real time though it demands presence of an infrastructure which does not exist yet in many regions. In methods 4–5, tomographic maps of N(h)-profiles were used for delay calculation. Results were obtained for the European zone and four stations: MAR6, GOPE, VILL, ANKR for several days of year 2002, and period 21 October–4 November 2003. By results of paper [30], it is possible to make Table 3 in which results are given in order of accuracy increase.

	Non-comp		Klobuchar		IRI2001		Forecast		MIDAS
	Mean	Max	Mean	Max	Mean	Max	Mean	Max	Mean	Max
MAR6	10	18	4	10	3	6	1.5	3	0.5	1.5
GOPE	11	20	3	9	3	6	1.5	3	0.5	1.5
VILL, ANKR	13	20	4	9	3	6	1.5	3	0.5	1.5

Table 3.

The positioning accuracy provided by various methods, by results of [31], in m.

Feature of reference [30] is the estimation of the positioning accuracy during the strongest geomagnetic perturbations which have paralyzed work of many satellite systems [33], however in [30] optimistic enough results are obtained at use of method MIDAS though conclusions have ambiguous character.

4. Use of a median of the equivalent slab thickness of the ionosphere τ for determination of NmF2

The presence of known advantages of TEC measurement (a great number of stations, continuous global monitoring) has made TEC attractive to calculation of NmF2 (the same foF2) in a global scale. For this purpose, it is necessary to know a constant of proportionality—an equivalent slab thickness τ of the ionosphere. Values of τ(IRI) are most often used [20, 34]. The surprising fact: There is a considerable quantity of publications in which morphological features of τ(obs) are described, but nobody has guessed to use it for calculation of NmF2. Probably, it was because practically nobody compared τ(IRI) and a median τ(med) of observational τ(obs). In Section 4.1, comparison of two types of τ is carried out and deviations of the calculated foF2 values from experimental magnitudes foF2(obs) are obtained. In Section 4.2, effectiveness coefficients Keff of use of a median τ(med) in comparison with τ(IRI) have introduced. Values of Keff will be presented as for separate stations of globe, and on a global scale, and it is shown that these coefficients for τ(med) are always higher than 1 unlike coefficients for τ(IRI). To use τ(med) on a global scale, it is necessary to have its model. The mention of a possibility of construction of the τ model practically does not meet in papers. Some variants are possible: (1) construction of superficial function of kriging using values of τ(med) in several points, (2) two-parameter model on the basis of hyperbolic approximation τ(hyp) = b0 + b1/NmF2, (3) the NGM model, (4) the IRI-Plas model. The doubts are stated in the paper [35] concerning the first variant, the model of the second variant is introduced in Section 4.3. Results of testing of the third and fourth models were given in Section 3.4 and in [27].

4.1. Comparison of model and observational values of τ

Assimilation of TEC into different models became one of the directions of ionospheric modeling. Results of TEC assimilation have a direct relation to use of models in real time. Use of observational TEC(obs) together with an equivalent slab thickness τ(IRI) to calculate foF2 values can be considered as the most simple procedure of assimilation. Magnitude of τ(IRI) is calculated from a relation τ(IRI) = TEC(IRI)/NmF2(IRI) where parameters TEC(IRI) and NmF2(IRI) are medians; therefore, τ(IRI) can be considered as a median. Using of values TEC(obs) provides values NmF2(calc) = TEC(obs)/τ(IRI) and foF2(τIRI) = 8.97 *SQRT(NmF2(calc)). In reference [36], it is offered to use a median τ(med) for calculation of foF2. The following expressions are used: τ(med) = med(TEC(obs)/NmF2(obs)), NmF2(calc) = TEC(obs)/τ(med), foF2(τmed) = 8.97 * SQRT(NmF2(calc)). Thus, differences of foF2 values calculated by two ways are defined by differences between τ(IRI) and τ(med). Though there is a considerable quantity of publications in which morphological features of τ(obs) are described [37, 38], practically, there are no papers in which values of τ(IRI) and τ(med) are compared. Especially, there are no papers comparing results of use of τ(IRI) and τ(med) together with observational TEC(obs) for foF2 calculation. In the given section, such comparison is carried out. For comparison of these values, effectiveness coefficients have introduced. Effectiveness coefficients are defined by means of deviations of calculated foF2 from the observational values. |ΔIRI| = |foF2(obs) − foF2(IRI)| is a difference between instantaneous values for the IRI model and experimental values. Monthly averages were calculated. This difference stays in numerators of effectiveness coefficients. The deviation |Δτ(IRI)| = |foF2(obs) − foF2(τIRI)| defines a difference between the values calculated with use τ(IRI) and experimental foF2(obs). The deviation |Δτ(med)| = |foF2(obs) − foF2(τmed)| defines a difference between the values calculated with use τ(med), and observational foF2(obs). Coefficient KτIRI = |ΔIRI|/|Δτ(IRI)| is the effectiveness coefficient for τ(IRI). Coefficient Keff = |ΔIRI|/|Δτ(med)| is the effectiveness coefficient for τ(med). Thus, the efficiency coefficients indicate in how many times increases consistency between the calculated and experimental values in these two cases. In reference [39], differences between τ(IRI) and τ(med) are illustrated for stations in various regions of globe: Juliusruh, Goosebay, Thule, Grahamstown, Ascension Island in a daily course for July and December of several years from 2002 to 2010. In Figure 10, illustration of differences is given on an example of July and December for reference station Juliusruh, map JPL and moderate level of solar activity (2004).

Figure 10:

Illustration of differences between model and experimental values of equivalent slab thicknesses for the middle-latitude station Juliusruh of the European region.

As it is known, observational values of TEC form the whole set of maps: JPL, CODE, UPC, ESA, La Plata, IONOLab TEC, RAL, and others. The corresponding values of τ are calculated for all these values. They can strongly differ. Differences between maps can be illustrated on an example of the data of reference [40]. Considering that τ does not depend on a latitude, on graphs of work [40], all values are given in a range of latitudes and longitudes of the European zone; therefore, it is possible to see an essential scatter of values on some graphs. These graphs are of interest for us, as they concern to period of low solar activity (2007–2010) and give the chance to compare experimental τ with τ(IRI). Calculations for all 12 cases of work [40] have shown good correspondence with map JPL. Figure 11 show results of comparison of τ(IRI) with τ(JPL) and τ (CODE) for station Juliusruh and July and December 2008. Period 2006–2009 was characterized by extremely low values of solar spots that have led to the increased errors of modeling [41]. Lack of latitudinal dependences were marked by other authors also for the European region; however, it is an essentially important point for calculation of foF2 using experimental values of TEC and medians of τ(med).

Figure 11.

Comparison of behavior of monthly medians of experimental and model values τ in a daily course on an example of the European region in low solar activity.

If latitudinal dependence of τ(med) did not exist, τ(med) of any ionospheric station could be used in all region for calculation of foF2, for example, by means of operative system Local Ionospheric Electron Density Reconstruction (LIEDR) which carries out monitoring of τ [35]. “Lack” of latitudinal dependence of τ is illustrated in Figure 12 for the stations lying in a range of latitudes used in [40] for January and July 2008 for maps JPL and CODE.

Figure 12.

An illustration of differences of τ for stations with various latitudes.

It is important to investigate what deviations of foF2 leads use of this τ with such various behaviors in a daily course to. Differences between experimental foF2 values and values calculated by means of medians τ for various maps are presented in Figure 13 for three stations of European region Tromso, Juliusruh, Athens.

Figure 13.

Deviations of calculated foF2 from experimental values for various maps.

Deviations for the IRI model are less 1 MHz. Deviations for medians of τ of global maps are 2–3 times less. It is necessary to note two important facts. For the high-latitude station Tromso, deviations for τ(CODE) exceed even deviations for the IRI model and they are maximum at night when TEC values are small. It can testify that the method of the CODE map can work insufficiently well at low TEC values. The second fact is connected rather with small differences between foF2 calculated by means of different maps and corresponding τ(med). It is a result of a good adjustment of τ(med) under TEC.

4.2. Coefficient efficiency of τ(med) usage

Since the efficiency coefficients of τ(med) are connected with the deviations, the results are given for the coefficients, and for deviations. Figure 14 shows the deviations and coefficients of efficiency for τ(IRI) and τ(med) for the Juliusruh station. The black dots on the figures of deviations concern to the IRI model, blue circles—to the usage of traditional τ(IRI), red dots refer to the usage of the median τ(med). In all cases, the new τ provides the smallest deviation, that is, most accurately determines the critical frequency. In the right-hand parts of figures, efficiency coefficients are given for the two cases. The black line shows the points K = 1. If the ratio is equal to 1, this indicates that the usage of the equivalent slab thickness and the experimental value of the TEC provide the same results as the model itself without the involvement of the TEC. If the ratio is greater than 1, then the use of TEC gives better results than the model. If the ratio is less than 1, the use of TEC worsens results compared with the model.

Figure 14.

Deviations and coefficients of efficiency for τ(IRI) and τ(med) for the station Juliusruh.

Figure 15.

Deviations and coefficients of efficiency for τ(IRI) and τ(med) for the Athens.

Figure 16.

Deviations and coefficients of efficiency for τ(IRI) and τ(med) for the Thule.

Figure 17.

The global picture of deviations and efficiency coefficients for April 2014 and March 2015.

The results are shown for all months, and it would be possible to see seasonal variation, but in this case, we are not interested in such details. More importantly, that there are too many cases for τ(IRI) when the ratio is less than 1, which means that the use of TEC worsens results. For the Athens station, this situation exists almost always (Figure 15). It is surprising, but the best results were obtained for the Thule station (Figure 16). Figure 17 give results on a global scale for April 2014 and March 2015.

These results lead to the following conclusions: (1) use of the TEC(obs) does not always improve coincidence between the calculated and experimental values of foF2 in comparison with the initial IRI model, (2) use of τ(med) leads to more exact values of foF2, (3) the coefficient Keff is always higher 1. Essential diurnal and seasonal variations are not visible. In the solar cycles including periods 2001–2011, 2002–2012, dependence of Keff on solar activity is characterized by maxima 2.5–3 in 2001–2002 and by values in a range 1.5–1.7 in remaining years.

4.3. About a global model of τ(med)

The mention of the possibility of constructing a model of τ practically does not occur in the articles, but in recent years articles on the use of TEC to determine NmF2 began to appear using equivalent thickness τ of the ionosphere. This shows the urgency of this task. In [42], the authors proposed the use of its two Neustrelitz models for the TEC and NmF2 [22, 23] to determine foF2, but without sufficient testing. These models can be named NGM (from the Neustrelitz Global Model). That is why, so much attention has been paid to comparison τ(NGM) with τ(IRI) and τ(med) in [27] and in Section 3. Authors [43] have reproached researchers that they are developing a model of the ionosphere but not a model of τ; however, authors [43] have done nothing. The latest step has been made in [44], where a model of the average values of τ was developed by using the Fourier series expansion according to the TEC and foF2 for 21 stations. Authors have taken monthly averages of the global map CODE for TEC, and monthly medians for foF2. To test the model, data from 13 stations are used in such a way to get results for multiple latitude zones (middle, low, equatorial). The results were obtained for quiet and disturbed conditions by comparison with the results of the IRI model, taking into account the STORM-factor. Formula (14) of their paper shows that the comparison is carried out not with respect to the observational values of foF2, but to this model. The assumptions made in constructing the model are as follows: (1) the linear dependence of the parameters of the TEC, foF2 and τ on the level of solar activity, (2) the lack of longitudinal dependence of these parameters at the same LT, (3) transition from a geographic to a geomagnetic coordinates does not affect the description of variations in the parameters of the ionosphere from the LT, (4) the constancy of τ in quiet and disturbed conditions. The results were obtained for the five magnetic storms of varying intensity in the period 2000–2014. They are described in detail for several stations during individual disturbances with the general conclusion that the new model provides improved compliance compared with the model IRI-STORM in middle and low latitudes and in equatorial latitudes worsens results in the quiet and in disturbed conditions. However, as shown in Table 4 of the paper [44], the deterioration takes place in quiet conditions for the midlatitude station Chilton, and the low-latitude station Ebre. Deterioration in quiet conditions is a surprising fact, since in this case, such a model should give better results than the IRI model. As is known, the model values of TEC(IRI) are very different from the observational ones. Since the model of the authors uses the observational values of TEC, it should always lead to improvement. Consider how the behavior of τ corresponds to the assumptions of the model. The behavior of τ depending on the level of solar activity can be obtained from [39], which shows the daily variations of τ(med) and τ(IRI) for July and December in different years in the range of 2002–2010 for the stations from different latitudinal zones of the globe from auroral to equatorial (Juliusruh, Goosebay, Thule, Grahamstown, Ascension Island). This behavior includes both nonlinear changes and the constancy of the values in the daytime. Dependence of foF2 and TEC on the level of solar activity does not play a significant role since the quotient is taken. The dependences of τ(med) from RZ12 for an etalon station Juliusruh are shown in Figure 18.

1	2	3	4	5	6	7	8	9	10
Station	b1, b0		IRI	rec	stat	reg2	reg4	Lat1	Lat2
Juliusruh	3295.5	full	0.73	0.41	0.43	0.68	0.67	1.03	0.57
reg2	273.2	dist	1.44	0.52	0.49	0.67	0.68	1.07	0.64
Athens	5929.3	full	0.91	0.36	0.46	0.56	0.48	0.52	0.58
reg2	253.2	dist	1.31	0.44	0.74	0.52	0.59	0.87	0.51
Grahams	3788.7	full	0.80	0.40	0.54	0.77	0.59	0.73	0.73
Lat2	293.2	dist	1.54	0.46	0.62	0.84	0.75	0.82	0.77
Longyear	4947.1	full	0.70	0.43	0.62	0.60	0.58	0.82	0.58
Lat2	244.2	dist	0.69	0.49	0.73	0.69	0.69	1.01	0.63
Thule	692.7	full	0.51	0.14	0.15	0.56	0.42	0.47	0.59
	437.6	dist	0.55	0.10	0.13	0.51	0.46	0.64	0.54
Millstone	4864.4	full	0.90	0.50	0.47	0.48	0.46	0.67	0.49
Lat1	265.4	dist	1.38	0.67	0.67	0.65	0.81	0.81	0.80
Bejing	5402.8	full	1.17	0.49	0.61	0.61	0.58	0.70	0.62
reg4	263.9	dist	1.99	0.42	0.64	0.45	0.51	0.84	0.45
Kokubunji	6176.7	full	1.29	0.47	0.65	0.61	0.69	0.85	0.62
reg4	228.4	dist	2.11	0.55	0.66	0.56	0.70	0.96	0.56
Niue	4874.7	full	1.85	1.15	1.36	1.35	1.28	1.43	1.29
reg4	285.0	dist	1.67	0.71	1.00	0.73	0.85	1.11	0.67
Cocos	5467.3	full	1.43	0.55	0.68	0.86	0.62	0.65	0.82
Lat2	267.8	dist	1.66	0.52	0.77	0.88	0.67	0.80	0.83
Mawson	1466.2	full	0.91	0.27	0.37	1.00	0.85	1.02	0.92
Lat2	386.8	dist	1.12	0.12	0.21	0.80	0.98	0.98	0.81

Table 4.

Deviations of frequencies, calculated by hyperbolic dependence, from the experimental values of March 2015.

For July, trend is visible in a linear relationship, but for transition from year to year, it cannot be. For December of moderate and low activity, there is a constancy of τ(obs) during daylight hours; in other periods, linearity is violated. With regard to the assumption 2, the authors themselves point out that the presence of longitudinal dependence may be the cause of the deterioration. Further illustration is shown in Figure 19 for stations in various zones during March 2015, which had the largest number of stations and which also contains moderate disturbance (min Dst = −223 nT). Figures are given for τ(med) and τ(IRI) in: (a) middle latitude zone, (b) lower latitudes, (c) equatorial areas. Latitudes of stations are very close. A couple Juliusruh–Novosibirsk belongs to the middle latitudes, couples Nicosia–Kokubunji and Perth–Grahamstown, respectively, to the low latitudes of the northern and southern hemispheres. A couple Ramey–Sanya lies in the area between the low and equatorial latitudes. A couple Cocos–Darwin is closer to the equatorial zone. A couple Sao Luis–Fortaleza is in the equatorial zone. Reference [44] does not apply to high-latitude and auroral zones and, however, as in [26] the possibility of using the IRI model in these areas was shown, the results for a couple Tromso–Amderma are given.

Figure 18.

Illustration of τ(obs) dependence on the level of solar activity on the example of the station Juliusruh.

Figure 19.

Effect of longitude dependence on the behavior of τ at the same LT.

We see a good agreement between the values of τ(IRI), however, large differences between τ(med). It is necessary to emphasize the differences between τ(IRI) and τ(med), which are precisely define the differences of ΔfoF2 using τ(IRI) and τ(med), reviewed in [26]. With regard to the assumption 3, if the transition is not affected, why is implemented it. Assumption 4 implies the use of the average value of τ. It goes without saying, but since the authors introduced the item, it should be noted that it is the difference between τ in quiet and disturbed conditions, especially differences from τ(IRI), are the main cause of discrepancies between the calculated and experimental values of foF2. Figure 20 shows a comparison of τ(obs) during the disturbances with a median τ(med) and the value of the model τ(IRI) for two moderate disturbances in July 2004 with a minimum Dst = −197 nT and in December 2006 with a minimum Dst = −147 nT.

Figure 20.

Illustration of differences of τ(obs) from τ(med) and τ(IRI) during the disturbances. Respective days are shown on the title of drawings.

These figures illustrate not only the difference between τ(IRI) and τ(med), but still big differences of τ(obs) from τ(med) and τ(IRI) during the disturbances. That is why, the use of τ(med) during the disturbances gives smaller deviation of foF2 than τ(IRI), but larger than the deviation in quiet conditions.

This paper also attempts to develop a global model of τ(med). In principle, there are several options: (1) the construction of a superficial function such as kriging of the values τ(med) at several points, (2) two-parameter model based on hyperbolic approximation τ(hyp) = b0 + b1/NmF2, (3) the NGM model which can be constructed on the basis of two empirical models for TEC [22] and NmF2 [23], (4) the IRI-Plas model [7, 11]. Regarding the first option in [35] were expressed some doubts. This section describes the model of the second option. Results of testing models of the third and fourth options were presented in [27] and Section 3.

Since the construction of the model using the values themselves is not possible because of the large variability of values (in particular, the pre-sunrise peak at some latitudes), we attempted to use a hyperbolic dependence on an example of March 2015 when there was the largest number of stations. For hyperbolic dependence, coefficients b0 and b1 from τ(med) = b0 + b1/NmF2 were modeled. The results are given for some of the most wide regions. The region 2 contains 8 stations of the European continent, the region 4 contains 9 stations of Far East area. Curves for a zone of latitudes Lat1 from −52° to +65° are constructed according to 13 stations, basically, of the American continent of northern and southern hemispheres. The area for a zone of latitudes Lat2 from −68° to +78° includes 20 stations of the European, Siberian, and Southeast regions. Behavior of coefficients b0 and b1 for these regions is shown in Figure 21.

Figure 21.

The behavior of the coefficients of a hyperbolic approximation for various regions.

Figure 22.

The behavior of the coefficients b0 and b1 of hyperbolic approximation for the Juliusruh station for April of several years.

The calculations use average values. They make up 250.62 km and 4757.36 m⁻² for region 2, 280.21 km and 4386.01 m⁻² for region 4, 282.92 km and 5581.81 m⁻² for zone Lat1, 257.63 km and 4276.64 m⁻² for zone Lat2. The results are shown in Table 4. This table includes the following data. Column 1 indicates the station name and the region which it belongs to. The second column shows the coefficients of the hyperbolic dependence of τ(obs) for the corresponding stations. The third column specifies the conditions which include two series of values. The top line shows average of all days of the month, at the bottom—the average for disturbed days (from 16 to 21 March). The fourth column shows the results for the initial IRI model, the fifth column—the absolute difference between the experimental values of foF2(obs) and the values calculated using τ(med) and TEC(obs). Column 6 contains the deviation of frequencies calculated using the coefficients b0 and b1 of hyperbolic approximation for a given station. Other columns give results using the coefficients of the regions indicated in the column heading. All of these values should be compared with the values for the IRI model selected in bold.

It is visible that all values are higher in disturbed days and distinctions are the greatest for initial IRI model. There is a certain possibility to use coefficients of one region for calculation of foF2 in another area. It testifies about a global character of τ(med) models. One of the important problems consists in dependence of coefficients on the level of solar activity. Figure 22 shows coefficients b0 and b1 for the various years arranged in decreasing order of solar activity.

Another method of constructing a global model of τ(med) would be to use the coefficients K(τ) = τ(obs)/τ(IRI). Definite advantage of this model may be the fact that in its denominator stays the value of τ(IRI), having a global nature, and a small change in K(τ) in regions with similar longitude.

5. Conclusion

The appearance of models of the total electron content of the ionosphere TEC shows the progress made in the modeling of this parameter. This allows us to compare and use these models to forecast of TEC for any level of solar activity and to estimate the positioning accuracy. The new result is their comparison. It is shown that the majority of them provide an adequate accuracy and reliability. However, it should be noted the impact of uncertainties of their determination. These inaccuracies can be compensated using relative values, but often absolute values are needed. For four global maps JPL, CODE, UPC, ESA, solution is to construct a weighted average IGS according to four maps [45] which is also available on the same site together with the values of the maps. The main application of the TEC, discussed in this chapter, is determination of NmF2 and the critical frequency foF2 of the ionosphere. Global and continuing measurement of TEC using navigation satellites allows us to pose the problem of determining foF2 in the global scale. To do this, we need to know the proportionality factor between the TEC and NmF2, that is, the equivalent slab thickness τ of the ionosphere. It is shown that the existing models of this parameter are not global and do not provide sufficient accuracy in determining foF2. It is proposed to use the median τ(med) of the experimental values of this parameter and an approach to build its global model is presented. The advantages of using τ(med) are: (1) obtaining instantaneous values of foF2, which are especially important for the disturbed conditions, (2) calibration of TEC values for any global map or any set of experimental TEC that mitigates the impact of the uncertainty of these values.

Acknowledgments

Authors thank Southern Federal University for support (the Grant No. 213.01-11/2014-22).