Ranging Sensor Fusion in LISA Data Processing: Treatment of Ambiguities, Noise, and On-Board Delays in LISA Ranging Observables

Interspacecraft ranging is crucial for the suppression of laser frequency noise via time-delay interferometry (TDI). So far, the effects of on-board delays and ambiguities on the LISA ranging observables were neglected in LISA modelling and data processing investigations. In reality, on-board delays cause offsets and timestamping delays in the LISA measurements, and pseudo-random noise (PRN) ranging is ambiguous, as it only determines the range up to an integer multiple of the PRN code length. In this article, we identify the four LISA ranging observables: PRN ranging, the sideband beatnotes at the interspacecraft interferometer, TDI ranging, and ground-based observations. We derive their observation equations in the presence of on-board delays, noise, and ambiguities. We then propose a three-stage ranging sensor fusion to combine these observables in order to gain accurate and precise ranging estimates. We propose to calibrate the on-board delays on ground and to compensate the associated offsets and timestamping delays in an initial data treatment (stage 1). We identify the ranging-related routines, which need to run continuously during operation (stage 2), and implement them numerically. Essentially, this involves the reduction of ranging noise, for which we develop a Kalman filter combining the PRN ranging and the sideband beatnotes. We further implement crosschecks for the PRN ranging ambiguities and offsets (stage 3). We show that both ground-based observations and TDI ranging can be used to resolve the PRN ranging ambiguities. Moreover, we apply TDI ranging to estimate the PRN ranging offsets.


I. INTRODUCTION
The Laser Interferometer Space Antenna (LISA), due for launch around the year 2035, is an ESA-led mission for space-based gravitational-wave detection in the frequency band between 0.1 mHz and 1 Hz [1].LISA consists of three satellites forming an approximate equilateral triangle with an armlength of 2.5 Gm, in a heliocentric orbit that trails or leads Earth by about 20 degrees.Six infrared laser links with a nominal wavelength of 1064 nm connect the three spacecraft (SC), whose relative motion necessitates the usage of heterodyne interferometry.Phasemeters are used to extract the phases of the corresponding beatnotes [2], in which gravitational-waves manifest in form of microcycle deviations equivalent to picometer variations in the interspacecraft ranges.
The phasemeter output, however, is obscured by various instrumental noise sources.They must be suppressed to fit in the LISA noise budget of 10 pm Hz −0.5 (single link) [3], otherwise they would bury the gravitationalwave signals.Dedicated data processing algorithms are being developed for each of these instrumental noise sources, their subsequent execution is referred to as initial noise reduction pipeline (INReP).The dominating noise source in LISA is by far the laser frequency noise, * janniklas.reinhardt@aei.mpg.de which must be reduced by more than 8 orders of magnitude.This is achieved by time-delay interferometry (TDI), which combines the various beatnotes with the correct delays to virtually form equal-optical-path-length interferometers, in which laser frequency noise naturally cancels [4,5].The exact definition of these delays depends on the location of TDI within the INReP (see fig. 1) [6], but wherever we place it, some kind of information about the absolute interspacecraft ranges is required.
Yet, absolute ranges are not a natural signal in a continuous-wave heterodyne laser interferometer such as LISA.Therefore, a ranging scheme based on pseudorandom noise (PRN) codes is implemented [7][8][9]: each SC houses a free-running ultra-stable oscillator (USO) as timing reference.It defines the spacecraft elapsed time (SCET).PRN codes generated according to the respective SCETs are imprinted onto the laser beams by phasemodulating the carrier.The comparison of a PRN code received from a distant SC, hence generated according to the distant SCET, with a local copy enables a measurement of the pseudorange: the pseudorange is commonly defined as the difference between the SCET of the receiving SC at the event of reception and the SCET of the emitting SC at the event of emission [10].It represents a combination of the true geometrical range (light travel time) with the offset between the two involved SCETs (see eq. A5).
In the baseline TDI topology (upper row in fig.1), TDI is performed after SCET synchronization to the barycen- In the baseline TDI topology (upper part) we perform TDI after clock synchronization to TCB, the delays are given by the light travel times.In the alternative TDI topology (lower part) we execute TDI without clock synchronization and apply the pseudoranges as delays [6].Both topologies rely on a ranging sensor fusion.
tric coordinate time (TCB), the light travel times are used as delays.The pseudoranges comprise information about both the light travel times and the SCET offsets required for synchronizing the clocks (see appendix A).A Kalman filter can be used to disentangle the pseudoranges in order to retrieve light travel times and SCET offsets [11].In the alternative TDI topology (lower row in fig.1), the pseudoranges are directly used as delays.In that topology, TDI is executed on the unsynchronized beatnotes sampled according to the respective SCETs [6].
However, PRN ranging (PRNR) does not directly provide the pseudoranges but requires three treatments.First, due to the finite PRN code length (we assume 400 km), PRNR measures the pseudoranges modulo an ambiguity [7].Secondly, on-board delays due to signal propagation and processing cause offsets and timestamping delays in the PRNR.Thirdly, PRNR is limited by white ranging noise with an rms amplitude of about 0.3 m, which is due to shot noise and PRN code interference [9,12].To overcome these difficulties, there are three additional pseudorange observables, which are actually designed for other purposes and serve that function secondarily: ground-based observations provide inaccurate but unambiguous pseudorange estimates; time-delay interferometric ranging (TDIR) turns TDI upside-down seeking a model for the delays that minimizes the laser frequency noise in the TDI combinations [13]; the sideband beatnotes, primarily designed for clock noise correction, provide a measurement of the pseudorange time derivatives [6].The combination of these four observables in order to form optimal pseudorange estimates is referred to as ranging sensor fusion in the course of this article.It is common to both TDI topologies (see fig. 1) and consequently a crucial stage of the INReP.
In section II, we specify the pseudorange definition in the context of on-board delays and identify the delays required for TDI.We then derive observation equations for the four pseudorange observables in section III.Here, we carefully consider the effects of on-board delays at  [14]).The SC are labeled clockwise.The MOSAs and their associated building blocks (lasers, interferometers, etc.) are labeled by 2 indices: the first one indicates the SC they are located at, the second one the SC they are oriented to.The measurements and related quantities (optical links, pseudoranges, etc.) share the indices of the MOSAs they are measured at.Below, we distinguish between left-handed MOSAs (12,23,31) and righthanded ones (13,32,21).
the interspacecraft interferometer. 1 In section IV, we introduce a three-stage ranging sensor fusion consisting of an initial data treatment, a core ranging processing, and crosschecks.In the initial data treatment, we propose to compensate for the offsets and timestamping delays caused by the on-board delays.We identify PRNR unwrapping and noise reduction as the ranging processing steps that need to run continuously during operation.In parallel to this core ranging processing, we propose crosschecks of the PRNR ambiguities and offsets.We implement the core ranging processing and the crosschecks numerically.In section V we discuss the performance of this implementation, and conclude in section VI.
Each SC houses an ultra-stable oscillator (USO) generating a signal at about 50 MHz.An 80 MHz clock signal, the phasemeter clock (PMC), is coherently derived from this USO.The PMC can be considered as the timing reference on board the SC (see fig. 4), its associated counter is referred to as spacecraft elapsed time (SCET): It is useful to consider the SCET as a continuous time scale, which we denote by τ .It differs from the barycentric coordinate time (TCB), denoted by t, due to instrumental clock drifts and jitters, and due to relativistic effects.Following the notation of [15], we use superscripts to indicate a quantity to be expressed as function of a certain time scale, e.g., τ t 1 denotes the SCET of SC 1 as function of TCB.Note that Each SC contains two movable optical sub-assemblies (MOSAs) connected by an optical fibre (see fig. 2 for the labeling conventions).Each MOSA has an associated laser and houses a telescope, a free-falling test mass marking the end of the corresponding optical link, and an optical bench with three interferometers: the interspacecraft interferometer (ISI), in which the gravitationalwave signals eventually appear, the reference interferometer (RFI) to compare local and adjacent lasers, and the test-mass interferometer (TMI) to sense the optical bench motion with respect to the free-falling test mass in direction of the optical link.The MHz beatnotes in these interferometers are detected with quadrant-photoreceivers (QPRs).They are digitized in analog-to-digital converters (ADCs) driven by the PMCs.Phasemeters extract the beatnote phases2 using digital phase-locked loops (DPLLs).These phases are then downsampled to 4 Hz in a multi-stage decimation procedure (DEC) and telemetered to Earth (see fig. 4).

B. The pseudorange and on-board delays
The pseudorange, denoted by R τi ij , is commonly defined as the difference between the SCET of the receiving SC at the event of reception and the SCET of the emitting SC at the event of emission [10].It represents a combination of the light travel time between the emission at SC j and the reception at SC i, and the differential SCET offset (see eq. A5).However, considering the complexity of the LISA metrology system, this definition appears to be rather vague: to what exactly do we relate the events of emission and reception?Two specifications are required here: we need to locate emission and reception, and we need to define the actual events.
It is convenient to consider emission and reception at the respective polarizing beam splitters (PBSs) in front of the telescopes (denoted PBS1 in [16]), and to treat the on-board signal propagation and processing on both SC as on-board delays (see fig. 3).Thus, we clearly separate the pseudorange from on-board delays.Note that this definition is not unique, the events of emission and reception could be located elsewhere, assuming that the on-board delays are defined accordingly.
The LISA optical links do not involve delta-pulse-like events.In order to define the actual events of emission and reception we, instead, use the instants when the light phase changes at the beginning of the first PRN code chip.At first glance, the PRN code might seem unfavorable for the pseudorange definition, as PRN and carrier phase are oppositely affected by the solar wind: the PRN phase is delayed by the group-delay, while the carrier phase is advanced by the phase delay.However, these effects are at the order of 10 pm (see appendix C), whereas our best pseudorange estimates are at 0.1 mm accuracy.Consequently, the solar wind dispersion can be neglected in the pseudorange definition.
When expressing the interferometric measurements according to this specified pseudorange definition, we need to consider the excluded on-board signal propagation and processing.For that purpose, we introduce two kinds of delay operators by their action on a function f τj .The on-board delay operator describes delays due to on-board signal propagation and processing and is defined on the same SCET as the function it is acting on: x is a place holder for any on-board delay, e.g., D qpr ← pbs denotes the optical path length from the PBS to the QPR and D dec the decimation filter group delay.The interspacecraft delay operator is defined on a different SCET than the function it is acting on and applies the pseudorange as delay: For on-board delays that differ between carrier, PRN, and sideband signals, we add the superscripts car, prn, and sb, respectively.To trace the full path of a signal from the distant SC, we need to combine the interspacecraft delay operator for the interspacecraft signal propagation and the SCET conversion (considered at the PBS of the receiving SC) with on-board delay operators on both SC.The application of a delay operator to another time-dependent delay operator results in nested delays: For a constant delay operator D x we can define the associated advancement operator A x acting as its inverse: Constituents of the pseudorange are marked purple.These are the light travel time between the PBSs (at the telescopes) and the transformation between the two SCETs.In section II C we identify the TDI delay that is required to cancel the noise of laser 12.For that purpose we mark the common path of the outgoing and the local beam in dark gray, the noncommon path of the local beam in light blue, and the noncommon path of the outgoing beam in green.
For advancement operators associated to propagation delays we write the subscript underlines that the advancement operator undoes the signal propagation.Below, we consider onboard delays as constant or slowly time varying so that their associated advancement operators are well-defined.

C. Delays for TDI
In [6] the pseudoranges are identified as the delays that need to be applied to cancel the laser noise in the alternative TDI topology (see fig. 1).Does this statement hold in the presence of on-board delays considering the refined pseudorange definition of section II B? To identify the delays D τi ij that are required in TDI combinations to suppress the laser noise, let us set up a simple TDI toy model: we consider the 2 MOSAs depicted in fig. 3 and the TDI combination, where we combine ISI τ2  21 with ISI τ1 12 delayed by D τ2 21 .We will identify the expression of D τ2 21 , which leads to a suppression of the noise from laser 12.
Let us first define the delays that are at play in this situation and highlighted in fig. 3. We denote the delay that is common to both the local and the outgoing beam by D τ1 A .It corresponds to the path from the laser source to the beam splitter (BS), which devides them (denoted BS2 in [16]).The delays associated to the paths only traversed by the local and the outgoing beam are called D τ1 B and D τ2 C .Note that they represent combinations of the delays from the BS to the recombining beam splitters at the respective interspacecraft interferometers (ISI BSs) (denoted BS12 in [16]), where the measurements are formed, with delays from the ISI BSs to the decimation filters, where the measurements are timestamped.Moreover, D τ2 C combines an interspacecraft delay operator with on-board delays on both SC.Accordingly, we define the delays D τ2 A , D τ2 B , and D τ1 C for laser 21.Now the TDI delay D τ2 21 can be derived as: Φ τ1 12 denotes the phase of laser 12.In the second step we factor out the common delay D τ1 A , which does not contribute to the TDI delay D τ2 21 .We ignore the phase Φ τ2 21 of laser 21, as we focus on canceling the noise of laser 12.For that purpose we need to choose the delay D τ2 21 , such that the first bracket vanishes, i.e.
The advancement operators are associated to the noncommon path of the local beam, the delay operators to the noncommon path of the outgoing beam.Hence, to cancel the noise of laser 12 we need to consider the difference between the delays applied to laser 12 in the ISI measurements on SC 1 and SC 2.
The TDI delays D τi ij represent combinations of interspacecraft delay operators (pseudoranges) with on-board delay and advancement operators on both SC.We propose to calibrate all required on-board delays on ground before mission start and to measure the pseudorange during operation.The next section covers the four pseudorange observables.Before, we close this section with a few remarks on the required on-board delays.D τ1 pbs ← bs and D τ2 qpr ← pbs are small optical path length delays of the outgoing beam before transmission and after reception at the distant SC.A τ1 bs ← qpr is a small optical path length advancement of the local beam.These optical path lengths are in the order of 10 cm to 1 m [16].

D τ2
dec ← qpr is the signal processing delay on the receiving SC.A τ1 qpr ← dec is the corresponding advancement on the local SC.It can be decomposed into The group delays of the quadrant-photo-receiver D car qpr and the analog backend electronics D car abee depend amongst others on the beatnote frequency [17].Hence, they change over time and differ between carrier, sideband, and PRN signals.Together with the cable delays D abee← qpr and D dpll← abee they can amount to 10 m.The DPLL delay D car dpll depends on the time-dependent beatnote amplitude [2,12].All time-dependent contributions should be calibrated for all combinations of the time-dependent parameters.Hence, during operation they can be constructed with the help of SC monitors providing the corresponding parameter values, e.g., beatnote frequency and amplitude.

A. PRN ranging (PRNR)
PRNR is the on-board ranging scheme in LISA.A set of 6 PRN sequences has been computed such that the cross-correlations and the auto-correlations for nonzero delays are minimized.These PRN codes are associated to the 6 optical links in the LISA constellation.The PRN codes are generated according to the respective PMCs and imprinted onto the laser beams by phase-modulating the carriers in electro-optical modulators (EOMs) (see fig. 4).In each phasemeter, DPLLs are applied to extract the beatnote phases.The PRN codes show up in the DPLL error signals since the DPLL bandwidth of 10 kHz to 100 kHz is lower than the PRN chipping rate of about 1 MHz.In a delay-locked loop (DLL), these error signals are correlated with PRN codes generated according to the local SCET.The local delay that maximizes the correlation yields a pseudorange measurement: the PRNR [7,8].
We now derive the PRNR observation equation carefully taking into account on-board delays.We model the path of the PRN code from the distant SC to the local DLL by applying delay operators to the distant SCET: The two on-board delays can be decomposed into D prn pbs ← pmc consists of the cable delays from the PMC to the EOM, the processing delay due to the PRN code generation, and the optical path length from the EOM to the PBS.All these delays are constant at the sensitive scale of PRNR, so that we do not have to consider delay nesting in D prn pbs ← pmc .We added the superscript prn because this path is different for the sideband signal (see fig. 4).D prn dll ← pbs is explained in the next paragraph as part of the PRNR timestamping delay.At the DLL, the received PRN codes are correlated with identical codes generated according to the local SCET.We model this correlation as the difference between the local SCET and the delayed distant SCET (eq.15), and we apply D prn dec to model the group delay of the decimation filters applicable to PRN ranging: To see how the on-board delays affect the PRNR we expand eq.18 applying eq.2:  We show the USO frequency distribution (follow the arrows after the USOs) and illustrate the on-board signal processing (follow the arrows after the QPRs).We highlight the timestamping delays and offsets in the pseudorange observables caused by on-board delays.In light blue, we mark the PRNR offset from the pseudorange.The PRNR timestamping delay is drawn pink, the sideband timestamping delay is marked green.
The on-board delays cause the PRNR timestamping delay D prn dec ← pbs and the PRNR offset O prn ij : The PRNR timestamping delay has similar constituents as the on-board delays appearing in the TDI delay D ij , they are marked pink in fig. 4.However, most of them are frequency or amplitude dependent.Therefore, they differ between carrier and PRN signals.As for the TDI delay, we propose to individually calibrate all constituents of the PRNR timestamping delay on ground before mission start.Hence, during operation D prn dec ← pbs can be compensated in an initial data treatment by application of its associated advancement operator A prn pbs ← dec .After that, the PRNR observation equation including ranging noise and PRN ambiguity can be written as: l denotes the finite PRN code length.We use 400 km as a placeholder, the final value has not been decided.The finite PRN code length leads to an ambiguity, a prn ij denote the associated ambiguity integers [8].N prn ij is the white ranging noise with an rms amplitude of about 0.3 m at the 4 Hz data rate, which is due to shot noise and PRN code interference [9,12].The latter refers to the interference between the PRN code modulated onto the received laser at the distant SC and the other code modulated onto the local laser at the local SC.The PRNR offset O prn ij involves contributions on the emitter and on the receiver side (see eq. 21), they are marked light blue in fig. 4. It can amount to 10 m and more [12,18].Similar to the PRNR timestamping delay, we propose to calibrate the PRNR offset on ground, so that it can be subtracted in an initial data treatment.

B. Sideband ranging (SBR)
For the purpose of in-band clock noise correction in the INReP, a clock noise transfer between the SC is implemented [9]: the 80 MHz PMC signals are upconverted to ν m l = 2.400 GHz and ν m r = 2.401 GHz for left and righthanded MOSAs, respectively (see fig. 2 for the definition of left and right-handed MOSAs).The EOMs phasemodulate the carriers with the upconverted PMC signals, thereby creating clock sidebands.We show below that the beatnotes between these clock sidebands constitute a pseudorange observable.
Considering on-board delays, the difference between carrier and sideband beatnotes can be written as D sb pbs ← pmc and D sb isi bs ← pmc are the delay operators associated to the paths from the PMC to the PBS and from the PMC to the ISI BS.They can be decomposed into: D up is the upconversion delay due to phase-locking a 2.40(1) GHz oscillator to the 80 MHz PMC signal, D eom ← pmc is the cable delay from the PMC to the EOM.ν m ij is the upconverted USO frequency associated to MOSA ij .Since eq.23 is expressed in the SCET, all clock imperfections are included into τ τi i (τ ).The modulation noise M τi ij contains any additional jitter collected on the path D sb (p|isi)bs ← pmc , e.g., due to the electrical frequency upconverters.The amplitude spectral densities (ASDs) of the modulation noise for left and right-handed MOSAs are specified to be below [6,19] The modulation noise on left-handed MOSAs is one order of magnitude lower, because the pilot tone for the ADC jitter correction (it is the ultimate phase reference) is derived from the 2.400 GHz clock signal.
To derive a pseudorange observation equation from the sideband beatnote we expand eq.23 using eq. 2. We ap-ply A sb pbs ← dec to avoid nested delays in the pseudorange: We subtract the 1 MHz ramp and then refer to eq. 27 as sideband ranging (SBR).Taking into account that the SBR phase is defined up to a cycle, the SBR can be written as a sb ij denote the SBR ambiguity integers.Expressed as length, the SBR ambiguity is 12.5 cm corresponding to the wavelength of the GHz sidebands.The SBR offset can be thought of as the differential phase accumulation of local and distant PMC signals on their paths to the respective PBSs.Similar to the PRNR offset the SBR offset could be measured on ground.N sb ij denotes the appearance of the modulation noise in the SBR: This is a combination of left and right-handed modulation noise, their rms amplitudes are 2.9 × 10 −5 m and 2.9 ×10 −4 m, respectively.As shown in [6], it is possible to combine carrier and sideband beatnotes from the RFI to form measurements of the dominating right-handed modulation noise, which can, thus, be subtracted from the SBRs (see appendix B).The advancement operator A sb pbs ← dec (see eq. 27) is associated to the delay operator D sb dec ← pbs , to which we refer as sideband timestamping delay.The sideband timestamping delay can be decomposed into: these constituents are marked green in fig. 4. As for the PRNR timestamping delay, we propose to individually calibrate all its constituents on ground.The sideband timestamping delay can then be compensated in an initial data treatment by application of its associated advancement operator (see eq. 27).
In reality, the beatnotes are expected to be delivered not in phase, but in frequency with occasional phase anchor points.Therefore, we consider the derivative of eq. 28, we refer to it as sideband range rate ( ṠBR): The sideband range rates are an offset-free and unambiguous measurement of the pseudorange time derivatives.Phase anchor points enable their integration, so that we recover eq.28.
C. Time-delay interferometric ranging (TDIR) TDI builds combinations of delayed ISI and RFI carrier beatnotes to virtually form equal-arm interferometers, in which laser frequency noise is suppressed.In the alternative TDI topology, the corresponding TDI delays D ij are given by the pseudoranges in combination with on-board delays (see eq. 13).Time-delay interferometric ranging (TDIR) turns the main scientific data streams themselves into an absolute ranging observable: it minimizes the power integral of the laser frequency noise in the TDI combinations by varying the delays that are applied to the beatnotes [13].When doing this before clock synchronization to TCB, i.e., with the beatnotes sampled according to the respective SCETs, the TDI delays D ij show up at the very minimum of that integral.Thus, TDIR constitutes an unbiased observable of the TDI delays D ij , which only requires the interferometric measurements.
Below, we consider TDI in frequency [14].Therefore, we introduce the Doppler-delay operator associated to the TDI delay: Here we assume the on-board delay constituents of D τi ij to be slowly time-varying, so that only the pseudorange time derivative appears in the Doppler factor.We use the shorthand notation to indicate chained Doppler-delay operators.In this paper we neglect on-board delays in the RFI beatnotes.We start our consideration of TDIR from the intermediary TDI variables η ij .These are combinations of the ISI and RFI carrier beatnotes to eliminate the laser frequency noise contributions of right-handed lasers.In terms of the η ij the second-generation TDI Michelson variables can be expressed as [20] Y τ2 2 (τ ) and Z τ3 2 (τ ) are obtained by cyclic permutation of the indices.For later reference, we also state the first generation TDI Michelson variables: In the framework of TDIR, the delays applied in TDI are parameterized by a model, e.g., by a polynomial model.We minimize the power integral of the TDI combinations by varying the model parameters.TDIR attempts to minimize the in-band laser frequency noise residual.Therefore, we apply a band-pass filter to first remove other contributions appearing out-of-band, i.e., slow drifts and contributions above 1Hz that are dominated by aliasing and interpolation errors.We use the TDIR estimator similar to the one proposed by [13]. 3Θ denotes the parameters of the delay model, the tilde indicates the filtered TDI combinations.
The TDIR accuracy, we denote it by σ tdir , increases with the integration time T (length of telemetry dataset).It is in the order of [13]: where d stands for day.Hence, we require about 1000s of data to reach meter accuracy.

D. Ground-observation based ranging (GOR)
The mission operation center (MOC) provides orbit determinations (ODs) via the ESA tracking stations and MOC time correlations (MOC-TCs).When combined properly, these two on-ground measurements form a pseudorange observable referred to as ground-observation based ranging (GOR).It has an uncertainty of about 50 km due to uncertainties in both the OD and the MOC-TC.Yet, it yields valuable information.It is unambiguous, hence it allows to resolve the PRNR ambiguities.
The OD yields information about the absolute positions and velocities of the three SC.New orbit determinations are published every few days.For the position and velocity measurements in the line of sight, radial (with respect to the sun) and cross-track direction conservative estimations by ESA state the uncertainties as 2 km and 4 mm s −1 , 10 km and 4 mm s −1 , 50 km and 5 cm s −1 , respectively [21].The MOC-TC is a measurement of the SCET desynchronization from TCB.It is determined during the telemetry contacts via a comparison of the SCET associated to the emission of a telemetry packet and the TCB of its reception on Earth taking into account the down link delay.We expect the accuracy of the MOC-TC to be better than 0.1 ms (corresponds to 30 km).This uncertainty is due to unexact knowledge of the SC-to-ground-station separation, as well as inaccuracies in the time tagging process on board and on ground.
As shown in appendix A, the pseudorange can be expressed in TCB as function of the reception time: d t ij denotes the light travel time from SC j to SC i, δτ t ij the offset between the involved SCETs, and δ τ t j the SCET drift of the emitting SC with respect to TCB.The light travel times can be expressed in terms of the ODs [22]: ⃗ r i denotes the position of the receiving SC, ⃗ r j and ⃗ v j are position and velocity of the emitting one.The O(c −3 ) terms contribute to the light travel time at the order of 10 m.They are negligible compared to the large uncertainties of the orbit determination.Combining these light travel time estimates with the MOC-TCs allows us to write the GOR as δτ t tc, i denotes the MOC-TC of SC i and N gor ∼ 50 km the GOR uncertainty.Note that the ODs and the MOC-TCs (hence, also the GOR) are given in TCB, while all other pseudorange observables are sampled in the respective SCETs.This desynchronization is negligible: the desynchronization can amount up to 10 s after the ten year mission time, the pseudoranges drift with 10 to 100 m s −1 (see central plot in fig.6).Hence, neglecting the desynchronization leads to errors of about 100 to 1000 m, which are negligible compared to the large GOR uncertainty.

IV. RANGING SENSOR FUSION
To combine the four pseudorange observables, we propose a three-stage ranging sensor fusion (RSF) consisting of an initial data treatment, a ranging processing, and crosschecks.The ranging processing (central part of fig.5) refers to the ranging-related routines, which need to run continuously during operation.These are the PRNR unwrapping, and the reduction of ranging and right-handed modulation noise.Simultaneously, the PRNR ambiguities and offsets are continuously crosschecked using TDIR and GOR (lower part of fig.5).Both ranging processing and crosschecks rely on a preceding initial data treatment (upper part of fig.5), in which the various delays and offsets are compensated for.Ranging processing and crosschecks can be categorized into four parts demonstrated below: PRNR ambiguity, noise reduction, PRNR offset, and SBR ambiguity.
The RSF delivers accurate and precise pseudorange estimates.To put the RSF into the context of LISA data processing we revisit fig.1: in the baseline TDI topology, the pseudorange estimates from the RSF are disentangled into light travel times and differential timer offsets (see eq. A5).The differential timer offsets are used to synchronize the interferometric measurements, the light travel times serve as delays in TDI.In the alternative topology TDI is executed on the unsynchronized beatnotes.Here the pseudorange estimates from the RSF are directly used as delays in TDI, after they have been combined with on-board delays according to eq. 13.

A. PRNR ambiguity
As part of the ranging processing, the PRNR needs to be steadily unwrapped: due to the finite PRN code length, the PRNR jumps back to 0 km when crossing 400 km and vice versa (see upper plot in fig.6).These jumps are unphysical but easy to identify and to remove.Apart from that, the PRNR ambiguities need to be crosschecked regularly.For that purpose we propose two independent methods below.
GOR represents an unambiguous pseudorange observable.Hence, the combination of PRNR and GOR enables an identification of the PRNR ambiguity integers a prn ij : 400 km is the value we assumed for the PRN code length.However, this procedure only succeeds if |N gor ij | does not exceed the PRN code's half length, i.e., 200 km.Otherwise, a wrong value for the associated PRN ambiguity integer is selected resulting in an estimation error of 400 km in the corresponding link.Note that GOR t ij (t) and PRNR τi ij (τ ) are sampled according to different time frames, but this desynchronization is negligible considering the low accuracy required here (see section III D).
TDIR constitutes an unambiguous pseudorange observable too.Hence, it can be applied as an independent crosscheck of the PRNR ambiguities.We linearly detrend the ISI, RFI, and TMI beatnotes.We then form  Raw datasets are drawn yellow, after the initial data treatment we add a green frame.In the central part we show the core ranging processing.Its output, the pseudoranges, are drawn with a blue frame.In the right box we show how the pseudoranges are combined with on-board delays to form the TDI delays.In the lower part we show simultaneous crosschecks of PRNR ambiguity, PRNR offset, and SBR ambiguity.Products of these crosschecks are drawn with a red frame.
the first-generation TDI Michelson variables (see eq. 36) assuming a constant model for the delays.The pseudoranges are actually drifting by 10 to 100 m s −1 mainly due to differential USO frequency offsets (see central plot in fig.6).Therefore, we choose a short integration time (we use 150 s), otherwise the constant delay model is not sufficient.We use the GOR estimates from above as initial delay values in the TDIR estimator.The TDIR estimator then estimates the six pseudoranges as constants, which can be used to resolve the PRNR ambiguity integers: It is not necessary to apply second-generation TDI, the first-generation already accomplishes the task (see fig. 9).

B. Noise reduction
For the ranging noise reduction in the ranging processing, we propose to combine PRNR and sideband range rates in a modified version of a linearized Kalman filter (KF).The conventional KF requires all measurements to be sampled according to one universal time grid.However, in LISA each SC involves its own SCET.We circumvent this difficulty by splitting up the system and build one KF per SC.Each KF only processes the measurements taken on its associated SC, so that the individual SCETs serve as time-grids.
The state vector of the KF belonging to SC 1 and its associated linear system model can be expressed as k being a discrete time index.Eq.48 describes the time evolution of the state vector from k to k + 1. w τ1 k denotes the process noise vector.In our implementation its covariance matrix W is set W = diag 0, 0, 0, 0, δ k, l denotes the Kronecker delta.The measurement vec-tor and the associated observation model are given by Eq. 52 relates the measurement vector to the state vector.v τ1 k denotes the measurement noise vector.In our implementation its covariance matrix V is set The diagonal entries denote the variances of the respective measurements.We assume the measurements to be uncorrelated, so that the off-diagonal terms are zero.The KFs for SC 2 and SC 3 are defined accordingly.Hence, we remove the ranging noise and obtain estimates for all the six pseudoranges and their time derivatives.These pseudorange estimates are dominated by the right-handed modulation noise, which is one order of magnitude higher than the left-handed one.As pointed out in [6], the right-handed modulation noise can be subtracted (see appendix B): we combine the RFI measurements to form the ∆M i , which are measurements of the right-handed modulation noise on SC i (see eq. B4).For right-handed MOSAs, the local right-handed modulation noise enters the sideband range rates and we just need to subtract the local ∆M i (see eq. B5b).For lefthanded MOSAs the Doppler-delayed right-handed modulation noise from the distant SC appears in the sideband range rates.Here we need to apply the Kalman filter estimates for the pseudoranges and their time derivatives to form the Doppler-delayed distant ∆M i , which then can be subtracted (see eq. B5a).We then process the three KFs again, this time with the corrected sideband range rates.Now they are limited by left-handed modulation noise, so that the respective noise levels are lower.Therefore, we need to adjust the measurement noise covariance matrix for the second run of the KFs: In this way we obtain estimates for the pseudoranges and their time derivatives, which are limited by the lefthanded modulation noise.

C. PRNR offset
The PRNR offset is calibrated on ground before mission start.During operation, it is constructed with the help of SC monitors and subtracted in the initial data treatment.
TDIR can be used as a crosscheck for residual PRNR offsets, as it is sensitive to offsets in the delays.To obtain optimal performance we choose the second-generation TDI Michelson variables (see eq. 35) to be ultimately limited by secondary noises.In-band clock noise is sufficiently suppressed, since we operate on beatnotes in total frequency and make use of the in-band ranging information provided by the preceding noise reduction step.Accordingly, the offset delay model is parameterized by Rτi ij denote the pseudorange estimates after noise reduction, O ij are the 6 offset parameters.As discussed in section III C, computing TDI in total frequency units generally results in a variable with residual trends.Those trends need to be removed prior to calculation of the TDIR integral to be sensitive to residual laser noise in band.This is achieved by an appropriate band-pass filter with a pass-band from 0.1 Hz to 1 Hz.The TDIR integral then reads where tilde indicates the filtered quantities.Ôij are the estimated offset parameters.

D. SBR ambiguity
Phase anchor points in combination with the pseudorange estimates after noise reduction enable the SBR ambiguity resolution (see eq. 28): SBR τi ij are the phase anchor points, Rτi ij the pseudorange estimates of the core ranging processing pipeline.Thus, we obtain estimates of the SBR ambiguity integers a sb ij .These can be used to compute unambiguous SBR pseudorange estimates associated to the phase anchor points, which serve as initial values for the integration of the sideband range rates (eq.32).This procedure is successful if the Rτi ij are more accurate than 6.25 cm (half the SBR ambiguity).Then, SBR constitutes a very accurate pseudorange observable, as not only its precision but also its accuracy are limited by the modulation noise, in contrast to the Rτi ij , whose accuracy is ultimately limited by the ranging noise.We assume here that the SBR offset (eq.29) is corrected in the initial data treatment.

V. RESULTS
In this section, we demonstrate the performance of our implementation of the core ranging processing and the crosschecks as proposed in section IV (central and lower part of fig.5).We did not implement the initial data treatment.Instead we assume that the PRNR and sideband timestamping delays are compensated beforehand.We further consider offset-free PRNR and apply TDIR as a crosscheck for residual offsets.
We use telemetry data simulated by LISA Instrument [23] and LISANode [24] based on orbits provided by ESA [21,25].We simulate phase anchor points for the SBR (see eq. 28).The SCET deviations from the respective proper times are modeled as the δτ i, 0 denote the initial SCET deviations set to 1 s, −1.2 s, and 0.6 s for SC 1, 2, and 3, respectively.The y i model the PMC frequency offsets corresponding to linear clock drifts.They are set to 10 −7 , −2 × 10 −7 , and 0.6 × 10 −7 for SC 1, 2, and 3, respectively.ẏi ∼ 10 −14 s −1 and ÿi ∼ 10 −23 s −2 are constants modeling the linear and quadratic PMC frequency drifts.The y ϵ i denote the stochastic clock noise in fractional frequency deviations, the associated ASD is given by S y ϵ (f ) = 6.32 × 10 −14 Hz −0.5 f Hz −0.5 .
We simulate laser frequency noise with an ASD of and ranging and modulation noise as specified in the sections III A and III B. Furthermore, we consider test-mass acceleration noise and readout noise where A = 6.35 × 10 −12 m Hz −0.5 for the ISI carrier and A = 1.25×10 −11 m Hz −0.5 for the ISI sideband beatnotes.
For the readout noise we set a saturation frequency of f sat = 0.1 mHz, below which we whiten.The orbit determinations are simulated by LISA Ground Tracking with the noise levels specified in section III D.

A. Ranging processing
Here we demonstrate the performance of our implementation of the core ranging processing for one day of telemetry data simulated by LISA Instrument [23].The first ranging processing step covers the PRNR unwrapping (see fig. 6).The upper plot shows the raw PRNR, which jumps back to 0 km when crossing 400 km and vice versa.These jumps are easy to identify and to remove.In our implementation we remove all PRNR jumps bigger than 200 km.The central plot shows the unwrapped but ambiguous PRNR.Here we see PRNR drifts of the order of 10 to 100 m s −1 , which are mainly due to differential USO frequency offsets.Inserting the PRNR ambiguity integers obtained from GOR and TDIR yields the unambiguous PRNR shown in the lower plot.
In the second step, we use the Kalman filter presented in section IV to reduce the ranging noise.Subsequently, we subtract the right-handed modulation noise applying the ∆M measurements constructed from the RFI beatnotes (see appendix B).After noise reduction, we resolve the SBR ambiguities combining the estimated pseudoranges with the simulated SBR phase anchor points (see eq. 58).We then integrate the sideband range rates, to obtain unambiguous SBR.
In fig. 7, we plot the ASDs of the residual pseudorange estimates (deviations of the estimates from the true pseudorange values in the simulation) for link 12 (upper plot) and link 21 (lower plot).Blue lines show the ASDs of the residual PRNR, which are essentially the ASDs of the white ranging noises.The residual pseudorange estimates after ranging noise reduction are plotted in orange.They are obtained by combining the PRNR with the sideband range rates.Therefore, they are limited by right-handed modulation noise (dashed black line).In green, we plot the residual pseudorange estimates after subtraction of right-handed modulation noise with the RFI beatnotes.Now the estimates are limited by lefthanded modulation noise (dash-dotted black line).The residual SBR are drawn red, they are limited by lefthanded modulation noise as well, but involve a smaller offset, since the SBR phase anchor points are more accurate than PRNR after ranging noise reduction (see fig. 8).In the case of left-handed MOSAs (see link 12) the RFI beatnotes need to be time shifted to form the delayed ∆M measurements.We apply the time shifting method of PyTDI [26], which consists in a Lagrange interpolation (we use order 5).The interpolation introduces noise in the high frequency band (see the bump at 2 Hz in the upper plot) but this is out of band.Fig. 8 shows the different residual pseudorange estimates as time series.The upper plot shows the 6 residual pseudorange estimates after ranging noise reduction, the second plot after subtraction of right-handed modulation noise.The third plot shows the SBR residuals.The subtraction of right-handed modulation noise reduces the noise floor, but it does not increase the accuracy of the pseudorange estimates.The accuracy can be increased by one order of magnitude through the resolution of the SBR ambiguities.After ambiguity resolution, SBR constitutes pseudorange estimates with sub-mm accuracy.

B. Crosschecks
Here we demonstrate the performance of our implementation of the crosschecks for PRNR ambiguity and PRNR offset.
The PRNR ambiguities can be resolved using either GOR (see eq. 45) or TDIR (see eq. 46).To evaluate the performance of both methods, we simulate 1000 short (150 s) telemetry datasets with LISA Instrument [23], and one set of ODs and MOC-TCs for each of them.We compute the GOR and TDIR pseudorange estimates for each of the 1000 datasets.Fig. 9 shows the GOR residuals (first row) and the TDIR residuals (second row) in km as histogram plots.We see that the GOR accuracy depends on the arm, because we obtain more accurate ODs for arms oriented in line of sight direction than for those oriented cross-track.The PRNR ambiguity resolution via GOR is successful for GOR deviations smaller than 200 km.In the case of the links 23, 31, 13, and 32 all PRNR ambiguity resolutions via GOR are successful.For each of the links 12 and 21, 2 out of the 1000 PRNR ambiguity resolutions fail.The GOR estimates are passed as initial values to TDIR, which then reduces the uncertainty by almost one order of magnitude (lower plot in fig.9), such that eventually all PRNR ambiguity resolutions are successful.The direction dependent offsets we observe in the TDIR estimates are due to the fact that we apply constant delays in the model used in the minimization of eq.57.In reality, these delays drift with time mainly due to the differential USO frequency offsets, and these drifts are opposite for counterpropagating links (central plot in fig.6).The offsets are not a problem, since the corresponding estimates are one order of magnitude more accurate than what we need.
TDIR can also be applied to estimate the PRNR offsets.Hence, it constitutes a cross-check of the on-ground PRNR offset calibration.We simulate year of teleme-try data using LISANode [24].We set the PRNR offsets to 160.3 m, −210.2 m, 137.3 m, −250.3 m, −188.8 m, and 105.1 m for the links 12, 23, 31, 13, 32, and 21, respectively.We divide the dataset into 1 day chunks (left plots in fig.10), 2 day chunks (central plots in fig.10), and 3 day chunks (right plots in fig.10).In each partition we apply the TDIR estimator presented in section IV C to each chunk in order to estimate the PRNR offsets.This computation was parallelized and executed on the AT-LAS cluster at the AEI Hannover.In the upper part of fig. 10 we show the offset estimation residuals for the three chunk sizes.The offset estimation accuracy increases with the chunk size in agreement with the order of magnitude estimate through eq.38.In the lower part of fig. 10 we plot the residual cumulative averages of the PRNR offset estimates for the different chunk sizes.Here, it can be seen that the TDIR estimator performs similarly for the different chunk sizes.With the 3 day chunk size we can estimate all PRNR offsets with an accuracy of better than 20 cm after 10 days.The dashed-black lines indicate 6.25 cm (half the SBR ambiguity).This is the required PRNR offset estimation accuracy for a successful SBR ambiguity resolution.All offset estimation residuals are below these 6.25 cm after 179 days.

VI. CONCLUSION
The on-board ranging system PRNR requires three treatments due to its ambiguity, offset, and noise.In this article we propose a ranging sensor fusion (RSF), which uses three further observables to solve these issues in order to obtain accurate and precise pseudorange estimates.We show that both GOR and TDIR enable the resolution of the PRNR ambiguity.We investigate how on-board delays affect the pseudorange observables, and propose an initial data treatment to compensate their effects (offsets, timestamping delays) based on an on-ground calibration of the on-board delays.We implement TDIR as a crosscheck for the PRNR offset calibration.We use a Kalman filter to remove the white ranging noise combining PRNR and sideband range rates and apply the RFI beatnotes to subtract the right-handed modulation noise.This leads to pseudorange estimates at sub-cm accuracy.We show that in combination with phase anchor points they allow to resolve the SBR ambiguities resulting in pseudorange estimates at sub-mm accuracy.Apart from that, we identify the delays that are to be applied in TDI in consideration of on-board delays.These are the pseudoranges in combination with on-board delays on both the emitter and the receiver side (see eq. 13).
TDI requires a ranging accuracy of about 10 m to suppress the laser frequency noise below the secondary noise levels [27].While the main building blocks of the here presented RSF are necessary, the reached sub-mm ranging accuracy does not translate into a higher sensitivity for gravitational waves in the final output.Nevertheless, a better ranging accuracy is beneficial: the sec-  ondary noise levels might decrease during the decade until launch; in the case that the cavities perform below expectations, i.e., the encountered laser frequency noise is higher than expected, TDI needs a better ranging accuracy; also in the context of next generation space-based gravitational wave missions a better ranging accuracy might be advantageous.
From the economic perspective, it could be argued whether the on-board ranging system could be dropped in favor of TDIR.Hence, it must be emphasized that PRNR is the most reliable LISA pseudorange observable.After ambiguity resolution and offset correction PRNR delivers a sub-m ranging accuracy directly with the first sample.The TDIR estimator needs more than 1000 s of telemetry data to reach a similar accuracy (see eq. 38).Hence, considering the possibility of data gaps, it is convenient to have PRNR as a direct pseudorange observable.Furthermore, PRNR is more robust than TDIR, which relies on the science data laser frequency noise being its signal.Not only the secondary noise sources but also e.g., gravitational wave signals, are noise to TDIR and degrade its performance.Finally, the PRN modulation has since long been a part of the LISA baseline architecture, as it serves another indisputable function: data transfer between the SC at 45 to 60 kbits s −1 .This requires some kind of digital phase modulation anyway.Implementing the ranging function on top imposes only uncritical extra constraints on that modulation, i.e., synchronization of the digital codes to the SCET.
While we carefully investigate the effects of on-board delays in the ISI, we neglect them in the RFI and in the TMI.However, any differential delay between the laser frequency noise terms appearing in the ISI and in the RFI will cause residual laser frequency noise.Therefore, a follow-up study should investigate the effects of onboard delays for the full LISA interferometeric system.
In the PRNR offset estimation via TDIR we assume these offsets to be constant.In reality, however, they are expected to be slowly time-varying.In a follow-up study, the PRNR offset estimation via TDIR should be extended to linearly time varying PRNR offsets.The delay model for the TDIR estimator would then become: The TDIR estimator would now have to fit the 6 constants O 0 ij and the 6 drifts O 1 ij .Tone-assisted TDIR [28] could be included in this study to reach a better accuracy.
Both, the Kalman filter algorithm presented in section IV B and the TDIR estimator presented in sec- tion IV C are pseudorange estimators.A follow-up study could investigate whether their combination in a single estimator leads to better pseudorange estimates [29].Furthermore, it would be interesting to include timevarying on-board delays and the associated SC monitors into the simulation.This would enable an inspection of the feasibility of the initial data treatment as proposed in section IV.
Finally, the RSF could be included into the different INReP topologies.Apart from that, the algorithms could be applied to real data as, e.g., produced by the hexagon experiment [30], [31].

Figure 1 .
Figure 1.In the baseline TDI topology (upper part) we perform TDI after clock synchronization to TCB, the delays are given by the light travel times.In the alternative TDI topology (lower part) we execute TDI without clock synchronization and apply the pseudoranges as delays[6].Both topologies rely on a ranging sensor fusion.

Figure 2 .
Figure2.LISA labeling conventions (from[14]).The SC are labeled clockwise.The MOSAs and their associated building blocks (lasers, interferometers, etc.) are labeled by 2 indices: the first one indicates the SC they are located at, the second one the SC they are oriented to.The measurements and related quantities (optical links, pseudoranges, etc.) share the indices of the MOSAs they are measured at.Below, we distinguish between left-handed MOSAs(12,23,31) and righthanded ones(13,32,21).

Figure 4 .
Figure 4. We revisit fig. 3. Again we trace the lasers 12 and 21 to the ISI BSs, where they interfere and form beatnotes.The carriers are phase-modulated with the GHz clock and the PRN signals (follow the arrows from the PMC and the PRN to the EOM).Blue arrows indicate analog signals, black arrows digital ones.We show the USO frequency distribution (follow the arrows after the USOs) and illustrate the on-board signal processing (follow the arrows after the QPRs).We highlight the timestamping delays and offsets in the pseudorange observables caused by on-board delays.In light blue, we mark the PRNR offset from the pseudorange.The PRNR timestamping delay is drawn pink, the sideband timestamping delay is marked green.

Figure 5 .
Figure 5.We illustrate the three-stage ranging sensor fusion.Processing elements are drawn with a black frame.In the upper part we show the initial data treatment.Products of the on-ground calibration (the various delays and offsets) are drawn green.Raw datasets are drawn yellow, after the initial data treatment we add a green frame.In the central part we show the core ranging processing.Its output, the pseudoranges, are drawn with a blue frame.In the right box we show how the pseudoranges are combined with on-board delays to form the TDI delays.In the lower part we show simultaneous crosschecks of PRNR ambiguity, PRNR offset, and SBR ambiguity.Products of these crosschecks are drawn with a red frame.

Figure 6 .
Figure 6.Upper plot: raw PRNR.The ambiguity jumps at 0 km and 400 km can be seen.Central plot: ambiguous PRNR, the jumps have been removed but the PRNR ambiguities have not been resolved yet.The large slopes are mainly due to USO frequency offsets.Lower plot: unambiguous PRNR.The large differences between the links are caused by differential SCET offsets.

Figure 7 .
Figure 7. ASDs of the residual pseudorange estimates for link 12 (upper plot) and link 21 (lower plot).In blue, residual PRNR.In orange, residual pseudorange estimates after ranging noise reduction.In green, residual pseudorange estimates after subtraction of right-handed modulation noise.In red, residual SBR.Dashed black lines, right-handed modulation noise model.Dash-dotted black lines, left-handed modulation noise model.

Figure 8 .
Figure 8. Upper plot: residual pseudorange estimates after ranging noise reduction.Second plot: residual pseudorange estimates after subtraction of right-handed modulation noise.Third plot: residual SBR estimates.

Figure 9 .
Figure 9. PRNR ambiguity resolution via GOR (upper plots) and TDIR (lower plots).The histogram plots show the residual GOR and TDIR pseudorange estimates (deviations from true pseudorange) for the different links.

Figure 10 .
Figure 10.We simulate one year of telemetry data with PRNR offsets in the order of 100 m.We divide this dataset into 1 day chunks (left plots), 2 day chunks (central plots), and 3 day chunks (right plots).We then we apply TDIR to each of these chunks in order to estimate the PRNR offsets.Upper plots: Residual offset estimates in m for the different chunk sizes.Lower plots: residual offset estimates after cumulative averaging in m for the different chunk sizes.Dashed-black lines: half the SBR ambiguity.