Real-Time Monitoring of Ionosphere VTEC Using Multi-GNSS Carrier-Phase Observations and B-Splines

Monitoring the Vertical Total Electron Content (VTEC) of the ionosphere is important for applications ranging from navigation to detection of space weather events. Therefore, emerging efforts have been made by several analysis centers to estimate the VTEC using different approaches in real‐time. Global Navigation Satellite Systems (GNSS) is a crucial technology for ionosphere modeling due to its worldwide distributed receivers, high temporal resolution, and low latency data dissemination capability. The selection of a convenient approach to extract ionosphere information from GNSS and the representation of VTEC by an appropriate mathematical model are essential factors for providing fast and accurate ionosphere products. Contrarily to the widespread phase‐leveling method, which uses noisy and erroneous code measurements, the modeling concept in this paper utilizes pure carrier‐phase measurements. Measurements acquired through the International GNSS Service (IGS) real‐time service in Radio Technical Commission for Maritime Services format are from GPS, GLONASS, and GALILEO. The measurement biases, including the ambiguity of carrier‐phase measurements, are simultaneously estimated along with VTEC model parameters. In our approach, VTEC is represented by B‐spline expansions embedded into a Kalman filter. Due to their localizing feature, B‐splines form a highly sparse structure in the filter measurement model. Thus, matrix operations for large‐scale problems can be performed fast using sparse matrix operations, as is done in this study. The differential slant total electron content (dSTEC) analysis and the comparison with Jason‐3 altimetry VTEC were performed for validation within selected periods in 2019. The dSTEC analysis shows that the quality of the generated real‐time VTEC maps slightly outperforms those provided by the other IGS analysis centers.

The accuracy of VTEC models varies according to the handling of carrier-phase ambiguities, satellite/receiver biases, spherical shell specifications, and mathematical approaches for the VTEC representation (Banville et al., 2011). The ionospheric information, that is, the slant total electron content (STEC), can be extracted from the geometry-free combination of GNSS signals at distinct frequencies. The pseudo-range (code) ionosphere combination generates unambiguous but rather noisy information. On the other hand, the carrier-phase ionosphere combination is about two orders of magnitude more precise, but it is biased by an unknown number of phase cycles (see, Hofmann-Wellenhof et al., 2008;Rovira-Garcia et al., 2016). The leveling technique is a widely used method to remove the ambiguity biases in carrier-phase observations with the aid of code observations (Mannucci et al., 1993(Mannucci et al., , 1998. Although this method is capable of eliminating the biases, the leveled observations are contaminated by significant errors up to several TEC units due to the high noise level, time-varying inter-frequency biases, and multi-path effects introduced by code observations (Ciraolo et al., 2007). Alternative approaches based solely on carrier-phase observations or precise point positioning (PPP) techniques have been used to avoid the leveling procedure, see Hernández-Pajares et al. (2011). The studies of Banville et al. (2011Banville et al. ( , 2012, Chen et al. (2018) and Xiang et al. (2019) demonstrated that the accuracy of STEC observations could significantly be improved by removing the phase ambiguities derived from a PPP solution compared to those obtained from the standard leveling approach. Plenty of methods relying on PPP techniques have been introduced to extract accurate ionospheric information from (un)differenced and (un)combined GNSS observations (see, e.g., T. Liu et al., 2018;Ren et al., 2020;Xiang et al., 2017;Zhang, 2016). The main drawback of many PPP-based ionosphere solutions is the requirement for precise auxiliary information, which usually includes, for example, satellite orbits, clock data, and even reference station coordinates in some studies. This requirement introduces another challenging task of obtaining the auxiliary data precisely and reliably in real-time.
As mentioned above, to avoid the leveling procedure, the phase biases can also be obtained along with other ionospheric target parameters by solely utilizing the ionosphere combination of carrier-phase observations. For instance, this approach was carried to estimate the global electron density distribution by employing a voxel-based tomographic model; see, Hernández-Pajares et al. (1999). Brunini and Azpilicueta (2009) analyzed the effect of the leveling technique and modeling errors on the calibration of STEC biases. The VTEC model used in their evaluations is based on the thin shell assumption and bi-quadratic expansions defined over a given receiver. They showed that the bias calibration using the phase leveling approach dominantly suffers from both the leveling and model errors. However, the use of the carrier-phase ionosphere combination alone is free of leveling error, and therefore, the modeling error is the primary source contributing to the total error budget. Another study for the STEC bias calibration based on carrier-phase observations was presented by Krypiak-Gregorczyk and Wielgosz (2018). This approach was employed to estimate the carrier-phase biases in a pre-processing step by the method of least squares, and in the following step, the calibrated STEC observations were utilized to generate the VTEC distribution over Europe using a thinplate-splines approximation (Krypiak-Gregorczyk et al., 2017).
In this study, we also solely utilize carrier-phase observations for real-time VTEC modeling. Instead of calibrating the observations in a pre-processing step, the presented approach simultaneously estimates the carrier-phase biases and the B-spline model coefficients by the Kalman filter running in real-time. Raw carrier-phase observations are from GPS, GLONASS, and GALILEO and were acquired in RTCM format from the IGS real-time service.
For VTEC and electron density modeling, B-splines yield a very powerful mathematical tool to deal with the heterogeneous data distribution and data gaps because of its localizing property . For example, Schmidt et al. (2011) showed that a large data gap at a region results in a dramatic global quality degradation in the estimated VTEC map due to artificial oscillating structures when a spherical harmonic (SH) expansion, which is of global support, is employed. Contrarily, in the same study, the artificial structures appear only in the close neighborhood of the data gap when the localizing B-spline expansion replaces the SH representation. The B-splines also support the construction of a multi-scale representation of signals, which paves the way for effective data compression techniques (Lyche & Schumaker, 2000;Schmidt, 2012). Moreover, the coefficients of the B-spline representation resemble global VTEC maps in terms of the spatio-temporal distribution and magnitude, which allow for inferring physical interpretations from the coefficients directly . B-splines were applied to model ionospheric variables in various studies, for instance, modeling regional VTEC (Durmaz & Karslioglu, 2015;Nohutcu et al., 2010;Schmidt et al., 2008;Zeilhofer et al., 2009) and global VTEC (Dettmering et al., 2014;Erdogan et al., 2017Erdogan et al., , 2020Goss et al., 2019;Schmidt et al., 2011), the computation of corrections to a climatological reference model such as IRI (Liang et al., 2015), the estimation of key parameters of the Chapman function for electron density modeling (Limberger et al., 2013), and the combination of different space geodetic techniques for VTEC modeling .
By taking advantage of the B-spline representation, Erdogan et al. (2017) presented an approach to sequentially estimate ultra-rapid global VTEC maps by a Kalman filter, which assimilates ionospheric observations extracted from the GPS and GLONASS via the conventional carrier-phase leveling technique. Later on, the sequential filtering approach was extended with adaptive methods to estimate the covariance matrices of the measurement and dynamics models in run-time  and then applied for regional VTEC modeling . In this study, we incorporate the adaptive approach presented by Erdogan et al. (2020) into the real-time modeling of the spatio-temporal VTEC distribution using B-spline functions. GNSS data is acquired from the IGS real-time service and sequentially processed over time by the filter to simultaneously estimate coefficients of the B-spline representation and the carrier-phase biases.
As stated above, the presented approach estimates the carrier-phase biases for each phase connected arc defined between a GNSS receiver and satellite. As a result, a large number of bias parameters needs to be estimated, which can degrade the stability of the Kalman filter implementation. In order to cope with this problem, the measurement model was extended to incorporate supplementary information derived from a nowcast model. The nowcast model aims to increase the filter stability, supports the model over regions with large data gaps, and provides VTEC maps in case of losing the connection to data providers. The nowcast model is based on the combination of a trigonometric series and an ARMA model. It uses the time series of high-resolution ultra-rapid products of DGFI-TUM (see, e.g., Goss et al., 2019) as input data generated according to Erdogan et al. (2020). This paper is outlined as follows. Section 2 explains the extraction of ionosphere information from real-time carrier-phase observations. In Section 3 the measurement model of the filter is presented. This section comprises the sub-parts: definition of the GNSS observation equations, B-spline representation of global VTEC, nowcasting of B-spline coefficients as supplementary information, and the overall measurement model. Section 4 shows the real-time estimation approach based on the adaptive Kalman filter implementation, the data editing procedures, and the overall procedures from an application point of view. Section 5 explains the validation techniques, and discusses the results. Finally, Section 6 provides the conclusion and future work.

GNSS Real-Time Ionosphere Measurements
For real-time applications, the GNSS data has been usually transmitted from GNSS receivers to users by casters through the Networked Transport of RTCM via Internet Protocol (NTRIP). NTRIP is a data transfer protocol that enables the streaming of GNSS data to stationary or mobile users via Internet (Weber et al., 2005). A list of GNSS NTRIP casters can be accessed, for example, from http://rtcm-ntrip.org/home.html. In this study, the data collected by IGS stations is downloaded through the IGS NTRIP caster (https://www.igs-ip. net/home) which is currently operated in support of the Real-Time IGS Working Group. The open-source software BNC provided by the Federal Agency for Cartography and Geodesy (BKG) as an NTRIP Client is executed to download data from the GNSS constellations (Weber et al., 2016). The downloaded data set includes carrier-phase observations from the L1 and L2 signals of GPS and GLONASS as well as the E1 and E5a signals of GALILEO. The challenging task in real-time modeling is to deal with high-rate raw GNSS data. Most receivers in the IGS caster provide data at a 1 Hz sampling rate. Since accomplishing the data processing and estimation steps in real-time in less than a second might not be possible in practice, a proper down-sampling is typically applied by considering the computational resources. In this study, the down-sampling rate is set to obtain raw data every 10 s.
Following the data download, the real-time data pre-processing step is carried out. First, the ionosphere combination   , , , The ionosphere combinations (Equation 1) for each receiver-satellite pair are collected in an arc container. A cycle-slip detection algorithm based on two frequency carrier-phase data is separately applied to each arc to detect jumps in the signals; see, for example, Subirana et al. (2013). The algorithm utilizes the first-order time difference of the carrier-phase ionosphere combination. The difference of the consecutive observations on the same arc amplifies the jumps and therefore increases the possibility to detect the jumps (Hofmann-Wellenhof et al., 2008). Hence, the difference equation reads  , the data point on the arc at 1 k E t is marked as a jump. Moreover, a data gap threshold of 120 s is also applied to mark time differences exceeding the threshold. It should be noted that a detected cycle-slip is not repaired; instead, a new arc is initiated. After the jump detection and the marking procedures, the ionosphere combinations , s r k E L are converted to TECU unit and then stored in a database.

Ionosphere Measurement Model
As mentioned before, the Kalman filter is carried out as a sequential estimator for real-time VTEC modeling. The following subsections explain the implementation details of the measurement model of the filter designed for real-time VTEC modeling.

Observation Equations
The geometry-free ionosphere combination (Equation 1) eliminates the non-dispersive effects, including the geometric line of sight distance between the satellite E s and the receiver E r as well as the receiver and satellite clock offsets, and can be written as Erdogan et al. (2017). The error term k E e comprises the measurement noise and the other effects such as multi-path and thermal disturbance.
For the sake of clarity, the ionosphere combination Equation 3 is re-written for each of the GNSS satellites as the observation equation

B-Spline Representation of Global VTEC
The global VTEC distribution in Equation 4 is represented by the tensor product expansion where , d t stands for the time-dependent B-spline coefficients . The basis functions along the latitude  E are defined by the polynomial B-spline functions The global VTEC representation (Equation 7) requires a properly handling of constraints to preserve spherical geometry. Two sets of constraint equations, namely the pole equality and the pole continuity, are introduced for the B-spline model representing a function defined on the sphere (Lyche & Schumaker, 2000). Accordingly, the overall constraint equation is a known matrix, MC,k E y equals to zero vector and the vector k E d comprises the B-spline coefficients.
Taking the dominant drivers into account for selecting an appropriate coordinate system is important to model ionosphere phenomena (Laundal & Richmond, 2016). The Earth's magnetic field has an essential role in forming the ionospheric electron distribution and the equatorial ionization anomaly. The Sun is the main driver for the spatio-temporal variation of the photo-ionization process. Therefore, a solar magnetic coordinate system is adapted for the representation of the global VTEC model (Equation 7); for further implementation details, we refer to Erdogan et al. (2020).

Nowcasted B-Spline Coefficients as Supplementary Information
Data gaps due to the heterogeneous GNSS data distribution and the large size of unknown parameters typically make estimators vulnerable to numerical problems, such as ill-conditioning and filter instability. The global real-time approach is extended by supplementary information to keep the Kalman filter numerically stable and to enhance the estimation quality at regions suffering from large data gaps. Furthermore, the supplementary information can be considered as a background ionosphere model providing homogeneous observations to support the real-time model over the oceans where the filter usually suffers from a lack of enough observations. Moreover, it feeds the estimator in case of interruptions in the real-time data streams, for example, due to loss of the connection to data providers.
The supplementary information is obtained in the form of B-spline coefficients from the ultra-rapid global VTEC product of DGFI-TUM with a latency of less than 3 hr (see, e.g., Erdogan et al., 2020;Goss et al., 2019). Because of its dissemination latency, the supplementary information is incorporated into real-time modeling via a nowcast model.
The nowcast model comprises a stochastic and a deterministic part. A linear trend (LT) model and a trigonometric series (TS) are chosen to represent the deterministic part in the time series of B-spline coefficients. An ARMA process consisting of an auto-regressive polynomial of order E p and a moving average polynomial of order E q is executed for the stochastic part. The overall nowcast model reads where the first term refers to LT model consisting of the coefficients 0 E c and 1 E c .  T  T  T  T  T  ,GPS, ,GLO, ,GAL, The corresponding design matrix k E X in Equation 11 is given as is a diagonal matrix, and, in practice, its diagonal components are usually set to very small values to avoid numerical problems (see, e.g., Erdogan et al., 2020;Koch, 1999;Teixeira et al., 2009 where " diag E " is the block-diagonal matrix operator, for example,

Sequential Filtering
The estimation of the ionospheric target parameters is carried out using the Kalman filter (Kalman, 1960) which provides an optimal estimator in terms of, for example, minimum variance estimation when the system of equations is linear and the model uncertainties have a Gaussian distribution; see, for example, Jazwinski (2007). The Kalman filter has been widely preferred for (near) real-time applications due to its recursive nature. It does not require the storage of past data and allows for an immediately updating of the unknown parameters upon the reception of new measurements. The Kalman filter extended by self-learning adaptive approaches can be carried out to estimate the covariance matrices or their parameters during the filter run-time, see Brown and Hwang (2012), Hide et al. (2004), Mehra (1972).
In this study, we employ the adaptive Kalman filter implementation developed by Erdogan et al. (2020) for VTEC modeling. Accordingly, the Kalman Filter consists of a prediction step (time update) and a correction step (measurement update), and these two steps are repeatedly executed over time. The prediction equations between the consecutive epochs k E t and 1  Gelb, 1974;Grewal & Andrews, 2008). Since the VTEC distribution is presented in the solar magnetic coordinate system, the effect of the Earth's diurnal motion in the temporal variation of the B-spline coefficients is mitigated due to the Sun alignment (see, e.g., Erdogan et al., 2017). The solar alignment allows the ionosphere to be treated as stationary (Hansen, 1998). Moreover, the carrier-phase biases, including the ambiguity biases and IFBs are stable over time. Therefore, the ionospheric variables in the time domain can be represented by a simple stochastic model, that is, the random walk process; see, for example, Mannucci et al. (1998). Accordingly, the state transition matrix k E F is set to an identity matrix, that is, are diagonal matrices, and they refer to the process noise covariance matrices of the B-spline coefficients and carrier-phase biases, respectively. Each diagonal component of , w k d E Σ varies in time with respect to the magnitude and spatial distribution of VTEC, which is computed in a self-adaptive manner during run-time of the filter as described by Erdogan et al. (2020). The matrix is set to an experimentally defined value to compensate the dynamic model errors and to keep the filter sensitive to observations.
The subsequent measurement update step of the Kalman Filter is executed following the time update step by considering the new information described by the measurement model (Equation 11) to correct the predicted state vector and covariance matrix obtained from the Equations 16 and 17. In this study, an alternative version of the update equations is employed as The validity of each carrier-phase bias is checked before the measurement update step at every filter cycle. The cleaning and initialization procedures are executed when broken carrier-phase biases are detected or new biases are introduced.

Overall Real-Time Modeling Algorithm
The flowchart in Figure 1 shows the different steps of the presented real-time modeling approach. The overall procedure can be classified into three main categories labeled by different colors in Figure 1; ultra-rapid VTEC product generation, nowcasting of the global VTEC, and real-time VTEC modeling.
First, hourly GNSS raw data is downloaded from the data centers in RINEX format. Ionosphere observations are extracted from the raw data in the ultra-rapid data pre-processing step. The observations are incorporated into the global VTEC model by the adaptive Kalman filter implementation to generate the ultra-rapid VTEC product labeled as"othg"; see Erdogan et al. (2020) for the computation and implementation details and Goss et al. (2019) for the naming of the product. Following the ultra-rapid product generation, after each hour, the coefficients of the nowcast model (Equation 9) are estimated for each of the global B-spline coefficients. The nowcast model is then executed to obtain the extrapolated B-spline coefficients at the present time (nowcasting). These nowcasted B-spline coefficients are used to feed the real-time VTEC model. In parallel to the ultra-rapid product generation and the nowcasting processes, real-time GNSS data are continuously downloaded in the RTCM format. As mentioned earlier, 10 s of data sampling are applied to download raw GNSS data. In the real-time data processing step, carrier-phase observations are extracted from the raw data set. Phase-continuous arcs between receivers and satellites are constructed by checking cycle slips and data gaps. Then the ionosphere combination of carrier-phase observations (Equation 3) are computed.
In the next step, the real-time VTEC modeling, the validity of each carrier-phase bias in the state vector is checked in the data editing process before incorporating the observations into the filter; outdated bias parameters are removed, and new biases are introduced into the filter state vector and its covariance matrix. After that, the ionospheric target parameters consisting of the B-spline coefficients and the carrier-phase biases are estimated sequentially using the real-time adaptive Kalman filter. It should be noted that the filter temporal step size is set to 30 s by taking computational resources into account, though the data download rate and the pre-processing step size are set to 10 s.
In the final step, the VTEC product generation in the desired form, for example, grid VTEC maps in IONEX format, are created from the estimated B-spline coefficients (Erdogan et al., 2017;Goss, Hernández-Pajares, et al., 2020).
In the ultra-rapid and real-time models, the criteria taken into account for selecting the resolution levels are essentially based on the distribution of the input data, the computational load, and the desired level of smoothness for the B-spline representation embedded into a Kalman filter (Erdogan et al., 2017. Accordingly, in this paper, the B-spline level values are set to  , leading to a total number of 816 time-dependent B-spline series coefficients updated at each cycle of the filter. The chosen values refer to maximum spherical harmonics degrees of 33 in latitude direction and 12 in longitude direction (Goss et al., 2019).
In the current implementation, one cycle of the sequential modeling approach involving the real-time data download, pre-processing, filtering, and product creation is accomplished in less than a minute.

Results and Discussion
The quality assessment and validation of the VTEC products derived by the proposed real-time approach are carried out using (a) the dSTEC analysis and (b) altimeter VTEC comparisons. The dSTEC analysis offers a precise method for validating VTEC over areas with installed GNSS receivers and has been widely used for quality assessments of VTEC products generated by various analysis centers; see, for example, Feltens et al. (2011), Hernández-Pajares et al. (2017, Rovira-Garcia et al. (2015). The abbreviation dSTEC refers to difference of a STEC observation with respect to a reference STEC with highest elevation angle on the same phase-continuous arc. VTEC observations acquired from the altimetry missions allow for validating and comparing VTEC maps over the oceans and have been in use for decades to assess VTEC products provided by the IGS analysis centers (Hernández-Pajares et al., 2009;Roma-Dollase et al., 2018).
The test data set comprises GNSS data from real-time streams collected between October 17 and 21, 2019 for the global VTEC modeling. Altimetry VTEC from the Jason-3 mission and additional GNSS data for the dSTEC analysis were also collected for the same period to validate and evaluate the results.

Internal Evaluation of VTEC Models
The growing number of receivers providing real-time and hourly data enhance the data distribution spatially. The GNSS receivers are usually operated as part of a network. The networks are designed to provide data covering the globe, such as the IGS network, or a region such as the EUREF network. The quality of the data distribution depends on the geographical distribution and density of receivers in a network.
The top and mid panels of Figure 2 show snapshot VTEC maps with the ionospheric pierce points of the measurements used for the ultra-rapid and real-time models, respectively. Their difference maps are also illustrated in the bottom panels. The maps on the left refer to 10:00 UTC at October 19, 2019, whereas the right side is for 18:00 UTC of the same day. The ultra-rapid approach uses hourly data collected from the IGS network, and its data distribution is enhanced with an additional hourly data set from manually selected receivers belonging to the UNAVCO and the EUREF networks; see the top panels of Figure 2. For example, additional data from the UNAVCO network mainly cover the Southern United States, Central America, Mexico, and Alaska. In the current implementation, the real-time model is only fed by the data from the receivers contributing to the IGS real-time service. However, the real-time model has more data over Brazil compared to the ultra-rapid model. The receivers in this region contribute to the IGS real-time service; but, no data was available in the hourly data set from IGS data centers.
As shown in the difference maps in Figure 2, the deviations between the ultra-rapid and real-time models are dependent on the VTEC magnitude and the data distribution. For instance, considerable deviations exist over Central America (for 18:00 UTC) and Indian (for 10:00 UTC) regions. A considerable deviation during a high geomagnetic activity between the two models is also expected over the Siberian region due to the data gap in the real-time model; however, it does not become significant for the test period. The real-time model performs well compared to the ultra-rapid model over Europe, where both approaches are fed with dense data, and the differences between the two models are close to zero.
It should also be mentioned that the amount of extracted ionospheric information from raw data also depends on the data pre-processing method. The ultra-rapid approach utilizes the phase-leveling method and requires a minimum of a 30-min length of consecutive code and carrier-phase measurements on the same arc between a receiver and a satellite. The measurements that are part of the short arcs are discarded in ultra-rapid modeling. However, the real-time data pre-processing, which is only based on carrier-phase measurements, does not require such a long arc length in the real-time data pre-processing step. The initialization of a phase-continuous arc is set to a minimum of 2-min, allowing for extracting more ionospheric measurements from GNSS satellites.

External Quality Assessment of the Real-Time VTEC Products
There are numerous ongoing experimental efforts for providing real-time products, such as real-time clock and orbit corrections for GPS and GLONASS as well as VTEC maps, provided by IGS; see, for example, Li et al. (2020) for combined VTEC maps. The real-time VTEC products of the Universitat Politècnica de Catalunya (UPC), Spain, and the Chinese Academy of Sciences (CAS), China, in the IONEX format are used in the comparisons to assess the relative performance of the VTEC maps generated in this study. The VTEC products of UPC and CAS will respectively be labeled "uadg" and "rt_casg", and the VTEC maps produced using the developed approach will be labeled "orhg" in the remainder of the paper. The ultra-rapid VTEC product of DGFI-TUM with the label "othg" and the rapid VTEC product of UPC with the label "uqrg", which provide both high spectral resolution, accuracy and temporal resolution, are also included in the analysis for the evaluation of the results. It should be noted that the real-time data dissemination format can change the quality of VTEC products. For example, a conversion of the VTEC product into a spherical harmonic representation complying with the RTCM format drastically decreases the accuracy of high-resolution products since the RTCM 3.x format allows VTEC representation by spherical harmonic expansions only up to a degree 16 (Goss, Hernández-Pajares, et al., 2020). Therefore, VTEC grids in the IONEX format were preferred to avoid such quality degradation for a fair evaluation of the results. Figure 3 shows the geographical distribution of the GNSS receivers used in the dSTEC analysis. The receivers are selected to reveal the performance of the products according to the VTEC activity showing significant variations with respect to geographic locations and time, such as the Equatorial Ionization Anomaly. The RMS deviations derived from the dSTEC analysis for each test receiver are depicted in Figure 4. The  . RMS values from the differential slant total electron content analysis at the Global Navigation Satellite System stations depicted in Figure 3. The results refer to the data sets covering the days between October 17 (DOY 290) and October 21 (DOY 294) of 2019. The label "orhg" stands for the presented real-time approach. The labels "uadg" and "rt_casg" refer to the real-time UPC and CAS products, "othg" and "uqrg" stand for the ultra-rapid DGFI-TUM product and the rapid UPC product.
legend in the right panel shows the average values of the RMS deviations for each analysis center. The average RMS deviations in Figure 4 are ranging from 0.76 to 1.22 TECU. The values reveal that the quality of the real-time VTEC product "orhg" is in line with those from the other analysis centers. The accuracy of "orhg" with RMS deviation of 1.05 TECU is slightly better than the real-time VTEC products of "uadg" and "rt_casg" during the test period. The biases, which reveal the systematic variations between the VTEC solutions, are not illustrated in the figure but range from −0.05 to 0.02 TECU for the test period. The real-time product "orhg" has an average bias of −0.01 TECU. The (ultra) rapid VTEC solutions "othg" and "uqrg" show superior performance than the real-time solutions by means of the RMS deviations. This performance can be associated with better data distribution. Noticeable quality degradation for all products is visible for the GNSS receivers located around the geomagnetic equator and is significantly larger for the real-time products than the rapid products; see the RMS deviations for the receivers NIUM, KOKB, BOGT, YKRO, IISC, and GUAM in Figure 4. This degradation is usually attributed to a lack of enough data coverage and VTEC activity at large magnitudes, which occurs at regions around the geomagnetic equator.
Furthermore, VTEC observations obtained from the dual-frequency altimeter of Jason-3 mission were used to evaluate the performance of generated VTEC maps over the oceans. Although VTEC data from satellite altimetry typically suffer from calibration biases that can reach up-to several TECU ) and pronounced noise (Hernández-Pajares et al., 2017), it provides a direct and independent way to evaluate VTEC maps (Roma-Dollase et al., 2018;Li et al., 2020). A pre-filtering is usually applied to smooth noisy raw altimetry observations (see, e.g., Li et al., 2020;Hernández-Pajares et al., 2017). In this study, a median filter with a window size of 20 s was carried out to smooth the raw data and discard outliers (see, Erdogan et al., 2017). Since the temporal sampling of Jason-3 data is higher than the VTEC products, a tri-linear spatio-temporal interpolator was performed in the Sun-fixed coordinate system to compute VTEC values from the maps for the corresponding time and location of the altimetry observations. The RMS of the differences between altimetry VTEC observations and the corresponding VTEC model values derived from products of each analysis center are used as a metric for the relative quality assessment. The RMS values of the differences computed for each hour for the test period are depicted in Figure 5. The results reveal that the average RMS values, given in the legend of the figure, for each product vary in a range between 2.5 and 3.2 TECU. The RMS value of the real-time solution "orhg" for the corresponding data sets is 3.2 TECU, and it is in compliance, considering the precision level of the altimetry observations, with those of the real-time products "uadg" and "rt_casg", which are 2.5 and 2.9 TECU, respectively. Additionally, among the real-time products, "uadg" has the smallest average bias of 0.5 TECU, whereas the "orhg" and "rt_casg" solutions have 0.8 TECU for the test period. In the current implementation, the presented approach focuses on continental areas, where the model accuracy is superior compared to the marine areas.
The slight under-performance of the B-spline approach over the oceans according to the RMS deviations can be attributed to the compact support of its basis functions compared to spherical harmonics, which provide an averaging effect over the globe. However, the accuracy over the marine areas can be improved by incorporating additional GNSS receivers established on islands and coasts from reliable networks. Furthermore, by adapting methods based on interpolations into ionosphere modeling, significant performance improvements over large data gaps can be achieved. For instance, the rapid "uqrg" product of UPC makes Figure 5. RMS values computed for each hour from Jason-3 altimetry comparisons. The results refer to the data sets covering the days between October 17 (DOY 290) and October 21 (DOY 294) of 2019. The label "orhg" stands for the presented real-time approach. The labels "uadg" and "rt_casg" refer to the real-time UPC and CAS products, "othg" and "uqrg" stand for the ultra-rapid DGFI-TUM product and the rapid UPC product. use of Kriging interpolation and the real-time product "uadg" employs a method called Atomic Decomposition Interpolator of GIMs (ADIGIM); see, for example, Yang et al. (2021).

Conclusions
This paper focuses (a) on embedding a B-spline representation of the global VTEC into a Kalman filter running in real-time, (b) on utilizing exclusively carrier-phase observations acquired from GPS, GLONASS, and GALILEO via the real-time data streams, and (c) on simultaneously estimating the unknown target parameters comprising the B-spline coefficients and the carrier-phase biases.
According to Erdogan et al. (2017Erdogan et al. ( , 2020, incorporating the global VTEC B-spline representation  into an adaptive Kalman filter results in an efficient recursive data ingestion framework for ultra-rapid modeling. This adaptive modeling approach was extended to real-time in the frame of this paper to exploit the efficient heterogeneous data handling capability and the localizing property of B-splines. Besides, the localizing feature of the B-spline functions results in highly sparse structures in the design matrix of the filter measurement model. By considering the large scale and sparsity of the problem, sparse matrix techniques were introduced for the matrix/vector operations in the filter measurement update step. This adaptation significantly decreased the filter computation time, which is a critical matter for real-time VTEC modeling. In the current implementation, one cycle of the real-time processes, including the filtering step, takes less than a minute. The measurement model of the presented approach is based on the geometry-free combination of carrier-phase observations, which is composed of an ionospheric signal term and a combined bias term. Instead of eliminating them in a data pre-processing step, all carrier-phase biases are estimated along with the B-spline coefficients in run-time. The advantage of the approach is that it does not require the use of code observations compared to the conventional carrier-phase leveling technique. Furthermore, the geometry-free combination does not require precise auxiliary data (such as precise orbit and clock products) contrary to methods for extracting ionosphere parameters from PPP-based solutions. Low latency and moderate-quality ultra-rapid orbits, which include predicted satellite positions and clock biases, are used to compute, for instance, ionosphere pierce-points and elevation angles of observations. On the other hand, these advantages come at the expense of estimating a large number of unknown carrier-phase biases, which are about 3,500 parameters in addition to the 816 B-spline coefficients at each epoch for the test period. However, the ultra-rapid modeling approach based on the phase leveling method requires the estimation of about 600 DCBs in addition to the B-spline coefficients. The large size of the state vector and large data gaps can result in an unstable filter implementation leading to quality degradation of VTEC products. Supplementary information derived from the nowcast model was incorporated into the filter to cope with such a problem.
The performance assessment of the implemented real-time VTEC modeling was carried out via the dSTEC analysis and Jason-3 altimetry VTEC comparisons. The results show that the quality of the real-time VTEC product "orhg" is slightly better in terms of average RMS of 1.05 TECU than the other real-time products, ranging from 1.12 to 1.22 TECU, according to the dSTEC analysis during the test period. Besides, the RMS of 3.2 TECU for the real-time "orhg" product is in close agreement with those, ranging from 2.5 to 2.9 TECU, used in the evaluations with respect to the altimetry VTEC analysis providing a noisy but reliable assessment method over the oceans. Nevertheless, as mentioned before, further performance improvement over the oceans for the product "orhg" is expected by incorporating new receivers from additional GNSS networks and implementing methods based on interpolation. Moreover, considering the rapid advances in technology, if data sets from other space geodetic techniques such as DORIS and satellite altimeter can be provided in real-time, these data can be included in VTEC modeling. This issue will be handled as a follow-up study.
The results encourage further research to extend the measurement model of the real-time Kalman filter to consider additional satellite constellations such as BeiDou. The current implementation of the approach utilizes data from the real-time IGS network. In order to obtain a better data distribution, complementary data sources from reliable real-time networks, for instance, the UNAVCO network, will be introduced. Rapidly enlarging regional networks providing a dense set of GNSS receivers, for instance, the real-time European network (EUREF), strongly encourage us to apply the presented approach to real-time regional VTEC modeling. The evaluations show that the presented approach has rather promising outcomes during the selected period (October 17-21, 2019). Nevertheless, further tests covering a longer time span and including solar events will be performed as part of a near-future study.

Data Availability Statement
Accessing the real-time GNSS resources from the IGS service (http://igs-ip.net/home) was possible by the open-source BKG Ntrip Client (BNC) software (see, Weber et al., 2016) after the completion of the registration procedure via the web page "https://register.rtcm-ntrip.org/cgi-bin/registration.cgi". The hourly GNSS data used for the dSTEC analysis can be obtained from the archives of IGS data centers, such as CDDIS (via the web portal "https://cddis.nasa.gov/"). The real-time and rapid VTEC GIMs from UPC and the real-time VTEC GIMs from CAS, used in the validation step, are publicly provided through the data repositories of UPC and CAS; see, Yang et al. (2021) and Li et al. (2020) for further information. The created VTEC GIMs based on the presented real-time approach can be accessed through the research data repository of TUM via the link "https://doi.org/10.14459/2021mp1621980". Altimetry data provided by NOAA is available from the website at "https://www.ncei.noaa.gov/data/oceans/jason3/".