1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17 package org.orekit.orbits;
18
19 import java.util.ArrayList;
20 import java.util.List;
21
22 import org.hipparchus.Field;
23 import org.hipparchus.geometry.euclidean.threed.FieldVector3D;
24 import org.hipparchus.util.Binary64;
25 import org.hipparchus.util.Binary64Field;
26 import org.hipparchus.util.FastMath;
27 import org.hipparchus.util.MathArrays;
28 import org.junit.jupiter.api.Assertions;
29 import org.junit.jupiter.api.BeforeEach;
30 import org.junit.jupiter.api.DisplayName;
31 import org.junit.jupiter.api.Test;
32 import org.orekit.Utils;
33 import org.orekit.frames.FramesFactory;
34 import org.orekit.propagation.analytical.FieldEcksteinHechlerPropagator;
35 import org.orekit.time.AbstractTimeInterpolator;
36 import org.orekit.time.FieldAbsoluteDate;
37 import org.orekit.time.FieldTimeInterpolator;
38 import org.orekit.utils.CartesianDerivativesFilter;
39 import org.orekit.utils.FieldPVCoordinates;
40
41 class FieldOrbitHermiteInterpolatorTest {
42
43 private final Field<Binary64> field = Binary64Field.getInstance();
44
45 @BeforeEach
46 public void setUp() {
47 Utils.setDataRoot("regular-data");
48 }
49
50 @Test
51 public void testCartesianInterpolationWithDerivatives() {
52
53 doTestCartesianInterpolation(true, CartesianDerivativesFilter.USE_P,
54 394, 0.1968, 3.21, 0.021630,
55 2474, 2707.6418, 6.6, 26.28);
56
57
58 doTestCartesianInterpolation(true, CartesianDerivativesFilter.USE_PV,
59 394, 1.085E-8, 3.21, 2.41518E-10,
60 2474, 0.06249, 6.6, 0.001170);
61
62
63
64 doTestCartesianInterpolation(true, CartesianDerivativesFilter.USE_PVA,
65 394, 1.92e-8, 3.21, 1.15e-9,
66 2474, 5227, 6.55, 142);
67 }
68
69 @Test
70 public void testCartesianInterpolationWithoutDerivatives() {
71
72 doTestCartesianInterpolation(false, CartesianDerivativesFilter.USE_P,
73 394, 0.1968, 3.21, 0.02163,
74 2474, 2707.6419, 6.55, 26.2826);
75
76
77 doTestCartesianInterpolation(false, CartesianDerivativesFilter.USE_PV,
78 394, 1.0857E-8, 3.21, 2.4151E-10,
79 2474, 0.06249, 6.55, 0.001170);
80
81
82
83 doTestCartesianInterpolation(false, CartesianDerivativesFilter.USE_PVA,
84 394, 2.61, 3.21, 0.154,
85 2474, 2.28e12, 6.55, 6.22e10);
86 }
87
88 @Test
89 public void testCircularInterpolationWithDerivatives() {
90 doTestCircularInterpolation(true,
91 397, 2.53e-8,
92 610, 4.72e-6,
93 4870, 110);
94 }
95
96 @Test
97 public void testCircularInterpolationWithoutDerivatives() {
98 doTestCircularInterpolation(false,
99 397, 0.0372,
100 610.0, 1.23,
101 4870, 8869);
102 }
103
104 @Test
105 public void testEquinoctialInterpolationWithDerivatives() {
106 doTestEquinoctialInterpolation(true,
107 397, 1.28e-8,
108 610, 4.35e-6,
109 4870, 114);
110 }
111
112 @Test
113 public void testEquinoctialInterpolationWithoutDerivatives() {
114 doTestEquinoctialInterpolation(false,
115 397, 0.0372,
116 610.0, 1.23,
117 4879, 8871);
118 }
119
120 @Test
121 public void testKeplerianInterpolationWithDerivatives() {
122 doTestKeplerianInterpolation(true,
123 397, 4.01, 4.75e-4, 1.28e-7,
124 2159, 1.05e7, 1.19e-3, 0.773);
125 }
126
127 @Test
128 public void testKeplerianInterpolationWithoutDerivatives() {
129 doTestKeplerianInterpolation(false,
130 397, 62.0, 4.75e-4, 2.87e-6,
131 2159, 79365, 1.19e-3, 3.89e-3);
132 }
133
134 private void doTestCartesianInterpolation(boolean useDerivatives, CartesianDerivativesFilter pvaFilter,
135 double shiftPositionErrorWithin, double interpolationPositionErrorWithin,
136 double shiftVelocityErrorWithin, double interpolationVelocityErrorWithin,
137 double shiftPositionErrorFarPast, double interpolationPositionErrorFarPast,
138 double shiftVelocityErrorFarPast, double interpolationVelocityErrorFarPast) {
139
140 Binary64 zero = field.getZero();
141 final Binary64 ehMu = zero.add(3.9860047e14);
142
143 final double ae = 6.378137e6;
144 final double c20 = -1.08263e-3;
145 final double c30 = 2.54e-6;
146 final double c40 = 1.62e-6;
147 final double c50 = 2.3e-7;
148 final double c60 = -5.5e-7;
149
150 final FieldAbsoluteDate<Binary64> date = FieldAbsoluteDate.getJ2000Epoch(field).shiftedBy(584.);
151 final FieldVector3D<Binary64> position =
152 new FieldVector3D<>(zero.add(3220103.), zero.add(69623.), zero.add(6449822.));
153 final FieldVector3D<Binary64> velocity =
154 new FieldVector3D<>(zero.add(6414.7), zero.add(-2006.), zero.add(-3180.));
155 final FieldCartesianOrbit<Binary64> initialOrbit =
156 new FieldCartesianOrbit<>(new FieldPVCoordinates<>(position, velocity),
157 FramesFactory.getEME2000(), date, ehMu);
158
159 FieldEcksteinHechlerPropagator<Binary64> propagator =
160 new FieldEcksteinHechlerPropagator<>(initialOrbit, ae, ehMu, c20, c30, c40, c50, c60);
161
162
163 List<FieldOrbit<Binary64>> sample = new ArrayList<>();
164 for (Binary64 dt = zero; dt.getReal() < 251.0; dt = dt.add(60.0)) {
165 FieldOrbit<Binary64> orbit = propagator.propagate(date.shiftedBy(dt)).getOrbit();
166 if (!useDerivatives) {
167
168 Binary64[] stateVector = MathArrays.buildArray(field, 6);
169 orbit.getType().mapOrbitToArray(orbit, PositionAngleType.TRUE, stateVector, null);
170 orbit = orbit.getType().mapArrayToOrbit(stateVector, null, PositionAngleType.TRUE,
171 orbit.getDate(), orbit.getMu(), orbit.getFrame());
172 }
173 sample.add(orbit);
174 }
175
176
177 final FieldTimeInterpolator<FieldOrbit<Binary64>, Binary64> interpolator =
178 new FieldOrbitHermiteInterpolator<>(sample.size(), 409, initialOrbit.getFrame(), pvaFilter);
179
180
181
182
183
184 double maxShiftPError = 0;
185 double maxInterpolationPError = 0;
186 double maxShiftVError = 0;
187 double maxInterpolationVError = 0;
188 for (Binary64 dt = zero; dt.getReal() < 240.0; dt = dt.add(1.0)) {
189 FieldAbsoluteDate<Binary64> t = initialOrbit.getDate().shiftedBy(dt);
190 FieldPVCoordinates<Binary64> propagated = propagator.propagate(t).getPVCoordinates();
191 FieldPVCoordinates<Binary64> shiftError = new FieldPVCoordinates<>(propagated,
192 initialOrbit.shiftedBy(dt.getReal())
193 .getPVCoordinates());
194 FieldPVCoordinates<Binary64> interpolationError = new FieldPVCoordinates<>(propagated,
195 interpolator.interpolate(t, sample)
196 .getPVCoordinates());
197 maxShiftPError = FastMath.max(maxShiftPError,
198 shiftError.getPosition().getNorm().getReal());
199 maxInterpolationPError = FastMath.max(maxInterpolationPError,
200 interpolationError.getPosition().getNorm().getReal());
201 maxShiftVError = FastMath.max(maxShiftVError,
202 shiftError.getVelocity().getNorm().getReal());
203 maxInterpolationVError = FastMath.max(maxInterpolationVError,
204 interpolationError.getVelocity().getNorm().getReal());
205 }
206 Assertions.assertEquals(shiftPositionErrorWithin, maxShiftPError, 0.01 * shiftPositionErrorWithin);
207 Assertions.assertEquals(interpolationPositionErrorWithin, maxInterpolationPError,
208 0.01 * interpolationPositionErrorWithin);
209 Assertions.assertEquals(shiftVelocityErrorWithin, maxShiftVError, 0.01 * shiftVelocityErrorWithin);
210 Assertions.assertEquals(interpolationVelocityErrorWithin, maxInterpolationVError,
211 0.01 * interpolationVelocityErrorWithin);
212
213
214 maxShiftPError = 0;
215 maxInterpolationPError = 0;
216 maxShiftVError = 0;
217 maxInterpolationVError = 0;
218 for (Binary64 dt = zero.add(500.0); dt.getReal() < 650.0; dt = dt.add(1.0)) {
219 FieldAbsoluteDate<Binary64> t = initialOrbit.getDate().shiftedBy(dt);
220 FieldPVCoordinates<Binary64> propagated = propagator.propagate(t).getPVCoordinates();
221 FieldPVCoordinates<Binary64> shiftError = new FieldPVCoordinates<>(propagated,
222 initialOrbit.shiftedBy(dt)
223 .getPVCoordinates());
224 FieldPVCoordinates<Binary64> interpolationError = new FieldPVCoordinates<>(propagated,
225 interpolator.interpolate(t, sample)
226 .getPVCoordinates());
227 maxShiftPError = FastMath.max(maxShiftPError,
228 shiftError.getPosition().getNorm().getReal());
229 maxInterpolationPError = FastMath.max(maxInterpolationPError,
230 interpolationError.getPosition().getNorm().getReal());
231 maxShiftVError = FastMath.max(maxShiftVError,
232 shiftError.getVelocity().getNorm().getReal());
233 maxInterpolationVError = FastMath.max(maxInterpolationVError,
234 interpolationError.getVelocity().getNorm().getReal());
235 }
236 Assertions.assertEquals(shiftPositionErrorFarPast, maxShiftPError, 0.01 * shiftPositionErrorFarPast);
237 Assertions.assertEquals(interpolationPositionErrorFarPast, maxInterpolationPError,
238 0.01 * interpolationPositionErrorFarPast);
239 Assertions.assertEquals(shiftVelocityErrorFarPast, maxShiftVError, 0.01 * shiftVelocityErrorFarPast);
240 Assertions.assertEquals(interpolationVelocityErrorFarPast, maxInterpolationVError,
241 0.01 * interpolationVelocityErrorFarPast);
242
243 }
244
245 private void doTestCircularInterpolation(boolean useDerivatives,
246 double shiftErrorWithin, double interpolationErrorWithin,
247 double shiftErrorSlightlyPast, double interpolationErrorSlightlyPast,
248 double shiftErrorFarPast, double interpolationErrorFarPast) {
249
250 Binary64 zero = field.getZero();
251 FieldAbsoluteDate<Binary64> date = new FieldAbsoluteDate<>(field);
252
253 final Binary64 ehMu = zero.add(3.9860047e14);
254 final double ae = 6.378137e6;
255 final double c20 = -1.08263e-3;
256 final double c30 = 2.54e-6;
257 final double c40 = 1.62e-6;
258 final double c50 = 2.3e-7;
259 final double c60 = -5.5e-7;
260
261 date = date.shiftedBy(584.);
262 final FieldVector3D<Binary64> position =
263 new FieldVector3D<>(zero.add(3220103.), zero.add(69623.), zero.add(6449822.));
264 final FieldVector3D<Binary64> velocity = new FieldVector3D<>(zero.add(6414.7), zero.add(-2006.), zero.add(-3180.));
265 final FieldCircularOrbit<Binary64> initialOrbit =
266 new FieldCircularOrbit<>(new FieldPVCoordinates<>(position, velocity),
267 FramesFactory.getEME2000(), date, ehMu);
268
269 FieldEcksteinHechlerPropagator<Binary64> propagator =
270 new FieldEcksteinHechlerPropagator<>(initialOrbit, ae, ehMu, c20, c30, c40, c50, c60);
271
272
273 List<FieldOrbit<Binary64>> sample = new ArrayList<>();
274 for (Binary64 dt = zero; dt.getReal() < 300.0; dt = dt.add(60.0)) {
275 FieldOrbit<Binary64> orbit = OrbitType.CIRCULAR.convertType(propagator.propagate(date.shiftedBy(dt)).getOrbit());
276 if (!useDerivatives) {
277
278 Binary64[] stateVector = MathArrays.buildArray(field, 6);
279 orbit.getType().mapOrbitToArray(orbit, PositionAngleType.TRUE, stateVector, null);
280 orbit = orbit.getType().mapArrayToOrbit(stateVector, null, PositionAngleType.TRUE,
281 orbit.getDate(), orbit.getMu(), orbit.getFrame());
282 }
283 sample.add(orbit);
284 }
285
286
287 final FieldTimeInterpolator<FieldOrbit<Binary64>, Binary64> interpolator =
288 new FieldOrbitHermiteInterpolator<>(sample.size(), 759, initialOrbit.getFrame(),
289 CartesianDerivativesFilter.USE_PVA);
290
291
292 double maxShiftError = 0.0;
293 double maxInterpolationError = 0.0;
294 for (Binary64 dt = zero; dt.getReal() < 241.0; dt = dt.add(1.0)) {
295 FieldAbsoluteDate<Binary64> Binary64 = initialOrbit.getDate().shiftedBy(dt);
296 FieldVector3D<Binary64> shifted = initialOrbit.shiftedBy(dt.getReal()).getPosition();
297 FieldVector3D<Binary64> interpolated = interpolator.interpolate(Binary64, sample).getPosition();
298 FieldVector3D<Binary64> propagated = propagator.propagate(Binary64).getPosition();
299 maxShiftError = FastMath.max(maxShiftError, shifted.subtract(propagated).getNorm().getReal());
300 maxInterpolationError =
301 FastMath.max(maxInterpolationError, interpolated.subtract(propagated).getNorm().getReal());
302 }
303 Assertions.assertEquals(shiftErrorWithin, maxShiftError, 0.01 * shiftErrorWithin);
304 Assertions.assertEquals(interpolationErrorWithin, maxInterpolationError, 0.01 * interpolationErrorWithin);
305
306
307 maxShiftError = 0;
308 maxInterpolationError = 0;
309 for (Binary64 dt = zero.add(240); dt.getReal() < 300.0; dt = dt.add(1.0)) {
310 FieldAbsoluteDate<Binary64> Binary64 = initialOrbit.getDate().shiftedBy(dt);
311 FieldVector3D<Binary64> shifted = initialOrbit.shiftedBy(dt).getPosition();
312 FieldVector3D<Binary64> interpolated = interpolator.interpolate(Binary64, sample).getPosition();
313 FieldVector3D<Binary64> propagated = propagator.propagate(Binary64).getPosition();
314 maxShiftError = FastMath.max(maxShiftError, shifted.subtract(propagated).getNorm().getReal());
315 maxInterpolationError =
316 FastMath.max(maxInterpolationError, interpolated.subtract(propagated).getNorm().getReal());
317 }
318 Assertions.assertEquals(shiftErrorSlightlyPast, maxShiftError, 0.01 * shiftErrorSlightlyPast);
319 Assertions.assertEquals(interpolationErrorSlightlyPast, maxInterpolationError,
320 0.01 * interpolationErrorSlightlyPast);
321
322
323
324 maxShiftError = 0;
325 maxInterpolationError = 0;
326 for (Binary64 dt = zero.add(300); dt.getReal() < 1000; dt = dt.add(1.0)) {
327 FieldAbsoluteDate<Binary64> t = initialOrbit.getDate().shiftedBy(dt);
328 FieldVector3D<Binary64> shifted = initialOrbit.shiftedBy(dt).getPosition();
329 FieldVector3D<Binary64> interpolated = interpolator.interpolate(t, sample).getPosition();
330 FieldVector3D<Binary64> propagated = propagator.propagate(t).getPosition();
331 maxShiftError = FastMath.max(maxShiftError, shifted.subtract(propagated).getNorm().getReal());
332 maxInterpolationError =
333 FastMath.max(maxInterpolationError, interpolated.subtract(propagated).getNorm().getReal());
334 }
335 Assertions.assertEquals(shiftErrorFarPast, maxShiftError, 0.01 * shiftErrorFarPast);
336 Assertions.assertEquals(interpolationErrorFarPast, maxInterpolationError, 0.01 * interpolationErrorFarPast);
337
338 }
339
340 private void doTestEquinoctialInterpolation(boolean useDerivatives,
341 double shiftErrorWithin, double interpolationErrorWithin,
342 double shiftErrorSlightlyPast, double interpolationErrorSlightlyPast,
343 double shiftErrorFarPast, double interpolationErrorFarPast) {
344
345 Binary64 zero = field.getZero();
346 FieldAbsoluteDate<Binary64> date = new FieldAbsoluteDate<>(field);
347
348 final double ehMu = 3.9860047e14;
349 final double ae = 6.378137e6;
350 final double c20 = -1.08263e-3;
351 final double c30 = 2.54e-6;
352 final double c40 = 1.62e-6;
353 final double c50 = 2.3e-7;
354 final double c60 = -5.5e-7;
355
356 date = date.shiftedBy(584.);
357 final FieldVector3D<Binary64> position =
358 new FieldVector3D<>(zero.add(3220103.), zero.add(69623.), zero.add(6449822.));
359 final FieldVector3D<Binary64> velocity = new FieldVector3D<>(zero.add(6414.7), zero.add(-2006.), zero.add(-3180.));
360 final FieldEquinoctialOrbit<Binary64> initialOrbit =
361 new FieldEquinoctialOrbit<>(new FieldPVCoordinates<>(position, velocity),
362 FramesFactory.getEME2000(), date, zero.add(ehMu));
363
364 FieldEcksteinHechlerPropagator<Binary64> propagator =
365 new FieldEcksteinHechlerPropagator<>(initialOrbit, ae, zero.add(ehMu), c20, c30, c40, c50, c60);
366
367
368 List<FieldOrbit<Binary64>> sample = new ArrayList<>();
369 for (double dt = 0; dt < 300.0; dt += 60.0) {
370 FieldOrbit<Binary64> orbit =
371 OrbitType.EQUINOCTIAL.convertType(propagator.propagate(date.shiftedBy(dt)).getOrbit());
372 if (!useDerivatives) {
373
374 Binary64[] stateVector = MathArrays.buildArray(field, 6);
375 orbit.getType().mapOrbitToArray(orbit, PositionAngleType.TRUE, stateVector, null);
376 orbit = orbit.getType().mapArrayToOrbit(stateVector, null, PositionAngleType.TRUE,
377 orbit.getDate(), orbit.getMu(), orbit.getFrame());
378 }
379 sample.add(orbit);
380 }
381
382
383 final FieldTimeInterpolator<FieldOrbit<Binary64>, Binary64> interpolator =
384 new FieldOrbitHermiteInterpolator<>(sample.size(), 759, initialOrbit.getFrame(),
385 CartesianDerivativesFilter.USE_PVA);
386
387
388 double maxShiftError = 0;
389 double maxInterpolationError = 0;
390 for (Binary64 dt = zero; dt.getReal() < 241.0; dt = dt.add(1.0)) {
391 FieldAbsoluteDate<Binary64> t = initialOrbit.getDate().shiftedBy(dt);
392 FieldVector3D<Binary64> shifted = initialOrbit.shiftedBy(dt.getReal()).getPosition();
393 FieldVector3D<Binary64> interpolated = interpolator.interpolate(t, sample).getPosition();
394 FieldVector3D<Binary64> propagated = propagator.propagate(t).getPosition();
395 maxShiftError = FastMath.max(maxShiftError, shifted.subtract(propagated).getNorm().getReal());
396 maxInterpolationError =
397 FastMath.max(maxInterpolationError, interpolated.subtract(propagated).getNorm().getReal());
398 }
399 Assertions.assertEquals(shiftErrorWithin, maxShiftError, 0.01 * shiftErrorWithin);
400 Assertions.assertEquals(interpolationErrorWithin, maxInterpolationError, 0.01 * interpolationErrorWithin);
401
402
403 maxShiftError = 0;
404 maxInterpolationError = 0;
405 for (Binary64 dt = zero.add(240); dt.getReal() < 300.0; dt = dt.add(1.0)) {
406 FieldAbsoluteDate<Binary64> t = initialOrbit.getDate().shiftedBy(dt);
407 FieldVector3D<Binary64> shifted = initialOrbit.shiftedBy(dt).getPosition();
408 FieldVector3D<Binary64> interpolated = interpolator.interpolate(t, sample).getPosition();
409 FieldVector3D<Binary64> propagated = propagator.propagate(t).getPosition();
410 maxShiftError = FastMath.max(maxShiftError, shifted.subtract(propagated).getNorm().getReal());
411 maxInterpolationError =
412 FastMath.max(maxInterpolationError, interpolated.subtract(propagated).getNorm().getReal());
413 }
414 Assertions.assertEquals(shiftErrorSlightlyPast, maxShiftError, 0.01 * shiftErrorSlightlyPast);
415 Assertions.assertEquals(interpolationErrorSlightlyPast, maxInterpolationError,
416 0.01 * interpolationErrorSlightlyPast);
417
418
419
420 maxShiftError = 0;
421 maxInterpolationError = 0;
422 for (Binary64 dt = zero.add(300); dt.getReal() < 1000; dt = dt.add(1.0)) {
423 FieldAbsoluteDate<Binary64> t = initialOrbit.getDate().shiftedBy(dt);
424 FieldVector3D<Binary64> shifted = initialOrbit.shiftedBy(dt).getPosition();
425 FieldVector3D<Binary64> interpolated = interpolator.interpolate(t, sample).getPosition();
426 FieldVector3D<Binary64> propagated = propagator.propagate(t).getPosition();
427 maxShiftError = FastMath.max(maxShiftError, shifted.subtract(propagated).getNorm().getReal());
428 maxInterpolationError =
429 FastMath.max(maxInterpolationError, interpolated.subtract(propagated).getNorm().getReal());
430 }
431 Assertions.assertEquals(shiftErrorFarPast, maxShiftError, 0.01 * shiftErrorFarPast);
432 Assertions.assertEquals(interpolationErrorFarPast, maxInterpolationError, 0.01 * interpolationErrorFarPast);
433
434 }
435
436 private void doTestKeplerianInterpolation(boolean useDerivatives,
437 double shiftPositionErrorWithin, double interpolationPositionErrorWithin,
438 double shiftEccentricityErrorWithin,
439 double interpolationEccentricityErrorWithin,
440 double shiftPositionErrorSlightlyPast,
441 double interpolationPositionErrorSlightlyPast,
442 double shiftEccentricityErrorSlightlyPast,
443 double interpolationEccentricityErrorSlightlyPast) {
444
445 final Binary64 zero = field.getZero();
446 final Binary64 ehMu = zero.add(3.9860047e14);
447 final double ae = 6.378137e6;
448 final double c20 = -1.08263e-3;
449 final double c30 = 2.54e-6;
450 final double c40 = 1.62e-6;
451 final double c50 = 2.3e-7;
452 final double c60 = -5.5e-7;
453
454 final FieldAbsoluteDate<Binary64> date = FieldAbsoluteDate.getJ2000Epoch(field).shiftedBy(584.);
455 final FieldVector3D<Binary64> position =
456 new FieldVector3D<>(field.getZero().add(3220103.), field.getZero().add(69623.),
457 field.getZero().add(6449822.));
458 final FieldVector3D<Binary64> velocity =
459 new FieldVector3D<>(field.getZero().add(6414.7), field.getZero().add(-2006.), field.getZero().add(-3180.));
460 final FieldKeplerianOrbit<Binary64> initialOrbit =
461 new FieldKeplerianOrbit<>(new FieldPVCoordinates<>(position, velocity),
462 FramesFactory.getEME2000(), date, ehMu);
463
464 FieldEcksteinHechlerPropagator<Binary64> propagator =
465 new FieldEcksteinHechlerPropagator<>(initialOrbit, ae, ehMu, c20, c30, c40, c50, c60);
466
467
468 List<FieldOrbit<Binary64>> sample = new ArrayList<>();
469 for (double dt = 0; dt < 300.0; dt += 60.0) {
470 FieldOrbit<Binary64> orbit =
471 OrbitType.KEPLERIAN.convertType(propagator.propagate(date.shiftedBy(dt)).getOrbit());
472 if (!useDerivatives) {
473
474 Binary64[] stateVector = MathArrays.buildArray(field, 6);
475 orbit.getType().mapOrbitToArray(orbit, PositionAngleType.TRUE, stateVector, null);
476 orbit = orbit.getType().mapArrayToOrbit(stateVector, null, PositionAngleType.TRUE,
477 orbit.getDate(), orbit.getMu(), orbit.getFrame());
478 }
479 sample.add(orbit);
480 }
481
482
483 final FieldTimeInterpolator<FieldOrbit<Binary64>, Binary64> interpolator =
484 new FieldOrbitHermiteInterpolator<>(sample.size(), 359, initialOrbit.getFrame(),
485 CartesianDerivativesFilter.USE_PVA);
486
487
488
489
490 double maxShiftPositionError = 0;
491 double maxInterpolationPositionError = 0;
492 double maxShiftEccentricityError = 0;
493 double maxInterpolationEccentricityError = 0;
494 for (double dt = 0; dt < 241.0; dt += 1.0) {
495 FieldAbsoluteDate<Binary64> t = initialOrbit.getDate().shiftedBy(dt);
496 FieldVector3D<Binary64> shiftedP = initialOrbit.shiftedBy(dt).getPosition();
497 FieldVector3D<Binary64> interpolatedP = interpolator.interpolate(t, sample).getPosition();
498 FieldVector3D<Binary64> propagatedP = propagator.propagate(t).getPosition();
499 Binary64 shiftedE = initialOrbit.shiftedBy(zero.add(dt)).getE();
500 Binary64 interpolatedE = interpolator.interpolate(t, sample).getE();
501 Binary64 propagatedE = propagator.propagate(t).getOrbit().getE();
502 maxShiftPositionError =
503 FastMath.max(maxShiftPositionError, shiftedP.subtract(propagatedP).getNorm().getReal());
504 maxInterpolationPositionError =
505 FastMath.max(maxInterpolationPositionError, interpolatedP.subtract(propagatedP).getNorm().getReal());
506 maxShiftEccentricityError =
507 FastMath.max(maxShiftEccentricityError, shiftedE.subtract(propagatedE).abs().getReal());
508 maxInterpolationEccentricityError =
509 FastMath.max(maxInterpolationEccentricityError, interpolatedE.subtract(propagatedE).abs().getReal());
510 }
511 Assertions.assertEquals(shiftPositionErrorWithin, maxShiftPositionError, 0.01 * shiftPositionErrorWithin);
512 Assertions.assertEquals(interpolationPositionErrorWithin, maxInterpolationPositionError,
513 0.01 * interpolationPositionErrorWithin);
514 Assertions.assertEquals(shiftEccentricityErrorWithin, maxShiftEccentricityError,
515 0.01 * shiftEccentricityErrorWithin);
516 Assertions.assertEquals(interpolationEccentricityErrorWithin, maxInterpolationEccentricityError,
517 0.01 * interpolationEccentricityErrorWithin);
518
519
520
521
522
523 maxShiftPositionError = 0;
524 maxInterpolationPositionError = 0;
525 maxShiftEccentricityError = 0;
526 maxInterpolationEccentricityError = 0;
527 for (double dt = 240; dt < 600; dt += 1.0) {
528 FieldAbsoluteDate<Binary64> t = initialOrbit.getDate().shiftedBy(dt);
529 FieldVector3D<Binary64> shiftedP = initialOrbit.shiftedBy(zero.add(dt)).getPosition();
530 FieldVector3D<Binary64> interpolatedP = interpolator.interpolate(t, sample).getPosition();
531 FieldVector3D<Binary64> propagatedP = propagator.propagate(t).getPosition();
532 Binary64 shiftedE = initialOrbit.shiftedBy(zero.add(dt)).getE();
533 Binary64 interpolatedE = interpolator.interpolate(t, sample).getE();
534 Binary64 propagatedE = propagator.propagate(t).getOrbit().getE();
535 maxShiftPositionError =
536 FastMath.max(maxShiftPositionError, shiftedP.subtract(propagatedP).getNorm().getReal());
537 maxInterpolationPositionError =
538 FastMath.max(maxInterpolationPositionError, interpolatedP.subtract(propagatedP).getNorm().getReal());
539 maxShiftEccentricityError =
540 FastMath.max(maxShiftEccentricityError, shiftedE.subtract(propagatedE).abs().getReal());
541 maxInterpolationEccentricityError =
542 FastMath.max(maxInterpolationEccentricityError, interpolatedE.subtract(propagatedE).abs().getReal());
543 }
544 Assertions.assertEquals(shiftPositionErrorSlightlyPast, maxShiftPositionError,
545 0.01 * shiftPositionErrorSlightlyPast);
546 Assertions.assertEquals(interpolationPositionErrorSlightlyPast, maxInterpolationPositionError,
547 0.01 * interpolationPositionErrorSlightlyPast);
548 Assertions.assertEquals(shiftEccentricityErrorSlightlyPast, maxShiftEccentricityError,
549 0.01 * shiftEccentricityErrorSlightlyPast);
550 Assertions.assertEquals(interpolationEccentricityErrorSlightlyPast, maxInterpolationEccentricityError,
551 0.01 * interpolationEccentricityErrorSlightlyPast);
552
553 }
554
555 @Test
556 @DisplayName("test default constructor")
557 void testDefaultConstructor() {
558
559 final int interpolationPoints = 2;
560
561
562 final FieldOrbitHermiteInterpolator<Binary64> interpolator =
563 new FieldOrbitHermiteInterpolator<>(FramesFactory.getGCRF());
564
565
566 Assertions.assertEquals(AbstractTimeInterpolator.DEFAULT_EXTRAPOLATION_THRESHOLD_SEC,
567 interpolator.getExtrapolationThreshold());
568 Assertions.assertEquals(interpolationPoints,
569 interpolator.getNbInterpolationPoints());
570 Assertions.assertEquals(CartesianDerivativesFilter.USE_PVA, interpolator.getPVAFilter());
571 }
572 }