GNSS-based precise orbit determination for maneuvering LEO satellites

Maneuverability is essential for low Earth orbit (LEO) satellites to fulfill various operational objectives. However, the precise orbit determination (POD) process might deteriorate due to imperfect satellite orbital dynamics modeling. This article develops a generic POD strategy with maneuver handling for LEO satellites equipped with high-performance spaceborne Global Navigation Satellite System (GNSS) receivers. Given the time span of an executed maneuver, a set of constant thrust accelerations in the satellite body-fixed reference frame is estimated without using a-priori maneuver accelerations. In addition, different numbers of velocity pulses are estimated at predefined epochs determined by the duration of a maneuver. POD experiments are done for the GRACE-FO and Sentinel-3 satellites, for which the orbit maneuvers vary significantly. The orbits are assessed via internal consistency checks and external orbit validations. Internally, in each direction, the agreement between the reduced-dynamic and kinematic orbits reaches a level of 1 cm, which is comparable with the reference day without maneuvers. Externally, comparisons with the official GRACE-FO products and orbits from the Sentinel-3 Copernicus POD Quality Working Group confirm the reliability of the new orbits with maneuver handling. Finally, satellite laser ranging and K-band ranging measurements indicate a 1-cm accuracy of the absolute orbits and a 2-mm accuracy of the GRACE-FO relative orbits. The maneuver handling strategy is tested in the Bernese GNSS Software, consistently developed at the Astronomical Institute of the University of Bern.


Introduction
Over the past decades, there has been an unprecedented increase in using low Earth orbit (LEO) satellites for various scientific research objectives, which usually necessitate the best-possible precise orbit determination (POD) products. Now an orbit accuracy of 1-2 cm is realistic in the routine operations of a few POD software packages using the highlow satellite-to-satellite tracking (hl-SST) data collected by the spaceborne Global Navigation Satellite System (GNSS) receivers (Fernández et al. 2022). Further improvements are possible by applying refined gravitational and non-gravitational force modeling, single-receiver integer ambiguity resolution (SR-IAR) and multi-GNSS observations (e.g., Mao et al. 2021;Montenbruck et al. 2021;Peter et al. 2022). Nevertheless, obtaining POD solutions for LEO satellites undergoing maneuvers is still challenging. In reality, orbit maneuvers are regularly executed to avoid threatening collisions, overcome orbit decaying, maintain predefined trajectories, etc.
If not properly handled, the executed orbit maneuvers might significantly deteriorate the (reduced-) dynamic POD solution since it relies heavily on the satellite orbital dynamics. A classical solution is introducing the a-priori maneuver accelerations and durations to the satellite equation-of-motion (EoM), which is numerically integrated to generate a continuous orbit. However, the planned maneuver information the flight dynamics and operations team provides is usually imperfect due to unpredictable causes such as in-flight software control, hardware delays and thruster performances. Further corrections are usually required when commanding a real maneuver but foreseeable errors w.r.t. the in-flight performances still exist (Yoon et al. 2006 The satellite radio telemetry, which in principle provides more accurate information about the performed maneuvers, is unfortunately not easily accessible for public use. The GNSS-based POD strategy thus needs to be improved based on the publicly available maneuver information. On the one hand, the orbit maneuver accelerations can be re-estimated or corrected as constant thrust accelerations over the entire maneuver duration in the local orbital reference frame (LORF), i.e., the classical radial, along-track and cross-track reference system Ju et al. 2017). Additional parameters such as maneuver execution delays can be estimated and adjusted for special cases, e.g., Allende-Alba et al. (2017) shift the maneuver start and/or end epochs of the GRACE twin satellites by about 0.04-0.9 s accordingly. Besides, satellite attitude control might be required to align the orbit control thrusters to a desired direction before executing a maneuver, e.g., large yaw-turns are necessary before and after an out-of-plane (cross-track) maneuver. However, the attitude control commands do not necessarily perform perfectly as planned and execution errors retain, e.g., up to 3 degrees in yaw control for the Sentinel-3 satellites (Taboada et al. 2018). Even during a strong maneuver, the thrust might perturb the attitude control stability, causing errors in the satellite EoM integration process.
On the other hand, an orbit maneuver can be modeled as instantaneous changes (the so-called impulsive maneuver) to the satellite velocities in LORF. The velocity pulses are treated as accumulated effects of the executed maneuvers. Jäggi et al. (2012) indicate that such a strategy works reliable for satellites experiencing short maneuvers (with durations up to a few seconds). For satellites performing long maneuvers, a delicate selection of epochs for estimating velocity pulses would be required (Springer 2009). Moreover, the estimates of velocity pulses are heavily determined by the quality and quantity of GNSS observations near the requested epochs. In case of poor tracking arcs or even data gaps, the POD performance will deteriorate. The estimates also need to be properly constrained by the predefined variances, which cannot universally apply to maneuver durations ranging from 1 s to 30 min. This strategy needs to be carefully tailored to prevent abnormal estimates of velocity pulses.
This article aims to develop a generic POD strategy for the different maneuvering LEO satellites using the Bernese GNSS Software (BSW), which has been continuously developed at the Astronomical Institute of the University of Bern (AIUB) (Dach et al. 2015 The remainder of this article will first discuss the generic maneuver handling strategy, shows the results of internal obit consistency checks and external orbit validations and draws conclusions at the end.

Maneuver handling strategy
The nominal reduced-dynamic POD (RPOD) strategy of LEO satellites implemented in BSW needs to be modified to consider satellite maneuvers. Representations of the satellite status in the Earth-centered inertial coordinate reference frame (ECIF) are given by, (1) r(t) = r t; a; e; i; Ω; ; u 0 ; Q 1,…,d ; P 1,…,s ; a M ; M where, r, a are the instantaneous position and acceleration vectors of the motion of the satellite's center of mass, respectively; a, e, i, Ω, , u 0 denote six Keplerian osculating orbital elements at a reference epoch t 0 ; Q represents d dynamic orbital parameters that are set as 3 constant accelerations over an orbit arc; P represents s groups of pseudo-stochastic parameters to compensate for modeling deficiencies of the orbital dynamics, e.g., is set as 240 groups of evenly spaced piecewise constant accelerations for a 24-h orbit arc (Jäggi 2007); GM is the gravity constant times the mass of the Earth; f represents the sum of gravitational and non-gravitational forces beyond the centered term of the Earth's gravity, and a satellite macro-model is used (Mao et al. 2021); a M and M are the vectors of maneuver accelerations and velocity pulses, respectively.
Instead of estimating maneuver accelerations in LORF (Allende- Ju et al. 2017), the proposed strategy selects the satellite body-fixed reference frame (SBRF), which is found to be crucial for satellites experiencing strong maneuvers, e.g., for S3A/B. An example out-of-plane maneuver is described in Fig. 1, which shows the yaw-turn to align the thrust vector to the designated direction during the maneuver. The thrust seems to perturb the attitude stability by roughly 0.5° after the execution epoch, which requires an attitude correction immediately. Estimating maneuver accelerations in LORF is thus not accurate in this case because the follow-up orbit integration will be influenced. As a comparison, similar attitude perturbation caused by the executed GFC/D out-of-plane maneuvers is less than 0.1°. Besides, the S3A/B orbit control thrust accelerations are larger by more than an order of magnitude than for GFC/D. In this article, maneuver accelerations are only estimated in SBRF, which is always aligned with the executed thrusts.
For GFC/D, the official spacecraft-event and sequence-ofevent files (TN-01a-SCE.txt and TN-01-SOE.txt from ftp:// isdcf tp. gfz-potsd am. de/ grace-fo/ DOCUM ENTS/ TECHN ICAL_ NOTES/) do not necessarily provide velocity changes for all the performed orbit maneuvers. For S3A/B, the planed maneuvers provided by GMV often display different accelerations at the start and end epochs of, especially the strong maneuvers. There will be further deviations w.r.t the real-commanded maneuvers in terms of both accelerations and durations (Jaime Fernández, private communication). Therefore, the proposed strategy does not introduce a-priori maneuver accelerations but instead estimates accelerations from scratch using the a-prior maneuver durations. Although the largest part of an executed maneuver can be estimated as constant thrust accelerations, empirical corrections are still necessary to overcome errors originating from maneuver planning, execution delays and/or thruster firing status. As a result, additional velocity pulses M are estimated. Depending on the duration of a maneuver, either 1 velocity pulse located at the end epoch of a maneuver (if shorter than 3 times of a nominal integration interval, i.e., 1 min in this article) or 3 velocity pulses located at the start, middle and end epochs of a maneuver (if longer), are estimated. Note that more velocity pulses do not necessarily further improve the orbit solution, particularly for days with short maneuvers, because the estimates are highly determined by the availability of GNSS observations near the requested epochs.
Assuming an a-priori orbit r 0 (t) is available, e.g., obtained from a dynamical fit of a coarse kinematic orbit, the true orbit r(t) can be expanded into the first order Taylor series w.r.t. the n unknown orbit parameters p: in BSW, orbit solutions are solved based on variational equations as explicitly explained by Jäggi (2007). The partial derivative of the sum of perturbing forces w.r.t. each parameter p i needs to be computed. Partial derivatives for the maneuver accelerations, which are only applied to the satellite dynamics during maneuvers, are written as, where, A denotes the attitude matrix between SBRF and ECIF. It is computed either from quaternions collected by star trackers or from a nominal attitude law only in case of missing quaternions. If multiple maneuvers are performed in one orbit arc (24 h), the same maneuver accelerations are estimated under the assumption that the thruster executions in SBRF are consistent in short time. In addition, the velocity pulses constrained by a-priori variances can be estimated at requested epochs in the parameter estimation process (Jäggi 2007). An estimated velocity pulse M,i is then applied as M,i • t − t i ⋅ e i (t) to the satellite EoM ( denotes the Dirac delta function and e contains the respective unit vector).
The orbit integrator uses a collocation method. Each integration interval is first assigned then subdivisions are determined by the order of the polynomial approximation of the satellite EoM (Beutler 2004). However, the start and end epochs of a maneuver usually do not exactly match with the predefined integration boundaries, resulting into mismodeling of maneuver accelerations (Fig. 2). As a solution, the two integration interval boundaries that are closest to the maneuver start and end epochs are shifted to fully cover the maneuver duration.
Three series of RPOD orbits (referred to as the solutions AF, BF and CF) are generated, allowing for direct assessments of the impact of the different implementations individually: the solution AF only models maneuvers as a set of constant thrust accelerations; the solution BF builds up on AF by estimating additional velocity pulses with a weak constraint of = 1m∕s ; the solution CF further applies SR-IAR and currently yields the best orbit. The observable-specific bias (OSB) product provided by the Center for Orbit Determination in Europe (CODE) is used (Schaer et al. 2021). Considering the fact that for S3A/B strong maneuvers may change the satellite velocity of up to 2.3 m/s (Taboada et al. 2018), in this article ambiguities crossing maneuvers are only estimated as float values to get rid of wrongly fixed integers near or during the maneuvers. Roughly 2% of all possible ambiguities are not fixed.

Results and discussion
This section shows internal orbit quality checks, orbit validations using external sources and further discusses modeling a maneuver using multiple different accelerations.

Internal consistency assessment
For assessing the RPOD orbits, the kinematic POD (KPOD) orbit solution, which heavily relies on the availability and quality of GNSS observations and orbit/clock products, is also generated. The GNSS receivers onboard the selected satellites all show stable tracking performance and hardly any data gaps occur, no further processing is required (Wang et al. 2020). The KPOD uses the same batch of GNSS carrier-phase observations screened in the RPOD strategy. The percentage of epochs with kinematic coordinates is first checked. GFC (99.50%), GFD (98.73%), S3A (99.96%) and S3B (99.79%) are all obtain percentages close to 100%. The performances for GFC/D are slightly worse than S3A/B due to more GNSS tracking gaps. When checking the days with strong maneuvers (i.e., > 300 s) for S3A (the solution CF), the percentage of maneuver durations is 1.03%, whereas the missing kinematic epochs only account for 0.08%, suggesting that the modeling of maneuvers does not affect GNSS data screening at a visible level. The kinematic positions during maneuvers thus allow for the evaluation of the RPOD orbits.
Comparisons between the KPOD and RPOD orbits are displayed in Table 1, and a time series analysis is shown in Fig. 3. Only estimating maneuvers as constant thrust accelerations (the solution AF) for GFC/D results in orbit consistencies of 1-2 cm (in each direction), which is similar to the POD performances for GRACE (Ju et al. 2017). The additional estimates of velocity pulses (the solution BF), which are at a level of 0.1 mm/s, marginally improve the orbits. Finally, SR-IAR (the solution CF), with success rates of 92.4/91.8% for GFC/D, leads to orbit consistencies of subcentimeter in particularly the along-track and cross-track directions (Montenbruck et al. 2018). Note that such orbit improvement starts from day 18/173 because the operational OSB products only became available from mid-2018. Figure 3 suggests that the solution AF for S3A is not satisfactory. For days (20/169, 20/337) with smaller maneuvers, the orbit consistency reaches a similar level for GFC/GFD, even a similar level for the reference day (20/245). However, for days (20/071, 20/246, 20/351) with stronger maneuvers, the adjacent hours around the maneuvers are heavily influenced such that the RPOD orbits are significantly degraded. It is interesting to note that the solution BF, with estimated velocity pulses at a level of around 1 cm/s, results in much better orbit consistencies for those problematic days. As anticipated, SR-IAR (the solution CF) further improves the orbit consistencies to 1.3/0.7/0.8 cm for S3A. Due to a higher orbit (GFC/D: 500 km, S3A/B: 800 km) that experiences less atmospheric disturbances, the SR-IAR success rates for S3A/B are higher (96.6/96.9%). The numbers are slightly smaller than the 99.0/99.2% reported by Mao et al. (2021), who analyzed a period of approximately 1.5 years but only for days without maneuvers. This is due to the fact that the GNSS carrier-phase ambiguities crossing maneuvers are not resolved to integer numbers.

Orbit cross-validations
The agreement with external orbits computed by other software packages is a widely acknowledged approach to evaluating the POD performances (CPOD team 2022). The KPOD orbits, which are independent of satellite dynamics, are in the following used to evaluate the different orbits. For GFC/D, the Jet Propulsion Laboratory (JPL) has been consistently generating RPOD orbits since the launch. Figure 4 indicates that the agreements between the reference KPOD orbit and the two RPOD orbits reach a level of 1-2 cm in each direction, no large influences caused by maneuvers are seen. Due to the fact of using the same batch of GNSS carrier-phase and integer ambiguities, the consistency between the two BSW orbits is slightly better than that of the JPL orbits. Note from now on the ultimate reduceddynamic orbits with SR-IAR (the solution CF) are used for orbit comparisons.
In the framework of the CPOD QWG, S3A/B orbits from different member institutes are routinely generated aiming to assess the official orbits (Fernández et al. 2022). Most of the QWG member institutes and external institutes, e.g., JPL, deliver orbits for the days with maneuvers. In addition, a combined product is routinely generated by averaging the 0.0 ± 1.9 0.0 ± 1.9 0.0 ± 1.6 − 0.4 ± 1.5 − 0.4 ± 1.5 − 0.0 ± 0.9 − 0.4 ± 1.4 − 0.2 ± 1.2 0.1 ± 0.6 S3A − 0.2 ± 7.0 − 0.0 ± 1.9 − 0.0 ± 1.3 − 0.0 ± 11.5 − 0.0 ± 1.6 − 0.0 ± 0.7 − 1.4 ± 15.6 − 0.2 ± 1.2 0.0 ± 0.8 S3B − 0.8 ± 9.5 − 0.0 ± 1.8 0.1 ± 1.1 0.2 ± 19.2 − 0.0 ± 1.3 − 0.1 ± 0.5 − 3.9 ± 23.1 − 0.1 ± 1.4 0.0 ± 0.7  Table 2 lists the orbits that are generated by different software packages, which model or calibrate satellite maneuvers as impulsive maneuvers (velocity pulses) and/or continuous accelerations. It must be noted that the delivered orbits are generated by using default processing strategies that are not necessarily fine-tuned to optimally cope with satellites experiencing strong maneuvers. Figure 5 shows orbit comparisons between the BSW KPOD orbit and other orbits. The agreements between the BSW KPOD and the CPOK orbits (a kinematic orbit with SR-IAR computed by NAPEOS) are at levels of 3.6/1.8/2.7 (S3A) and 1.7/0.8/6.9 (S3B) cm under the fact that no outlier screening is applied, showing reliable quality of the two KPOD orbits. In general, all reduced-dynamic orbits show good agreements with the BSW KPOD orbit when containing small maneuvers. Nevertheless, a few orbits are degraded by up to a few meters when containing strong maneuvers. The impact may last a few hours near the maneuvers. Some orbits are not influenced for hours but instead show large discrepancies close to maneuvers. Although occasionally appearing, those discrepancies are found to be even larger than 10 m. Finally, the BSW RPOD orbit seems to show the best orbit agreements at a level of 0.5-1.3 cm in each direction. On the one hand, the BSW RPOD orbit benefits from the proposed maneuver handling strategy, showing less discrepancy near maneuvers; on the other hand, using the same screened GNSS observations leads to better internal consistency.

Agreement with a-priori maneuvers
As discussed above, the accurate satellite radio telemetry information is not publicly available. However, the flight dynamics and operations team will plan a scheme of maneuver operations depending on the satellite status and orbit prediction (Kahle et al. 2012). It can be compared with the estimated maneuver accelerations as a cross-validation, as displayed in Fig. 6. For the days with long GFC/D maneuvers (> 100 s), the relative difference is within 11.0%. When checking the agreement for days with shorter maneuvers, a larger discrepancy is shown. For instance, on day 20/253 GFD only performed a small maneuver of 5 s but a difference of 51.7% occurs. Note that the GRACE-FO sequenceof-event file declares all the maneuver execution starting epochs with an offset of 0.6 s to integer seconds, and they  always last for integer numbers of seconds. In fact, these epochs will not be perfectly followed considering the diverse delays (Allende-Alba et al. 2017). Therefore, it is not surprising to see larger discrepancies when handling shorter maneuvers. The same conclusions can be drawn for S3A/B. For day 20/211 when S3B experienced a 1.9-s maneuver, a difference of 76.8% occurred. Otherwise, it is interesting to see that for the days with S3A/B maneuvers longer than 15 s, the relative difference is only less than 0.6%. Although it is impossible to decide which set of accelerations is more accurate without checking the radio telemetry, the orbit agreements show a reliable quality (Fig. 5).

Satellite laser ranging
The SLR observations offer an opportunity to assess orbit accuracy independently. Note that given the days with maneuver plans, some SLR stations will lower the observation priority of a satellite. However, analysis in this article shows that S3B had been well tracked during all of its maneuver days in 2020 (no SLR tracking exists exactly during the maneuvers). A list of ten high-performance SLR stations ) are used in the validation, and an elevation cut off angle of 10° is applied. Figure 7 illustrates that the STDs of SLR residuals are 9.8/7.6 mm for the BSW KPOD/RPOD orbits, showing a similar orbit quality for days without maneuvers (Mao et al. 2021). Again, orbits from some institutes deteriorate due to maneuver handling that influences orbit parameter estimates even during the maneuver-free parts of an orbit arc. Besides, a special SLR tracking pass crossing a 774.3-s maneuver is found for S3B on day 21/034, when S3B was flying above the Tahiti SLR station and 15 normal points were successfully processed centering at the maneuver execution epoch. Figure 8 shows that the BSW orbits are hardly affected by the maneuver as the SLR residuals remain stable at levels of 1 cm.

K-band radar ranging
The GRACE-FO mission offers a unique possibility to compare the baseline derived from the two independent GNSSbased orbits with the directly measured distance through KBR, as illustrated in Fig. 9. Note that KBR measurements are unavailable for the early days since launch and for orbits near large maneuvers, e.g., out-of-plane shifting. POD for day 20/245, without any maneuvers (as used for S3A), is performed as reference. The official orbits, with maneuver handling and SR-IAR using the JPL GNSS orbit/clock products, obtain a consistency level of 3.3 mm. The consistency for the BSW RPOD orbit is at a level of 2.0 mm, probably because of the different methodologies employed for maneuver handling and SR-IAR. Please note days 20/245 and 20/253 show slightly larger discrepancy because the GNSS receivers seem to suffer from the activation of the flex power mode IV on the GPS Block-IIR-M and -IIF satellites since Feb 24, 2020 (Steigenberger et al. 2020

Test on estimating multiple accelerations per maneuver
The current maneuver handling strategy estimates only one set of constant thrust accelerations to model the maneuver accelerations, whereas the modeling deficiency is compensated by the additional velocity pulses. It is worth checking if an increased set of estimated acceleration parameters will be better due to a larger degree of freedom for maneuver modeling (Oliver Montenbruck, private communication). An S3A POD test for day 20/246 was done by Oliver Montenbruck using GHOST, where the three subarcs and maneuver execution delays were intently tuned given the a-prior maneuver duration. The new DLR orbit agrees with the same BSW KPOD orbit at a level of only 1.0/0.4/0.9 cm, which is significantly better than the current operational DLR orbit as reported in Fig. 5. A similar strategy, which defines two identical subarcs around the middle epoch of a maneuver and estimates additional velocity pulses, is tested in BSW. For day 20/246, the orbit consistency reaches a level of 0.9/0.4/0.5 cm (the solution DF in Fig. 10), which is better than the 1.0/0.4/0.8 cm as obtained by the solution CF. The improvement is mostly due to orbit parts near the executed maneuver. Unfortunately, such an extended maneuver handling requires a case-to-case fine-tuning. Besides, it cannot be flexibly extended to process days with multiple maneuvers and/or short maneuvers. The current BSW maneuver handling strategy stays with the solution CF, which in principle can be more easily conducted to generically fit into different situations. Fig. 6 Comparisons between the a-priori (for S3A/B, take mean acceleration; for GFC/D, the velocity changes are transformed into maneuver accelerations after being divided by the maneuver dura-tions) and the estimated maneuver accelerations for the selected satellites. The differences are described as percentages listed below the plots Fig. 7 Residuals of SLR validations (933 measurements) of the different orbits for S3B. AIKN is the BSW KPOD orbits Fig. 8 Residuals of a single SLR tracking pass collected by the Tahiti station for S3B on day 21/034. Note that one orbit is outside of the figure and the ESOC orbit is unavailable

Conclusions
A POD strategy for maneuvering LEO satellites is developed and assessed using the Bernese GNSS Software (BSW). Given the a-priori maneuver time, a set of constant thrust accelerations defined in the satellite body-fixed reference frame is estimated over the maneuver durations.
To compensate for maneuver modeling deficiency, one or three sets of additional velocity pulses are estimated at requested epochs determined by the duration of a maneuver. The orbit integration steps should also be adjusted to fully cover a maneuver. POD investigations are done with the Sentinel-3 and GRACE-FO satellites, which experience significantly different maneuvers and carry a set of independent tools for POD and orbit validations. Based on the same batch of GNSS carrier-phase observations screened with maneuver handling, kinematic and reduceddynamic orbits are generated and compared. GRACE-FO POD tests show that a satisfactory orbit can be already obtained when only modeling maneuvers as a set of constant thrust accelerations, mostly because GRACE-FO maneuver accelerations are much smaller than those of the Sentinel-3 satellites. The latter will benefit significantly from the estimated velocity pulses, particularly for days with strong maneuvers. Therefore, the generic maneuver handling strategy is designed as a combination of estimating constant thrust accelerations and velocity pulses simultaneously. With a further step of single-receiver integer ambiguity resolution at a success rate of over 91% for all satellites, the obtained reduced-dynamic orbits show consistencies of 1-2 cm w.r.t. the kinematic orbits in each direction.
Orbit cross-validations are made by comparing the BSW orbits with external orbits, which are generated using different operational software packages that might suffer from imperfect modeling of maneuvers. When referring to the same BSW kinematic orbit, the consistencies for the BSW reduced-dynamic orbit seem to be the best. In addition, the estimated maneuver accelerations are compared with the publicly available maneuver information. For days with strong maneuvers, the agreements reach levels of less than 0.6% and 11% for the Sentinel-3 and GRACE-FO satellites, respectively. However, larger discrepancies are seen for days with small maneuvers. Validations by satellite laser ranging and K-band radar ranging measurements reveal that a 1-cm level absolute orbit accuracy is achieved for the Sentinel-3B satellite and a 2-mm level relative orbit accuracy is achieved for the GRACE-FO satellites. Finally, a preliminary test for the Sentinel-3A satellite shows that further orbit improvement is possible with the estimation of multiple accelerations per maneuver. Instead of case-to-case fine-tuning, future work can be done to make this strategy more robustly implemented. More satellite maneuver detection and validation techniques can be explored, e.g., checking the real telemetry data and the SLR full-rate data for a specified tracking pass (Dilssner et al. 2020).
Acknowledgements This work was partially funded by the European Organisation for the Exploitation of Meteorological Satellites (EUMETSAT) through contract EUM/CO/21/4600002516/CMa. We acknowledge Dr. Christian Marquardt, Dr. Francisco Sancho and Dr. Yago Andres for their input and feedback on this project. The generous support of the anonymous reviewers and editors is also gratefully acknowledged.  Consistencies between the same BSW KPOD reference orbit and the two RPOD orbits for S3A on day 20/246. Only orbit parts near maneuver are shown. Statistics is for the entire orbit arc Data availability The GRACE-FO data are obtained from the German Research Centre for Geosciences (ftp:// isdcf tp. gfz-potsd am. de/ grace-fo/) and the Jet Propulsion Laboratory (https:// podaac. jpl. nasa. gov/ GRACE-FO). The Sentinel-3 data are provided by the ESA's Copernicus Precise Orbit Determination service (https:// cpode xtftp. gmv. com/ COPPOD/ DATA). GNSS ephemeris, clock and OSB products are obtained from the Center for Orbit Determination in Europe (http:// ftp. aiub. unibe. ch/ CODE/). Last accessed date to all the mentioned weblinks is Sep 1, 2022.

Conflict of interest
The authors declare no competing interests.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http:// creat iveco mmons. org/ licen ses/ by/4. 0/.