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.bodies;
18
19
20 import org.hipparchus.CalculusFieldElement;
21 import org.hipparchus.Field;
22 import org.hipparchus.geometry.euclidean.threed.FieldRotation;
23 import org.hipparchus.geometry.euclidean.threed.FieldVector3D;
24 import org.hipparchus.geometry.euclidean.threed.Rotation;
25 import org.hipparchus.geometry.euclidean.threed.RotationConvention;
26 import org.hipparchus.geometry.euclidean.threed.Vector3D;
27 import org.hipparchus.util.Precision;
28 import org.orekit.bodies.JPLEphemeridesLoader.EphemerisType;
29 import org.orekit.bodies.JPLEphemeridesLoader.ZeroRawPVProvider;
30 import org.orekit.frames.FieldKinematicTransform;
31 import org.orekit.frames.FieldStaticTransform;
32 import org.orekit.frames.FieldTransform;
33 import org.orekit.frames.Frame;
34 import org.orekit.frames.KinematicTransform;
35 import org.orekit.frames.OriginTransformProvider;
36 import org.orekit.frames.Predefined;
37 import org.orekit.frames.StaticTransform;
38 import org.orekit.frames.Transform;
39 import org.orekit.frames.TransformProvider;
40 import org.orekit.time.AbsoluteDate;
41 import org.orekit.time.FieldAbsoluteDate;
42 import org.orekit.time.TimeOffset;
43 import org.orekit.utils.FieldPVCoordinates;
44 import org.orekit.utils.PVCoordinates;
45 import org.orekit.utils.TimeStampedFieldPVCoordinates;
46 import org.orekit.utils.TimeStampedPVCoordinates;
47
48 import java.util.concurrent.TimeUnit;
49
50 /** Implementation of the {@link CelestialBody} interface using JPL or INPOP ephemerides.
51 * @author Luc Maisonobe
52 */
53 class JPLCelestialBody implements CelestialBody {
54
55 /** Name of the body. */
56 private final String name;
57
58 /** Regular expression for supported files names. */
59 private final String supportedNames;
60
61 /** Ephemeris type to generate. */
62 private final JPLEphemeridesLoader.EphemerisType generateType;
63
64 /** Raw position-velocity provider. */
65 private final JPLEphemeridesLoader.RawPVProvider rawPVProvider;
66
67 /** Attraction coefficient of the body (m³/s²). */
68 private final double gm;
69
70 /** Scaling factor for position-velocity. */
71 private final double scale;
72
73 /** IAU pole. */
74 private final IAUPole iauPole;
75
76 /** Body's PV coordinates are defined in this frame. */
77 private final Frame definingFrameAlignedWithIcrf;
78
79 /** Body centered frame aligned with ICRF. */
80 private final Frame icrfAlignedFrame;
81
82 /** Inertially oriented, body-centered frame. */
83 private final Frame inertialFrame;
84
85 /** Body oriented, body-centered frame. */
86 private final Frame bodyFrame;
87
88 /** Build an instance and the underlying frame.
89 * @param name name of the body
90 * @param supportedNames regular expression for supported files names
91 * @param generateType ephemeris type to generate
92 * @param rawPVProvider raw position-velocity provider
93 * @param gm attraction coefficient (in m³/s²)
94 * @param scale scaling factor for position-velocity
95 * @param iauPole IAU pole implementation
96 * @param definingFrameAlignedWithICRF frame in which celestial body coordinates are defined,
97 * this frame <strong>must</strong> be aligned with ICRF
98 * @since 14.0
99 */
100 JPLCelestialBody(final String name, final String supportedNames,
101 final JPLEphemeridesLoader.EphemerisType generateType,
102 final JPLEphemeridesLoader.RawPVProvider rawPVProvider,
103 final double gm, final double scale,
104 final IAUPole iauPole,
105 final Frame definingFrameAlignedWithICRF) {
106 this.name = name;
107 this.gm = gm;
108 this.scale = scale;
109 this.supportedNames = supportedNames;
110 this.generateType = generateType;
111 this.rawPVProvider = rawPVProvider;
112 this.iauPole = iauPole;
113 this.definingFrameAlignedWithIcrf = definingFrameAlignedWithICRF;
114 if (rawPVProvider instanceof ZeroRawPVProvider) {
115 // no translation or rotation needed, use directly
116 // might be better to have a method instead of using "instanceof"
117 // but the classes are tightly coupled and package private
118 this.icrfAlignedFrame = definingFrameAlignedWithICRF;
119 } else {
120 // translation needed
121 final String icrfName;
122 if (EphemerisType.SOLAR_SYSTEM_BARYCENTER == generateType) {
123 // in Orekit FramesFactory.getICRF() is implemented by
124 // CelestialBodyFactor.getSsb().getInertiallyOrientedFrame()
125 // so have to match Predefined.ICRF
126 icrfName = Predefined.ICRF.getName();
127 } else {
128 icrfName = name + "/ICRF";
129 }
130 this.icrfAlignedFrame = new Frame(
131 definingFrameAlignedWithICRF,
132 new OriginTransformProvider(this, definingFrameAlignedWithICRF),
133 icrfName,
134 true);
135 }
136 if (iauPole == null || iauPole.isGcrfAligned()) {
137 // Body "fixed" and inertial frames are GCRF aligned.
138 this.inertialFrame = icrfAlignedFrame;
139 this.bodyFrame = icrfAlignedFrame;
140 } else {
141 this.inertialFrame = new InertiallyOriented(icrfAlignedFrame);
142 this.bodyFrame = new BodyOriented();
143 }
144 }
145
146 /** {@inheritDoc} */
147 @Override
148 public TimeStampedPVCoordinates getPVCoordinates(final AbsoluteDate date, final Frame frame) {
149
150 // apply the scale factor to raw position-velocity
151 final PVCoordinates rawPV = rawPVProvider.getRawPV(date);
152 final TimeStampedPVCoordinates scaledPV = new TimeStampedPVCoordinates(date, scale, rawPV);
153
154 // the raw PV are relative to the parent of the body centered inertially oriented frame
155 final Transform transform = definingFrameAlignedWithIcrf.getTransformTo(frame, date);
156
157 // convert to requested frame
158 return transform.transformPVCoordinates(scaledPV);
159
160 }
161
162 /** Get the {@link FieldPVCoordinates} of the body in the selected frame.
163 * @param date current date
164 * @param frame the frame where to define the position
165 * @param <T> type of the field elements
166 * @return time-stamped position/velocity of the body (m and m/s)
167 */
168 @Override
169 public <T extends CalculusFieldElement<T>> TimeStampedFieldPVCoordinates<T> getPVCoordinates(final FieldAbsoluteDate<T> date,
170 final Frame frame) {
171
172 // apply the scale factor to raw position-velocity
173 final FieldPVCoordinates<T> rawPV = rawPVProvider.getRawPV(date);
174 final TimeStampedFieldPVCoordinates<T> scaledPV = new TimeStampedFieldPVCoordinates<>(date, scale, rawPV);
175
176 // the raw PV are relative to the parent of the body centered inertially oriented frame
177 final FieldTransform<T> transform = definingFrameAlignedWithIcrf.getTransformTo(frame, date);
178
179 // convert to requested frame
180 return transform.transformPVCoordinates(scaledPV);
181
182 }
183
184 /** {@inheritDoc} */
185 @Override
186 public Vector3D getVelocity(final AbsoluteDate date, final Frame frame) {
187
188 // apply the scale factor to raw position-velocity
189 final PVCoordinates rawPV = rawPVProvider.getRawPV(date);
190 final TimeStampedPVCoordinates scaledPV = new TimeStampedPVCoordinates(date, scale, rawPV);
191
192 // the raw PV are relative to the parent of the body centered inertially oriented frame
193 final KinematicTransform transform = definingFrameAlignedWithIcrf.getKinematicTransformTo(frame, date);
194
195 // convert to requested frame
196 return transform.transformOnlyPV(scaledPV).getVelocity();
197
198 }
199
200 /** {@inheritDoc} */
201 @Override
202 public Vector3D getPosition(final AbsoluteDate date, final Frame frame) {
203
204 // apply the scale factor to raw position
205 final Vector3D rawPosition = rawPVProvider.getRawPosition(date);
206 final Vector3D scaledPosition = rawPosition.scalarMultiply(scale);
207
208 // the raw position is relative to the parent of the body centered inertially oriented frame
209 final StaticTransform transform = definingFrameAlignedWithIcrf.getStaticTransformTo(frame, date);
210
211 // convert to requested frame
212 return transform.transformPosition(scaledPosition);
213 }
214
215 /** {@inheritDoc} */
216 @Override
217 public <T extends CalculusFieldElement<T>> FieldVector3D<T> getPosition(final FieldAbsoluteDate<T> date, final Frame frame) {
218
219 // apply the scale factor to raw position
220 final FieldVector3D<T> rawPosition = rawPVProvider.getRawPosition(date);
221 final FieldVector3D<T> scaledPosition = rawPosition.scalarMultiply(scale);
222
223 // the raw position is relative to the parent of the body centered inertially oriented frame
224 final FieldStaticTransform<T> transform = definingFrameAlignedWithIcrf.getStaticTransformTo(frame, date);
225
226 // convert to requested frame
227 return transform.transformPosition(scaledPosition);
228 }
229
230 /** {@inheritDoc} */
231 @Override
232 public <T extends CalculusFieldElement<T>> FieldVector3D<T> getVelocity(final FieldAbsoluteDate<T> date, final Frame frame) {
233 // apply the scale factor to raw position-velocity
234 final FieldPVCoordinates<T> rawPV = rawPVProvider.getRawPV(date);
235 final TimeStampedFieldPVCoordinates<T> scaledPV = new TimeStampedFieldPVCoordinates<>(date, scale, rawPV);
236
237 // the raw PV are relative to the parent of the body centered inertially oriented frame
238 final FieldKinematicTransform<T> transform = definingFrameAlignedWithIcrf.getKinematicTransformTo(frame, date);
239
240 // convert to requested frame
241 return transform.transformOnlyPV(scaledPV).getVelocity();
242 }
243
244 /** {@inheritDoc} */
245 public String getName() {
246 return name;
247 }
248
249 /** {@inheritDoc} */
250 public double getGM() {
251 return gm;
252 }
253
254 @Override
255 public Frame getIcrfAlignedFrame() {
256 return icrfAlignedFrame;
257 }
258
259 /** {@inheritDoc} */
260 public Frame getInertiallyOrientedFrame() {
261 return inertialFrame;
262 }
263
264 /** {@inheritDoc} */
265 public Frame getBodyOrientedFrame() {
266 return bodyFrame;
267 }
268
269 /** Inertially oriented body centered frame. */
270 private class InertiallyOriented extends Frame {
271
272 /** Suffix for inertial frame name. */
273 private static final String INERTIAL_FRAME_SUFFIX = "/inertial";
274
275 /**
276 * Simple constructor.
277 *
278 * @param bodyIcrf body centered ICRF aligned frame.
279 * @since 14.0
280 */
281 InertiallyOriented(final Frame bodyIcrf) {
282 super(bodyIcrf, new TransformProvider() {
283
284 /** {@inheritDoc} */
285 public Transform getTransform(final AbsoluteDate date) {
286
287 // compute rotation from ICRF frame to self,
288 // as per the "Report of the IAU/IAG Working Group on Cartographic
289 // Coordinates and Rotational Elements of the Planets and Satellites"
290 // These definitions are common for all recent versions of this report
291 // published every three years, the precise values of pole direction
292 // and W angle coefficients may vary from publication year as models are
293 // adjusted. These coefficients are not in this class, they are in the
294 // specialized classes that do implement the getPole and getPrimeMeridianAngle
295 // methods
296 final Vector3D pole = iauPole.getPole(date);
297 final Vector3D qNode = iauPole.getNode(date);
298 final Transform rotation =
299 new Transform(date, new Rotation(pole, qNode, Vector3D.PLUS_K, Vector3D.PLUS_I));
300
301 // update transform from parent to self
302 return rotation;
303
304 }
305
306 @Override
307 public StaticTransform getStaticTransform(final AbsoluteDate date) {
308 // compute rotation from ICRF frame to self,
309 // as per the "Report of the IAU/IAG Working Group on Cartographic
310 // Coordinates and Rotational Elements of the Planets and Satellites"
311 // These definitions are common for all recent versions of this report
312 // published every three years, the precise values of pole direction
313 // and W angle coefficients may vary from publication year as models are
314 // adjusted. These coefficients are not in this class, they are in the
315 // specialized classes that do implement the getPole and getPrimeMeridianAngle
316 // methods
317 final Vector3D pole = iauPole.getPole(date);
318 final Vector3D qNode = iauPole.getNode(date);
319 final Rotation rotation =
320 new Rotation(pole, qNode, Vector3D.PLUS_K, Vector3D.PLUS_I);
321
322 // update transform from parent to self
323 return StaticTransform.of(date, rotation);
324 }
325
326 /** {@inheritDoc} */
327 public <T extends CalculusFieldElement<T>> FieldTransform<T> getTransform(final FieldAbsoluteDate<T> date) {
328
329 // compute rotation from ICRF frame to self,
330 // as per the "Report of the IAU/IAG Working Group on Cartographic
331 // Coordinates and Rotational Elements of the Planets and Satellites"
332 // These definitions are common for all recent versions of this report
333 // published every three years, the precise values of pole direction
334 // and W angle coefficients may vary from publication year as models are
335 // adjusted. These coefficients are not in this class, they are in the
336 // specialized classes that do implement the getPole and getPrimeMeridianAngle
337 // methods
338 final FieldVector3D<T> pole = iauPole.getPole(date);
339 FieldVector3D<T> qNode = FieldVector3D.crossProduct(Vector3D.PLUS_K, pole);
340 if (qNode.getNorm2Sq().getReal() < Precision.SAFE_MIN) {
341 qNode = FieldVector3D.getPlusI(date.getField());
342 }
343 final FieldTransform<T> rotation =
344 new FieldTransform<>(date,
345 new FieldRotation<>(pole,
346 qNode,
347 FieldVector3D.getPlusK(date.getField()),
348 FieldVector3D.getPlusI(date.getField())));
349
350 // update transform from parent to self
351 return rotation;
352
353 }
354
355 @Override
356 public <T extends CalculusFieldElement<T>> FieldStaticTransform<T> getStaticTransform(final FieldAbsoluteDate<T> date) {
357 // field
358 final Field<T> field = date.getField();
359
360 // compute rotation from ICRF frame to self,
361 // as per the "Report of the IAU/IAG Working Group on Cartographic
362 // Coordinates and Rotational Elements of the Planets and Satellites"
363 // These definitions are common for all recent versions of this report
364 // published every three years, the precise values of pole direction
365 // and W angle coefficients may vary from publication year as models are
366 // adjusted. These coefficients are not in this class, they are in the
367 // specialized classes that do implement the getPole and getPrimeMeridianAngle
368 // methods
369 final FieldVector3D<T> pole = iauPole.getPole(date);
370 final FieldVector3D<T> qNode = iauPole.getNode(date);
371 final FieldRotation<T> rotation =
372 new FieldRotation<>(pole, qNode, FieldVector3D.getPlusK(field), FieldVector3D.getPlusI(field));
373
374 // update transform from parent to self
375 return FieldStaticTransform.of(date, rotation);
376 }
377
378 }, name + INERTIAL_FRAME_SUFFIX, true);
379 }
380
381 }
382
383 /** Body oriented body centered frame. */
384 private class BodyOriented extends Frame {
385
386 /**
387 * Suffix for body frame name.
388 */
389 private static final String BODY_FRAME_SUFFIX = "/rotating";
390
391 /**
392 * Simple constructor.
393 *
394 * @since 14.0
395 */
396 BodyOriented() {
397 super(inertialFrame, new TransformProvider() {
398
399 /** {@inheritDoc} */
400 public Transform getTransform(final AbsoluteDate date) {
401 final TimeOffset dt = new TimeOffset(10, TimeUnit.SECONDS);
402 final double w0 = iauPole.getPrimeMeridianAngle(date);
403 final double w1 = iauPole.getPrimeMeridianAngle(date.shiftedBy(dt));
404 return new Transform(date,
405 new Rotation(Vector3D.PLUS_K, w0, RotationConvention.FRAME_TRANSFORM),
406 new Vector3D((w1 - w0) / dt.toDouble(), Vector3D.PLUS_K));
407 }
408
409 /** {@inheritDoc} */
410 public <T extends CalculusFieldElement<T>> FieldTransform<T> getTransform(final FieldAbsoluteDate<T> date) {
411 final TimeOffset dt = new TimeOffset(10, TimeUnit.SECONDS);
412 final T w0 = iauPole.getPrimeMeridianAngle(date);
413 final T w1 = iauPole.getPrimeMeridianAngle(date.shiftedBy(dt));
414 return new FieldTransform<>(date,
415 new FieldRotation<>(FieldVector3D.getPlusK(date.getField()), w0,
416 RotationConvention.FRAME_TRANSFORM),
417 new FieldVector3D<>(
418 w1.subtract(w0).divide(dt.toDouble()),
419 Vector3D.PLUS_K));
420 }
421
422 }, name + BODY_FRAME_SUFFIX, false);
423 }
424 }
425 }