Assessing the Performance of GPS Precise Point Positioning Under Different Geomagnetic Storm Conditions during Solar Cycle 24

The geomagnetic storm, which is an abnormal space weather phenomenon, can sometimes severely affect GPS signal propagation, thereby impacting the performance of GPS precise point positioning (PPP). However, the investigation of GPS PPP accuracy over the global scale under different geomagnetic storm conditions is very limited. This paper for the first time presents the performance of GPS dual-frequency (DF) and single-frequency (SF) PPP under moderate, intense, and super storms conditions during solar cycle 24 using a large data set collected from about 500 international GNSS services (IGS) stations. The global root mean square (RMS) maps of GPS PPP results show that stations with degraded performance are mainly distributed at high-latitude, and the degradation level generally depends on the storm intensity. The three-dimensional (3D) RMS of GPS DF PPP for high-latitude during moderate, intense, and super storms are 0.393 m, 0.680 m and 1.051 m, respectively, with respect to only 0.163 m on quiet day. RMS errors of mid- and low-latitudes show less dependence on the storm intensities, with values less than 0.320 m, compared to 0.153 m on quiet day. Compared with DF PPP, the performance of GPS SF PPP is inferior regardless of quiet or disturbed conditions. The degraded performance of GPS positioning during geomagnetic storms is attributed to the increased ionospheric disturbances, which have been confirmed by our global rate of TEC index (ROTI) maps. Ionospheric disturbances not only lead to the deteriorated ionospheric correction but also to the frequent cycle-slip occurrence. Statistical results show that, compared with that on quiet day, the increased cycle-slip occurrence are 13.04%, 56.52%, and 69.57% under moderate, intense, and super storms conditions, respectively.


Introduction
Geomagnetic storms, i.e., large-scale disturbances in the Earth's near-space environment, are caused by the enhanced solar wind and its interaction with the magnetosphere-ionosphere-thermosphere system. According to the Dst index derived from near-equatorial geomagnetic measurement, the principal feature of a geomagnetic storm represents an obvious decrease of the horizontal intensity of Earth's magnetic field followed by a recovery [1]. The disturbance that drives geomagnetic storm can be due to the occurrence

Data Sources
The Dst index is a good quantitative measure of the intensity of the geomagnetic storm [20]. In this work, we focus only on the geomagnetic storm with minimum Dst (Dst min ) ≤ −50 nT, which is commonly adopted as the threshold by many past studies [21,22]. Further, if a period of high activity showed multiple Dst ≤ −50 nT, we arbitrarily treat them as a single storm event if the minima were separated by less than 24 h, rather than define each minimum as a separate storm [23]. Depending on Dst index, the storms are classified as moderate (−100 nT < Dst min ≤ −50 nT), intense (−200 nT < Dst min ≤ −100 nT), and super storm (Dst min ≤ −200 nT). describes the data and methodology in details. Section 3 presents GPS DF and SF PPP solutions under three typical geomagnetic storms conditions, i.e., moderate, intense, and super storms, during solar cycle 24 from 1 January 2011 to 31 December 2017. Furthermore, the corresponding analyses for the experiments are also shown in Section 3. Finally, the discussion and conclusion are given in Sections 4 and 5, respectively.

Data and Methodology
This section describes the data sources, including geomagnetic activity, solar activity, global TEC maps, and GPS measurements, used in our study. In addition, the data processing models and strategies used to perform GPS DF and SF PPP are presented in this section.

Data Sources
The Dst index is a good quantitative measure of the intensity of the geomagnetic storm [20]. In this work, we focus only on the geomagnetic storm with minimum Dst (Dstmin) ≤ −50 nT, which is commonly adopted as the threshold by many past studies [21,22]. Further, if a period of high activity showed multiple Dst ≤ −50 nT, we arbitrarily treat them as a single storm event if the minima were separated by less than 24 h, rather than define each minimum as a separate storm [23]. Depending on Dst index, the storms are classified as moderate (−100 nT < Dstmin ≤ −50 nT), intense (−200 nT < Dstmin ≤ −100 nT), and super storm (Dstmin ≤ −200 nT).    In this work, the GPS dual-frequency data provided by the IGS are used to investigate the effects of storms with different intensities on GPS PPP solutions (ftp://cddis.gsfc.nasa.gov/). There are about 500 continuous operation tracking stations distributed as shown in Figure 2. Besides, the global TEC maps with a time resolution of 5 min derived from the worldwide GNSS stations are used to analyze the ionospheric response to various geomagnetic storms. The TEC map data can be freely accessible from the CEDAR Madrigal database (http://cedar.openmadrigal.org/ftp/). In this work, the GPS dual-frequency data provided by the IGS are used to investigate the effects of storms with different intensities on GPS PPP solutions (ftp://cddis.gsfc.nasa.gov/). There are about 500 continuous operation tracking stations distributed as shown in Figure 2. Besides, the global TEC maps with a time resolution of 5 min derived from the worldwide GNSS stations are used to analyze the ionospheric response to various geomagnetic storms. The TEC map data can be freely accessible from the CEDAR Madrigal database (http://cedar.openmadrigal.org/ftp/).

GPS PPP Model
The raw observables of dual-frequency GPS pseudorange and carrier phase between a receiver "r" and a satellite "s" are generally expressed as below [24]: where , and Φ , are the pseudorange and carrier phase on frequency in meters; is the geometric range in meters; and are the receiver and satellite clock offset; is the zenith tropospheric delay; is the mapping function; and are the frequency dependent signal delay for receiver and satellite; is the line-of-sight TEC with the frequency dependent factor = 40.3/ ; is the wavelength; is the float ambiguity and is the phase windup error in cycle; and and are the measurement noise of pseudorange and carrier phase, respectively. In this study, we adopted undifferenced PPP model with raw GPS measurements instead of ionospheric-free PPP model. Compared with the latter one, in the undifferenced model, individual signals of each frequency are treated as independent observables, thus avoiding noise amplification in the linear combinations. For the ionospheric delay in DF and SF PPP models, a priori ionospheric model with proper constraints is utilized as [25]: where is the vertical TEC of the ionospheric pierce point (IPP); denotes the ionosphere mapping function; (i = 0, 1, 2, 3, 4) are the coefficients used to describe the deterministic behavior of the ionospheric delay; and are the longitude and latitude difference between IPP and the approximate location of station; the scalar field represents the stochastic component from a second-order stationary process; and is the vertical ionospheric delay correction interpolated from global ionospheric map (GIM) with corresponding noise . In the data processing, all GPS observations have a sampling interval of 30 s, and the satellite elevation mask angle is set to 10°. GPS precise satellite orbit and clock products provided by IGS are adopted for PPP processing. The data processing models and strategies for GPS DF and SF PPP are summarized in Table 1. Note that in this study the experiments and results of DF PPP are introduced

GPS PPP Model
The raw observables of dual-frequency GPS pseudorange and carrier phase between a receiver "r" and a satellite "s" are generally expressed as below [24]: where P s r, f and Φ s r, f are the pseudorange and carrier phase on frequency f in meters; ρ is the geometric range in meters; t r and t s are the receiver and satellite clock offset; T z is the zenith tropospheric delay; α is the mapping function; b r and b s are the frequency dependent signal delay for receiver and satellite; I is the line-of-sight TEC with the frequency dependent factor β f = 40.3/ f 2 ; λ is the wavelength; N is the float ambiguity and ϕ is the phase windup error in cycle; and ε P and ε Φ are the measurement noise of pseudorange and carrier phase, respectively.
In this study, we adopted undifferenced PPP model with raw GPS measurements instead of ionospheric-free PPP model. Compared with the latter one, in the undifferenced model, individual signals of each frequency are treated as independent observables, thus avoiding noise amplification in the linear combinations. For the ionospheric delay in DF and SF PPP models, a priori ionospheric model with proper constraints is utilized as [25]: I(z) s r = a 0 + a 1 dL + a 2 dL 2 + a 3 dB + a 4 dB 2 + r s r (4) where I(z) s r is the vertical TEC of the ionospheric pierce point (IPP); γ denotes the ionosphere mapping function; a i (i = 0, 1, 2, 3, 4) are the coefficients used to describe the deterministic behavior of the ionospheric delay; dL and dB are the longitude and latitude difference between IPP and the approximate location of station; the scalar field r s r represents the stochastic component from a second-order stationary process; and I(z) s r is the vertical ionospheric delay correction interpolated from global ionospheric map (GIM) with corresponding noise ε I(z) s r .
In the data processing, all GPS observations have a sampling interval of 30 s, and the satellite elevation mask angle is set to 10 • . GPS precise satellite orbit and clock products provided by IGS are adopted for PPP processing. The data processing models and strategies for GPS DF and SF PPP are summarized in Table 1. Note that in this study the experiments and results of DF PPP are introduced first and those of SF PPP are shown next since the DF receivers are more common used in real surveying applications. In addition, the following experiments only focus on the kinematic PPP solutions since geomagnetic storms have no obvious influence on static PPP with a daily processing.

Results and Analyses
The variations of solar wind velocity (Vsw), IMF Bz component, Kp index, and Dst index during three typical storms are shown in Figure 3. For comparison, the observations on the quiet day of 17 March 2017 are also given in the figure. Obviously, the time series of Vsw and Bz show a smooth variation on 17 March 2017 due to the quiet solar activity. For the moderate storm, the Vsw is found to gradually increase on 27 March 2017, and the Bz shows rapidly variations but the magnitude is small within ±8 nT. The geomagnetic index Dst shows a decrease until 14:00 coordinated universal time (UTC) when arriving at the minimum of −74 nT. For the intense and super storms, the Vsw shows quick increases at the beginning of the storm sudden commencement (SSC) time. The IMF Bz components for the two storms are also found being rapid fluctuations when the solar winds are very dynamic. Most importantly, the duration time of negative Bz for the intense and super storms can reach several hours, which are much longer than that of moderate storm; moreover, the larger magnitude of negative Bz around −15 nT can be clearly seen in those two great storms. That is why those storms are much stronger than the moderate storm of 27 March 2017. From the Kp index panels, we can find that the values of Kp in the main phase of intense and super storms are in the range of 6 to 8, which are generally larger than those of moderate storm as the range of 4 to 6. EIA can reach around ±22° MLAT, spreading toward mid-latitude. In addition, their TEC during 18:00-22:00 UTC can exceed 40 TECU. For the super storm on 17 March 2015, TEC values of more than 40 TECU are easily observed around ±30° MLAT on the whole day. In Figure 4, we can find that with the increase of geomagnetic storms intensity, the EIA also expands far from their regular position of ±15° MLAT and the TEC values also increase to a larger magnitude as 40-50 TECU.   To analyze the ionospheric response to geomagnetic storms with different intensities, the corresponding magnetic latitude (MLAT) versus UTC of TEC maps are presented in Figure 4. The TEC measurements are averaged over 5 min in UTC and 2.5 • MLAT over all longitudes. Note that no data of TEC map during 23:17-24:00 UTC on 17 and 27 March 2017, so the corresponding regions are blank in Figure 4a,b. On the quiet day, large values of TEC are found in low-latitude, and the peaks are found in the equatorial ionization anomaly (EIA) regions around ±15 • MLAT during 18:00-23:17 UTC. During the moderate geomagnetic storm, the crests of EIA further extend to about ±17 • MLAT and the larger magnitude TEC are observed during 15:00-23:17 UTC when the western hemisphere was on dayside, which is the same as on the quiet day. It is easy to understand that the daytime TEC is larger than that on the nightside. During the intense geomagnetic disturbed condition, the crests of EIA can reach around ±22 • MLAT, spreading toward mid-latitude. In addition, their TEC during 18:00-22:00 UTC can exceed 40 TECU. For the super storm on 17 March 2015, TEC values of more than 40 TECU are easily observed around ±30 • MLAT on the whole day. In Figure 4, we can find that with the increase of geomagnetic storms intensity, the EIA also expands far from their regular position of ±15 • MLAT and the TEC values also increase to a larger magnitude as 40-50 TECU.  Table 2. In Figure 5 and Table 2, we can find that stronger intensity geomagnetic storms generally cause more significant degradation of GPS PPP.    Table 2. In Figure 5 and Table 2, we can find that stronger intensity geomagnetic storms generally cause more significant degradation of GPS PPP.
From the ROTI panels, it is observed that under geomagnetic storms conditions the time series of ROTI fluctuate significantly during the main phase as well as the recovery phase period, and the values can increase to around 1.0 TECU/min. Conversely, most of them in the quiet ionospheric condition show smooth variation, being lower than 0.05 TECU/min. Generally, ROTI > 0.5 TECU/min indicates the presence of ionospheric irregularities at scale length of a few kilometers [27]. Note that there are few abnormal ROTIs being larger than 0.3 TECU/min under quiet condition. The abnormal results can be caused by many factors such as multipath. When satellite signals pass through the ionospheric irregularities, they would encounter signal attenuation, even sometimes complete loss of lock [28]. Hence, the degradation of PPP during the geomagnetic storm period should be attributed to ionospheric irregularities. The detailed relationship between ionospheric irregularities and degradation of PPP over the global scale are shown below.  From the ROTI panels, it is observed that under geomagnetic storms conditions the time series of ROTI fluctuate significantly during the main phase as well as the recovery phase period, and the values can increase to around 1.0 TECU/min. Conversely, most of them in the quiet ionospheric condition show smooth variation, being lower than 0.05 TECU/min. Generally, ROTI > 0.5 TECU/min indicates the presence of ionospheric irregularities at scale length of a few kilometers [27]. Note that there are few abnormal ROTIs being larger than 0.3 TECU/min under quiet condition. The abnormal results can be caused by many factors such as multipath. When satellite signals pass through the ionospheric irregularities, they would encounter signal attenuation, even sometimes complete loss of lock [28]. Hence, the degradation of PPP during the geomagnetic storm period should be attributed to ionospheric irregularities. The detailed relationship between ionospheric irregularities and degradation of PPP over the global scale are shown below.   Table 3. During geomagnetic storms, it can be clearly seen that more serious degradation of DF PPP is found in the high-latitude region compared with other two latitude regions. The degradations of mid-and low-latitude regions are comparable. For high-latitude region, the 3D RMSs during moderate, intense, and super storms are 0.393 m, 0.680 m, and 1.051 m, respectively, while under quiet condition is only 0.163 m (see Table 3). The performances of DF PPP for mid-and low-latitude regions are not very influenced by the moderate and intense storms, but few stations show major degradation during the super storm period. some stations located in the high-latitude, e.g. North America region, experience the positioning quality deterioration as their RMS can reach around 0.5 m. With the increase of storm intensity, more and more stations are affected by the ionospheric disturbances, which can be obviously seen in the third and fourth panels. In the vertical component, stations with RMS > 0.5 m account for 1.7%, 7.0%, and 14.6% during moderate, intense, and super storms period, respectively. In addition, there are 19 stations with RMS > 1.0 m during the super storm of 17 March 2015.  The corresponding results of GPS SF PPP are shown in Figures 8 and 9. Similar to DF PPP, the degree of deterioration for SF PPP mainly depends on the intensity of geomagnetic storm. Due to the lack of L2 observations, the performance of SF PPP is inferior to that of DF PPP regardless of quiet or disturbed condition. Under quiet condition, the 3D RMS of SF PPP are 0.341 m, 0.361 m, and 0.456 m in high-, mid-, and low-latitude regions (see Table 3), which is consistent with the previous study [19]. During geomagnetic storms periods, it can be found that stations experiencing degraded positioning of SF PPP are mainly distributed in high-latitude region. Their 3D RMSs during moderate, intense, and super storms are 0.963 m, 0.911 m, and 1.624 m (see Table 3). Note that, during a super storm period, the degradation of SF PPP can also be clearly seen in mid-and low-latitude regions. Nava et al. [29] pointed out that, during the super storm period on 17 March 2015, the extreme ionospheric conditions have been observed in mid-and low-latitudes. Therefore, the deteriorated ionospheric correction caused by large gradients of ionosphere under super storm condition should be treated as the main reason for the degradation of positioning in mid-and low-latitude regions. It is easy to understand that, under severe ionospheric condition, ionospheric correction is a more serious challenge in SF solution compared with that in DF solution. The 90% confidence level of time series of positioning errors of DF PPP over all stations located at high-latitude (60°-90° N/S), mid-latitude (30°-60° N/S), and low-latitude (0°-30° N/S) regions during three typical storms are presented in Figure 7. The background gray scatters are the positioning results for all IGS stations.For clarity, the first two hours results in the state of convergence are excluded in Figure 7. The corresponding RMS statistics in the north, east, and up components are also given in the figure. The three-dimensional (3D) RMS statistics for all stations are presented in Table 3. During geomagnetic storms, it can be clearly seen that more serious degradation of DF PPP is found in the high-latitude region compared with other two latitude regions. The degradations of mid-and low-latitude regions are comparable. For high-latitude region, the 3D RMSs during moderate, intense, and super storms are 0.393 m, 0.680 m, and 1.051 m, respectively, while under quiet condition is only 0.163 m (see Table 3). The performances of DF PPP for mid-and low-latitude regions are not very influenced by the moderate and intense storms, but few stations show major degradation during the super storm period.   From the above analysis, it is found that the degraded performance of GPS DF and SF PPP should be attributed to the complex ionospheric conditions caused by geomagnetic storms. Therefore, it is necessary to show the specific ROTI data to further analyze the ionospheric disturbances during different intensities storms. Figure 10 presents the global ROTI maps derived from entire IGS network in the three perturbation periods (i.e., initial phase, storm development, and deep main phase) of each geomagnetic storm. The ROTI values are averaged in grid of 2.5 • latitude and 5 • longitude with a sliding window for each 5 min interval. As the GPS phase fluctuation index, ROTI can be used to monitor global activity of ionospheric irregularities [27,30]. From the topside panels, we can clearly observe that most ROTI values are smaller than 0.3 TECU/min, indicating no obvious irregularity occurrences on the quiet day of 17 March 2017. Consequently, the corresponding results of GPS PPP are in agreement with the expectations. In the case of the geomagnetic storms, the ROTI values would be different. During the initial phase period, as shown in the left three panels, the irregularities appeared at high-latitude of northern America sector as well as Antarctic sector where some ROTI values are larger than 0.5 TECU/min. Since this period is in the initial phase of the storms, the irregularities are not commonly seen in other sectors. Afterwards, the storms went into the main phase. The significant irregularities can be easily seen in high-latitude regions including North America, northern Europe as well as Antarctic region, and some ROTI values even can reach 1 TECU/min. Note that the stronger intensity of geomagnetic storms, the higher ROTI values and larger scale of irregularities are found in general. Consequently, more stations are affected by the irregularities (see Figures 6 and 8). In addition, during the main phase period around 18:00 UTC of the super storm, the irregularities can also be found in the equatorial regions such as the India sector. During the deep main phase close to recovery phase, shown in the right three panels, the intensity of irregularities decreases gradually although the auroral activities were also active especially in the Antarctic and South America regions during the super storm. The corresponding results of GPS SF PPP are shown in Figures 8 and 9. Similar to DF PPP, the degree of deterioration for SF PPP mainly depends on the intensity of geomagnetic storm. Due to the lack of L2 observations, the performance of SF PPP is inferior to that of DF PPP regardless of quiet or disturbed condition. Under quiet condition, the 3D RMS of SF PPP are 0.341 m, 0.361 m, and 0.456 m in high-, mid-, and low-latitude regions (see Table 3), which is consistent with the previous study [19]. During geomagnetic storms periods, it can be found that stations experiencing degraded positioning of SF PPP are mainly distributed in high-latitude region. Their 3D RMSs during moderate, intense, and super storms are 0.963 m, 0.911 m, and 1.624 m (see Table 3). Note that, during a super storm period, the degradation of SF PPP can also be clearly seen in mid-and low-latitude regions. Nava et al. [29] pointed out that, during the super storm period on 17 March 2015, the extreme ionospheric conditions have been observed in mid-and low-latitudes. Therefore, the deteriorated ionospheric correction caused by large gradients of ionosphere under super storm condition should be treated as the main reason for the degradation of positioning in mid-and low-latitude regions. It is easy to understand that, under severe ionospheric condition, ionospheric correction is a more serious challenge in SF solution compared with that in DF solution.   From the above analysis, it is found that the degraded performance of GPS DF and SF PPP should be attributed to the complex ionospheric conditions caused by geomagnetic storms. Therefore, it is necessary to show the specific ROTI data to further analyze the ionospheric disturbances during different intensities storms. Figure 10 presents the global ROTI maps derived from entire IGS network in the three perturbation periods (i.e., initial phase, storm development, and deep main phase) of each geomagnetic storm. The ROTI values are averaged in grid of 2.5° latitude and 5° longitude with a sliding window for each 5 min interval. As the GPS phase fluctuation index, ROTI can be used to monitor global activity of ionospheric irregularities [27,30]. From the topside panels, we can clearly observe that most ROTI values are smaller than 0.3 TECU/min, indicating no obvious irregularity occurrences on the quiet day of 17 March 2017. Consequently, the corresponding results of GPS PPP are in agreement with the expectations. In the case of the geomagnetic storms, the ROTI values would be different. During the initial phase period, as shown in the left three panels, the irregularities appeared at high-latitude of northern America sector as well as Antarctic sector where some ROTI values are larger than 0.5 TECU/min. Since this period is in the initial phase of the storms, the irregularities are not commonly seen in other sectors. Afterwards, the storms went into the main phase. The significant irregularities can be easily seen in high-latitude regions including North America, northern Europe as well as Antarctic region, and some ROTI values even can reach 1 TECU/min. Note that the stronger intensity of geomagnetic storms, the higher ROTI values and larger scale of irregularities are found in general. Consequently, more stations are affected by the irregularities (see Figures 6 and 8). In addition, during the main phase period around 18:00 UTC of the super storm, the irregularities can also be found in the As discussed before, the ionospheric irregularities caused by geomagnetic storms possibly lead to the occurrence of cycle-slip and even loss of signal tracking [31]. The frequent cycle-slip occurrence can result in degradation of GNSS positioning. To further analyze the characteristics of cycle-slip occurrence under different intensities storms, Figure 11 gives the scatters of cycle-slip occurrence rate against geographical latitude derived from all IGS stations data during the quiet, moderate, intense, and super storm periods. In this study, the Melbourne-Wübbena wide-land combination and rate of geometry-free combination proposed by Luo et al. [9] are jointly used to detect cycle-slip. Note that cycle-slip occurrence rate is calculated as the total number of cycle-slips divided by total available measurements. equatorial regions such as the India sector. During the deep main phase close to recovery phase, shown in the right three panels, the intensity of irregularities decreases gradually although the auroral activities were also active especially in the Antarctic and South America regions during the super storm. As discussed before, the ionospheric irregularities caused by geomagnetic storms possibly lead to the occurrence of cycle-slip and even loss of signal tracking [31]. The frequent cycle-slip occurrence can result in degradation of GNSS positioning. To further analyze the characteristics of cycle-slip occurrence under different intensities storms, Figure 11 gives the scatters of cycle-slip occurrence rate against geographical latitude derived from all IGS stations data during the quiet, moderate, intense, and super storm periods. In this study, the Melbourne-Wübbena wide-land combination and rate of geometry-free combination proposed by Luo et al. [9] are jointly used to detect cycle-slip. Note that cycle-slip occurrence rate is calculated as the total number of cycle-slips divided by total available measurements.
From the topside two panels, we can see that most scatters are within the range of 0-0.5%. Compared with the South Pole, more cycle-slips are detected in the North Pole, which is consistent with the results reported by Astafyeva et al. [32]. From the bottom two panels, we can find that the scatters are more discrete mainly ranging from 0% to 1%. That means many cycle-slips occurred From the topside two panels, we can see that most scatters are within the range of 0-0.5%. Compared with the South Pole, more cycle-slips are detected in the North Pole, which is consistent with the results reported by Astafyeva et al. [32]. From the bottom two panels, we can find that the scatters are more discrete mainly ranging from 0% to 1%. That means many cycle-slips occurred under intense and super storms conditions. Meanwhile, although there are fewer stations located in the South Pole, the cycle-slip occurrence rate can reach around 2% under super storm condition. That implies larger scale irregularities induced by super storm resulted in cycle-slip occurrence in the two poles. The statistics show that the mean values of cycle-slip occurrence rate are 0.23%, 0.26%, 0.36%, and 0.39% for quiet, moderate, intense, and super storm days, respectively. Compared with quiet days, the increased cycle-slips occurrence percentages are 13.04%, 56.52%, and 69.57% during moderate, intense, and super storm days. From the above analysis, we conclude that the degraded performance of GPS PPP under geomagnetic storm conditions should be attributed to two factors: (1) the deteriorated ionospheric correction under severe ionospheric conditions; and (2) the contaminated satellite measurements, such as serious signal attenuation and frequent cycle-slip occurrence, when satellites signals pass through the ionospheric irregularities. quiet days, the increased cycle-slips occurrence percentages are 13.04%, 56.52%, and 69.57% during moderate, intense, and super storm days. From the above analysis, we conclude that the degraded performance of GPS PPP under geomagnetic storm conditions should be attributed to two factors: (1) the deteriorated ionospheric correction under severe ionospheric conditions; and (2) the contaminated satellite measurements, such as serious signal attenuation and frequent cycle-slip occurrence, when satellites signals pass through the ionospheric irregularities.

Discussion
This study presents a comprehensive investigation of the performance of GPS DF and SF PPP under three typical geomagnetic storms, using a large data set collected from about 500 IGS stations, Under geomagnetic storms conditions, the time series of positioning results derived from high-latitude stations data for GPS DF and SF PPP vary significantly and the range can reach several meters in the horizontal and vertical components. This is consistent with the previous investigations based on RTK [13][14][15][16] and PPP techniques [17]. With the increase of storm intensity, more and more stations experienced the degradation of positioning with most located at high-latitude regions. It is well known that the high-latitude region has strong coupling among the interplanetary medium, the magnetosphere, the thermosphere, and the ionosphere [29]. During geomagnetic storms, abundant energetic particles precipitate into high-latitude region, and heat the upper atmosphere in the form of Joule heating and Lorentz forces, resulting in the steep ionospheric density gradients and irregularities [5,33]. Our global ROTI maps clearly show that the ionospheric irregularities are mainly distributed in high-latitude regions including North America, northern Europe and Antarctic region during the main phase period of geomagnetic storms (see Figure 10).

Discussion
This study presents a comprehensive investigation of the performance of GPS DF and SF PPP under three typical geomagnetic storms, using a large data set collected from about 500 IGS stations, Under geomagnetic storms conditions, the time series of positioning results derived from high-latitude stations data for GPS DF and SF PPP vary significantly and the range can reach several meters in the horizontal and vertical components. This is consistent with the previous investigations based on RTK [13][14][15][16] and PPP techniques [17]. With the increase of storm intensity, more and more stations experienced the degradation of positioning with most located at high-latitude regions. It is well known that the high-latitude region has strong coupling among the interplanetary medium, the magnetosphere, the thermosphere, and the ionosphere [29]. During geomagnetic storms, abundant energetic particles precipitate into high-latitude region, and heat the upper atmosphere in the form of Joule heating and Lorentz forces, resulting in the steep ionospheric density gradients and irregularities [5,33]. Our global ROTI maps clearly show that the ionospheric irregularities are mainly distributed in high-latitude regions including North America, northern Europe and Antarctic region during the main phase period of geomagnetic storms (see Figure 10). As a result, most stations located at high-latitude region experienced the degradation of positioning.
To GNSS users, adverse effects caused by the geomagnetic storm should be carefully mitigated to ensure the reliable positioning requirement. The mitigation strategies can be summarized as three parts: improvement of the success rate of cycle-slip detection and correction during the data processing [34]; consideration of the second-or high-order ionospheric delays in ionospheric corrections [15,35]; and combination of other GNSS system satellites to improve the weak satellite geometries [36]. The frequent cycle-slip occurrence during geomagnetic storms has been demonstrated in Figure 11. Since most approaches of cycle-slip detection and correction lack proper handling of ionospheric delay variations, a geometry-based approach with rigorous handling of the ionosphere was presented by Banville and Langley [34]. After applying this method to GNSS data collected in northern Canada, some kinematic PPP results showed a dm-level improvement. However, the increased measurement noise associated with an active ionosphere still makes cycle-slip detection and correction difficulty. Under geomagnetic storm conditions, the first-order ionospheric delays can reach about 30 m and the second-order delays are very small but with about 3 mm magnitude [35,37]. Therefore, consideration of high-order ionospheric delays can only contribute to an improvement of mm-level in GNSS PPP. The integration of multi-constellation GNSS can significantly increase the number of satellites, thereby increasing the measurement redundancy and improving the satellite geometry. However, the positioning errors of GPS/GLONASS PPP in the height component can also reach around 1 m during the geomagnetic storm period on 13 May 2015 [36]. From the above analysis, we can conclude that aforementioned strategies may also be influenced by severe geomagnetic storm conditions. More reliable mitigation strategies require further investigation.

Conclusions
A comprehensive investigation of geomagnetic storms effects on GPS PPP has been presented in our study. Results show that under normal condition without geomagnetic storms occurrence GPS DF PPP accuracy of most IGS stations is better than 0.15 m and 0.2 m in the horizontal and vertical components, respectively. The degradation performance of DF PPP can be easily found in geomagnetic storms conditions especially at high-latitude region. The 3D RMS at this region based on all stations' results during moderate, intense, and super storms are 0.393 m, 0.680 m, and 1.051 m, respectively, while that value is only 0.163 m under quiet condition. In addition, our results show that stronger intensity storms generally cause more significant degradation of positioning. The statistics indicate that the number of stations with RMS > 0.5 m accounts for 1.7%, 7.0%, and 14.6% during moderate, intense, and super storm periods, respectively. Besides, there are 19 stations with RMS > 1.0 m during the super storm of 17 March 2015. Compared with DF PPP, the performance of GPS SF PPP is inferior regardless of quiet or disturbed condition. Under quiet condition, GPS SF PPP accuracy of most IGS stations is around 0.35 m and 0.6 m in the horizontal and vertical components. Under disturbed condition, the 3D RMS of SF PPP during moderate, intense, and super storms can reach 0.963 m, 0.910 m, and 1.624 m at high-latitude region, respectively. Different from DF PPP, more serious degradation of SF PPP can also be clearly found at mid-and low-latitude regions, especially for the super storm condition.
During geomagnetic storm periods, the degraded performance of GPS positioning regardless of DF or SF PPP should be attributed to the increased ionosphere disturbances. The strong response between positioning errors and ionospheric irregularities index ROTI has been confirmed in this study based on single station data collected at high-latitude. Moreover, the global ROTI maps further indicated the strong correlation between the degradation positioning and ionospheric irregularities. The ionospheric disturbances can not only lead to the deteriorated ionospheric correction but also to frequent cycle-slip occurrence. Statistical results show that the increased cycle-slip occurrence percentages are 13.04%, 56.52%, and 69.57% under moderate, intense, and super storms conditions, respectively, compared with those of quiet day. Note that this study only focuses on GPS PPP performance analysis under different geomagnetic storm conditions. The corresponding multi-GNSS PPP analysis requires further investigation.
Author Contributions: X.L. and S.G. designed this study, performed the analysis, and wrote the manuscript. Y.L., C.X., B.C. and X.J. provided guidance and helped polish the manuscript. All authors reviewed the manuscript.