1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17 package org.orekit.models.earth.ionosphere.nequick;
18
19 import org.hipparchus.CalculusFieldElement;
20 import org.hipparchus.Field;
21 import org.hipparchus.geometry.euclidean.threed.FieldVector3D;
22 import org.hipparchus.geometry.euclidean.threed.Vector3D;
23 import org.hipparchus.util.Binary64Field;
24 import org.hipparchus.util.FastMath;
25 import org.junit.jupiter.api.Assertions;
26 import org.junit.jupiter.api.BeforeEach;
27 import org.junit.jupiter.api.Test;
28 import org.orekit.Utils;
29 import org.orekit.bodies.FieldGeodeticPoint;
30 import org.orekit.bodies.GeodeticPoint;
31 import org.orekit.bodies.OneAxisEllipsoid;
32 import org.orekit.frames.FramesFactory;
33 import org.orekit.frames.TopocentricFrame;
34 import org.orekit.gnss.PredefinedGnssSignal;
35 import org.orekit.orbits.CartesianOrbit;
36 import org.orekit.orbits.FieldCartesianOrbit;
37 import org.orekit.orbits.FieldOrbit;
38 import org.orekit.orbits.Orbit;
39 import org.orekit.propagation.FieldSpacecraftState;
40 import org.orekit.propagation.SpacecraftState;
41 import org.orekit.time.AbsoluteDate;
42 import org.orekit.time.FieldAbsoluteDate;
43 import org.orekit.time.TimeScalesFactory;
44 import org.orekit.utils.Constants;
45 import org.orekit.utils.FieldPVCoordinates;
46 import org.orekit.utils.IERSConventions;
47 import org.orekit.utils.PVCoordinates;
48
49
50
51
52
53 public class NeQuickItuTest {
54
55 @BeforeEach
56 public void setUp() {
57 Utils.setDataRoot("regular-data");
58 }
59
60 @Test
61 public void testMediumSolarActivityItu() {
62
63
64 final NeQuickItu model = new NeQuickItu(137.568737, TimeScalesFactory.getUTC());
65
66
67 Assertions.assertEquals(137.568737, model.getF107(), 1.0e-6);
68 Assertions.assertEquals("UTC", model.getUtc().getName());
69
70
71 final GeodeticPoint recP = new GeodeticPoint(FastMath.toRadians(-31.80), FastMath.toRadians(115.89), 12.78);
72 final GeodeticPoint satP = new GeodeticPoint(FastMath.toRadians(-14.31), FastMath.toRadians(124.09), 20100697.90);
73
74
75 final AbsoluteDate date = new AbsoluteDate(2018, 1, 2, 16, 0, 0, TimeScalesFactory.getUTC());
76
77
78 final double stec = model.stec(date, recP, satP);
79 Assertions.assertEquals(14.607, stec, 1.0e-3);
80 }
81
82 @Test
83 public void testFieldMediumSolarActivityItu() {
84 doTestFieldMediumSolarActivityItu(Binary64Field.getInstance());
85 }
86
87 private <T extends CalculusFieldElement<T>> void doTestFieldMediumSolarActivityItu(final Field<T> field) {
88
89
90 final T zero = field.getZero();
91
92
93 final NeQuickItu model = new NeQuickItu(137.568737, TimeScalesFactory.getUTC());
94
95
96 final FieldGeodeticPoint<T> recP = new FieldGeodeticPoint<>(zero.newInstance(FastMath.toRadians(-31.80)),
97 zero.newInstance(FastMath.toRadians(115.89)),
98 zero.newInstance(12.78));
99 final FieldGeodeticPoint<T> satP = new FieldGeodeticPoint<>(zero.newInstance(FastMath.toRadians(-14.31)),
100 zero.newInstance(FastMath.toRadians(124.09)),
101 zero.newInstance(20100697.90));
102
103
104 final FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field, 2018, 1, 2, 16, 0, 0, TimeScalesFactory.getUTC());
105
106
107 final T stec = model.stec(date, recP, satP);
108 Assertions.assertEquals(14.607, stec.getReal(), 1.0e-3);
109 }
110
111
112 @Test
113 public void testValidationItu63() {
114 doTestValidationItu(63, 4.36189);
115 }
116
117
118 @Test
119 public void testValidationItu128() {
120 doTestValidationItu(128, 11.83090);
121 }
122
123
124 @Test
125 public void testValidationItu193() {
126 doTestValidationItu(193, 20.85232);
127 }
128
129 private void doTestValidationItu(final double flux, final double expected) {
130
131
132 final NeQuickItu model = new NeQuickItu(flux, TimeScalesFactory.getUTC());
133
134
135 final GeodeticPoint recP = new GeodeticPoint(FastMath.toRadians(82.494293510492980204),
136 FastMath.toRadians(-62.340460202287275138),
137 78.107446296427042398);
138 final GeodeticPoint satP = new GeodeticPoint(FastMath.toRadians(54.445029415916160076),
139 FastMath.toRadians(-118.47006897550868132),
140 20370730.845002099872);
141
142
143 final AbsoluteDate date = new AbsoluteDate(2007, 4, 1, TimeScalesFactory.getUTC());
144
145
146 final double stec = model.stec(date, recP, satP);
147 Assertions.assertEquals(expected, stec, 1.0e-3);
148 }
149
150
151 @Test
152 public void testValidationItu63Field() {
153 doTestValidationItu(Binary64Field.getInstance(), 63, 4.36189);
154 }
155
156
157 @Test
158 public void testValidationItu128Field() {
159 doTestValidationItu(Binary64Field.getInstance(), 128, 11.83090);
160 }
161
162
163 @Test
164 public void testValidationItu193Field() {
165 doTestValidationItu(Binary64Field.getInstance(), 193, 20.85232);
166 }
167
168 private <T extends CalculusFieldElement<T>> void doTestValidationItu(final Field<T> field,
169 final double flux, final double expected) {
170
171
172 final NeQuickItu model = new NeQuickItu(flux, TimeScalesFactory.getUTC());
173
174
175 final FieldGeodeticPoint<T> recP =
176 new FieldGeodeticPoint<>(field,
177 new GeodeticPoint(FastMath.toRadians(82.494293510492980204),
178 FastMath.toRadians(-62.340460202287275138),
179 78.1074462964270423981));
180 final FieldGeodeticPoint<T> satP =
181 new FieldGeodeticPoint<>(field,
182 new GeodeticPoint(FastMath.toRadians(54.4450294159161600763),
183 FastMath.toRadians(-118.47006897550868132),
184 20370730.8450020998725));
185
186
187 final FieldAbsoluteDate<T> date =
188 new FieldAbsoluteDate<>(field, new AbsoluteDate(2007, 4, 1, TimeScalesFactory.getUTC()));
189
190
191 final T stec = model.stec(date, recP, satP);
192 Assertions.assertEquals(expected, stec.getReal(), 1.0e-3);
193 }
194
195 @Test
196 public void testHeights1() {
197 doTestHeights(30.0, 205000000.0, 10.730);
198 }
199
200 @Test
201 public void testHeights2() {
202 doTestHeights(1300000.0, 205000000.0, 0.943);
203 }
204
205 @Test
206 public void testHeights3() {
207 doTestHeights(2300000.0, 205000000.0, 0.580);
208 }
209
210 @Test
211 public void testHeights4() {
212 doTestHeights(30.0, 1950000.0, 29.670);
213 }
214
215 @Test
216 public void testHeights5() {
217 doTestHeights(1030000.0, 1950000.0, 2.594);
218 }
219
220 private void doTestHeights(final double hRec, final double hSat, final double expected) {
221
222
223 final NeQuickItu model = new NeQuickItu(128.0, TimeScalesFactory.getUTC());
224
225
226 final GeodeticPoint recP = new GeodeticPoint(FastMath.toRadians( 82.494),
227 FastMath.toRadians(-62.340),
228 hRec);
229 final GeodeticPoint satP = new GeodeticPoint(FastMath.toRadians(54.445),
230 FastMath.toRadians(-118.470),
231 hSat);
232
233
234 final AbsoluteDate date = new AbsoluteDate(2007, 4, 1, TimeScalesFactory.getUTC());
235
236
237 final double stec = model.stec(date, recP, satP);
238 Assertions.assertEquals(expected, stec, 1.0e-3);
239
240 }
241
242 @Test
243 public void testHeights1Field() {
244 doTestHeights(Binary64Field.getInstance(), 30.0, 205000000.0, 10.730);
245 }
246
247 @Test
248 public void testHeights2Field() {
249 doTestHeights(Binary64Field.getInstance(), 1300000.0, 205000000.0, 0.943);
250 }
251
252 @Test
253 public void testHeights3Field() {
254 doTestHeights(Binary64Field.getInstance(), 2300000.0, 205000000.0, 0.580);
255 }
256
257 @Test
258 public void testHeights4Field() {
259 doTestHeights(Binary64Field.getInstance(), 30.0, 1950000.0, 29.670);
260 }
261
262 @Test
263 public void testHeights5Field() {
264 doTestHeights(Binary64Field.getInstance(), 1030000.0, 1950000.0, 2.594);
265 }
266
267 private <T extends CalculusFieldElement<T>> void doTestHeights(final Field<T> field,
268 final double hRec, final double hSat, final double expected) {
269
270 final T zero = field.getZero();
271
272
273 final NeQuickItu model = new NeQuickItu(128.0, TimeScalesFactory.getUTC());
274
275
276 final FieldGeodeticPoint<T> recP = new FieldGeodeticPoint<>(FastMath.toRadians(zero.newInstance( 82.494)),
277 FastMath.toRadians(zero.newInstance(-62.340)),
278 zero.newInstance(hRec));
279 final FieldGeodeticPoint<T> satP = new FieldGeodeticPoint<>(FastMath.toRadians(zero.newInstance(54.445)),
280 FastMath.toRadians(zero.newInstance(-118.470)),
281 zero.newInstance(hSat));
282
283
284 final FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field, 2007, 4, 1, TimeScalesFactory.getUTC());
285
286
287 final T stec = model.stec(date, recP, satP);
288 Assertions.assertEquals(expected, stec.getReal(), 1.0e-3);
289
290 }
291
292 @Test
293 public void testDelay() {
294
295
296 final NeQuickItu model = new NeQuickItu(128.0, TimeScalesFactory.getUTC());
297
298
299 final GeodeticPoint recP = new GeodeticPoint(FastMath.toRadians(-31.80), FastMath.toRadians(115.89), 12.78);
300 final GeodeticPoint satP = new GeodeticPoint(FastMath.toRadians(-14.31), FastMath.toRadians(124.09), 20100697.90);
301
302
303 final AbsoluteDate date = new AbsoluteDate(2018, 4, 2, 16, 0, 0, TimeScalesFactory.getUTC());
304
305
306 final OneAxisEllipsoid ellipsoid = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
307 Constants.WGS84_EARTH_FLATTENING,
308 FramesFactory.getITRF(IERSConventions.IERS_2010, true));
309
310 final Vector3D satPosInITRF = ellipsoid.transform(satP);
311 final Vector3D satPosInEME2000 = ellipsoid.getBodyFrame().getStaticTransformTo(FramesFactory.getEME2000(), date).transformPosition(satPosInITRF);
312
313
314 final PVCoordinates pv = new PVCoordinates(satPosInEME2000, new Vector3D(1.0, 1.0, 1.0));
315 final Orbit orbit = new CartesianOrbit(pv, FramesFactory.getEME2000(), date, Constants.WGS84_EARTH_MU);
316 final SpacecraftState state = new SpacecraftState(orbit);
317
318 final double delay = model.pathDelay(state, new TopocentricFrame(ellipsoid, recP, null),
319 PredefinedGnssSignal.G01.getFrequency(), model.getParameters());
320
321
322 Assertions.assertEquals(1.32, delay, 0.01);
323 }
324
325 @Test
326 public void testFieldDelay() {
327 doTestFieldDelay(Binary64Field.getInstance());
328 }
329
330 private <T extends CalculusFieldElement<T>> void doTestFieldDelay(final Field<T> field) {
331
332
333 final T zero = field.getZero();
334 final T one = field.getOne();
335
336
337 final NeQuickItu model = new NeQuickItu(128.0, TimeScalesFactory.getUTC());
338
339
340 final double recLat = FastMath.toRadians(-31.80);
341 final double recLon = FastMath.toRadians(115.89);
342 final double recAlt = 12.78;
343 final GeodeticPoint recP = new GeodeticPoint(recLat, recLon, recAlt);
344 final FieldGeodeticPoint<T> satP = new FieldGeodeticPoint<>(zero.newInstance(FastMath.toRadians(-14.31)),
345 zero.newInstance(FastMath.toRadians(124.09)), zero.newInstance(20100697.90));
346
347
348 final FieldAbsoluteDate<T> date = new FieldAbsoluteDate<>(field, 2018, 4, 2, 16, 0, 0, TimeScalesFactory.getUTC());
349
350
351 final OneAxisEllipsoid ellipsoid = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
352 Constants.WGS84_EARTH_FLATTENING,
353 FramesFactory.getITRF(IERSConventions.IERS_2010, true));
354
355 final FieldVector3D<T> satPosInITRF = ellipsoid.transform(satP);
356 final FieldVector3D<T> satPosInEME2000 = ellipsoid.getBodyFrame().getStaticTransformTo(FramesFactory.getEME2000(), date).transformPosition(satPosInITRF);
357
358
359 final FieldPVCoordinates<T> pv = new FieldPVCoordinates<>(satPosInEME2000, new FieldVector3D<>(one, one, one));
360 final FieldOrbit<T> orbit = new FieldCartesianOrbit<>(pv, FramesFactory.getEME2000(), date, zero.newInstance(Constants.WGS84_EARTH_MU));
361 final FieldSpacecraftState<T> state = new FieldSpacecraftState<>(orbit);
362
363 final T delay = model.pathDelay(state, new TopocentricFrame(ellipsoid, recP, null),
364 PredefinedGnssSignal.G01.getFrequency(), model.getParameters(field));
365
366
367 Assertions.assertEquals(1.32, delay.getReal(), 0.01);
368 }
369
370 @Test
371 public void testAntiMeridian() {
372
373
374 final NeQuickItu model = new NeQuickItu(128.0, TimeScalesFactory.getUTC());
375
376
377 final AbsoluteDate date = new AbsoluteDate(2018, 11, 2, 16, 0, 0, TimeScalesFactory.getUTC());
378
379
380 final GeodeticPoint recP = new GeodeticPoint(FastMath.toRadians(-31.80), FastMath.toRadians(-179.99), 12.78);
381 final GeodeticPoint satP = new GeodeticPoint(FastMath.toRadians(-14.31), FastMath.toRadians(-177.43), 20100697.90);
382 double stec = model.stec(date, recP, satP);
383 Assertions.assertEquals(13.269, stec, 0.001);
384
385 }
386
387 @Test
388 public void testFieldAntiMeridian() {
389 doTestFieldAntiMeridian(Binary64Field.getInstance());
390 }
391
392 private <T extends CalculusFieldElement<T>> void doTestFieldAntiMeridian(final Field<T> field) {
393
394 final T zero = field.getZero();
395
396
397 final NeQuickItu model = new NeQuickItu(128.0, TimeScalesFactory.getUTC());
398
399
400 final FieldAbsoluteDate<T> date =
401 new FieldAbsoluteDate<>(field,
402 new AbsoluteDate(2018, 11, 2, 16, 0, 0, TimeScalesFactory.getUTC()));
403
404
405 final FieldGeodeticPoint<T> recP = new FieldGeodeticPoint<>(FastMath.toRadians(zero.newInstance(-31.80)),
406 FastMath.toRadians(zero.newInstance(-179.99)),
407 zero.newInstance(12.78));
408 final FieldGeodeticPoint<T> satP = new FieldGeodeticPoint<>(FastMath.toRadians(zero.newInstance(-14.31)),
409 FastMath.toRadians(zero.newInstance(-177.43)),
410 zero.newInstance(20100697.90));
411 T stec = model.stec(date, recP, satP);
412 Assertions.assertEquals(13.269, stec.getReal(), 0.001);
413
414 }
415
416 @Test
417 public void testPole() {
418
419
420 final NeQuickItu model = new NeQuickItu(128.0, TimeScalesFactory.getUTC());
421
422
423 final AbsoluteDate date = new AbsoluteDate(2018, 11, 2, 16, 0, 0, TimeScalesFactory.getUTC());
424
425
426 final GeodeticPoint recP = new GeodeticPoint(FastMath.toRadians(90), FastMath.toRadians(0), 100);
427 final GeodeticPoint satP = new GeodeticPoint(FastMath.toRadians(-14.31), FastMath.toRadians(-177.43), 20100697.90);
428 double stec = model.stec(date, recP, satP);
429 Assertions.assertEquals(70.372, stec, 0.001);
430
431 }
432
433 @Test
434 public void testFieldPole() {
435 doTestFieldPole(Binary64Field.getInstance());
436 }
437
438 private <T extends CalculusFieldElement<T>> void doTestFieldPole(final Field<T> field) {
439
440 final T zero = field.getZero();
441
442
443 final NeQuickItu model = new NeQuickItu(128.0, TimeScalesFactory.getUTC());
444
445
446 final FieldAbsoluteDate<T> date =
447 new FieldAbsoluteDate<>(field,
448 new AbsoluteDate(2018, 11, 2, 16, 0, 0, TimeScalesFactory.getUTC()));
449
450
451 final FieldGeodeticPoint<T> recP = new FieldGeodeticPoint<>(FastMath.toRadians(zero.newInstance(90)),
452 FastMath.toRadians(zero.newInstance(0)),
453 zero.newInstance(12.78));
454 final FieldGeodeticPoint<T> satP = new FieldGeodeticPoint<>(FastMath.toRadians(zero.newInstance(-14.31)),
455 FastMath.toRadians(zero.newInstance(-177.43)),
456 zero.newInstance(20100697.90));
457 T stec = model.stec(date, recP, satP);
458 Assertions.assertEquals(70.373, stec.getReal(), 0.001);
459
460 }
461
462 @Test
463 public void testZenith() {
464
465
466 final NeQuickItu model = new NeQuickItu(128.0, TimeScalesFactory.getUTC());
467
468
469 final AbsoluteDate date = new AbsoluteDate(2018, 4, 2, 16, 0, 0, TimeScalesFactory.getUTC());
470
471
472 final GeodeticPoint recP = new GeodeticPoint(FastMath.toRadians(51.678), FastMath.toRadians(-9.448), 0.0);
473 final GeodeticPoint satP = new GeodeticPoint(FastMath.toRadians(51.678), FastMath.toRadians(-9.448), 600000.0);
474 double stec = model.stec(date, recP, satP);
475 Assertions.assertEquals(16.618, stec, 0.001);
476
477 }
478
479 @Test
480 public void testFieldZenith() {
481 doTestFieldZenith(Binary64Field.getInstance());
482 }
483
484 private <T extends CalculusFieldElement<T>> void doTestFieldZenith(final Field<T> field) {
485
486 final T zero = field.getZero();
487
488
489 final NeQuickItu model = new NeQuickItu(128.0, TimeScalesFactory.getUTC());
490
491
492 final FieldAbsoluteDate<T> date =
493 new FieldAbsoluteDate<>(field,
494 new AbsoluteDate(2018, 4, 2, 16, 0, 0, TimeScalesFactory.getUTC()));
495
496
497 final FieldGeodeticPoint<T> recP = new FieldGeodeticPoint<>(FastMath.toRadians(zero.newInstance(51.678)),
498 FastMath.toRadians(zero.newInstance(-9.448)),
499 zero.newInstance(0.0));
500 final FieldGeodeticPoint<T> satP = new FieldGeodeticPoint<>(FastMath.toRadians(zero.newInstance(51.678)),
501 FastMath.toRadians(zero.newInstance(-9.448)),
502 zero.newInstance(600000.0));
503 T stec = model.stec(date, recP, satP);
504 Assertions.assertEquals(16.618, stec.getReal(), 0.001);
505
506 }
507
508 }