# GipsyX/RTGx, a new tool set for space geodetic operations and research

Willy Bertiger, Yoaz Bar-Sever, Angie Dorsey, Bruce Haines, Nate Harvey, Dan Hemberger, Michael Heflin, Wenwen Lu, Mark Miller, Angelyn W. Moore, et al.

► To cite this version:

Willy Bertiger, Yoaz Bar-Sever, Angie Dorsey, Bruce Haines, Nate Harvey, et al. GipsyX/RTGx, a new tool set for space geodetic operations and research. Advances in Space Research, 2020, 66, pp.469-489. 10.1016/j.asr.2020.04.015. insu-03584748

HAL Id: insu-03584748

https://insu.hal.science/insu-03584748v1

Submitted on 16 May 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# GipsyX/RTGx, a new tool set for space geodetic operations and research

Willy Bertiger<sup>a,*</sup>, Yoaz Bar-Sever<sup>a</sup>, Angie Dorsey<sup>a</sup>, Bruce Haines<sup>a</sup>, Nate Harvey<sup>a</sup>, Dan Hemberger<sup>a</sup>, Michael Heflin<sup>a</sup>, Wenwen Lu<sup>a</sup>, Mark Miller<sup>a</sup>, Angelyn W. Moore<sup>a</sup>, Dave Murphy<sup>a</sup>, Paul Ries<sup>a</sup>, Larry Romans<sup>a</sup>, Aurore Sibois<sup>a</sup>, Ant Sibthorpe<sup>a</sup>, Bela Szilagyi<sup>a</sup>, Michele Vallisneri<sup>a</sup>, Pascal Willis<sup>b,c</sup>

Received 8 August 2019; received in revised form 8 April 2020; accepted 13 April 2020

Available online 22 April 2020

## Abstract

GipsyX/RTGx is the Jet Propulsion Laboratory’s (JPL) next generation software package for positioning, navigation, timing, and Earth science using measurements from three geodetic techniques: Global Navigation Satellite Systems (GNSS), Satellite Laser Ranging (SLR), and Doppler Orbitography and Radiopositioning Integrated by Satellite (DORIS); with Very Long Baseline Interferometry (VLBI) under development. The software facilitates combined estimation of geodetic and geophysical parameters using a Kalman filter approach on real or simulated data in both post-processing and in real-time. The estimated parameters include station coordinates and velocities, satellite orbits and clocks, Earth orientation, ionospheric and tropospheric delays. The software is also capable of full realization of a dynamic terrestrial reference through analysis and combination of time series of ground station coordinates.

Applying lessons learned from its predecessors, GIPSY-OASIS and Real Time GIPSY (RTG), GipsyX/RTGx was re-designed from the ground up to offer improved precision, accuracy, usability, and operational flexibility. We present some key aspects of its new architecture, and describe some of its major applications, including Real-time orbit determination and ephemeris predictions in the U.S. Air Force Next Generation GPS Operational Control Segment (OCX), as well as in JPL’s Global Differential GPS (GDGPS) System, supporting User Range Error (URE) of <5 cm RMS; precision post-processing GNSS orbit determination, including JPL’s contributions to the International GNSS Service (IGS) with URE in the 2 cm RMS range; Precise point positioning (PPP) with ambiguity resolution, both statically and kinematically, for geodetic applications with 2 mm horizontal, and 6.5 mm vertical repeatability for static positioning; Operational orbit and clock determination for Low Earth Orbiting (LEO) satellites, such as NASA’s Gravity Recovery and Climate Experiment (GRACE) mission with GRACE relative clock alignment at the 20 ps level; calibration of radio occultation data from LEO satellites for weather forecasting and climate studies; Satellite Laser Ranging (SLR) to GNSS and LEO satellites, DORIS-based and multi-technique orbit determination for LEO; production of terrestrial reference frames and Earth rotation parameters in support of JPL’s contribution to the International Terrestrial Reference Frame (ITRF).

© 2020 COSPAR. Published by Elsevier Ltd. This is an open access article under the CC BY license (http://creativecommons.org/licenses/by/4.0/).

Keywords: Multi-technique space geodesy; Orbit determination and satellite clocks; Precise Point Positioning (PPP); GNSS; SLR; DORIS

MSC 2010: 00-01; 99-00

<footnote><sup>a</sup> Jet Propulsion Laboratory, California Institute of Technology, 4800 Oak Grove Drive, Pasadena, CA 91109, United States</footnote>

<footnote><sup>b</sup> Université de Paris, Institut de physique du globe de Paris, CNRS, IGN, F-75005 Paris, France</footnote>

<footnote><sup>c</sup> ENSG-Géomatique, IGN, F-77455 Marne-la-Vallée, France</footnote>

<footnote>* Corresponding author.</footnote>

<footnote>E-mail address: william.bertiger@jpl.nasa.gov (W. Bertiger).</footnote>
# 1. Introduction

Over the past thirty years, Global Navigation Satellite Systems (GNSS), starting with GPS, have become ubiquitous in our daily lives. A growing number of scientific, industrial, and public safety applications, spanning geodesy, geophysics, seismology, atmospheric sciences, weather forecast and climatology, time and frequency transfer, have come to depend on precision modeling of all aspects of GNSS at the cm and mm accuracy. Very few software packages exist with the capabilities to obtain these highest accuracy results, keeping up with the latest GNSS models and data processing techniques. Fewer still are capable of also modeling and processing data from other space geodetic techniques such as SLR (Satellite laser ranging), DORIS (Doppler orbit determination and radio-positioning integrated by satellite), or VLBI (Very long baseline interferometry). Recently, JPL has embarked upon a major overhaul of its widely-used space geodetic data analysis software packages, GIPSY-OASIS and RTG, to derive GipsyX/RTGx, with many advanced new capabilities and quality features. This paper is intended for the users and developers of space geodesy tools.

We present the design and some implementation notes for a new, but already widely used software set, and describe its performance in a broad range of applications spanning orbit determination and positioning operations, geodesy, and remote sensing. We first give some historical context, which influences any large software and research project design. We then present the overall design of the software, and sample use cases with tests of accuracy and precision.

The Jet Propulsion Laboratory (JPL) has a long history in space geodesy and precision orbit determination (POD) for Earth orbiting satellites, and GipsyX/RTGx is the 4th major redesign of JPL's GNSS data analysis software. For an overview of POD and measurement processing see the books (Tapley et al., 2004; Bierman, 1977). Born from Very Long Baseline Interferometry (VLBI) data analysis expertise and software (Sovers and Fanselow, 1987), the first JPL GPS data models were coded during 1984–1985 in the GPSOMC FORTRAN software (Sovers and Border, 1990). Immediately thereafter we began development of a more comprehensive software set, the GPS Inferred Positioning System and Orbit Analysis Simulation Software package (GIPSY-OASIS, forever nicknamed GIPSY) (Wu et al., 1986), also coded in FORTRAN, to support the early campaign-style GPS orbit determination (GiG'91 for example (Melbourne et al., 1993)) and positioning experiments (Blewitt, 1989). That very early version of GIPSY should be considered a prototype for the more formally designed, enhanced, and documented version (Wu et al., 1990), written to support the first precision GPS flight experiment on the TOPEX/Poseidon altimetry mission (Bertiger et al., 1994) as well as the continuous GPS orbit determination operations at JPL that are the cornerstone for all precision GIPSY applications.

With the maturation of the GPS constellation, which became fully operational in 1994, real-time operational applications, such as differential positioning, started to emerge. The state-space formulation of real-time differential systems, today a universally accepted practice also known as wide area differential GPS (WADGPS) or global differential GPS (GDGPS), was conceived and patented at JPL (Yunck et al., 1989), and led to the creation of the Real Time GIPSY (RTG) software. RTG ported most of GIPSY's precise satellite and signal models from FORTRAN into C, but redesigned the software architecture, and in particular the estimation filter, for efficient real-time operations. Successful early demonstration of RTG (Whitehead et al., 1998) has led the United States Federal Aviation Authority's (FAA) to adopt JPL's state-space algorithms for its Wide Area Augmentation System (WAAS), with RTG as the prototype software for this critical aviation infrastructure, which has since been adopted by several other countries as well (Bertiger et al., 1997). In 2000, RTG became core software for JPL's GDGPS system, as well as for other commercial systems.

As with any JPL-developed software, GIPSY and RTG are owned by the California Institute of Technology (Caltech). Licensing requests are available from [https://gipsy-oasis.jpl.nasa.gov]. With hundreds of issued licenses to academia (which are free of license fees), industry, and government, GIPSY and RTG are among JPL and NASA's most licensed software, and the recipient of multiple industry and NASA awards, including the 2000 NASA Software of the Year Award for RTG, and 2004 Space Technology Hall of Fame induction for GIPSY and RTG (as "Precision GPS Software System") [https://spinoff.nasa.gov/Spinoff2008/award_winners.html].

In 2005 Raytheon (as lead), ITT (now L3Harris Corporation), and JPL teamed up to develop an innovative navigation concept of operation for the U.S. Air Force next generation GPS operational control segment, known as OCX. The OCX project is a complete overhaul of the GPS operational control segment, including new architecture, infrastructure, hardware, and software required to comply with a set of demanding performance and quality specifications. The team's proposed navigation concept was patterned after JPL's GDGPS System and built upon JPL's proven GIPSY and RTG navigation software technology. In 2010 the team was competitively awarded the $1B + project to develop and deploy OCX [http://insidegnsss.com/raytheon-wins-1-5-billion-gps-ocx-contract/], (Bertiger et al., 2010), and we set out to create JPL's 4th generation GPS data processing software, entitled RTGx. Later, the NASA Space Geodesy Project (SGP) provided additional support to enhance post-processing for space geodesy (Merkowitz et al., 2018). All programs synergistically benefited.

This was an opportunity to modernize the entire geodetic data processing software set at JPL and unify the agile real-time processing of RTG with the post-processing flexibility of GIPSY, all within a new architecture that can sup-
port all of JPL's, NASA's and the Air Force's OCX current and anticipated needs. This unified multi-capability software set was named GipsyX/RTGx (or RTGx/GipsyX, depending on the main use case. The choice of capitalization is intentional).

The development of GipsyX/RTGx was mostly completed by 2016. By 2014 it had replaced RTG as the orbit determination engine for most GDGPS real-time GNSS operations, and by 2017 it had replaced GIPSY as the orbit determination engine of all post-processed GPS orbit determination operations at JPL, including the reprocessing of the entire record of GPS ground tracking data back to the early 1990s (in order to establish a consistent long-term set of orbit and clock solutions in ITRF2014). However, with the ever-changing GNSS landscape and the on-going quest to combine the various geodetic measurement techniques, GipsyX/RTGx is constantly being updated and upgraded.

GipsyX/RTGx, as was GIPSY-OASIS, is a member of a small 'club' of available precision space geodetic data analysis packages, with diverse capabilities, and with a track record of contributing high quality GNSS data analysis products. These include Bernese (Dach et al., 2015), EPOS.P8 (Uhlemann et al., 2015), GAMIT (Herring et al., 2018), GEODYN (Pavlis et al., 2017), GINS (Marty et al., 2011), NAPEOS (Springer et al., 2011), and PANDA (Shi et al., 2010).

Since 2018, GIPSY and RTG are no longer supported by JPL. The large GIPSY licensee community is transitioning to GipsyX/RTGx with the help of JPL-provided classes, online documentation, tutorials and licensing information on (https://gipsy-oasis.jpl.nasa.gov) with over 230 licenses of which 80% are free academic licenses.

## 2. Software design, overview

GIPSY's core computational engines were written in FORTRAN with glue code and automation in C-shell, Bourne-Shell, and Perl. RTG was written in C. Both GIPSY and RTG implemented precise sets of models consistent with the 2010 Conventions of the International Earth Rotation Services (IERS) (Petit and Luzum, 2010). RTG was designed for lean real-time processing, and did not have post-processing capabilities required by many of NASA's science applications and by the geodetic user community. While FORTRAN has been modernized since our original GIPSY writing, it is not currently as widely used as C or C++, and is hardly taught at universities. The U.S. Air Force contract for OCX forbade outright the use of FORTRAN. C does not support object-oriented programming, which we deemed necessary for ease of development and maintenance. We chose, therefore, to write the main computational parts of the new software set in C++. We wanted to avoid the hodgepodge of scripting languages used in GIPSY and chose Python3, which has a large user base and community support in the sciences and engineering, as the sole scripting language. For small problems and identical use cases (for instance precise point positioning

with GPS or low Earth orbit determination with GPS) which are implemented in GIPSY and GipsyX/RTG, run times are nearly identical. Larger problems, that can be handled by GIPSY and GipsyX/RTGx can be faster with the new software, due to optimization and use of parallel computation. Table 1, shows almost a factor of three increase in performance for a moderate size problem, the problem is too small for realizing large increases in speed with multiple cores/CPUs. In addition to new use cases, which are not implemented in GIPSY, there are very large problems that can only be handled in GipsyX/RTGx, simultaneous adjustment of a large gravity field, GNSS orbits, GRACE orbits and additional low Earth orbiters. For use cases that are common to GIPSY and GipsyX/RTGx, we expect almost identical precision and accuracy, but a much improved user interface. It is the preservation of the accuracy and precision, along with the additional use cases that adds new value to space geodesy.

Among the key new capabilities of GipsyX/RTGx are:

* Complete support for most Linux distributions and Mac OS X; tailored installation on embedded flight hardware

* Full GNSS satellite and signal modeling, including GPS, GLONASS, BeiDou, Galileo, and QZSS

* Seamless processing of data from shared memory or files supports both real-time and post-processing operations

* High fidelity data simulator capable of simulating all the data types it can process

* Efficient square root information filter (SRIF) formulation featuring multi-threading and parallelization with Message Passing Interface (MPI)

* Hot start capability with archived filter state and covariance

* Integer GNSS phase ambiguity resolution in post-processing and in real-time

* Inertial measurement modeling to optimally combine GNSS and inertial data in real-time and in post-processing

* Hyperlinked documentation using Doxygen and Sphinx

* User input interface featuring inheritance

* Reference frame determination and time series analysis of Earth station positions.

Table 2 shows the size of the entire GipsyX/RTGx distribution in terms of the number of lines of code. Raw lines are just the total number of lines in all the source files. Physical lines are lines that are not comments or blank. The majority of the C++ code occurs in the single executable module,

<table><caption>Table 1 CPU Time comparison, GIPSY versus GipsyX/RTGx for a moderate size problem, 45 stations, 32 GPS satellites with bias fixing. Unloaded Intel® Xeon® Gold 6130 CPU at 2.1 GHz. GIPSY does not support multiple cores.</caption><thead><tr><th></th><th>CPU time (min)</th></tr></thead><tbody><tr><td>GIPSY</td><td>45</td></tr><tr><td>GipsyX/RTGx 1-core</td><td>16</td></tr><tr><td>GipsyX/RTGx 2-core</td><td>11</td></tr></tbody></table>
<table><caption>Table 2
GipsyX source lines, 1k = 1,000.</caption><thead><tr><td>Language</td><td>Physical Lines</td><td>Raw Lines</td><td>Comment Lines</td></tr></thead><tbody><tr><td>C++</td><td>202.7k</td><td>353.7k</td><td>103.3k</td></tr><tr><td>Python3</td><td>35.6k</td><td>60.4k</td><td>14.6k</td></tr></tbody></table>

“rtgx”, which models the GNSS data and other tracking data types, and estimates parameters such as station and satellite position. We describe the structure of “rtgx” below. There are approximately 140 individual executable programs in GipsyX/RTGx. Included are many small utilities, such as time format manipulation, time series and data manipulation, trend analysis, data archive access and plotting utilities. Most users will only use a small set of these. In addition to the tracking data analysis tools, there is a set of Python3 executables dedicated to time series analysis of positions on the Earth and reference frame determination as well as many other small utility tools (for example, to convert between file formats or timing formats).

We have attempted to apply in this development the many lessons learned from the previous iterations of GIPSY and RTG. GIPSY was written by many people under a fairly loose set of coding rules. The main computational pieces of GIPSY were divided into several independent FORTRAN modules (orbit integration, signal/Earth model, filter, smoother...). There was considerable inconsistency in the input to these pieces, which were mostly FORTRAN name-lists. Common information across these FORTRAN executables was not consistently named or configured. To mitigate the complex configuration of the various GIPSY modules, we wrote a Perl script (gd2p.pl) to hide the complexity and act as the interface for most of our external users who were doing precise point positioning (PPP) or precise orbit determination (POD) of low Earth Orbiters using GPS (Bertiger et al., 2010). GipsyX/RTGx improves upon GIPSY in two ways here. One, there is a uniform input format, referred to as a ‘tree’ input defined in detail below, which has a general structure similar to YAML. We did not use YAML because, at the time of our development, it did not have some object-oriented inheritance properties or arbitrary script execution that would allow human readable, compact input for GNSS constellations where many of the satellites have identical properties. Second, instead of writing many individual executable modules as in GIPSY, the object-oriented C++ allowed us to write a single main executable, rtgx.

For version control, we are using `svn`, subversion (Pilato et al., 2008) chosen before distributed version control systems such as `git` (Chacon and Straub, 2014) were more popular. We believe that documentation is best maintained inside the source code; thus, we chose to use `Doxygen` for C++ code, and `Sphinx` for Python. `Doxygen` and `Sphinx` both output `html` documentation and it is easy to link the two sets of documentation. The mathematical description for the radio-metric signals, solid tides, Earth orbiting satellite forces, ordinary differential equation integrator,

etc., are best written in `LATEX`. Our html documentation includes a pdf containing the mathematical description, generated with `LATEX`, with links to the source code documentation.

Having a single main executable, nearly all use cases can be constructed as a very thin wrapper around a single command line:

```
rtgx myInput.tree
```

To automate the most common tasks carried out by our user community, namely precise static point positioning with single receiver ambiguity resolution, we have provided a Python3 wrapper, where the only required input is a `RINEX2` (Receiver Independent Exchange Format) or `RINEX3` GNSS data file.

```
gd2e.py -rnx my.rnx
```

This is the analog of `GIPSY`'s `gd2p.pl` automation script, and it contains many enhancements over its predecessor in terms of ease of use and documentation.

### 3. User interface - the input tree

Several of our executables, especially those requiring more complicated and finely controllable user inputs, are configured by a single interface called a tree. A tree is a text file with a hierarchical structure identified using Python-like indentations that consist of roots (the highest level), branches, and leaves (the lowest level; the leaves of the tree are the parts of a tree with no branches and are generally the ones that contain specific data). See electronic supplement 1 for details.

### 4. Main C++ Software Modules/Classes

The computationally intensive software breaks into several broad categories:

* Data editing

* Orbit Integration and force models

* Signal Models

* Earth Models

* Filter, optimized fit to the linearized model forward in time

* Smoother, optimized fit to the linearized model over all time.

All of these items are implemented as C++ classes or functions in a single executable module, `rtgx`, except for data editing, which is provided by the `gde` module.

#### 4.1. gde Data Editor

GNSS data must be edited for phase breaks and gross outliers before it is processed to fit model parameters.
The GNSS Data Editor (gde), is a C++ code that implements two distinct algorithms: turbo-edit (Blewitt, 1990), and simple continuity checks of linear combinations of phase data based on low-degree polynomial fits. Turbo-edit looks at averages of pseudorange – phase measurements to detect jumps in the phase measurements. It is well suited to low rate data, such as the large number of RINEX files typically recorded at a 30-s rate. When data rates are high, for example sampled every second as often seen in real-time applications, it is possible to detect phase breaks simply by monitoring the change in time of differences in phase observations on two frequencies. Such a difference cancels out receiver and transmitter clock jumps. Removing a low order time polynomial from these differences effectively removes the slow-changing ionospheric signal in this difference. Single differences of dual-frequency combinations across two satellites are also available as an option.

## 4.2. Orbit Integration and force models

To model the path of a satellite in orbit about the Earth, we implement precise force models that include gravitational forces and surface forces on the satellite. These forces may contain model parameters. Given initial values for the satellite epoch state (position and velocity at the initial epoch time) and parameter values, the second order differential system of equations for the time evolution of the satellite state and the partial derivatives of the satellite epoch state with respect to the model parameters are integrated numerically forward in time. Force models represented in GipsyX/RTGx include the Earth’s gravity field as a spherical harmonic expansion to arbitrary degree and order including time variability of these coefficients, gravitational effects of Earth solid tides, pole tides and ocean tides, first order General Relativistic acceleration and atmospheric drag (DTM-2000, (Bruinsma et al., 2003)). Conical models for solar eclipse by the Earth and Moon are implemented. Generic Earth albedo (Knocke, 1989; Knocke et al., 1988) and solar radiation pressure (Milani et al., 1987) models are available through a flexible tree input interface which enables the specification of arbitrary sets of panels along with their orientation and reflective properties. More specific satellite radiation pressure models, including some empirical GPS solar radiation pressure models specially developed at JPL (Bar-Sever and Kuang, 2005; Sibthorpe et al., 2011), may also be input either as arbitrary Fourier series, or bi-linearly interpolated tables dependent on Sun azimuth and elevation angles in the satellite body-fixed frame.

In order to accurately compute surface force models, spacecraft orientation is required, and may be input in a variety of ways, including time-series of attitude quaternions and, especially useful for simulations, dynamic pointing specified via mathematical operations on vectors from the satellite to the Sun and Earth, and the satellite’s velocity with a defined input calculus. Attitude may also be

selected from a number of specific models including those for GPS (Bar-Sever, 1996; Bar-Sever et al., 1996), GNONASS (Dilssner et al., 2011), Galileo (https://www.gseuropa.eu/support-to-developers/galileo-iov-satellite-metadata) and configurable yaw steering and orbit normal modes for BeiDou geostationary and inclined geostationary satellites (Hugentobler and Montenbruck, 2017). For other orbiting satellites including altimetry missions, attitude may be selected from models or quaternions.

We have implemented a general-purpose ordinary differential equation (ODE) solver which we use to obtain numerical solutions of the time evolution of the satellite state. This ODE solver includes multistep methods (Hairer et al., 1993) as well as high-order embedded Runge-Kutta methods (Prince and Dormand, 1981). The solver employs typical local-truncation-error-estimate adaptive step size techniques. In order to allow for discontinuities in the integrated source functions (e.g., the solar radiation pressure transition between shadow regimes, and satellite thruster firings have discontinuous time derivatives) and yet still employ high-order integration techniques, the ODE solver employs a model-specific discontinuity detection logic, following the “gstop-function” approach specialized at JPL by Fred Krogh for integrating spacecraft trajectories (Krogh, 1972) https://mathalacarte.com/fkrogh/pub/dxrkmeth.pdf.

## 4.3. Signal model

The basic signal model is the time of propagation from a transmitter to a receiver. For GNSS range measurements recorded at receiver time, $\tilde{t}_r$, determined by the receiver’s clock and transmitted at time $\tilde{t}_t$, we model the measured range as

$$R = c(\tilde{t}_r - \tilde{t}_t) + d_{\text{trop}} + d_{\text{iono}} \quad (1)$$

where c is the speed of light, $d_{\text{trop}}$ is a delay due to troposphere and $d_{\text{iono}}$ is the delay due to the ionosphere. For many GNSS measurements the first order delay of the ionosphere is removed by a dual-frequency combination of the data (Parkinson et al., 1996) and a model for higher order effects (Kedar et al., 2003). The troposphere is modeled as a delay at zenith plus gradient parameters (Bar-Sever et al., 1998; Böhm and Schuh, 2004; Böhm et al., 2006; Böhm et al., 2015). The difference in the time of reception, $\tilde{t}_r$, and the time of transmission, $\tilde{t}_t$, in Eq. (1) can be modeled as

$$\tilde{t}_r - \tilde{t}_t = \tilde{t}_r - \tilde{t}_r + \tilde{t}_r - t_r + t_r - t_t + t_t - \tilde{t}_t + \tilde{t}_t - \tilde{t}_t \quad (2)$$

where $\tilde{t}_r - \tilde{t}_r$ is the difference between the time on the receiver clock and proper time(time clock would have read with no errors) and $\tilde{t}_r - t_r$ is the difference between proper time
and coordinate time (General Relativistic effects, (Moyer, 1971; Thomas, 1975; Moyer, 1981)). The last two lines of Eq. (2) contain similar differences for the transmitter. The middle term, difference in coordinate time, $t_r - t_t$, is modeled as the geodesic distance (General Relativity with the only mass being the Earth’s point mass) between the transmitter and receiver phase centers. The fit to measurements of the range is optimized by linearizing about the nominal parameter values in Eqs. (1), (2), first analytically taking partial derivatives with respect to the parameters and then in the case of satellites using the numerically calculated partials in Section 4.2.

For multi-GNSS, Eq. (1), must be modified due to different delays of ranging codes through the receiver for different constellations. The code adds a bias parameter for each receiver by constellation (Odijk and Teunissen, 2013). One constellation must be held as reference to prevent singularity. Typically we hold GPS as reference, but it is arbitrary.

The model for GNSS phase data is identical to range, except two terms must be added to Eq. (1),

$$ \phi = c(\tilde{t}_r - \tilde{t}_t) + d_{\text{trop}} + d_{\text{iono}} + \omega(t_r) + B_r^t \quad (3) $$

where $\phi$ is the modeled phase, $B_r^t$, is an arbitrary phase bias between the transmitter and receiver over time periods where the receiver has not lost lock, and $\omega(t_r)$ is the phase windup (Wu et al., 1993). The $B_r^t$ may further be modeled as integer number of wavelengths and fractional hardware delays in the given receiver and transmitter; see (Blewitt, 1989) for treatment of network ambiguity resolution and (Bertiger et al., 2010) for treatment of the integer ambiguities with a single receiver.

## 4.4. Earth models

For receivers or transmitters that are not orbiting the Earth and are either attached to or associated with the Earth’s crust, we must model the deformation of the Earth’s crust and the orientation of the Earth in inertial space since satellite positions are integrated in inertial space. GipsyX/RTGx implements the IERS standards (Petit and Luzum, 2010). Adjustable parameters include the Earth orientation parameters, polar motion and hour angle.

## 4.5. Filter, smoother, ambiguity resolution

Except for single receiver use cases, most of the computation time is dominated by filtering, smoothing and ambiguity resolution functions. If we treated all the parameters as constant in time, the filter and smoother parts of the code would just be a standard least squares; but many of the parameters do vary in time, for instance the zenith troposphere delay or the receiver and transmitter clock errors. Many force parameters, including empirical accelerations affecting satellite dynamics, are also best treated as stochastic variables in a process called reduced-dynamic orbit

determination (Yunck et al., 1990; Wu et al., 1991). We chose to handle these time variations as a first order Markov process with a forward-in-time filter implemented as a Square Root Information Filter (SRIF) detailed in Bierman (Bierman, 1977). The square root implementation allows for better conditioned matrices, and the use of Householder (orthogonal transformation in multidimensional space) transformations in our implementation of the Bierman algorithms lends itself naturally to more efficient use of processor specific pipelining (Jeong et al., 2012), when available, along with the Message Passing Interface (MPI (Gropp et al., 1999)) for processing across multiple CPU cores and computer clusters (Beowulf (Sterling et al., 2002)).

The forward filter is the best fit of all the past data up to the current epoch. In order to have the parameters fit the future as well as the past data optimally, one needs to smooth the information back in time. We again follow Bierman’s algorithm (Bierman, 1977) for a SRIF smoother, which uses a series of Givens (two-dimensional orthogonal rotation) and Householder transformations. We again optimize efficiency using machine-specific-low-level routines and MPI where appropriate. After editing for outliers in the smoother, integer ambiguities may be constrained for GNSS data. For single receiver ambiguity resolution, we follow the procedure in (Bertiger et al., 2010) extended to multi-GNSS, applying constraints to single differences. For multiple GNSS receivers, we can resolve double difference integer ambiguities following (Blewitt, 1989) extended to multi-GNSS to constrain the double-differenced phase ambiguities, as well as adding single-receiver constraints when appropriate external information is available. The single receiver methods are related but not identical to methods developed in (Laurichesse et al., 2009; Laurichesse, 2011; Loyer et al., 2012). Our ambiguity resolution algorithms, do not currently extend to GLONASS Frequency Division Multiple Access (FDMA) signals.

# 5. Main executable, rtgx

A detailed description of the main loop of the primary executable is contained in electronic supplement 2, rtgx main loop.

# 6. Sample use cases, accuracy, precision

## 6.1. Precision post-processed GNSS orbit and clock determination

### 6.1.1. JPL GPS post-processed operational products

JPL delivers to the geodetic community, and for combination by the International GNSS Service (IGS) (Johnston et al., 2017), two sets of GPS products characterized by different latencies, precision and accuracy levels, but all based on filtering and smoothing of ground tracking data in post-processing. “Rapid” products are generated and distributed daily. “Final” products are generated weekly and typically
have a 14-day latency. We also distribute “Ultra-rapid” products, generated every hour. While ultra-rapid and rapid products are tightly tied to the IGS realization of the ITRF by fixing the coordinates of a subnetwork of stations in the POD process to their ITRF coordinates, Final products come in three flavors. Fiducial-free products, identified by the suffix ‘nf’ (non-fiducial) in their names, are not tied to any specific terrestrial frame; the positions of the stations constituting the network are left floating with a 1 km a priori sigma, so that the frame of these solutions changes from day-to-day. No-net-rotation (‘nnr’ suffix) products are such that 3 no-net-rotation (relative to the ITRF) constraints are enforced in their generation process. Final products with no suffix are tied to the ITRF by means of 7 constraints (3 rotations, 3 translations and 1 scale parameter) and are referred to as NNRTS (no-net-rotation, translation, or scale) or ‘fiducial’ products, though their fiducialization differs from Rapid and Ultra-Rapid. Every set of products contain files with consistent orbital state estimates, transmitter clock estimates, spacecraft attitude information, Earth rotation parameter (ERP) estimates, and widelane phase biases information. In addition, the no-net-rotation and fiducial-free solutions are associated with coordinate transformation files referred to as x-files. These x-files contain the set of 7 Helmert parameters (3 small rotations, 3 translations, and 1 scale parameter) needed to transform from the frame of the day to the ITRF. The format of these different products is detailed in the GipsyX software documentation. Following extensive testing and validation, the JPL GNSS analysis center transitioned seamlessly from products generated using the GIPSY software package to products created using GipsyX on January 29, 2017. More recently, the entire span of 1994 through 2018 were reprocessed using GipsyX, consistent with IGS repro2 standards. All products are available at https://sideshow.jpl.nasa.gov/pub/JPL_GNSS_Products.

Table 3 displays the characteristics of each type of post-processed product. Latencies are indicated as well as commonly used performance metrics. JPL orbit and clock products cover 30 h centered at noon, so that each daily solution overlaps with the next-day solution by 6 h. After removing 30-min tails at both ends of this overlap period, RMS statistics on orbit and clock differences are computed in the 5-h overlap period as an internal measure of the precision of the product. The long-term medians of the RMS clock and orbit differences in the overlap is shown as “Precision” in Table 3. We measure accuracy in Table 3 by comparing JPL’s orbit and clock solutions with the combined final IGS orbit and clock solutions for the same days. Although we have listed this as “Accuracy” in Table 3, we note that the IGS Final product is a weighted average of all the contributing analysis centers including JPL.

As mentioned above, each set of products includes a file containing estimates of the Earth rotation parameters. The ERPs consist of two biases and two rates to describe the motion of the Earth’s pole with respect to the Earth’s crust:

Xp, Yp and their respective rates. The excess length-of-day (LOD) is defined as the time derivative of UT1-UTC, where UT1 is the time scale associated to the total rotational phase angle of the Earth and UTC is the Universal Time Coordinated. All 6 parameters are directly observable by GNSS, with the exception of UT1-UTC. As a result, the polar motion coordinates and rates along with LOD are freely estimated in the GNSS precise orbit determination process whereas UT1-UTC is tightly constrained to a nominal value reported in the IERS combined ERP solution, the Bulletin A file (IERS, 2019). Operational processes at JPL rely on the predicted portion of Bulletin A for a priori values for Earth orientation parameters, while reprocessing campaigns typically use the Bulletin A final values of UT1-UTC. All JPL estimates are routinely contributed to the IGS for their combined ERP products. The IGS, in turn, contribute these combined GNSS ERP products to the IERS for their final multi-technique combined products. Figs. 1–3 show the differences between the JPL Final no-net-rotation solutions (contributed to the combined IGS products) relative to the final combined IGS solution for Earth Rotation Parameters. The noise in the differences in the polar motion values amounts to about 20 µas for both coordinates, which corresponds to 0.6 mm at the surface of the Earth. The scatter in the differences of the LOD estimates is higher at 16 µs/day, translating into 1.2 mm at the surface of the Earth over one day. These numbers are consistent, even slightly better for the pole coordinates, than the statistics reported by (Rebischung et al., 2016): they found that the weighted RMS of the IGS analysis centers (AC) pole coordinate residual time series range between 25 and 40 µas, while the weighted RMS of the AC LOD residual time series range from 8 to 20 µs/day. Over the timespan extending from Jan. 29, 2017 to Apr. 20, 2019 (since the transition from JPL’s legacy software to GipsyX), the scatter in the differences between the JPL Xp, Yp, and LOD solutions relative to 6 other IGS analysis centers using different software packages ranges from 27 to 49 µas, 26 to 44 µas, and 13 to 18 µs/day, respectively. Again, these statistics are in full agreement with the numbers cited in (Rebischung et al., 2016). Orbit modeling differences, e.g., in solar radiation pressure modeling, are known to impact the accuracy of GNSS-based ERP determination. In particular, they could explain the bias visible in Fig. 3. The JPL versus IGS ERP differences may be compared directly to the published uncertainties of the final combined ERP values from the IERS Bulletin A relative to the IGS combined series (Bizouard et al., 2017). Table 4 of (Bizouard et al., 2017) shows Xp and Yp scatter of IGS combined at 31 and 27 µas respectively and 10 µs for LOD for ITRF 2014 from 2010–2015.

### 6.1.2. *Four constellation GNSS orbits and clocks*

A greater number of satellites may provide significant benefits to users, including additional coverage in urban canyons, and perhaps better estimates of parameters such as troposphere, Earth orientation, and geocenter.
<caption>Table 3</caption>

JPL post-processed operational products. Units for all precision and accuracy metrics are centimeters. Statistics were computed from time of switch (Jan. 29, 2017) to Apr. 20, 2019. Each clock metric is the median value of the daily root-mean-square values of the overlaps/differences across the GPS satellites after removing a linear trend from the entire GPS constellation to account for reference clock differences. Each orbit metric is the median value of the daily median of the 3D RSS of RMS positions across the GPS constellation.

<table><thead><tr><th>Product</th><th>Latency</th><th>Orbit Precision (cm)</th><th>Clock Precision (cm)</th><th>Orbit Accuracy (cm)</th><th>Clock Accuracy (cm)</th></tr></thead><tbody><tr><td>Ultra</td><td>&lt;2 h</td><td>3.4</td><td>3.0</td><td>2.2</td><td>4.3</td></tr><tr><td>Rapid</td><td>next day by 16:00 UTC</td><td>2.6</td><td>2.4</td><td>2.0</td><td>4.0</td></tr><tr><td>Final</td><td>&lt;14 days</td><td>2.3</td><td>2.3</td><td>1.9</td><td>3.8</td></tr></tbody></table>

<figure><img src="image_2.png" /><figcaption>Fig. 1. Differences between the Earth's pole coordinates (Xp, Yp) determined by JPL's Final no-net-rotation solutions and values from the IGS final combined solution. The start date of the period covered corresponds to the date of the operational transition from GIPSY to GipsyX at the JPL IGS Analysis Center. The scatter for both coordinates, once converted to equivalent linear distance, is of the order of 0.6 mm at the surface of the Earth.</figcaption></figure>

<figure><img src="image_3.png" /><figcaption>Fig. 3. Differences between Length-of-Day recovered by JPL's Final no-net-rotation solutions and LOD final combined solution distributed by the IGS. The start date of the period covered corresponds to the date of the operational transition from GIPSY to GipsyX at the JPL IGS Analysis Center. The scatter around the mean is of the order of 16 µs per day, equivalent to 1.2 mm at the surface of the Earth over a day. The visible –35 µs/d bias, equivalent to 2.6 mm at the surface of the Earth over a day, is attributed to orbit modeling errors and in particular solar radiation pressure modeling errors in GNSS data processing, as discussed in (Sibthorpe et al., 2011; Ray, 1996.) for instance.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>Fig. 2. Differences between the Earth's pole coordinate (Xp, Yp) rates determined by JPL's Final no-net-rotation solutions and the IGS final combined ERP product. The start date of the period covered corresponds to the date of the operational transition from GIPSY to GipsyX at the JPL IGS Analysis Center. The scatter for both rates correspond to an equivalent linear distance of about 2.8 mm/day at the surface of the Earth.</figcaption></figure>

<caption>Table 4</caption>

Daily overlaps of JPL four-constellation products (GPS, BeiDou, Galileo and GLONASS). GPS (JPL Rapid), provided for comparison, is computed using JPL's standard GPS-only Rapid products described in Section 6.1.1. Units for all precision metrics are centimeters. All statistics were computed from Aug. 1 to Aug. 31, 2019. Each clock metric is the median value of the daily root-mean-square values of the overlaps/differences across the GNSS satellites after removing a linear trend from the entire GPS constellation to account for reference clock differences, a subsequent linear trend from each of Galileo and BeiDou to account for constellation bias reference differences, and a linear trend from each GLONASS satellite due to the presence of range biases – *likely making the GLONASS clock overlaps spuriously small. Each orbit metric is the median value of the daily median of the 3D RSS of RMS positions across the GNSS constellations.*

<table><thead><tr><th>Product</th><th>Orbit Precision (cm)</th><th>Clock Precision (cm)</th></tr></thead><tbody><tr><td>GPS (JPL Rapid)</td><td>1.9</td><td>1.9</td></tr><tr><td>GPS (4-constellation)</td><td>2.0</td><td>1.5</td></tr><tr><td>BeiDou GEO</td><td>67.9</td><td>15.3</td></tr><tr><td>BeiDou IGSO</td><td>10.2</td><td>7.8</td></tr><tr><td>BeiDou MEO</td><td>5.3</td><td>5.4</td></tr><tr><td>Galileo</td><td>2.8</td><td>1.4</td></tr><tr><td>GLONASS</td><td>5.5</td><td>1.1*</td></tr></tbody></table>
To this end, many IGS Multi-GNSS Experiment (MGEX, (Montenbruck et al., 2017), http://mgex.igs.org/) analysis centers are now routinely processing multiple constellations, particularly GPS, BeiDou, Galileo and GLONASS, although as yet no official combined product is available. It is worth considering however, that processing new satellites might also have unintended impacts on, e.g. GPS quality compared to existing GPS-only products, because such new satellites may require individual treatment to allow for accurate modeling of their specific observation types (e.g. GLONASS frequency dependent code biases) and the spacecraft themselves (e.g. solar radiation pressure, attitude and antenna offsets (Montenbruck et al., 2015)). Although JPL has been processing these constellations for some years now in real-time (see Section 6.2), we are still refining our station-data selection procedures for fully automated post-processed products, as well as our force model choices, alongside internal quality testing, before making an operational version corresponding to one or more of our GPS-only product lines (see Section 6.1.1). Such internal testing showed a dramatic improvement in quality as of GPS week 2056 (starting June 2, 2019) with the release of a new IGS ANTEX file (https://kb.igs.org/hc/en-us/articles/216104678-ANTEX-format-description) which, for the first time, includes antenna corrections for all satellites from the four major constellations. With our most recent fully automated software, we have produced daily, Rapid-like 4-constellation products for the month of August, 2019, using 120 stations for GPS and 80 for the remaining constellations. Accuracy is more difficult to measure in this case owing to the lack of a provably more accurate multi-GNSS product with which to compare, whether from IGS MGEX or elsewhere, and therefore we leave off such an assessment until a future date. Individual MGEX analysis center contributions may continue to exhibit significant differences until consistent antenna information, only relatively recently available, begins to see widespread use. Our four-constellation products contain an average of 102 satellites per day, approximately comprised of 31 GPS, 28 Beidou, 22 Galileo and 21 GLONASS. While this dataset is limited, it gives a preliminary indication of the sort of performance that might be expected from multi-GNSS products processed with GipsyX/RTGx, which should also improve with time as models and procedures are further enhanced.

Table 4 displays the overlap precision statistics of each constellation/sub-constellation along with statistics for the JPL GPS Rapid products described in Section 6.1.1 over the same time period for comparison. As with Table 3, each daily solution overlaps with the next-day solution by 6 h. After removing 30-min tails at both ends of this overlap period to avoid edge effects, statistics on orbit and clock differences are computed in the 5-h overlap period as an internal measure of the precision of the product. The long-term medians of the clock and orbit differences in the overlap is shown as “Precision” in Table 4.

## 6.2. Real-time GNSS orbit and clock solutions, GNSS differential corrections

The real-time configuration of GipsyX/RTGx is typically called just RTGx. RTGx has been driving real-time orbit determination processes in the GDGPS System since 2014 (RTG was the engine prior to that, with some overlap). These include GPS-only orbit determination filters as well as GPS plus various other GNSS constellations. A typical GNSS orbit determination filter may include 1 Hz phase and pseudorange measurements from 150 ground sites tracking 32 GPS satellites, 24 GLONASS satellites, 35 BDS satellites, and 24 Galileo satellites (as of June 2019). Such a filter produces orbital states every 60 s and clock solutions at 1 Hz, with a latency(difference in the time corrections are available to a user – epoch of the data used for the clock correction) that never exceeds 6 s relative to the measurement epochs. GNSS accuracy is typically quoted in terms of User Range Error (URE), Eq. (4), where h is the radial orbit error, c is the cross-track error, l is the along-track error, and clk is the clock error:

$$URE = \sqrt{(h - clk)^2 + (l^2 + c^2)/50} \quad (4)$$

For the full derivation, see the appendix in (Zumberge and Bertiger, 1996). Here we have used the standard approximations of the coefficients. The real-time URE relative to post-processed solutions is typically 5 cm RMS (Fig. 4).

The specialized requirements of real-time operations have driven the design of some unique and powerful RTGx features. One such feature is the ability to reconfigure a real-time running filter into partitions such that some parameters associated with specific satellites or specific stations, are estimated but have no influence on any other parameters. At the GDGPS System these “decoupled” partitions are used to accommodate unhealthy satellites in order to keep monitoring them while protecting the rest of the estimated constellations from the potential mismodeling of unhealthy satellites, for example due to maneuvers. Our decoupling algorithm is equivalent to infinitely weighting the data from the associated satellite or station following a white-noise re-set of the parameters associated with the satellite or station. For earthquake monitoring the GDGPS System operates orbit determination filters with multiple decoupled partitions where the position of hundreds of ground sites (a site per partition) are estimated kinematically at 1 Hz (Larson et al., 2007), but these sites are prevented from impacting the satellites and stations in other partitions. Only the main partition is allowed to influence other partitions. The incorporation of kinematic ground site positioning within an orbit determination filter, albeit in a decoupled partition, is more accurate than forward filtering point positioning due to the availability within the real-time filter of the full temporal correlations
Forward Filter GPS-only 5-minute POD, May 1 - 13, 2018

(reference: JPL Rapid)

<figure><img src="image_5.png" /></figure>

Fig. 4. RMS URE (black curves) and worst URE (red curves) across all GPS satellites as a function of tracking network size for real-time orbit determination by RTGx within the GDGPS System during May 1–13, 2018. Significant improvement due to integer phase ambiguity resolution (dashed curves) is evident relative to the solutions with float ambiguity resolution (solid curves). The reference solutions are provided by JPL's daily post-processed solutions (see Table 3 above for post-processing accuracy). The RMS URE approaches 5 cm with ambiguity resolution and 90 or more real-time tracking sites.

among all filter parameters. Few centimeter positioning accuracy (from high-quality, continuous tracking sites) is accomplished at 1 Hz (Fig. 5). Real-time time series of hundreds of ground sites are produced by the GDGPS System and made available for natural hazard monitoring on a global scale (via https://ga.gdgps.net/).

A cold start of a real-time GNSS orbit determination may take from hours to days to fully converge (i.e. to achieve the accuracy level of a long running continuous filter), depending on the quality of the initial states. In certain operational scenarios it is desirable to restart the filter with some modified configuration, while achieving instant convergence. This capability is enabled in RTGx through the commanded (or scheduled) saving of a “snapshot” of the entire filter state and model state, and the ability to read this snapshot upon restarting to achieve a hot-start capability. By instant convergence we mean that if the identical data are fed to the filter after the snapshot epoch, the solution will be identical to the continuous running filter. This feature allows the operator to rapidly correct for errors in the past and catch back up to real-time.

Other specialized features designed to support robust real-time operations include special handling for clock references to ensure smooth transition from one reference clock to another, on-the-fly decoupling of satellites or sites, pairwise filter operations combining a low cadence orbit filter (e.g., every 60 s) with a high cadence (e.g. every second) clock filter for optimal throughput.

The GDGPS System produces a variety of real-time differential corrections for the GNSS broadcast ephemerides,

Horizontal (1D) PPP Error, September 2017

(average / standard deviation of 24-hour RMS error)

<figure><img src="image_6.png" /></figure>

Vertical PPP Error, September 2017

(average / standard deviation of 24-hour RMS error)

<figure><img src="image_7.png" /></figure>

Fig. 5. Horizontal (top) and vertical (bottom) real-time positioning accuracy of real-time kinematic 5-min point-positioning of 125 GDGPS tracking sites with RTGx as decoupled partition of the GDGPS GPS orbit determination filter during September 2017. Each point represents a monthly average of daily RMS positioning error. The error bars depict the standard deviation of the daily RMS values over the entire month. Integer phase ambiguity resolution was used to generate the black curves, and float phase ambiguity resolution was used to generate the red curves. Stations are ordered according to their integer ambiguity-resolved positioning accuracy. The reference solutions are provided by long-term static point-positioning.

based on the real-time orbit and clock solutions, including formats that enable real-time ambiguity resolutions such as the state-space representation standards 10403.3 of the Radio Technical Commission For Maritime Services (RTCM) http://www.rtcm.org/differential-global-naviga-tion-satellite-dgnss-standards.html.

## 6.3. Precise point positioning with GNSS

### 6.3.1. Realization of the terrestrial reference frame with precise point positioning

In (Bertiger et al., 2010), it was demonstrated that the IGS implementation of the current ITRF frame at the time,
ITRF2005, could be realized with 24-h of GPS data to an accuracy of 2 mm in the horizontal and 6 mm in the vertical processing 6-months of data from a 106 reference frame sites. Here, we do a similar experiment with 59 stations in ITRF2014. The stations selected are available on at least 90% of the days from 2002–2012, are part of IGS14 (IGS implementation of ITRF2014) (Rebischung et al., 2016), and are selected to be well-distributed geographically. Fig. 6 shows the geographic distribution of sites used for this test.

First, each site was point positioned using a tree very similar to the default tree included with GipsyX-1.0. The only substantive difference with the default tree is the use of Ionosphere Exchange (IONEX) files (Schaer et al., 1998) for second-order ionospheric corrections (except for one no-IONEX case). A GPT2-based (Böhm et al., 2015) nominal troposphere and mapping function was used. Each site was positioned using JPL's repro3.0 IGS14 Final products, with separate PPP runs for the non-fiducial (NF), no-net rotation (NNR), and no-net translation, rotation, or scale constraint (NNRTS also referred to as “fiducial”) products from 2008–06-01 through 2008–11-30. In the case of the NNR and NF products, a 7-parameter Helmert transformation is performed using the appropriate product x-file to convert the resulting position into the reference frame. Positions are compared with those given in ITRF2014. For each station, a 5-sigma edit is performed on total position differences, then the root-mean-square (RMS) is calculated in east, north and vertical (E, N, V) components. We then report the median of each of these components in Table 5. Fig. 7 shows a histogram of the frame repeatability for all 59 stations including all outliers. All three methods (NF, NNR, NNRTS) produce nearly identical results. When compared to the GIPSY-OASIS results in (Bertiger et al., 2010), E and N are almost identical but the results here are a bit worse in V (6.5 vs. 6.0 mm). We attribute the difference in vertical to station selection differences (i.e. 59 IGS14 stations vs. 106 IGS08), as vertical repeatability is highly site-dependent.

While using IONEX files does not improve repeatability in this analysis, we still recommend its use to be consistent with the models used in the Final POD process, and because other analyses show that omitting it can introduce an offset in the z-component (Ries et al., 2018). In this test, the offset is an average of 0.4 mm. The ENV frame repeatability of about 2, 2, and 6 mm in Table 5 is significantly better than the 24-h RMS repeatability of 3, 6, and 9 mm shown in (Soycan, 2011) using Bernese 5.0 software. We note that significant improvements are seen in the east component with bias fixing (Bertiger et al., 2010) with the GipsyX/RTGx bias fixing methods as well as the result with PANDA software (Ge et al., 2007) with bias fixing based on fractional satellite delays. In (Ge et al., 2007), they compute repeatability relative to IGS weekly solutions after removing a 7-parameter Helmert transformation and get mean RMS E, N, and V repeatability of 2.4, 2.7, and 5.3 mm with bias fixing. These are not exactly the same statistics as shown in Table 5. The weekly solutions will take out some vertical signal from seasonal loading not contained in the ITRF2014 frame.

Fig. 8 shows the seasonality of PPP residuals relative to the ITRF solution for ALGO (Algonquin Park, Canada) using IONEX corrections and NF products for 16 years. Most of the vertical residual can be explained by atmospheric and hydrological loading (Tregoning and Watson, 2009), models of which are also plotted. The horizontal displacements also show seasonality due to other effects such as snow cover (Larson and Nievinski, 2013). This example illustrates one complexity of repeatability calculations, and how results can be biased by short time periods and/or inhomogenous distributions of stations.

### 6.3.2. Kinematic positioning

For many applications the receiver is moving, and it is the accuracy of the kinematic position of the GNSS receiver that matters. Applications, among many, include synthetic aperture radar (SAR) from aircraft, ocean floating buoys used to determine sea surface height, co-seismic

<figure><img src="image_8.png" /><figcaption>Fig. 6. Geographic distribution of the 59 sites used for PPP testing.</figcaption></figure>
<table><caption>Table 5
Median daily RMS repeatability relative to ITRF2014 solution.</caption><thead><tr><td>PPP type</td><td>East (mm)</td><td>North (mm)</td><td>Vertical (mm)</td></tr></thead><tbody><tr><td>NF no-IONEX</td><td>1.90</td><td>1.98</td><td>6.47</td></tr><tr><td>NF</td><td>1.89</td><td>2.07</td><td>6.49</td></tr><tr><td>NNR</td><td>1.92</td><td>2.11</td><td>6.53</td></tr><tr><td>NNRTS</td><td>1.96</td><td>2.11</td><td>6.45</td></tr></tbody></table>

<figure><img src="image_9.png" /><figcaption>Fig. 7. Histogram: East, North, and Vertical ITRF repeatability for 59 sites used for static PPP testing using NF products and IONEX corrections. Each count represents a single day of processing for a single station.</figcaption></figure>

deformation, and ice sheet movement, see (Hensley et al., 2008; Bonnefond et al., 2003; Born et al., 1994; Simons et al., 2011; Doake et al., 2002) for example. To give an

idea of the expected performance from GipsyX, we kinematically positioned a set of 45 ITRF2014 frame stations every 5 min for the month of June 2008 using GPS data. Table 6 shows the RMS differences with the ITRF2014 frame positions where the position has been adjusted as a loose random walk, 1 m/√s, while the troposphere is a bit more constrained than in the static point positioning. Instead of adjusting the zenith delay and tropospheric gradient parameters, only a more constrained zenith delay is adjusted to help break the correlation with the random walk vertical position.

Some kinematic applications may have better information about the tropospheric delay. For instance, high flying aircraft typically have a small wet delay and the dry zenith delay may be computed using a pressure sensor. To demonstrate the expected performance of kinematic positioning with better known troposphere, we repeated the experiment from Table 6, but fixed the zenith delay and the gradient parameters to the values from static 24-h position solutions. Comparing Table 6 with Table 7, one can see a significant and larger improvement in the vertical.

## 6.4. Precise orbit and clock determination of low Earth orbiters

Both GIPSY and RTG supported POD of satellites in Low Earth Orbit (LEO), starting with the launch of TOPEX/Poseidon (Bertiger et al., 1994). As a superset of these two software, GipsyX/RTGx supports precise post-processing, reduced-dynamic orbit determination, as well

<figure><img src="image_10.png" /><figcaption>Fig. 8. Daily of PPP ALGO over 13 years vs ITRF frame solution compared with displacements from http://loading.u-strasbg.fr/displ_all.php (Gegout et al., 2010) using the sum of atmospheric loading (ATMMO) and hydrological loading (GLDAS2, Rodell et al., n.d.) models in E (upper-left), N (upper-right), and V (lower-left) components. A higher temporal-resolution plot of V is provided in the lower-right.</figcaption></figure>
<table><caption>**Table 6**
Kinematic positioning 45 ITRF2014 frame stations, adjusting troposphere zenith delay. For each station RMS difference relative to the station frame position is computed for the month of June 2008. The mean and median of the station RMS values are listed.</caption><thead><tr><th rowspan="2"></th><th colspan="3">RMS(mm)</th></tr><tr><th>East</th><th>North</th><th>Vertical</th></tr></thead><tbody><tr><th>Mean</th><td>6.3</td><td>7.2</td><td>20.3</td></tr><tr><th>Median</th><td>6.0</td><td>7.2</td><td>20.0</td></tr></tbody></table>

<table><caption>**Table 7**
Kinematic positioning 45 ITRF2014 frame stations, tropospheric delay from static positioning, statistics as in 6.</caption><thead><tr><th rowspan="2"></th><th colspan="3">RMS(mm)</th></tr><tr><th>East</th><th>North</th><th>Vertical</th></tr></thead><tbody><tr><th>Mean</th><td>5.6</td><td>6.3</td><td>15.4</td></tr><tr><th>Median</th><td>5.4</td><td>6.3</td><td>14.9</td></tr></tbody></table>

as embedded real-time POD. For post-processed LEO POD using GPS data, the GipsyX/RTGx provides identical capability to GIPSY, but with a simpler user interface along with additional capability for the newer GNSS constellations, DORIS, and SLR. Tests were performed with GRACE (Tapley et al., 2004), the follow-on GRACE mission (GRACE-FO), and Jason-2 (Lambin et al., 2010; Bertiger et al., 2010; Cerri et al., 2010). Since the GPS determined orbits for Jason-2 are sub-centimeter accurate (Bertiger et al., 2010), we can use those orbits to measure our accuracy for Jason-2 orbits determined with other data types. Mixing GNSS, SLR, and DORIS data is also possible in GipsyX/RTGx, but has not had extensive testing yet. In the subsections below, we detail the expected accuracy and software capabilities for LEO POD using GPS, SLR, and DORIS data.

### 6.4.1. SLR Orbit determination

GipsyX has many capabilities for SLR processing, including producing residuals to orbits determined from other techniques and independent orbit determination. For one experiment, We used 30-h arcs of SLR data from 15 well-distributed, low-bias SLR sites, centered on noon of each day in 2015 to determine Jason 2's orbit with a median radial RMS difference to GPS-determined orbits of 1.9 cm. We also validated our SLR capabilities by performing orbit determination on the LAGEOS 1 and 2 satellites using 7-day arcs for all of 2017. We then compared our orbits to the ILRSA combined orbits. Our LAGEOS 1 and 2 orbits have an overall RMS difference of 5 mm in radial, 16 mm in cross-track and 18 mm in along-track to LAGEOS, which is comparable in magnitude to the inter-center differences measured by the ILRSA combination to other centers over the course of the year (averages of 5 mm, 24 mm, and 25 mm in the weekly ilrsa weekly sum files across all centers for both LAGEOS satellites).

### 6.4.2. GPS — GRACE, GRACE-FO POD and time synchronization

GipsyX/RTGx is the operational software used for GRACE-FO to determine the orbit and clock of the GRACE-FO spacecraft.

GRACE and GRACE-FO are among the best low Earth orbiting satellites to test POD with GPS since they have a measurement of the inter-satellite range (up to a bias) that is good to the micron-level. This is the dual-way range measurement made at K and Ka band, K-Band Range (KBR) (Dunn et al., 2003; Thomas, 1999). It is the KBR's precision along with an accelerometer measuring the non-gravitational forces that allows the precise recovery of the Earth's time varying gravity field. Although GipsyX/RTGx, can process KBR range data and use accelerometer data, neither is used in the operational POD for GRACE or GRACE-FO discussed here. Fig. 9 shows the daily standard deviation of the difference of the range determined by GPS POD and the KBR (clocks aligned with the independent operational code). The processing used GPS data sampled every 10-s, antenna calibrations based on one month of data, August 2010, reduced dynamics selected to perform well over the full time period, and bias fixing. The average daily baseline standard deviation was 1.6 mm. This is significantly better than the 10 mm found in Table 8 of (Kang et al., 2006) using the Center for Space Research's MSODP software which does not include bias fixing. It improves over the results in (Bertiger et al., 2010) mostly due to the use of higher rate data. (Kroes et al., 2005; Jäggi et al., 2007) show improved baseline determination of a little under a mm if biases are fixed between the two spacecraft with a simultaneous POD of both spacecraft. In this processing, both spacecraft are processed independently.

Although the KBR measurement is not strongly dependent on the GPS POD/Clock solution for GRACE at the 100 micron level, it is necessary to have the error of the rel-

<figure><img src="image_11.png" /><figcaption>Fig. 9. GRACE KBR-GPS Range, Daily GPS POD Baseline Accuracy, Four 1-Month Samples: Aug. 2008, Aug. 2010, June 2016, April 2017.</figcaption></figure>
<caption>Table 8</caption>

Comparison of GipsyX frame with ITRF2014 on January 1, 2005. Three translations, one scale, and three rotations are given in mm while the corresponding rates are given in mm/yr. Rotations and scale are multiplied by the Earth radius to get mm at the surface.

<table><thead><tr><th>Parameter</th><th>TX</th><th>TY</th><th>TZ</th><th>S</th><th>RX</th><th>RY</th><th>RZ</th><th>Unit</th></tr></thead><tbody><tr><td>Offset</td><td>0.508</td><td>-0.192</td><td>-0.144</td><td>2.903</td><td>0.543</td><td>-0.053</td><td>-0.354</td><td>mm</td></tr><tr><td>Offset s.d.</td><td>0.102</td><td>0.095</td><td>0.200</td><td>0.521</td><td>0.117</td><td>0.141</td><td>0.204</td><td>mm</td></tr><tr><td>Rate</td><td>-0.184</td><td>0.044</td><td>0.498</td><td>0.153</td><td>0.001</td><td>0.076</td><td>0.156</td><td>mm/yr</td></tr><tr><td>Rate s.d.</td><td>0.013</td><td>0.012</td><td>0.027</td><td>0.036</td><td>0.018</td><td>0.021</td><td>0.026</td><td>mm/yr</td></tr></tbody></table>

ative clocks between the two GRACE spacecraft less than 160 ps, 1-sigma, after removing a bias (see (Thomas, 1999) Eq. 3.19) to have the errors in the KBR range due to the clock alignment less than 0.5 microns. Fig. 10 shows a measure of the precision of the relative clock between the two GRACE spacecraft. Since we are solving for the GRACE clock as a white noise process, relative to a fixed GPS reference with 30-h of data centered on noon of each day, there is a six hour overlap from one day to the next. If our reference clock were the same, and our solution were perfect, the difference of the clock solutions between these two days would be zero. In Eq. (5), the first term is the difference of the two GRACEA clocks, the second GRACEB. Taking the difference of the A and B differences removes any reference clock, thus the RMS of this quantity over the central 5 h of the 6-h is a good measure of the relative clock precision. Fig. 10 shows that we are significantly under the 160 ps requirement with a median value of 20.3 ps.

$$dd_{clk} = (A_{clk}^{day1} - A_{clk}^{day2}) - (B_{clk}^{day1} - B_{clk}^{day2}) \quad (5)$$

### 6.4.3. Doppler System POD

In addition to GNSS phase and range data, GipsyX/RTGx can process almost any radio-metric phase and range data and is easily extended to data types that are linear combinations of radio-metric range/phase data; for instance, Satellite Laser Ranging or Doppler Orbitography

<figure><img src="image_12.png" /><figcaption>Fig. 10. RMS Double Difference GRACE Clock Overlaps, Four 1-Month Samples: Aug. 2008, Aug. 2010, June 2016, April 2017.</figcaption></figure>

and Radiopositioning Integrated by Satellite (DORIS) (Willis et al., 2010). Current DORIS data is given in a RINEX like format (Auriol and Tourain, 2010) with dual-frequency phase measurements and range measurements similar to GNSS. Traditionally this is processed as Doppler

$$D = \frac{\phi(\tau_1) - \phi(\tau_2)}{\tau_1 - \tau_2} \quad (6)$$

where $\phi(\tau)$ is the measured phase on an orbiter from a ground transmitter at its local time $\tau$. Prior to June 20, 2008 the Doppler measurements were given directly with corrections for errors in the satellite clock and phase center offsets. The count time, $\tau_1 - \tau_2$, is 10 s of local satellite time in Eq. (6). When forming DORIS Doppler from the phase on the RINEX DORIS data, one must model the satellite receiver clock errors for the computed Doppler model, C,

$$ \begin{aligned} C &= \frac{\frac{\psi(\tau_1) - \psi(\tau_2) + c(E(\tau_1) - E(\tau_2))}{t_1 - t_2 + E(\tau_1) - E(\tau_2)}} \\ &= \frac{\frac{\psi(\tau_1) - \psi(\tau_2) + c(E(\tau_1) - E(\tau_2))}{t_1 - t_2} + \frac{E(\tau_1) - E(\tau_2)}{1 + \frac{E(\tau_1) - E(\tau_2)}{t_1 - t_2}}}{1 + \frac{E(\tau_1) - E(\tau_2)}{t_1 - t_2}} \end{aligned} \quad (7) $$

The phase measurement, with units of length, is written as $\phi(\tau) = \psi(\tau) + cE(\tau)$, where $E(\tau)$ is the difference between coordinate time ($t$) and the local time of the receiver in seconds and $c$ is the speed of light. For the DORIS system, the receiver is driven by a stable clock (Allan deviation on the order of $10^{-13}$ at a day) and the rate of the clock errors including relativistic effects, $\frac{E(\tau_1) - E(\tau_2)}{t_1 - t_2}$ is small, typically on the order of $10^{-7}$ s/s, we can expand Eq. (6)

$$ \begin{aligned} C &\approx \frac{\psi(\tau_1) - \psi(\tau_2)}{t_1 - t_2} + c \frac{E(\tau_1) - E(\tau_2)}{t_1 - t_2} \\ &\quad - \frac{E(\tau_1) - E(\tau_2)}{t_1 - t_2} \frac{\psi(\tau_1) - \psi(\tau_2)}{t_1 - t_2} \\ &\quad - c \left( \frac{E(\tau_1) - E(\tau_2)}{t_1 - t_2} \right)^2 \end{aligned} \quad (8) $$

The computed measurement model in Eq. (8) is very close to the difference of phase measurements and fairly easy to code in terms of the basic phase measurement used for GNSS. This model also may be used for the legacy DORIS doppler data since the clock error terms are essentially zero, due to pre-calibration. On the DORIS RINEX file, the nominal receiver clock values are given as a time series derived from the range data. A quadratic fit over time periods on the order of a day maybe used for $E(\tau)$.

We have tested this with data from the Jason-2 spacecraft (Couhert et al., 2015) using station coordinates from
the latest DPOD (DORIS Precise Orbit Determination) solution, DPOD2014, aligned on ITRF2014 (Moreaux et al., 2019). To take into account the effect of the South Atlantic Anomaly (SAA) on the onboard oscillator (Willis et al., 2016), we did not use any correction model but removed a few stations in the South America region. In Fig. 11, we compare the Jason-2 radial orbit position determined with the DORIS RINEX data from Feb. 14, 2014 through August 23, 2014 with a GPS-determined orbit whose radial accuracy is about 5 mm (Desai et al., 2018) using ITRF2014 coordinates for tracking stations and typical RMS differences with the GipsyX/RTGx DORIS determined orbit are at the 9.9 mm level. Errors in radial position go directly into the mission's prime objective to measure sea surface height and are the most important metric for ocean altimetry missions.

To better evaluate the SAA effect on Jason2, we will soon investigate which stations are the most affected and use the new available correction model from (Belli and Exertier, 2018), based on T2L2 measurements. In the future, we plan to analyze the Sentinel-3A DORIS RINEX data. The fact that satellite clock for Sentinel-3A is the same for the GPS and the DORIS on-board instruments should allow us to compare the GPS clock solution and the DORIS clock. The GPS solution for the on-board clock could be used to properly correct the SAA effects for this satellite, (Jalabert and Mercier, 2018).

## 6.5. Multi-technique geodesy

As an example of multi-technique use, we processed GPS and SLR data simultaneously using GPS and SLR data from low Earth orbiters to tie the systems together. For the development of reference frames, the ties are traditionally made through local ground surveys to measure the vector between the GPS and SLR instruments. With only three days of data, and no ground survey ties, we show agreement at the cm level.

<figure><img src="image_13.png" /><figcaption>Fig. 11. RMS Radial DORIS Doppler radial orbit errors relative to GPS JPL RLSE18.</figcaption></figure>

Depicted in Fig. 12 are summaries (by satellite) of tracking data residuals from a multi-day GipsyX/RTGx solution in which SLR and GPS are combined at the observation level. Represented in the solution are 31 GPS satellites, 5 dedicated SLR targets in space (LAGEOS 1 and 2, Starlette, Stella and Ajisai), and 3 low-Earth orbiters (LEO) with both GPS receivers and corner reflectors (Jason-2 and GRACE A and B). The latter three LEO missions provide a means of accurately linking the SLR and GPS systems without the benefit of ground survey ties. The uniform and precise (cm-level) tracking residuals across both data types (SLR and GPS) testify to the coherence of these diverse space-geodetic observations (radio-metric and optical) in this grand GipsyX/RTGx network solution.

The solution also features 45 and 18 GPS and SLR ground stations respectively. To compensate for know systematic SLR errors, range biases (per satellite) were estimated for a minority (40%) of the stations. Satellite-specific estimates of the optical range variations as a function of the line-of-sight could further improve the results, but were not applied in this case. The positions of all satellites and stations were estimated simultaneously in the GipsyX/RTGx run (Ref.), along with Earth orientation. Using this 'fiducial-free' strategy (Ref.), the network of satellites and stations will be subject to rotation, but the scale, origin, and relative positions will be well determined. Two of the SLR stations in the run (Matera, Italy and Tahiti, French Polynesia) also have participating GPS stations linked by ground survey ties. With only 3 days of tracking data, the GipsyX/RTGx solution reproduced these independent ground ties to 1 cm (3D). The tie between SLR and GPS in the GipsyX/RTGx solution comes solely from the space ties on Jason-1 and the GRACE tandem, testifying not only to the benefits of combining data at the observation level but also to the fidelity of the GipsyX/RTGx modeling.

## 6.6. Building terrestrial reference frames and fitting ground site time series

GipsyX includes station coordinate processing software which can combine solutions from several techniques to build reference frames, simulate future performance, compute transformation parameters between reference frames, and fit individual time series to estimate positions, velocities, seasonals, and possible discontinuities in the linear and periodic fit, which we refer to as breaks. Each tool has command line help available. On-line training provides detailed instructions and example data sets for the most common use cases.

A terrestrial reference frame is typically defined by a set of reference positions with a model of time evolution, the simplest example being a table of positions and velocities. Reference frames are built by first combining daily or weekly solutions for each individual technique and then combining all techniques with ties into a single frame. A
<figure><img src="image_14.png" /><figcaption>Fig. 12. GPS dual-freq. phase residuals and SLR residuals combining SLR and GPS Data; Jason-2, GRACE A, B = Sat. Number 32, 33, 34; LAGEOS 1, 2, Starlette, Stella, Ajisai = 35, 36, 37, 38, 39.</figcaption></figure>

terrestrial reference frame built using GipsyX with inputs from all four geodetic techniques and ties is compared to ITRF2014 (Altamimi et al., 2016) in Table 8 and DTRF2014 (Seitz et al., 2016) in Table 9. See electronic supplement 3 for details.

GipsyX includes a network simulation tool which can generate output files at a given frequency for a specified time span, for example daily files spanning one year. The user provides a model consisting of positions, velocities, breaks, and/or seasonal terms. Various forms of white noise can then be injected. Coordinate noise can be added which is independent from site to site. Geocenter, scale, and/or rotational noise can be injected which is correlated and impacts all sites. Simulated solutions can be processed just like real ones to test software or predict performance of future networks.

The transformation between two reference frames at a particular time can be described by three translations, three rotations, and one scale. Daily transformation parameters known as x-files are computed using an input reference frame to predict coordinates at a particular epoch and then using those predictions to estimate transformation parameters to GNSS positions observed at that epoch. Daily x-files are publicly available as described in Section 6.1.1. A plot of x-file parameters between IGS14 and GipsyX free-network solutions is shown in Fig. 13. The TX parameter has a mean value of −0.5 mm and a standard deviation of 4.7 mm. The TY parameter has a mean value of −0.7 mm and a standard deviation of 5.6 mm. The TZ parameter has a mean value of −5.3 mm and a standard deviation of 6.8 mm. The Scale parameter has mean value of 0.0 mm and a standard deviation of 1.7 mm. These daily corrections indicate how far GNSS only solutions are from a standard reference frame such as IGS14. The GNSS only geocenter values have approached IGS14 more closely over

<figure><img src="image_15.png" /><figcaption>Fig. 13. X-file geocenter and scale parameters relative to IGS14.</figcaption></figure>

time as modeling and data quality have improved and recent research has shown that adding LEO satellites such as GRACE and Jason allow GNSS to be competitive with SLR for both linear and annual geocenter estimation (Haines et al., 2015).

<caption>Table 9</caption>

Comparison of GipsyX frame with DTRF2014 on January 1, 2005. Three translations, one scale, and three rotations are given in mm while the corresponding rates are given in mm/yr. Rotations and scale are multiplied by the Earth radius to get mm at the surface.

<table><thead><tr><th>Parameter</th><th>TX</th><th>TY</th><th>TZ</th><th>S</th><th>RX</th><th>RY</th><th>RZ</th><th>Unit</th></tr></thead><tbody><tr><td>Offset</td><td>0.434</td><td>-0.032</td><td>-0.861</td><td>0.049</td><td>0.693</td><td>1.221</td><td>-1.023</td><td>mm</td></tr><tr><td>Offset s.d.</td><td>0.078</td><td>0.077</td><td>0.174</td><td>0.079</td><td>0.048</td><td>0.067</td><td>0.079</td><td>mm</td></tr><tr><td>Rate</td><td>-0.112</td><td>0.077</td><td>0.629</td><td>0.057</td><td>-0.098</td><td>-0.121</td><td>0.329</td><td>mm/yr</td></tr><tr><td>Rate s.d.</td><td>0.009</td><td>0.010</td><td>0.021</td><td>0.012</td><td>0.010</td><td>0.012</td><td>0.014</td><td>mm/yr</td></tr></tbody></table>
<figure><img src="image_16.png" /><figcaption>Fig. 14. Time series and model fit for ALGO, Algonquin Park, Canada.</figcaption></figure>

Time series for an individual site in a particular reference frame can be fit to estimate positions, velocities, breaks, and seasonals. GipsyX includes tools for automatic break detection and outlier removal. Break detection carries out an exhaustive search based on all possible break locations. An F-test is carried out for the location which minimizes chi squared per degree of freedom. If the F-test is passed that break is kept and the search continues at all remaining possible break locations. The search stops when no remaining break location satisfies the F-test. Outlier removal is based on a maximum value for either the sigma of each point or the residual/sigma ratio for each point. Estimation is carried out using the same time variable square root information filter used for reference frame combinations.

Fitting tools have been optimized to compute large velocity fields quickly. Processing can be done in parallel because point positions at different sites are independent. Given 2000 global sites with daily measurements spanning 20 years; break detection, outlier removal, and parameter estimation can be carried out in about 2 h on a cluster containing thirty-two Intel(R) Xeon(R) X5450 @ 3.00 GHz processors with eight core nodes each. GipsyX is used to

update parameter estimates and edited time series for thousands of globally distributed GNSS sites every week and they are made available via our web service https://side-show.jpl.nasa.gov/post/series.html. Results for the site ALGO are shown in Fig. 14. The plot shows daily measurements with error bars in black, the model fit in red, and break times in green. Velocity estimates with error bars are provided along with residual repeatability/WRMS values for each component.

## 7. Summary

The new GipsyX/RTGx software is a robust and powerful tool for geodetic data analysis and simulations, incorporating decades of expertise and lessons learned from the design and operations of previous software generations, and their applications to the most challenging positioning, navigation, timing, and science applications. GipsyX/RTGx now underlie all GNSS orbit determination operations at JPL, and with hundreds of academic and research licenses, powers geodetic analyses and science operations across the globe.
## Acknowledgements

We are grateful to the following individuals for contributions to GipsyX/RTGx: Jason Gross, Miquel Garcia Fernandez, Jan Weiss.

The research was carried out at the Jet Propulsion Laboratory, California Institute of Technology, under a contract with the National Aeronautics and Space Administration.

We acknowledge support from the U.S Air Force GPS OCX program, from the JPL Global Differential GPS (GDGPS) System, and from NASA's Space Geodesy Program.

Part of this work was performed supported by the Centre National d'Etudes Spatiales (CNES) and is based on DORIS satellite observations.

©2018 California Institute of Technology. Government sponsorship acknowledged.

## Appendix A. Supplementary material

Supplementary data associated with this article can be found, in the online version, at https://doi.org/10.1016/j.asr.2020.04.015.

## References

Altamimi, Z., Rebischung, P., Métivier, L., Collilieux, X., 2016. ITRF2014: a new release of the international terrestrial reference frame modeling nonlinear station motions. J. Geophys. Res.: Solid Earth 121 (8), 6109–6131. https://doi.org/10.1002/2016jb013098.

Auriol, A., Tourain, C., 2010. DORIS system: the new age. Adv. Space Res. 46 (12), 1484–1496. https://doi.org/10.1016/j.asr.2010.05.015.

Bar-Sever, Y., Kroger, P., Borjesson, J., 1998. Estimating horizontal gradients of tropospheric path delay with a single GPS receiver. J. Geophys. Res. 103 (B3), 5019–5035. https://doi.org/10.1029/97jb03534.

Bar-Sever, Y., Kuang, D., 2005. New empirically derived solar radiation pressure model for global positioning system satellites during eclipse seasons. IPN Progress Report, Jet Propulsion Laboratory, Pasadena CA, 42–160.

Bar-Sever, Y.E., 1996. A new model for GPS yaw attitude. J. Geodesy 70 (11), 714–723. https://doi.org/10.1007/s00190050060.

Bar-Sever, Y.E., Bertiger, W.I., Davis, E.A.S., Anselmi, J.A., 1996. Fixing the GPS bad attitude: modeling GPS satellite yaw during eclipse seasons. Navigation 43 (1), 25–40. https://doi.org/10.1002/j.2161-4296.1996.tb01915.x.

Bell, A., Exertier, P., 2018. Long-term behavior of the DORIS oscillator under radiation: The Jason-2 case. IEEE Trans. Ultrason., Ferroelectr., Freq. Control 65 (10), 1965–1976. https://doi.org/10.1109/TUFFC.2018.2855085.

Bertiger, W., Bar-Sever, Y., Christensen, E., Davis, E., Guinn, J., Haines, B., Ibanez-Meier, R., Jee, J., Lichten, S., Melbourne, W., et al., 1994. GPS precise tracking of TOPEX/Poseidon: results and implications. J. Geophys. Res.: Oceans 99 (C12), 24449–24464. https://doi.org/10.1029/94jc01171.

Bertiger, W., Bar-Sever, Y., Haines, B., Iijima, B., Lichten, S., Lindqwister, U., Mannucci, A., Mueller Schoen, R., Munson, T., Moore, A., et al., 1997. A real-time wide area differential GPS system. Navigation 44 (4), 433–447. https://doi.org/10.1002/j.2161-4296.1997.tb02359.x.

Bertiger, W., Bar-Sever, Y., Harvey, N., Miller, K., Romans, L., Weiss, J., Doyle, L., Solorzano, T., Petzinger, J., Stell, A., 2010a. Next

generation gps ground control segment (ocx) navigation design. In: Proceedings of the 23rd International Technical Meeting of the Satellite Division of The Institute of Navigation (ION GNSS 2010), Portland OR, 964–977. https://www.ion.org/publications/abstract.cfm?articleID=9213>.

Bertiger, W., Desai, S.D., Dorsey, A., Haines, B.J., Harvey, N., Kuang, D., Sibthorpe, A., Weiss, J.P., 2010. Sub-centimeter precision orbit determination with GPS for ocean altimetry. Mar. Geodesy 33 (S1), 363–378. https://doi.org/10.1080/01490419.2010.487800.

Bertiger, W., Desai, S.D., Haines, B., Harvey, N., Moore, A.W., Owen, S., Weiss, J.P., 2010. Single receiver phase ambiguity resolution with GPS data. J. Geodesy 84 (5), 327–337. https://doi.org/10.1007/s00190-010-0371-9.

Bierman, G., 1977. Factorization methods for discrete sequential estimation, vol. 128. Academic Press, Cambridge. https://doi.org/10.1016/s0076-5392(08)x6052-7.

Bizouard, C., Lambert, S., Becker, O., Richard, J.Y., 2017. Combined solution C04 for Earth rotation parameters consistent with International Terrestrial Reference Frame 2014, DOI: 10.1007/978-3-642-00860-3_41, https://doi.org/10.1007/978-3-642-00860-3_41, <https://hpiers.obspm.fr/iers/eop/eopc04/C04 guide.pdf>.

Blewitt, G., 1989. Carrier phase ambiguity resolution for the Global Positioning System applied to geodetic baselines up to 2000 km. J. Geophys. Res.: Solid Earth 94 (B8), 10187–10203. https://doi.org/10.1029/jb094ib08p10187.

Blewitt, G., 1989. Carrier phase ambiguity resolution for the global positioning system applied to geodetic baselines up to 2000 km. J. Geophys. Res. 94 (B8), 10187–10203.

Blewitt, G., 1990. An automatic editing algorithm for GPS data. Geophys. Res. Lett. 17 (3), 199–202. https://doi.org/10.1029/gl017i003p00199.

Böhm, J., Möller, G., Schindelegger, M., Pain, G., Weber, R., 2015. Development of an improved empirical model for slant delays in the troposphere (gpt2w). GPS Solutions 19 (3), 433–441. https://doi.org/10.1007/s10291-014-0403-7.

Böhm, J., Niell, A., Tregoning, P., Schuh, H., 2006. Global mapping function (GMF): a new empirical mapping function based on numerical weather model data. Geophys. Res. Lett. 33 (7), L07304. https://doi.org/10.1029/2005gl025546.

Böhm, J., Schuh, H., 2004. Vienna mapping functions in VLBI analyses. Geophys. Res. Lett. 31 (1), L01603. https://doi.org/10.1029/2003gl018984.

Bonnefond, P., Exertier, P., Laurain, O., Ménard, Y., Orsoni, A., Jeansou, E., Haines, B.J., Kubitschek, D.G., Born, G., 2003. Leveling the sea surface using a GPS-catamaran, special issue: Jason-1 calibration/validation. Mar. Geodesy 26 (3–4), 319–334. https://doi.org/10.1080/714044524.

Born, G.H., Parke, M.E., Axelrad, P., Gold, K.L., Johnson, J., Key, K. W., Kubitschek, D.G., Christensen, E.J., 1994. Calibration of the TOPEX altimeter using a GPS buoy. J. Geophys. Res.: Oceans 99 (C12), 24517–24526. https://doi.org/10.1029/94jc00920.

Bruinsma, S., Thuillier, G., Barlier, F., 2003. The DTM-2000 empirical thermosphere model with new data assimilation and constraints at lower boundary: accuracy and properties. J. Atmos. Solar-terrestrial Phys. 65 (9), 1053–1070. https://doi.org/10.1016/s1364-6826(03)00137-8.

Cerri, L., Berthias, J., Bertiger, W., Haines, B., Lemoine, F., Mercier, F., Ries, J., Willis, P., Zelensky, N., Ziebart, M., 2010. Precision orbit determination standards for the Jason series of altimeter missions. Mar. Geodesy 33 (S1), 379–418. https://doi.org/10.1080/01490419.2010.488966.

Chacon, S., Straub, B., 2014. Distributed git. Pro Git. Apress, 97–130. https://doi.org/10.1007/978-1-4842-0076-6_5, https://git-scm.com/book/en/v2.

Couhert, A., Cerri, L., Legeais, J.-F., Ablain, M., Zelensky, N.P., Haines, B.J., Lemoine, F.G., Bertiger, W.I., Desai, S.D., Otten, M., 2015. Towards the 1 mm/y stability of the radial orbit error at regional scales. Adv. Space Res. 55 (1), 2–23. https://doi.org/10.1016/j.asr.2014.06.041.
Dach, R., et al., 2015. Bernese GNSS Software Version 5.2. User manual, Astronomical Institute, University of Bern. Tech. rep. https://doi.org/10.7892/boris.72297, ISBN: 978-3-906813-05-9.

Desai, S., Bertiger, W., Haines, B., Sibois, A., September 2018. GPS-Based Jason-2 and Jason-3 Precision Orbit Determination Solutions in the IGS14 Reference Frame. Presentation at the Ocean Science Topography Team Meeting, OSTST 2018 Meeting (Ponta Delgada, Portugal).

Dilssner, F., Springer, T., Gienger, G., Dow, J., 2011. The GLONASS-m satellite yaw-attitude model. Adv. Space Res. 47 (1), 160–171. https://doi.org/10.1016/j.asr.2010.09.007.

Doake, C., Corr, H., Nicholls, K., Gaffikin, A., Jenkins, A., Bertiger, W., King, M., 2002. Tide-induced lateral movement of brunt ice shelf, antarctica. Geophys. Res. Lett., vol. 29, 8, pp. 67–1, 1226, doi: https://doi.org/10.1029/2001gl014606.

Dunn, C., Bertiger, W., Bar-Sever, Y., Desai, S., Haines, B., Kuang, D., Franklin, G., Harris, I., Kruizinga, G., Meehan, T., 2003. Instrument of GRACE GPS augments gravity measurements. GPS World 14 (2), 16–29.

Ge, M., Gendt, G., Rothacher, M., Shi, C., Liu, J., 2007. Resolution of GPS carrier-phase ambiguities in precise point positioning (PPP) with daily observations. J. Geodesy 82 (7), 389–399. https://doi.org/10.1007/s00190-007-0187-4.

Gegout, P., Boy, J.P., Hinderer, J., Ferhat, G., 2010. Modeling and observation of loading contribution to time-variable gps sites positions. In: Mertikas, S.P. (Ed.), Gravity, Geoid and Earth Observation. Springer Berlin Heidelberg, Berlin, Heidelberg, pp. 651–659.

Gropp, W.D., Lusk, E., Skjellum, A., 1999. Using MPI: portable parallel programming with the message-passing interface, Vol. 1. MIT press, Cambridge. https://doi.org/10.7551/mitpress/7056.001.0001.

Haines, B.J., Bar-Sever, Y.E., Bertiger, W.I., Desai, S.D., Harvey, N., Sibois, A.E., Weiss, J.P., 2015. Realizing a terrestrial reference frame using the Global Positioning System. J. Geophys. Res.: Solid Earth 120 (8), 5911–5939. https://doi.org/10.1002/2015jb012225.

Hairer, E., Norsett, S., Wanner, G., 1993. Solving ordinary differential equations I. nonstiff problems. In: Springer Series in Comput. Mathematics, 2nd Edition. Vol. 8. Springer, Manhattan.

Hensley, S., Wheeler, K., Sadowy, G., Jones, C., Shaffer, S., Zebker, H., Miller, T., Heavey, B., Chuang, E., Chao, R., et al., 2008. The UAVSAR instrument: Description and first results. In: 2008 IEEE Radar Conference. IEEE, pp. 1–6. https://doi.org/10.1109/radar.2008.4720722.

Herring, T.A. and King, R., Floyd, M.A. and McClusky, S., 2018. GAMIT Reference Manual: GPS Analysis at MIT, Release 10.7. Department of Earth. Tech. rep., Massachusetts Institute of Technology, Cambridge, Mass. URL <http://geoweb.mit.edu/gg/Intro_GG.pdf>.

Hugentobler, U., Montenbruck, O., 2017. Satellite orbits and attitude. In: Teunissen, P.J., Montenbruck, O. (Eds.), Springer Handbook of Global Navigation Satellite Systems. Springer International Publishing, Cham, pp. 59–90. https://doi.org/10.1007/978-3-319-42928-1, Ch. 3.

IERS, 2019. IERS Bulletin A. <https://www.iers.org/IERS/EN/Publications/Bulletins/bulletins.html>.

Jalabert, E., Mercier, F., 2018. Analysis of south atlantic anomaly perturbations on sentinel-3a ultra stable oscillator. impact on DORIS phase measurement and doris station positioning. Adv. Space Res. 62 (1), 174–190. https://doi.org/10.1016/j.asr.2018.04.005.

Jeong, H., Kim, S., Lee, W., Myung, S.-H., 2012. Performance of SSE and AVX instruction sets. arXiv preprint arXiv:1211.0820.

Jäggi, A., Hugentobler, U., Bock, H., Beutler, G., 2007. Precise orbit determination for GRACE using undifferenced or doubly differenced GPS data. Adv. Space Res. 39 (10), 1612–1619. https://doi.org/10.1016/j.asr.2007.03.012.

Johnston, G., Riddell, A., Hausler, G., 2017. The international GNSS service. In: Springer Handbook of Global Navigation Satellite Systems. Springer International Publishing, Cham, pp. 967–982, doi: https://doi.org/10.1007/978-3-319-42928-1_33.

Kang, Z., Tapley, B., Bettadpur, S., Ries, J., Nagel, P., Pastor, R., 2006. Precise orbit determination for the GRACE mission using only GPS data. J. Geodesy 80 (6), 322–331. https://doi.org/10.1007/s00190-006-0073-5.

Kedar, S., Hajj, G.A., Wilson, B.D., Heflin, M.B., 2003. The effect of the second order GPS ionospheric correction on receiver positions. Geophys. Res. Lett. 30 (16), 1829. https://doi.org/10.1029/2003gl017639.

Knocke, P., 1989. Earth radiation pressure effects on satellites. Ph.D. thesis, The University of Texas at Austin, Center for Space Research (CSR-89-1).

Knocke, P., Ries, J., Tapley, B., 1988. Earth radiation pressure effects on satellites. In: AIAA/AAS Astrodynamics Conference, Minneapolis, MN., pp. 577–587. https://doi.org/10.2514/6.1988-4292.

Kroes, R., Montenbruck, O., Bertiger, W., Visser, P., 2005. Precise GRACE baseline determination using GPS. GPS Solutions 9 (1), 21–31. https://doi.org/10.1007/s10291-004-0123-5.

Krogh, F.T., 1972. Changing stepsize in the integration of differential equations using modified divided differences. In: Proceedings of the Conference on the Numerical Solution of Ordinary Differential Equations. Springer, Manhattan, pp. 22–71, doi: 10.1007/BFb0066582, https://doi.org/10.1007/BFb0066582.

Lambin, J., Morrow, R., Fu, L.-L., Willis, J.K., Bonekamp, H., Lillibridge, J., Perbos, J., Zaouche, G., Vaze, P., Bannoura, W., et al., 2010. The OSTM/Jason-2 mission. Mar. Geodesy 33 (S1), 4–25. https://doi.org/10.1080/01490419.2010.491030.

Larson, K.M., Bilich, A., Axelrad, P., 2007. Improving the precision of high-rate GPS. J. Geophys. Res.: Solid Earth 112. https://doi.org/10.1029/2006jb004367.

Larson, K.M., Nievinski, F.G., 2013. GPS snow sensing: results from the EarthScope Plate Boundary Observatory. GPS Solutions 17 (1), 41–52. https://doi.org/10.1007/s10291-012-0259-7.

Laurichesse, D., 2011. The CNES Real-time PPP with undifferenced integer ambiguity resolution demonstrator. In: Proceedings of the ION GNSS, pp. 654–662.

Laurichesse, D., Mercier, F., BERTHIAS, J.-P., Broca, P., Cerri, L., 2009. Integer ambiguity resolution on undifferenced GPS phase measurements and its application to PPP and satellite precise orbit determination. Navigation 56 (2), 135–149, doi: 10.1002/j.2161-4296.2009.tb01750.x, https://doi.org/10.1002/j.2161-4296.2009.tb01750.x.

Loyer, S., Perosanz, F., Mercier, F., Capdeville, H., Marty, J.-C., 2012. Zero-difference GPS ambiguity resolution at CNES-CLS IGS analysis center. J. Geodesy 86 (11), 991–1003. https://doi.org/10.1007/s00190-012-0559-2.

Marty, J.C., L., S., Perosanz, F., et al., 2011. GINS: the CNES/GRGS GNSS scientific software. In: 3rd International Colloquium Scientific and Fundamental Aspects of the Galileo Programme, ESA Proceedings WPP326. <http://ids-doris.org/images/documents/report/publications/GINS_GRGSsoftware-Marty-2011.pdf>, Copenhagen, Denmark, 31 August - 2 September.

Melbourne, W., Fisher, S., Neilan, R., Yunck, T., Engen, B., Reigber, C., Tatevjan, S., 1993. The first GPS IERS and geodynamics experiment’1991. In: Permanent satellite tracking networks for geodesy and geodynamics. Springer, Manhattan, pp. 65–80, doi: https://doi.org/10.1007/978-3-642-77726-4_6.

Merkowitz, S.M., Bolotin, S., Elosegui, P., Esper, J., Gipson, J., Hilliard, L., Himwich, E., Hoffman, E.D., Lakins, D.D., Lamb, R.C., Lemoine, F.G., Long, J.L., McGarry, J.F., MacMillan, D.S., Michael, B.P., Noll, C., Pavlis, E.C., Pearlman, M.R., Ruszczyk, C., Shappirio, M. D., Stowers, D.A., 2018. Modernizing and expanding the NASA Space Geodesy Network to meet future geodetic requirements. J. Geodesy 93 (11), 2263–2273. https://doi.org/10.1007/s00190-018-1204-5.

Milani, A., Nobili, A.M., Farinella, P., 1987. Non-gravitational perturbations and satellite geodesy. Non-gravitational perturbations and satellite geodesy. A. Milani, AM Nobili, P. Farinella. Adam Hilger Ltd., Bristol, UK. 125 pp. (1987). ISBN 0-85274-538-9.

Montenbruck, O., Schmid, R., Mercier, F., Steigenberger, P., Noll, C., Fatkulin, R., Kogure, S., Ganeshan, A., 2015. GNSS satellite
geometry and attitude models. Adv. Space Res. 56 (6), 1015–1029.
https://doi.org/10.1016/j.asr.2015.06.019.

Montenbruck, O., Steigenberger, P., Prange, L., Deng, Z., Zhao, Q., Perosanz, F., Romero, I., Noll, C., Stürze, A., Weber, G., et al., 2017. The multi-GNSS experiment (MGEX) of the International GNSS Service (IGS)—achievements, prospects and challenges. Adv. Space Res. 59 (7), 1671–1697. https://doi.org/10.1016/j.asr.2017.01.011.

Moreaux, G., Willis, P., Lemoine, F.G., Zelensky, N.P., Couhert, A., Lakbir, H.A., Ferrage, P., 2019. DPOD2014: a new DORIS extension of ITRF2014 for precise orbit determination. Adv. Space Res. 63 (1), 118–138. https://doi.org/10.1016/j.asr.2018.08.043.

Moyer, T.D., 1971. Mathematical formulation of the double-precision orbit determination program. Tech. rep., technical Report 32–1527. Jet Propulsion Laboratory, Pasadena, CA, USA.

Moyer, T.D., 1981. Transformation from proper time on earth to coordinate time in solar system barycentric space-time frame of reference. Celest. Mech. 23 (1), 57–68. https://doi.org/10.1007/bf01228544.

Odijk, D., Teunissen, P.J., 2013. Characterization of between-receiver galileo inter-system biases and their effect on mixed ambiguity resolution. GPS Solutions 17 (4), 521–533. https://doi.org/10.1007/s10291-012-0298-0.

Parkinson, B.W., Enge, P., Axelrad, P., Spilker Jr, J.J., 1996. Global positioning system: Theory and applications, Volume I. American Institute of Aeronautics and Astronautics, Reston, doi: 10.2514/4.866388, https://doi.org/10.2514/4.866388.

Pavlis, D., Nicholas, J., Rowlands, D., 2017. GEODYN II system description (vols. 1-5) (contractor report). Tech. rep., SGT Inc., Greenbelt, Maryland.

Petit, G., Luzum, B., 2010. IERS conventions (2010). Tech. rep., Bureau International Des Poids Et Mesures Sevres (France), doi: 10.1017/s1743921314005535, doi: 10.1017/s1743921314005535.

Pilato, C.M., Collins-Sussman, B., Fitzpatrick, B.W., 2008. Version Control with Subversion: Next Generation Open Source Version Control. O’Reilly Media Inc, Boston, http://svnbook.red-bean.com/en/1.7/svn-book.pdf.

Prince, P., Dormand, J., 1981. High order embedded Runge-Kutta formulae. J. Comput. Appl. Math. 7 (1), 67–75. https://doi.org/10.1016/0771-050x(81)90010-3.

Ray, J.R., 1996. Measurements of length of day using the Global Positioning System. J. Geophys. Res.: Solid Earth 101 (B9), 2014–2019. https://doi.org/10.1029/96jb01889.

Rebischung, P., Altamimi, Z., Ray, J., Garayt, B., 2016. The IGS contribution to ITRF2014. J. Geodesy 90 (7), 611–630. https://doi.org/10.1007/s00190-016-0897-6.

Ries, P., Dietrich, A., Heflin, M., Hemberger, D., Sibois, A., Sibthorpe, A., Murphy, D., October 2018. Status of IGS14 reprocessing at the JPL IGS analysis center. Presentation at International GNSS Service Workshop, IGS 2018 Workshop (Wuhan, China).

Rodell, M., Houser, P.R., Jambor, U., Gottschalck, J., Mitchell, K., Meng, C.J., Arsenault, K., Cosgrove, B., Radakovich, J., Bosilovich, M., Entin, J.K., Walker, J.P., Lohmann, D., Toll, D., Mar 2004. The Global Land Data Assimilation System. Bull. Am. Meteorol. Soc., vol. 85, 3, pp. 381–394.

Schaer, S., Gurtner, W., Feltens, J., 1998. Ionex: The ionosphere map exchange format version 1. In: Proceedings of the IGS AC workshop, Darmstadt, Germany, Vol. 9.

Seitz, M., Bloßfeld, M., Angermann, D., Schmid, R., Gerstl, M., Seitz, F., 2016. The new DGFI-TUM realization of the ITRS: DTRF2014 (data). Deutsches Geodäisches Forschungsinstitut, Munich.

Shi, C., Zhao, Q., Ge, M., Liu, J., Lou, Y., Zhang, H., Geng, J., Li, M., 2010. Introduction to PANDA software and the latest development for high precision GNSS data processing and application. Presentation at the IGS Analysis Center Workshop, Newcastle upon Tyne. ftp://stella.ncl.ac.uk/pub/IGSposters/Shi.pdf.

Sibthorpe, A., Bertiger, W., Desai, S.D., Haines, B., Harvey, N., Weiss, J. P., 2011. An evaluation of solar radiation pressure strategies for the

GPS constellation. J. Geodesy 85 (8), 505–517. https://doi.org/10.1007/s00190-011-0450-6.

Simons, M., Minson, S.E., Sladen, A., Ortega, F., Jiang, J., Owen, S.E., Meng, L., Ampuero, J.-P., Wei, S., Chu, R., et al., 2011. The 2011 magnitude 9.0 Tohoku-Oki earthquake: Mosaicking the megathrust from seconds to centuries. Science 332 (6036), 1421–1425 https://doi.org/10.1126/science.1206731.

Sovers, O., Border, J., 1990. Observation model and parameter partials for the JPL geodetic GPS modeling software GPSOMC. JPL publication 87–21, Rev 2. Jet Propulsion Laboratory, Pasadena, California.

Sovers, O., Fanselow, J., 1987. Observation model and parameter partials for the JPL VLBI parameter estimation software MASTERFIT. JPL publication 83–39, Rev 3. Jet Propulsion Laboratory, Pasadena, California.

Soycan, M., 2011. A quality evaluation of precise point positioning within the bernese GPS software version 5.0. Arab. J. Sci. Eng. 37 (1), 147–162. https://doi.org/10.1007/s13369-011-0162-5.

Springer, T., Dilssner, F., Escobar, D., Flohrer, C., Otten, O., Svehla, D., Zandbergen, R., 2011. NAPEOS: The ESA/ESOC tool for space geodesy. Eur. Geosci. Union, https://meetingorganizer.copernicus.org/EGU2011/EGU2011-8287.pdf.

Sterling, T.L., Gropp, W., Lusk, E., 2002. Beowulf cluster computing with Linux. MIT press, Cambridge.

Tapley, B.D., Bettadpur, S., Ries, J.C., Thompson, P.F., Watkins, M.M., 2004. GRACE measurements of mass variability in the Earth system. Science 305 (5683), 503–505. https://doi.org/10.1126/science.1099192.

Tapley, B.D., Schutz, B.E., Born, G.H., 2004. Statistical Orbit Determination. Elsevier. https://doi.org/10.1016/b978-0-12-683630-1.x5019-x.

Thomas, J., 1999. An analysis of gravity-field estimation based on intersatellite dual-1-way biased ranging. JPL publication 98–15. Jet Propulsion Laboratory, Pasadena, California.

Thomas, J.B., 1975. Reformulation of the relativistic conversion between coordinate time and atomic time. Astron. J. 80, 405–411. https://doi.org/10.1086/111756.

Tregoning, P., Watson, C., 2009. Atmospheric effects and spurious signals in GPS analyses. J. Geophys. Res.: Solid Earth 114 (B9), B09403. https://doi.org/10.1029/2009jb006344.

Uhlemann, M., Gendt, G., Ramatschi, M., Deng, Z., 2015. GFZ global multi-GNSS network and data processing results. In: International Association of Geodesy Symposia. Springer International Publishing, Cham, pp. 673–679, doi: 10.1007/1345_2015_120, https://doi.org/10.1007/1345_2015_120.

Whitehead, M., Penno, G., Feller, W., Messinger, I., Bertiger, W., Muellerschoen, R., Iijima, B., 1998. A close look at Satloc’s real-time WADGPS system GPS solutions. GPS Solutions 2 (2), 46–63. https://doi.org/10.1007/pl000000036.

Willis, P., Fagard, H., Ferrage, P., Lemoine, F.G., Noll, C.E., Noomen, R., Otten, M., Ries, J.C., Rothacher, M., Soudarin, L., et al., 2010. The international DORIS service (IDS): toward maturity. Adv. Space Res. 45 (12), 1408–1420. https://doi.org/10.1016/j.asr.2009.11.018.

Willis, P., Heflin, M.B., Haines, B.J., Bar-Sever, Y.E., Bertiger, W.I., Mandea, M., 2016. Is the Jason-2 DORIS oscillator also affected by the South Atlantic Anomaly? Adv. Space Res. 58 (12), 2617–2627. https://doi.org/10.1016/j.asr.2016.09.015.

Wu, J., Wu, S., Hajj, G., Bertiger, W., Lichten, S., 1993. Effects of antenna orientation on GPS carrier phase. Manuscripta Geodaetica 18, 91–98.

Wu, S., Bar-Sever, Y., Bertiger, W., Bassiri, S., Hajj, G., Lichten, S., Malla, R., Trinkle, R., Wu, J., 1990. Global positioning system (GPS) precision orbit determination (POD) software design. JPL publication D-725, TOPEX/Poseidon Project. Jet Propulsion Laboratory, Pasadena, California.

Wu, S., Bertiger, W., Border, J., Lichten, S., Sunseri, R., Williams, B., Wolff, P., Wu, J., April 1986. OASIS mathematical description, V. 1.0. JPL publication D-3139, Jet Propulsion Laboratory, Pasadena, California.

Wu, S.-C., Yunck, T.P., Thornton, C.L., 1991. Reduced-dynamic technique for precise orbit determination of low Earth satellites. J. Guid., Control, Dynam. 14 (1), 24–30. https://doi.org/10.2514/3.20600.
Yunck, T., Bertiger, W., Lichten, S., Mannucci, A., Muellerschoen, R., Wu, S., 1989. Robust real-time wide-area differential GPS navigation. United States patent 5,828,336. Jet Propulsion Laboratory, Pasadena, California.

Yunck, T.P., Wu, S.-C., Wu, J.-T., Thornton, C.L., 1990. Precise tracking of remote sensing satellites with the global positioning system. IEEE Trans. Geosci. Remote Sens. 28 (1), 108–116. https://doi.org/10.1109/36.45753.

Zumberge, J., Bertiger, W., 1996. Ephemeris and clock message accuracy. In: Parkinson, B.W., Enge, P., Axelrad, P., Spilker Jr, J. J. (Eds.), Global positioning system: Theory and applications, Volume I. American Institute of Aeronautics and Astronautics, Reston, Ch. 16, pp. 585–599, doi: 10.2514/4.866388, https://doi.org/10.2514/4.866388.