1   /* Copyright 2002-2025 CS GROUP
2    * Licensed to CS GROUP (CS) under one or more
3    * contributor license agreements.  See the NOTICE file distributed with
4    * this work for additional information regarding copyright ownership.
5    * CS licenses this file to You under the Apache License, Version 2.0
6    * (the "License"); you may not use this file except in compliance with
7    * the License.  You may obtain a copy of the License at
8    *
9    *   http://www.apache.org/licenses/LICENSE-2.0
10   *
11   * Unless required by applicable law or agreed to in writing, software
12   * distributed under the License is distributed on an "AS IS" BASIS,
13   * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14   * See the License for the specific language governing permissions and
15   * limitations under the License.
16   */
17  
18  package org.orekit.estimation.measurements.modifiers;
19  
20  import org.hipparchus.util.FastMath;
21  import org.orekit.estimation.measurements.EstimatedMeasurementBase;
22  import org.orekit.frames.Frame;
23  import org.orekit.orbits.KeplerianOrbit;
24  import org.orekit.propagation.SpacecraftState;
25  import org.orekit.utils.Constants;
26  import org.orekit.utils.TimeStampedPVCoordinates;
27  
28  /**
29   * Class modifying theoretical measurements with relativistic J2 clock correction.
30   * <p>
31   * Relativistic clock correction of the effects caused by the oblateness of Earth on
32   * the gravity potential.
33   * </p>
34   * <p>
35   * The time delay caused by this effect is computed based on the orbital parameters of the
36   * emitter's orbit.
37   * </p>
38   *
39   * @author Louis Aucouturier
40   * @since 11.2
41   *
42   * @see "Teunissen, Peter, and Oliver Montenbruck, eds. Springer handbook of global navigation
43   * satellite systems. Chapter 19.2. Equation 19.18 Springer, 2017."
44   */
45  public class AbstractRelativisticJ2ClockModifier {
46  
47      /**
48       * Relativistic J2 effect constant.
49       */
50      private final double cJ2;
51  
52      /** Central attraction coefficient. */
53      private final double gm;
54  
55      /**
56       * Constructor for the Relativistic J2 Clock modifier.
57       * @param gm Earth gravitational constant (mu) in m³/s².
58       * @param c20 Earth un-normalized second zonal coefficient (Signed J2 constant, is negative) (Typical value -1.0826e-3).
59       * @param equatorialRadius Earth equatorial radius in m.
60       */
61      public AbstractRelativisticJ2ClockModifier(final double gm,
62                                                 final double c20,
63                                                 final double equatorialRadius) {
64          this.cJ2 = 1.5 * c20 * equatorialRadius * equatorialRadius /
65                  (Constants.SPEED_OF_LIGHT * Constants.SPEED_OF_LIGHT);
66          this.gm = gm;
67      }
68  
69      /** Get the name of the effect modifying the measurement.
70       * @return name of the effect modifying the measurement
71       * @since 13.0
72       */
73      public String getEffectName() {
74          return "J₂ clock relativity";
75      }
76  
77      /**
78       * Computes the relativistic J2 clock time delay correction.
79       *
80       * @param estimated EstimatedMeasurements on which to calculate the correction
81       * @return dt_relJ2clk Time delay due to the relativistic J2 clock effect in seconds
82       */
83      protected double relativisticJ2Correction(final EstimatedMeasurementBase<?> estimated) {
84  
85          // Extracting the state of the receiver to determine the frame and mu
86          //
87          // The satellite states are stored at the creation of the estimated measurements
88          // and can contain up to 2 elements. In most cases, only the receiver's state and
89          // therefore frame is stored, with the emitter's frame corresponding to the receiver's.// Still, in the InterSatellites case, the states of the 2 spacecrafts are stored, // and can contain different frames. This case is treated by looking at the length
90          // of SpacecraftState stored in the Estimated Measurements, with the only length 2
91          // case is the InterSatellites case.
92          //
93          final SpacecraftState[] states = estimated.getStates();
94          final SpacecraftState state =  (states.length < 2) ? states[0] : states[1];
95  
96          final Frame remoteFrame = state.getFrame();
97  
98          // Getting Participants to extract the remote PV
99          final TimeStampedPVCoordinates[] pvs = estimated.getParticipants();
100 
101         // Checking if the correction is applied on a two-way GNSS problem
102         // In that case the emitter is at index 1, else index 0
103         final TimeStampedPVCoordinates pvRemote = (pvs.length < 3) ? pvs[0] : pvs[1];
104 
105         // Define a Keplerian orbit to extract the orbital parameters needed to compute the correction
106         final KeplerianOrbit remoteOrbit = new KeplerianOrbit(pvRemote, remoteFrame, gm);
107         final double orbitInclination = remoteOrbit.getI();
108 
109         // u = perigee argument + true anomaly
110         final double orbitU = remoteOrbit.getTrueAnomaly() + remoteOrbit.getPerigeeArgument();
111         final double n = remoteOrbit.getKeplerianMeanMotion();
112 
113         // Returning the value of the time delay
114         return cJ2 * n * FastMath.sin(2 * orbitU) * FastMath.sin(orbitInclination) * FastMath.sin(orbitInclination);
115     }
116 
117 }