© M. Limberger et al., Published by EDP Sciences 2015

1. Introduction

The application of electron density profiles (EDPs) derived from radio occultation (RO) measurements for ionospheric electron density (N_e) modeling is widely accepted, helping to gain a global data coverage and increased sensitivity for the vertical structure of the N_e distribution. The ionospheric F2 layer peak, defined by the maximum electron density NmF2 and its corresponding peak altitude hmF2, can be directly accessed from the profiles’ shape. The first successful proof-of-concept on limb sounding the Earth atmosphere was carried out in 1993 at the University Corporation for Atmospheric Research (UCAR) by means of the Global Positioning System/Meteorology (GPS/MET) experiment with the Microlab-1 satellite (Ware et al. 1996). The studies with GPS/MET were mainly focused on numerical weather predictions but besides, the RO technique has been exploited for profiling the ionosphere and deriving EDPs (Hajj & Romans 1998). The number of low earth orbit (LEO) satellites applicable for ionospheric sounding continuously rose in the last decade and radio occultations have become increasingly important. Among them, the Formosat-3/COSMIC (F-3/C) constellation consisting of six satellites located on six individual orbital planes. F-3/C was launched on 15 April 2006 and contributes significantly towards a global coverage of ionospheric sounding data (Tsai et al. 2009).

One popular and probably the most widespread algorithm to derive EDPs from ionospheric RO measurements is the classical Abel transform inversion that can be driven by single-frequency bending angle data or directly by dual-frequency GNSS carrier-phase measurements. Hajj & Romans (1998) considered, for instance, the GPS single-frequency signal bending in the ionosphere and used the Abel inversion to derive EDPs to perform a validation between GPS/MET profiles with ionosonde measurements during low solar cycle conditions with a bending angle <0.01°. It was concluded that the accuracy of NmF2 yields around 20% of its absolute level, which was later confirmed by Tsai et al. (2001). Furthermore, to estimate the absolute precision of the classical Abel inversion, Schreiner et al. (2007) studied the accuracy of N_e by co-located F-3/C EDPs and obtained a root-mean-square (RMS) difference of approximately 10⁹ el/m³, a value which indicates already, that rough discrepancies cannot be ruled out by using the classical approach.

The use of GPS L1 frequency measurements for bending angle calculations is one possible method for the retrieval of EDPs but requires precise orbit and clock information. The processing of GNSS dual-frequency measurements on the contrary enables the elimination of geometry-dependent terms through forming the ionosphere (or geometry-free) linear combination L₄, respectively. Various facilities provide EDPs routinely as preprocessed products mostly by exploitation of the carrier-phase excess. F-3/C EDPs are, for instance, processed and officially released by the COSMIC Data Analysis and Archival Center (CDAAC, US)¹ and the Taiwan Analysis Center for COSMIC (TACC, Taiwan)² as level 2 products. After registration, the data is freely available via http servers as daily tarballs. For the automatic real-data retrieval, CDAAC considers the classical Abel transform, most probably due to its simplicity and relatively low computational burden.

It is well known that the spherical symmetry assumption is one of the most significant error sources of the classical Abel inversion transform that degrades the N_e accuracy in particular in the presence of the equatorial ionization anomaly (EIA). This deficiency led to the situation that various facilities developed their own modified retrieval algorithms with different concepts to account for horizontal electron density variations in the data processing. Research facilities, such as the Universitat Politécnica de Catalunya (UPC, Spain) or the Center for Space and Remote Sensing Research (CSRSR, Taiwan), perform an independent retrieval based on improved Abel inversion techniques under consideration of horizontal ionospheric gradients. The CSRSR approach is based on the interpolation of information from nearby EDPs with a cubic spline fit (Tsai & Tsai 2004) while the UPC method considers the Vertical Total Electron Content (VTEC) derived from terrestrial GNSS observations (Hernández-Pajares et al. 2000). The latter approach is capable of improving the retrieval algorithm with relatively low implementation effort under support of global ionosphere maps (GIMs) as shown, for instance, by Hernández-Pajares et al. (2000), Garcia-Fernandez et al. (2003), and Aragon-Angel (2010). Therefore, the UPC approach has been applied in this work involving GIMs (Hernández-Pajares et al. 2009) provided by the International GNSS Service (IGS; Dow et al. 2009) in the Ionosphere Map Exchange (IONEX) format. Few potential drawbacks, however, remain. (1) Since the processing of IGS GIMs is driven by terrestrial GNSS data, it is a matter of fact that the improvement is limited over the oceans. (2) A single layer mapping is taken into account by most of the IGS data analysis centers, i.e., the accuracy of the GIMs partly depends on the quality of the associated mapping function. The accuracy improvement for the F2 layer peak achieved through the inclusion of horizontal VTEC gradients has for instance been studied by Yue et al. (2012) who obtained an enhancement of around 7% for the determination of NmF2 and a slightly weaker improvement for hmF2. From the comparison of EDPs against co-located Jicamarca ionosonde measurements during 2007, Aragon-Angel et al. (2011) obtained an error improvement of 10% by means of the UPC retrieval method.

Besides the evaluation of the GIM-aided improvement of the EDP retrieval, several authors published works focusing on the accuracy assessment of the F2 layer peak by using diverse evaluation methods with study periods of different lengths. For instance, Wu et al. (2009) compared NmF2 against ionosonde observations during the period of 13 July 2006 and 31 December 2007 and Lei et al. (2007) evaluated correlations of maximum electron densities between F-3/C EDPs and ionosondes at the early stage of the constellation in July 2006. Krankowski et al. (2011) validated F2 peak parameters of F-3/C measurements in the European sector during 2008 and Yue et al. (2010) studied the retrieval accuracy by means of synthetic data driven by background information of the NeQuick model to determine absolute and relative accuracies for NmF2 and hmF2.

In contrast to the aforementioned publications, the study presented here deals with the global evaluation of NmF2 and hmF2 by exploitation of the long-term availability of F-3/C occultation measurements for a period of more than half a solar cycle. As a key innovation, the paper focuses on the monitoring and global accuracy assessment of the F2 layer peak against ionosonde measurements between 10 June 2006 and 29 April 2014, i.e., during a time frame of almost 8 years. The study interval thus includes a solar minimum and maximum period. Such a validation is based on an immense amount of data where the manual selection of reliable ionosonde data is no longer possible. Therefore, electron density peaks observed by the global SPIDR ionosonde network are considered as a reference. The SPIDR database provides automatically scaled measurements to be used as one of only few data sources for global long-term assessments. The evaluation presented here concentrates on the long-term aspect and less on the impact of diurnal and small-scale variations. For this purpose, a resolution of 30 days and consequently around one day per month has been chosen to account for the annual and seasonal trends of the ionosphere. As required for the electron density retrieval following the UPC strategy (Garcia-Fernandez et al. 2003; Aragon-Angel et al. 2011), horizontal VTEC gradients provided by the ionospheric IGS analysis facility CODE (Center of Orbit Determination in Europe) are considered in this study.

Starting with Section 2, the separability algorithms according to Hernández-Pajares et al. (2000), Garcia-Fernandez et al. (2003), and Aragon-Angel (2010) are briefly recapitulated as they constitute the fundamental methodology including the observation equation for the establishment of the EDP database containing the F2 layer peaks. The process chain as realized in this study is described in Section 3 together with the relevant data screening steps for the database generation. In Section 4, the global trends of NmF2 and hmF2 are investigated and compared with the solar activity. The assessment of ionospheric F2 layer peaks from the comparison between F-3/C and ionosonde observations for evaluating the accuracy level is presented in Section 5 in terms of a global evaluation and in Section 6 by distinction of local times and different magnetic latitude regions. Section 7 closes the paper with a summary and final remarks.

2. Electron density retrieval

The propagation of microwave signals is affected by the ionosphere in terms of a group delay and phase advance (or excess of phase). The ionosphere is a dispersive medium and the higher the signal frequency, the smaller the delay impact. Changes in the phase data are proportional to the Slant Total Electron Content (STEC) and can be exploited in dual-frequency measurements to extract information about the current N_e distribution. In terms of occultations, i.e., based on the condition that the signal transmitted by the GPS satellite reaches the LEO satellite from below the local horizon, N_e profiles below the LEO orbit can be derived. An occultation pass is defined by a sequence of consecutive measurements in occultation geometry referring to the epochs i ∈ {1, …, M}. There are rising and setting occultations but due to the retrieval concept that will be described, all occultation passes in this work are treated as setting occultations starting from the highest to the lowest occultation. The so-called tangent point (TP) is defined as the closest point to Earth along the line-of-sight (LOS) between GPS and LEO satellite at a specific occultation epoch and defines the anchor point for the retrieved N_e. Consequently, each EDP is composed of a TP arc with associated N_e values at decreasing altitudes and also changing horizontal coordinates where each N_e observation corresponds to an individual occultation epoch. It is important to have in mind, that the whole occultation event can last up to several minutes and the distance between the first and last TP can differ in the magnitude of 20°–30° (≈ 2000–3000 km) in latitude (φ) and longitude (λ).

The UPC retrieval approach is of particular importance for the following validation studies and the main steps for deriving the fundamental observation equation will therefore be recapitulated in accordance with Hernández-Pajares et al. (2000). The retrieval technique is based on the Abel transform inversion where each occultation is processed independently. From each dual-frequency occultation measurement, STEC can be derived as(1)following the definition of the N _e integral along the signal path s from GPS to LEO satellite. Under the assumption of spherical symmetry as implied by the classical Abel retrieval, the integration can alternatively be performed between TP and LEO. The parameter STEC^TOP has been introduced in Eq. (1) to take the ionospheric contribution above the orbit into account. STEC can be obtained from phase geometry-free linear combinations L₄ as(2)where ϕ₁ and ϕ₂ are the dual-frequency carrier-phase measurements in cycles, λ₁ and λ₂ are the wavelengths of both carrier frequencies, b₄ is the carrier-phase bias on L₄, and α is a frequency-dependent factor. ε denotes the unmodeled random measurement noise.

To solve Eq. (1) by means of integral discretization, an onion shell structure can be established with the heights of successive TPs as mean radii. This implies of course that the shell thicknesses depend on the given measurement sampling. The full set of layers j ∈ {1, …, N} is consequently defined after collecting all measurements of a specific occultation pass. l_ii is the distance between the TP_i at epoch i and the adjacent layer. l_ij are the pass lengths of the GPS-LEO LOS at epoch i in shell j. A solution of Eq. (1) at epoch i can then be determined from(3)

Equation (3) is valid under the assumption of spherical symmetry, i.e., a constant density N_e in each shell, parameterized only by the radius r. As a consequence, no horizontal gradients are considered in N_e. To account for this drawback, the separability hypothesis is introduced with(4)

It describes the N_e distribution as a function of VTEC scaled by a shape function parameter SF(h). The unknown parameters to be estimated are SF values located in the TP positions. VTEC is assumed to be given and can, for instance, be extracted from IGS GIMs. The consideration of Eq. (4) in Eq. (3) leads to the observation equation(5)where differences to a pass-specific reference observation STEC_R (or L_4,R are taken into account. For simplicity, it is assumed that the b₄ term only contains the phase ambiguities whereas effects of secondary importance such as phase windup or antenna phase center variations, have been neglected and b₄ consequently cancels out. In the present work, the reference observation must satisfy the condition of being observed within the elevation interval [0°, −5°]. For the reason that the F-3/C satellites orbit in altitudes of approximately 800 km, the difference between STEC^TOP and STEC_R is expected to be negligibly small. This assumption should be considered with caution since it may be valid for F-3/C but typically not for LEO satellites in lower orbits. Regarding Eq. (5), it shall further be noticed that SF_i and SF_j indicate the affiliation of the shape function to the corresponding TP_i and TP_j, respectively. The VTEC values are specified by different superscripts (TP, F3C, and GPS) denoting the position to which the VTEC refers to. For instance, VTEC^TP refers to the VTEC observed at the TP position while VTEC^F3C and VTEC^GPS are the VTEC values at the layer intersection points in LEO and GPS direction. This concept is depicted in Figure 1 including the separation of the LOS into a TP-LEO and TP-GPS directed part.

Fig. 1. Improved Abel inversion under consideration of the separability hypothesis with integral discretization by means of the onion shell structure. Colored arrows illustrate the discretized integration steps between adjacent layers. This draft shows the composition with measurements at two epochs i1 and i2 yielding the associated tangent points TP1 and TP2 that are located at the mean radii of the layers j1 and j2.

After collecting all data of an occultation pass, Eq. (5) can be solved recursively starting from the outermost toward the innermost shell as a triangular linear equations system x = A⁻¹ y that is established for each occultation pass. x contains the set of shape function parameters, y is the vector of STEC observations, and A is the triangular coefficient matrix(6)where rows represent the observations i ∈ {1, …, M} with decreasing TP altitudes and the columns contain partial derivatives of the shape function parameters starting from the outer toward the innermost shell j ∈ {1, …, N}.

3. Database of electron density profiles

The paper structure orients on the flowchart in Figure 2. Processing steps are distinguished by different colors: Green boxes are representing input data, yellow boxes indicate the retrieval computation steps, and blue boxes highlight the peak evaluation. To perform the retrieval, various data resources are incorporated. According to the topmost segment of Figure 2, observation data in the Receiver Independent Exchange (RINEX) Format and F-3/C orbits in the Standard Product 3 Orbit Format (SP3) are downloaded from the TACC or CDAAC data servers. Both data types are labeled with unique identifiers, podObs and leoOrb, as declared by CDAAC. For the podObs data, we used the re-processed products which are available since 22 April 2006, i.e., one week after the launch of F-3/C, until 20 April 2014 covering almost 8 years of measurements. GPS orbits³ as well as 2 h sampled GIMs⁴ are acquired from the network servers of the IGS Analysis Center CODE.

Fig. 2. Flowchart for the validation procedure providing an overview of data sources (green boxes), the retrieval (yellow boxes), and assessment steps (blue boxes).

This study has been performed between 10 June 2006 and 29 April 2014 starting around 2 months after the launch of F-3/C due to missing data and orbits in the initial period of the constellation. A homogeneous sampling of 30 days has been chosen resulting in totally 97 days under investigation. For each day, the data sets depicted in Figure 2 are downloaded to perform the electron density retrieval as described in Section 2. Among others, the computation of geometry-free linear combinations from occultation measurements is performed in a pre-processing step. Occultation passes are composed of sequences of uninterrupted measurements between a specific GPS transmitter and LEO receiver in occultation geometry where various conditions have to be fulfilled. For instance, maximum thresholds for time gaps between adjacent measurements and for the STEC acceleration have been considered. Furthermore, a minimum number of observations per pass is mandatory. After performing the pre-processing, each retrieved profile is screened in a post-processing step including the

1. identification of hmF2 within a defined altitude interval,

2. verification of the profile extending between given bottom and top height limits, and

3. exclusion of profiles with large electron density jumps between adjacent observations.

For each day included in this study, a database consisting of NmF2 and hmF2 observations derived from RO is finally established and checked for rough outliers by testing against a 3σ threshold. Those values, that exceed the interval of ±3σ around the global average, are iteratively excluded. The daily numbers of profiles that passed the pre-processing and screening procedure are depicted in Figure 3.

Fig. 3. Number of retrieved F2 layer peaks at each investigated day during the study period (10 June 2006–29 April 2014).

In total, around 85,000 suitable EDPs have been stored for the validation. The minimum amount of 97 profiles has been obtained on the first investigation day, i.e., on 10 June 2006, where only FM-1 (or COSMIC-1) observations and orbit information between 00:00 UT and 05:59 UT were available. The maximum is reached on 5 June 2007 with 2301 profiles observed by all six F-3/C satellites. The distribution of the retrieved electron density peaks for both days, 10 June 2006 on the left and 5 June 2007 on the right, is shown in Figure 4. Different colors indicate the different F-3/C satellites (FM-1 = red, FM-2 = green, FM-3 = yellow, FM-4 = magenta, FM-5 = cyan, FM-6 = white).

Fig. 4. Daily distribution of the electron density peaks derived from EDPs of FM-1 (left, 109 peaks) at 10 June 2006 vs. all F-3/C satellites FM-1 to FM-6 (right, 2447 peaks) at 5 June 2007.

4. Global trend of the F2 layer peak

The global trends of NmF2 and hmF2 have been evaluated in a first step. Therefore, all retrieved electron density profiles that passed the pre- and post screening steps as described in Section 3 are taken into account for the derivation of F2 layer electron density peaks and computation of global, daily averages. The time series of NmF2 and hmF2 are provided by Figure 5. As can be expected, an increase of the daily averages from low toward high solar activity becomes visible for both parameters, indicated by the dark blue line. Light blue bands represent the corresponding daily standard deviations, σ_NmF2 as well as σ_hmF2, and characterize the variability of the parameters in the global domain. It is obvious, that σ_NmF2 increases with the solar activity while the magnitude of σ_hmF2 remains at a rather constant level and almost independent of the solar conditions. The effect may be interpreted as an increase of the global variations for NmF2 in correlation with the solar conditions while hmF2 likewise increases but with a consistent fluctuation magnitude.

Fig. 5. Trend of the global NmF2 (left) and hmF2 (right) as derived from global averages and standard deviations of the retrieved F2 layer peak computed for the days under investigation. The data is plotted for the study period between 10 June 2006 and 29 April 2014.

Taking into account that the global mean may be a good representative for the solar conditions, comparisons with the daily sunspot number R and the solar radio flux parameter F10.7 as solar physical quantities have been performed. Both quantities are commonly accepted as indices for reflecting the solar condition, typically on a day-to-day basis, where R is determined from optical surveying and F10.7 is a measure of the solar radio flux per unit frequency at a wavelength of 10.7 cm. The data of R and F10.7 are downloaded via the NASA GSFC/SPDF OMNIWeb⁵ interface (King & Papitashvili 2005). With e being a N × 1 unit vector and v_F10.7 as well as v_F2 being N × 1 data vectors containing daily values of F10.7 and global averages of the F2 peak parameters observed by F-3/C, then the time series of NmF2 and hmF2 can be scaled by(7)to the units of F10.7 to enable a graphical comparison of the relative trends. This scaling procedure was accordingly carried out to scale NmF2 and hmF2 with respect to R.

The resulting time series are shown in Figure 6 where red squares represent the scaled averages of NmF2, green squares are the hmF2 samples, and the blue squares identify sunspot numbers (top) and F10.7 indices (bottom). According to the color of the y-axis labeling, its units thus refer to the blue squares. Figure 6 indicates that the retrieved averages, both for NmF2 and hmF2, characterize the solar conditions very well with a slightly better performance of NmF2 as can be expected due to its higher sensitivity with respect to solar variations. The individual correlation coefficients are provided in Table 1. A high correlation coefficient of 0.92 between NmF2 and hmF2 indicates the proportional increase of both parameters toward higher solar activity. R and F10.7 as established indicators for reflecting the solar conditions are as expected highly correlated with an index of 0.94. Furthermore, correlations between NmF2 (hmF2) and R as well as F10.7 are computed with 0.83 (0.87) and 0.89 (0.88), respectively. The equivalent calculations have been performed to obtain also correlations with F2 parameters observed by ionosondes of the SPIDR network (not depicted here). In this case, the correlations yield 0.83 (NmF2) and 0.80 (hmF2) as compared to F10.7 as well as 0.81 (NmF2) and 0.80 (hmF2) with respect to R. In this scenario, the global averages of electron density peaks derived from occultations seem to represent the purely physical quantities R and F10.7 slightly better than the averages of those peaks observed by ionosondes. It can be assumed, that the better representation follows mainly from the improved global coverage of the occultation data. Global daily averages of RO-based F2 layer peaks can thus be identified as very suitable quantities to characterize the solar conditions.

Fig. 6. Comparison of scaled daily NmF2 (red squares) and hmF2 (green squares) averages with the daily sunspot number R (top, blue squares) and daily solar flux F10.7 index (bottom, blue squares) during the study period (10 June 2006–29 April 2014). The unit of the y-axis units refers to the sunspot number and F10.7 and the NmF2 as well as hmF2 samples have been scaled.

5. Assessment of the global F2 layer peak

To assess the accuracy of the F2 layer peak retrieved from RO measurements, ionosonde data tracked by stations of the SPIDR network have been taken into account. The data is available via the National Geophysical Data Center (NGDC) website⁶ as a public service by the US Department of Commerce, National Oceanic and Atmospheric Administration (NOAA), National Environmental Satellite, Data and Information Service (NESDIS).

Initial investigations showed, that the ionosonde data set contained anomalous observations mostly located in critical regions with high ionospheric gradients, for instance while the EIA passes, or due to irregularities such as scintillations, sporadic E-layer or spread-F. To eliminate these outliers, additional ionosonde specific data screening steps have been implemented. Thresholds related to minimum and maximum peak heights as previously considered for the profile retrieval should of course also hold for the ionosonde measurements. Additional ionosonde-specific preprocessing algorithms reject single measurements that appear isolated with time gaps of more than 30 min difference to the previous and next observation and analyze the parameter variability for the rejection of observations that exceed a relative change of 30% (NmF2) and 20% (hmF2) within a time interval of 30 min.

After cleaning the peaks observed by the ionosondes, the vicinity of each EDP peak within ±8° and 20 min in space and time is scanned for co-located ionosonde measurements. The occultation peak together with the closest ionosonde observation is then considered as a value pair introduced into the assessment procedure. In total, 7780 (NmF2) and 5902 (hmF2) samples in co-location with ionosonde observations have been detected.

The magnitudes of NmF2 and hmF2 as determined from the EDPs are provided by daily RMS values in Figure 7. Two prominent peaks are contained in the NmF2 time series that may be caused by a large amount of profiles located in ionospheric active regions for the specific days. It is clearly visible how the trend is driven by the solar cycle. Total averages with variations of 4.34 ± 2.19 × 10¹¹ el/m³ and 270.27 ± 22.18 km are obtained for NmF2 and hmF2, respectively. Minimum RMS values during the solar minimum phase (≈ between 2006 and 2010) are determined with around 2 × 10¹¹ el/m³ and 240 km. The maximum RMS values are analogously observed during the current solar maximum with up to 8 × 10¹¹ el/m³ for NmF2 under negligence of the peaks and 330 km for hmF2.

Fig. 7. Daily RMS of the retrieved peak density NmF2 (left) and corresponding height hmF2 (right) derived from radio occultations between 10 June 2006 and 29 April 2014.

Relative differences to the nearby ionosonde measurements have been determined as(8)and accordingly for ΔhmF2_% where subscripts IS and F3C indicate the affiliation of the measurement to the ionosonde or RO and the subscript % identifies the relative differences given in percentage units. The results are depicted in Figure 8 for ΔNmF2_% (left) and ΔhmF2_% (right). Total averages with standard deviations of avg(ΔNmF2_%) = −3.13 ± 5.67% and avg(ΔhmF2_%) = 0.87 ± 3.07% indicate no significance for systematic biases. Few abnormal peaks are contained in the series of ΔhmF2_% reaching 5%–14% which can be assigned to the middle of the years 2006–2009, i.e., during low activity. The largest differences for these days are originating from significantly lower ionosonde hmF2 observations and clearly higher EDP peaks. Further investigation is required to locate the source of these discrepancies.

Fig. 8. Daily averaged relative differences between global NmF2 (left) and hmF2 (right) observations as measured by ionosondes and F-3/C radio occultations between 10 June 2006 and 29 April 2014.

The trend of both series ΔNmF2_% and ΔhmF2_% features a slight increase correlated with the solar conditions over the years. The variations of ΔNmF2_% and ΔhmF2_% oscillate roughly in the range of 10% and 5%, respectively, under negligence of outliers. Compared to ΔhmF2_%, an increased variability is obtained for ΔNmF2_%.

In summary, the global evaluation showed a high agreement of the NmF2 and hmF2 trend with the sunspot number R and the F10.7 index. A rough approximate of the overall global accuracy of the F2 layer peak derived from RO measurements can be estimated from the global RMS together with the uncertainty computed from the relative differences to ionosondes. For NmF2, a standard deviation of around 5% together with a RMS varying with the solar activity between 2 × 10¹¹ el/m³ and 8 × 10¹¹ el/m³ yields an uncertainty of the retrieved peak density in the range of 1 × 10¹⁰ el/m³ and 4 × 10¹⁰ el/m³. Accordingly, the uncertainty of hmF2 with 3% and RMS between 240 km and 320 km leads to an approximation of the absolute difference in the range of 7 km and 10 km.

6. Peak evaluation under consideration of magnetic latitude and local time

With the outcomes of the Sections 4 and 5, initial estimates for the overall accuracy of NmF2 and hmF2 have been determined. In order to assess the accuracies under exclusion of strongly varying ionospheric conditions, different magnetic latitude sectors and local time intervals have been selected to account for daily ionospheric variations as well as the orientation of the structure (in particular the EIA) along the magnetic equator.

The time interval around the daily peak at 14:00 LT ± 2 h has been selected to study the quality of the retrieval around the ionospheric day maximum and shifted by 12 h at 02:00 LT ± 2 h to consider quiet conditions at night. Additionally, a timespan of moderate conditions at dawn between 8:00 LT ± 2 h was investigated. Different magnetic latitude sectors Δφ have been defined as Δφ_EQU = ±[0°, 20°], Δφ_MID = ±]20°, 60°], and Δφ_POL = ±]60°, 90°] to account for equatorial, mid-latitude, and polar regions. Due to the rather symmetric formation of the EIA along both sides of the magnetic equator, the latitude bands are defined with equal boundaries for the northern and southern hemisphere. The magnetic equator has been computed from yearly dependent North Pole coordinates as published by NOAA⁷ with annual Pole coordinates derived from the 10th generation of the International Geomagnetic Reference Field (IGRF). The variation of the magnetic equator in geographic coordinates as a consequence of the annual magnetic North Pole movement is depicted in Figure 9. The magnetic North Pole converged slowly from 2006 to 2014 toward the geographic Pole as can be seen from the polar plot in the top-left corner. The Pole latitudes vary approximately between φ_NP = 83.5° in 2006 and φ_NP = 86° in 2014. Consequently, the innermost and outermost curves of the main plot belong to 2014 and 2006, respectively. The magnetic equator varies within 2.5° latitude and exhibits a difference of up to 6.5° in 2006 with respect to the geographic equator that is depicted as a gray, dashed line.

Fig. 9. Variation of the magnetic equator in geographic coordinates during 2006 (outermost blue curve) and 2014 (innermost green curve).

Absolute differences ΔNmF2 = NmF2_F3C − NmF2_IS and ΔhmF2 = hmF2_F3C − hmF2_IS of both F2 peak parameters have been computed for each time interval and magnetic latitude sector. Exemplarily, the results of the daytime comparison in the interval 14:00 LT ± 2 h are depicted in Figure 10 (NmF2) and Figure 11 (hmF2). Three plots on the right column include the data distribution with red circles for the peaks retrieved from F-3/C occultations and blue triangles as co-located ionosondes for the three magnetic latitude regions. It should be noticed that all electron density profiles within a specific latitude band are detected first and afterwards the search algorithm screens for co-located ionosondes. Thus, the ionosonde location must not mandatory lie within the same latitude band but may be located slightly outside of the bounds. The data has been collected and plotted for all days under investigation and therefore each ionosonde that contributes data at least once during the study time frame is contained in the plot. Furthermore, the magnetic latitude bands computed from the North Pole coordinates of the latest year 2014 are exemplary depicted. The parameter differences for each latitude interval are provided in the left columns of Figures 10 and 11. For each interval, the absolute averaged difference and standard deviation have been computed and included in the legends above the plots.

Fig. 10. Left column: Differences ΔhmF2 between ionosonde observations and F2 layer peaks retrieved from F-3/C occultations collected during the days under investigation between 10 June 2006 and 29 April 2014 within 12:00 and 16:00 LT and magnetic latitude sectors ΔφEQU = ±[0°, 20°] (top), ΔφMID = ±]20°,60°] (center) and ΔφPOL = ±]60°, 90°] (bottom). Right column: Retrieved NmF2 positions (red circles) and co-located ionosondes (blue triangles) within the defined magnetic latitude bands. The latitude intervals are exemplary depicted here under consideration of the magnetic North Pole from 2014.

Fig. 11. Left column: Differences ΔhmF2 between ionosonde observations and F2 layer peaks retrieved from F-3/C occultations. The data have been collected during the days under investigation between 10 June 2006 and 29 April 2014 within 12:00 and 16:00 LT and magnetic latitude sectors ΔφEQU = ±[0°, 20°] (top), ΔφMID = ±]20°,60°] (center) and ΔφPOL = ±]60°, 90°] (bottom). Right column: Retrieved hmF2 positions (red circles) and co-located ionosondes (blue triangles) within the defined magnetic latitude bands. The latitude intervals are exemplary depicted here under consideration of the magnetic North Pole from 2014.

The increased presence of ionosondes in mid-latitudes in populated and urbanized regions is clearly visible from the figures. Consequently, the number of value pairs varies strongly with a majority in the mid-latitudes. For ΔNmF2, averages and standard deviations of 1.22 ± 1.91 × 10¹¹ el/m³ in the equatorial region, 0.10 ± 1.43 × 10¹¹ el/m³ for the mid-latitude band and −0.27 ± 1.16 × 10¹¹ el/m³ for the polar region, are determined. Since the differences are obtained during the daily ionospheric maximum, the standard deviation increase from the polar toward the equatorial band (or quiet toward active regions, respectively) can be related to the signal strength.

Accordingly, the averages and standard deviations of ΔhmF2 are determined with 7.50 ± 30.01 km in the equatorial region, 5.43 ± 23.81 km for the mid-latitudes and 7.45 ± 29.91 km for the polar caps. Maximum differences of approximately 100 km for only few comparisons appear. All three plots related to ΔhmF2 exhibit a positively biased mean, i.e., a possible overdetermination of the F2 layer peak height by the occultations. However, the significance of a systematic offset cannot be concluded as the variations of ΔhmF2 in the range of 24 km–30 km are fairly high.

Furthermore, the correlations of the absolute measurements for NmF2 and hmF2 referring to the daytime interval during 14:00 ± 2 h (according to the Figures 10 and 11) are depicted in Figure 12 for all detected pairs of F2 layer peaks. Those parameters derived from F-3/C are given along the x-axis and the associated ionosonde values refer to the y-axis. NmF2 is depicted in the left column and hmF2 on the right. Similar to the previous figures, the top subplots represent the equatorial band, the center plots refer to the mid-latitudes, and the bottom plots show the results obtained for the polar areas. x- and y-axes of each plot are restricted to the same limits with edges connected by a black line indicating a 100% correlation.

Fig. 12. Correlations between F2 layer peaks derived from F-3/C (x-axis) and ionosonde observations (y-axis) for NmF2 (left column) and hmF2 (right column). According to Figures 10 and 11, the subplots from top to bottom refer to the equatorial, mid-latitude, and polar regions.

Very good agreements have been obtained for NmF2 with correlation coefficients of 0.94 and 0.92 for the equatorial belt and mid-latitudes, respectively. A slightly weaker correlation coefficient of 0.84 was found at the polar caps. In the case of hmF2, the correlations are generally weaker compared to NmF2 but nevertheless reach 0.76 along the Equator and 0.75 in mid-latitudes. A moderate correlation of only 0.56 is obtained for the polar area. Results for the time periods at dawn and night are not shown here, but included in the Tables 2 and 3, summarizing the outcomes of Section 6.

Absolute and relative averaged differences and standard deviations of NmF2 and hmF2 for all local time intervals and magnetic latitude sectors as well as correlation coefficients between ionosonde and EDP peaks are compared. Additionally, the number of samples (#) considered in the evaluation is provided.

For ΔNmF2, mean differences in the magnitude of −0.2 × 10¹¹ el/m³ to 1.2 × 10¹¹ el/m³ have been obtained together with standard deviations between 0.6 × 10¹¹ el/m³ and 1.9 × 10¹¹ el/m³. It can be noticed for the night interval, that the NmF2 measurements by ionosondes are consistently lower compared to the F-3/C observations. The strongest variation can be found at the equator for all three time intervals where relative differences of 15%, 2.89%, and −6.95% are found during daytime, dawn, and night. The uncertainties in the range of 22%–29% are determined during day and dawn and 30%–35% at night.

High NmF2 correlations are found for all times and latitudes with marginally reduced conformity at the Poles and generally during night. In the equatorial and mid-latitude bands at day and dawn the correlations even reach values larger than 0.9. The strong correlation between ionosondes with NmF2 supports the outcome of Lei et al. (2007) where 276 co-located value pairs between 31 globally distributed ionosondes and F-3/C measurement are analyzed and a correlation coefficient of 0.85 was found. Furthermore, a relative accuracy level of around 20% for NmF2 was determined by Hajj & Romans (1998), Tsai et al. (2001), and Wu et al. (2009). The slightly worse accuracy level resulting from this study is probably caused by the use of automatically detected ionogram peaks contained in the SPIDR database. Schreiner et al. (1999) compared GPS/MET NmF2 with ionosondes and determined variations at the magnitude of 26%, i.e., similar to the results obtained here. Schreiner et al. (1999) additionally found an insignificant bias of few percentages which was later confirmed by Wu et al. (2009) and various other authors.

In the case of hmF2, mean differences in the range of −7 km and 7.5 km with an outlier of −22 km at the Pole during night are determined. Standard deviations of the order of 24 km–30 km discard the existence of systematic biases. Again, an outlier of 48 km standard deviation is calculated at the Pole during night. Significant differences of the variations with respect to the analyzed time intervals or latitude bands cannot be observed. Furthermore, maximum relative averages between −6% and 5% with relative standard deviations in the range of 10%–16% are determined for ΔhmF2. It can be noticed, that the largest variations are determined within the Polar Sector. Moderate correlations of around 0.7 have been found at the Equator and mid-latitudes during daytime and dawn with weaker correlations of 0.5 at the Poles. The reliability is clearly reduced during night where correlations of 0.6 are obtained for the equatorial and mid-latitude sectors and 0.5 for the Poles.

7. Summary

The use of RO measurements and derived electron density profiles for describing the ionospheric electron density distribution helps to gain a global data distribution and to access information about the vertical N _e structure. The Formosat-3/COSMIC mission collects occultation measurements since 2006 and has been exploited for establishing a database of occultation data between 10 June 2006 and 29 April 2014 using a sampling of 30 days. This data set is compared to ionosonde measurements provided by the global SPIDR network to obtain information about the reliability of retrieved F2 peak parameters during changing solar conditions.

In a first step, the global trend of NmF2 and hmF2 has been expressed by daily averages and compared with sunspot numbers and solar flux F10.7 indices where high correlations between 0.83 and 0.89 have been found for both parameters. Afterwards, a database of globally distributed value pairs including EDP peaks and co-located ionosonde observations has been established to determine an overall accuracy level. For this purpose, RMS values of the retrieved peak together with relative differences between ionosonde and EDP peaks have been computed. As a preliminary assessment, the overall accuracy level was determined in the range of 1 × 10¹¹ el/m³ and 4 × 10¹¹ el/m³ for NmF2 and 7 km–10 km for hmF2.

A more detailed assessment was provided in Section 6 where the data was separated by different latitude sectors and local time intervals. The results are summarized in the Tables 2 and 3 for NmF2 and hmF2. High correlations between ionosonde and EDP data have been found for NmF2. A significant systematic bias was not observed. The corresponding relative variations are determined in the range 22%–35% depending on the latitude sector and time interval. Particularly high correlations above 0.9 are found at the Equator and mid-latitudes during day and dawn. A marginally weaker agreement is observed for measurements at the Poles and generally during night.

As for NmF2, the absence of systematic biases can be also concluded for hmF2. Standard variations are found in the range of 24 km and 30 km as well as 10% and 16% referring to absolute and relative units, respectively. Good correlations of 0.7 have been obtained during day and dawn for the equatorial and mid-latitude regions. However, the correlations at the Poles and during night are significantly reduced.

The accuracies obtained in this work show quiet a good agreement with other publications that concentrated on the assessment by only few selected ionosonde stations. Generally high correlations between the F2 layer electron density peaks as derived from F-3/C occultations and ionosonde measurements have been obtained in this study. However, an extensive pre-processing and data screening was necessary for both the RO and ionosonde data to eliminate disturbances in the data and improve the data reliability. Many authors like Krankowski et al. (2011), Hu et al. (2014), Wu et al. (2009), Yue et al. (2010), and Lei et al. (2007) evaluated the accuracy of electron density peaks, retrieved from F-3/C data, either on a regional scale by means of selected ionosonde stations allowing for a manual ionogram scaling or even used synthetic data to validate the accuracy under different conditions, e.g., under the influence of varying F10.7 indices. In this work, SPIDR network data based on automatically scaled measurements are taken into account. Such observations enable the global long-term assessment but hold the potential risk of a reduced data quality. The processing algorithms of most modern ionosondes are based on the ARTIST-5 software (Galkin et al. 2008; Reinisch et al. 2009) that provide quality flags and uncertainty boundaries. Such quality indicators were not contained in the SPIDR data used here but could be of high relevance for further studies to verify the reliability of ionosondes in such validations.

In this paper, the accuracy of NmF2 and hmF2 derived from F-3/C occultations has been assessed from a global long-term study. Based on an improved retrieval technique, the potential of radio occultations for ionospheric sounding has been demonstrated. Depending on the number of LEO satellites and the geometrical constellation, a main advantage of the RO technique compared to ionosonde measurements can be found in the global availability, in particular by providing data over the oceans. The outstanding role of occultations for multi-dimensional ionosphere modeling will be additionally strengthened with the launch of the F-3/C follow-on mission Formosat-7/COSMIC-2 scheduled for 2016.

Acknowledgments

The authors would like to thank the following services and institutions for providing the input data for this study: CDAAC and TACC (Formosat-3/COSMIC observations and orbits), IGS and CODE (GPS orbits and IONEX data) as well as NOAA/NESDIS (SPIDR ionosonde measurements and magnetic North Pole coordinates). Special thanks go to the IonSAT research group (UPC, Barcelona, Spain) for their support of any kind. Moreover, we are grateful to the editor and the reviewers for their constructive feedback and suggestions to improve the paper.

The editor thanks Ljiljana R. Cander and an anonymous referee for their assistance in evaluating this paper.

This work was supported by the German Research Foundation (DFG) and the Technische Universität München within the funding programme Open Access Publishing.

References

All Tables

All Figures

Fig. 1. Improved Abel inversion under consideration of the separability hypothesis with integral discretization by means of the onion shell structure. Colored arrows illustrate the discretized integration steps between adjacent layers. This draft shows the composition with measurements at two epochs i1 and i2 yielding the associated tangent points TP1 and TP2 that are located at the mean radii of the layers j1 and j2.

In the text

	Fig. 2. Flowchart for the validation procedure providing an overview of data sources (green boxes), the retrieval (yellow boxes), and assessment steps (blue boxes).
In the text

	Fig. 3. Number of retrieved F2 layer peaks at each investigated day during the study period (10 June 2006–29 April 2014).
In the text

	Fig. 4. Daily distribution of the electron density peaks derived from EDPs of FM-1 (left, 109 peaks) at 10 June 2006 vs. all F-3/C satellites FM-1 to FM-6 (right, 2447 peaks) at 5 June 2007.
In the text

	Fig. 5. Trend of the global NmF2 (left) and hmF2 (right) as derived from global averages and standard deviations of the retrieved F2 layer peak computed for the days under investigation. The data is plotted for the study period between 10 June 2006 and 29 April 2014.
In the text

	Fig. 6. Comparison of scaled daily NmF2 (red squares) and hmF2 (green squares) averages with the daily sunspot number R (top, blue squares) and daily solar flux F10.7 index (bottom, blue squares) during the study period (10 June 2006–29 April 2014). The unit of the y-axis units refers to the sunspot number and F10.7 and the NmF2 as well as hmF2 samples have been scaled.
In the text

	Fig. 7. Daily RMS of the retrieved peak density NmF2 (left) and corresponding height hmF2 (right) derived from radio occultations between 10 June 2006 and 29 April 2014.
In the text

	Fig. 8. Daily averaged relative differences between global NmF2 (left) and hmF2 (right) observations as measured by ionosondes and F-3/C radio occultations between 10 June 2006 and 29 April 2014.
In the text

	Fig. 9. Variation of the magnetic equator in geographic coordinates during 2006 (outermost blue curve) and 2014 (innermost green curve).
In the text

Fig. 10. Left column: Differences ΔhmF2 between ionosonde observations and F2 layer peaks retrieved from F-3/C occultations collected during the days under investigation between 10 June 2006 and 29 April 2014 within 12:00 and 16:00 LT and magnetic latitude sectors ΔφEQU = ±[0°, 20°] (top), ΔφMID = ±]20°,60°] (center) and ΔφPOL = ±]60°, 90°] (bottom). Right column: Retrieved NmF2 positions (red circles) and co-located ionosondes (blue triangles) within the defined magnetic latitude bands. The latitude intervals are exemplary depicted here under consideration of the magnetic North Pole from 2014.

In the text

Fig. 11. Left column: Differences ΔhmF2 between ionosonde observations and F2 layer peaks retrieved from F-3/C occultations. The data have been collected during the days under investigation between 10 June 2006 and 29 April 2014 within 12:00 and 16:00 LT and magnetic latitude sectors ΔφEQU = ±[0°, 20°] (top), ΔφMID = ±]20°,60°] (center) and ΔφPOL = ±]60°, 90°] (bottom). Right column: Retrieved hmF2 positions (red circles) and co-located ionosondes (blue triangles) within the defined magnetic latitude bands. The latitude intervals are exemplary depicted here under consideration of the magnetic North Pole from 2014.

In the text

	Fig. 12. Correlations between F2 layer peaks derived from F-3/C (x-axis) and ionosonde observations (y-axis) for NmF2 (left column) and hmF2 (right column). According to Figures 10 and 11, the subplots from top to bottom refer to the equatorial, mid-latitude, and polar regions.
In the text