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