###### Abstract

The GRACE Follow-On satellite mission measures distance variations between the
two satellites in order to derive monthly gravity field maps, indicating mass
variability on Earth on a few 100 km scale due to hydrology, seismology,
climatology and others. This mission hosts two ranging instruments, a
conventional microwave system based on K(a)-band ranging (KBR) and a novel
laser ranging instrument (LRI), both relying on interferometric phase readout.
In this paper we show how the phase measurements can be converted into range
data using a time-dependent carrier frequency (or wavelength) that takes
potential intraday variability in the microwave or laser frequency into
account. Moreover, we analyze the KBR-LRI residuals and discuss which error
and noise contributors limit the residuals at high and low Fourier
frequencies. It turns out that the agreement between KBR and LRI biased range
observations can be slightly improved by considering intraday carrier
frequency variations in the processing. Although the effect is probably small
enough to have little relevance for gravity field determination at the current
precision level, the analysis is of relevance for detailed instrument
characterization and potentially for future more precise missions.

###### keywords:

GRACE Follow-On; instrumentation; laser ranging; microwave ranging

1 1 0 https://doi.org/ Comparing GRACE-FO KBR and LRI Ranging Data with Focus
on Carrier Frequency Variations Comparing GRACE-FO KBR and LRI ranging data
with focus on carrier frequency variations Vitali Müller 1*, Markus Hauk1,2,3,
Malte Misfeldt 1, Laura Müller 1, Henry Wegener 1, Yihao Yan 1, and Gerhard
Heinzel 1* Müller, V.; Hauk, M.; Misfeldt, M.; Müller, L.;Wegener, H.;Yan, Y.;
Heinzel, G.; Müller, V.; Hauk, M.; Misfeldt, M.; Müller, L.;Wegener, H.;Yan,
Y.; Heinzel, G.; Correspondence<EMAIL_ADDRESS>(V.M.) and
<EMAIL_ADDRESS>(G.H.)

##  1 Introduction

The Gravity Recovery and Climate Experiment (GRACE) satellites have measured
temporal variations in Earth’s gravity field from 2002 until 2017. The
successor mission GRACE Follow-On (GRACE-FO), launched in 2018, continues the
valuable dataset of monthly gravity field maps Landerer et al. (2020), which
are used for instance in the fields of climate research (Tapley et al., 2019),
hydrology Landerer and Swenson (2012) and seismology Han et al. (2006).
Prominent GRACE(-FO) results have quantified the ice mass loss in regions such
as Greenland and Antarctic over the past two decades Ciracì et al. (2020) and
allowed to attribute the effect of mass influx into the oceans in the observed
global mean sea level rise Chen et al. (2020). Some changes in the gravity
field are likely caused by groundwater depletion Chen et al. (2014). The
GRACE(-FO) results lead to more than 2600 publications, and many of them are
frequently cited in reports of the Intergovernmental Panel on Climate Change
(IPCC).

The primary measurement that contains the gravity field information in a
GRACE-like mission is the biased inter-satellite range, which was measured to
micrometer precision in GRACE using a microwave K-band ranging system (KBR).
GRACE-FO hosts two ranging instruments, the KBR and a novel and more precise
laser ranging instrument (LRI, Kornfeld et al. (2019); Ghobadi-Far et al.
(2020)). The ranging data of either KBR or LRI, or a combination of both, is
processed in the course of gravity field retrieval together with various other
measurements, such as attitude and orbit data as well as with measurements of
non-gravitational accelerations from accelerometers.

LRI and KBR measure the inter-satellite distance variations in GRACE-FO in
parallel. Future missions such as the European-lead Next Generation Geodesy
Mission (NGGM, Nicklaus et al. (2020)) and the US-lead Mass Change Mission
Wiese et al. (2020), with counterpart GRACE-I(carus) Flechtner (2020) in
Germany, are currently being studied or developed. Since these new missions
will solely rely on laser-based ranging, it is important to understand what is
currently limiting the ranging data.

The GRACE-FO data discussed here is publicly available in terms of raw
observations (level1a), in processed observations (level1b) and in terms of
gravity field products (level2). The unique setup of two almost independent
ranging systems allows us to study the differential instrument behavior and
small instrument errors at frequency regions which are usually dominated by
signal several orders of magnitude larger than the expected noise.

In this paper, we compare the LRI and KBR ranging data in the spectral domain,
revisit the main error and noise contributors in both instruments and attempt
to explain the current level of KBR-LRI residuals. We pay special attention to
carrier frequency variations in both instruments and on the approach to
account for them when converting the observed phase (level1a) to a range
(level1b).

##  2 Phase Observable

The phase of microwave or optical radiation can be described in a general-
relativistic context as

$\displaystyle\Phi=\nu_{0}\cdot\tau_{\mathrm{osc}}(t)+q_{0}=\nu_{0}\cdot[\tau(t)+\tau_{\mathrm{FV}}(\tau(t))]+q_{0},$
(1)

where $t$ is the coordinate time or GPS time, $\nu_{0}$ the nominal frequency,
and $q_{0}$ is a constant that cannot be resolved in general using
interferometry, because phase observations exhibit an integer-ambiguity.
$\tau_{\mathrm{osc}}$ is a virtual oscillator time that is given by the sum of
the proper time $\tau$ of the satellite or emitting device (e.g. laser) and a
contributor $\tau_{\mathrm{FV}}$ describing frequency variations or, in other
words, deviations of the oscillator time from proper time due to
imperfections, noise, etc. The oscillator time for the KBR is realized inside
the instrument processing unit (IPU) of the microwave ranging instrument (MWI)
assembly and can be called the ultra-stable oscillator (USO) time
$\tau^{\textrm{USO}}$ or the IPU receiver time. The difference between USO
time and GPS time is the clock error, which is estimated during precise orbit
determination (POD) together with the satellite’s position and velocity state.
The clock error is reported as time series in the CLK1B data product as a
eps_Time ($\epsilon_{\textrm{time}}$), i.e.

$\displaystyle\epsilon_{\textrm{time}}=t-\tau^{\textrm{USO}}(t).$ (2)

$\epsilon_{\textrm{time}}$ contains the proper time, clock variations and a
frequency offset, as can be seen from fig. 1. The proper time for a satellite
can easily be computed by numerically integrating Burt et al. (2021)

$\displaystyle 1-\frac{\mathrm{d}\tau}{\mathrm{d}t}$
$\displaystyle=\frac{GM}{rc_{0}^{2}}-\frac{GM}{rc_{0}^{2}}\cdot
J_{2}\cdot\Bigl{(}\frac{a_{e}}{r}\Bigr{)}^{2}\cdot\frac{3z^{2}-r^{2}}{2r^{2}}+\frac{v^{2}}{2c_{0}^{2}}+\frac{\phi_{0}}{c_{0}^{2}},$
(3)

where $c_{0}$ is speed of light, $GM$ is Earth’s gravitational parameter,
$a_{e}$ is Earth’s mean readius, $J_{2}$ is Earths oblateness coefficient ,
$r$ is the distance between geocenter and satellite, $z$ is the z-component of
the position vector, $v$ is the satellite velocity and $\phi_{0}$ is a
potential constant.

The frequency offsets of the USOs w.r.t. its nominal frequency yield y-offsets
in the left plot of fig. 1. As apparent from the legend, they have values of
$7.4\text{\,}\mathrm{p}\mathrm{p}\mathrm{b}$ on satellite 1 (GF-1) and
$6.6\text{\,}\mathrm{p}\mathrm{p}\mathrm{b}$ on satellite 2 (GF-2) on that
particular day. The oscillator time of the laser is not accessible and is
usually expressed in terms of the laser frequency. The instantaneous frequency
reads

$\displaystyle\nu(\tau):=\frac{\mathrm{d}\Phi}{\mathrm{d}\tau}(\tau)=\nu_{0}\cdot\left(1+\frac{\mathrm{d}\tau_{\mathrm{FV}}}{\mathrm{d}\tau}\right),$
(4)

where $\nu_{0}\cdot\mathrm{d}\tau_{\mathrm{FV}}/\mathrm{d}\tau$ describes the
frequency variations or noise with respect to a nominal frequency $\nu_{0}$.

Figure 1: (Left:) Typical time-domain plot of the proper time $\tau$ (red)
determined from the satellite orbit product and the CLK1B
$\epsilon_{\textrm{time}}$(light and dark blue), both in terms of the first
derivative. The proper time is clearly visible in $\epsilon_{\textrm{time}}$.
(Right): Amplitude spectral density (ASD) of the time-series shown on the left
side together with an estimate of the USO stability based on Thomas (1999).
The traces labelled refined USO model and POD/Clock noise model are described
in sec. 6.4. See fig. 6 for a representation as Allan deviations. The
Equivalent Noise Bandwidth (ENBW) given here and in other spectra can be used
to convert the spectral density units to spectrum (rms) units Heinzel et al.
(2002).

The phase $\Phi$ is invariant under Lorentz transformations Shiozawa (2004),
which means it does not change its value when transformed from the local
system of the satellite to the geocentric celestial reference system (GCRS),
e.g.

$\displaystyle\Phi(\tau_{\textrm{A}},\vec{r}_{\textrm{A}})=\Phi(t_{\textrm{A}},\vec{r}_{\textrm{A}}^{~{}\prime}),$
(5)

where $(\tau_{\textrm{A}},\vec{r}_{\textrm{A}})$ is the four-valued event in
the local satellite system, while
$(t_{\textrm{A}},\vec{r}_{\textrm{A}}^{~{}\prime})$ denotes the same event in
the GCRS.

The frequency is not relativistically invariant, since it appears in the GCRS
as

$\displaystyle\nu(t):=\frac{\mathrm{d}\Phi}{\mathrm{d}\tau}(\tau(t))\cdot\frac{\mathrm{d}\tau}{\mathrm{d}t}(t)=\nu_{0}\dot{\tau}(t)\cdot\left(1+\frac{\mathrm{d}\tau_{\mathrm{FV}}}{\mathrm{d}\tau}(\tau(t))\right),$
(6)

where the short form for the derivative
$\dot{\tau}=\mathrm{d}\tau/\mathrm{d}t$ was used. We note that even if a
perfect laser system would produce a constant freqeuncy in the satellite
frame, the frequency is varying in other systems like the geocentric one used
for gravity field determination.

The optical phase is the time-integral of the optical frequency, i.e.

$\displaystyle\Phi=\int_{\tau_{0}}^{\tau}\nu(\tau^{\prime})~{}\mathrm{d}\tau^{\prime}=\int_{t(\tau_{0})}^{t(\tau)}\nu(t^{\prime})~{}\mathrm{d}t^{\prime}.$
(7)

##  3 LRI Ranging Phase

For the LRI, the ranging phase $\varphi_{\textrm{TWR}}$ is called a two-way
ranging (TWR) quantity in the GRACE-FO project and it is determined from the
phase difference between transponder $\varphi_{\textrm{T}}$ and reference
$\varphi_{\textrm{R}}$ satellite Wen et al. (2019). The phase measurements on
both satellites are performed with phase meters. The transponder satellite
uses a phase-locked loop to lock its laser to the incoming light using a high
gain and high bandwidth control loop Abich et al. (2019). This means the
transponder phase is a rather trivial measurement given by a phase ramp with
constant slope (Müller, 2017, eq. (2.207))

$\displaystyle\varphi_{\textrm{T}}=$10\text{\,}\mathrm{M}\mathrm{H}\mathrm{z}$\cdot\tau^{\textrm{USO}}_{\textrm{T}},\quad[\varphi_{\textrm{T}}]=\textrm{cycles}.$
(8)

On the reference satellite, the laser is stabilized using an optical cavity
and the phase meter output $\varphi_{\textrm{R}}$ contains both the ranging
information and the transponder phase ramp. Since the same phase ramp is
present on both satellites, it cancels out in the differential measurement
$\varphi_{\textrm{TWR}}$ and one can write the TWR phase as Yan et al. (2021)

$\displaystyle\varphi_{\textrm{TWR}}=\varphi_{\textrm{T}}(t-\Delta
t_{\textrm{TR}})-\varphi_{\textrm{R}}(t)=\Phi_{\textrm{R}}(t)-\Phi_{\textrm{R}}(t-\Delta
t_{\textrm{RTR}})+q_{\textrm{R}},$ (9)

where $\Phi_{\textrm{R}}$ is the optical phase of the reference laser in units
of cycles and $\Delta t_{\textrm{RTR}}$ is the round-trip light propagation
time from reference to transponder and back to the reference satellite. The
round-trip propagation time $\Delta t_{\textrm{RTR}}$ contains the ranging
information, while the propagation time from transponder to reference
satellite $\Delta t_{\textrm{TR}}$ can be derived from GPS orbit data and is
needed to account for the propagation of the phase value from transponder to
reference. The constant $q_{\textrm{R}}$ arises from the inability of
interferometers to measure the absolute phase or distance.

Using the definition of the optical phase from eq. (7), one obtains

$\displaystyle\varphi_{\textrm{TWR}}$
$\displaystyle=\int^{\tau}_{\tau-\Delta\tau_{\textrm{RTR}}}\nu_{\textrm{R}}(\tau^{\prime})~{}\mathrm{d}\tau^{\prime}+q_{\textrm{R}}$
(10)
$\displaystyle=\nu_{\textrm{R}}\cdot[\tau-\tau+\Delta\tau_{\textrm{RTR}}(\tau)+\tau_{\mathrm{FV}}(\tau)-\tau_{\mathrm{FV}}(\tau-\Delta\tau_{\textrm{RTR}})]+q_{\textrm{R}}$
(11)
$\displaystyle\approx\nu_{\textrm{R}}\cdot\Delta\tau_{\textrm{RTR}}+\nu_{\textrm{R}}\cdot\frac{\mathrm{d}\tau_{\mathrm{FV}}}{\mathrm{d}\tau}\cdot\Delta\tau_{\textrm{RTR}}+q_{\textrm{R}}$
(12)
$\displaystyle=\nu_{\textrm{R}}(\tau)\cdot\Delta\tau_{\textrm{RTR}}(\tau)+q_{\textrm{R}}=\nu_{\textrm{R}}(t)\cdot\Delta
t_{\textrm{RTR}}(t)+q_{\textrm{R}},$ (13)

where the last line shows the phase observable in terms of proper or local
quantities on the satellite ($\nu(\tau)\cdot\Delta\tau_{\textrm{RTR}}(\tau)$)
and in terms of GCRS quantities ($\nu(t)\cdot\Delta t_{\textrm{RTR}}(t)$). We
use $\varphi$ to indicate phase quantities that have a slope in the radio-
frequency domain (kHz or MHz), i.e. which are tracked and recorded on GRACE-
FO, while the phase quantities $\Phi$ have phase slopes in the THz or GHz
domain, which cannot be recorded.

The LRI level1a (LRI1A) data of GRACE-FO contains the raw phases of reference
$\varphi_{\textrm{R}}$ and transponder $\varphi_{\textrm{T}}$ satellite, while
the LRI1B product contains the raw biased range $\rho_{\textrm{TWR}}$ and the
light-time correction $\rho_{\textrm{LTC}}$. Adding the raw range and the LTC
yields the corrected biased range

$\displaystyle\rho_{\textrm{LRI}}$
$\displaystyle=\rho_{\textrm{TWR}}(\varphi_{\textrm{TWR}})+\rho_{\textrm{LTC,TWR}}+\rho_{\textrm{TTL}}+\rho_{\textrm{media}}$
(14)
$\displaystyle\approx\rho_{\textrm{TWR}}(\varphi_{\textrm{TWR}})+\rho_{\textrm{LTC,TWR}}=|\vec{r}_{\textrm{A}}(t)-\vec{r}_{\textrm{B}}(t)|+\mathrm{bias}$
(15)

which corresponds to the biased instantaneous Euclidean distance between
satellites in the GCRS and that serves as input for the gravity field
recovery. The effects from tilt-to-length coupling $\rho_{\textrm{TTL}}$ and
from the light propagation in a medium $\rho_{\textrm{media}}$ are negligible
for LRI, as will be shown in sec. 6. The light-time correction is computed as

$\displaystyle\rho_{\textrm{LTC,TWR}}:=|\vec{r}_{\textrm{A}}^{\textrm{~{}GPS}}(t)-\vec{r}_{\textrm{B}}^{\textrm{~{}GPS}}(t)|-\frac{c_{0}\Delta
t_{\textrm{RTR}}^{\textrm{GPS}}(t)}{2}$ (16)

from GPS orbit data of the satellites and using formulas for the light time
propagation time ($\Delta t_{\textrm{RTR}}^{\textrm{GPS}}$) that depend on the
orbit data Yan et al. (2021). The magnitude of $\rho_{\textrm{LTC,TWR}}$ is
below a millimeter such that the GPS accuracy of a few millimeter over the
inter-satellite baseline of approx. $200\text{\,}\mathrm{k}\mathrm{m}$ is
sufficient in order to not limit the data quality of $\rho_{\textrm{LRI}}$.

The derivation of the range $\rho$ from the (round-trip) phase $\varphi$ is
trivial, if one assumes a constant frequency $\nu$

$\displaystyle\rho(t)=\frac{\lambda}{2}\cdot\varphi(t)=c_{0}\cdot\frac{\varphi(t)}{2\nu}.$
(17)

The wavelength $\lambda$ or frequency $\nu$ serves as a conversion factor from
phase with units of cycles to a half round-trip range in units of meters. The
currently available (v04) level1b data products of LRI and KBR assume a day-
wise constant frequency such that the equation is applicable and used to
convert data from individual days.

However, such an approach has the disadvantage that the range data might have
some offsets and slope changes at day bounds due to the daily changes of the
frequency. Moreover, we anticipate that intraday frequency variations $\nu(t)$
are a noticeable error term in both laser and microwave interferometers. Thus,
we investigated methods to convert the full round-trip phase
$\varphi_{\textrm{TWR}}$ to a half round-trip range $\rho_{\textrm{TWR}}$ when
frequency variations are perfectly known. This means we are searching for
functions $\rho_{\textrm{TWR}}(\varphi_{\textrm{TWR}})$ that, for the LRI,
yield

$\displaystyle\rho_{\textrm{TWR}}(\varphi_{\textrm{TWR}})=c_{0}\frac{\Delta
t_{\textrm{RTR}}(t)}{2}+\textrm{const}.$ (18)

Once the error-free equations are understood, we address errors in the
knowledge of $\nu$ in sec. 3.4.

###  3.1 Formula 1: Phase-Frequency Ratio

The direct naive approach using eq. (17) with a time-dependent frequency and
approximation in eq. (13) yields

$\displaystyle\rho_{\textrm{TWR,1}}(t):=c_{0}\cdot\frac{\varphi_{\textrm{TWR}}(t)}{2\nu_{\textrm{R}}(t)}\approx
c_{0}\frac{\Delta
t_{\textrm{RTR}}(t)}{2}+c_{0}\frac{q_{\textrm{R}}}{2\nu_{\textrm{R}}(t)}.$
(19)

The divisor is the apparent frequency in the GCRS $\nu(t)$, since $\nu(\tau)$
would result in a contracted range $c_{0}\Delta\tau_{\textrm{RTR}}$. The naive
formula $\rho_{\textrm{TWR,1}}$ yields the correct result of eq. (18), if
$q_{\textrm{R}}=0$, i.e. if the phase measurement $\varphi_{\textrm{RTR}}$ is
an absolute measurement without integer ambiguity due to indistinguishable
cycles. In actual interferometry, the phase time-series
$\varphi_{\textrm{RTR}}(t)$ exhibits an unknown arbitrary bias and it is
debiased in order to avoid large numerical values, loss of numerical precision
due to floating-point arithmetic or overflow in finite-sized registers of a
computer. We assume without loss of generality that phase time-series start at
0 cycles at an initial epoch, i.e. $\varphi_{\textrm{RTR}}(t=0)=0$ in a
stretch of data that is being converted from phase to range. This implies that

$\displaystyle q_{\textrm{R}}$
$\displaystyle=-\nu_{\textrm{R}}(t=0)\cdot\Delta t_{\textrm{RTR}}(t=0)$ (20)
$\displaystyle\rho_{\textrm{TWR,1}}(t)$ $\displaystyle\approx
c_{0}\frac{\Delta t_{\textrm{RTR}}(t)}{2}-c_{0}\Delta
t_{\textrm{RTR}}(0)\cdot\frac{\nu_{\textrm{R}}(t=0)}{2\nu_{\textrm{R}}(t)}$
(21)

where the first summand on the right-hand side in line (21) is the desired
ranging signal and the second term describes laser frequency variations that
couple via the baseline length $c_{0}\Delta
t_{\textrm{RTR}}(0)/2\approx$220\text{\,}\mathrm{k}\mathrm{m}$$. Since eq.
(21) shows the deviation of $\rho_{\textrm{TWR,1}}$ from the correct result,
we can directly define an improved formula as

$\displaystyle\rho_{\textrm{TWR,1corr}}(t)$
$\displaystyle:=c_{0}\cdot\frac{\varphi_{\textrm{TWR}}(t)}{2\nu_{\textrm{R}}(t)}+c_{0}\Delta
t^{\textrm{GPS}}_{\textrm{RTR}}(0)\cdot\left(\frac{\nu_{\textrm{R}}(0)}{2\nu_{\textrm{R}}(t)}-\frac{1}{2}\right),$
(22) $\displaystyle\approx
c_{0}\cdot\frac{\varphi_{\textrm{TWR}}(t)}{2\nu_{\textrm{R}}(t)}-\frac{c_{0}\Delta
t^{\textrm{GPS}}_{\textrm{RTR}}(0)}{2}\cdot\left(\frac{\nu_{\textrm{R}}(t)-\nu_{\textrm{R}}(0)}{\nu_{\textrm{R}}(0)}\right),$
(23)

where the second correction term on the right-hand side in line (22) is
usually much smaller than the first ranging term, as the expression in the
bracket is close to zero due to the introduced constant $-1/2$. The round-trip
propagation time at the initial epoch received the superscript GPS
($c_{0}\Delta\tau_{\textrm{RTR}}^{\mathrm{GPS}}(0)$) in order to indicate that
the value can be derived from GPS orbit data and algorithms to compute the
light propagation time Yan et al. (2021). The millimeter to centimeter
precision of GPS is sufficient to derive the absolute distance $c_{0}\Delta
t_{\textrm{RTR}}^{\mathrm{GPS}}(0)\approx$440\text{\,}\mathrm{k}\mathrm{m}$$
without spoiling the high precision from LRI in $\rho_{\textrm{TWR,1corr}}$,
since the correction term is small due to
$\nu_{\textrm{R}}(0)/\nu_{\textrm{R}}(t)\approx 1$. The approximation in line
(23) was based on $1/(1+x)\approx 1-x$, with
$x=\nu_{\textrm{R}}(t)/\nu_{\textrm{R}}(0)-1$, and it illustrates more clearly
that the correction term is the product of inter-satellite distance and
fractional frequency deviations.

###  3.2 Formula 2: Integral of Differentiated Phase

An alternative way to compute the half round-trip range from the round-trip
ranging phase $\varphi_{\textrm{TWR}}$ has been suggested as (Yan et al.,
2021, eq. (60))

$\displaystyle\rho_{\textrm{TWR,2}}\approx
c_{0}\int^{t}_{0}\frac{\mathrm{d}\varphi_{\textrm{TWR}}(t^{\prime})/\mathrm{d}t^{\prime}}{2\nu_{\textrm{R}}(t^{\prime})}~{}\mathrm{d}t^{\prime}$
(24)

though it is easy to see that the equation is not exact, since

$\displaystyle\frac{\mathrm{d}\varphi_{\textrm{TWR}}}{\mathrm{d}t}$
$\displaystyle=\frac{\mathrm{d}missing}{\mathrm{d}t}\int_{t-\Delta
t}^{t}\nu_{\textrm{R}}(t^{\prime})~{}\mathrm{d}t^{\prime}=\nu_{\textrm{R}}(t)-\nu_{\textrm{R}}(t-\Delta
t)\cdot(1-\mathrm{d}\Delta t/\mathrm{d}t)$ (25)
$\displaystyle=\nu_{\textrm{R}}(t-\Delta t)\cdot\mathrm{d}\Delta
t/\mathrm{d}t+\nu_{\textrm{R}}(t)-\nu_{\textrm{R}}(t-\Delta t)$ (26)
$\displaystyle\approx\nu_{\textrm{R}}(t)\cdot\mathrm{d}\Delta
t/\mathrm{d}t+(1-\mathrm{d}\Delta
t/\mathrm{d}t)\cdot\dot{\nu}_{\textrm{R}}(t)\cdot\Delta t$ (27)

The first term in the line (26) contains the ranging information that is
multiplied with the optical frequency at the emission event and not at the
reception event as the divisor in eq. (24) suggests. Moreover, the second and
third terms indicate needed correction terms that are missing in eq. (24).
Thus, the exact formula to convert the round-trip phase
$\varphi_{\textrm{TWR}}$ to a half round-trip range reads

$\displaystyle\rho_{\textrm{TWR,2exact}}(t)$
$\displaystyle:=\frac{c_{0}}{2}\int^{t}_{0}\frac{\mathrm{d}\varphi_{\textrm{TWR}}(t^{\prime})/\mathrm{d}t^{\prime}}{\nu_{\textrm{R}}\left(t^{\prime}-\Delta
t_{\textrm{RTR}}(t^{\prime})\right)}-\left(\frac{\nu_{\textrm{R}}(t^{\prime})}{\nu_{\textrm{R}}\left(t^{\prime}-\Delta
t_{\textrm{RTR}}(t^{\prime})\right)}-1\right)~{}\mathrm{d}t^{\prime}$ (28)
$\displaystyle=\frac{c_{0}\Delta t_{\textrm{RTR}}(t)}{2}-\frac{c_{0}\Delta
t_{\textrm{RTR}}(0)}{2}\ .$ (29)

The term

$\displaystyle\frac{c_{0}}{2}\int_{0}^{t}\frac{\nu_{\textrm{R}}(t^{\prime})}{\nu_{\textrm{R}}(t^{\prime}-\Delta
t_{\mathrm{RTR}}(t^{\prime}))}-1\,\mathrm{d}t^{\prime}$
$\displaystyle\quad\approx\frac{c_{0}}{2}\int_{0}^{t}\frac{\dot{\nu}_{\textrm{R}}(t^{\prime})}{\nu_{\textrm{R}}(t^{\prime})}\cdot\Delta
t_{\mathrm{RTR}}(t^{\prime})\,\mathrm{d}t^{\prime}\approx\frac{c_{0}}{2}\langle\Delta
t_{\mathrm{RTR}}\rangle\cdot\frac{\log(\nu_{\textrm{R}}(t))}{\log(\nu_{\textrm{R}}(0))}\approx\frac{c_{0}}{2}\langle\Delta
t_{\mathrm{RTR}}\rangle\cdot\frac{\nu_{\textrm{R}}(t)-\nu_{\textrm{R}}(0)}{\nu_{\textrm{R}}(0)}$
(30)

approximately describes the range correction due to laser phase changes on the
reference satellite while the laser light propagates forth and back between
satellites with average propagation time $\langle\Delta
t_{\mathrm{RTR}}\rangle$. The last approximation of this expression resembles
the correction term shown in eq. (23).

A handy approximation of the exact form in eq. (28) with all quantities being
evaluated at the same time is given as

$\displaystyle\rho_{\textrm{TWR,2approx}}=c_{0}\int^{t}_{0}\frac{\mathrm{d}\varphi_{\textrm{TWR}}(t^{\prime})/\mathrm{d}t^{\prime}}{2\nu_{\textrm{R}}\left(t^{\prime}\right)}-\left(1-\frac{\mathrm{d}\Delta
t_{\textrm{RTR}}^{\textrm{GPS}}}{\mathrm{d}t}\right)\cdot\frac{\dot{\nu}_{\textrm{R}}(t^{\prime})}{2\nu_{\textrm{R}}\left(t^{\prime}\right)}\Delta
t_{\textrm{RTR}}^{\textrm{GPS}}(t^{\prime})~{}\mathrm{d}t^{\prime},$ (31)

which is based on the approximation shown in eq. (27). The light propagation
time in the correction term received a superscript GPS to indicate that this
quantity can be obtained from the orbit product, as in eq. (22).

###  3.3 Accuracy of Approximations

In order to assess the accuracy of the formulas given in the two previous
subsections, we have used an analytical model for the optical frequency
$\nu_{\textrm{R}}$ and the true inter-satellite distance, computed the optical
phase $\Phi_{\textrm{R}}$ and, retrieved from that, a simulated range
$\rho_{\textrm{TWR}}$. The models and numerical values are given in Appendix
A.1 and emulate conditions present in GRACE-FO. Based on figure 2, one can
state that the simplistic formula without frequency correction term (eq. (21))
yields a linearly drifting error of approx. 68 µm per day (dashed red trace)
if the laser frequency is linearly drifting by 87 kHz per day
($\nu_{d}/\nu_{0}=3.6\cdot 10^{-15}$ 1/sec with $\nu_{0}=$ 282 THz). An
oscillating frequency with a fractional peak amplitude of
$\nu_{1}/\nu_{0}=4\cdot 10^{-12}$, which is approximately the modulation
amplitude in the laser (and microwave) frequency from the proper time, yields
errors in the order of 1 µm at the 1/rev oscillation frequency (solid trace).

The $\rho_{\textrm{TWR,2approx}}$ expression has no error when the frequency
$\nu_{\textrm{R}}$ is drifting (dashed blue trace not visible, since zero not
present in log-scale plot), and slightly lower errors than
$\rho_{\textrm{TWR,1corr}}$ when $\nu_{\textrm{R}}$ is oscillating (cf. solid
blue vs. solid green trace) . However, the accuracy of the
$\rho_{\textrm{TWR,1corr}}$ is already sufficient for GRACE-FO like
applications, since the error is at the picometer level for a one day long
time-series in the drifting and oscillating case.

The error for $\rho_{\textrm{TWR,2exact}}$ from eq. (28) cannot be shown in
the log-scale figure because the equation is exact with zero error. The exact
solution is our recommended way to transform the phase to range.

Figure 2: Evaluation of accuracy for different formulas and approximations
$\rho_{\textrm{TWR,...}}(t)$ shown in red for eq. (21), in green for eq. (22)
and in blue for eq. (27). The simulation parameters for a time-series of
$86\,400\text{\,}\mathrm{s}\mathrm{e}\mathrm{c}$ length are given in Appendix
A.1 and Table 2 therein. The error is defined as
$|\rho_{\textrm{TWR,...}}(t)-(L(t)-L_{0})|$, where $(L(t)-L_{0})$ denotes the
true (debiased) range (cf. Appendix A.1). The dashed traces assume a drifting
frequency $\nu_{\textrm{R}}$, while the solid lines assume an oscillating
frequency $\nu_{\textrm{R}}$. The dashed-blue trace is not visible because the
error is zero. $\nu_{\textrm{R}}$ is assumed to be known without any error.

###  3.4 Model Errors: Scale Factor and Light Propagation Time

The formulas used for converting phase to range considered so far error-free
knowledge of the laser frequency ($\nu_{\textrm{R}}$) and light travel time
($\Delta t_{\textrm{RTR}}^{\textrm{GPS}}$ = $\Delta t_{\textrm{RTR}}$). In
order to assess the effect of inaccuracies in these models, one can replace

$\displaystyle\nu_{\textrm{R}}(t)\longrightarrow\nu_{\textrm{R}}(t)+\delta\nu_{\textrm{R}}(t)$
(32) $\displaystyle\Delta
t_{\textrm{RTR}}^{\textrm{GPS}}(t)\longrightarrow\Delta
t_{\textrm{RTR}}^{\textrm{GPS}}(t)+\Delta
t_{\textrm{RTR}}^{\textrm{GPS,e}}(t)$ (33)

such that $\delta\nu_{\textrm{R}}$ and $\Delta
t_{\textrm{RTR}}^{\textrm{GPS,e}}$ account for errors in the absolute laser
frequency and propagation time, respectively. There are two ways to apply the
laser frequency errors. One way is to consider them in the phase domain (cf.
eq. (9))

$\displaystyle\varphi_{\textrm{TWR}}^{\textrm{Noisy}}(t)-\varphi_{\textrm{TWR}}^{\textrm{NoiseFree}}(t)$
$\displaystyle=\Phi_{\textrm{R,Noise}}(t)-\Phi_{\textrm{R,Noise}}(t-\Delta
t_{\textrm{RTR}})=\int_{t-\Delta
t_{\textrm{RTR}}}^{t}\delta\tilde{\nu}_{\textrm{R}}(t^{\prime})~{}\mathrm{d}t^{\prime}$
(34) $\displaystyle\approx\delta\tilde{\nu}_{\textrm{R}}(t)\cdot\Delta
t_{\textrm{RTR}}(t),$ (35)

where we just renamed $\delta{\nu}$ into $\delta\tilde{\nu}_{\textrm{R}}$. Eq.
(34) illustrates that the (cavity) phase variations $\Phi_{\textrm{R,Noise}}$
are suppressed in the TWR combination due to the self-comparison with a
delayed instance, which can be expressed in the frequency domain with a
multiplier $(1-e^{-2\pi if\Delta t_{\textrm{RTR}}})$. Similar equations apply
to the KBR, where this multiplier is called a _ranging filter_ (Kim, 2000, Eq.
(4.20)) (Thomas, 1999, Eq. (3.38)). The second way of applying the errors in
line (32) and (33) is to consider them in the quantities used to convert the
phase to range. When using both ways at the same time, and comparing the range
with model errors ($\rho^{\textrm{mod}}_{\textrm{...}}$) to the original
formula, we obtain after expansion to leading order111For the sake of
simplicity, we additionally neglected terms that contain $\Delta
t_{\textrm{RTR}}^{\textrm{GPS,e}}$ in eq. (38).

$\displaystyle\rho_{\textrm{TWR,1corr}}^{\textrm{mod}}-\rho_{\textrm{TWR,1corr}}$
$\displaystyle\approx+\frac{c_{0}}{2}\tilde{\epsilon}_{\textrm{SCF}}(t)\cdot\Delta
t_{\textrm{RTR}}(t)-\frac{c_{0}}{2}\epsilon_{\textrm{SCF}}(t)\cdot\left(\Delta
t_{\textrm{RTR}}(t)+\Delta
t_{\textrm{RTR}}(0)\cdot\left(\frac{\nu_{\textrm{R}}(0)}{\nu_{\textrm{R}}(t)}-1\right)\right)$
$\displaystyle\qquad+\frac{c_{0}}{2}\cdot\Delta
t_{\textrm{RTR}}^{\textrm{GPS,e}}(0)\cdot\left(\frac{\nu_{\textrm{R}}(0)}{\nu_{\textrm{R}}(t)}-1\right)+\frac{c_{0}}{2}\cdot\epsilon_{\textrm{SCF}}(0)\cdot\Delta
t_{\textrm{RTR}}(0)\cdot\frac{\nu_{\textrm{R}}(0)}{\nu_{\textrm{R}}(t)}+\textrm{const.}$
(36)
$\displaystyle\approx+\frac{c_{0}}{2}\tilde{\epsilon}_{\textrm{SCF}}(t)\cdot\Delta
t_{\textrm{RTR}}(t)-\frac{c_{0}}{2}\epsilon_{\textrm{SCF}}(t)\cdot\Delta
t_{\textrm{RTR}}(t)+\textrm{const.}$ (37)
$\displaystyle\rho_{\textrm{TWR,2approx}}^{\textrm{mod}}-\rho_{\textrm{TWR,2approx}}$
$\displaystyle\approx-\frac{c_{0}}{2}\Delta
t_{\textrm{RTR}}\cdot(\epsilon_{\textrm{SCF}}-\tilde{\epsilon}_{\textrm{SCF}})+\int_{0}^{t}\frac{c_{0}\Delta
t_{\textrm{RTR}}}{2}\cdot\left(\Delta\dot{t}_{\textrm{RTR}}~{}\frac{\delta\dot{\nu}_{\textrm{R}}}{\nu_{\textrm{R}}}-\epsilon_{\textrm{SCF}}\frac{\dot{\nu}_{\textrm{R}}}{\nu_{\textrm{R}}}\right)\mathrm{d}t^{\prime}+\textrm{const.}$
(38)
$\displaystyle\approx+\frac{c_{0}}{2}\tilde{\epsilon}_{\textrm{SCF}}(t)\cdot\Delta
t_{\textrm{RTR}}(t)-\frac{c_{0}}{2}\epsilon_{\textrm{SCF}}(t)\cdot\Delta
t_{\textrm{RTR}}(t)+\textrm{const.},$ (39)

where the scale factor (SCF) error was introduced in two versions as

$\displaystyle\epsilon_{\textrm{SCF}}(t)=\frac{\delta\nu_{\textrm{R}}(t)}{\nu_{\textrm{R}}(t)},\qquad\tilde{\epsilon}_{\textrm{SCF}}(t)=\frac{\delta\tilde{\nu}_{\textrm{R}}(t)}{\nu_{\textrm{R}}(t)}.$
(40)

Both formulas ($\rho_{\textrm{TWR,1corr}}$ and $\rho_{\textrm{TWR,2approx}}$)
and both ways to apply the error yield to leading order the same result: the
scale factor error ($\epsilon_{\textrm{SCF}}$ and
$\tilde{\epsilon}_{\textrm{SCF}}$) arising from the limited knowledge of the
frequency $\nu_{\textrm{R}}$ produces ranging errors proportional to the
satellite separation $c_{0}\Delta t_{\textrm{RTR}}(t)/2\approx L(t)$. This is
the expected result and the usual coupling of frequency noise into the range
measurement (Sheard et al., 2012, eq. (11)) at least to leading order. The
exact error coupling can become rather complicated.

In order to assess the magnitude of scale factor errors, we once more assume
an oscillatory and a drifting behavior as in the previous subsection. For the
drift, we reuse the value of $3.6\cdot 10^{-15}$ 1/sec (cf. table 2), which
results in a scale factor change $\Delta\epsilon_{\textrm{SCF}}=3.1\cdot
10^{-10}$ per day. Such a scale factor error is present if the cavity
resonance frequency in the LRI drifts by 87 kHz per day without knowledge
about the drift, which we regard to be potentially possible. The resulting
ranging error per day is then

$\displaystyle\langle
L\rangle\cdot\Delta\epsilon_{\textrm{SCF}}=220\,\textrm{km}\cdot 3.1\cdot
10^{-10}=68\,\textrm{\textmu m}.$ (41)

A realistic error for the oscillatory behavior arises when the proper time
induced modulation is neglected, as currently in the regular LRI data
processing. This means the sinusoidal modulation with peak amplitude
$\nu_{\textrm{R}}(t)/\langle\nu_{\textrm{R}}\rangle\approx 4\cdot 10^{-12}$ at
1/rev frequency is neglected, which implies an oscillating scale factor error
$\epsilon_{\textrm{SCF}}\approx 4\cdot 10^{-12}$ at the orbital frequency. The
resulting ranging error is

$\displaystyle\langle
L\rangle\cdot\epsilon_{\textrm{SCF}}\approx$220\text{\,}\mathrm{k}\mathrm{m}$\cdot
4\cdot 10^{-12}\approx 0.9\,\textrm{\textmu m}_{\textrm{peak}}\textrm{~{}at
1/rev frequency}.$ (42)

The calculated changes of $68\,\textrm{\textmu m}$ per day and
$0.9\,\textrm{\textmu m}_{\textrm{peak}}$ at 1/rev frequency are actually
apparent from the red traces in fig. 2 for the error from eq. (17). This
equation exhibits a time-dependent frequency $\nu_{\textrm{R}}(t)$ only in the
divisor, which has very little effect on the range, as will be shown later in
fig. 7. Thus, eq. (17) is approximately equivalent to a formula with static or
mean frequency and the dominant error in eq. (17) is the product
$L(t=0)\cdot\epsilon_{\textrm{SCF}}(t)$, which is very similar to the scale
factor error $\langle L\rangle\cdot\epsilon_{\textrm{SCF}}(t)$ discussed here.

Any static error in the scale factor, e.g.
$\langle\epsilon_{\textrm{SCF}}\rangle=10^{-8}$ corresponding to 2.81 MHz
frequency error of the LRI cavity, produces an error proportional to the
separation $L(t)$. If one assumes a drift of 0.01 m/s or 864 m per day and a
1/rev amplitude of 400 m (peak), the scale factor errors correspond to ranging
errors of 8.6 µm per day and 4 µm at 1/rev frequency.

Although the derivation of the scale factor error was based on the LRI, we
anticipate that the KBR instrument features a similar error coupling. The
sinusoidal modulation amplitude of approx. $4\cdot 10^{-12}$ from the proper
time is a fractional quantity and appears with the same amplitude in the
microwave (cf. sec. 2). Thus, the $0.9\,\textrm{\textmu m}_{\textrm{peak}}$
error would be applicable to the KBR as well. However, day-to-day changes in
the USO frequency are usually a factor 10 or more lower than $3.1\cdot
10^{-10}$ (cf. left plot of fig. 1). Thus, the ranging error from linear drift
should be at the level of a few micron level instead of $68\,\textrm{\textmu
m}$ per day. In addition, static scale factor errors are expected to be much
smaller in KBR, as the microwave frequency can be determined with high
accuracy from GPS data in post-processing during orbit determination.

##  4 KBR Dual-One Way Ranging

In the KBR, the microwave radiation with carrier frequencies
$\nu_{A/B}^{\textrm{K/Ka}}$ on satellite A and B is produced by upconversion
of the USO base frequencies $f_{\textrm{A/B,USO}}$ using fixed integer
multipliers (see table 1). Two phase measurements $\varphi^{\textrm{K/Ka}}$
are performed per satellite, one at the K- and another at the Ka-band. These
are combined into

$\displaystyle\varphi_{\textrm{DOWR}}^{\textrm{K/Ka}}(t)$
$\displaystyle:=\varphi_{\textrm{B}}^{\textrm{K/Ka}}(t)-\varphi_{\textrm{A}}^{\textrm{K/Ka}}(t)$
(43)
$\displaystyle=\Phi_{\textrm{B}}^{\textrm{K/Ka}}(t)-\Phi_{\textrm{B}}^{\textrm{K/Ka}}(t-\Delta
t_{\textrm{BA}}^{\textrm{K/Ka}})+\Phi_{\textrm{A}}^{\textrm{K/Ka}}(t)-\Phi_{\textrm{A}}^{\textrm{K/Ka}}(t-\Delta
t_{\textrm{AB}}^{\textrm{K/Ka}})$ (44)
$\displaystyle\approx\nu_{\textrm{B}}^{\textrm{K/Ka}}(t)\cdot\Delta
t_{\textrm{BA}}^{\textrm{K/Ka}}(t)+q_{\textrm{B}}^{\textrm{K/Ka}}+\nu_{\textrm{A}}^{\textrm{K/Ka}}(t)\cdot\Delta
t_{\textrm{AB}}^{\textrm{K/Ka}}(t)+q_{\textrm{A}}^{\textrm{K/Ka}},$ (45)

where $\Delta t_{\textrm{AB}}^{\textrm{K/Ka}}(t)$ describes the light
propagation time from satellite A to B at the respective K- or Ka frequency
band at a reception time $t$. The propagation time can be expressed as

$\displaystyle\Delta t_{\textrm{AB}}^{\textrm{K/Ka}}(t)=\Delta
t_{\textrm{media,AB}}^{\textrm{K/Ka}}(t)+\Delta t_{\textrm{AB}}(t)=\Delta
t_{\textrm{media,AB}}^{\textrm{K/Ka}}(t)+\Delta t_{\textrm{LTC,AB}}(t)+\Delta
t_{\textrm{inst}}(t),$ (46)

where the first term mainly describes the effect of the ionosphere on the
propagation time, which is proportional to the squared frequency ($\Delta
t_{\textrm{media,AB}}^{\textrm{K/Ka}}\propto 1/\nu_{\textrm{A}}^{2}$). The
second term on the right-hand side is the light-time correction Yan et al.
(2021)

$\displaystyle\Delta t_{\textrm{LTC,AB}}(t):=\Delta t_{\textrm{AB}}(t)-\Delta
t_{\textrm{inst}}(t)$ (47)

and the third term is the desired instantaneous range rescaled to light
propagation time

$\displaystyle\Delta
t_{\textrm{inst}}(t):=|\vec{r}_{\textrm{A}}(t)-\vec{r}_{\textrm{B}}(t)|/c_{0}.$
(48)

One can easily show with the previous definitions that the linear combination
Kim (2000)

$\displaystyle\rho_{\textrm{DOWR}}(t):=c_{0}\cdot
a^{\textrm{K}}\cdot\frac{\varphi_{\textrm{DOWR}}^{\textrm{K}}(t)}{\nu_{\textrm{A}}^{\textrm{K}}+\nu_{\textrm{B}}^{K}}+c_{0}\cdot
a^{\textrm{Ka}}\cdot\frac{\varphi_{\textrm{DOWR}}^{\textrm{Ka}}(t)}{\nu_{\textrm{A}}^{\textrm{Ka}}+\nu_{\textrm{B}}^{\textrm{Ka}}}$
(49)

with $a^{\textrm{K}}=-9/7$ and $a^{\textrm{Ka}}=16/7$ is a so-called
ionosphere-free combination where the $\Delta t_{\textrm{media}}$ terms
vanish. This equation is usually employed with constant frequencies
$\nu_{\textrm{A/B}}^{\textrm{K/Ka}}$ in the denominator (cf. eq. (3.1)) and on
daily segments in the regular KBR processing. It is reasonable to assume that
the four phase measurements $\varphi_{\textrm{A/B}}^{\textrm{K/Ka}}(t)$ are
debiased such that the phase is zero at the initial epoch $t=0$ of the time-
series in order to remove the arbitrary constant. This implies that the
$q$-constants in eq. (45) are

$\displaystyle
q_{\textrm{A/B}}^{\textrm{K/Ka}}=-\nu_{\textrm{A/B}}^{\textrm{K/Ka}}(t=0)\cdot\Delta
t_{\textrm{AB/BA}}^{\textrm{K/Ka}}(t=0).$ (50)

Finally, by employing $q$-constants and time-dependent frequencies in eq.
(49), inserting eq. (46) into eq. (49), and considering the relation of
carrier frequencies to USO frequencies from table 1, one obtains

$\displaystyle\frac{a^{\textrm{K}}\cdot\varphi^{\textrm{K}}_{\textrm{DOWR}}(t)}{\nu^{\textrm{K}}_{\textrm{A}}(t)+\nu^{\textrm{K}}_{\textrm{B}}(t)}+\frac{a^{\textrm{Ka}}\cdot\varphi^{\textrm{Ka}}_{\textrm{DOWR}}(t)}{\nu^{\textrm{Ka}}_{\textrm{A}}(t)+\nu^{\textrm{Ka}}_{\textrm{B}}(t)}$
$\displaystyle=\Delta t_{\textrm{inst}}(t)-\Delta
t_{\textrm{DOWR,FV}}(t)-\Delta t_{\textrm{LTC,DOWR}}(t)+\textrm{const},$ (51)

where the first term $\Delta t_{\textrm{inst}}$ on the right-hand side is the
desired ranging signal expressed as propagation time, the second term
describes the frequency variations that couple with the satellite separation
($\Delta t_{\textrm{AB}}\approx\Delta t_{\textrm{BA}}\approx L/c_{0}$) as

$\displaystyle\Delta t_{\textrm{FV,DOWR}}(t):=\frac{\Delta
t_{\textrm{AB}}(0)\cdot f_{\textrm{A,USO}}(0)+\Delta t_{\textrm{BA}}(0)\cdot
f_{\textrm{B,USO}}(0)}{f_{\textrm{A,USO}}(t)+f_{\textrm{B,USO}}(t)}-\frac{\Delta
t_{\textrm{AB}}(0)+\Delta t_{\textrm{BA}}(0)}{2}$ (52)

and the third term is the light-time correction Yan et al. (2021)

$\displaystyle\Delta t_{\textrm{LTC,DOWR}}(t):=\Delta
t_{\textrm{LTC,AB}}(t)\cdot b_{\textrm{AB}}(t)+\Delta
t_{\textrm{LTC,BA}}(t)\cdot b_{\textrm{BA}}(t),$ (53)

where the coefficients $b_{\textrm{AB}}$ and $b_{\textrm{BA}}$ are given in
table 1. The second and third term on the right-hand side of eq. (51) show the
needed correction terms in order to obtain the range independent of USO
frequency variations, i.e.

$\displaystyle\rho_{\textrm{KBR}}$
$\displaystyle:=c_{0}\cdot\left(\frac{a^{\textrm{K}}\cdot\varphi^{\textrm{K}}_{\textrm{DOWR}}}{\nu^{\textrm{K}}_{\textrm{A}}+\nu^{\textrm{K}}_{\textrm{B}}}+\frac{a^{\textrm{Ka}}\cdot\varphi^{\textrm{Ka}}_{\textrm{DOWR}}}{\nu^{\textrm{Ka}}_{\textrm{A}}+\nu^{\textrm{Ka}}_{\textrm{B}}}+\Delta
t_{\textrm{FV,DOWR}}+\Delta t_{\textrm{LTC,DOWR}}+\Delta
t_{\textrm{AOC}}\right),$ (54)

where only $c_{0}$ and $a_{\textrm{K/Ka}}$ are time-independent.

The instantaneous biased range $\rho_{\textrm{KBR}}$ can be directly employed
in gravity field recovery. The so-called antenna offset correction
$\rho_{\textrm{AOC}}$ was added to be consistent with the official SDS data
processing and is discussed in sec 6.2.

The light-propagation times $\Delta t_{\textrm{AB}}$ and $\Delta
t_{\textrm{BA}}$ needed to compute the frequency correction term in eq. (52)
can be obtained from the GPS orbit products and are also used for the light-
time correction computation. The quantities differ by $c_{0}|\Delta
t_{\textrm{AB}}-\Delta t_{\textrm{BA}}|\approx$10\text{\,}\mathrm{m}$$ in
GRACE-FO due to the finite speed of light and the relative velocity between
satellites Yan et al. (2021). Note that the USO frequencies
$f_{\textrm{A/B,USO}}$ in eq. (52) can be replaced with the carrier frequency
of the K- or the Ka-band, since

$\displaystyle\frac{f_{\textrm{A/B,USO}}(0)}{f_{\textrm{A/B,USO}}(t)}=\frac{\nu_{\textrm{A/B}}^{\textrm{K}}(0)}{\nu_{\textrm{A/B}}^{\textrm{K}}(t)}=\frac{\nu_{\textrm{A/B}}^{\textrm{Ka}}(0)}{\nu_{\textrm{A/B}}^{\textrm{Ka}}(t)}.$
(55)

Compared to previous derivations and the official SDS processing described in
Case et al. (2010), the novel aspect of eq. (54) is that it considers a time-
variable frequency when converting the phase to range, and it contains a
correction term $\Delta t_{\textrm{FV,DOWR}}$ for the frequency variations. In
the following sections, we show that this term might be of relevance when LRI
and KBR data is compared at low frequencies.

Table 1: Numerical values for frequencies and coefficients used to describe dual one-way ranging of the KBR. $\dagger$: assumes $\mathrm{d}\tau^{\textrm{USO}}_{\textrm{A}}/\mathrm{d}t=\mathrm{d}\tau^{\textrm{USO}}_{\textrm{B}}/\mathrm{d}t=1$.  Name |  Value/Formula |  Comment
---|---|---
$\hat{f}_{\textrm{A,USO}}$ |  $4.832\,000\text{\,}\mathrm{M}\mathrm{H}\mathrm{z}$ |  exact
$\hat{f}_{\textrm{B,USO}}$ |  $4.832\,099\text{\,}\mathrm{M}\mathrm{H}\mathrm{z}$ |  exact
${f}_{\textrm{A,USO}}(t)$ |  $$4.832\,000\text{\,}\mathrm{M}\mathrm{H}\mathrm{z}$\cdot\mathrm{d}\tau^{\textrm{USO}}_{\textrm{A}}/\mathrm{d}t$ |
${f}_{\textrm{B,USO}}(t)$ |  $$4.832\,099\text{\,}\mathrm{M}\mathrm{H}\mathrm{z}$\cdot\mathrm{d}\tau^{\textrm{USO}}_{\textrm{B}}/\mathrm{d}t$ |
$\nu_{\textrm{A}}^{\textrm{K}}(t)$ |  $f_{\textrm{A,USO}}(t)\cdot 5076$ |  $$24.527\,232\text{\,}\mathrm{G}\mathrm{H}\mathrm{z}$^{\dagger}$
$\nu_{\textrm{A}}^{\textrm{Ka}}(t)$ |  $f_{\textrm{A,USO}}(t)\cdot 6768$ |  $$32.702\,976\text{\,}\mathrm{G}\mathrm{H}\mathrm{z}$^{\dagger}$
$\nu_{\textrm{B}}^{\textrm{K}}(t)$ |  $f_{\textrm{B,USO}}(t)\cdot 5076$ |  $$24.527\,734\,524\text{\,}\mathrm{G}\mathrm{H}\mathrm{z}$^{\dagger}$
$\nu_{\textrm{B}}^{\textrm{Ka}}(t)$ |  $f_{\textrm{B,USO}}(t)\cdot 6768$ |  $$32.703\,646\,032\text{\,}\mathrm{G}\mathrm{H}\mathrm{z}$^{\dagger}$
$a^{\textrm{K}}$ |  $\frac{-\nu_{\textrm{A}}^{\textrm{K}}\cdot\nu_{\textrm{B}}^{\textrm{K}}}{(\nu_{\textrm{A}}^{\textrm{Ka}}\cdot\nu_{\textrm{B}}^{\textrm{Ka}}-\nu_{\textrm{A}}^{\textrm{K}}\cdot\nu_{\textrm{B}}^{\textrm{K}})}=-9/7$ |  exact
$a^{\textrm{Ka}}$ |  $\frac{\nu_{\textrm{A}}^{\textrm{Ka}}\cdot\nu_{\textrm{B}}^{\textrm{Ka}}}{(\nu_{\textrm{A}}^{\textrm{Ka}}\cdot\nu_{\textrm{B}}^{\textrm{Ka}}-\nu_{\textrm{A}}^{\textrm{K}}\cdot\nu_{\textrm{B}}^{\textrm{K}})}=16/7$ |  exact
$b_{\textrm{AB}}^{\textrm{K}}(t)$ |  $\frac{(\nu_{\textrm{A}}^{\textrm{K}})^{2}\cdot\nu_{\textrm{B}}^{\textrm{K}}}{(\nu_{\textrm{A}}^{\textrm{K}}+\nu_{\textrm{B}}^{\textrm{K}})(\nu_{\textrm{A}}^{\textrm{K}}\nu_{\textrm{B}}^{\textrm{K}}-\nu_{\textrm{A}}^{\textrm{Ka}}\nu_{\textrm{B}}^{\textrm{Ka}})}=\frac{-9{f}_{\textrm{A,USO}}(t)}{7\cdot({f}_{\textrm{A,USO}}(t)+{f}_{\textrm{B,USO}}(t))}$ |  $\frac{-43488000}{67648693}\approx-0.642851^{\dagger}$
$b_{\textrm{AB}}^{\textrm{Ka}}(t)$ |  $-\frac{(\nu_{\textrm{A}}^{\textrm{Ka}})^{2}\cdot\nu_{\textrm{B}}^{\textrm{Ka}}}{(\nu_{\textrm{A}}^{\textrm{Ka}}+\nu_{\textrm{B}}^{\textrm{Ka}})(\nu_{\textrm{A}}^{\textrm{K}}\nu_{\textrm{B}}^{\textrm{K}}-\nu_{\textrm{A}}^{\textrm{Ka}}\nu_{\textrm{B}}^{\textrm{Ka}})}=\frac{16{f}_{\textrm{A,USO}}(t)}{7\cdot({f}_{\textrm{A,USO}}(t)+{f}_{\textrm{B,USO}}(t))}$ |  $\frac{77312000}{67648693}\approx 1.1428454^{\dagger}$
$b_{\textrm{BA}}^{\textrm{K}}(t)$ |  $\frac{\nu_{\textrm{A}}^{\textrm{K}}\cdot(\nu_{\textrm{B}}^{\textrm{K}})^{2}}{(\nu_{\textrm{A}}^{\textrm{K}}+\nu_{\textrm{B}}^{\textrm{K}})(\nu_{\textrm{A}}^{\textrm{K}}\nu_{\textrm{B}}^{\textrm{K}}-\nu_{\textrm{A}}^{\textrm{Ka}}\nu_{\textrm{B}}^{\textrm{Ka}})}=\frac{-9{f}_{\textrm{B,USO}}(t)}{7\cdot({f}_{\textrm{A,USO}}(t)+{f}_{\textrm{B,USO}}(t))}$ |  $\frac{-43488891}{67648693}\approx-0.642864^{\dagger}$
$b_{\textrm{BA}}^{\textrm{Ka}}(t)$ |  $-\frac{\nu_{\textrm{A}}^{\textrm{Ka}}\cdot(\nu_{\textrm{B}}^{\textrm{Ka}})^{2}}{(\nu_{\textrm{A}}^{\textrm{Ka}}+\nu_{\textrm{B}}^{\textrm{Ka}})(\nu_{\textrm{A}}^{\textrm{K}}\nu_{\textrm{B}}^{\textrm{K}}-\nu_{\textrm{A}}^{\textrm{Ka}}\nu_{\textrm{B}}^{\textrm{Ka}})}=\frac{16{f}_{\textrm{B,USO}}(t)}{7\cdot({f}_{\textrm{A,USO}}(t)+{f}_{\textrm{B,USO}}(t))}$ |  $\frac{77313584}{67648693}\approx 1.142869^{\dagger}$
$b_{\textrm{AB}}(t)$ |  $b_{\textrm{AeBr}}^{\textrm{K}}+b_{\textrm{AeBr}}^{\textrm{Ka}}=\frac{{f}_{\textrm{A,USO}}(t)}{{f}_{\textrm{A,USO}}(t)+{f}_{\textrm{B,USO}}(t)}$ |  $\approx 0.499995^{\dagger}$
$b_{\textrm{BA}}(t)$ |  $b_{\textrm{BeAr}}^{\textrm{K}}+b_{\textrm{BeAr}}^{\textrm{Ka}}=\frac{{f}_{\textrm{B,USO}}(t)}{{f}_{\textrm{A,USO}}(t)+{f}_{\textrm{B,USO}}(t)}$ |  $\approx 0.500005^{\dagger}$

##  5 Range Differences between LRI and KBR

The LRI and KBR ranging data can be directly compared at the instantaneous
range level, meaning that light-time corrections of KBR and LRI and KBR
antenna offset correction are applied. The light-time correction is different
for KBR DOWR and LRI TWR. The ionospheric effect is already removed by default
in KBR1B datasets, but the correction is provided separately in case one wants
to derive the variations in electron content of the ionosphere. Since KRB1B
data is provided at a rate of $0.2\text{\,}\mathrm{H}\mathrm{z}$ and LRI1B
data at a rate of $0.5\text{\,}\mathrm{H}\mathrm{z}$, the data has to be
resampled, ideally through low-pass filtering and decimation in order to avoid
aliasing.

The LRI data as provided in the official v04 LRI1B dataset has already been
deglitched, which is a process that removes steps and glitches in the ranging
data that are correlated with some thruster activations Abich et al. (2019).
Moreover, the LRI ranging data is rescaled and time-shifted such that KBR-LRI
residuals are reduced. This rescaling is necessary due to the LRI absolute
laser frequency $\nu_{\textrm{R}}$, defining the conversion factor from phase
to range, only being known with low accuracy in GRACE-FO. In contrast to the
KBR frequencies $\nu_{\textrm{A/B}}^{\textrm{K/Ka}}$, the laser frequency
$\nu_{\textrm{R}}$ cannot be directly measured in-flight on GRACE-FO, although
attempts to derive the frequency from laser telemetry are being made Misfeldt
et al. (2021). However, future missions will likely host a LRI scale factor
unit that enable measurements of a time-resolved frequency Rees et al. (2021).

The deglitching and rescaling of LRI data was sub-optimal in v04 LRI1B data
before 2020-06-30, but improved significantly thereafter due to changes in
processing Hui Ying Wen (2019).

We previously derived an alternative LRI1B product from level1a data, which we
label as v50 in order to distinguish it from the official v04 product Müller
et al. (2020). A comparison of both datasets is shown for January 2019 in
figure 3. Other time-periods exhibit the same characteristics. The range
signal in GRACE-FO has dominant 1/rev (and 2/rev) oscillations of several
$100\text{\,}\mathrm{m}$ amplitude and a continuous spectrum containing the
interesting gravity field information and non-gravitational accelerations up
to approx. $35\text{\,}\mathrm{m}\mathrm{H}\mathrm{z}$. At high frequencies,
the LRI and KBR range data is dominated by the respective instrument noise.

The traces showing the difference between KBR and LRI reveal that these
residuals are dominated by KBR noise at frequencies above
$35\text{\,}\mathrm{m}\mathrm{H}\mathrm{z}$ and the residual level increases
towards lower frequencies. The residuals show pronounced 2/rev oscillations in
LRI1Bv04 and LRI1Bv50 with a few micron amplitude. However, LRI1Bv04 has an
additional 1/rev oscillation. When we re-estimate a scale and time-shift of
LRI1Bv04, this 1/rev peak vanishes in the residuals and the results get closer
to LRI1Bv50. This is an indication of the scale and time-shift being slightly
inaccurate in the official v04 data. Therefore, we will often show the
original LRI1Bv04 data as well as our rescaled and time-shifted version
without the 1/rev oscillations in the following.

Figure 3: Typical LRI and KBR ranging data shown in the spectral domain (top)
and in the time-domain (bottom). The left plots use the official LRI1Bv04 data
product, while the right plots show an in-house derived alternative product
(LRI1Bv50).

The day 2019-01-16 is visible as an outlier in the residuals (upper left plot
of fig. 3), because the deglitching did not properly remove a glitch in the
LRI1Bv04 data, which then disturbed the time-shift estimation.

##  6 Instrument Noises and Corrections

In order to understand which noises or error sources are limiting the KBR-LRI
residuals, we have plotted many of the known effects in KBR and LRI in fig. 4.

We have analyzed the light-time correction extensively in Yan et al. (2021)
and skip it therefore in our analysis given here, as we found no indication
that the accuracy of this correction could be limiting KBR-LRI residuals in
our previous work.

Figure 4: Noise contributors and corrections in LRI and KBR in terms of
amplitude spectral densities (ASDs). The grey KBR-LRI residuals are described
in sec. 5 and replotted from fig. 3. The KBR tilt-to-length correction
(Antenna-Offset-Correction) in orange is provided in the KBR1B dataset and
already removed in the KBR range, while the corresponding effect in the LRI
has a lower magnitude and was derived based on Wegener et al. (2020) \- see
sec. 6.2 for more details. The ionospheric correction (red) is relevant for
KBR, but negligible in LRI due to the scaling with wavelength. The smaller
wavelength of the laser scales the phase readout noise (magenta) to a
negligible level as well, while in KBR it is dominant at high frequencies (cf.
sec. 6.3). Raw frequency variations (dashed blue) assume a static frequency
(within one day) and are related to the stability of the optical cavity and
proper time effect in the LRI, whereas they are related to the USO stability
in KBR. These raw frequency variations can be corrected by using a time-
dependent frequency (cf. sec. 6.5). The effect of timetag errors is being
corrected in LRI and KBR by employing CLK1B $\epsilon_{\textrm{time}}$
information (dashed vs. solid black traces). For some legend entries, multiple
days were plotted in order to illustrate the stationarity and variability of
the traces. The mean satellite separation is
$L=$194\text{\,}\mathrm{k}\mathrm{m}$$ for January 2019. The ENBW for both
plots is 52 µHz.

###  6.1 Ionospheric correction

The propagation of electromagnetic waves like microwaves or laser light is
affected by the non-ideal vacuum, in particular by charged particles forming
the ionosphere. The KBR system uses the dual-band observations to measure and
correct the effect. The magnitude of this correction, as provided in the
KBR1Bv04 data, is shown as red trace in the right plot of fig. 4. The
correction is well above the KBR-LRI residuals and the accuracy of the
correction could in principle cause some error. However, we assessed the
second-order ionospheric effect Hernández-Pajares et al. (2007) and concluded
that the ionosphere is likely not limiting KBR-LRI residuals, since the
second-order effect is direction-dependent and highly rejected in the DOWR
combination, in which the contribution of the path from satellite A to B is
counteracted by the contribution of the path from B to A.

The ionospheric effect as measured by the KBR can be rescaled to the LRI
optical wavelength of $1064.5\text{\,}\mathrm{n}\mathrm{m}$, see red trace on
left plot of fig. 4. For such near-infrared wavelengths, the ionospheric delay
is a negligible error contributor Sheard et al. (2012).

###  6.2 Tilt-To-Length Coupling

Satellite rotations cause pathlength changes when the reference points (RPs)
of the LRI or the antenna phase centers (APCs) of KBR are not co-located with
the respective center-of-mass of satellites, which are the pivot points for
rotations in space (Müller, 2017, sec. 2.6.3). The RPs and APCs are
effectively the fiducial points used to determine the biased inter-satellite
range from interferometric phase readout and are fixed in the satellite body
frame, i.e. moving in inertial space when the satellites rotate.

The main effect of the coupling can be expressed in terms of a vector
describing the offset between center-of-mass and the RP or APC that is
projected onto the line-of-sight as the measurement axis. To first order,
lateral components perpendicular to the line-of-sight ($y$ & $z$) couple
linearly with yaw and pitch, i.e. 100 µm offset produce a yaw and pitch
coupling of 100 µm/rad per satellite into the biased range. The $x$-offset
couples quadratically with yaw and pitch as 100 µm/$\textrm{rad}^{2}$.

In the KBR context, the tilt-to-length coupling is called antenna offset
correction Wen et al. (2019); Case et al. (2010) and it is provided in the
KRB1B data set in terms of the DOWR-range correction, including effects from
both spacecraft (cf. $c_{0}\Delta t_{\textrm{AOC}}$ in eq. (54)). Since the
KBR antenna is at the front panel of the satellite, the APCs are offset by
approx. $1.47\text{\,}\mathrm{m}$ in $x$ direction from the center-of-mass.
This means that the correction contains a significant quadratic contribution
from the pointing angles. The coupling factors can be measured in-orbit using
so-called KBR calibration maneuver during which the satellites are
deliberately rotated.

In the LRI, the RPs on both satellites were co-located to the center-of-mass
to the level of the mechanical integration accuracy in the order of 100 µm in
$y$ and $z$, which yields the dominant linear contribution. The LRI tilt-to-
length coupling factors can be measured using the center-of-mass calibration
maneuvers Wegener et al. (2020) that are usually repeated every six months. As
of version 04 of LRI1B data, the tilt-to-length correction is not derived or
applied to LRI data, but this may change in future releases. More information
in regard to LRI tilt-to-length coupling can be found in Wegener (2022).

The magnitude of the tilt-to-length effect on the KBR and LRI range
measurement is shown in orange in fig. 4. In the LRI case, the effect is
approximately one order of magnitude smaller compared to the KBR. Based on
these results, we can conclude that the tilt-to-length effect does not limit
the KBR-LRI residuals.

###  6.3 Phase Readout Noise

Figure 5: CNR-limited phase readout noise in LRI (left) and KBR (right). In
the LRI, the phase readout noise is well below the actual ranging signal,
which is limited by laser frequency noise (magenta, left plot). In KBR, the
phase readout noise is likely limiting the KBR range at frequencies above
$10\text{\,}\mathrm{m}\mathrm{H}\mathrm{z}$. The measured noise level in the
K-band is slightly higher than the expected one derived from the CNR-values.
This also yields a small discrepancy in the final ionosphere-free DOWR
combination (magenta, dashed $0.4\,\textrm{\textmu m}/\sqrt{\textrm{Hz}}$ vs.
solid $\approx 0.6\,\textrm{\textmu m}/\sqrt{\textrm{Hz}}$).

The KBR and LRI both track the phase of an interference beatnote that is
present as an electrical signal in the corresponding processing unit. Phase
changes are proportional to changes in the light or microwave propagation time
up to first order and, thus, to distance changes between spacecraft. The
readout noise of the phase is usually determined by the carrier-to-noise ratio
(CNR) entering the phase-tracking loop. which can be expressed in terms of a
ranging error as Müller (2017)

$\displaystyle\textrm{ASD}[\rho_{\textrm{CNR}}](f)=\frac{\lambda\cdot$1\text{\,}\mathrm{r}\mathrm{a}\mathrm{d}$}{2\cdot
2\pi\cdot\sqrt{\textrm{CNR}}},\quad\textrm{CNR}=\frac{I_{\textrm{rms}}^{2}}{\textrm{PSD}[\mathcal{N}](f_{\textrm{beat}})},\quad[\textrm{ASD}[\rho_{\textrm{CNR}}]]=$\text{\,}\mathrm{m}\mathrm{/}\sqrt{Hz}$,$
(56)

where the half-wavelength $\lambda/2$ is used to rescale the white phase noise
level to units of half-round-trip range.

The CNR is defined as the ratio of root-mean-square (rms) carrier power
$I_{\textrm{rms}}^{2}$, i.e. squared rms amplitude of the beatnote
photocurrent, to the noise power spectral density
($\textrm{PSD}[\mathcal{N}]$) of the photocurrent at the beatnote frequency.
Typical values for the LRI are above
$$80\text{\,}\mathrm{d}\mathrm{B}\mathrm{H}\mathrm{z}$=10^{8}\,$\text{\,}\mathrm{H}\mathrm{z}$$
for both the reference and transponder satellite, which is well above the
minimum requirement of
$$70\text{\,}\mathrm{d}\mathrm{B}\mathrm{H}\mathrm{z}$=10^{7}\,$\text{\,}\mathrm{H}\mathrm{z}$$
of the instrument. The transponder phase measurement is used as a sensor in a
feedback control loop to control the transponder laser frequency with high
gain and bandwidth, e.g. the phase variations in terms of range on the
transponder side are below
$10\text{\,}\mathrm{p}\mathrm{m}\mathrm{/}\sqrt{\textrm{Hz}}$ at
$1\text{\,}\mathrm{H}\mathrm{z}$ and even smaller for lower frequencies. Thus,
it is an in-loop measurement close to zero. Any phase variations arising on
the transponder, e.g. from laser noise or readout noise, are imprinted onto
the laser beam and transmitted to the reference satellite. Hence, the
reference satellite measures the phase readout noise of transponder and
reference satellite

$\displaystyle\textrm{ASD}[\rho_{\textrm{LRI,CNR}}](f)=\frac{\lambda_{\textrm{R}}\cdot$1\text{\,}\mathrm{r}\mathrm{a}\mathrm{d}$}{2\cdot
2\pi}\cdot\sqrt{\frac{1}{\textrm{CNR}_{\textrm{T}}}+\frac{1}{\textrm{CNR}_{\textrm{M}}}}.$
(57)

This formula is actually independent of the loop gain of the frequency-
controller, i.e. even with low gain and significant signal in the transponder
phase, the ranging signal $\rho_{\textrm{LRI}}$ is usually formed by the
combination of both satellites in order to remove the common phase ramp (cf.
eq. (9)) and to recover the complete ranging signal. Thus, it will exhibit the
phase readout noise of both satellites.

The CNR-limited phase readout noise $\textrm{ASD}[\rho_{\textrm{CNR}}]$ is,
strictly speaking, non-stationary, as the CNR is time-dependent. Therefore, we
plot short segments of data in order to prevent potential artifacts due to
non-stationarity. The left plot of fig. 5 shows the calculated phase readout
noise compared to the ranging signal of the LRI (magenta trace). A CNR value
of
$$77.5\text{\,}\mathrm{d}\mathrm{B}\mathrm{H}\mathrm{z}$=10^{7.75}\,$\text{\,}\mathrm{H}\mathrm{z}$$
from a single satellite yields approx.
$11\text{\,}\mathrm{p}\mathrm{m}\mathrm{/}\sqrt{Hz}$ of white phase readout
noise, or $\sqrt{2}\cdot$11\text{\,}\mathrm{p}\mathrm{m}\mathrm{/}\sqrt{Hz}$$
from two satellites. Such low values are not observable in the ranging phase
of the LRI (mangenta trace in left plot of fig. 5) due to other noise sources,
mainly laser frequency variations from the cavity.

However, the LRI system measures the interference beatnote on a satellite with
four channels in order to make use of the Differential Wavefront Sensing (DWS)
technique Sheard et al. (2012). The ranging information is given by the
average of the four phase measurements, which is the linear combination of
four channels with multipliers $+1/4$. The CNR values reported by the LRP of
the LRI refer to the coherent sum (average) of four channels. However, if the
noise among different channels is uncorrelated and the phase difference among
different channels is small (cf. (Müller, 2017, sec. 2.6.9) for details), the
same CNR values would apply to other linear combinations as well. These other
three linear combinations with multipliers $\pm 1/4$ can be called DWSyaw,
DWSpitch and DWSx. DWSyaw and DWSpitch represent the angular misalignment
between the interfering laser beams, and are actually zeroed by a feedback
control loop using a fine-steering mirror as actuator. The remaining DWSx
combination is not zeroed and exhibits white noise behavior for frequencies
above $0.1\text{\,}\mathrm{H}\mathrm{z}$ (cf. left plot in fig. 5) with a
magnitude consistent with the expected CNR-limited phase readout noise. It is
noteworthy that DWSx is a local quantity that needs to be compared against the
single spacecraft CNR-limit (eq. (56)). If noise is transferred between
satellites via the laser light, it will affect all channels on the receiver,
thus appearing in the average phase of the four channels, but not in the DWS
combinations of the receiver.

The consistency between DWSx noise level and analytical CNR-limited phase
readout noise level is at least an indicator that the CNR-limited phase
readout noise is at the expected level in the ranging signal (magenta trace),
though it is not directly measurable.

The phase readout noise in the context of the KBR is usually called system
noise Thomas (1999); Kim (2000). Based on eq. (56) and the DOWR combination
discussed in sec. 4, it is straightforward to express the CNR-limited phase
readout noise in the ionosphere-free DOWR range as

$\displaystyle\textrm{ASD}[\rho_{\textrm{KBR,CNR}}](f)=\frac{$1\text{\,}\mathrm{r}\mathrm{a}\mathrm{d}$}{2\cdot
2\pi}\cdot\sqrt{\frac{(a^{\textrm{K}}\lambda_{\textrm{A}}^{\textrm{K}})^{2}}{\textrm{CNR}_{\textrm{A}}^{\textrm{K}}}+\frac{(a^{\textrm{K}}\lambda_{\textrm{B}}^{\textrm{K}})^{2}}{\textrm{CNR}_{\textrm{B}}^{\textrm{K}}}+\frac{(a^{\textrm{Ka}}\lambda_{\textrm{A}}^{\textrm{Ka}})^{2}}{\textrm{CNR}_{\textrm{A}}^{\textrm{Ka}}}+\frac{(a^{\textrm{Ka}}\lambda_{\textrm{B}}^{\textrm{Ka}})^{2}}{\textrm{CNR}_{\textrm{B}}^{\textrm{Ka}}}},$
(58)

where coefficients $a^{\textrm{K}}$ and $a^{\textrm{Ka}}$ as well as
wavelength $\lambda=c_{0}/\nu$ can be obtained from table 1. This expression
is consistent with the derivation given in (Thomas, 1999, sec. B-1) and shown
as magenta-dashed trace in the right plot of fig. 5. Unfortunately, the actual
DOWR ranging noise is slightly higher than predicted based on the CNR values
reported in the KBR1B data product. When plotting the individual range
measurements at K- and Ka-band, we can identify that the discrepancy is caused
by the K-band range, which exhibits slightly higher noise than predicted by
the CNR. Though the reason of this discrepancy is unclear to the authors of
this paper, the deviation is rather small. Based on fig. 4, we conclude that
the phase readout noise is negligible for the LRI, while it likely limits the
KBR and therefore the KBR-LRI residuals at high frequencies to a white noise
level of approx. 0.6 µm/$\sqrt{\textrm{Hz}}$.

###  6.4 Timetag Errors

The telemetry of KBR and LRI contains recorded instrument timetags, which are
equally sampled at their respective nominal rate. Timetag errors arise when
the timetags in the instrument timeframe are converted to another time system,
usually to GPS time.

The KBR and LRI are driven by the same ultra-stable oscillator (USO). Thus,
the KBR time, which is often called IPU receiver time, and the LRI time are
both realizations of the USO time, but KBR/IPU and LRI time have a quasi-
static offset relative to each other. It is quasi-static, because it changes
whenever one of the instruments reboots. And the offset is usually less than a
second, since the instruments synchronize their respective time to the integer
GPS second after the reboot occurred.

Since KBR and GPS measurements are performed within the IPU, the KBR timetags
can be directly converted to GPS time using the $\epsilon_{\textrm{time}}$
value from the CLK1B product (cf. eq. (2)) that is determined during precise
orbit determination for each satellite.

For the LRI, the conversion requires addtional steps. First, datation reports
as given in the LHK1A/B data product provide the timetag difference between
LRI time and OBC time. Second, the TIM1B data product contains the offset
between OBC time and KBR/IPU time, which is virtually zero during nominal
operation, as the OBC time is steered towards the IPU time. In the end, the
$\epsilon_{\textrm{time}}$ from the CLK1B data product yields the final
relation needed to convert the timetag to GPS time.

In fig. 4, the curves labeled _timetag raw_ shown in dashed black indicate a
hypothetical noise calculated under the assumption of using a daily mean value
of $\epsilon_{\textrm{time}}$, e.g. no $\epsilon_{\textrm{time}}$ information
is employed for frequencies above $\approx 1/$86\,400\text{\,}\mathrm{s}$$.
Thereby, the error in the timetag conversion from instrument to GPS time would
be given by the variability of $\epsilon_{\textrm{time}}$. The in-flight
variability was plotted as spectral density in the right plot of fig. 1 for
GF-1 (blue) and GF-2 (cyan blue). The curves increase below
$3\text{\,}\mathrm{m}\mathrm{H}\mathrm{z}$, mainly due to the USO instability
and proper time modulations at 1/rev and 2/rev frequencies, while the increase
above $3\text{\,}\mathrm{m}\mathrm{H}\mathrm{z}$ is likely caused by
limitations of precise orbit determination and clock solution precision, i.e.
the measurement precision of $\epsilon_{\textrm{time}}$.

The analytical formula for the USO stability from (Thomas, 1999, B-1), shown
in solid magenta in fig. 1, is in poor agreement with the GRACE-FO in-flight
data, justifying a refined USO model that we derive ad-hoc as

$\displaystyle\textrm{ASD}[\delta
t_{\textrm{USO,mod}}](f)=$1\text{\,}\mathrm{s}$\cdot\sqrt{2.556\cdot
10^{-33}\,$\text{\,}\mathrm{H}\mathrm{z}^{4}$/f^{5}+3.325\cdot
10^{-26}\,$\text{\,}\mathrm{H}\mathrm{z}$/f^{2}}.$ (59)

The first time-derivative of this USO model is shown by the magenta dashed
trace in the right plot of fig. 1. It is only a rough model which only
approximately describes the variability in $\epsilon_{\textrm{time}}$ in the
spectral domain, because the spectra differ slightly between satellites and
the $\epsilon_{\textrm{time}}$ spectra change slightly from day to day (not
shown). Furthermore, we take an educated guess for the measurement precision
for eps-time due to POD/Clock errors of the GPS measurement system. The model
reads

$\displaystyle\textrm{ASD}[\delta
t_{\epsilon}](f)=\frac{$1.5\text{\,}\mathrm{m}\mathrm{m}$}{c_{0}\cdot\sqrt{f}}+\frac{10\,\textrm{\textmu
m}\cdot\textrm{Hz}}{c_{0}\cdot\sqrt{f^{2.5}}},$ (60)

and is shown as a black trace in the spectral domain in fig. 1 and as an Allan
deviation in fig. 6. The model has a rms of
$$11\text{\,}\mathrm{m}\mathrm{m}$/c_{0}$ in the band
$$0.1\text{\,}\mathrm{m}\mathrm{H}\mathrm{z}$-$0.05\text{\,}\mathrm{H}\mathrm{z}$$
and (roughly) agrees with the GPS POD/clock measurement noise in the Deep
Space Atomic Clock (DSAC) experiment conducted in a low-Earth orbit Burt et
al. (2021) with respect to Allan deviations. The noise level of DSAC is
consistent with the GRACE-FO USO Allan deviations in fig. 6 at lowest
integration times (i.e. highest frequencies), though it bears mentioning that
the GRACE-FO traces are for a single day and exhibit some day-to-day
variability (not shown).

Figure 6: Allan deviations of the most of the traces shown in the right plot
of fig. 1. The additional green trace shows the in-orbit stability of the
trapped-ion atomic clock Burt et al. (2021), which is dominated by GPS
measurement noise for the shown values of $\tau$.

The noise models for the timetag error in eq. (59) and (60) can produce an
error in the phase measurement of KBR and LRI. In order to assess its
magnitude, we replace the argument of the phase measurement from correct time
$t$ with an erroneous time $t+\delta t$:

$\displaystyle\varphi(t)\longrightarrow\varphi(t+\delta
t)\approx\varphi(t)+\dot{\varphi}\cdot\delta t.$ (61)

The effect of the timetag error on the LRI ranging phase can be derived from
eq. (9), whereas its effect on the LRI range can be determined from eq. (17)
as

$\displaystyle\textrm{ASD}[\delta\varphi_{\textrm{TWR}}]$
$\displaystyle=\textrm{ASD}[\dot{\varphi}_{\textrm{T}}\cdot\delta
t_{\textrm{T}}-\dot{\varphi}_{\textrm{R}}\cdot\delta
t_{\textrm{R}}]\approx\dot{\varphi}\cdot\sqrt{2}\cdot\textrm{ASD}[\delta
t]\approx$10\text{\,}\mathrm{M}\mathrm{H}\mathrm{z}$\cdot\sqrt{2}\cdot\textrm{ASD}[\delta
t]$ (62) $\displaystyle\textrm{ASD}[\delta\rho_{\textrm{TWR}}]$
$\displaystyle\approx\textrm{ASD}[\varphi_{\textrm{TWR}}\cdot\lambda_{\textrm{R}}/2]\approx$5\text{\,}\mathrm{m}\mathrm{/}\mathrm{s}$\cdot\sqrt{2}\cdot\textrm{ASD}[\delta
t],$ (63)

using the fact that the LRI phases measured on transponder and reference
satellites have a common and dominant slope of
$10\text{\,}\mathrm{M}\mathrm{H}\mathrm{z}$ (cf. eq. (8)), and we assumed that
the timetag noise on both spacecraft is uncorrelated and has the same
stochastic characteristics, such that the differential combination yields
$\sqrt{2}\delta t$ in the spectral domain.

The same approach for the KBR at K- and Ka-band yields a timetag induced noise
of

$\displaystyle\textrm{ASD}[\delta\varphi_{\textrm{DOWR}}^{\textrm{K/Ka}}]$
$\displaystyle=\textrm{ASD}[\dot{\varphi}_{\textrm{B}}^{\textrm{K/Ka}}\cdot\delta
t_{\textrm{B}}-\dot{\varphi}_{\textrm{A}}^{\textrm{K/Ka}}\cdot\delta
t_{\textrm{A}}]$ (64)
$\displaystyle\approx\dot{\varphi}\cdot\sqrt{2}\cdot\textrm{ASD}[\delta
t]\approx\sqrt{2}\cdot\textrm{ASD}[\delta
t]\cdot\begin{cases}$500\text{\,}\mathrm{k}\mathrm{H}\mathrm{z}$,&\text{\quad
for K-band}\\\ $670\text{\,}\mathrm{k}\mathrm{H}\mathrm{z}$,&\text{\quad for
Ka-band}\end{cases}$ (65)
$\displaystyle\textrm{ASD}[\delta\rho_{\textrm{DOWR}}^{\textrm{K/Ka}}]$
$\displaystyle\approx\textrm{ASD}[\varphi_{\textrm{DOWR}}^{\textrm{K/Ka}}\cdot
c_{0}/(\nu_{\textrm{A}}^{\textrm{K/Ka}}+\nu_{\textrm{B}}^{\textrm{K/Ka}})]\approx$3071\text{\,}\mathrm{m}\mathrm{/}\mathrm{s}$\cdot\sqrt{2}\cdot\textrm{ASD}[\delta
t].$ (66)

for the phase $\varphi$ and range $\rho$. Though the timetag error is
different for the K- and Ka-phase due to the difference in beatnote
frequencies ($500\text{\,}\mathrm{k}\mathrm{H}\mathrm{z}$ vs.
$670\text{\,}\mathrm{k}\mathrm{H}\mathrm{z}$), it is the same in terms of
range. When forming the final ionosphere-free DOWR combination, one has to
consider that the timetag error is correlated between the K- and Ka-band, such
that both contributors add linearly as

$\displaystyle\textrm{ASD}[\delta\rho_{\textrm{DOWR}}]=a^{\textrm{K}}\cdot\textrm{ASD}[\delta\rho_{\textrm{DOWR}}^{\textrm{K}}]+a^{\textrm{Ka}}\cdot\textrm{ASD}[\delta\rho_{\textrm{DOWR}}^{\textrm{Ka}}]\approx$3071\text{\,}\mathrm{m}\mathrm{/}\mathrm{s}$\cdot\sqrt{2}\cdot\textrm{ASD}[\delta
t],$ (67)

with the coefficients $a^{\textrm{K/Ka}}$ given in table 1.

Eq. (63) and (67) depend on $\delta t$ and have been used to produce the black
timetag curves of LRI and KBR, respectively, in fig. 4. Traces labeled with
_timetag raw_ are based on $\delta t=\delta t_{\textrm{USO,mod}}$ given in eq.
(59), while _timetag corrected_ use the $\delta t=\delta t_{\epsilon}$ noise
from eq. (60). As apparent from fig. 4, the timetag noise in KBR likely is a
significant contributor to the increase of KBR-LRI residuals below
$10\text{\,}\mathrm{m}\mathrm{H}\mathrm{z}$.

However, we note that the timetag error depends on the differential clock
jitter between both satellites and our assumption of uncorrelated noise
between satellites ($\textrm{ASD}[\delta t_{\epsilon}]\cdot\sqrt{2}$) neglects
potential common-mode rejection of errors that might appear between close-by
GPS receiver. Moreover, for the sake of simplicity, we considered just a
stochastic USO noise model and omitted potential sinusoidal variations due to
the proper time and secular drifts due to uncertainties in the USO frequency
for the calculation of the trace labeled _timetag raw_. The proper time effect
is not severe, because the differential proper time
($\tau_{\textrm{A}}-\tau_{\textrm{B}}$) between satellites is much smaller
than the individual contribution ($\tau_{\textrm{A}}$ or $\tau_{\textrm{B}}$).
We regard these approximations as acceptable, as the _timetag raw_ trace is
anyway rather hypothetical due to the general availability of GPS data to
correct the raw values.

###  6.5 Carrier Frequency Variations

Carrier frequencies are the laser frequency $\nu_{\textrm{R}}$ in LRI and the
K- and Ka-band frequencies $\nu_{\textrm{A/B}}^{\textrm{K/Ka}}$ in KBR.
Variations in these frequencies correspond to changes in the conversion factor
from phase to range, as discussed in sec. 3. If the frequency variations are
measured and known, they can be accounted for in the conversion and do not
falsify the range measurement, when employing eq. (28). For the application in
GRACE-FO, the approximation formulas eq. (22), (23) and (31) usually are
sufficiently accurate.

For the LRI, we assume the following model for the laser frequency at the
reference satellite and for a time-series of one day

$\displaystyle\nu_{\textrm{R}}(t)\approx\langle\nu_{\textrm{R}}(t)\rangle+\langle\nu_{\textrm{R}}(t)\rangle\cdot\left(\frac{\mathrm{d}\tau_{\textrm{R}}(t)}{\mathrm{d}t}-1\right)+\delta\nu_{\textrm{cav}}(t),$
(68)

where $\langle\nu_{\textrm{R}}(t)\rangle$ is the daily mean frequency
determined from correlating LRI phase to KBR range (cf. sec. 5), the second
term contains the modulation of the laser frequency due to the proper time
$\tau_{\textrm{R}}$ (cf. eq. (6)) and the last term accounts for the
instability of the cavity, which was modelled in sec. 2 as
$\delta\nu_{\textrm{cav}}(t)=\mathrm{d}\tau_{\textrm{FV}}(t)/\mathrm{d}t\cdot\langle\nu_{\textrm{R}}(t)\rangle$.
The daily mean frequency $\langle\nu_{\textrm{R}}(t)\rangle$ is used as the
conversion factor from phase to range (cf. eq. (17)) in the LRI1Bv04
processing. The error arising from neglecting the intraday time-variability of
the frequency represents a scale factor error as discussed in sec. 3.4. Such
errors couple into the range to first order as the product of fractional
frequency variation and inter-satellite distance $L(t)$. The ranging error due
to $\delta\nu_{\textrm{cav}}$ is usually called laser frequency noise in the
LRI and caused by the limited stability of the cavity resonance frequency,
e.g. due to thermal fluctuations. It can be written as

$\displaystyle\quad\textrm{ASD}[\delta\rho_{\textrm{FV,corr,LRI}}](f)=\frac{\textrm{ASD}[\delta\nu_{\textrm{cav}}](f)}{\langle\nu_{\textrm{R}}(t)\rangle}\cdot
L(t)\quad\textrm{with}\quad\textrm{ASD}[\delta\nu_{\textrm{cav}}](f)=\frac{10^{-15}}{\sqrt{f}}\cdot\langle\nu_{\textrm{R}}(t)\rangle,$
(69)

which is based on an extrapolation of the LRI ranging noise from high to low
frequencies Spero (2021) and wherein
$L\approx$220\text{\,}\mathrm{k}\mathrm{m}$\pm$50\text{\,}\mathrm{k}\mathrm{m}$$
is the satellite separation. This noise level is in good agreement with pre-
flight measurements of the cavity stability Abich et al. (2019).
$\textrm{ASD}[\delta\rho_{\textrm{FV,corr,LRI}}]$ is shown as a solid blue
trace in the left plots of fig. 4. It is labeled as _frequency variations
corrected_ because this trace assumes that the second proper time term in eq.
(68) was removed. This removal can be accomplished using

$\displaystyle\rho_{\textrm{LRI1Bv04}}^{\textrm{corr}}(t)=\rho_{\textrm{LRI1Bv04}}(t)+\left(\frac{\mathrm{d}\tau_{\textrm{R}}(t)}{\mathrm{d}t}-1\right)\cdot
L(t).$ (70)

The proper time of the reference satellite $\tau_{\textrm{R}}$ and absolute
distance $L$ can be computed from GRACE-FO orbit data (GNI1B data product).
The correction is added to the biased range of LRI1B in order to obtain a
corrected biased range, which accounts for the variations of proper time in
the reference laser frequency. The spectrum of the raw laser frequency noise
can be approximated as the sum of the proper-time term and the cavity
stability as

$\displaystyle\textrm{PSD}[\delta\rho_{\textrm{FV,raw,LRI}}](f)\approx\textrm{PSD}\left[\frac{\mathrm{d}\tau_{\textrm{R}}}{\mathrm{d}t}\right]\cdot
L^{2}+\textrm{PSD}[\delta\rho_{\textrm{FV,corr,LRI}}],$ (71)

which is shown as a dashed blue trace in the left plots of fig. 4.

For the KBR, the carrier frequency variations at low frequencies can be
computed directly from $\epsilon_{\textrm{time}}$ in the CLK1B data product.
Currently, the KBR1B data processing uses daily mean values of the carrier
frequencies, which neglects intraday variations. In order to assess the
ranging error, we use the USO model introduced in eq. (59) and shown as a
dashed magenta trace in the right plot of fig. 1, and then compute the
contribution as

$\displaystyle\textrm{PSD}[\delta\rho_{\textrm{FV,raw,KBR}}]=\frac{\textrm{PSD}[\delta
t_{\textrm{USO,mod}}](f)\cdot(2\pi f)^{2}}{2}\cdot
L^{2}+\textrm{PSD}\left[\frac{\mathrm{d}(\tau_{\textrm{A}}+\tau_{\textrm{B}})}{2\mathrm{d}t}\right]\cdot
L^{2},$ (72)

which is depicted as a dashed blue line in the right plot of fig. 4. The first
term on the right side is again is the product of fractional frequency
variations and inter-satellite separation $L$. The factor of 2 in the
denominator (cf. (Thomas, 1999, B-7)) arises from the fact that two
independent USOs are used in the KBR and we assume that both USOs exhibit the
same amount of uncorrelated variations. The second summand accounts for the
effect of the proper time.

These frequency variations $\textrm{ASD}[\delta\rho_{\textrm{FV,raw,KBR}}]$
can be corrected by using time-dependent carrier frequencies as defined in
table 1, which are related to the CLK1B $\epsilon_{\textrm{time}}$ by eq. (2),
i.e.

$\displaystyle\nu_{\textrm{A/B}}^{\textrm{K}}(t)=\hat{f}_{\textrm{A/B,USO}}\cdot
5076\cdot\left(1-\frac{\mathrm{d}\epsilon_{\textrm{time,A/B}}}{\mathrm{d}t}\right)$
(73)
$\displaystyle\nu_{\textrm{A/B}}^{\textrm{Ka}}(t)=\hat{f}_{\textrm{A/B,USO}}\cdot
6768\cdot\left(1-\frac{\mathrm{d}\epsilon_{\textrm{time,A/B}}}{\mathrm{d}t}\right)$
(74)

The time-dependent carrier frequencies should be low-pass filtered in order to
remove high frequency fluctuations above
$3\text{\,}\mathrm{m}\mathrm{H}\mathrm{z}$, which limit the measurement
precision of $\epsilon_{\textrm{time}}$ and do not stem from real variations
in the carrier frequency (cf. right plot in fig. 1).

The corrected KBR range should ideally be computed in a modified level1A to
level1B processing scheme using eq. (54), which assumes that the re-scaling
from phase to range is performed with time-dependent carrier frequencies and
also takes the coupling of the phase-bias ($\Delta t_{\textrm{FV,DOWR}}$)
given in eq. (52) into account.

However, we point out that a direct correction of the KBR1Bv04 range data

$\displaystyle\rho_{\textrm{KBR1Bv04}}^{\textrm{corr}}(t)\approx\rho_{\textrm{KBR1Bv04}}(t)+c_{0}\cdot\Delta
t_{\textrm{FV,DOWR}}(t),$ (75)

might be performed as well, as the magnitude of $c_{0}\Delta
t_{\textrm{FV,DOWR}}$ has a much larger effect than dividing the phase by a
time-dependent frequency compared to division by a daily mean frequency. This
is illustrated in fig. 7, which also confirms that the difference between both
correction approaches is below the KBR noise level ($\approx$ 0.6
µm/$\sqrt{\textrm{Hz}}$, black trace).

Figure 7: Difference between the official KBR1B biased range (v04) and a
corrected KBR1B data product, from which USO frequency variations were
removed. The blue curve shows the difference after frequency variations have
been fully removed, while the black trace shows the difference calculated if
only the division of the phase with a time-dependent frequency is considered
without $\Delta t_{\textrm{FV,DOWR}}$-term. The red trace serves as a sanity
check to show that our self-derived KBR1B product can reproduce the official
KBR1Bv04 results up to some negligible deviations, which probably stem from
different interpolation methods. The spectrum was computed with GRACE-FO data
from 2019-01-01.

In order to derive the residual effect of USO frequency fluctuations in the
KBR ranging data, i.e. after applying the correction using the modified
level1a to level1b conversion, we use the model for the precision of CLK1B
$\epsilon_{\textrm{time}}$ data derived in sec. 6.4 as $\textrm{ASD}[\delta
t_{\epsilon}]$. The solid blue trace shows the residual KBR frequency
variations in the right plot of fig. 4 based on

$\displaystyle\textrm{ASD}[\delta\rho_{\textrm{FV,corr,KBR}}]=\frac{\textrm{ASD}[\delta
t_{\epsilon}](f)\cdot 2\pi f}{\sqrt{2}}\cdot L.$ (76)

We can conclude that correcting for KBR frequency variations might slightly
reduce the KBR-LRI residuals at frequencies between
$0.1\text{\,}\mathrm{m}\mathrm{H}\mathrm{z}$ and
$1\text{\,}\mathrm{m}\mathrm{H}\mathrm{z}$ (dashed blue vs solid blue trace).
However, the change is expected to be rather small as the timetag error (solid
black trace) becomes limiting.

##  7 Impact of Frequency Variations on KBR-LRI residuals

We have corrected the relativistic (proper-time) frequency variations in LRI
and the KBR frequency variations as discussed in previous sub-section 6.5 on a
daily basis for most of the days between December 2018 and March 2022 in the
GRACE-FO v04 data. Days where LRI1B, KBR1B or CLK1B data was incomplete, e.g.
due to gaps from instrument reboots or diagnostic events, were skipped. For
each day, we evaluated the spectrum of KBR-LRI residuals as shown exemplarily
in fig. 8. On day 2019-01-10, the spectrum of KBR-LRI residuals was reduced
when applying the correction, while day 2021-10-28 shows an increase in the
residual level (compare red and blue traces).

Figure 8: KBR-LRI residuals at low frequencies without applying the frequency
corrections (blue) and with applied LRI and KBR frequency correction (red).
The corrections are shown as yellow and purple traces.

As a metric, we use the root-mean-square (rms) value of KBR-LRI residuals at
low frequencies, which is obtained by integrating the power spectral density
(squared ASD) for frequencies below
$1\text{\,}\mathrm{m}\mathrm{H}\mathrm{z}$. The plots titled _With Rescale_
(right column) indicate that we re-estimated the differential scale and time-
shift between LRI and KBR together with a linear trend. This is in particular
necessary for days before 2019-08-12, as the LRI1B data is sometimes not
properly scaled and time-shifted in this period.

The change in daily rms value of LRI-LBR residuals at the low frequency due to
LRI and KBR frequency correction is shown in fig. 9. For January 2019, the rms
value could be lowered noticeably, but, the values oscillate with the beta
angle of the orbit, i.e. the angle between orbital plane and Sun direction.
The absolute rms values are depicted in fig. 10 and exhibit the oscillatory
pattern as well. Both plots indicate only a marginal reduction in rms due to
the frequency correction, e.g. a change in mean rms from 3.97 µm to 3.90 µm.

Figure 9: Daily change in the low frequency rms value
($<$1\text{\,}\mathrm{m}\mathrm{H}\mathrm{z}$$) of the KBR-LRI residuals when
applying the KBR and LRI frequency variation correction. Both plots show the
difference in rms value before and after the correction. The left plot uses
original LRI1Bv04 data without re-estimating the scale and time-shift between
LRI and KBR, while the right plot assumes re-estimation (rescale) of scale,
time-shift and linear trend on a daily basis. Values above the cyan dashed
line indicate a reduction in rms value due to the correction. On days with a
high absolute value of beta angle, i.e. close to or in full-sun phase of the
satellites, the frequency correction yields the largest rms reductions, e.g.
in January 2019. Figure 10: Daily absolute value of low frequency
($<$1\text{\,}\mathrm{m}\mathrm{H}\mathrm{z}$$) rms of the KBR-LRI residuals.
The rms value of residuals is low when the beta angle is large, i.e. when
satellites are close to or in full-sun phase and not entering shadow. The left
plot uses original LRI1Bv04 data with its scale and time-shift, while the
right plot assumes re-estimation (rescale) of scale, time-shift and linear
trend on a daily basis.

We repeated the analysis with more empirical parameters being estimated in the
daily rescale step, namely, the 1/rev, 2/rev and 2/day oscillation amplitudes
in terms of sine and cosine components and a quadratic trend, in addition to
the differential scale, time-shift and linear trend already in use. The 1/rev
and 2/rev frequencies in the ranging data are likely to have significant
errors, as environmental parameters such as temperature, inertial orientation
and magnetic field, predominately vary at these frequencies. The 2/day
parameter was added due to the CLK1B $\epsilon_{\textrm{time}}$ exhibiting a
pronounced 2/day feature on several days (cf. fig. 1) which could be related
to Earth’s rotation period or to the orbit period of satellites in the GPS
constellation.

We then estimated the KBR-LRI residuals together with the many empirical
parameters again, once without the frequency correction and once with the
correction for individual days. The results are shown in fig. 11. The
additional empirical parameters absorbed the oscillatory beta angle-related
behavior in the rms value to a large extent and lowered the rms value from
approx. 4 µm to approx. 1 µm, as had been expected. A change in rms value due
to the frequency correction thereby became visible more clearly, i.e. from
mean rms of 1.14 µm to 0.99 µm (cf. legend in figure). This 0.15 µm change in
rms in the 1 mHz bandwidth, which would imply a spectral density improvement
of 0.15 µm/$\sqrt{\textrm{1\,mHz}}=$ 4.7 µm/$\sqrt{\textrm{Hz}}$ assuming
white noise behavior, excludes effects at 1/rev, 2/rev and 2/day, because KBR-
LRI residuals at these frequencies are highly reduced due the empirical
parameters. Therefore, this change should be understood as a lower bound as it
neglects the frequencies where the correction has most of its signal.

The reduction in rms demonstrates that applying the frequency corrections,
i.e. considering intraday variability of the carrier frequencies, slightly
improves the agreement of LRI and KBR ranging data. We also emphasize that the
relativistic effect due to the proper time, i.e. 0.9 µm-peak at 1/rev (cf. eq.
(42)), is omitted in both instruments, if no frequency correction is applied.
If the corresponding frequency correction is applied to both instruments, the
0.9 µm-peak is present in both instruments. Thus, it never appears in the KBR-
LRI residuals. However, neglecting the effect would introduce an error
compared to a true error-free range.

The effect of the frequency corrections on LRI and KBR is small, i.e. at the
micron scale at very low frequencies, which means it should have negligible
impact on gravity fields, as typical pre-fit and post-fit residuals of gravity
fields are much larger than the KBR-LRI residuals. However, these corrections
reduce the errors in the instruments and are a necessary step towards a better
understanding of the KBR-LRI residuals at low frequencies, e.g. when studying
the tone errors of instruments at 1/rev and 2/rev frequencies.

Figure 11: Results as in fig. 9 and 10, but now using more empirical
parameters for the daily estimates: 1/rev, 2/rev, 2/day, scale, time-shift and
with linear and quadratic trends.

##  8 Summary and Conclusions

GRACE Follow-On offers the unique opportunity to compare two almost
independent ranging instruments, the KBR and the LRI. Both instruments rely on
readout of phase changes present in the beat note of interfering microwave or
near-infrared laser radiation. The measured phase values can be converted to a
range with units of meter using the respective carrier frequencies
$\nu_{\textrm{A/B}}^{\textrm{K/Ka}}$ or $\nu_{\textrm{R}}$.

In this paper, we have derived and analyzed different formulas for converting
the phase to a range when the time-dependency of the carrier frequency is
taken into account. It turned out that the previously suggested approach of
differentiating the phase, rescaling the derivative using the frequency and
performing an integration afterwards is not exact, since a summand was
neglected. In this work, we presented the exact formula, which also is likely
to be relevant for future laser ranging systems. Besides the analytically
exact formula, several approximate formulas can be used that yield results
with picometer precision, which is usually sufficient in a GRACE-like context.

We stressed that the time-variability of the carrier frequency $\nu$, if it is
neglected in the phase-to-range conversion, yields an error that is
proportional to the satellite separation $L$, commonly called laser frequency
noise in the LRI. The time-variability is, on the one hand, caused by
relativistic effects from the proper time, which is common to both LRI and
KBR, and on the other hand by the instability of the optical cavity for LRI
and by the instability of the USO for KBR. In the case of the KBR, the
microwave carrier frequencies can be determined in post-processing from GPS
precise orbit and clock error determination.

In the second part of this paper, the difference between measured range of LRI
and KBR in the official v04 data products was compared to some potential noise
and error contributors. It was pointed out that the KBR-LRI residuals are
bound at frequencies above a few mHz to a level of approx. 0.6
µm/$\sqrt{\textrm{Hz}}$ due to the CNR-limited phase readout noise from KBR,
though the analytical model predicts a noise level of approx. 0.4
µm/$\sqrt{\textrm{Hz}}$. This small discrepancy has not been understood so far
and is caused by the K-band observations, while Ka-band measurements seem to
be consistent with the model.

According to our analysis, two effects are limiting the KBR-LRI residuals at
low frequencies. Residual timetag errors in the KBR system after applying the
GPS-derived clock correction, which cannot be further improved with post-
processing, and the intraday carrier frequency variations in the KBR and LRI,
which can be corrected using the formulas introduced in the first half of this
paper. The LRI frequency variations can only be corrected for the proper time
effect using GNI1B orbit data, which contain dominant 1/rev and 2/rev
components. In GRACE-FO, there is no direct measurement of the absolute laser
frequency, and thus, the laser frequency is assumed to be static on a daily
basis and its daily values are derived from correlating LRI and KBR
observations. For KBR, the microwave frequencies are derived from the USO and
potential frequency variations can be determined with respect to GPS data. The
corresponding KBR frequency correction accounting for intraday variability can
be computed directly from the CLK1B data product and exhibits a continuous
spectrum as well as a 1/rev peak due to the proper time effect.

In the end, we applied the KBR and LRI frequency corrections and were able to
demonstrate a slight improvement in the agreement between KBR and LRI range,
i.e. the rms value of KBR-LRI residuals below 1 mHz was reduced by 0.15 µm.
However, that figure excludes the 1/rev, 2/rev and 2/day frequencies as well
as linear and quadratic drifts. Hence, the actual reduction might be larger.

Moreover, applying the frequency correction to KBR and LRI should reduce the
error at the 1/rev frequency by approx. 0.9 µm-peak due to the relativistic
proper time effect, as the current regular processing that generates v04 data
is not accounting for this effect in LRI and KBR data.

We expect that the micron-scale effects at low frequencies discussed here have
negligible impact on the current gravity field determination, as typical pre-
fit and post-fit residuals of gravity field determination are much larger than
direct KBR-LRI residuals. However, the corrections are of interest, if one
attempts to understand the low-frequency behavior of the instruments, e.g. the
tone errors at 1/rev and 2/rev frequency. These corrections might also be
important for future missions employing laser ranging that aim for higher
precision at sub-mHz frequencies.

Conceptualization, Vitali Müller; Funding acquisition, Gerhard Heinzel;
Investigation, Markus Hauk, Malte Misfeldt, Laura Müller and Henry Wegener;
Project administration, Gerhard Heinzel; Writing – original draft, Vitali
Müller and Yihao Yan; Writing – review & editing, Markus Hauk, Malte Misfeldt,
Laura Müller, Henry Wegener and Gerhard Heinzel.

Data can be made available upon reasonable request.

This work has been supported by: The Deutsche Forschungsgemeinschaft (DFG,
German Research Foundation, Project-ID 434617780, SFB 1464). Clusters of
Excellence “QuantumFrontiers: Light and Matter at the Quantum Frontier:
Foundations and Applications in Metrology” (EXC-2123, project number:
390837967); the European Space Agency in the framework of Next Generation
Geodesy Mission development and ESA’s third-party mission support for GRACE-
FO; the Chinese Academy of Sciences (CAS) and the Max Planck Society (MPG) in
the framework of the LEGACY cooperation on low-frequency gravitational-wave
astronomy (M.IF.A.QOP18098).

The authors declare no conflict of interest.

no

## Appendix A

### A.1

In GRACE-like missions, the inter-satellite distance $L(t)$ can be modelled
with the following dominant components

$\displaystyle L(t)=L_{0}+L_{1}\sin(2\pi f_{\textrm{orb}}t)+L_{d}\cdot t$ (77)

such that the propagation time can be described as

$\displaystyle\Delta t(t)\approx L(t)/c_{0}.$ (78)

We assume the following frequency model for the LRI laser with constant,
oscillatory and linear drift components:

$\displaystyle\nu_{\textrm{R}}(t)=\nu_{0}+\nu_{1}\sin(2\pi f_{\nu
1}t)+\nu_{d}\cdot t.$ (79)

Finally, the phase of the laser can be written as

$\displaystyle\Phi_{\textrm{R}}(t)=\int_{0}^{t}\nu(t^{\prime})~{}\mathrm{d}t^{\prime}=\nu_{0}t-\frac{\nu_{1}\cos(2\pi
f_{\nu 1}t)}{2\pi f_{\nu 1}}+\frac{\nu_{d}\cdot t^{2}}{2}$ (80)

or as

$\displaystyle\varphi_{\textrm{TWR}}(t)=\Phi_{\textrm{R}}(t)-\Phi_{\textrm{R}}(t-\Delta
t)$ (81)

in terms of the ranging phase $\varphi_{\textrm{TWR}}$. These analytical
models of the range and phase, with the numerical values given in table 2, can
be used to validate the expressions derived in sec. 3.

Table 2: Numerical values for a GRACE-FO scenario.  Name | Formula/Value | Comment
---|---|---
$\nu_{0}$ | $282\text{\,}\mathrm{T}\mathrm{H}\mathrm{z}$ | LRI optical frequency (1064.5 nm)
$\nu_{1}$ | $4\cdot 10^{-12}\cdot\nu_{0}$ | proper time 1/rev variation (cf. fig. 1)
$\nu_{d}$ | $3.6\cdot 10^{-15}\cdot 1/s\cdot\nu_{0}$ | assumption, $87\text{\,}\mathrm{k}\mathrm{H}\mathrm{z}\mathrm{/}\mathrm{d}\mathrm{a}\mathrm{y}$
$L_{0}$ | $220\text{\,}\mathrm{k}\mathrm{m}$ | satellite separation DC
$L_{1}$ | $400\text{\,}\mathrm{m}$ | 1/rev range variation
$L_{d}$ | $0.01\text{\,}\mathrm{m}\mathrm{/}\mathrm{s}$ | linear range drift
$f_{\textrm{orb}}$ | $0.176\text{\,}\mathrm{m}\mathrm{H}\mathrm{z}$ | orbit frequency, 1/rev

References

## References

  * Landerer et al. (2020) Landerer, F.W.; Flechtner, F.M.; Save, H.; Webb, F.H.; Bandikova, T.; Bertiger, W.I.; Bettadpur, S.V.; Byun, S.H.; Dahle, C.; Dobslaw, H.; et al.  Extending the global mass change data record: GRACE Follow-On instrument and science data performance.  Geophysical Research Letters 2020, 47, e2020GL088306.  https://doi.org/10.1029/2020GL088306.
  * Tapley et al. (2019) Tapley, B.D.; Watkins, M.M.; Flechtner, F.; Reigber, C.; Bettadpur, S.; Rodell, M.; Sasgen, I.; Famiglietti, J.S.; Landerer, F.W.; Chambers, D.P.; et al.  Contributions of GRACE to understanding climate change.  Nature climate change 2019, 9, 358–369.  https://doi.org/10.1038/s41558-019-0456-2.
  * Landerer and Swenson (2012) Landerer, F.W.; Swenson, S.  Accuracy of scaled GRACE terrestrial water storage estimates.  Water resources research 2012, 48.  https://doi.org/10.1029/2011WR011453.
  * Han et al. (2006) Han, S.C.; Shum, C.K.; Bevis, M.; Ji, C.; Kuo, C.Y.  Crustal dilatation observed by GRACE after the 2004 Sumatra-Andaman earthquake.  Science 2006, 313, 658–662.  https://doi.org/10.1126/science.1128661.
  * Ciracì et al. (2020) Ciracì, E.; Velicogna, I.; Swenson, S.  Continuity of the mass loss of the world’s glaciers and ice caps from the GRACE and GRACE Follow-On missions.  Geophysical Research Letters 2020, 47, e2019GL086926.  https://doi.org/10.1029/2019GL086926.
  * Chen et al. (2020) Chen, J.; Tapley, B.; Wilson, C.; Cazenave, A.; Seo, K.W.; Kim, J.S.  Global ocean mass change from GRACE and GRACE Follow-On and altimeter and Argo measurements.  Geophysical Research Letters 2020, 47, e2020GL090656.  https://doi.org/10.1029/2020GL090656.
  * Chen et al. (2014) Chen, J.; Li, J.; Zhang, Z.; Ni, S.  Long-term groundwater variations in Northwest India from satellite gravity measurements.  Global and Planetary Change 2014, 116, 130–138.  https://doi.org/https://doi.org/10.1016/j.gloplacha.2014.02.007.
  * Kornfeld et al. (2019) Kornfeld, R.P.; Arnold, B.W.; Gross, M.A.; Dahya, N.T.; Klipstein, W.M.; Gath, P.F.; Bettadpur, S.  GRACE-FO: The Gravity Recovery and Climate Experiment Follow-On Mission.  Journal of Spacecraft and Rockets 2019, 56, 931–951.  https://doi.org/10.2514/1.A34326.
  * Ghobadi-Far et al. (2020) Ghobadi-Far, K.; Han, S.C.; McCullough, C.M.; Wiese, D.N.; Yuan, D.N.; Landerer, F.W.; Sauber, J.; Watkins, M.M.  GRACE Follow-On laser ranging interferometer measurements uniquely distinguish short-wavelength gravitational perturbations.  Geophysical Research Letters 2020, 47, e2020GL089445.  https://doi.org/10.1029/2020GL089445.
  * Nicklaus et al. (2020) Nicklaus, K.; Cesare, S.; Massotti, L.; Bonino, L.; Mottini, S.; Pisani, M.; Silvestrin, P.  Laser metrology concept consolidation for NGGM.  CEAS Space Journal 2020, 12, 313–330.  https://doi.org/10.1007/s12567-020-00324-6.
  * Wiese et al. (2020) Wiese, D.; Boening, C.; Zlotnicki, V.; Luthcke, S.; Loomis, B.; Rodell, M.; Sauber, J.; Bearden, D.; Chrone, J.; Horner, S.; et al.  The NASA Mass Change Designated Observable Study: Overview, Progress, and Future Plans.  In Proceedings of the IGARSS 2020-2020 IEEE International Geoscience and Remote Sensing Symposium. IEEE, 2020, pp. 6121–6124.  https://doi.org/10.1109/IGARSS39084.2020.9323225.
  * Flechtner (2020) Flechtner, F.  Realization of a satellite mission “GRACE-I” for parallel observation of changing global water resources and biodiversity.  Technical report, Copernicus Meetings, 2020.
  * Burt et al. (2021) Burt, E.; Prestage, J.; Tjoelker, R.; Enzer, D.; Kuang, D.; Murphy, D.; Robison, D.; Seubert, J.; Wang, R.; Ely, T.  Demonstration of a trapped-ion atomic clock in space.  Nature 2021, 595, 43–47.  https://doi.org/10.1038/s41586-021-03571-7.
  * Thomas (1999) Thomas, J.  An analysis of gravity-field estimation based on intersatellite dual-1-way biased ranging (No. JPL-98-15.).  Technical report, Jet Propulsion Laboratory (JPL), Pasadena, CA, 1999\.  ftp://isdcftp.gfz-potsdam.de/grace/DOCUMENTS/ANCILLARY/.
  * Heinzel et al. (2002) Heinzel, G.; Rüdiger, A.; Schilling, R.  Spectrum and spectral density estimation by the Discrete Fourier transform (DFT), including a comprehensive list of window functions and some new at-top windows.  Technical report, Albert-Einstein-Institute Hannover, 2002.  http://hdl.handle.net/11858/00-001M-0000-0013-557A-5.
  * Shiozawa (2004) Shiozawa, T.  Classical relativistic electrodynamics: theory of light emission and application to free electron lasers; Springer Science & Business Media, 2004\.  https://doi.org/10.10071978-3-662-06261-6.
  * Wen et al. (2019) Wen, H.; Gerhard, K.; William, B.; Meegyeong, P.; Carly, S.; Felix, L.  Gravity Recovery and Climate Experiment (GRACE) Follow-On (GRACE-FO) Level-1 Data Product User Handbook.  JPL D-56935, 2019.  https://podaac-tools.jpl.nasa.gov/drive/files/allData/gracefo/docs/GRACE-FO_L1_Handbook.pdf.
  * Abich et al. (2019) Abich, K.; Abramovici, A.; Amparan, B.; Baatzsch, A.; Okihiro, B.B.; Barr, D.C.; Bize, M.P.; Bogan, C.; Braxmaier, C.; Burke, M.J.; et al.  In-Orbit Performance of the GRACE Follow-on Laser Ranging Interferometer.  Physical review letters 2019, 123, 031101.  https://doi.org/10.1103/physrevlett.123.031101.
  * Müller (2017) Müller, V.  Design Considerations for Future Geodesy Missions and for Space Laser Interferometry.  PhD thesis, Gottfried Wilhelm Leibniz Universität Hannover, Germany, 2017.  https://doi.org/10.15488/9029.
  * Yan et al. (2021) Yan, Y.; Müller, V.; Heinzel, G.; Zhong, M.  Revisiting the light time correction in gravimetric missions like GRACE and GRACE follow-on.  Journal of Geodesy 2021, 95, 1–19.  https://doi.org/10.1007/s00190-021-01498-5.
  * Kim (2000) Kim, J.  Simulation study of a low-low satellite-to-satellite tracking mission.  PhD thesis, University of Texas, Austin, USA, 2000.  ftp://isdcftp.gfz-potsdam.de/grace/DOCUMENTS/ANCILLARY/.
  * Sheard et al. (2012) Sheard, B.; Heinzel, G.; Danzmann, K.; Shaddock, D.; Klipstein, W.; Folkner, W.  Intersatellite laser ranging instrument for the GRACE follow-on mission.  Journal of Geodesy 2012, 86, 1083–1095.  https://doi.org/10.1007/s00190-012-0566-3.
  * Case et al. (2010) Case, K.; Gerhard, K.; Wu, S.C.  Gravity Recovery and Climate Experiment (GRACE) Level-1B Data Product User Handbook.  JPL D-22027, 2010.  https://opendap.jpl.nasa.gov/opendap/allData/grace/docs/Handbook_1B_v1.3.pdf.
  * Misfeldt et al. (2021) Misfeldt, M.; Müller, V.; Müller, L.; Wegener, H.; Heinzel, G.  Thermal Influence on the LRI Scale Factor.  In Proceedings of the EGU General Assembly Conference Abstracts, 2021, pp. EGU21–1242.  https://doi.org/10.5194/egusphere-egu21-1242.
  * Rees et al. (2021) Rees, E.R.; Wade, A.R.; Sutton, A.J.; Spero, R.E.; Shaddock, D.A.; Mckenzie, K.  Absolute frequency readout derived from ULE cavity for next generation geodesy missions.  Opt. Express 2021, 29, 26014–26027.  https://doi.org/10.1364/OE.434483.
  * Hui Ying Wen (2019) Hui Ying Wen, S.  Release Notes of GRACE-FO Level1B data.  Technical report, PODAAC, 2019.  https://podaac.jpl.nasa.gov/gravity/gracefo-documentation.
  * Müller et al. (2020) Müller, L.; Müller, V.; Misfeldt, M.; Wegener, H.; Heinzel, G.  An alternative derivation of the GRACE Follow-On LRI1B product: current status.  Technical report, Copernicus Meetings, GRACE-FO Science Team Meeting, 2020\.  https://doi.org/10.5194/gstm2020-27.
  * Wegener et al. (2020) Wegener, H.; Müller, V.; Heinzel, G.; Misfeldt, M.  Tilt-to-Length Coupling in the GRACE Follow-On Laser Ranging Interferometer.  Journal of Spacecraft and Rockets 2020, 57, 1362–1372.  https://doi.org/10.2514/1.A34790.
  * Hernández-Pajares et al. (2007) Hernández-Pajares, M.; Juan, J.; Sanz, J.; Orús, R.  Second-order ionospheric term in GPS: Implementation and impact on geodetic estimates.  Journal of geophysical research: solid earth 2007, 112.  https://doi.org/10.1029/2006JB004707.
  * Wegener (2022) Wegener, H.  Analysis of Tilt-to-Length Coupling in the GRACE Follow-On Laser Ranging Interferometer.  PhD thesis, Gottfried Wilhelm Leibniz Universität Hannover, Germany, 2022.  https://doi.org/10.15488/11984.
  * Spero (2021) Spero, R.  Point-mass sensitivity of gravimetric satellites.  Advances in Space Research 2021, 67, 1656–1664.  https://doi.org/10.1016/j.asr.2020.12.019.