TEC and Instrumental Bias Estimation of GAGAN Station Using Kalman Filter and SCORE Algorithm ()
Received 2 February 2016; accepted 22 February 2016; published 25 February 2016
1. Introduction
The Global Positioning System (GPS) is a satellite-based navigation system capable of providing three dimensional position, velocity and timing information to users anywhere on or above the surface of the earth. GPS has been in use for a wide variety of applications. These include during flight in oceanic routes, enroute over the domestic airspace, and in crowded metropolitan airspaces. Aircraft on the final approach to airports, demands the greatest safety and reliability. To use GPS for precision approach (PA) and landings of civil aviation, the navigation system has to meet the Required Navigation Performance (RNP) parameters. These include accuracy, integrity, availability and continuity of service. Standalone GPS does not meet these precision approach requirements. The horizontal and vertical accuracy required for Category-I PA is 16 m and 4 m (95%) respectively [1] .
As a part of the global Communications, Navigation, and Surveillance/Air Traffic Management (CNS/ATM) plan adopted by the International Civil Aviation Organisation (ICAO), the Airports Authority of India (AAI) and Indian Space Research Organisation (ISRO) have jointly undertaken a programme for implementation of the Indian Satellite Based Augmentation System (SBAS) known as GPS Aided Geo Augmented Navigation (GAGAN) [2] .
Figure 1 shows the GAGAN system architecture and various corrections broadcasted to users through a geostationary (GEO) satellite. The GAGAN system consists of a network of ground reference stations called Indian Reference Stations (INRES) to provide necessary augmentations to the GPS navigation signal. Each INRES is equipped with a high quality dual frequency GPS receiver. These stations are strategically positioned across the country to collect data from GPS and geostationary (GEO) satellites. The measurement data from each of the INRES are collected by the Indian Mission Control Centre (INMCC) in real time. The chief functions of the INMCC are: network management; integrity monitoring; estimation of ephemeris, clock and ionospheric corrections; and issue of necessary commands. These wide area corrections and integrity information are packed into a GAGAN message and sent to the ground earth station known as the Indian Land Uplink Station (INLUS). The INLUS uplinks the GAGAN message to the geostationary satellite using C-band (6455.2 MHz), for broadcast to the user aircraft [3] . The GAGAN system will also provide benefits beyond aviation to all modes of transportation including maritime, highways, and railroads.
2. Literature Review
Several error sources that affect the positional accuracy of GPS are ionosphere, troposphere, satellite and receiver clock offsets, instrumental biases of the receiver and satellite, receiver measurement noise and multipath. Among all of these, the ionospheric delay is the most predominant error in GPS precise positioning and navigation. The measurement of ionospheric delay involves estimation of ionospheric total electron content (TEC). TEC is defined as number of free electrons contained in a tube of 1 m2 cross-section extending from the satellite to the receiver. The dual frequency GPS observations can be used to estimate the TEC, taking advantage of the
![]()
Figure 1. Indian SBAS “GAGAN” architecture.
dispersive nature of the ionosphere. However, line-of-sight ionosphere measurements derived from dual frequency GPS data are corrupted by instrumental biases in both the receiver and GPS satellite transmitters. The instrumental biases occur due to the frequency dependent delays of analog hardware within the GPS satellite and receiver [4] . An instrumental delay (bias) exists in the signals of each of the two GPS frequencies. The difference of the instrumental biases of individual frequencies is known as the differential instrumental bias [5] . Possible sources of interfrequency biases include the antennas, cables, and filters used in the GPS satellites and receivers. The biases for the individual frequencies do not affect the ionospheric delay measurement; only the differential instrumental biases affect the measurement. Therefore, the differential instrumental bias estimation and correction is necessary in the dual frequency GPS receivers, installed in the GAGAN network to obtain accurate ionospheric corrections that are transmitted to users as part of GAGAN.
3. Prominent Methods for TEC and Bias Estimation
Two prominent techniques based on the Kalman filter and the Self Calibration Of pseudoRange Error (SCORE) algorithm are applied to the low latitude GAGAN station data for estimation of TEC and instrumental bias. The Kalman filter is a computational algorithm that can optimally estimate the states of a system from a multidimensional signal contaminated by noise. Both these techniques use dual frequency GPS data of a single station to estimate the instrumental bias.
3.1. Instrumental Bias Estimation Using Kalman Filter
The TEC measurements from a dual frequency GPS receiver are affected by the thermal noise as well as the differential instrumental biases within the satellite and receiver hardware. A two step method is proposed to optimally combine both code and carrier phase observables for improving the TEC estimation accuracy. The ionospheric delay measurements derived from the code observables are unambiguous but are affected by the measurement noise and multipath errors. In the first step, line-of-sight ionospheric delay derived from the code observables is smoothed using the carrier phase derived ionosphere measurements. In the second step, a single layer ionosphere model is used to estimate the vertical TEC and instrumental biases from the smoothed line-of- sight TEC, using a five state Kalman filter.
The dual frequency GPS code and carrier phase measurements in metres can be described by the following equations (subscript i = 1, 2, refers to GPS frequencies, f1 and f2) [12] ,
(1)
(2)
where r is the true geometric range (m);
c is the speed of light (m/s);
dtu, dts are the receiver and satellite clock offsets, respectively (s); TD is the tropospheric delay (m);
Ii is the ionospheric delay at frequency fi (m);
SBPi and RBPi are the satellite and receiver instrumental group delay biases at frequency fi, respectively (m);
SBLi and RBLi are the satellite and receiver instrumental phase delay biases at frequency fi, respectively (m);
λi is the carrier wavelength at frequency fi (m);
Ni is the carrier phase integer ambiguity (cycles);
e(.) includes measurement noise and multipath error (m).
The ionospheric delay at frequency fi can be expressed as,
(3)
where 
and I is the differential ionospheric delay (m).
Using the Equations (1), (2) and (3), the differential ionospheric delay, I can be obtained as,
(4)
(5)
where SBP = SBP1 − SBP2, RBP = RBP1 − RBP2, SBL = SBL1 − SBL2, and RBL = RBL1 − RBL2. SBP and RBP are referred to as the satellite and receiver differential instrumental biases respectively.
The smoothing procedure used in this paper is described in [13] , which is usually called the Hatch filter. The Hatch filter is based on the concept that the change in code range between observations at different epochs approximately equals the change in carrier range. Smoothing reduces the noise level in the ionospheric measurements. In order to further improve the accuracy of smoothed ionospheric measurements, instrumental biases are estimated using the Kalman filter. In this approach, the smoothed line-of-sight differential delay is modeled as the sum of a receiver bias, a satellite bias, and the actual line-of-sight TEC.
The vertical TEC and instrumental biases are estimated using a five state Kalman filter. The biased smoothed line-of-sight TEC is represented using a single layer ionosphere model as [14] ,
(6)
where A1, A2, and A3 are the parameters for the spatial linear approximation of vertical TEC;
S is the slant function and e is elevation angle;
dlIP is the difference between the geomagnetic longitude of the ionospheric pierce point (IPP) and the mean longitude of the Sun;
dfIP is the difference between the geomagnetic latitude of the ionospheric pierce point and that of the receiver.
Equation (6) forms the measurement model of the Kalman filter. The parameters A1, A2, and A3 describing the vertical TEC on the IPP are estimated for each time t along with the instrumental biases using a five state Kalman filter. In this investigation, the instrumental bias obtained using Fitted Receiver Bias (FRB) method is used as the initial state of the receiver bias in the Kalman filter. However, only nighttime data are used in the FRB method for achieving better results. The details of the FRB method are discussed in a later section.
3.2. Instrumental Bias Estimation Using Score Algorithm
The Self-Calibration Of pseudo Range Error (SCORE) technique can be used to improve the accuracy of TEC measurement from the GPS observations. A dual frequency receiver contains several components such as antenna, low noise amplifier (LNA), cables, and filters in the RF and IF sections, which contribute to instrumental bias errors. SCORE technique can be used to calibrate such a dual frequency GPS receiver system. With this ability to calibrate and monitor the integrity of pseudo range measurements, SCORE algorithm can be used for ionospheric error measurement for GAGAN. The SCORE technique does not require use of any hardware calibrators or ionospheric models.
The SCORE technique can be applied to estimate the combined satellite plus receiver instrumental biases for each satellite. It uses a self-consistency constraint on the receiver’s measurements of ionospheric delay to derive the instrumental bias errors. This self-consistency can be understood by considering a conjunction occurring between two satellite paths, i.e. the two satellite paths arrive at the same moment at an ionospheric pierce point. In such an event, the same ionospheric pseudo range error (TEC value) should be seen on each satellite [15] . The algorithm operates by assigning bias values to each satellite, so as to minimize the difference in equivalent vertical TEC derived from two satellite observations having the same IPP latitude and local time. By considering multiple such pairs of satellites for various observations, a consistent set of bias corrections can be obtained. For each satellite, a correction term is produced, which represents the sum of the satellite plus receiver bias. The absolute satellite biases are provided by many agencies including the International GPS Service for Geodynamics (IGS). To obtain an absolute receiver system calibration, the IGS derived satellite biases may be subtracted from the respective combined satellite plus receiver bias corrections of obtained using the SCORE technique.
The mathematical quantity E describes the equivalent vertical TEC difference for multiple observations [16] ,
(7)
where, α: satellite PRN number; i: sample number;
= weighting factor between samples αi and βj;
Tγk = calculated equivalent vertical TEC for sample γk, using the appropriate local zenith angle and satellite bias,
(8)
for Sγk = slant TEC for the data sample;
Bγ = combined satellite plus receiver bias for PRN γ, in TECu;
εγk = elevation angle for satellite sample, at observing site;
μ = altitude scale factor for conversion to IPP zenith angle,
(9)
where Re = 6378 km (earth radius);
HIPP = 350 km (altitude of ionospheric pierce point);
The Gaussian function is chosen as an appropriate weighting factor,
(10)
for θk = Latitude for sample k (degrees).
λk = MJD + LT/24
where MJD is Modified Julian Day (day) and LT is local time (fraction of day), for sample k;
θ0 = reference latitude difference, for scaling (degrees);
λ0 = reference local time difference, for scaling (days).
For selection of IPP crossover points, an IPP latitude band of 4˚, centered 3˚ north of the receiver is considered. For Hyderabad station (Longitude: 78.47˚E, Latitude: 17.45˚N), this corresponds to 20.47˚N ± 2˚. In the computations, latitude scale parameter (θ0) is chosen as 16˚. This is much larger than the IPP latitude band range (4˚). The local time scale parameter (λ0) corresponds to a smaller spatial domain than θ0, and is comparable to the IPP latitude band range [16] . In the computations, λ0 is chosen as 4.5˚ in longitude.
3.3. Fitted Receiver Bias Method
A method for estimating the receiver bias of a single receiver is described by [7] . In most of the regions in the world, at any point of time, signals from 5 or more GPS satellites are received simultaneously by a GPS receiver. The slant TECs due to various satellites at any observation epoch are different due to different delays experienced by the individual satellites. The assumption made in the fitted receiver bias (FRB) method is that the vertical TEC values of all the satellites is identical, if the satellite and receiver biases are correctly removed. When the ionosphere has a horizontal gradient and vertical structure as in the actual case, the scattering of vertical TECs is assumed to be the smallest, if the instrumental biases are correctly removed. The FRB method is based on the minimization of standard deviation of vertical TEC computed from different satellites. The satellite biases are computed by several receivers data in the IGS network and are available in the public domain by many agencies such as the Centre for Orbit Determination (CODE), Europe [17] . Using the known satellite biases, the receiver bias is estimated by changing the value of receiver bias in a particular range and finding out a bias value that gives the minimum deviation of VTEC values to their mean. Here, a trial receiver bias b(i) is chosen, and the standard deviation of VTEC values due to all the visible satellites, with respect to the mean value is calculated at each observation epoch. Then the sum of standard deviations of the VTECs is computed for the whole day. The receiver bias value b(io) for which the standard deviation sum is minimum is considered as the correct receiver bias.
The standard deviation of VTEC data over the measurement period is given as,
(11)
where
(12)
where Mt denotes the total number of satellites and Nt is duration of the desired measurement time interval in samples. In equation (11), the total standard deviation is obtained by summing the standard deviation values of each measurement sample where Nt is chosen equal to the number of measurement samples present in 24 hours
of GPS data. ![]()
is average of all VTEC from Mt satellites.
4. Results and Discussion
In this analysis, dual frequency GPS data of Hyderabad GAGAN station is used. The data is provided by the Space Applications Centre (SAC), Indian Space Research Organisation, Ahmedabad, India. The raw data obtained from the NovAtel GPS card is converted into the desired Receiver Independent Exchange (RINEX) format using the “Convert” software. The RINEX navigation and observation data files are used in the processing.
4.1. TEC and Bias Estimation Results Using Kalman Filter
From the navigation data, position of all the visible satellites in the Earth Centered Earth Fixed (ECEF) reference system is computed [18] . From the observation data, dual frequency code and carrier phase observables of all the visible satellites is extracted. These are used to compute the biased phase smoothed slant TEC using the Hatch filter technique. Four satellites giving lowest PDOP are selected at each epoch. Using the selected satellites’ position and the corresponding pseudorange information, receiver coordinates in ECEF reference system are computed at each epoch using the Bancroft algorithm [19] . The geomagnetic latitude of the receiver is computed at each epoch for later use. Using the satellite and receiver position information, elevation and azimuth angle are computed at each epoch. Further, the slant factor, geomagnetic latitude and longitude of IPP, mean longitude of Sun are calculated at each epoch.
Figure 2 shows the phase smoothed slant TEC obtained using the Hatch filter. The smoothing process reduces the code measurement noise and multipath errors to a reasonable extent.
The computed biased phase smoothed slant TEC, slant factor, geomagnetic latitude and longitude of IPP, geomagnetic latitude of receiver, and the mean longitude of Sun form the inputs to the Kalman filter. The satellite instrumental biases determined by the Centre for Orbit Determination (CODE), Europe are used as initial state of the biases of various satellites, in the proposed Kalman filter model.
Normally, FRB method is used for rough estimation of instrumental bias of GPS receiver considering the 24-hour data of several days. The FRB method is applied to GPS data of Hyderabad GAGAN station (0 - 24 hours LT, March 4, 2005). In this method, the slant TEC is converted into equivalent vertical TEC using the
![]()
Figure 2. Biased phase smoothed slant TEC estimated using hatch filter.
slant function. The sampling period of the data is 60 s. Hence, Nt is chosen as 1440. Literature suggests that the receiver bias can be as large as about ±15 ns [9] . In the computations, the trial receiver bias values are varied from −30 ns to 30 ns in the steps of 0.1. For each receiver bias value, total standard deviation of VTEC at each observation time is calculated. The variation of the total standard deviation of VTEC with respect to the receiver bias is shown in Figure 3. From the figure, the receiver bias is estimated as −0.6 ns corresponding to the minimum total standard deviation. As we know that the TEC is diurnal and random in nature, we considered only nighttime data which are almost constant and smooth. The receiver bias value obtained using the FRB method applied to nighttime data (20:00-08:00 hours local time) is 4.1 ns. This value is taken as initial state value for the Kalman filter for further improving the prediction.
The estimated vertical TEC after removal of the instrumental bias error is shown in Figure 4.
The estimated maximum vertical TEC of two satellites (PRN 4 and PRN 7) visible during mid-day, after correcting for the instrumental biases are 47.12 and 52.87 TECU. The mean value of the receiver bias estimated using the Kalman filter is obtained as 3.2 ns (1 ns of differential delay = 2.852 TECU).
4.2. TEC and Bias Estimation Results Using Score Algorithm
The ionospheric TEC and the combined instrumental bias of satellite and receiver are estimated using SCORE technique. The dual frequency GPS data of Hyderabad GAGAN station is used in the estimation. The slant TEC is estimated using carrier phase measurements for all the visible satellites over a 24 hour period (4 March 2005). Knowing the satellite and receiver position information, elevation and azimuth angles due to various satellites are computed. Further, the IPP coordinates, and mapping (slant) function is computed for later use. The slant TEC values are converted into equivalent vertical TEC estimates by using the mapping function.
![]()
Figure 3. Estimated receiver bias using FRB method.
![]()
Figure 4. Estimated vertical TEC using Kalman filter.
For input to the SCORE algorithm, vertical TEC values and the IPP latitudes over a narrow band of 4˚, centered 3˚ north of the receiver (i.e. 20.5˚N ± 2˚) is considered. The vertical TEC and IPP latitude variation over this limited band are shown in Figure 5. The maximum TEC before bias correction is observed to be 36.34 TECU at 15:48 Hrs LT for SV 15. The conjunction points of IPP latitudes of various satellite pairs called as crossover points are identified within the band of 4˚ spacing.
The vertical TEC values of each satellite pair (corresponding to the crossover points) in the chosen latitude band are used in computing E. The combined satellite plus receiver biases due to various satellites are estimated using a minimization procedure for E. The estimated biases due to various satellites are shown in Table 1.
The absolute values of the combined biases are found to range from 0.66 to 20 TECU over a 24 hour observation period. The system calibration parameter (SCP) is computed by taking the average of bias corrections of all the visible satellites. This can be used to calibrate the GPS receiver system. The SCP is estimated as −8.29 TECU. The vertical TEC profile after application of the SCORE algorithm is shown in Figure 6.
![]()
Figure 5. Vertical TEC and IPP latitude variations (18.5˚ - 22.5˚N), before applying score algorithm.
![]()
Figure 6. Vertical TEC and IPP latitude variations (18.5˚ - 22.5˚N), after applying score algorithm.
![]()
Table 1. Combined satellite plus receiver instrumental biases estimated using score algorithm.
5. Conclusion
In this paper, two prominent techniques based on Kalman filter and SCORE algorithms are used for estimation of TEC and instrumental biases. In order to reduce the noise level in the GPS pseudorange data, ionospheric delay measurements are smoothed using precise carrier phase data. The phase smoothed slant TEC measurements obtained using the Hatch filter technique show considerable improvement over the code derived slant TEC. To further improve the accuracy of TEC estimation, a five state Kalman filter is developed for estimating the differential instrumental bias. An FRB method is used for estimating the initial state of the receiver bias. However, only nighttime data are considered in the estimation of receiver bias using FRB method, as TEC variations are relatively small during nighttime. The biases are found to be consistent over the observation period and agree with other reported values in open literature. Using the SCORE technique, the combined satellite and receiver biases are estimated. The SCORE approach is distinct in that it allows the dual frequency GPS receivers to autonomously maintain its pseudorange accuracy without use of any hardware calibrators or ionospheric models. The ionospheric delay corrections obtained after removal of instrumental bias would enable the user aircraft having a SBAS enabled GPS receiver to determine its position accurately for Category-I Precision Approach applications.
Acknowledgements
The author would like acknowledge the Director, Space Applications Centre, Indian Space Research Organisation, Ahmedabad, India for providing the data.