Positioning performance of the NTCM model driven by GPS Klobuchar model parameters

This is anOp Abstract – Users of the Global Positioning System (GPS) utilize the Ionospheric Correction Algorithm (ICA) also known as Klobuchar model for correcting ionospheric signal delay or range error. Recently, we developed an ionosphere correction algorithm called NTCM-Klobpar model for single frequency GNSS applications. The model is driven by a parameter computed from GPS Klobuchar model and consecutively can be used instead of the GPS Klobuchar model for ionospheric corrections. In the presented work we compare the positioning solutions obtained using NTCM-Klobpar with those using the Klobuchar model. Our investigation using worldwide ground GPS data from a quiet and a perturbed ionospheric and geomagnetic activity period of 17 days each shows that the 24-hour prediction performance of the NTCMKlobpar is better than the GPS Klobuchar model in global average. The root mean squared deviation of the 3D position errors are found to be about 0.24 and 0.45m less for the NTCM-Klobpar compared to the GPS Klobuchar model during quiet and perturbed condition, respectively. The presented algorithm has the potential to continuously improve the accuracy of GPS single frequency mass market devices with only little software modification.


Introduction
The ionospheric delay is considered as one of the biggest errors for single frequency use of space based Global Navigation and Satellite System (GNSS).At GNSS operating frequencies in the order of 1-2 GHz the ionospheric delay may cause link related range errors of up to 100 m.Thus, GNSS single frequency operations need ionospheric delay information for mitigating ionospheric propagation errors.The ionospheric propagation delay is inversely proportional to the square of the signal frequency and directly proportional to the integral of the electron density along the ray path called total electron content (TEC) which is commonly expressed in units of TEC termed TECU (1 TECU = 10 16 electrons/m 2 ).Therefore, single frequency GNSS positioning needs either TEC or equivalent ionospheric delay information for mitigating ionospheric propagation errors.
Global Positioning System (GPS) utilizes the Ionospheric Correction Algorithm (ICA), also known as Klobuchar model, for correcting ionospheric signal delay or range error (Klobuchar, 1987;IS-GPS-200G, 2012).In order to do this, GPS transmits 8 ionospheric correction coefficients in the navigation message on a daily basis as driving parameter set for the Klobuchar model.The Klobuchar model gives a representation of the mean vertical delay at the GPS L1 frequency as a half-cosine function with varying amplitude and period.The peak of the cosine is fixed at 14 hour local time (LT) and during night-time hours the vertical ionospheric delay is fixed at a constant value of 5 ns or 9.24 TECU.The amplitude and period of the cosine function are modelled by 3rd order polynomials whose coefficients are broadcasted in the GPS navigation message.
The European satellite navigation system Galileo uses a three dimensional time dependent ionospheric electron density model called NeQuick (Nava et al., 2008) for single frequency ionospheric correction.The original NeQuick model uses a monthly averaged solar radio flux index F10.7 as a proxy measure of the solar activity level.However, to achieve higher accuracy, the Galileo ionospheric correction model uses an effective ionization level called Az as a primary input parameter for the NeQuick model.The Az approach whose polynomial coefficients are derived from dual frequency measurements at selected ground stations takes implicitly into account the daily variation of the solar activity and the user's local geomagnetic conditions.The Galileo satellites broadcast Az coefficients via the navigation message for computing Az at user level all over the globe.
In our former study (Jakowski et al., 2011a) we developed an empirical global TEC model called Neustrelitz TEC Model (NTCM) for estimating trans-ionospheric radio wave propagation errors.The NTCM approach explicitly describes the TEC dependencies on local time, geographic/geomagnetic location and solar irradiance and activity using only 12 model coefficients.The model coefficients were computed by nonlinear fitting of global TEC data covering one decade, i. e. about one solar cycle to the model approach in least squares sense.So we used TEC data obtained from Center for Orbit Determination System (CODE) at the University of Bern (http://aiuws.unibe.ch/ionosphere/)for the years 1998-2007.The only driving parameter of the model is the daily solar radio flux proxy F10.7.Our subsequent study (Hoque et al., 2015) shows that NTCM model can successfully be used in GNSS applications for ionospheric correction estimation.
In a recent study (Hoque & Jakowski, 2015) we simplified the original NTCM model in order to use it as a broadcast ionosphere model for future satellite navigation systems and called it NTCM-BC.In contrast to the original NTCM, the NTCM-BC coefficients are not constant because they are adapted to the current ionospheric conditions on a daily base.Being valid for a period of typically 24 h worldwide, NTCM-BC reacts to global ionospheric dynamics and therefore can achieve a higher accuracy than NTCM which uses fixed coefficients.
In a more recent study (Hoque et al., 2017) we recomputed the NTCM model coefficients using more recent (2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014)(2015) TEC data from CODE.Additionally, we substituted the daily solar radio flux index F10.7 that drives NTCM by providing a proxy of the solar activity level, by a new quantity called Klobpar.The Klobpar is computed from the GPS Klobuchar model driven by the broadcasted coefficients as pointed out in section 2 in detail.Since F10.7 is not available to the GPS receiver without additional data link, the user cannot use F10.7-driven-NTCM model in operational purposes.Keeping this in mind we proposed a new NTCM model driven by a parameter computed from the GPS Klobuchar model.Therefore, GPS users can easily compute the value of Klobpar using broadcasted coefficients and then use it as a primary input parameter for the NTCM-Klobpar model.Our research, using post processed reference TEC data from more than one solar cycle, showed that on average the RMS modelled TEC errors are up to 40% less for the proposed NTCM model compared to the Klobuchar model during high solar activity period, and about 10% less during low solar activity period (Hoque et al., 2017).
In the work presented here, we evaluate for the first time the performance of the NTCM-Klobpar model in the position domain by comparing position estimates for selected quiet and perturbed ionospheric conditions with those obtained by using the mother Klobuchar model.Additionally, we compare position solution estimates for original F10.7 driven NTCM and NTCM-BC models.We calculated numerous test user positions using a standard Single Point Positioning (SPP) approach (IS-GPS-200F, 2012) in which the ionospheric correction is provided either by the GPS Klobuchar or NTCM-Klobpar or NTCM or NTCM-BC.The actual 3D user positions are known and subtracted from each SPP solution to obtain associated position errors.Then, the models are compared in terms of estimated position errors.In addition, the actual 3D positions are compared with ionosphere uncorrected SPP solutions.

Data sources and data processing
To compare the different models we selected a period of 17 days from 15-31 January, 2011 (day of year or DOY 15-31) during which the solar activity level proxy F10.7 lies between 77-84 solar flux units (sfu, 1 sfu = 10 À22 Wm À2 Hz À1 ) and the geomagnetic activity proxy Kp does not exceed the 3.0 level (see Fig. 1a) indicating quiet ionospheric and geomagnetic conditions.Again, we selected a period of 17 days from 23rd October-8th November, 2011 (DOY 296-312) during which F10.7 lies between 121 and 178 sfu and Kp exceeds the 4.0 threshold for two times (see Fig. 1a) indicating ionospheric perturbed conditions.
The daily GPS Klobuchar coefficients which are driving parameters for the Klobuchar model were downloaded from the NOAA's (National Oceanic and Atmospheric Administration) National Geodetic Survey (NGS) archive ftp://www.ngs.noaa.gov/cors/rinex/The only driving parameter of NTCM which is the daily F10.7 is obtained from the Space Physics Interactive Data Resource (SPIDR).The NTCM-BC coefficients are computed each day using dual-frequency GPS data from globally distributed International GNSS Service (IGS, Dow et al., 2009) stations downloaded from the NASA's Crustal Dynamics Data Information System (CDDIS) archive (Noll, 2010; ftp://cddis.gsfc.nasa.gov/gnss/data/).The procedures used for inter-frequency satellite and receiver bias estimation and TEC calibration, NTCM-BC optimization are discussed by Jakowski et al. (2011a), Hoque & Jakowski (2015) and Hoque et al. (2015).
The only driving parameter of NTCM-Klobpar is the Klobpar which is defined as the sum of TEC obtained by the GPS Klobuchar model at two specific geodetic locations A (latitude f = 10 °N, longitude l = 90 °W) and B (f = 10 °S, l = 90 °W) at 14 Universal Time (UT).The GPS Klobuchar model is used to compute the GPS L1 ionospheric delay T iono at points A and B in seconds.The vertical delays are then converted to TEC to compute the daily Klobpar parameter in TECU as where speed of light c = 299792458 m/s and factor ¼ 6:1587 (for details we refer to Hoque et al., 2017).
As discussed in Hoque & Jakowski (2015), the NTCM-BC can be optimized on a daily basis in the GNSS control segment using TEC data of the previous day from monitor stations like the GPS Klobuchar or Galileo NeQuick.To accomplish the model comparisons as close as possible to conditions applicable to an operational GNSS, a set of globally distributed IGS stations is selected as sensor/monitor stations for NTCM-BC optimization and a set of test stations is chosen from remaining IGS stations for the analysis of modelling errors.In the present study, we followed the same procedure.The number of monitor and test stations selected each day during quiet and perturbed period is plotted in Figure 1b.The selected test stations are used to perform a global analysis of model performance in position solution estimates.The considered models are Klobuchar, NTCM, NTCM-BC and NTCM-Klobpar.
As an example, on 16th January (DOY 16) and 24th October (DOY 297) 2011, we selected 35 and 31 test stations worldwide, respectively.The location of test stations (marked with red crosses) as well as monitor stations (marked with green dots) on global map is shown in Figure 2a, b for DOY 16 and DOY 297, respectively.We used the GPS P1 pseudorange measurements from RINEX (Receiver Independent Exchange) observation files for SPP solution computation.The RINEX observation files from worldwide IGS stations are downloaded from the archive ftp://cddis.gsfc.nasa.gov/gnss/data/for the selected quiet and perturbed period.As reference user position we used IGS weekly solutions of station coordinate given in the SINEX (Solution Independent Exchange) file.The SINEX data are downloaded from the same archive mentioned above.

Single-point positioning approach
We implemented a single-point positioning SPP algorithm based on DOD SPS ( 2008) and IS-GPS-200F (2012) to validate the performance of different ionospheric corrections.We used the GPS P1 pseudorange measurements recorded at 30-second interval as input to the SPP module.The pseudorange is a measure of the distance between the satellite and the receiver which is different from the geometric distance between them due to the errors of the both clocks and the influences of the signal propagation mediums.Taking the satellite and receiver clock errors dt s and dt r , the ionospheric effects d ion , tropospheric effects d tro into account, the pseudorange R s r can be written as where c is the speed of light, r s r is the geometric distance between the satellite and receiver, t e and t r are the GPS signal emission and reception time of the satellite and receiver, respectively.The term e represents errors due to effects that are not modelled such as multipath, receiver noise, etc.
During the signal transmission, the receiver rotates with the Earth, therefore the so-called Sagnac corrections need to be considered.At the equator, the Earth rotation effect is equivalent to about 31 meters position displacement (Xu, 2007).To mitigate the Sagnac effect we corrected the original satellite coordinates using the method described by Seeber (2003).The GPS satellite broadcasts navigation messages for satellite clock correction coefficients and orbit parameters.Using these parameters we corrected the satellite clock error as described in IS- GPS-200F (2012).The correction accounts for clock error characteristics of bias, drift and aging, as well as for the satellite implementation characteristics of group delay bias and mean differential group delay.

Tropospheric correction
The troposphere is a non-ionized medium extending from the Earth's surface up to about 10-13 km of altitude.The tropospheric delay on GNSS signal can be separated into a dry or hydrostatic and a wet component.The hydrostatic component in zenith direction (ZHD) amounts to about 2.3 m whereas zenith wet delay (ZWD) is only in the range of about 0.15 m in global average.The dry (hydrostatic) component contributes about 90% of the total tropospheric error while the wet component contributes about 10% of the error.Typical tropospheric correction models can calculate the zenith total delay (ZTD = ZHD þ ZWD) from measured, predicted or longtime mean values of pressure, temperature and humidity on the earth surface and derive the delay in other elevation angles by multiplying with so called mapping functions.
For tropospheric correction we used the ESA blind model (Martellucci & Blarzino, 2003) which is an extension to the approach suggested by the RTCA-MOPS (2001) model by using input parameters derived from Numerical Weather Predictions (NWP) spatial fields.The use of spatial fields produced by the European Centre for Meteorological Weather Forecast (ECMWF) permits to increase the spatial resolution and the temporal resolution of the meteorological parameters to be used as input to the propagation model of tropospheric delay (Galileo, 2004).The ZTD is calculated by summing up the ZWD and ZHD, and the Niell mapping function (Neill, 1996) is used for zenith to slant delay conversion.

Ionospheric correction
The variability of the Earth's ionosphere is much larger than that of the troposphere, and it is more difficult to model.The ionospheric range error can vary from a few meters to a few tens of meters at the zenith, whereas the tropospheric range error at the zenith is generally between two to three meters.The range error of the ionosphere frequently changes by at least one order of magnitude during the course of each day.
As already mentioned, the ionospheric correction is provided either by the GPS Klobuchar, NTCM, NTCM-BC or NTCM-Klobpar.All these models are 2 dimensional vertical delay models indicating that the ionospheric delay at real slant ray paths has to be derived by means of an appropriate mapping function.
The GPS Klobuchar model uses an elevation dependent mapping function to convert vertical to slant delay at user level.The used vertical delay corresponds at the ionospheric pierce point (IPP) location of the satellite-receiver ray path at the height of 350 km.In case of NTCM models, a thin-shell ionosphere mapping function is used (e.g., Jakowski et al., 2011b).As already mentioned, NTCM and NTCM-Klobpar model coefficients are computed based on CODE TEC maps which represent vertical delays at an IPP height of 450 km.Therefore, the vertical delay corresponds at the IPP height of 450 km is used for both versions of NTCM driven by F10.7 and Klobpar.However, NTCM-BC model coefficients are computed based on GPS data from worldwide IGS stations and in the ionosphere estimation procedure an IPP height of 400 km is used.So for the NTCM-BC we used the IPP height of 400 km for mapping function.
We computed receiver position by first linearizing the range observation equations and then using ordinary leastsquares principle.We estimated residuals, i.e., the difference between the actual observations and the new estimated model for the observations.Several iterations are done before obtaining the final solution.The basic approach is given in associated matlab routines published along with the paper by Borre (2003).

Evaluation of models and discussion
In order to analyze and discuss the ionospheric delay correction effectiveness of the models described in the previous section, the following analysis method has been applied to compare the SPP results obtained by each of the models.
The approximate 3D user positions (X, Y, Z coordinates in Earth Centered Earth Fixed system) are given in RINEX observation files.However, their accuracy is not known.Therefore, we looked for station coordinates provided by the International GNSS Service (IGS) community.We found that the IGS routinely generates a number of weekly station coordinates among other products by combining independent estimates from at least seven IGS Analysis Centers (ACs) and distributes in SINEX file format.The combined solutions are aimed to align to an IGS realization (IGS05) of the International Terrestrial Reference Frame (ITRF, 2005).A measure of the internal coordinate consistency is given by Ferland & Piraszewski (2009) by analyzing the residual standard deviations between the ACs and the IGS weekly combination solution.They found the station coordinate consistency during the GPS weeks 1400-1480 as about 2-3 mm for the horizontal components and about 7 mm for the vertical component which are representative of the coordinate's accuracy.The SINEX provided station coordinates are assumed as reference values and subtracted from each SPP solution to calculate positioning errors.For simplicity the weekly station coordinates are assumed to be constant during the given GPS week.
As already mentioned the reference station coordinates are taken from SINEX files.However, for few stations the SINEX station coordinates were not available.In such cases, we used approximate station coordinates from RINEX files.As for example in Figure 3 for "kir0", SINEX station coordinates were not available for the quiet period.However, during the perturbed period coordinates were available and we compared the results obtained using both SINEX and RINEX station coordinates and found very similar results.So for "kir0" station, we presented the results obtained using RINEX station coordinates during the quiet period.Similarly, for "mtwa" station, SINEX station coordinates were not available for both periods.However, we found similar performance results for stations "tah2", "suth" and "chat" when using SINEX and RINEX station coordinates and therefore, for "mtwa" station we presented the results obtained using RINEX station coordinates.In the title of each sub plot (see Figures 3-5) the used reference station coordinate source is mentioned as either "sinex" or "rinex".
Figure 3 shows significant improvement in the position solution when using NTCM-Klobpar instead of the mother Klobuchar model.During the ionospheric perturbed period the improvement is more evident (see right panel plots).We found that all three NTCM models show similar performances during both periods.We found that at kir0 station (see Fig. 3a) the hourly mean 3D position estimates obtained by the Klobuchar model are worse than those without using ionosphere correction model especially at early and evening hours during quiet period.At mdvj station (see Fig. 3e) we found similar trend during early hours.Left panel plots (b), (d) and (f) show that at high latitude stations the peak of diurnal variation is shifted to about 18 LT for Klobuchar model whereas the peak for No Corr case is at about 14 LT.The exact reason for this deviation is not known.However, the maximum of mismodelling at 18 LT may be related to the perturbations.During the perturbation period we see this behavior also at other latitudes and also in our NTCM models.It disappears in the analysis for the quiet period.
Figure 4 shows positioning results obtained at middle and low latitude stations during the selected quiet and perturbed periods.It is evident that the NTCM-Klobpar model performs better than the Klobuchar model at middle latitudes.As before, NTCM and NTCM-BC models perform very similar to the NTCM-Klobpar model.Like high latitudes, at middle and low latitudes the Klobuchar model performs worse than No Corr solution during night time hours especially during quiet period (see Figure 4a, c and e).This is may be due to the reason that during night-time hours the vertical ionospheric delay of the Klobuchar model is fixed at a constant value of 5 ns or 9.24 TECU.However, during high solar activity conditions, it seems that the period of the cosine function of the Klobuchar model is better modelled and we found better results compared to the No Corr solution.
Figure 5 compares positioning results at several southern hemisphere stations.We found that all models perform very similar; that means we don't have significant improvement using NTCM models over the Klobuchar model.One reason may be that the lack of sufficient number of IGS stations in southern hemisphere which are used in CODE TEC map generation as well as NTCM-BC coefficients updates.Since the CODE TEC maps from previous and current solar cycles are used in NTCM and NTCM-Klobpar coefficients generation (Jakowski et al., 2011a;Hoque et al., 2017), they inherit inaccuracy due to disperse data distribution in the southern hemisphere from CODE TEC maps.
For global analysis of model comparisons, we computed the mean, standard deviation (STD), Root Mean Squared (RMS) of 3D position errors as well as 65 and 95 percentiles of position errors.A percentile is a measure of the position error below which a given percentage of observations in the data set fall.The statistical estimates are computed over the data from test stations during quiet and perturbed periods separately.However, to obtain statistically representative results, the data are arranged into the following groups: i) geographic latitude range 0 ≥ f ≥ |90°| and local time LT range 0-24; ii) 0 ≥ f ≥ |90°| and LT 6-18 h; iii) 0 ≥ f ≥ |90°| and LT 0-6 and 18-24 h; iv) 0 ≥ f ≥ |30°| and LT 0-24 h; v) 30>f ≥ |60°| and LT 0-24 h; vi) 60>f ≥ |90°| and LT 0-24 h.
We performed statistical analysis for each case and the results are given in Table 1 and 2  ii) daytime condition: we see similar performances as mentioned in case i; iii) nighttime condition: we see that for all models the values are less compared to daytime condition (i.e., case for both quiet and perturbed periods; iv) low latitude region: The RMS error is the same about 3.1 m during quiet period whereas about 0.4 m less for the NTCM-Klobpar compared to the Klobuchar model during perturbed period.We noted that the 65 percentile error is slightly higher for NTCM-Klobpar during perturbed period than those for Klobuchar model.Bar plots in Figure 6a-d show the mean, STD, RMS, 65 and 95 percentile of 3D position errors considering all samples (case i) during quiet and perturbed period separately.We found that NTCM-Klobpar error estimates are less compared to those of the Klobuchar model.It is also evident that the F10.7driven-NTCM and NTCM-BC also perform better than the Klobuchar model.Among all three NTCM versions, the NTCM-BC performs the best.

Conclusion
Summarizing the evaluation results, it can be stated that all three different NTCM model approaches NTCM, NTCM-BC and NTCM-Klobpar achieve a better performance in the position domain than the Klobuchar model regularly provided by GPS.Focusing on estimating the accuracy of the NTCM-Klobpar approach we found an improvement in the order of 0.24 m and 0.45 m in global average for unperturbed low solar activity and perturbed medium solar activity conditions, respectively.It can be expected that the improvement is even more pronounced at high solar conditions and during ionospheric perturbations.This will be shown in further studies that include severe storm and high solar activity conditions.Since the NTCM-Klobpar approach uses the Klobuchar coefficients, regularly provided in the GPS navigation message, for estimating the current solar activity level, this model could further improve single frequency positioning performed by mass market devices.To reach this goal NTCM-Klobpar must be implemented in mass market GPS receivers.Since the model approach is very compact, the required technology modification is rather easy to handle.

Fig. 1 .
Fig. 1.(a) F10.7 and Kp variation during selected quiet and perturbed days and (b) Number of monitor and test stations used in NTCM-BC optimization and all models validation, respectively.

Fig. 2 .
Fig. 2. Location of monitor stations marked as green dots and test stations marked as red crosses.The red line indicates the geomagnetic equator whereas magenta lines at ± 20 of the red line bound equatorial anomaly regions.

Fig. 4 .
Fig. 4. Hourly mean 3D position error at middle and low latitude stations wtzz, rabt and kokb during quiet (left panel-plots a, c, e) and perturbed (right panel-plots b, d, f) period.
However, the RMS, mean as well as STD and 95 percentile errors are less for the NTCM-Klobpar.v) mid latitude region: the RMS error is about 0.3 and 0.5 m less for the NTCM-Klobpar compared to the Klobuchar model during quiet and perturbed condition, respectively.The mean, STD, 65 and 95 percentile values are less for the NTCM-Klobpar model.We found that NTCM and NTCM-BC perform better than the Klobuchar model.vi) high latitude region: the RMS error is about 0.6 and 0.8 m less for the NTCM-Klobpar compared to the Klobuchar model during quiet and perturbed condition.The mean, 5. Hourly mean 3D position error at southern hemisphere stations mtwa, tah2, suth and chat during quiet (left panel-plots a, c, e, g) and perturbed (right panel-plots b, d, f, h) period.STD, 65 and 95 percentile values are significantly less for the NTCM-Klobpar model.NTCM and NTCM-BC perform better than the Klobuchar model during both periods.Like middle latitudes, NTCM and NTCM-BC perform better than Klobuchar model at high latitudes.
. The reference station coordinate values are taken from SINEX files.So we only used stations for which SINEX data were available.The number of test stations exceeds 30 for each day during quiet and perturbed periods.It should be noted the statistical estimates may change if another set of test stations is considered.

Table 1 .
Statistical estimates of 3D position errors for the Klobuchar model and NTCM driven by Klobpar during quiet and perturbed ionospheric conditions.

Table 2 .
Statistical estimates of 3D position errors for NTCM and NTCM-BC during quiet and perturbed ionospheric conditions.