1   /* Copyright 2002-2026 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  package org.orekit.estimation;
18  
19  import org.hipparchus.geometry.euclidean.threed.Vector3D;
20  import org.hipparchus.util.FastMath;
21  import org.hipparchus.util.Pair;
22  import org.orekit.bodies.CelestialBody;
23  import org.orekit.bodies.GeodeticPoint;
24  import org.orekit.bodies.OneAxisEllipsoid;
25  import org.orekit.estimation.measurements.GroundStation;
26  import org.orekit.estimation.measurements.ObserverSatellite;
27  import org.orekit.forces.drag.DragSensitive;
28  import org.orekit.forces.gravity.potential.NormalizedSphericalHarmonicsProvider;
29  import org.orekit.forces.gravity.potential.UnnormalizedSphericalHarmonicsProvider;
30  import org.orekit.forces.radiation.RadiationSensitive;
31  import org.orekit.frames.TopocentricFrame;
32  import org.orekit.models.earth.displacement.StationDisplacement;
33  import org.orekit.orbits.*;
34  import org.orekit.propagation.PropagationType;
35  import org.orekit.propagation.analytical.BrouwerLyddanePropagator;
36  import org.orekit.propagation.analytical.tle.TLE;
37  import org.orekit.propagation.analytical.tle.generation.FixedPointTleGenerationAlgorithm;
38  import org.orekit.propagation.conversion.*;
39  import org.orekit.time.TimeScale;
40  import org.orekit.time.UT1Scale;
41  import org.orekit.utils.IERSConventions;
42  import org.orekit.utils.PVCoordinates;
43  
44  import java.util.List;
45  import java.util.Map;
46  
47  public class Context implements StationDataProvider {
48  
49      /** IERS conventions. */
50      public IERSConventions conventions;
51  
52      /** Central body shape. */
53      public OneAxisEllipsoid earth;
54  
55      /** Sun model. */
56      public CelestialBody sun;
57  
58      /** Moon model. */
59      public CelestialBody moon;
60  
61      /** Radiation pressure model. */
62      public RadiationSensitive radiationSensitive;
63  
64      /** Drag model. */
65      public DragSensitive dragSensitive;
66  
67      /** Normalized spherical harmonics coefficients. */
68      public NormalizedSphericalHarmonicsProvider normalizedProvider;
69  
70      /** Unnormalized spherical harmonics coefficients. */
71      public UnnormalizedSphericalHarmonicsProvider unnormalizedProvider;
72  
73      /** UTC. */
74      public TimeScale utc;
75  
76      /** UT1. */
77      public UT1Scale ut1;
78  
79      /** Initial orbit. */
80      public Orbit initialOrbit;
81  
82      /** Initial TLE for SGP4 orbit determination. */
83      public TLE initialTLE;
84  
85      /** Reference points displacement models. */
86      public StationDisplacement[] displacements;
87  
88      /** Ground stations for orbit determination. */
89      public List<GroundStation> stations;
90  
91      /** Stations for turn-around range.
92       * Map entry = primary station
93       * Map value = secondary station associated
94       */
95      public List<ObserverSatellite> satellites;
96  
97      /** Stations for turn-around range.
98       * Map entry = primary station
99       * Map value = secondary station associated
100      */
101     public Map<GroundStation, GroundStation> TARstations;
102 
103     /** Stations for bi-static range rate.
104      * Map entry = primary station
105      * Map value = secondary station associated
106      */
107     public Pair<GroundStation, GroundStation> BRRstations;
108 
109     /** Stations for TDOA.
110      * Map entry = primary station
111      * Map value = secondary station associated
112      */
113     public Pair<GroundStation, GroundStation> TDOAstations;
114 
115     /** Stations for FDOA.
116      * Map entry = primary station
117      * Map value = secondary station associated
118      */
119     public Pair<GroundStation, GroundStation> FDOAstations;
120 
121     /**
122      * Creates a DSST propagator builder.
123      * <p>
124      *  By default propagation type and initial state type are defined as {@link PropagationType#MEAN MEAN}.
125      * </p>
126      * @param perfectStart if false, orbit estimation will start from a wrong point
127      * @param minStep the minimum step size for numerical integration
128      * @param maxStep the maximum step size for numerical integration
129      * @param dP position error
130      * @param forces the set of force models to include in the numerical propagator
131      * @return a configured DSSTPropagatorBuilder instance
132      */
133     public DSSTPropagatorBuilder createDsst(final boolean perfectStart,
134                                             final double minStep, final double maxStep, final double dP,
135                                             final DSSTForce... forces) {
136         return createDsst(PropagationType.MEAN, PropagationType.MEAN, perfectStart, minStep, maxStep, dP, forces);
137     }
138 
139     /**
140      * Creates a DSST propagator builder.
141      *
142      * @param propagationType type of the orbit used for the propagation (mean or osculating)
143      * @param stateType  type of the elements used to define the initial orbital state (mean or osculating)
144      * @param perfectStart if false, orbit estimation will start from a wrong point
145      * @param minStep the minimum step size for numerical integration
146      * @param maxStep the maximum step size for numerical integration
147      * @param dP position error
148      * @param forces the set of force models to include in the numerical propagator
149      * @return a configured DSSTPropagatorBuilder instance
150      */
151     public DSSTPropagatorBuilder createDsst(final PropagationType propagationType,
152                                             final PropagationType stateType, final boolean perfectStart,
153                                             final double minStep, final double maxStep, final double dP,
154                                             final DSSTForce... forces) {
155         // Initialize builder
156         final DSSTPropagatorBuilder propagatorBuilder =
157                 new DSSTPropagatorBuilder(createInitialOrbit(perfectStart),
158                                           new DormandPrince853IntegratorBuilder(minStep, maxStep, dP),
159                                           dP,  propagationType, stateType);
160         // Add force models
161         for (DSSTForce force : forces) {
162             propagatorBuilder.addForceModel(force.getForceModel(this));
163         }
164         // Return the configured builder
165         return propagatorBuilder;
166     }
167 
168     /**
169      * Creates a numerical propagator builder.
170      *
171      * @param orbitType the type of orbit to be used by the builder
172      * @param positionAngleType the position angle type to be used by the builder
173      * @param perfectStart if false, orbit estimation will start from a wrong point
174      * @param minStep the minimum step size for numerical integration
175      * @param maxStep the maximum step size for numerical integration
176      * @param dP position error
177      * @param forces the set of force models to include in the numerical propagator
178      * @return a configured NumericalPropagatorBuilder instance
179      */
180     public NumericalPropagatorBuilder createNumerical(final OrbitType orbitType, final PositionAngleType positionAngleType,
181                                                       final boolean perfectStart,
182                                                       final double minStep, final double maxStep, final double dP,
183                                                       final Force... forces) {
184         // Initialize builder
185         final NumericalPropagatorBuilder propagatorBuilder =
186                 new NumericalPropagatorBuilder(orbitType.convertType(createInitialOrbit(perfectStart)),
187                                                new DormandPrince853IntegratorBuilder(minStep, maxStep, dP),
188                                                positionAngleType, dP);
189         // Add force models
190         for (Force force : forces) {
191             propagatorBuilder.addForceModel(force.getForceModel(this));
192         }
193         // Return the configured builder
194         return propagatorBuilder;
195     }
196 
197     /**
198      * Creates a keplerian propagator builder.
199      *
200      * @param angleType the position angle type to be used by the builder
201      * @param perfectStart if false, orbit estimation will start from a wrong point
202      * @param dP scaling factor used for orbital parameters normalization
203      * @return a configured KeplerianPropagatorBuilder instance
204      */
205     public KeplerianPropagatorBuilder createKeplerian(final PositionAngleType angleType, final boolean perfectStart, final double dP) {
206         // Return the configured builder
207         return new KeplerianPropagatorBuilder(createInitialOrbit(perfectStart), angleType, dP);
208     }
209 
210     /**
211      * Creates a Eckstein-Hechler propagator builder.
212      *
213      * @param angleType the position angle type to be used by the builder
214      * @param perfectStart if false, orbit estimation will start from a wrong point
215      * @param dP scaling factor used for orbital parameters normalization
216      * @return a configured EcksteinHechlerPropagatorBuilder instance
217      */
218     public EcksteinHechlerPropagatorBuilder createEcksteinHechler(final PositionAngleType angleType, final boolean perfectStart, final double dP) {
219         // Return the configured builder
220         return new EcksteinHechlerPropagatorBuilder(createInitialOrbit(perfectStart), unnormalizedProvider, angleType, dP);
221     }
222 
223     /**
224      * Creates a Brouwer-Lyddane propagator builder.
225      *
226      * @param angleType the position angle type to be used by the builder
227      * @param perfectStart if false, orbit estimation will start from a wrong point
228      * @param dP scaling factor used for orbital parameters normalization
229      * @return a configured BrouwerLyddanePropagatorBuilder instance
230      */
231     public BrouwerLyddanePropagatorBuilder createBrouwerLyddane(final PositionAngleType angleType, final boolean perfectStart, final double dP) {
232         // Return the configured builder
233         return new BrouwerLyddanePropagatorBuilder(createInitialOrbit(perfectStart), unnormalizedProvider, angleType, dP, BrouwerLyddanePropagator.M2);
234     }
235 
236     /**
237      * Creates a TLE propagator builder.
238      *
239      * @param dP scaling factor used for orbital parameters normalization
240      * @return a configured TLEPropagatorBuilder instance
241      */
242     public TLEPropagatorBuilder createTleBuilder(final double dP) {
243         return new TLEPropagatorBuilder(initialTLE, PositionAngleType.MEAN, dP, new FixedPointTleGenerationAlgorithm());
244     }
245 
246     /**
247      * Creates the initial orbit.
248      *
249      * @param perfectStart if false, orbit estimation will start from a wrong point
250      * @return the initial orbit.
251      */
252     public Orbit createInitialOrbit(boolean perfectStart) {
253         final Orbit startOrbit;
254         if (perfectStart) {
255             // orbit estimation will start from a perfect orbit
256             startOrbit = initialOrbit;
257         } else {
258             // orbit estimation will start from a wrong point
259             final Vector3D initialPosition = initialOrbit.getPosition();
260             final Vector3D initialVelocity = initialOrbit.getVelocity();
261             final Vector3D wrongPosition   = initialPosition.add(new Vector3D(1000.0, 0, 0));
262             final Vector3D wrongVelocity   = initialVelocity.add(new Vector3D(0, 0, 0.01));
263             startOrbit                     = new CartesianOrbit(new PVCoordinates(wrongPosition, wrongVelocity),
264                                                                 initialOrbit.getFrame(),
265                                                                 initialOrbit.getDate(),
266                                                                 initialOrbit.getMu());
267         }
268         return startOrbit;
269     }
270 
271     /**
272      * Creates a ground station given its geographic coordinates and name.
273      *
274      * @param latitudeInDegrees  the latitude of the station in degrees
275      * @param longitudeInDegrees the longitude of the station in degrees
276      * @param altitude           the altitude of the station in meters
277      * @param name               the name of the station
278      * @return a new instance of {@code GroundStation} created with the provided parameters
279      */
280     public GroundStation createStation(double latitudeInDegrees, double longitudeInDegrees,
281                                 double altitude, String name) {
282         final GeodeticPoint gp = new GeodeticPoint(FastMath.toRadians(latitudeInDegrees),
283                                                    FastMath.toRadians(longitudeInDegrees),
284                                                    altitude);
285         return new GroundStation(new TopocentricFrame(earth, gp, name), ut1.getEOPHistory(), displacements);
286     }
287     
288     ObserverSatellite createObserverSatellite(final Vector3D deltaPosition, final Vector3D deltaVelocity,
289             final String name, final Force... forces) {
290 
291         final Vector3D initialPosition = initialOrbit.getPosition();
292         final Vector3D initialVelocity = initialOrbit.getVelocity();
293         final Vector3D finalPosition = initialPosition.add(deltaPosition);
294         final Vector3D finalVelocity = initialVelocity.add(deltaVelocity);
295 
296         final Orbit orbit = new CartesianOrbit(new PVCoordinates(finalPosition, finalVelocity),
297                 initialOrbit.getFrame(),
298                 initialOrbit.getDate(),
299                 initialOrbit.getMu());
300 
301         final NumericalPropagatorBuilder propagatorBuilder = new NumericalPropagatorBuilder(orbit,
302                 new DormandPrince853IntegratorBuilder(1.0, 5.0, 10.0),
303                 PositionAngleType.TRUE, 10.0);
304         for (Force force : forces) {
305             propagatorBuilder.addForceModel(force.getForceModel(this));
306         }
307 
308         return new ObserverSatellite(name, propagatorBuilder.buildPropagator());
309     }
310 
311     @Override
312     public List<GroundStation> getStations() {
313         return stations;
314     }
315     
316     @Override
317     public List<ObserverSatellite> getSatellites() {
318         return satellites;
319     }
320 
321     @Override
322     public OneAxisEllipsoid getEarth() {
323         return earth;
324     }
325 
326 }