Error Characteristics of GNSS Derived TEC

The Global Navigation Satellite System (GNSS) allows for the cost-effective estimation of the ionospheric total electron content (TEC). However, research on error characteristics of the derived TEC is scarce, which provides insights into the quality of the GNSS ionospheric observation. We investigate characteristics of errors in the derived TEC with data from ~260 GNSS dual-frequency receivers of the Crustal Movement Observation Network of China (CMONOC). The slant TEC is calculated from carrier phase measurements and the vertical TEC over China is fitted with a spatial resolution of 1° by 1° in latitude and longitude in four seasons of 2014. It is found that the errors of both the slant TEC and the derived TEC follow Laplace distribution rather than Gaussian distribution in all seasons. The errors of the slant TEC have sharper peaks than those of the derived TEC. The Mean Absolute Error (MAE) and Root Mean Square Error (RMSE) of the slant TEC are typically 0.04 TECU and 0.2 TECU, while the MAE and RMSE of the fitting residuals for the derived TEC are typically 1 TECU and under 2 TECU, respectively. Both MAEs and RMSEs of the derived TEC have the largest value in spring and the smallest value in summer, while the seasonal dependence is only observed in RMSE of the slant TEC.


Introduction
In the last three decades, beacon signals transmitted from global navigation satellite systems (GNSS) have widely been utilized for ionosphere studies. The primary parameter obtained from the satellite beacon signals is the total electron content (TEC), derived from measurements of differential propagation delay and carrier phase between two frequencies, f 1 and f 2 , transmitted from GNSS satellites and propagating through the dispersive ionosphere. A basic problem in TEC derivation is the signal delays in the electronics of satellites and receivers being different for the two frequencies. It is a type of systematic error, which is often referred to as inter-frequency instrumental bias. To remove the instrumental biases, a mathematical model should be established to describe the relations among TEC, biases, and measurements, in which data characteristics and computer resources should also be considered, and then data fitting should be carried out. This makes the TEC derivation a challenging task and efforts never cease to pursue goals, such as better accuracy and faster processing time. Various algorithms have been developed for a single GNSS receiver, small-scale GNSS networks, regional, and global GPS networks. The ionosphere is generally modeled as a thin shell with a spatially smooth distribution of condensed electrons. With functional expressions of TECs in time, latitude, and longitude, equations are built to relate TEC, instrumental biases, and measurements. Then, the equations are solved with numerical computation [1][2][3][4][5][6][7][8][9].
As to the random errors in the derived GNSS TEC, the term "fitting error" is often used to specify the statistical characteristics, since the TEC derivation is essentially the result of data fitting. Most of the available algorithms are based on pseudo-range measurement, which is smoothed with carrier phase for better precision. The fitting error has been reported within~5 TECU [10][11][12][13]. Recently, ionospheric TEC with an uncertainty lower than 1.0 TECU was achieved by using the carrier phase measurements from the EUREF Permanent Network (EPN) and European Position Determination System (EUPOS) [14]. A root mean square error of 0.40 TECU in TEC derivation was obtained with carrier phase from the dense GNSS Earth Observation Network (GEONET) of Japan [15].
The current fitting errors for GNSS TEC are quantified simply by methods which implicitly assume the random errors follow a normal or Gaussian distribution. However, the error characteristics of GNSS TEC have rarely been addressed or even not verified to be the Gaussian distribution. In construction of the global empirical background TEC model from GNSS TEC Mukhtarov et al. [16] noticed that the probability function of the fitting residuals resembles more of the Laplace distribution than Gaussian distribution. Alternatively, the process of GNSS TEC derivation involves using data from different satellites and receivers. As a matter of fact, it is widely recognized that different equipment can produce widely divergent measurements. Hsu [17] demonstrates that the probability distribution of a mixture of non-uniform Gaussians in air and sea navigation tend to be Laplace with a more pronounced peak and two longer tails. Deviation from Gaussian tends to occur in situations when heterogeneity exists in a large number of observations collected under comparable conditions [18]. Therefore, the ideal assumption of the uniform Gaussian distribution used in GNSS TEC derivation is questionable and needs to be evaluated.
We investigate characteristics of errors in both GNSS measurement and the derived TEC with data from the Crustal Movement Observation Network of China (CMONOC). Section 2 describes the GNSS observations of CMONOC, the grid-based algorithm for deriving TEC from the differential carrier phases between two frequencies of satellite signals, and the error analysis method. Section 3 presents the result from the application of the methods to CMONOC data in 2014. We show the derived TEC with its spatiotemporal variations. Error analysis is made in detail by investigating not only the error distributions, but also the statistical parameters. Summary of the work is given in the last section.

GNSS Observation
The CMONOC includes~260 continuous operating GNSS receivers. As shown in

Algorithm
As a dispersive medium, the ionosphere causes different group delays and phase advances for radio waves in different frequencies. The group delay and the phase advance are proportional to the TEC along the line of sight. A GNSS receiver measures the

Algorithm
As a dispersive medium, the ionosphere causes different group delays and phase advances for radio waves in different frequencies. The group delay and the phase advance are proportional to the TEC along the line of sight. A GNSS receiver measures the pseudoranges from satellites to the receiver and the carrier phases of the satellite signals, which provide a mean to obtain TEC.

STEC Measurement and Arc Bias
Having higher precision than that of the pseudo-range, the carrier phase observables L 1 and L 2 are used to obtain the slant TEC (stec) from a satellite to a receiver.
where k, related to the ionosphere refraction, is 80.62 (m 3 /s 2 ); λ 1 and λ 2 are the wavelengths of the two signals; and stec is a relative value, since the instrumental biases and 2π ambiguity are included in the carrier phase measurements. Generally, the instrumental bias can be assumed to be constant for a period of time; the ambiguity does not change in a consecutive measurement (referred to as an arc). Here we define an arc bias, b arc , to represent the inter-frequency hardware bias mixed with integer ambiguity. Then the absolute slant TEC, sTEC, can be obtained by taking out the arc bias from the relative stec measurement.
When calculating stec, a data screening step is needed, especially to check cycle-slips in either L 1 or L 2 phase measurements [19]. If cycle-slips are identified for a certain arc, it would be separated to different arcs, since the 2π ambiguity changes with cycle slip. The arc biases are hence different even though the measurements are taken by the same pair of satellite and receiver.

The Ionosphere Model
The ionosphere is modeled as a thin shell of a hollow sphere at height h (400 km) and co-centered with the Earth at O, as shown in Figure 2. The line of sight for a satellite and receiver pair crosses the ionosphere at point P, which is called an ionospheric piercing point (IPP). The latitude (ϕ), longitude (θ), and zenith (χ) of the IPP are related with the satellite and receiver positions as the following, φ = arcsin(sin φ Rx cos β + cos φ Rx sin β cos γ) (3) where, ϕ Rx and θ Rx are the latitude and longitude of the receiver, α and γ are the elevation and azimuth of the satellite, and R e is the Earth radius. With a smoothly varying assumption of the ionosphere, from the slant TEC given by (2) the vertical TEC at the IPP can be converted as To suppress the effect of multipath, the cutoff angle for the satellite elevation is set to be 30 • when extracting the slant TEC stec and calculating the zenith angle χ.
where, φRx and θRx are the latitude and longitude of the receiver, α and γ are the elevation and azimuth of the satellite, and Re is the Earth radius. With a smoothly varying assumption of the ionosphere, from the slant TEC given by (2) the vertical TEC at the IPP can be converted as To suppress the effect of multipath, the cutoff angle for the satellite elevation is set to be 30° when extracting the slant TEC stec and calculating the zenith angle χ.

Matrix Equations Construction
For regional GNSS networks, their receivers tend to be densely distributed in certain areas. The IPPs are accordingly accumulated. With a smoothly varying behavior, the ionospheric TECs at adjacent IPPs would have approximately the same value, which enables functional relations built among slant TECs and arc biases. We divide the ionosphere into grids of the same size over the areas concerned. Further, we assume that the vertical TEC is identical at any point within one grid block at a given time, and the vertical TECs can differ in different grid blocks at the time. Next, mathematical relations among the vertical TEC, arc biases, and the slant TEC can be constructed, as long as two or more than two lines of sight pass through the same grid.
For the signals from satellite j(m) to receiver k(n) we use subscript jk(mn) to denote their consecutive measurement or arc. If their lines of sight propagate through the same grid l at time t of epoch i, the vertical TECs, VTECijk from satellite (j) to receiver (k), and VTECimn from satellite (m) to receiver (n) would be identical. Then using (7), we have the following equations

Matrix Equations Construction
For regional GNSS networks, their receivers tend to be densely distributed in certain areas. The IPPs are accordingly accumulated. With a smoothly varying behavior, the ionospheric TECs at adjacent IPPs would have approximately the same value, which enables functional relations built among slant TECs and arc biases. We divide the ionosphere into grids of the same size over the areas concerned. Further, we assume that the vertical TEC is identical at any point within one grid block at a given time, and the vertical TECs can differ in different grid blocks at the time. Next, mathematical relations among the vertical TEC, arc biases, and the slant TEC can be constructed, as long as two or more than two lines of sight pass through the same grid.
For the signals from satellite j(m) to receiver k(n) we use subscript jk(mn) to denote their consecutive measurement or arc. If their lines of sight propagate through the same grid l at time t of epoch i, the vertical TECs, VTEC ijk from satellite (j) to receiver (k), and VTEC imn from satellite (m) to receiver (n) would be identical. Then using (7), we have the following equations b arc jk cos(χ ijk ) − b arc mn cos(χ imn ) = stec ijk cos(χ ijk ) − stec imn cos(χ imn ) (8) where, stec ijk (stec imn ) is the slant TEC from satellite j(m) to receiver k(n) at time t of epoch i, b arc jk (b arc mn ) is the arc bias for the consecutive measurement of satellite j(m) and receiver k(n). A set of equations can be obtained by using a data set from GNSS network, which is expressed in the form of matrix equation as the following: The coefficient matrix H has only two non-zero elements in every row. The non-zero elements on a column of H come from the same arc. It can be speculated that the columns of H are linearly independent.

Solution of the Matrix Equations
For GNSS observation with 1/30 Hz sampling rate, there are~20,000 arcs for a network with~260 receivers in one day. The arcs are increased if cycle slips are detected. Consequently, the number of unknowns is about 20,000, and the number of equations can be over one million with 30 • cutoff angle of satellite elevation. Equation (9) is an overdetermined system. It has no solution in the ordinary sense, but rather has a least-squares solution.
There are several ways to solve the least-squares problem. The normal equation is the most popular equation, which is obtained by pre-multiplying (9) with the transpose of the coefficient matrix, H T . However, the multiplication amplifies the error in the data, and the computations also accumulate round-off errors. A theoretically equivalent is the singular value decomposition (SVD) approach. The SVD factorizes the coefficient matrix H into two unitary matrices and a diagonal matrix composed of singular values. It is widely used for rank deficient and ill-conditioned matrices [20]. Another solution is determined by the QR factorization of the matrix H. It is more reliable than that of the normal equation if the columns of H are linearly independent [21].
The SVD fixes the round-off problem and produces a best approximation solution in the least squares sense [20]. However, the SVD is computationally expensive. Additionally, memory problems generally occur when creating and handling large matrices directly in MATLAB on desktop computers. Note that the coefficient matrix H in (9) has only two non-zero elements in every row; it can be stored in sparse form and then handled. MATLAB has efficient solving tools based on QR decomposition, which is applicable to a large sparse matrix.
With a QR factorization of H, we have: Then from (9) and (10) we have the arc biases as given by: which is the least squares solution of (9). It can be obtained by using matrix left division function of MATLAB.

Program Flowchart
The algorithm is carried out in MATLAB. The corresponding program is named Grid_TEC and developed in MATLAB 2016a (operational in higher versions). The flowchart of Grid_TEC is shown in Figure 3. The inputs are GNSS observations recorded in RINEX format and precise ephemeris of satellites in sp3 format. Carrier phase measurements are read out and stored for each receiver. Possible cycle slip is identified to separate the arc, since it leads 2π ambiguities change in the phase measurements. Next, slant TEC is extracted and stored for each arc. Satellites positions are read out from the sp3 files and interpolated for observation times when the precise ephemeris is not provided. The  The matrix equation is built with slant TECs, their IPP positions, and grids information. Next, Linear Least Squares solution of ambiguity combined biases is obtained. VTECs are calculated from slant TEC and stored with grids according to their IPP positions.

Error Analysis Method
We analyze two kinds of measurement dealt with in the Grid_TEC algorithm. The first concerns the input, which is stec from (1). The second is the output VTEC in (7), which is related with the input stec and the fitting result of arc bias.

Error Estimation
An error is the difference between the observed or calculated value of a quantity and its true value. Since a true value for any measurable quantity is not known, an average of many repeated independent measurements is used to replace true value. However, it is impossible to repeat measurements in many practical observations where the measured quantity varies with time. A moving average or running mean is often used.
The measurement error of stec is estimated by examining the dispersion of the stec to its running mean, For the GNSS observation with a 30-s resolution, Δsteci is calculated by taking the standard deviation of stec for 5-min interval (11 consecutive epochs).
Concerning the error of the derived VTEC, we can take the residual defined by the difference between the left-hand-side and the right-hand-side of (8), which assumes the identity of two VTECs within a grid.
The smaller VTEC Δ is, the better the assumption and the fitting are.

Error Analysis Method
We analyze two kinds of measurement dealt with in the Grid_TEC algorithm. The first concerns the input, which is stec from (1). The second is the output VTEC in (7), which is related with the input stec and the fitting result of arc bias.

Error Estimation
An error is the difference between the observed or calculated value of a quantity and its true value. Since a true value for any measurable quantity is not known, an average of many repeated independent measurements is used to replace true value. However, it is impossible to repeat measurements in many practical observations where the measured quantity varies with time. A moving average or running mean is often used.
The measurement error of stec is estimated by examining the dispersion of the stec to its running mean, For the GNSS observation with a 30-s resolution, ∆stec i is calculated by taking the standard deviation of stec for 5-min interval (11 consecutive epochs).
Concerning the error of the derived VTEC, we can take the residual defined by the difference between the left-hand-side and the right-hand-side of (8), which assumes the identity of two VTECs within a grid.
The smaller ∆VTEC is the better the assumption and the fitting are.

Statistical Characteristics
The statistical parameters, the median of the error (µ), mean absolute error (MAE, λ), and root mean square error (RMSE, σ), are used to qualify the observation and fitting result. With ∆ i representing the general error, the statistical parameters are calculated by where, N is the amount of the observations or number of the equations in (9). The MAE result will always be smaller or equal to the RMSE. Their difference demonstrates the variance in the individual errors. The ratio of the MAE to the RMSE is a fixed value for a certain dataset. Both RMSE and MAE are important to error analysis [22]. It should be noted that the MAE is mostly suitable to describe uniformly distributed errors. When using the RMSE a normal (Gaussian) distribution is generally an underlying assumption. Generally, the errors are assumed to follow Gaussian distributions, since it is the most important probability distribution in statistics, which fits many natural phenomena. The MAE is related with the RMSE by λ = σ √ 2/π (λ ≈ 0.8σ) for Gaussian distribution. However, practical data and its processing result would not follow an ideal Gaussian distribution. We also use histograms to graphically show the distributions of the errors. Proper distributional model can be estimated with probability plot.

Results
Following the flowchart of Grid_TEC algorithm, data processing is made on a daily base. The number of established equations is about 1.5 million and the number of unknowns is about 15 thousand. It takes about 25 min to run the program, from reading the GNSS data to solving the Equation (11), by a desktop computer (PC) using an Intel ® Core I i7-700 CPU.

Derived VTEC
With the bias given by (11), VTEC is calculated from (7). Concerning the diurnal variation, the variations of VTEC both in line and circle shows a high degree of conformity. Note that the local time (LT) at 110 • E is~7.3 h advanced. VTEC reaches its peak at~12:30-15:30 LT in China sector. The peak appears at an earlier time with the increase of latitude. The value of the peak at low latitude is much larger than that at middle latitude, implying a large gradient from 25 • N to 35 • N. VTEC gets to a valley at~21:00-22:00 LT, with insignificant latitudinal difference. Seasonal variation is remarkable. VTEC is larger in equinoxes than in summer. The peak value of VTEC in spring can be close to 130 TECU, while it is around 60 TECU in summer. It is noteworthy that VTECs in daytime are larger in winter than in autumn. This agrees with previous studies, which attribute the result to greater contribution from plasmasphere.
The TEC map shows more details of the spatial variation of VTEC. The longitudinal range from 70 to 140 • E at 7:00 UT corresponds to 11:40~16:20 LT. Latitudinal variation or VTEC latitude gradient is conspicuous, due to the equatorial ionization anomaly (EIA), which centered around 20 • N. The magnitude and latitudinal range of the EIA are the largest in spring. They are the smallest in summer, and larger in winter than autumn. If converting the longitude into LT, the largest latitude EIA reaches~25 • N at the time between 12:30 and 13:00 LT in spring. The location of EIA crest does not change significantly in other seasons.  Concerning the diurnal variation, the variations of VTEC both in line and circle shows a high degree of conformity. Note that the local time (LT) at 110° E is ~7.3 h advanced. VTEC reaches its peak at ~12:30-15:30 LT in China sector. The peak appears at an earlier time with the increase of latitude. The value of the peak at low latitude is much larger than that at middle latitude, implying a large gradient from 25° N to 35° N. VTEC gets to a valley at ~21:00-22:00 LT, with insignificant latitudinal difference. Seasonal variation is remarkable. VTEC is larger in equinoxes than in summer. The peak value of VTEC in spring can be close to 130 TECU, while it is around 60 TECU in summer. It is noteworthy that VTECs in daytime are larger in winter than in autumn. This agrees with previous studies, which attribute the result to greater contribution from plasmasphere.
The TEC map shows more details of the spatial variation of VTEC. The longitudinal range from 70 to 140° E at 7:00 UT corresponds to 11:40~16:20 LT. Latitudinal variation or VTEC latitude gradient is conspicuous, due to the equatorial ionization anomaly (EIA), which centered around 20° N. The magnitude and latitudinal range of the EIA are the largest in spring. They are the smallest in summer, and larger in winter than autumn. If converting the longitude into LT, the largest latitude EIA reaches ~25° N at the time between 12:30 and 13:00 LT in spring. The location of EIA crest does not change significantly in other seasons.

Error Analysis
Error analysis is carried out by using the definitions and procedures described in the above error analysis method.

Measurement Error
The measurement error Δsteci is estimated with (12) over a 5-min running window for all receivers. Table 1 gives the statistical parameters based on three days for each season in 2014. The biases of the errors are 0 for all the seasons. The MAE is 0.04 TECU,

Error Analysis
Error analysis is carried out by using the definitions and procedures described in the above error analysis method.

Measurement Error
The measurement error ∆stec i is estimated with (12) over a 5-min running window for all receivers. Table 1 gives the statistical parameters based on three days for each season in 2014. The biases of the errors are 0 for all the seasons. The MAE is 0.04 TECU, except for autumn when the MAE is 0.03 TECU. The RMSE is typically~0.2 TECU, from 0.12 TECU in summer to 0.22 TECU in spring. It is 0.19 TECU in winter, larger than that of 0.15 TECU in autumn. Therefore, the RMSE is different among different seasons though the MAE is the same as or equal to one another. Moreover, it is noteworthy that λ is much smaller than σ √ 2/π. The typical ratio of λ/σ is 0.2. This indicates that the measurement errors do not follow the normal distribution. To examine the error characteristics, Figure 5 presents the normalized frequency distribution of stec measurement errors with a bin of 0.025 TECU for the four seasons of 2014. Nearly symmetrically, each histogram peaks sharply around the center, then decreases steeply to a low level and spreads noticeably afterwards to both sides. The error distribution is found to follow a double exponential or Laplace distribution rather than a Gaussian distribution. The probability density function of Gaussian distribution is associated with the median and RMSE parameters.
The Laplace distribution can be expressed with the median and MAE.
We use N(µ, σ 2 ) and La(µ, λ) to express the Gaussian and Laplace distributions with statistical parameters, respectively. As shown in Figure 5, the stec error distributions in the four seasons are fitted by curves of Laplace distribution (in orange) quite closely. They differ from Gaussian distribution (in red) greatly. To examine the error characteristics, Figure 5 presents the normalized frequency distribution of stec measurement errors with a bin of 0.025 TECU for the four seasons of 2014. Nearly symmetrically, each histogram peaks sharply around the center, then decreases steeply to a low level and spreads noticeably afterwards to both sides. The error distribution is found to follow a double exponential or Laplace distribution rather than a Gaussian distribution. The probability density function of Gaussian distribution is associated with the median and RMSE parameters.
The Laplace distribution can be expressed with the median and MAE.
We use La μ λ to express the Gaussian and Laplace distributions with statistical parameters, respectively. As shown in Figure 5, the stec error distributions in the four seasons are fitted by curves of Laplace distribution (in orange) quite closely. They differ from Gaussian distribution (in red) greatly.   Table 2 lists the cumulative probability (%) of the measurement errors for the typical ranges in four seasons of 2014. Theoretically, the cumulative probability for Laplace distribution is 63.2%, 86.5%, and 95.0% for ±λ, ±2λ, and ±3λ ranges. It can be seen that the cumulative probability of ∆stec is larger than that of theory within ±λ and ±2λ. It agrees with that of the theory for ±3λ. For an ideal Laplace distribution, we have λ/σ = 0.5. The λ/σ in Table 1 ranges from 0.19 to 0.33, which are smaller than 0.5. The measurement error has a sharper peak and shorter tails than the ideal Laplace distribution. It is widely recognized that different equipment can produce widely divergent measurements. Hsu [17] demonstrates that the probability distribution of a mixture of nonuniform Gaussians in air and sea navigation tend to be Laplace with a more pronounced peak and two longer tails. Deviation from Gaussian tends to occur in situations when heterogeneity exists in the population [18]. This can explain the Laplace distribution of the errors in slant TEC measurement, since the noise level of different satellite-receiver pair is affected by both the hardware and space environment, which are different from one another.

Error of the Derived TEC
The error of VTEC is calculated with (13). Table 3 lists the statistical parameters of the residuals in the four seasons of 2014. The bias of the fitting errors, µ, is 0 for all seasons. The MAE, λ, is in 1 TECU level. In spring/summer it is 1.12/0.70 TECU, which is the largest/smallest in the four seasons. λ is 0.87/0.89 TECU in autumn/winter. The RMSE, σ, has the same dependence on seasons as λ. It is under 2 TECU for all seasons. It is 1.94/1.12 TECU in spring/summer, which is the largest/smallest. σ is 1.46 TECU and 1.60 TECU in autumn and winter, respectively. For residuals λ is not equal to σ √ 2/π either. λ/σ is about 0.6, which is smaller than that of Gaussian distribution and larger than that of Laplace distribution. Therefore, the fitting errors do not follow the normal distribution even though they are unbiased.  Figure 6 presents the normalized frequency distribution with a bin of 0.5 TECU for the residuals in the four seasons of 2014. Overlaid are curves for probability density functions of Gaussian distribution (N(µ, σ 2 )) in red, and Laplace distribution ( La(µ, λ)) in orange, respectively. The fitting residuals resemble Laplace distribution more than Gaussian distributions.    Table 4 lists the cumulative probability (%) of the fitting errors for the typical ranges in four seasons of 2014. The cumulative probability tends to be the similar level for different seasons. They generally agree with those of ideal Laplace distribution, especially for ±2λ and ±3λ ranges. With RMSE, 80.9%, 95.1%, and 97.8% of errors fall within ±σ, ±2σ, and ±3σ, respectively. Noting that for Gaussian distribution it is 68.3%, 95.5%, and 99.7%, respectively, there are~13% more errors fall in ±1σ, there are~2% more cases outside ±3σ. The result indicates that most of the derived VTEC have an error less than 1σ, however, the outliers of 3σ tend to be 2%. It can be seen that both MAE and RMSE are important parameters to characterize the errors deviate from Gaussian distribution. A smaller value of λ/σ means larger deviation from Gaussian distribution and less outliers outside ±3σ. The error of the VTEC follows the Laplace distribution as that of stec does since the VTEC is derived from stec. In evaluation of the global TEC maps from GNSS derived TEC, Mukhtarov et al. reported that the probability function of the residuals follows the Laplace distribution rather than Gaussian distribution [16]. However, the peak of the distribution for VTEC is not as sharp as that of stec. Both of the MAE and RMSE for VTEC are much larger than those for stec. It is reasonable to speculate that the assumption used in the TEC derivation should be responsible. Simulation should be needed to clarify how the error characteristics of the derived VTEC vary with the ionospheric model and the measurement error.

Summary
We develop the Grid_TEC algorithm for deriving TEC with carrier phase observables from regional GNSS networks. The ionosphere over the area of the network is divided into grids with the same size. A set of equations is established to relate the vertical TEC in each grid with the observations. The equations are solved with QR decomposition. Written in MATLAB, the algorithm can be run in a desktop computer, conveniently. The performance of the Grid_TEC is evaluated by applying it to data of 2014 from CMONOC Network of China and making error analysis of the measurement and the fitting result. The TEC over China is derived with a spatial resolution of 1 • by 1 • in latitude and longitude. Daily and seasonal variations of TEC are observed. The north crest of the equatorial anomaly is captured with TEC map, which also display seasonal variations.
The error characteristics of both the measured stec and the derived VTEC are investigated. While the MAE and RMSE of the slant TEC error are typically 0.04 TECU and 0.2 TECU, the MAE and RMSE of the fitting residuals for VTEC are typically 1 TECU and under 2 TECU, respectively. Both the error distributions of the GNSS measurement and the derived TEC follow Laplace distribution rather than Gaussian distribution. The error distribution for the slant TEC has a sharper peak and shorter tails than the ideal Laplace distribution. The distribution of the fitting residuals has a sharper peak and longer tails than Gaussian distribution. Both MAE and RMSE are useful to represent the error characteristics. A smaller ratio of MAE to RMSE, λ/σ, leads to a sharper peak and two shorter tails of the distribution. The ratio of MAE to RMSE can also be used to indicate how much the error distribution deviates from Gaussian. The smaller the value of λ/σ is, the lager deviation of the distribution from Gaussian. Moreover, there tends to have less outliers outside ±3σ for a smaller value of λ/σ.