1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17 package org.orekit.models.earth.atmosphere;
18
19 import org.hipparchus.CalculusFieldElement;
20 import org.hipparchus.Field;
21 import org.hipparchus.analysis.differentiation.DSFactory;
22 import org.hipparchus.analysis.differentiation.DerivativeStructure;
23 import org.hipparchus.analysis.differentiation.FiniteDifferencesDifferentiator;
24 import org.hipparchus.exception.LocalizedCoreFormats;
25 import org.hipparchus.geometry.euclidean.threed.FieldVector3D;
26 import org.hipparchus.geometry.euclidean.threed.Vector3D;
27 import org.hipparchus.random.RandomGenerator;
28 import org.hipparchus.random.Well19937a;
29 import org.hipparchus.util.Binary64;
30 import org.hipparchus.util.Binary64Field;
31 import org.hipparchus.util.FastMath;
32 import org.hipparchus.util.Precision;
33 import org.junit.jupiter.api.Assertions;
34 import org.junit.jupiter.api.BeforeEach;
35 import org.junit.jupiter.api.Test;
36 import org.orekit.Utils;
37 import org.orekit.bodies.CelestialBody;
38 import org.orekit.bodies.CelestialBodyFactory;
39 import org.orekit.bodies.GeodeticPoint;
40 import org.orekit.bodies.OneAxisEllipsoid;
41 import org.orekit.errors.OrekitException;
42 import org.orekit.errors.OrekitMessages;
43 import org.orekit.frames.Frame;
44 import org.orekit.frames.FramesFactory;
45 import org.orekit.models.earth.ReferenceEllipsoid;
46 import org.orekit.time.AbsoluteDate;
47 import org.orekit.time.DateComponents;
48 import org.orekit.time.FieldAbsoluteDate;
49 import org.orekit.time.TimeComponents;
50 import org.orekit.time.TimeScalesFactory;
51 import org.orekit.utils.Constants;
52 import org.orekit.utils.IERSConventions;
53 import org.orekit.utils.PVCoordinatesProvider;
54 import org.orekit.utils.TimeStampedFieldPVCoordinates;
55 import org.orekit.utils.TimeStampedPVCoordinates;
56
57 import java.lang.reflect.Constructor;
58 import java.lang.reflect.InvocationTargetException;
59 import java.lang.reflect.Method;
60
61
62 class NRLMSISE00Test {
63
64 @Test
65 void testLegacy() throws
66 NoSuchMethodException, SecurityException, InstantiationException,
67 IllegalAccessException, IllegalArgumentException, InvocationTargetException {
68
69
70 Frame itrf = FramesFactory.getITRF(IERSConventions.IERS_2010, true);
71 OneAxisEllipsoid earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
72 Constants.WGS84_EARTH_FLATTENING, itrf);
73 NRLMSISE00 atm = new NRLMSISE00(null, null, earth);
74
75
76 final int doy = 172;
77 final double sec = 29000.;
78 final double alt = 400.;
79 final double lat = 60.;
80 final double lon = -70.;
81 final double hl = 16.;
82 final double f107a = 150.;
83 final double f107 = 150.;
84 double[] ap = {4., 100., 100., 100., 100., 100., 100.};
85 final boolean print = false;
86
87 Method gtd7 = getOutputClass().getDeclaredMethod("gtd7", Double.TYPE);
88 gtd7.setAccessible(true);
89
90
91 final Object out1 = createOutput(atm, doy, sec, lat, lon, hl, f107a, f107, ap);
92 gtd7.invoke(out1, alt);
93 checkLegacy(1, out1, print);
94
95
96 final int doy2 = 81;
97 final Object out2 = createOutput(atm, doy2, sec, lat, lon, hl, f107a, f107, ap);
98 gtd7.invoke(out2, alt);
99 checkLegacy(2, out2, print);
100
101
102 final double sec3 = 75000.;
103 final double alt3 = 1000.;
104 final Object out3 = createOutput(atm, doy, sec3, lat, lon, hl, f107a, f107, ap);
105 gtd7.invoke(out3, alt3);
106 checkLegacy(3, out3, print);
107
108
109 final double alt4 = 100.;
110 final Object out4 = createOutput(atm, doy, sec, lat, lon, hl, f107a, f107, ap);
111 gtd7.invoke(out4, alt4);
112 checkLegacy(4, out4, print);
113
114
115 final double lat5 = 0.;
116 final Object out5 = createOutput(atm, doy, sec, lat5, lon, hl, f107a, f107, ap);
117 gtd7.invoke(out5, alt);
118 checkLegacy(5, out5, print);
119
120
121 final double lon6 = 0.;
122 final Object out6 = createOutput(atm, doy, sec, lat, lon6, hl, f107a, f107, ap);
123 gtd7.invoke(out6, alt);
124 checkLegacy(6, out6, print);
125
126
127 final double hl7 = 4.;
128 final Object out7 = createOutput(atm, doy, sec, lat, lon, hl7, f107a, f107, ap);
129 gtd7.invoke(out7, alt);
130 checkLegacy(7, out7, print);
131
132
133 final double f107a8 = 70.;
134 final Object out8 = createOutput(atm, doy, sec, lat, lon, hl, f107a8, f107, ap);
135 gtd7.invoke(out8, alt);
136 checkLegacy(8, out8, print);
137
138
139 final double f1079 = 180.;
140 final Object out9 = createOutput(atm, doy, sec, lat, lon, hl, f107a, f1079, ap);
141 gtd7.invoke(out9, alt);
142 checkLegacy(9, out9, print);
143
144
145 ap[0] = 40.;
146 final Object out10 = createOutput(atm, doy, sec, lat, lon, hl, f107a, f107, ap);
147 gtd7.invoke(out10, alt);
148 checkLegacy(10, out10, print);
149 ap[0] = 4.;
150
151
152 final double alt11 = 0.;
153 final Object out11 = createOutput(atm, doy, sec, lat, lon, hl, f107a, f107, ap);
154 gtd7.invoke(out11, alt11);
155 checkLegacy(11, out11, print);
156
157
158 final double alt12 = 10.;
159 final Object out12 = createOutput(atm, doy, sec, lat, lon, hl, f107a, f107, ap);
160 gtd7.invoke(out12, alt12);
161 checkLegacy(12, out12, print);
162
163
164 final double alt13 = 30.;
165 final Object out13 = createOutput(atm, doy, sec, lat, lon, hl, f107a, f107, ap);
166 gtd7.invoke(out13, alt13);
167 checkLegacy(13, out13, print);
168
169
170 final double alt14 = 50.;
171 final Object out14 = createOutput(atm, doy, sec, lat, lon, hl, f107a, f107, ap);
172 gtd7.invoke(out14, alt14);
173 checkLegacy(14, out14, print);
174
175
176 final double alt15 = 70.;
177 final Object out15 = createOutput(atm, doy, sec, lat, lon, hl, f107a, f107, ap);
178 gtd7.invoke(out15, alt15);
179 checkLegacy(15, out15, print);
180
181
182 NRLMSISE00 otherAtm = atm.withSwitch(9, -1);
183 final Object out16 = createOutput(otherAtm, doy, sec, lat, lon, hl, f107a, f107, ap);
184 gtd7.invoke(out16, alt);
185 checkLegacy(16, out16, print);
186
187
188 final double alt17 = 100.;
189 final Object out17 = createOutput(otherAtm, doy, sec, lat, lon, hl, f107a, f107, ap);
190 gtd7.invoke(out17, alt17);
191 checkLegacy(17, out17, print);
192
193 }
194
195 @Test
196 void testDensity() throws
197 InstantiationException, IllegalAccessException,
198 IllegalArgumentException, InvocationTargetException,
199 NoSuchMethodException, SecurityException {
200
201 final InputParams ip = new InputParams();
202
203 final CelestialBody sun = CelestialBodyFactory.getSun();
204
205 final Frame itrf = FramesFactory.getITRF(IERSConventions.IERS_2010, true);
206 final OneAxisEllipsoid earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
207 Constants.WGS84_EARTH_FLATTENING, itrf);
208
209 final NRLMSISE00 atm = new NRLMSISE00(ip, sun, earth).withSwitch(9, -1);
210
211 final AbsoluteDate date = new AbsoluteDate(new DateComponents(2003, 172),
212 new TimeComponents(29000.),
213 TimeScalesFactory.getUT1(IERSConventions.IERS_2010, true));
214
215 final double alt = 400.;
216 final double lat = 60.;
217 final double lon = -70.;
218 final GeodeticPoint point = new GeodeticPoint(FastMath.toRadians(lat),
219 FastMath.toRadians(lon),
220 alt * 1000.);
221 final Vector3D pos = earth.transform(point);
222
223
224 try {
225 atm.getDensity(date.shiftedBy(2 * Constants.JULIAN_YEAR), pos, itrf);
226 Assertions.fail("an exception should have been thrown");
227 } catch (OrekitException oe) {
228 Assertions.assertEquals(OrekitMessages.NO_SOLAR_ACTIVITY_AT_DATE, oe.getSpecifier());
229 }
230
231 final double rho = atm.getDensity(date, pos, itrf);
232 final double lst = 29000. / 3600. - 70. / 15.;
233 final double[] ap = {4., 100., 100., 100., 100., 100., 100.};
234
235 Class<?> outputClass = getOutputClass();
236 Constructor<?> cons = outputClass.getDeclaredConstructor(NRLMSISE00.class,
237 Integer.TYPE,
238 Double.TYPE,
239 Double.TYPE,
240 Double.TYPE,
241 Double.TYPE,
242 Double.TYPE,
243 Double.TYPE,
244 double[].class);
245 cons.setAccessible(true);
246 Method gtd7d = outputClass.getDeclaredMethod("gtd7d", Double.TYPE);
247 gtd7d.setAccessible(true);
248 Method getDensity = outputClass.getDeclaredMethod("getDensity", Integer.TYPE);
249 getDensity.setAccessible(true);
250
251 final Object out = createOutput(atm, 172, 29000., 60., -70, lst, 150., 150., ap);
252 gtd7d.invoke(out, 400.0);
253 Assertions.assertEquals(rho, ((Double) getDensity.invoke(out, 5)).doubleValue(), rho * 1.e-3);
254
255 }
256
257 @Test
258 void testDensityField() {
259
260 final InputParams ip = new InputParams();
261
262 final CelestialBody sun = CelestialBodyFactory.getSun();
263
264 final Frame itrf = FramesFactory.getITRF(IERSConventions.IERS_2010, true);
265 final OneAxisEllipsoid earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
266 Constants.WGS84_EARTH_FLATTENING, itrf);
267
268 final NRLMSISE00 atm = new NRLMSISE00(ip, sun, earth);
269
270 final AbsoluteDate date = new AbsoluteDate(new DateComponents(2003, 172),
271 new TimeComponents(29000.),
272 TimeScalesFactory.getUT1(IERSConventions.IERS_2010, true));
273
274 final double alt = 400.;
275 final double lat = 60.;
276 final double lon = -70.;
277 final GeodeticPoint point = new GeodeticPoint(FastMath.toRadians(lat),
278 FastMath.toRadians(lon),
279 alt * 1000.);
280 final Vector3D pos = earth.transform(point);
281 Field<Binary64> field = Binary64Field.getInstance();
282
283
284 try {
285 atm.getDensity(new FieldAbsoluteDate<>(field, date).shiftedBy(2 * Constants.JULIAN_YEAR),
286 new FieldVector3D<>(field.getOne(), pos),
287 itrf);
288 Assertions.fail("an exception should have been thrown");
289 } catch (OrekitException oe) {
290 Assertions.assertEquals(OrekitMessages.NO_SOLAR_ACTIVITY_AT_DATE, oe.getSpecifier());
291 }
292
293 final double rho = atm.getDensity(date, pos, itrf);
294 final Binary64 rho64 = atm.getDensity(new FieldAbsoluteDate<>(field, date),
295 new FieldVector3D<>(field, pos),
296 itrf);
297
298 Assertions.assertEquals(rho, rho64.getReal(), rho * 2.0e-11);
299
300 }
301
302 @Test
303 void testDensityGradient() {
304
305 final InputParams ip = new InputParams();
306
307 final CelestialBody sun = CelestialBodyFactory.getSun();
308
309 final Frame itrf = FramesFactory.getITRF(IERSConventions.IERS_2010, true);
310 final OneAxisEllipsoid earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
311 Constants.WGS84_EARTH_FLATTENING, itrf);
312
313 final NRLMSISE00 atm = new NRLMSISE00(ip, sun, earth);
314
315 final AbsoluteDate date = new AbsoluteDate(new DateComponents(2003, 172),
316 new TimeComponents(29000.),
317 TimeScalesFactory.getUT1(IERSConventions.IERS_2010, true));
318
319 final double alt = 400.;
320 final double lat = 60.;
321 final double lon = -70.;
322 final GeodeticPoint point = new GeodeticPoint(FastMath.toRadians(lat),
323 FastMath.toRadians(lon),
324 alt * 1000.);
325 final Vector3D pos = earth.transform(point);
326
327
328 DerivativeStructure zero = new DSFactory(1, 1).variable(0, 0.0);
329 FiniteDifferencesDifferentiator differentiator = new FiniteDifferencesDifferentiator(5, 10.0);
330 DerivativeStructure rhoX = differentiator.
331 differentiate((double x) -> {
332 try {
333 return atm.getDensity(date, new Vector3D(1, pos, x, Vector3D.PLUS_I), itrf);
334 } catch (OrekitException oe) {
335 return Double.NaN;
336 }
337 }). value(zero);
338 DerivativeStructure rhoY = differentiator.
339 differentiate((double y) -> {
340 try {
341 return atm.getDensity(date, new Vector3D(1, pos, y, Vector3D.PLUS_J), itrf);
342 } catch (OrekitException oe) {
343 return Double.NaN;
344 }
345 }). value(zero);
346 DerivativeStructure rhoZ = differentiator.
347 differentiate((double z) -> {
348 try {
349 return atm.getDensity(date, new Vector3D(1, pos, z, Vector3D.PLUS_K), itrf);
350 } catch (OrekitException oe) {
351 return Double.NaN;
352 }
353 }). value(zero);
354
355 DSFactory factory3 = new DSFactory(3, 1);
356 Field<DerivativeStructure> field = factory3.getDerivativeField();
357 final DerivativeStructure rhoDS = atm.getDensity(new FieldAbsoluteDate<>(field, date),
358 new FieldVector3D<>(factory3.variable(0, pos.getX()),
359 factory3.variable(1, pos.getY()),
360 factory3.variable(2, pos.getZ())),
361 itrf);
362
363 Assertions.assertEquals(rhoX.getValue(), rhoDS.getReal(), rhoX.getValue() * 2.0e-13);
364 Assertions.assertEquals(rhoY.getValue(), rhoDS.getReal(), rhoY.getValue() * 2.0e-13);
365 Assertions.assertEquals(rhoZ.getValue(), rhoDS.getReal(), rhoZ.getValue() * 2.0e-13);
366 Assertions.assertEquals(rhoX.getPartialDerivative(1),
367 rhoDS.getPartialDerivative(1, 0, 0),
368 FastMath.abs(2.0e-10 * rhoX.getPartialDerivative(1)));
369 Assertions.assertEquals(rhoY.getPartialDerivative(1),
370 rhoDS.getPartialDerivative(0, 1, 0),
371 FastMath.abs(2.0e-10 * rhoY.getPartialDerivative(1)));
372 Assertions.assertEquals(rhoZ.getPartialDerivative(1),
373 rhoDS.getPartialDerivative(0, 0, 1),
374 FastMath.abs(2.0e-10 * rhoY.getPartialDerivative(1)));
375
376 }
377
378 @Test
379 void testWrongNumberLow() {
380 try {
381 new NRLMSISE00(null, null, null).withSwitch(0, 17);
382 Assertions.fail("an exception should have been thrown");
383 } catch (OrekitException oe) {
384 Assertions.assertEquals(LocalizedCoreFormats.OUT_OF_RANGE_SIMPLE, oe.getSpecifier());
385 Assertions.assertEquals( 0, oe.getParts()[0]);
386 Assertions.assertEquals( 1, oe.getParts()[1]);
387 Assertions.assertEquals(23, oe.getParts()[2]);
388 }
389 }
390
391 @Test
392 void testWrongNumberHigh() {
393 try {
394 new NRLMSISE00(null, null, null).withSwitch(24, 17);
395 Assertions.fail("an exception should have been thrown");
396 } catch (OrekitException oe) {
397 Assertions.assertEquals(LocalizedCoreFormats.OUT_OF_RANGE_SIMPLE, oe.getSpecifier());
398 Assertions.assertEquals(24, oe.getParts()[0]);
399 Assertions.assertEquals( 1, oe.getParts()[1]);
400 Assertions.assertEquals(23, oe.getParts()[2]);
401 }
402 }
403
404 @Test
405 void testGlobe7SwitchesOn() {
406 RandomGenerator random = new Well19937a(0xb9d06451353d23cbl);
407 NRLMSISE00 atm = new NRLMSISE00(null, null, null);
408 for (int i = 1; i <= 23; ++i) {
409 atm = atm.withSwitch(i, 1);
410 }
411 doTestDoubleMethod(atm, random, "globe7", 2.0e-14, 2.2e-16);
412 }
413
414 @Test
415 void testGlobe7SwitchesOff() {
416 RandomGenerator random = new Well19937a(0x778b486a40464b8fl);
417 NRLMSISE00 atm = new NRLMSISE00(null, null, null);
418 for (int i = 1; i <= 23; ++i) {
419 atm = atm.withSwitch(i, 0);
420 }
421 doTestDoubleMethod(atm, random, "globe7", 1.0e-50, 1.0e-50);
422 }
423
424 @Test
425 void testGlobe7SwitchesRandom() {
426 RandomGenerator random = new Well19937a(0xe20a69235cc9583dl);
427 NRLMSISE00 atm = new NRLMSISE00(null, null, null);
428 for (int i = 1; i <= 23; ++i) {
429 atm = atm.withSwitch(i, random.nextInt(3) - 2);
430 }
431 doTestDoubleMethod(atm, random, "globe7", 3.0e-14, 4.0e-16);
432 }
433
434 @Test
435 void testGlob7sSwitchesOn() {
436 RandomGenerator random = new Well19937a(0xc7c218fabec5e98cl);
437 NRLMSISE00 atm = new NRLMSISE00(null, null, null);
438 for (int i = 1; i <= 23; ++i) {
439 atm = atm.withSwitch(i, 1);
440 }
441 doTestDoubleMethod(atm, random, "glob7s", 4.0e-15, 9.0e-16);
442 }
443
444 @Test
445 void testGlob7sSwitchesOff() {
446 RandomGenerator random = new Well19937a(0x141f7aa933299a83l);
447 NRLMSISE00 atm = new NRLMSISE00(null, null, null);
448 for (int i = 1; i <= 23; ++i) {
449 atm = atm.withSwitch(i, 0);
450 }
451 doTestDoubleMethod(atm, random, "glob7s", 1.0e-50, 1.0e-50);
452 }
453
454 @Test
455 void testGlob7sSwitchesRandom() {
456 RandomGenerator random = new Well19937a(0x3671893ce741fc5cl);
457 NRLMSISE00 atm = new NRLMSISE00(null, null, null);
458 for (int i = 1; i <= 23; ++i) {
459 atm = atm.withSwitch(i, random.nextInt(3) - 2);
460 }
461 doTestDoubleMethod(atm, random, "glob7s", 1.0e-50, 1.0e-50);
462 }
463
464 @Test
465 void testgts7SwitchesOn() {
466 RandomGenerator random = new Well19937a(0xb6dcf73ed5e5d985l);
467 NRLMSISE00 atm = new NRLMSISE00(null, null, null);
468 for (int i = 1; i <= 23; ++i) {
469 atm = atm.withSwitch(i, 1);
470 }
471 doTestVoidMethod(atm, random, "gts7", 1.0e-50, 6.0e-14);
472 }
473
474 @Test
475 void testgts7SwitchesOff() {
476 RandomGenerator random = new Well19937a(0x0c953641bea0f6d2l);
477 NRLMSISE00 atm = new NRLMSISE00(null, null, null);
478 for (int i = 1; i <= 23; ++i) {
479 atm = atm.withSwitch(i, 0);
480 }
481 doTestVoidMethod(atm, random, "gts7", 1.0e-50, 6.0e-14);
482 }
483
484 @Test
485 void testgts7SwitchesRandom() {
486 RandomGenerator random = new Well19937a(0x7347cacb946cb93bl);
487 NRLMSISE00 atm = new NRLMSISE00(null, null, null);
488 for (int i = 1; i <= 23; ++i) {
489 atm = atm.withSwitch(i, random.nextInt(3) - 2);
490 }
491 doTestVoidMethod(atm, random, "gts7", 1.0e-50, 6.0e-14);
492 }
493
494 @Test
495 void testgtd7SwitchesOn() {
496 RandomGenerator random = new Well19937a(0x3439206bdd4dff5dl);
497 NRLMSISE00 atm = new NRLMSISE00(null, null, null);
498 for (int i = 1; i <= 23; ++i) {
499 atm = atm.withSwitch(i, 1);
500 }
501 doTestVoidMethod(atm, random, "gtd7", 5.0e-16, 3.0e-14);
502 }
503
504 @Test
505 void testgtd7SwitchesOff() {
506 RandomGenerator random = new Well19937a(0x3dc1f824e1033d1bl);
507 NRLMSISE00 atm = new NRLMSISE00(null, null, null);
508 for (int i = 1; i <= 23; ++i) {
509 atm = atm.withSwitch(i, 0);
510 }
511 doTestVoidMethod(atm, random, "gtd7", 1.0e-50, 3.0e-14);
512 }
513
514 @Test
515 void testgtd7SwitchesRandom() {
516 RandomGenerator random = new Well19937a(0xa12175ef0b689b04l);
517 NRLMSISE00 atm = new NRLMSISE00(null, null, null);
518 for (int i = 1; i <= 23; ++i) {
519 atm = atm.withSwitch(i, random.nextInt(3) - 2);
520 }
521 doTestVoidMethod(atm, random, "gtd7", 1.0e-50, 3.0e-14);
522 }
523
524 @Test
525 void testgtd7dSwitchesOn() {
526 RandomGenerator random = new Well19937a(0x4bbb424422a1b909l);
527 NRLMSISE00 atm = new NRLMSISE00(null, null, null);
528 for (int i = 1; i <= 23; ++i) {
529 atm = atm.withSwitch(i, 1);
530 }
531 doTestVoidMethod(atm, random, "gtd7d", 3.0e-16, 3.0e-14);
532 }
533
534 @Test
535 void testgtd7dSwitchesOff() {
536 RandomGenerator random = new Well19937a(0x7f6da37655e30103l);
537 NRLMSISE00 atm = new NRLMSISE00(null, null, null);
538 for (int i = 1; i <= 23; ++i) {
539 atm = atm.withSwitch(i, 0);
540 }
541 doTestVoidMethod(atm, random, "gtd7d", 1.0e-50, 2.0e-14);
542 }
543
544 @Test
545 void testgtd7dSwitchesRandom() {
546 RandomGenerator random = new Well19937a(0x4a75e29ddf23ccd7l);
547 NRLMSISE00 atm = new NRLMSISE00(null, null, null);
548 for (int i = 1; i <= 23; ++i) {
549 atm = atm.withSwitch(i, random.nextInt(3) - 2);
550 }
551 doTestVoidMethod(atm, random, "gtd7d", 1.0e-50, 3.0e-14);
552 }
553
554
555 @Test
556 void testIssue1365AvoidNan1() {
557
558
559
560
561
562 final NRLMSISE00InputParameters ip = new InputParamsIssue1365AvoidNan1();
563
564
565 final Field<Binary64> field = Binary64Field.getInstance();
566
567
568 final AbsoluteDate date = new AbsoluteDate("2005-09-10T00:26:47.232", TimeScalesFactory.getUTC());
569 final FieldAbsoluteDate<Binary64> fieldDate = new FieldAbsoluteDate<>(field, date);
570
571
572 final Frame j2000 = FramesFactory.getEME2000();
573 final Vector3D sunPosition = new Vector3D(-1.469767604504155E11, 3.030095780108449E10, 1.3136383992886505E10);
574
575
576 final CelestialBody sunAtDate = new SunPositionIssue1365AvoidNan(sunPosition, j2000);
577
578
579 final Frame itrf = FramesFactory.getITRF(IERSConventions.IERS_2010, true);
580 final OneAxisEllipsoid earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
581 Constants.WGS84_EARTH_FLATTENING, itrf);
582
583 final NRLMSISE00 atm = new NRLMSISE00(ip, sunAtDate, earth);
584
585
586 final Vector3D position = new Vector3D(-2519211.5855839024, -135107.58355852086, -6238233.867025304);
587 final FieldVector3D<Binary64> fieldPosition = new FieldVector3D<>(field, position);
588
589
590
591
592 final double density = atm.getDensity(date, position, j2000);
593 final Binary64 fieldDensity = atm.getDensity(fieldDate, fieldPosition, j2000);
594
595
596
597
598
599 Assertions.assertFalse(Double.isNaN(density));
600 Assertions.assertFalse(Double.isNaN(fieldDensity.getReal()));
601 }
602
603
604
605
606 @Test
607 void testIssue1365AvoidNan2() {
608
609
610
611
612
613 final NRLMSISE00InputParameters ip = new InputParamsIssue1365AvoidNan2();
614
615
616 final Field<Binary64> field = Binary64Field.getInstance();
617
618
619 final AbsoluteDate date = new AbsoluteDate("2005-09-10T00:00:50.000", TimeScalesFactory.getUTC());
620 final FieldAbsoluteDate<Binary64> fieldDate = new FieldAbsoluteDate<>(field, date);
621
622
623 final Frame j2000 = FramesFactory.getEME2000();
624 final Vector3D sunPosition = new Vector3D(-1.469673551020242E11, 3.0342331223223766E10, 1.3154322212769707E10);
625
626
627 final CelestialBody sunAtDate = new SunPositionIssue1365AvoidNan(sunPosition, j2000);
628
629
630 final Frame itrf = FramesFactory.getITRF(IERSConventions.IERS_2010, true);
631 final OneAxisEllipsoid earth = ReferenceEllipsoid.getWgs84(itrf);
632
633 final NRLMSISE00 atm = new NRLMSISE00(ip, sunAtDate, earth);
634
635
636 final Vector3D position = new Vector3D(-2094968.638528762, 58037.02211604678, -6396195.230946325);
637 final FieldVector3D<Binary64> fieldPosition = new FieldVector3D<>(field, position);
638
639
640
641
642
643 final double density = atm.getDensity(date, position, j2000);
644 final Binary64 fieldDensity = atm.getDensity(fieldDate, fieldPosition, j2000);
645
646
647
648
649
650 Assertions.assertFalse(Double.isNaN(density));
651 Assertions.assertFalse(Double.isNaN(fieldDensity.getReal()));
652 }
653
654
655
656
657 @Test
658 void testIssue1365AvoidNanExceptionRaised() {
659
660
661
662
663
664 final NRLMSISE00InputParameters ip = new InputParamsIssue1365AvoidNanExceptionRaised();
665
666
667 final Field<Binary64> field = Binary64Field.getInstance();
668
669
670 final AbsoluteDate date = new AbsoluteDate("2005-10-30T06:49:10.000", TimeScalesFactory.getUTC());
671 final FieldAbsoluteDate<Binary64> fieldDate = new FieldAbsoluteDate<>(field, date);
672
673
674 final Frame j2000 = FramesFactory.getEME2000();
675 final Vector3D sunPosition = new Vector3D(-1.1883503376589482E11, -8.178200712683229E10, -3.545541568807778E10);
676
677
678 final CelestialBody sunAtDate = new SunPositionIssue1365AvoidNan(sunPosition, j2000);
679
680
681 final Frame itrf = FramesFactory.getITRF(IERSConventions.IERS_2010, true);
682 final OneAxisEllipsoid earth = ReferenceEllipsoid.getWgs84(itrf);
683
684
685 final NRLMSISE00 atm = new NRLMSISE00(ip, sunAtDate, earth);
686
687
688 final Vector3D position = new Vector3D(1.e49, 1.e49, 1.e49);
689 final FieldVector3D<Binary64> fieldPosition = new FieldVector3D<>(field, position);
690
691
692
693
694 try {
695 atm.getDensity(date, position, j2000);
696 Assertions.fail("an exception should have been thrown");
697 } catch (OrekitException oe) {
698 Assertions.assertEquals(OrekitMessages.INFINITE_NRLMSISE00_DENSITY, oe.getSpecifier());
699 }
700
701 try {
702 atm.getDensity(fieldDate, fieldPosition, j2000);
703 Assertions.fail("an exception should have been thrown");
704 } catch (OrekitException oe) {
705 Assertions.assertEquals(OrekitMessages.INFINITE_NRLMSISE00_DENSITY, oe.getSpecifier());
706 }
707 }
708
709 private void doTestDoubleMethod(NRLMSISE00 atm, RandomGenerator random, String methodName,
710 double absTolerance, double relTolerance)
711 {
712 try {
713
714 final int doy = 172;
715 final double sec = 29000.;
716 final double lat = 60.;
717 final double lon = -70.;
718 final double hl = 16.;
719 final double f107a = 149.;
720 final double f107 = 150.;
721 double[] ap = {4., 100., 100., 100., 100., 100., 100.};
722
723 Method methodD = getOutputClass().getDeclaredMethod(methodName, double[].class);
724 methodD.setAccessible(true);
725
726 Method methodF = getFieldOutputClass().getDeclaredMethod(methodName, double[].class);
727 methodF.setAccessible(true);
728
729 double[] p = new double[150];
730
731 Object output = createOutput(atm, doy, sec, lat, lon, hl, f107a, f107, ap);
732 Object fieldOutput = createFieldOutput(Binary64Field.getInstance(),
733 atm, doy, sec, lat, lon, hl, f107a, f107, ap);
734 double maxAbsoluteError = 0;
735 double maxRelativeError = 0;
736 for (int i = 0; i < 100; ++i) {
737 for (int k = 0; k < p.length; ++k) {
738 p[k] = random.nextDouble();
739 }
740 double resDouble = ((Double) methodD.invoke(output, p)).doubleValue();
741 double resField = ((CalculusFieldElement<?>) methodF.invoke(fieldOutput, p)).getReal();
742 maxAbsoluteError = FastMath.max(maxAbsoluteError, FastMath.abs(resDouble - resField));
743 maxRelativeError = FastMath.max(maxRelativeError, FastMath.abs((resDouble - resField) / resDouble));
744 }
745 Assertions.assertEquals(0.0, maxAbsoluteError, absTolerance);
746 if (maxAbsoluteError != 0.0) {
747 Assertions.assertEquals(0.0, maxRelativeError, relTolerance);
748 }
749
750 } catch (NoSuchMethodException | SecurityException |
751 IllegalAccessException | IllegalArgumentException | InvocationTargetException e) {
752 Assertions.fail(e.getLocalizedMessage());
753 }
754 }
755
756 private void doTestVoidMethod(NRLMSISE00 atm, RandomGenerator random, String methodName,
757 double temperatureRelativeTolerance, double densityRelativeTolerance)
758 {
759 try {
760
761
762 final int doy = 172;
763 final double sec = 29000.;
764 final double lat = 60.;
765 final double lon = -70.;
766 final double hl = 16.;
767 final double f107a = 149.;
768 final double f107 = 150.;
769 double[] ap = {4., 100., 100., 100., 100., 100., 100.};
770
771 Method methodD = getOutputClass().getDeclaredMethod(methodName, Double.TYPE);
772 methodD.setAccessible(true);
773
774 Method methodF = getFieldOutputClass().getDeclaredMethod(methodName, CalculusFieldElement.class);
775 methodF.setAccessible(true);
776
777 Object output = createOutput(atm, doy, sec, lat, lon, hl, f107a, f107, ap);
778 Object fieldOutput = createFieldOutput(Binary64Field.getInstance(),
779 atm, doy, sec, lat, lon, hl, f107a, f107, ap);
780 double maxTemperatureError = 0;
781 double maxDensityError = 0;
782 for (int i = 0; i < 100; ++i) {
783 double alt = 500.0 * random.nextDouble();
784 methodD.invoke(output, alt);
785 methodF.invoke(fieldOutput, new Binary64(alt));
786 for (int index = 0; index < 2; ++index) {
787 double tD = getOutputTemperature(output, index);
788 double tF = getFieldOutputTemperature(fieldOutput, index);
789 maxTemperatureError = FastMath.max(maxTemperatureError, FastMath.abs((tD - tF) / tF));
790 }
791 for (int index = 0; index < 9; ++index) {
792 double dD = getOutputDensity(output, index);
793 double dF = getFieldOutputDensity(fieldOutput, index);
794 if (Double.isNaN(dD)) {
795
796
797
798 Assertions.assertTrue(Double.isNaN(dF));
799 } else if (dD == 0) {
800
801 Assertions.assertEquals(dD, dF, Precision.SAFE_MIN);
802 } else {
803 maxDensityError = FastMath.max(maxDensityError, FastMath.abs((dD - dF) / dD));
804 }
805 }
806 }
807 Assertions.assertEquals(0.0, maxTemperatureError, temperatureRelativeTolerance);
808 Assertions.assertEquals(0.0, maxDensityError, densityRelativeTolerance);
809
810 } catch (NoSuchMethodException | SecurityException |
811 IllegalAccessException | IllegalArgumentException | InvocationTargetException e) {
812 Assertions.fail(e.getLocalizedMessage());
813 }
814 }
815
816 private Class<?> getOutputClass() {
817 for (final Class<?> c : NRLMSISE00.class.getDeclaredClasses()) {
818 if (c.getName().endsWith("$Output")) {
819 return c;
820 }
821 }
822 return null;
823 }
824
825 private Object createOutput(final NRLMSISE00 atm,
826 final int doy, final double sec,
827 final double lat, final double lon, final double hl,
828 final double f107a, final double f107, final double[] ap) {
829 try {
830 Class<?> outputClass = getOutputClass();
831 Constructor<?> cons = outputClass.getDeclaredConstructor(NRLMSISE00.class,
832 Integer.TYPE,
833 Double.TYPE,
834 Double.TYPE,
835 Double.TYPE,
836 Double.TYPE,
837 Double.TYPE,
838 Double.TYPE,
839 double[].class);
840 cons.setAccessible(true);
841
842 return cons.newInstance(atm, doy, sec, lat, lon, hl, f107a, f107, ap);
843 } catch (NoSuchMethodException | SecurityException | InstantiationException |
844 IllegalAccessException | IllegalArgumentException | InvocationTargetException e) {
845 Assertions.fail(e.getLocalizedMessage());
846 return null;
847 }
848
849 }
850
851 private double getOutputDensity(Object o, int index) {
852 try {
853 Method getDensity = getOutputClass().
854 getDeclaredMethod("getDensity", Integer.TYPE);
855 getDensity.setAccessible(true);
856 return ((Double) getDensity.invoke(o, index)).doubleValue();
857 } catch (NoSuchMethodException | SecurityException |
858 IllegalAccessException | IllegalArgumentException |
859 InvocationTargetException e) {
860 Assertions.fail(e.getLocalizedMessage());
861 return Double.NaN;
862 }
863 }
864
865 private double getOutputTemperature(Object o, int index) {
866 try {
867 java.lang.reflect.Field temperaturesField = getOutputClass().getDeclaredField("temperatures");
868 temperaturesField.setAccessible(true);
869 return ((double[]) temperaturesField.get(o))[index];
870 } catch (NoSuchFieldException | SecurityException |
871 IllegalAccessException | IllegalArgumentException e) {
872 Assertions.fail(e.getLocalizedMessage());
873 return Double.NaN;
874 }
875
876 }
877
878 private Class<?> getFieldOutputClass() {
879 for (final Class<?> c : NRLMSISE00.class.getDeclaredClasses()) {
880 if (c.getName().endsWith("$FieldOutput")) {
881 return c;
882 }
883 }
884 return null;
885 }
886
887 private <T extends CalculusFieldElement<T>> Object createFieldOutput(Field<T>field,
888 final NRLMSISE00 atm,
889 final int doy, final double sec,
890 final double lat, final double lon,
891 final double hl,
892 final double f107a, final double f107,
893 final double[] ap) {
894 try {
895 Class<?> fieldOutputClass = getFieldOutputClass();
896 Constructor<?> cons = fieldOutputClass.getDeclaredConstructor(NRLMSISE00.class,
897 Integer.TYPE,
898 CalculusFieldElement.class,
899 CalculusFieldElement.class,
900 CalculusFieldElement.class,
901 CalculusFieldElement.class,
902 Double.TYPE,
903 Double.TYPE,
904 double[].class);
905 cons.setAccessible(true);
906 return cons.newInstance(atm, doy,
907 field.getZero().add(sec),
908 field.getZero().add(lat),
909 field.getZero().add(lon),
910 field.getZero().add(hl),
911 f107a, f107, ap);
912 } catch (NoSuchMethodException | SecurityException | InstantiationException |
913 IllegalAccessException | IllegalArgumentException | InvocationTargetException e) {
914 Assertions.fail(e.getLocalizedMessage());
915 return null;
916 }
917
918 }
919
920 private double getFieldOutputDensity(Object o, int index) {
921 try {
922 Method getDensity = getFieldOutputClass().
923 getDeclaredMethod("getDensity", Integer.TYPE);
924 getDensity.setAccessible(true);
925 return ((CalculusFieldElement<?>) getDensity.invoke(o, index)).getReal();
926 } catch (NoSuchMethodException | SecurityException |
927 IllegalAccessException | IllegalArgumentException |
928 InvocationTargetException e) {
929 Assertions.fail(e.getLocalizedMessage());
930 return Double.NaN;
931 }
932 }
933
934 private double getFieldOutputTemperature(Object o, int index) {
935 try {
936 java.lang.reflect.Field temperaturesField = getFieldOutputClass().getDeclaredField("temperatures");
937 temperaturesField.setAccessible(true);
938 return ((CalculusFieldElement[]) temperaturesField.get(o))[index].getReal();
939 } catch (NoSuchFieldException | SecurityException |
940 IllegalAccessException | IllegalArgumentException e) {
941 Assertions.fail(e.getLocalizedMessage());
942 return Double.NaN;
943 }
944
945 }
946
947 private void checkLegacy(final int nb, final Object out, final boolean print)
948 throws IllegalAccessException, IllegalArgumentException, InvocationTargetException {
949 final double[] tInfRef = {1.250540E+03, 1.166754E+03, 1.239892E+03, 1.027318E+03,
950 1.212396E+03, 1.220146E+03, 1.116385E+03, 1.031247E+03,
951 1.306052E+03, 1.361868E+03, 1.027318E+03, 1.027318E+03,
952 1.027318E+03, 1.027318E+03, 1.027318E+03, 1.426412E+03,
953 1.027318E+03};
954 final double[] tAltRef = {1.241416E+03, 1.161710E+03, 1.239891E+03, 2.068878E+02,
955 1.208135E+03, 1.212712E+03, 1.112999E+03, 1.024848E+03,
956 1.293374E+03, 1.347389E+03, 2.814648E+02, 2.274180E+02,
957 2.374389E+02, 2.795551E+02, 2.190732E+02, 1.408608E+03,
958 1.934071E+02};
959 final double[] dHeRef = {6.665177E+05, 3.407293E+06, 1.123767E+05, 5.411554E+07,
960 1.851122E+06, 8.673095E+05, 5.776251E+05, 3.740304E+05,
961 6.748339E+05, 5.528601E+05, 1.375488E+14, 4.427443E+13,
962 2.127829E+12, 1.412184E+11, 1.254884E+10, 5.196477E+05,
963 4.260860E+07};
964 final double[] dORef = {1.138806E+08, 1.586333E+08, 6.934130E+04, 1.918893E+11,
965 1.476555E+08, 1.278862E+08, 6.979139E+07, 4.782720E+07,
966 1.245315E+08, 1.198041E+08, 0.000000E+00, 0.000000E+00,
967 0.000000E+00, 0.000000E+00, 0.000000E+00, 1.274494E+08,
968 1.241342E+11};
969 final double[] dN2Ref = {1.998211E+07, 1.391117E+07, 4.247105E+01, 6.115826E+12,
970 1.579356E+07, 1.822577E+07, 1.236814E+07, 5.240380E+06,
971 2.369010E+07, 3.495798E+07, 2.049687E+19, 6.597567E+18,
972 3.170791E+17, 2.104370E+16, 1.874533E+15, 4.850450E+07,
973 4.929562E+12};
974 final double[] dO2Ref = {4.022764E+05, 3.262560E+05, 1.322750E-01, 1.225201E+12,
975 2.633795E+05, 2.922214E+05, 2.492868E+05, 1.759875E+05,
976 4.911583E+05, 9.339618E+05, 5.498695E+18, 1.769929E+18,
977 8.506280E+16, 5.645392E+15, 4.923051E+14, 1.720838E+06,
978 1.048407E+12};
979 final double[] dARRef = {3.557465E+03, 1.559618E+03, 2.618848E-05, 6.023212E+10,
980 1.588781E+03, 2.402962E+03, 1.405739E+03, 5.501649E+02,
981 4.578781E+03, 1.096255E+04, 2.451733E+17, 7.891680E+16,
982 3.792741E+15, 2.517142E+14, 2.239685E+13, 2.354487E+04,
983 4.993465E+10};
984 final double[] dHRef = {3.475312E+04, 4.854208E+04, 2.016750E+04, 1.059880E+07,
985 5.816167E+04, 3.686389E+04, 5.291986E+04, 8.896776E+04,
986 3.244595E+04, 2.686428E+04, 0.000000E+00, 0.000000E+00,
987 0.000000E+00, 0.000000E+00, 0.000000E+00, 2.500078E+04,
988 8.831229E+06};
989 final double[] dNRef = {4.095913E+06, 4.380967E+06, 5.741256E+03, 2.615737E+05,
990 5.478984E+06, 3.897276E+06, 1.069814E+06, 1.979741E+06,
991 5.370833E+06, 4.889974E+06, 0.000000E+00, 0.000000E+00,
992 0.000000E+00, 0.000000E+00, 0.000000E+00, 6.279210E+06,
993 2.252516E+05};
994 final double[] dAnORef = {2.667273E+04, 6.956682E+03, 2.374394E+04, 2.819879E-42,
995 1.264446E+03, 2.667273E+04, 2.667273E+04, 9.121815E+03,
996 2.667273E+04, 2.805445E+04, 0.000000E+00, 0.000000E+00,
997 0.000000E+00, 0.000000E+00, 0.000000E+00, 2.667273E+04,
998 2.415246E-42};
999 final double[] rhoRef = {4.074714E-15, 5.001846E-15, 2.756772E-18, 3.584426E-10,
1000 4.809630E-15, 4.355866E-15, 2.470651E-15, 1.571889E-15,
1001 4.564420E-15, 4.974543E-15, 1.261066E-03, 4.059139E-04,
1002 1.950822E-05, 1.294709E-06, 1.147668E-07, 5.881940E-15,
1003 2.914304E-10};
1004 final double deltaT = 1.e-2;
1005 final double deltaD = 5.e-7;
1006 final int id = nb - 1;
1007 if (print) {
1008 System.out.printf("Case #%d\n", nb);
1009 System.out.printf("Tinf: %E %E\n", tInfRef[id], getOutputTemperature(out, 0));
1010 System.out.printf("Talt: %E %E\n", tAltRef[id], getOutputTemperature(out, 1));
1011 System.out.printf("He: %E %E\n", dHeRef[id], getOutputDensity(out, 0) * 1e-6);
1012 System.out.printf("O: %E %E\n", dORef[id], getOutputDensity(out, 1) * 1e-6);
1013 System.out.printf("N2: %E %E\n", dN2Ref[id], getOutputDensity(out, 2) * 1e-6);
1014 System.out.printf("O2: %E %E\n", dO2Ref[id], getOutputDensity(out, 3) * 1e-6);
1015 System.out.printf("Ar: %E %E\n", dARRef[id], getOutputDensity(out, 4) * 1e-6);
1016 System.out.printf("H: %E %E\n", dHRef[id], getOutputDensity(out, 6) * 1e-6);
1017 System.out.printf("N: %E %E\n", dNRef[id], getOutputDensity(out, 7) * 1e-6);
1018 System.out.printf("AnO: %E %E\n", dAnORef[id], getOutputDensity(out, 8) * 1e-6);
1019 System.out.printf("Rho: %E %E\n\n", rhoRef[id], getOutputDensity(out, 5) * 1e-3);
1020 } else {
1021 Assertions.assertEquals(tInfRef[id], getOutputTemperature(out, 0), deltaT);
1022 Assertions.assertEquals(tAltRef[id], getOutputTemperature(out, 1), deltaT);
1023 Assertions.assertEquals(dHeRef[id], getOutputDensity(out, 0) * 1e-6, dHeRef[id] * deltaD);
1024 Assertions.assertEquals(dORef[id], getOutputDensity(out, 1) * 1e-6, dORef[id] * deltaD);
1025 Assertions.assertEquals(dN2Ref[id], getOutputDensity(out, 2) * 1e-6, dN2Ref[id] * deltaD);
1026 Assertions.assertEquals(dO2Ref[id], getOutputDensity(out, 3) * 1e-6, dO2Ref[id] * deltaD);
1027 Assertions.assertEquals(dARRef[id], getOutputDensity(out, 4) * 1e-6, dARRef[id] * deltaD);
1028 Assertions.assertEquals(dHRef[id], getOutputDensity(out, 6) * 1e-6, dHRef[id] * deltaD);
1029 Assertions.assertEquals(dNRef[id], getOutputDensity(out, 7) * 1e-6, dNRef[id] * deltaD);
1030 Assertions.assertEquals(dAnORef[id], getOutputDensity(out, 8) * 1e-6, dAnORef[id] * deltaD);
1031 Assertions.assertEquals(rhoRef[id], getOutputDensity(out, 5) * 1e-3, rhoRef[id] * deltaD);
1032 }
1033
1034 }
1035
1036 @BeforeEach
1037 public void setUp() {
1038 Utils.setDataRoot("regular-data");
1039 }
1040
1041 private static class InputParams implements NRLMSISE00InputParameters {
1042
1043
1044 private static final long serialVersionUID = 1L;
1045
1046
1047 public InputParams() {
1048
1049 }
1050
1051 @Override
1052 public AbsoluteDate getMinDate() {
1053 return new AbsoluteDate(2003, 1, 1, TimeScalesFactory.getUTC());
1054 }
1055
1056 @Override
1057 public AbsoluteDate getMaxDate() {
1058 return new AbsoluteDate(2003, 12, 31, TimeScalesFactory.getUTC());
1059 }
1060
1061 @Override
1062 public double getDailyFlux(AbsoluteDate date) {
1063 return 150.;
1064 }
1065
1066 @Override
1067 public double getAverageFlux(AbsoluteDate date) {
1068 return 150.;
1069 }
1070
1071 @Override
1072 public double[] getAp(AbsoluteDate date) {
1073 return new double[] {4., 100., 100., 100., 100., 100., 100.};
1074 }
1075 }
1076
1077
1078 @SuppressWarnings("serial")
1079 private static class InputParamsIssue1365AvoidNan1 implements NRLMSISE00InputParameters {
1080
1081 @Override
1082 public AbsoluteDate getMinDate() { return new AbsoluteDate(2005, 9, 10, TimeScalesFactory.getUTC()); }
1083 @Override
1084 public AbsoluteDate getMaxDate() { return new AbsoluteDate(2005, 9, 11, TimeScalesFactory.getUTC()); }
1085 @Override
1086 public double getDailyFlux(AbsoluteDate date) { return 707.6; }
1087 @Override
1088 public double getAverageFlux(AbsoluteDate date) { return 98.8; }
1089 @Override
1090 public double[] getAp(AbsoluteDate date) { return new double[] {33.0, 9.0, 18.0, 32.0, 32.0, 8.875, 7.25}; }
1091 }
1092
1093
1094 @SuppressWarnings("serial")
1095 private static class InputParamsIssue1365AvoidNan2 implements NRLMSISE00InputParameters {
1096
1097 @Override
1098 public AbsoluteDate getMinDate() { return new AbsoluteDate(2005, 9, 10, TimeScalesFactory.getUTC()); }
1099 @Override
1100 public AbsoluteDate getMaxDate() { return new AbsoluteDate(2005, 9, 11, TimeScalesFactory.getUTC()); }
1101 @Override
1102 public double getDailyFlux(AbsoluteDate date) { return 707.6; }
1103 @Override
1104 public double getAverageFlux(AbsoluteDate date) { return 98.8; }
1105 @Override
1106 public double[] getAp(AbsoluteDate date) { return new double[] {33.0, 9.0, 18.0, 32.0, 32.0, 8.875, 7.25}; }
1107 }
1108
1109
1110 @SuppressWarnings("serial")
1111 private static class InputParamsIssue1365AvoidNanExceptionRaised implements NRLMSISE00InputParameters {
1112
1113 @Override
1114 public AbsoluteDate getMinDate() { return new AbsoluteDate(2005, 10, 30, TimeScalesFactory.getUTC()); }
1115 @Override
1116 public AbsoluteDate getMaxDate() { return new AbsoluteDate(2005, 10, 31, TimeScalesFactory.getUTC()); }
1117 @Override
1118 public double getDailyFlux(AbsoluteDate date) { return 707.6; }
1119 @Override
1120 public double getAverageFlux(AbsoluteDate date) { return 98.8; }
1121 @Override
1122 public double[] getAp(AbsoluteDate date) { return new double[] {33.0, 9.0, 18.0, 32.0, 32.0, 8.875, 7.25}; }
1123 }
1124
1125
1126 @SuppressWarnings("serial")
1127 private static class SunPositionIssue1365AvoidNan implements CelestialBody {
1128
1129 private final Vector3D sunPosition;
1130 private final Frame j2000;
1131
1132
1133
1134
1135
1136 private SunPositionIssue1365AvoidNan(Vector3D sunPosition,
1137 Frame j2000) {
1138 this.sunPosition = sunPosition;
1139 this.j2000 = j2000;
1140 }
1141 @Override
1142 public Vector3D getPosition(AbsoluteDate date, Frame frame) {
1143 return sunPosition;
1144 }
1145 @Override
1146 public <T extends CalculusFieldElement<T>> FieldVector3D<T> getPosition(FieldAbsoluteDate<T> date, Frame frame) {
1147 return new FieldVector3D<>(date.getField(), getPosition(date.toAbsoluteDate(), frame));
1148 }
1149 @Override
1150 public String getName() { return "SUN"; }
1151 @Override
1152 public Frame getInertiallyOrientedFrame() { return j2000; }
1153 @Override
1154 public double getGM() { return Constants.JPL_SSD_SUN_GM; }
1155 @Override
1156 public Frame getBodyOrientedFrame() { return j2000; }
1157 }
1158 }