Introduction

Sporadic E (Es) layers are ionization enhancements in the ionospheric E region at altitudes usually between 90 and 120 km (Whitehead 1970, 1989; Kelley 2009; Haldoupis 2011). A characteristic feature of Es layers is that they are thin layers with 0.6–2 km thicknesses and 10–1000 km horizontal extension. The Es layer could also appear as a non-uniform wave layer, multiple layers occurring simultaneously and separated at 6–10 km, or as composition of irregular elongated clouds of intense ionization within the lower E region (Zeng and Sokolovskiy 2010). Therefore, this gradient plasma instability can give rise to irregularities with scale sizes from a few tens of meters to a few kilometers (Shume et al. 2015), which can produce VHF/UHF radio wave scatters and cause signal scintillations on satellite radio communication and navigation system links (Ogawa et al. 1989). Over the past decades, Es layers have been extensively studied using a variety of instruments, including ionosondes, incoherent scatter radars and coherent scatter VHF radars operating from the ground and rocket payloads with in situ measurements (Smith and Matsushita 1962; Special issues on Radio Sci. 1972, 1975; Kelley 2009; Haldoupis 2012, and references therein). Several excellent reviews on theories and observations of ionospheric Es layers have been published (Whitehead 1970, 1989; Mathews 1998; Haldoupis 2011). In summary, they concluded that mid- and low-latitude Es layers are generally caused by the plasma convergence effect of the neutral wind shear through ion-neutral collision process and geomagnetic Lorentz forcing, i.e., the wind shear theory. Equatorial Es layers arise from the gradient instability and depend on the electro-jet current. Additionally, high-latitude and auroral Es layers are formed in appropriate large-scale convective electric field structures with the wind system playing a lesser role.

From 1995 to 1997, the proof-of-concept Global Positioning System/Meteorology (GPS/MET) experiment successfully demonstrated active limb sounding of the earth’s neutral atmosphere and ionosphere via GPS RO observations from a Low Earth Orbit (LEO) satellite (Ware et al. 1996; Kursinski et al. 1997). The GPS RO technique provides an innovative view on the ionosphere including Es layer using horizontal limb sounding. After GPS/MET, there have been many successful GPS-LEO RO missions, such as the German CHAMP (CHAllenging Mini-satellite Payload), the U.S–German GRACE (Gravity Recovery and Climate Experiment), the GPS-IOX (Ionospheric Occultation Experiment), the Argentinean SAC-C (Satellite de Aplicaciones Cientificas-C), the Taiwan-U.S. FS3/COSMIC, or the GRAS (GNSS Receiver for Atmospheric Sounding) experiment on the European Metop-A satellite (Wickert et al. 2001, 2009; Tapley et al. 2004; Schmidt et al. 2006; Anthes et al. 2008; von Engeln et al. 2011). Based on the data from these missions, several Es-layer studies were presented (Hocke et al. 2001; Wu et al. 2005; Arras et al. 2008, 2009; Zeng and Sokolovskiy 2010; Chu et al. 2014). Hocke et al. (2001) observed small-scale amplitude and phase fluctuations from GPS/MET data and confirmed that the Es irregularities occur mainly at heights 90–110 km and their seasonal occurrences have a maximum in the summer hemisphere. Wu et al. (2005) obtained the amplitude and phase scintillations from CHAMP data and studied the Es-layer morphology, including their global coverage, diurnal and seasonal variations, and their activities depending on geomagnetic dip angle. Using the RO amplitude measurements during the CHAMP, GRACE, and COSMIC missions, Arras et al. (2008) also confirmed strong seasonal variations with highest Es occurrence rates (up to 50%) during summer at mid-latitudes and showed the effect of the Southern Atlantic Anomaly (SAA) in geomagnetic field intensity on Es occurrence reduction. Arras et al. (2009) evaluated semidiurnal tidal signature in Es occurrence rate at higher mid-latitudes. Zeng and Sokolovskiy (2010) modeled RO signals through Es clouds and compared the resulting amplitude structures with observations from the FS3/COSMIC. Using the 1 Hz RO amplitude and excess phase data from the COSMIC mission, Chu et al. (2014) obtained similar Es morphology as the previous investigations by Wu et al. (2005) and Arras et al. (2008). Chu et al. (2014) also presented that the Es occurrence is very likely attributed to the metallic ion flux convergence induced by the zonal wind shear (simulated by the Horizontal Wind Model 07) in the E region.

In the following section, we introduce the GPS-LEO RO amplitude and excess phase data sets from the FS3/COSMIC mission and describe our analysis method to derive amplitude standard deviation and relative electron density profiles for E- and/or Es-layer observations. We also describe and discuss the criteria to identify Es events. In the subsequent section, we analyze and organize different types of Es occurrence distributions into typical geomagnetic equatorial, mid-latitude, and auroral zone morphologies. During the last sunspot cycle from mid-2006 to the end of 2017, the Es-event climatology, namely its variation with identified zone (referring geomagnetic latitude), altitude, season, and local time have been documented. In conclusion, we summarize some typical RO Es-event characteristics and point out several problems to be addressed in future studies.

FS3/COSMIC data analysis of Es-layer observation

The FS3/COSMIC is a joint Taiwan-U.S. constellation, consisting of six identical LEO micro-satellites, that were launched in mid-April, 2006, and were still operating at the end of 2017 (Anthes et al. 2008). The main payload of each spacecraft is a GPS receiver, developed by the Jet Propulsion Laboratory. It utilizes four GPS antennas: two limb viewing occultation antennas for atmospheric/lower ionospheric remote sensing at 50 Hz, and two slant observing antennas for Precise Orbit Determination (POD) and also ionospheric remote sensing at 1 Hz (Schreiner et al. 2007). The use of both limb viewing and POD antennas in parallel could produce simultaneous occultation observations, which consist of two separate sets of limb-viewing links with different altitude resolutions of about 0.064 and 3.2 km at E-region altitudes and different vertical sounding ranges starting from about 125 km altitude and the LEO altitude, respectively, down to the earth’s surface. Because a typical Es layer could have a thickness of less than 1 km as introduced in the first section, we propose to observe and then detect Es layer using the 50 Hz data. One disadvantage of these data is that its upper limits of the GPS-LEO tangent point altitude are irregular distributed between 100 and 135 km and some of the profiles reach not high enough to detect high-altitude Es layers. The intermediate RO measurements include satellite orbits, carrier excess phase and carrier signal-to-noise ratio (SNR) amplitude for dual GPS L-band signals (f1 = 1575.42 MHz and f2 = 1227.60 MHz), i.e., the “atmPhs” data from the COSMIC Data Analysis and Archive Center (CDAAC, http://cdaac-www.cosmic.ucar.edu/) and the Taiwan Analysis Center for COSMIC (TACC, http://tacc.cwb.gov.tw/cdaac/).

For the statistical analysis of Es-layer observation, we calculated mean SNR amplitudes using the 50 Hz “atmPhs” data at both L1 and L2 bands with a sliding altitude window of 20 km. Then, we derived normalized amplitude standard deviations separately as functions SDL1() and SDL2() of GPS-LEO tangent point altitude (r) with another sliding window of 1 km. We also identified an Es event if the following conditions are satisfied: (1) the peak SDL1 value (SDL1max) is higher than 0.2, (2) the SDL1 value at the top “atmPhs” RO link, i.e., SDL1(rtop), is less than 0.2, and (3) the altitude with peak SDL1 is higher than 80 km. We note that the peak SDL1 criteria of 0.2 to identify the Es event is the same as those used in former investigations by Arras et al. (2008) and Arras and Wickert (2017). However, the sliding window width of 1 km is narrower than 2.0 and 2.5 km (used by the above investigations) and is chosen in accordance with the first Fresnel zone, where radio waves occur diffracting and scattering if the irregularity scale sizes in the medium are comparable or smaller than the first Fresnel zone (Yeh and Liu 1982). Additionally, to include multiple Es-layer events, we did not limit our Es-layer thickness but exclude the SDL1(rtop) > 0.2 observations which could be induced by F-region irregularities but not by the Es layer. When ionospheric F-region irregularities occur along the limb-viewing GPS-LEO rays, the observed amplitude profile could have strong fluctuations at the altitudes of the F and also E region. We also note that some SDL2 values are not reliable, and therefore, cannot be used to identify and characterize Es events because of weak and/or insensitive L2-band signals.

Figure 1 shows two FS3/COSMIC “atmPhs” RO observations: one single Es-layer event (Fig. 1 top) and the other multiple Es-layer events (Fig. 1 bottom). The figures show the limb-viewing SNR amplitude profiles in GPS-LEO tangent point altitude (from 60 km to the top of each observation) in black and gray colors for the L1 and L2 bands, respectively, and the resulting normalized amplitude standard deviation (SDL1 and SDL2) profiles in blue and green. As observed, there are significant amplitude fluctuations from L1-band signals but not from L2-band signals, which are much weaker and do not have enough sensitivity to derive reliable SDL2 values. In Fig. 1 the retrieved electron density (Ne) profiles are also shown in red color. The top panel displays a thin Ne-enhanced layer at 112 km altitude and with about 3 km thickness where SDL1s > 0.2, and the bottom panel displays multiple Es layers occurring simultaneously from 115 to 95 km altitude. These Ne profiles were inverted using the Abel integral transform (Tricomi 1985) through relative GPS-LEO Total Electron Contents (TECs) under assumptions that include locally spherical symmetry on ionospheric Ne, straight-line ray propagation, static ionosphere, known initial Ne at the top RO observation altitude, and an earth’s spherical shape over the occultation duration. An Ne profile can thus be derived from top to bottom in altitude by applying the following integral equation to the relative path TEC values (TEC′), which are derived using L1-band excess phase measurements (S1) in meters,

Fig. 1
figure 1

Two examples of the FS3/COSMIC RO observations are displayed in the top and bottom panels separately for single and multiple Es-layer events. Both images show limb-viewing L1- and L2-band SNR amplitudes versus GPS-LEO tangent point altitudes in black and gray, the corresponding normalized amplitude standard deviation profiles, i.e., SDL1 and SDL2 profiles, in blue and green, and the retrieved and biased Ne profile in red. The retrieved Ne values are shifted by 106 el/cm3 for better visualization

$$\begin{aligned} {N_{\text{e}}}({r_t}) & ={N_{\text{e}}}({r_{{\text{top}}}}) - \frac{1}{\pi }\int_{{\;{r_t}}}^{{{r_{{\text{top}}}}}} {\frac{{{\raise0.7ex\hbox{${{\text{dTE}}{{\text{C}}^\prime }(r)}$} \!\mathord{\left/ {\vphantom {{{\text{dTE}}{{\text{C}}^\prime }(r)} {{\text{d}}r}}}\right.\kern-0pt}\!\lower0.7ex\hbox{${{\text{d}}r}$}}}}{{\sqrt {{r^2} - r_{t}^{2}} }}{\text{d}}r} \\ {\text{TE}}{{\text{C}}^\prime } & = - \frac{{{f_1}^{2}{S_1}}}{{40.3}}+B \\ \end{aligned}$$
(1)

where rt is a tangent GPS-LEO path, Ne(rtop) is an assumed Ne at top RO observation altitude, and B is the TEC bias during an occultation observation. The standard Abel-inversion technique applying on path TECs were independently developed by Hajj and Romans (1998) and Schreiner et al. (1999). The Ne value at a tangent point’s radial distance rt can be computed in a recursive way starting from the outer, i.e., upper, rays with a known and/or assumed initial Ne. In general, the assumption of spherical symmetry of local ionosphere used in the Abel integral transform is not realistic. It was found that the retrieval Ne errors were caused mostly by the spherical symmetry assumption imposed on the ionosphere and could be accumulated larger from the top to bottom of an RO observation. For the GPS RO observations to remotely sense the ionosphere, several investigations (Hajj and Romans 1998; Schreiner et al. 1999; Lei et al. 2007; Tsai et al. 2009) examined and estimated the accuracy of retrieved Abel-inversion Ne profiles with measurements from incoherent scatter radars and globally distributed ionosondes. Statistical comparisons of peak ionospheric Ne indicated a mean difference of a few percent and a root-mean-square difference of about 20%.

We note that the normal peak Ne at the E region is more or less 10% of the peak ionospheric Ne, i.e., the retrieved E-region Ne values using the standard Abel inversion are not reliable because they could be even less than the retrieved Ne error. Nicolls et al. (2009) developed an alternative approach for estimate E-region Ne using the Abel inversion while the F-region specification is obtained from the ionospheric data assimilation four-dimensional (IDA4D) algorithm. Nevertheless, this study does not attempt to retrieve the absolute Ne values but the Ne differences at E- or Es-layer peaks, derived by the difference in relative peak Ne value and the corresponding relative medium Ne value using a window of ± 10 km. As shown in Fig. 1, the Ne enhancement ranges are consistent with the altitudes where SDL1 values are larger, and both SDL1(rtop) values are much less than 0.2. It can be assumed that strong amplitude fluctuations are induced by the E-region Ne enhancement, i.e., the Es layer. For a statistical distribution study, we could reasonably position each Es event nearby at the GPS-LEO tangent point location having peak SDL1 value. We focus here on the altitude range from 60 to 130 km, i.e., the E-layer region, to study Es morphology and use L1-band data only hereafter.

The FS3/COSMIC GPS RO observation locations cover the whole globe except for two small areas around the North and South Poles. Vertical profile information can be obtained, for instance, over oceans or over other depopulated areas. Compared with other ground-based ionospheric-sounding techniques and/or in situ measurements from rocket payloads, the GPS occultation technique has the great advantage of all-altitude capability and global coverage. Because the FS3/COSMIC program consists of six identical LEO micro-satellites and each satellite was equipped with two occultation antennas viewing forth and back separately, the FS3/COSMIC provides much more RO observations than other missions. Figure 2 shows five-day averaged numbers of the FS3/COSMIC 50 Hz RO observations, i.e., atmPhs data, collected during the years of the last sunspot cycle from mid-2006 to the end of 2017. As shown by the blue profile in Fig. 2, the FS3/COSMIC performed more than 2000 “complete” observations per day from 2007 to 2009, and then the numbers of complete observations declined to 1000 in 2015 and became less, year by year, due to FS3/COSMIC spacecraft bus and GPS receiver/antenna payload degradation. In this study, we defined a “complete” atmPhs RO observation where its SDL1 and Ne profiles can be successfully retrieved at the tangent point altitudes from 60 km to the RO observation top, which must be higher than 115 km include most Es-layer altitudes. These complete observations constitute more than 70% of the total atmPhs RO measurements (shown by the black profile in Fig. 2) which also include failed observations mostly because the obtained top RO observation altitudes are lower than a ceiling of 115 km, and/or the recorded average L1-band SNR amplitudes are less than 3000. Meanwhile, the 5-day averaged numbers of Es events and the corresponding percentages of complete observations are also shown by the red and green profiles, respectively, in Fig. 2. We note that Es events generally occurred more or less 31% of complete atmPhs RO observations, and two diffuse seasonal maximums of Es occurrence rate occurred at the middle and end of each year. More seasonal Es occurrence behavior will be represented and discussed in the later sections.

Fig. 2
figure 2

Five-day averaged numbers (left y-axis) of the FS3/COSMIC RO observations are shown in red for observations with Es events, in blue for complete observations retrieving SDL1 and Ne profiles successfully. The black profile shows the numbers of observations including failed ones because the obtained top RO observation altitudes are lower than a ceiling of 115 km, mostly, and/or the recorded average L1-band SNR amplitudes are less than 3000. The corresponding percentages (right y-axis) of Es events rated with complete observations are also shown by the green profile

Figure 3 shows a scatter plot of retrieved Ne difference (dNe) at E- or Es-layer peak in a log scale versus peak normalized L1-band amplitude standard deviation (SDL1max) in a linear scale and its least-square fitting line on Es events (SDL1max > 0.2) based on the complete FS3/COSMIC atmPhs RO observations from the day of year (DOY) 161–190, 2008. From the correlative analysis of logarithm dNe versus SDL1max of Es events, as shown in the figure, the result has a high linear correlation coefficient of 0.70 and a fitting line slope of 1.29. We found that Es-event SDL1max values are approximately proportional to the logarithms of corresponding dNe values and the obtained dNe values have a minimum of about 104 el/cm3 and are mostly less than 3 × 105 el/cm3. The other correlative analyses using different FS3/COSMIC atmPhs RO observations have similar results as shown in the figure. The figure also shows the numbers of complete RO observations and their corresponding accumulated percentages as functions of SDL1max. We note that the Es occurrence rate based on the FS3/COSMIC RO observations around mid-2008 is about 40%, which is larger than the average rate (about 31%).

Fig. 3
figure 3

Scatter plot of retrieved Ne difference (dNe, left y-axis) at E- or Es-layer peak versus peak SDL1 (SDL1max) and its least-squares fitting line (yellow) on Es events (SDL1max > 0.2) are shown based on the complete FS3/COSMIC atmPhs RO observations from DOY 161–190, 2008. The right y-axis refers to the numbers of complete RO observations (red line) and their corresponding accumulated percentages (green line), as functions of SDL1max in 0.01 bins

Figure 4 shows another scatter plot of retrieved peak SDL1 (SDL1max) versus Es-layer thickness (EsW) based on the same complete FS3/COSMIC atmPhs RO observations used in Fig. 3. The image also shows the numbers of Es events and their corresponding accumulated percentages as functions of EsW in 0.1 km bins. We found that the detected Es events are mostly thick from 0.1 to 10 km and have a peak distribution at 1.2 km. There are about 80% Es events having an Es layer thickness less than 6 km and possibly being not observed from 1 Hz limb-viewing data, which has an altitude resolution of about 3 km. We also note that the use of a narrower sliding altitude window, i.e., 1 km in this study, on deriving SDL1 profiles could narrow the Es-layer thickness and increase the peak SDL1 value, and thus increase the Es occurrence rate. Additionally, as shown in Fig. 4, there are about 7% Es events having an Es-layer thickness larger than 10 km, which was used as an identifying criterion in the former investigations by Arras et al. (2008) and Arras and Wickert (2017). This Es occurrence distribution in E-layer thickness can explain that the obtained maximum Es occurrence rate (up to 80% and presented in the later section) is much higher than the obtained rate (about 15%) from Chu et al. (2014) and is also higher than the observed rate (up to 50%) from Arras et al. (2008) Arras and Wickert (2017).

Fig. 4
figure 4

Scatter plot of derived peak SDL1 (SDL1max referred to the left y-axis) versus Es-layer thickness (EsW) based on the same complete FS3/COSMIC atmPhs RO observations in Fig. 3. The image also shows the numbers of Es events and their corresponding accumulated percentages (red and green, respectively, and referred to the right y-axis) as functions of EsW in 0.1 km bins

Results and discussions on global Es-event morphology

Based on the FS3/COSMIC atmPhs, i.e., 50 Hz, RO observations described in the previous section, Fig. 5 shows 12 monthly global Es occurrence distributions in 2008. There are 229,003 Es events from 736,507 complete E-layer RO observations, i.e., 31.09 percent on average. The coded colors represent the Es occurrence rates from 0 to 80% within every 5° by 5° in the geographic bin. The white bins around the North and South Poles are not taken into account and indicate less than five complete E-layer RO observations located within each bin area. The monthly global Es occurrence distributions of other years are similar to Fig. 5 in outline. The monthly Es occurrence distribution pattern depends highly on the geomagnetic and/or geographical location. We obtained five major zones in which the seasonal and diurnal Es occurrence patterns are markedly different. The five zones include geomagnetic equatorial zone (Zone A: − 5° < ML < 5°), two extended geomagnetic mid-latitude zones (Zone B: 15° < ML < 55°, and Zone C: − 55° < ML < − 15°), and two auroral zones (Zone D: 70° < ML, and Zone E: ML < − 70°), which are shown and enclosed by black lines in each occurrence distribution map of Fig. 5.

Fig. 5
figure 5

Twelve monthly global occurrence distribution of Es events (SDL1max > 0.2) in 2008. The coded color represents the Es occurrence rate from 0 to 80% within every 5° by 5° in the geographic bin. Five typical zones enclosed by black lines are chosen and identified based on the occurrence statistics, and the curves of geomagnetic latitudes at 70°, 55°, 15°, 5°, − 5°, − 55°, and − 70° are also shown and labeled

We found that Es events are strong and dense in the two extended geomagnetic mid-latitude zones, weak in the geomagnetic equatorial zone, and weak and minimum in the two auroral zones. We also obtained that maximum Es activities in two mid-latitude zones occurred during hemispheric summer seasons. In the northern extended geomagnetic mid-latitude zone (Zone B), Es activity had a maximum occurrence rate of > 80% and an occurrence depression around the American area (the longitude sector of 70°–120° W), where the Es occurrence rates were 20% lower (at most) than anywhere else along the zone bands. In the southern extended geomagnetic mid-latitude zone (Zone C), the Es activity had a maximum occurrence rate (approaching 80%) less than the maximum rate in Zone B. Meanwhile, there is a stronger depression of occurrence around the South Atlantic and African areas (the longitude sector of 40°W–50°E), where the Es occurrence rates were 50% lower (at most) than elsewhere in Zone C. The two occurrence depression areas along the mid-latitude zone band of summer hemisphere were also obtained in the former investigation from Arras et al. (2008). This Es occurrence depression result can be explained by the wind shear theory (Whitehead 1961) from a connection with the horizontal component of the earth’s magnetic field which also shows distinctive depressions in the two areas. We found that two Es occurrence depression areas in geomagnetic mid-latitude zones shown in Fig. 5 are typical and coincident with the horizontal component of the International Geomagnetic Reference Field (IGRF) (Finlay et al. 2010; http://www.ngdc.noaa.gov/ geomagmodels/IGRFWMM.jsp), which is not shown here. During hemispheric winter seasons, on the other hand, two geomagnetic mid-latitude zones had fairly general low-Es activities and displayed lower occurrence rates where the geomagnetic latitudes are higher. More validation and discussion will be presented in later sections.

To study more detailed Es occurrence zone characteristics related to season and local time, Fig. 6 depicts the temporal profiles of 5-day averaged Es occurrence rates and Es occurring local times for the five identified zones separately during the last sunspot cycle from mid-2006 to the end of 2017. For each left panel, the yellow background denotes the “yellow” season from April to September annually, and the white background denotes the “white” season from October to March. Note that the “yellow” seasons include the northern hemispheric summers, and, on the other hand, the “white” seasons include the southern hemispheric summers. The five right panels in Fig. 6 depict the corresponding zone Es occurrence rates as functions of local time in red and blue lines for “yellow” and “white” seasons, respectively. We note that the duration from 2007 to 2009 is a solar minimum period and the other duration from 2012 to 2014 is a solar maximum period. Nevertheless, from our studies, there are no typical correlations for the five identified zone Es occurrences with solar activity. The implications of the characteristic differences in zone Es occurrences are discussed below.

Fig. 6
figure 6

5-day averaged Es occurrence rates (red lines referring to the left y-axes) and Es occurring local times (blue dots referring to the right y-axes) for the five identified zones, i.e., Zones A, B, C, D, and E, during the time from mid-2006 to the end of 2017. For each left panels, the yellow background denotes the “yellow” season from April to September annually, and the white background denotes the “white” season from October to March. The five right panels also show the Es occurrence rates for the corresponding zone as functions of local time in red lines, and blue lines for both “yellow” and “white” seasons, respectively

Geomagnetic equatorial zone (Zone A)

From Fig. 6 (top row, left panel), the 5-day averaged temporal profile of Es occurrence rate in the geomagnetic equatorial zone (Zone A) shows more or less 30% occurred from mid-2006 to the end of 2017. There is no clear seasonal variation but it shows small diffuse maxima approaching 40% in the middle of each year from 2010 to 2015, which include the solar maximum period of the last cycle. Meanwhile, we note that, as shown in the top row (right panel), showing the Es occurrence distributions in local time for Zone A, the diurnal occurrence behavior is prevalent and similar during both “yellow” and “white” seasons, and the Es occurrences are slightly more active during the day than night. Es activity has a minimum occurrence rate of about 7% occurring around sunrise (i.e., 6 LT), increases slowly and steadily after that, and reaches a maximum of about 45% around sunset, i.e., 18 LT. After sunset, Es activity decreases rapidly to 30% within a couple of hours, stays more or less at 30% from 20 to 24 LT, and then decreases again slowly to reach a minimum occurrence rate of 7% at sunrise.

Two extended geomagnetic mid-latitude zones (Zones B and C)

From Fig. 6 (second and third rows from the top, left panels), it is generally observed that the five-day averaged temporal profile of Es occurrence rate in each geomagnetic mid-latitude zone is often similar by years, and the Es activity is dense in hemispheric summers but relatively weak in hemispheric winters. In the northern (southern) geomagnetic mid-latitude zone, Es activity increases very rapidly during the first two months of “yellow” (“white”) seasons, i.e., from April to May (from October to November), the maximum activity approaches a rate of up to 70% (60%) and occurs around the middle (the end) of each year, and then Es activity decreases during the last 2 months of “yellow” (“white”) seasons, i.e., from August to September (from February to March). It also appears that, in the northern (southern) geomagnetic mid-latitude zone, the Es activity maintains occurrence rates of more or less 25% during “white” (“yellow”) seasons, and that during the 1st month and around the middles of “white” (“yellow”) seasons there are two minor maxima of activity. We note that the Es activity is stronger in northern hemispheric summers than in southern hemispheric summers but, on average, has a similar occurrence rate of about 25% in both northern and southern hemispheric winters. Meanwhile, the local time variations of Es occurrence rate are shown in the respective right panels. It is obvious that semidiurnal behavior of Es occurrence is dominant in both northern and southern hemispheric summers, but the diurnal behavior is prevalent in both northern and southern hemispheric winters. As shown by the red (blue) curve in the second (third) right panel for the local time Es occurrence rates of northern (southern) hemispheric summers, two maxima are presented at about 10 and 18 LT, where the former is slightly stronger. Meanwhile, a minimum occurrence rate of 25% happened before 6 LT at about sunrise. On the other hand, as shown by the blue and red curve in the second and third right panels, respectively, diurnal occurrence behavior is prevalent and similar in both northern and southern hemispheric winters. Es activity is more active during the day than night. Furthermore, Es activity has a minimum occurrence rate of about 7% occurred around sunrise (i.e., 7 LT), increases after that and reaches a maximum of about 40% around sunset, i.e., 17 LT. After sunset, Es activity decreases rapidly to 25% in a couple of hours, stays more or less at 25% from 20 to 24 LT, and then decreases again slowly to reach a minimum occurrence rate of 7% at sunrise. We will show altitude versus local time distribution of Es events later, and that the vertical propagation of tidal waves with different periods is responsible for the patterns as shown in the 2nd and 3rd right panels.

Two auroral zones (Zones D and E)

From Fig. 6 (fourth and fifth rows from the top, left panels), it is generally observed that the 5-day averaged temporal profiles of Es occurrence rates in both auroral zones are slightly less than 20% on average and have similar seasonal variations. The seasonal occurrence behavior depicts two small diffuse maxima of more or less 30% around the middles of “yellow” and “white” seasons, i.e., summer and winter seasons, with the stronger one in summer. It also has a minimum of < 10% that happened during March (September), i.e., the last month of “white” (“yellow”) season in the northern (southern) auroral zone, except for the years of 2016 and 2017, when the RO observations are much less and not uniformly distributed in season and local time due to spacecraft bus and/or payload degradation in several FS3/COSMIC satellites. Meanwhile, as shown in the respective right panels for zone Es occurrence distributions in local time, it is obvious that diurnal Es occurrence patterns are dominated and similar in both auroral zones and during both “yellow” and “white” seasons. Es activity is most prevalent during the night, has a maximum occurrence rate of about 35% at 22 LT, i.e., 2 h before midnight, stays at low occurrence rate more or less at 7% from 6 to 15 LT, and has a minimum occurrence rate at about 9 LT.

As discussed, the Es activity is geomagnetically controlled and is often similar by year. To study more detailed Es-layer characteristics related to geomagnetic latitude and season, Fig. 7 (top to bottom) depicts the geomagnetic latitude versus day of year variations in Es occurrence rate, mean SDL1max (peak normalized amplitude standard deviation) of Es layers, and mean EsW (Es-layer thickness) separately, in which averaged values in the years from 2006 to 2017 are presented in bins of 2° geomagnetic latitude and 2 days. Based on the same complete FS3/COSMIC atmPhs RO observations as in Fig. 2, there are 1806,841 Es events from 5828,710 complete E-layer RO observations totally, i.e., 31.1 percent on average. The white bins in Fig. 7 are not considered and indicate the Es events located within each bin being less than five.

Fig. 7
figure 7

Geomagnetic latitude versus day of year variations in Es occurrence rate (top), mean Es-layer SDL1max (middle), and mean EsW (bottom). The coded color represents the averaged occurrence rate, the mean SDL1max, and the mean EsW values of Es layers over the years 2006–2017 using bins of 2° geomagnetic latitude and 2 days

The top panel depicts the same seasonal zone Es occurrence variation results as discussed in Fig. 6. In the geomagnetic equatorial zone, a small diffuse maximum occurs around the middle of the year. In the northern (southern) geomagnetic mid-latitude zone Es activity increases rapidly from April to May (from October to November), approaches a maximum around the middle (the end) of year, decreases from August to September (from February to March), and has two minor maxima in October (April) and around the end (middle) of year. In both auroral zones, Es activity has two small diffuse maxima around the middles of hemispheric summer and winter with the stronger one in summer. In the northern (southern) auroral zone the minimum seasonal occurrence rate happened around March (September). Additionally, as shown in the top panel, in the northern (southern) hemisphere and geomagnetic low- and mid-latitude regions Es occurrences distribute over higher latitudes from April (October) to the middle (end) of year, then decline back to lower latitudes until September (March); Es activity stays low in other seasons and also has lower occurrence rate where the geomagnetic latitude is higher except for near the geomagnetic equator. The peak Es occurrence rates move in geomagnetic latitude from about 20° (− 20°) in April (October) to about 35° (− 35°) at the middle (end) of year, back to about 20° (− 20°) in September (March), and staying 15°–10° (− 15° to − 10°) in the other seasons. As shown in the middle and bottom panels, larger mean Es-layer SDL1max and mean EsW values occur where the Es occurrence rates are higher in the magnetic equatorial and mid-latitude zones. Nevertheless, the mean Es-layer SDL1max and EsW values in two auroral zones and during hemispheric summers are comparatively larger than those at low- and mid-latitudes.

Figure 8 shows the Es-event altitude versus local time (0–24 LT) distributions summed in a targeted month, e.g., January, February, etc., during the years from 2006 to 2017 and for five identified zones. Because there is no reference on the characterized altitude of GPS RO observations except for the peaks of Es layers, we show the Es-event altitude distributions not in occurrence rate but occurrence number. Nevertheless, because of much few Es events in the geomagnetic equatorial and auroral zones, we present the summed numbers of Es events that occurred in January, February, and so on, from 2006 to 2017. Consequently, it is easier to recognize the number differences on logarithm than on linearity. In general, most Es-layer peaks are located at 90–120 km altitudes in both geomagnetic mid-latitude zones and at 85–115 km altitudes in the geomagnetic equatorial zone and the two auroral zones.

Fig. 8
figure 8

Es-event altitude versus local time (0–24 LT) distributions summed in a targeted month, e.g., January, February etc., during the years from 2006 to 2017, and for five identified zones. The coded color represents the corresponding logarithm of summed number of Es events in 2 km bins of altitude and 0.5 h bins of local time

From the 2nd and 3rd panels in Fig. 8, counting from the top, we find that for a local season, the Es-event altitude versus local time distributions are often similar in both extended geomagnetic mid-latitude zones (Zones B and C). During hemispheric summers, the Es activity is dominated primarily by a semidiurnal feature, which is generally believed to be induced by east–west zonal winds in terms of semidiurnal tides (Whitehead 1989; Arras et al. 2009; Chu et al. 2014). The semidiurnal tides generally start around 6 and 14 LT, continue for 14 h, and then fade out around 20 and 4 LT separately. As shown, the semidiurnal feature of the Es layers are distributed within a relatively narrow altitude range of about 10 km and gradually descend from 120 to 90 km at an average rate of about 2.1 km/h. We note that, from our analyzed results, the starting time of 14 LT for the afternoon tide is 4 h earlier than the result (18 LT) obtained in Chu et al. (2014), and the descending rate of the Es-layer altitude is close to the result (1.5–2.5 km/h) from the former investigation by Arras et al. (2009) but is faster than the result (0.9–1.6 km/h) from Chu et al. (2014). As discussed for Fig. 6, two Es occurrence rate maxima happened at about 10 and 18 LT, i.e., 4 h after the starts of the semidiurnal tides, with altitudes of about 105 km and 110 km, respectively. On the other hand, during mid-latitude hemispheric winters, a diurnal feature is prevalent as shown in the 2nd and 3rd panels of Fig. 8. The diurnal-feature Es layers are distributed within a wider altitude range of more than 20 km, but it is difficult to derive the corresponding vertical movement.

In the geomagnetic equatorial zone (Zone A), the Es activity is dominated by a diurnal behavior and is more active during afternoons, which is similar to the diurnal feature of Es occurrence in mid-latitude hemispheric winters. Note that the obtained numbers of Es events in Zone A are less than in the hemispheric winters of Zones B and C because of a smaller zone latitude range. However, the averaged Es occurrence rate in Zone A is about the same as those in the hemispheric winters of Zones B and C. The Es activity in the two auroral zones is also dominated by diurnal behavior but is more active at night.

Conclusions

We investigated the global morphology of ionospheric Es events using the 50 Hz amplitude and excess phase measurements from the FS3/COSMIC mission. A number of Es morphology properties are presented, and, as discussed, some of those were also obtained from earlier investigations using GPS RO observations. These are dense Es occurrences at geomagnetic mid-latitudes but weak over the geomagnetic equator, strong seasonal Es occurrence variations at mid-latitudes with highest rates during hemispheric summers, Es occurrence depletion at mid-latitudes around the American area and around the South Atlantic and African area during northern and southern hemispheric summers separately, semidiurnal behavior in mid-latitude hemispheric summers with a descending movement.

Nevertheless, several global Es morphology results from GPS RO observations are newly obtained and/or identified from this study. (1) The maximum normalized amplitude standard deviations for GPS-LEO limb-viewing observations on Es layers are approximately proportional to the logarithms of peak Ne differences. (2) The detected Es events had thickness mostly less than 6 km (about 80%) and had a peak distribution at 1.2 km. (3) Diurnal Es occurrence behavior is prevalent for all seasons in geomagnetic equatorial zone and has highest rates around 18 LT and lowest rates around 6 LT. (4) In two geomagnetic mid-latitude zones, two minor maxima of seasonal Es activity occurred in the 1st month, i.e., October for the north and April for the south, and around mid of hemispheric winters. (5) The Es semidiurnal tide behavior in mid-latitude hemispheric summers generally starts around 6 and 14 LT separately with layer altitudes descending from 120 to 90 km at an average rate of about 2.1 km/h. (6) The diurnal Es occurrence behavior is prevalent and similar in both northern and southern mid-latitude hemispheric winters although with lower amplitudes compared to the semidiurnal features in hemispheric summers. (7) The diurnal Es occurrence variation in two auroral zones is first observed using GPS RO measurements, and Es activity in two auroral zones is also dominated by diurnal behavior but is more active at night with the highest rate around 22 LT. In future investigations, a follow-on program called FS7/COSMIC2 is in progress with proposed satellite launching in the 2nd half of 2018. Similar to the FS3/COSMIC, it is a six-satellite constellation mission orbiting at 24° inclination and 550 km altitude and enhanced GPS receiver capability to accommodate signals from GPS and GLONASS satellites. The FS7/COSMIC2 constellation mission is designed for providing more than 5000 GPS/GLONASS RO observations per day and focusing on the regions between − 50° and 50° latitude. It is expected that denser RO observations could be used to structure and/or model ionospheric F- and Es-layer irregularities in the equatorial, low- and mid-latitude zones.