1 /*
2 * Copyright (C) 2020 Alberto Irurueta Carro (alberto@irurueta.com)
3 *
4 * Licensed under the Apache License, Version 2.0 (the "License");
5 * you may not use this file except in compliance with the License.
6 * You may obtain a copy of the License at
7 *
8 * http://www.apache.org/licenses/LICENSE-2.0
9 *
10 * Unless required by applicable law or agreed to in writing, software
11 * distributed under the License is distributed on an "AS IS" BASIS,
12 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13 * See the License for the specific language governing permissions and
14 * limitations under the License.
15 */
16 package com.irurueta.navigation.inertial.wmm;
17
18 import com.irurueta.navigation.frames.NEDPosition;
19 import com.irurueta.navigation.geodesic.Constants;
20 import com.irurueta.units.Angle;
21 import com.irurueta.units.AngleConverter;
22 import com.irurueta.units.AngleUnit;
23 import com.irurueta.units.Distance;
24 import com.irurueta.units.DistanceConverter;
25 import com.irurueta.units.DistanceUnit;
26
27 import java.io.IOException;
28 import java.util.Calendar;
29 import java.util.Date;
30 import java.util.GregorianCalendar;
31
32 /**
33 * Estimates Earth magnetic flux density resolved around NED frame at
34 * a given Earth location.
35 */
36 public class WMMEarthMagneticFluxDensityEstimator {
37
38 /**
39 * Guaranteed accuracy of estimated angles by the WMM during
40 * the valid timespan of a model.
41 * Accuracy is expressed in degrees and refers to estimated
42 * declination and dip angles.
43 */
44 public static final double ANGLE_ACCURACY_DEGREES = 5e-3;
45
46 /**
47 * Guaranteed accuracy of estimated angles by the WMM during
48 * the valid timespan of a model.
49 * Accuracy is expressed in radians and refers to estimated
50 * declination and dip angles.
51 */
52 public static final double ANGLE_ACCURACY_RADIANS = Math.toRadians(ANGLE_ACCURACY_DEGREES);
53
54 /**
55 * Guaranteed accuracy of estimated magnetic flux density
56 * values by the WMM during the valid timespan of a model.
57 * This value refers to intensity norm, vertical intensity,
58 * horizontal intensity, north intensity and east intensity.
59 */
60 public static final double INTENSITY_ACCURACY = 5e-2;
61
62 /**
63 * Number of coefficients.
64 */
65 private static final int N = WorldMagneticModel.N;
66
67 /**
68 * Coefficients file.
69 */
70 private static final String COEFFICIENTS_FILE = "wmm.cof";
71
72 /**
73 * Mean radius of IAU-66 ellipsoid expressed in Km.
74 */
75 private static final double RE_KM = 6371.2;
76
77 /**
78 * Converts to nanos.
79 */
80 private static final double FROM_NANO = 1e-9;
81
82 /**
83 * Time value used in the previous calculation.
84 * This is used to save on calculation time if some
85 * inputs don't change.
86 * Old time is expressed in decimal years.
87 */
88 private Double oldTime;
89
90 /**
91 * Geodetic height (a.k.a. altitude) value used in the previous
92 * calculation.
93 * This is used to save on calculation time if some
94 * inputs don't change.
95 * Old height is expressed in Kilometers (Km).
96 */
97 private Double oldHeight;
98
99 /**
100 * Old geodetic latitude value used in previous
101 * calculation.
102 * This is used to save on calculation time if some
103 * inputs don't change.
104 * Old latitude is expressed in degrees (deg).
105 */
106 private Double oldLatitude;
107
108 /**
109 * Old geodetic longitude value used in previous
110 * calculation.
111 * This is used to save on calculation time if some
112 * inputs don't change.
113 * Old longitude is expressed in degrees (deg).
114 */
115 private Double oldLongitude;
116
117 /**
118 * Geomagnetic declination in degrees.
119 * East is positive, West is negative.
120 * (The negative of variation).
121 */
122 private double dec;
123
124 /**
125 * Geomagnetic inclination in degrees.
126 * Down is positive, up is negative.
127 */
128 private double dip;
129
130 /**
131 * Geomagnetic total intensity expressed in nano Teslas (nT)
132 */
133 private double ti;
134
135 /**
136 * A World Magnetic Model containing all required coefficients.
137 */
138 private final WorldMagneticModel model;
139
140 /**
141 * The time-adjusted geomagnetic gauss coefficients (nt).
142 */
143 private final double[][] tc = new double[N][N];
144
145 /**
146 * The theta derivative of p(n,m) (un-normalized).
147 */
148 private final double[][] dp = new double[N][N];
149
150 /**
151 * The sine of (m*spherical coord. longitude).
152 */
153 private final double[] sp = new double[N];
154
155 /**
156 * The cosine of (m*spherical coord. longitude).
157 */
158 private final double[] cp = new double[N];
159
160 /**
161 * The associated Legendre polynomials for m=1 (un-normalized).
162 */
163 private final double[] pp = new double[N];
164
165 /**
166 * The north-south field intensity expressed in nano Teslas (nT).
167 */
168 private double bx;
169
170 /**
171 * The east-west field intensity expressed in nano Teslas (nT).
172 */
173 private double by;
174
175 /**
176 * The vertical field intensity positive downward expressed
177 * in nano Teslas (nT)
178 */
179 private double bz;
180
181 /**
182 * The horizontal field intensity expressed in nano Teslas (nT)
183 */
184 private double bh;
185
186 /**
187 * Semi-major axis of WGS-84 ellipsoid, in Km, squared.
188 */
189 private final double a2;
190
191 /**
192 * Semi-minor axis of WGS-84 ellipsoid, in Km, squared.
193 */
194 private final double b2;
195
196 /**
197 * The difference between the squared semi-axes.
198 * c2 = a2 - b2
199 */
200 private final double c2;
201
202 /**
203 * {@link #a2} squared.
204 */
205 private final double a4;
206
207 /**
208 * The difference between a4 and b4
209 * c4 = a4 - b4
210 */
211 private final double c4;
212
213 // Below: internal values being reused. These values are only
214 // recalculated if the height (altitude) changes.
215
216 /**
217 * This is the magnetic field strength at the magnetic pole.
218 */
219 private double r;
220
221 /**
222 * This is the cosine of the magnetic inclination at the magnetic pole.
223 */
224 private double ca;
225
226 /**
227 * This is the sine of the magnetic inclination at the magnetic pole.
228 */
229 private double sa;
230
231 /**
232 * This is the cosine of the magnetic declination at the magnetic pole.
233 */
234 private double ct;
235
236 /**
237 * This is the sine of the magnetic declination at the magnetic pole.
238 */
239 private double st;
240
241 /**
242 * Constructor.
243 *
244 * @throws IOException if an I/O error occurs while loading
245 * model coefficients.
246 */
247 public WMMEarthMagneticFluxDensityEstimator() throws IOException {
248 this(WMMLoader.loadFromResource(COEFFICIENTS_FILE));
249 }
250
251 /**
252 * Constructor.
253 *
254 * @param model a World Magnetic Model.
255 * @throws NullPointerException if provided model is null.
256 */
257 public WMMEarthMagneticFluxDensityEstimator(final WorldMagneticModel model) {
258 if (model == null) {
259 throw new NullPointerException();
260 }
261 this.model = model;
262 sp[0] = 0.0;
263 cp[0] = 1.0;
264 pp[0] = 1.0;
265
266 // semi-major axis of WGS-84 ellipsoid, in Km (6378.137 Km).
267 final var a = DistanceConverter.convert(Constants.EARTH_EQUATORIAL_RADIUS_WGS84, DistanceUnit.METER,
268 DistanceUnit.KILOMETER);
269
270 // semi-minor axis of WGS-84 ellipsoid, in Km (6356.7523142 Km).
271 final var b = DistanceConverter.convert(Constants.EARTH_POLAR_RADIUS_WGS84, DistanceUnit.METER,
272 DistanceUnit.KILOMETER);
273 a2 = a * a;
274 b2 = b * b;
275 c2 = a2 - b2;
276 a4 = a2 * a2;
277
278 final var b4 = b2 * b2;
279 c4 = a4 - b4;
280 }
281
282 /**
283 * Gets World Magnetic Model containing all required coefficients.
284 *
285 * @return World Magnetic Model.
286 */
287 public WorldMagneticModel getModel() {
288 return model;
289 }
290
291 /**
292 * Returns the declination from the Department of Defense geomagnetic
293 * model and data, in radians.
294 * The magnetic heading + declination is the true heading of a device
295 * in terms of geographical north pole.
296 * This method uses default altitude (0.0 - mean sea level) and time
297 * (half way through the valid 5 year period of the model).
298 *
299 * @param latitude latitude expressed in radians.
300 * @param longitude longitude expressed in radians.
301 * @return magnetic field declination expressed in radians.
302 */
303 public double getDeclination(final double latitude, final double longitude) {
304 final var defaultTime = model.epoch + WorldMagneticModel.LIFESPAN / 2.0;
305 calcGeoMag(Math.toDegrees(latitude), Math.toDegrees(longitude), 0.0, defaultTime);
306 return Math.toRadians(dec);
307 }
308
309 /**
310 * Returns the declination from the Department of Defense geomagnetic
311 * model and data, in radians.
312 * The magnetic heading + declination is the true heading of a device
313 * in terms of geographical north pole.
314 * This method uses default altitude (0.0 - mean sea level) and time
315 * (half way through the valid 5 year period of the model).
316 *
317 * @param latitude latitude angle.
318 * @param longitude longitude angle.
319 * @return magnetic field declination expressed in radians.
320 */
321 public double getDeclination(final Angle latitude, final Angle longitude) {
322 return getDeclination(convertAngle(latitude), convertAngle(longitude));
323 }
324
325 /**
326 * Returns the declination from the Department of Defense geomagnetic
327 * model and data, in radians.
328 * The magnetic heading + declination is the true heading of a device
329 * in terms of geographical north pole.
330 * This method uses default altitude (0.0 - mean sea level) and time
331 * (half way through the valid 5 year period of the model).
332 *
333 * @param latitude latitude expressed in radians.
334 * @param longitude longitude expressed in radians.
335 * @param result instance where magnetic field declination will be
336 * stored.
337 */
338 public void getDeclinationAsAngle(final double latitude, final double longitude, final Angle result) {
339 result.setUnit(AngleUnit.RADIANS);
340 result.setValue(getDeclination(latitude, longitude));
341 }
342
343 /**
344 * Returns the declination from the Department of Defense geomagnetic
345 * model and data, in radians.
346 * The magnetic heading + declination is the true heading of a device
347 * in terms of geographical north pole.
348 * This method uses default altitude (0.0 - mean sea level) and time
349 * (half way through the valid 5 year period of the model).
350 *
351 * @param latitude latitude angle.
352 * @param longitude longitude angle.
353 * @return magnetic field declination.
354 */
355 public Angle getDeclinationAsAngle(final double latitude, final double longitude) {
356 return new Angle(getDeclination(latitude, longitude), AngleUnit.RADIANS);
357 }
358
359 /**
360 * Returns the declination from the Department of Defense geomagnetic
361 * model and data, in radians.
362 * The magnetic heading + declination is the true heading of a device
363 * in terms of geographical north pole.
364 * This method uses default altitude (0.0 - mean sea level) and time
365 * (half way through the valid 5 year period of the model).
366 *
367 * @param latitude latitude angle.
368 * @param longitude longitude angle.
369 * @param result instance where magnetic field declination will be
370 * stored.
371 */
372 public void getDeclinationAsAngle(final Angle latitude, final Angle longitude, final Angle result) {
373 result.setUnit(AngleUnit.RADIANS);
374 result.setValue(getDeclination(latitude, longitude));
375 }
376
377 /**
378 * Returns the declination from the Department of Defense geomagnetic
379 * model and data, in radians.
380 * The magnetic heading + declination is the true heading of a device
381 * in terms of geographical north pole.
382 * This method uses default altitude (0.0 - mean sea level) and time
383 * (half way through the valid 5 year period of the model).
384 *
385 * @param latitude latitude angle.
386 * @param longitude longitude angle.
387 * @return magnetic field declination.
388 */
389 public Angle getDeclinationAsAngle(final Angle latitude, final Angle longitude) {
390 return new Angle(getDeclination(latitude, longitude), AngleUnit.RADIANS);
391 }
392
393 /**
394 * Returns the declination from the Department of Defense geomagnetic
395 * model and data, in radians.
396 * The magnetic heading declination is the true heading of a device
397 * in terms of geographical north pole.
398 *
399 * @param latitude latitude expressed in radians.
400 * @param longitude longitude expressed in radians.
401 * @param height height expressed in meters.
402 * @param year year expressed in decimal years.
403 * @return magnetic field declination expressed in radians.
404 */
405 public double getDeclination(
406 final double latitude, final double longitude, final double height, final double year) {
407 final var heightKm = DistanceConverter.convert(height, DistanceUnit.METER, DistanceUnit.KILOMETER);
408 calcGeoMag(Math.toDegrees(latitude), Math.toDegrees(longitude), heightKm, year);
409 return Math.toRadians(dec);
410 }
411
412 /**
413 * Returns the declination from the Department of Defense geomagnetic
414 * model and data, in radians.
415 * <p>
416 * The magnetic heading + declination is the true heading of a device
417 * in terms of geographical north pole.
418 *
419 * @param latitude latitude expressed in radians.
420 * @param longitude longitude expressed in radians.
421 * @param height height expressed in meters.
422 * @param calendar a calendar containing a specific instant.
423 * @return magnetic field declination expressed in radians.
424 */
425 public double getDeclination(
426 final double latitude, final double longitude, final double height, final GregorianCalendar calendar) {
427 return getDeclination(latitude, longitude, height, convertTime(calendar));
428 }
429
430 /**
431 * Returns the declination from the Department of Defense geomagnetic
432 * model and data, in radians.
433 * <p>
434 * The magnetic heading + declination is the true heading of a device
435 * in terms of geographical north pole.
436 *
437 * @param latitude latitude expressed in radians.
438 * @param longitude longitude expressed in radians.
439 * @param height height expressed in meters.
440 * @param time a specific time instant.
441 * @return magnetic field declination expressed in radians.
442 */
443 public double getDeclination(
444 final double latitude, final double longitude, final double height, final Date time) {
445 final var calendar = new GregorianCalendar();
446 calendar.setTime(time);
447 return getDeclination(latitude, longitude, height, calendar);
448 }
449
450 /**
451 * Returns the declination from the Department of Defense geomagnetic
452 * model and data, in radians.
453 * The magnetic heading + declination is the true heading of a device
454 * in terms of geographical north pole.
455 *
456 * @param latitude latitude angle.
457 * @param longitude longitude angle.
458 * @param height height.
459 * @param year year expressed in decimal years.
460 * @return magnetic field declination expressed in radians.
461 */
462 public double getDeclination(
463 final Angle latitude, final Angle longitude, final Distance height, final double year) {
464 return getDeclination(convertAngle(latitude), convertAngle(longitude), convertDistance(height), year);
465 }
466
467 /**
468 * Returns the declination from the Department of Defense geomagnetic
469 * model and data, in radians.
470 * <p>
471 * The magnetic heading + declination is the true heading of a device
472 * in terms of geographical north pole.
473 *
474 * @param latitude latitude angle.
475 * @param longitude longitude angle.
476 * @param height height.
477 * @param calendar a calendar containing a specific instant.
478 * @return magnetic field declination expressed in radians.
479 */
480 public double getDeclination(
481 final Angle latitude, final Angle longitude, final Distance height, final GregorianCalendar calendar) {
482 return getDeclination(convertAngle(latitude), convertAngle(longitude), convertDistance(height), calendar);
483 }
484
485 /**
486 * Returns the declination from the Department of Defense geomagnetic
487 * model and data, in radians.
488 * <p>
489 * The magnetic heading + declination is the true heading of a device
490 * in terms of geographical north pole.
491 *
492 * @param latitude latitude angle.
493 * @param longitude longitude angle.
494 * @param height height.
495 * @param time a specific time instant.
496 * @return magnetic field declination expressed in radians.
497 */
498 public double getDeclination(
499 final Angle latitude, final Angle longitude, final Distance height, final Date time) {
500 return getDeclination(convertAngle(latitude), convertAngle(longitude), convertDistance(height), time);
501 }
502
503 /**
504 * Returns the declination from the Department of Defense geomagnetic
505 * model and data, in radians.
506 * The magnetic heading + declination is the true heading of a device
507 * in terms of geographical north pole.
508 *
509 * @param position a position expressed in geodetic coordinates.
510 * @param year year expressed in decimal years.
511 * @return magnetic field declination expressed in radians.
512 */
513 public double getDeclination(final NEDPosition position, final double year) {
514 return getDeclination(position.getLatitude(), position.getLongitude(), position.getHeight(), year);
515 }
516
517 /**
518 * Returns the declination from the Department of Defense geomagnetic
519 * model and data, in radians.
520 * The magnetic heading + declination is the true heading of a device
521 * in terms of geographical north pole.
522 *
523 * @param position a position expressed in geodetic coordinates.
524 * @param calendar a calendar containing a specific instant.
525 * @return magnetic field declination expressed in radians.
526 */
527 public double getDeclination(final NEDPosition position, final GregorianCalendar calendar) {
528 return getDeclination(position.getLatitude(), position.getLongitude(), position.getHeight(), calendar);
529 }
530
531 /**
532 * Returns the declination from the Department of Defense geomagnetic
533 * model and data, in radians.
534 * The magnetic heading + declination is the true heading of a device
535 * in terms of geographical north pole.
536 *
537 * @param position a position expressed in geodetic coordinates.
538 * @param time a specific time instant.
539 * @return magnetic field declination expressed in radians.
540 */
541 public double getDeclination(final NEDPosition position, final Date time) {
542 return getDeclination(position.getLatitude(), position.getLongitude(), position.getHeight(), time);
543 }
544
545 /**
546 * Returns the declination from the Department of Defense geomagnetic
547 * model and data, in radians.
548 * The magnetic heading + declination is the true heading of a device
549 * in terms of geographical north pole.
550 *
551 * @param latitude latitude expressed in radians.
552 * @param longitude longitude expressed in radians.
553 * @param height height expressed in meters.
554 * @param year year expressed in decimal years.
555 * @param result instance where magnetic field declination will be
556 * stored.
557 */
558 public void getDeclinationAsAngle(final double latitude, final double longitude, final double height,
559 final double year, final Angle result) {
560 result.setUnit(AngleUnit.RADIANS);
561 result.setValue(getDeclination(latitude, longitude, height, year));
562 }
563
564 /**
565 * Returns the declination from the Department of Defense geomagnetic
566 * model and data, in radians.
567 * The magnetic heading + declination is the true heading of a device
568 * in terms of geographical north pole.
569 *
570 * @param latitude latitude expressed in radians.
571 * @param longitude longitude expressed in radians.
572 * @param height height expressed in meters.
573 * @param year year expressed in decimal years.
574 * @return magnetic field declination.
575 */
576 public Angle getDeclinationAsAngle(
577 final double latitude, final double longitude, final double height, final double year) {
578 return new Angle(getDeclination(latitude, longitude, height, year), AngleUnit.RADIANS);
579 }
580
581 /**
582 * Returns the declination from the Department of Defense geomagnetic
583 * model and data, in radians.
584 * <p>
585 * The magnetic heading + declination is the true heading of a device
586 * in terms of geographical north pole.
587 *
588 * @param latitude latitude expressed in radians.
589 * @param longitude longitude expressed in radians.
590 * @param height height expressed in meters.
591 * @param calendar a calendar containing a specific instant.
592 * @param result instance where magnetic field declination will be
593 * stored.
594 */
595 public void getDeclinationAsAngle(
596 final double latitude, final double longitude, final double height, final GregorianCalendar calendar,
597 final Angle result) {
598 result.setUnit(AngleUnit.RADIANS);
599 result.setValue(getDeclination(latitude, longitude, height, calendar));
600 }
601
602 /**
603 * Returns the declination from the Department of Defense geomagnetic
604 * model and data, in radians.
605 * <p>
606 * The magnetic heading + declination is the true heading of a device
607 * in terms of geographical north pole.
608 *
609 * @param latitude latitude expressed in radians.
610 * @param longitude longitude expressed in radians.
611 * @param height height expressed in meters.
612 * @param calendar a calendar containing a specific instant.
613 * @return magnetic field declination.
614 */
615 public Angle getDeclinationAsAngle(
616 final double latitude, final double longitude, final double height, final GregorianCalendar calendar) {
617 return new Angle(getDeclination(latitude, longitude, height, calendar), AngleUnit.RADIANS);
618 }
619
620 /**
621 * Returns the declination from the Department of Defense geomagnetic
622 * model and data, in radians.
623 * <p>
624 * The magnetic heading + declination is the true heading of a device
625 * in terms of geographical north pole.
626 *
627 * @param latitude latitude expressed in radians.
628 * @param longitude longitude expressed in radians.
629 * @param height height expressed in meters.
630 * @param time a specific time instant.
631 * @param result instance where magnetic field declination will be
632 * stored.
633 */
634 public void getDeclinationAsAngle(
635 final double latitude, final double longitude, final double height, final Date time, final Angle result) {
636 result.setUnit(AngleUnit.RADIANS);
637 result.setValue(getDeclination(latitude, longitude, height, time));
638 }
639
640 /**
641 * Returns the declination from the Department of Defense geomagnetic
642 * model and data, in radians.
643 * <p>
644 * The magnetic heading + declination is the true heading of a device
645 * in terms of geographical north pole.
646 *
647 * @param latitude latitude expressed in radians.
648 * @param longitude longitude expressed in radians.
649 * @param height height expressed in meters.
650 * @param time a specific time instant.
651 * @return magnetic field declination.
652 */
653 public Angle getDeclinationAsAngle(
654 final double latitude, final double longitude, final double height, final Date time) {
655 return new Angle(getDeclination(latitude, longitude, height, time), AngleUnit.RADIANS);
656 }
657
658 /**
659 * Returns the declination from the Department of Defense geomagnetic
660 * model and data, in radians.
661 * The magnetic heading + declination is the true heading of a device
662 * in terms of geographical north pole.
663 *
664 * @param latitude latitude angle.
665 * @param longitude longitude angle.
666 * @param height height.
667 * @param year year expressed in decimal years.
668 * @param result instance where magnetic field declination will be
669 * stored.
670 */
671 public void getDeclinationAsAngle(
672 final Angle latitude, final Angle longitude, final Distance height, final double year, final Angle result) {
673 result.setUnit(AngleUnit.RADIANS);
674 result.setValue(getDeclination(latitude, longitude, height, year));
675 }
676
677 /**
678 * Returns the declination from the Department of Defense geomagnetic
679 * model and data, in radians.
680 * The magnetic heading + declination is the true heading of a device
681 * in terms of geographical north pole.
682 *
683 * @param latitude latitude angle.
684 * @param longitude longitude angle.
685 * @param height height.
686 * @param year year expressed in decimal years.
687 * @return magnetic field declination.
688 */
689 public Angle getDeclinationAsAngle(
690 final Angle latitude, final Angle longitude, final Distance height, final double year) {
691 return new Angle(getDeclination(latitude, longitude, height, year), AngleUnit.RADIANS);
692 }
693
694 /**
695 * Returns the declination from the Department of Defense geomagnetic
696 * model and data, in radians.
697 * <p>
698 * The magnetic heading + declination is the true heading of a device
699 * in terms of geographical north pole.
700 *
701 * @param latitude latitude angle.
702 * @param longitude longitude angle.
703 * @param height height.
704 * @param calendar a calendar containing a specific instant.
705 * @param result instance where magnetic field declination will be
706 * stored.
707 */
708 public void getDeclinationAsAngle(
709 final Angle latitude, final Angle longitude, final Distance height, final GregorianCalendar calendar,
710 final Angle result) {
711 result.setUnit(AngleUnit.RADIANS);
712 result.setValue(getDeclination(latitude, longitude, height, calendar));
713 }
714
715 /**
716 * Returns the declination from the Department of Defense geomagnetic
717 * model and data, in radians.
718 * <p>
719 * The magnetic heading + declination is the true heading of a device
720 * in terms of geographical north pole.
721 *
722 * @param latitude latitude angle.
723 * @param longitude longitude angle.
724 * @param height height.
725 * @param calendar a calendar containing a specific instant.
726 * @return magnetic field declination.
727 */
728 public Angle getDeclinationAsAngle(
729 final Angle latitude, final Angle longitude, final Distance height, final GregorianCalendar calendar) {
730 return new Angle(getDeclination(latitude, longitude, height, calendar), AngleUnit.RADIANS);
731 }
732
733 /**
734 * Returns the declination from the Department of Defense geomagnetic
735 * model and data, in radians.
736 * <p>
737 * The magnetic heading + declination is the true heading of a device
738 * in terms of geographical north pole.
739 *
740 * @param latitude latitude angle.
741 * @param longitude longitude angle.
742 * @param height height.
743 * @param time a specific time instant.
744 * @param result instance where magnetic field declination will be
745 * stored.
746 */
747 public void getDeclinationAsAngle(
748 final Angle latitude, final Angle longitude, final Distance height, final Date time, final Angle result) {
749 result.setUnit(AngleUnit.RADIANS);
750 result.setValue(getDeclination(latitude, longitude, height, time));
751 }
752
753 /**
754 * Returns the declination from the Department of Defense geomagnetic
755 * model and data, in radians.
756 * <p>
757 * The magnetic heading + declination is the true heading of a device
758 * in terms of geographical north pole.
759 *
760 * @param latitude latitude angle.
761 * @param longitude longitude angle.
762 * @param height height.
763 * @param time a specific time instant.
764 * @return magnetic field declination.
765 */
766 public Angle getDeclinationAsAngle(
767 final Angle latitude, final Angle longitude, final Distance height, final Date time) {
768 return new Angle(getDeclination(latitude, longitude, height, time), AngleUnit.RADIANS);
769 }
770
771 /**
772 * Returns the declination from the Department of Defense geomagnetic
773 * model and data, in radians.
774 * The magnetic heading + declination is the true heading of a device
775 * in terms of geographical north pole.
776 *
777 * @param position a position expressed in geodetic coordinates.
778 * @param year year expressed in decimal years.
779 * @param result instance where magnetic field declination will be
780 * stored.
781 */
782 public void getDeclinationAsAngle(final NEDPosition position, final double year, final Angle result) {
783 result.setUnit(AngleUnit.RADIANS);
784 result.setValue(getDeclination(position, year));
785 }
786
787 /**
788 * Returns the declination from the Department of Defense geomagnetic
789 * model and data, in radians.
790 * The magnetic heading + declination is the true heading of a device
791 * in terms of geographical north pole.
792 *
793 * @param position a position expressed in geodetic coordinates.
794 * @param year year expressed in decimal years.
795 * @return magnetic field declination.
796 */
797 public Angle getDeclinationAsAngle(final NEDPosition position, final double year) {
798 return new Angle(getDeclination(position, year), AngleUnit.RADIANS);
799 }
800
801 /**
802 * Returns the declination from the Department of Defense geomagnetic
803 * model and data, in radians.
804 * The magnetic heading + declination is the true heading of a device
805 * in terms of geographical north pole.
806 *
807 * @param position a position expressed in geodetic coordinates.
808 * @param calendar a calendar containing a specific instant.
809 * @param result instance where magnetic field declination will be
810 * stored.
811 */
812 public void getDeclinationAsAngle(
813 final NEDPosition position, final GregorianCalendar calendar, final Angle result) {
814 result.setUnit(AngleUnit.RADIANS);
815 result.setValue(getDeclination(position, calendar));
816 }
817
818 /**
819 * Returns the declination from the Department of Defense geomagnetic
820 * model and data, in radians.
821 * The magnetic heading + declination is the true heading of a device
822 * in terms of geographical north pole.
823 *
824 * @param position a position expressed in geodetic coordinates.
825 * @param calendar a calendar containing a specific instant.
826 * @return magnetic field declination.
827 */
828 public Angle getDeclinationAsAngle(
829 final NEDPosition position, final GregorianCalendar calendar) {
830 return new Angle(getDeclination(position, calendar), AngleUnit.RADIANS);
831 }
832
833 /**
834 * Returns the declination from the Department of Defense geomagnetic
835 * model and data, in radians.
836 * The magnetic heading + declination is the true heading of a device
837 * in terms of geographical north pole.
838 *
839 * @param position a position expressed in geodetic coordinates.
840 * @param time a specific time instant.
841 * @param result instance where magnetic field declination will be
842 * stored.
843 */
844 public void getDeclinationAsAngle(final NEDPosition position, final Date time, final Angle result) {
845 result.setUnit(AngleUnit.RADIANS);
846 result.setValue(getDeclination(position, time));
847 }
848
849 /**
850 * Returns the declination from the Department of Defense geomagnetic
851 * model and data, in radians.
852 * The magnetic heading + declination is the true heading of a device
853 * in terms of geographical north pole.
854 *
855 * @param position a position expressed in geodetic coordinates.
856 * @param time a specific time instant.
857 * @return magnetic field declination.
858 */
859 public Angle getDeclinationAsAngle(final NEDPosition position, final Date time) {
860 return new Angle(getDeclination(position, time), AngleUnit.RADIANS);
861 }
862
863 /**
864 * Returns the magnetic field dip angle from the Department of
865 * Defense geomagnetic model and data, in radians.
866 * This method uses default altitude (0.0 - mean sea level) and time
867 * (half way through the valid 5 year period of the model).
868 *
869 * @param latitude latitude expressed in radians.
870 * @param longitude longitude expressed in radians.
871 * @return the magnetic field dip angle expressed in radians.
872 */
873 public double getDip(final double latitude, final double longitude) {
874 final var defaultTime = model.epoch + WorldMagneticModel.LIFESPAN / 2.0;
875 calcGeoMag(Math.toDegrees(latitude), Math.toDegrees(longitude), 0.0, defaultTime);
876 return Math.toRadians(dip);
877 }
878
879 /**
880 * Returns the magnetic field dip angle from the Department of
881 * Defense geomagnetic model and data, in radians.
882 * This method uses default altitude (0.0 - mean sea level) and time
883 * (half way through the valid 5 year period of the model).
884 *
885 * @param latitude latitude angle.
886 * @param longitude longitude angle.
887 * @return the magnetic field dip angle expressed in radians.
888 */
889 public double getDip(final Angle latitude, final Angle longitude) {
890 return getDip(convertAngle(latitude), convertAngle(longitude));
891 }
892
893 /**
894 * Returns the magnetic field dip angle from the Department of
895 * Defense geomagnetic model and data, in radians.
896 * This method uses default altitude (0.0 - mean sea level) and time
897 * (half way through the valid 5 year period of the model).
898 *
899 * @param latitude latitude expressed in radians.
900 * @param longitude longitude expressed in radians.
901 * @param result instance where magnetic field dip angle will be
902 * stored.
903 */
904 public void getDipAsAngle(final double latitude, final double longitude, final Angle result) {
905 result.setUnit(AngleUnit.RADIANS);
906 result.setValue(getDip(latitude, longitude));
907 }
908
909 /**
910 * Returns the magnetic field dip angle from the Department of
911 * Defense geomagnetic model and data, in radians.
912 * This method uses default altitude (0.0 - mean sea level) and time
913 * (half way through the valid 5 year period of the model).
914 *
915 * @param latitude latitude angle.
916 * @param longitude longitude angle.
917 * @return the magnetic field dip angle.
918 */
919 public Angle getDipAsAngle(final double latitude, final double longitude) {
920 return new Angle(getDip(latitude, longitude), AngleUnit.RADIANS);
921 }
922
923 /**
924 * Returns the magnetic field dip angle from the Department of
925 * Defense geomagnetic model and data, in radians.
926 * This method uses default altitude (0.0 - mean sea level) and time
927 * (half way through the valid 5 year period of the model).
928 *
929 * @param latitude latitude angle.
930 * @param longitude longitude angle.
931 * @param result instance where magnetic field dip angle will be
932 * stored.
933 */
934 public void getDipAsAngle(final Angle latitude, final Angle longitude, final Angle result) {
935 result.setUnit(AngleUnit.RADIANS);
936 result.setValue(getDip(latitude, longitude));
937 }
938
939 /**
940 * Returns the magnetic field dip angle from the Department of
941 * Defense geomagnetic model and data, in radians.
942 * This method uses default altitude (0.0 - mean sea level) and time
943 * (half way through the valid 5 year period of the model).
944 *
945 * @param latitude latitude angle.
946 * @param longitude longitude angle.
947 * @return the magnetic field dip angle.
948 */
949 public Angle getDipAsAngle(final Angle latitude, final Angle longitude) {
950 return new Angle(getDip(latitude, longitude), AngleUnit.RADIANS);
951 }
952
953 /**
954 * Returns the magnetic field dip angle from the Department of
955 * Defense geomagnetic model and data, in radians.
956 *
957 * @param latitude latitude expressed in radians.
958 * @param longitude longitude expressed in radians.
959 * @param height height expressed in meters.
960 * @param year year expressed in decimal years.
961 * @return the magnetic field dip angle expressed in radians.
962 */
963 public double getDip(final double latitude, final double longitude, final double height, final double year) {
964 final var heightKm = DistanceConverter.convert(height, DistanceUnit.METER, DistanceUnit.KILOMETER);
965 calcGeoMag(Math.toDegrees(latitude), Math.toDegrees(longitude), heightKm, year);
966 return Math.toRadians(dip);
967 }
968
969 /**
970 * Returns the magnetic field dip angle from the Department of
971 * Defense geomagnetic model and data, in radians.
972 *
973 * @param latitude latitude expressed in radians.
974 * @param longitude longitude expressed in radians.
975 * @param height height expressed in meters.
976 * @param calendar a calendar containing a specific instant.
977 * @return the magnetic field dip angle expressed in radians.
978 */
979 public double getDip(
980 final double latitude, final double longitude, final double height, final GregorianCalendar calendar) {
981 return getDip(latitude, longitude, height, convertTime(calendar));
982 }
983
984 /**
985 * Returns the magnetic field dip angle from the Department of
986 * Defense geomagnetic model and data, in radians.
987 *
988 * @param latitude latitude expressed in radians.
989 * @param longitude longitude expressed in radians.
990 * @param height height expressed in meters.
991 * @param time a specific time instant.
992 * @return the magnetic field dip angle expressed in radians.
993 */
994 public double getDip(final double latitude, final double longitude, final double height, final Date time) {
995 final var calendar = new GregorianCalendar();
996 calendar.setTime(time);
997 return getDip(latitude, longitude, height, calendar);
998 }
999
1000 /**
1001 * Returns the magnetic field dip angle from the Department of
1002 * Defense geomagnetic model and data, in radians.
1003 *
1004 * @param latitude latitude angle.
1005 * @param longitude longitude angle.
1006 * @param height height.
1007 * @param year year expressed in decimal years.
1008 * @return the magnetic field dip angle expressed in radians.
1009 */
1010 public double getDip(final Angle latitude, final Angle longitude, final Distance height, final double year) {
1011 return getDip(convertAngle(latitude), convertAngle(longitude), convertDistance(height), year);
1012 }
1013
1014 /**
1015 * Returns the magnetic field dip angle from the Department of
1016 * Defense geomagnetic model and data, in radians.
1017 *
1018 * @param latitude latitude angle.
1019 * @param longitude longitude angle.
1020 * @param height height.
1021 * @param calendar a calendar containing a specific instant.
1022 * @return the magnetic field dip angle expressed in radians.
1023 */
1024 public double getDip(final Angle latitude, final Angle longitude, final Distance height,
1025 final GregorianCalendar calendar) {
1026 return getDip(convertAngle(latitude), convertAngle(longitude), convertDistance(height), calendar);
1027 }
1028
1029 /**
1030 * Returns the magnetic field dip angle from the Department of
1031 * Defense geomagnetic model and data, in radians.
1032 *
1033 * @param latitude latitude angle.
1034 * @param longitude longitude angle.
1035 * @param height height.
1036 * @param time a specific time instant.
1037 * @return the magnetic field dip angle expressed in radians.
1038 */
1039 public double getDip(final Angle latitude, final Angle longitude, final Distance height, final Date time) {
1040 return getDip(convertAngle(latitude), convertAngle(longitude), convertDistance(height), time);
1041 }
1042
1043 /**
1044 * Returns the magnetic field dip angle from the Department of
1045 * Defense geomagnetic model and data, in radians.
1046 *
1047 * @param position a position expressed in geodetic coordinates.
1048 * @param year year expressed in decimal years.
1049 * @return the magnetic field dip angle expressed in radians.
1050 */
1051 public double getDip(final NEDPosition position, final double year) {
1052 return getDip(position.getLatitude(), position.getLongitude(), position.getHeight(), year);
1053 }
1054
1055 /**
1056 * Returns the magnetic field dip angle from the Department of
1057 * Defense geomagnetic model and data, in radians.
1058 *
1059 * @param position a position expressed in geodetic coordinates.
1060 * @param calendar a calendar containing a specific instant.
1061 * @return the magnetic field dip angle expressed in radians.
1062 */
1063 public double getDip(final NEDPosition position, final GregorianCalendar calendar) {
1064 return getDip(position.getLatitude(), position.getLongitude(), position.getHeight(), calendar);
1065 }
1066
1067 /**
1068 * Returns the magnetic field dip angle from the Department of
1069 * Defense geomagnetic model and data, in radians.
1070 *
1071 * @param position a position expressed in geodetic coordinates.
1072 * @param time a specific time instant.
1073 * @return the magnetic field dip angle expressed in radians.
1074 */
1075 public double getDip(final NEDPosition position, final Date time) {
1076 return getDip(position.getLatitude(), position.getLongitude(), position.getHeight(), time);
1077 }
1078
1079 /**
1080 * Returns the magnetic field dip angle from the Department of
1081 * Defense geomagnetic model and data, in radians.
1082 *
1083 * @param latitude latitude expressed in radians.
1084 * @param longitude longitude expressed in radians.
1085 * @param height height expressed in meters.
1086 * @param year year expressed in decimal years.
1087 * @param result instance where magnetic field dip angle will be
1088 * stored.
1089 */
1090 public void getDipAsAngle(
1091 final double latitude, final double longitude, final double height, final double year, final Angle result) {
1092 result.setUnit(AngleUnit.RADIANS);
1093 result.setValue(getDip(latitude, longitude, height, year));
1094 }
1095
1096 /**
1097 * Returns the magnetic field dip angle from the Department of
1098 * Defense geomagnetic model and data, in radians.
1099 *
1100 * @param latitude latitude expressed in radians.
1101 * @param longitude longitude expressed in radians.
1102 * @param height height expressed in meters.
1103 * @param year year expressed in decimal years.
1104 * @return the magnetic field dip angle.
1105 */
1106 public Angle getDipAsAngle(final double latitude, final double longitude, final double height, final double year) {
1107 return new Angle(getDip(latitude, longitude, height, year), AngleUnit.RADIANS);
1108 }
1109
1110 /**
1111 * Returns the magnetic field dip angle from the Department of
1112 * Defense geomagnetic model and data, in radians.
1113 *
1114 * @param latitude latitude expressed in radians.
1115 * @param longitude longitude expressed in radians.
1116 * @param height height expressed in meters.
1117 * @param calendar a calendar containing a specific instant.
1118 * @param result instance where magnetic field dip angle will be
1119 * stored.
1120 */
1121 public void getDipAsAngle(
1122 final double latitude, final double longitude, final double height, final GregorianCalendar calendar,
1123 final Angle result) {
1124 result.setUnit(AngleUnit.RADIANS);
1125 result.setValue(getDip(latitude, longitude, height, calendar));
1126 }
1127
1128 /**
1129 * Returns the magnetic field dip angle from the Department of
1130 * Defense geomagnetic model and data, in radians.
1131 *
1132 * @param latitude latitude expressed in radians.
1133 * @param longitude longitude expressed in radians.
1134 * @param height height expressed in meters.
1135 * @param calendar a calendar containing a specific instant.
1136 * @return the magnetic field dip angle.
1137 */
1138 public Angle getDipAsAngle(
1139 final double latitude, final double longitude, final double height, final GregorianCalendar calendar) {
1140 return new Angle(getDip(latitude, longitude, height, calendar), AngleUnit.RADIANS);
1141 }
1142
1143 /**
1144 * Returns the magnetic field dip angle from the Department of
1145 * Defense geomagnetic model and data, in radians.
1146 *
1147 * @param latitude latitude expressed in radians.
1148 * @param longitude longitude expressed in radians.
1149 * @param height height expressed in meters.
1150 * @param time a specific time instant.
1151 * @param result instance where magnetic field dip angle will be
1152 * stored.
1153 */
1154 public void getDipAsAngle(
1155 final double latitude, final double longitude, final double height, final Date time, final Angle result) {
1156 result.setUnit(AngleUnit.RADIANS);
1157 result.setValue(getDip(latitude, longitude, height, time));
1158 }
1159
1160 /**
1161 * Returns the magnetic field dip angle from the Department of
1162 * Defense geomagnetic model and data, in radians.
1163 *
1164 * @param latitude latitude expressed in radians.
1165 * @param longitude longitude expressed in radians.
1166 * @param height height expressed in meters.
1167 * @param time a specific time instant.
1168 * @return the magnetic field dip angle.
1169 */
1170 public Angle getDipAsAngle(
1171 final double latitude, final double longitude, final double height, final Date time) {
1172 return new Angle(getDip(latitude, longitude, height, time), AngleUnit.RADIANS);
1173 }
1174
1175 /**
1176 * Returns the magnetic field dip angle from the Department of
1177 * Defense geomagnetic model and data, in radians.
1178 *
1179 * @param latitude latitude angle.
1180 * @param longitude longitude angle.
1181 * @param height height.
1182 * @param year year expressed in decimal years.
1183 * @param result instance where magnetic field dip angle will be
1184 * stored.
1185 */
1186 public void getDipAsAngle(
1187 final Angle latitude, final Angle longitude, final Distance height, final double year, final Angle result) {
1188 result.setUnit(AngleUnit.RADIANS);
1189 result.setValue(getDip(latitude, longitude, height, year));
1190 }
1191
1192 /**
1193 * Returns the magnetic field dip angle from the Department of
1194 * Defense geomagnetic model and data, in radians.
1195 *
1196 * @param latitude latitude angle.
1197 * @param longitude longitude angle.
1198 * @param height height.
1199 * @param year year expressed in decimal years.
1200 * @return the magnetic field dip angle.
1201 */
1202 public Angle getDipAsAngle(final Angle latitude, final Angle longitude, final Distance height, final double year) {
1203 return new Angle(getDip(latitude, longitude, height, year), AngleUnit.RADIANS);
1204 }
1205
1206 /**
1207 * Returns the magnetic field dip angle from the Department of
1208 * Defense geomagnetic model and data, in radians.
1209 *
1210 * @param latitude latitude angle.
1211 * @param longitude longitude angle.
1212 * @param height height.
1213 * @param calendar a calendar containing a specific instant.
1214 * @param result instance where magnetic field dip angle will be
1215 * stored.
1216 */
1217 public void getDipAsAngle(final Angle latitude, final Angle longitude, final Distance height,
1218 final GregorianCalendar calendar, final Angle result) {
1219 result.setUnit(AngleUnit.RADIANS);
1220 result.setValue(getDip(latitude, longitude, height, calendar));
1221 }
1222
1223 /**
1224 * Returns the magnetic field dip angle from the Department of
1225 * Defense geomagnetic model and data, in radians.
1226 *
1227 * @param latitude latitude angle.
1228 * @param longitude longitude angle.
1229 * @param height height.
1230 * @param calendar a calendar containing a specific instant.
1231 * @return the magnetic field dip angle.
1232 */
1233 public Angle getDipAsAngle(
1234 final Angle latitude, final Angle longitude, final Distance height, final GregorianCalendar calendar) {
1235 return new Angle(getDip(latitude, longitude, height, calendar), AngleUnit.RADIANS);
1236 }
1237
1238 /**
1239 * Returns the magnetic field dip angle from the Department of
1240 * Defense geomagnetic model and data, in radians.
1241 *
1242 * @param latitude latitude angle.
1243 * @param longitude longitude angle.
1244 * @param height height.
1245 * @param time a specific time instant.
1246 * @param result instance where magnetic field dip angle will be
1247 * stored.
1248 */
1249 public void getDipAsAngle(
1250 final Angle latitude, final Angle longitude, final Distance height, final Date time, final Angle result) {
1251 result.setUnit(AngleUnit.RADIANS);
1252 result.setValue(getDip(latitude, longitude, height, time));
1253 }
1254
1255 /**
1256 * Returns the magnetic field dip angle from the Department of
1257 * Defense geomagnetic model and data, in radians.
1258 *
1259 * @param latitude latitude angle.
1260 * @param longitude longitude angle.
1261 * @param height height.
1262 * @param time a specific time instant.
1263 * @return the magnetic field dip angle.
1264 */
1265 public Angle getDipAsAngle(
1266 final Angle latitude, final Angle longitude, final Distance height, final Date time) {
1267 return new Angle(getDip(latitude, longitude, height, time), AngleUnit.RADIANS);
1268 }
1269
1270 /**
1271 * Returns the magnetic field dip angle from the Department of
1272 * Defense geomagnetic model and data, in radians.
1273 *
1274 * @param position a position expressed in geodetic coordinates.
1275 * @param year year expressed in decimal years.
1276 * @param result instance where magnetic field dip angle will be
1277 * stored.
1278 */
1279 public void getDipAsAngle(final NEDPosition position, final double year, final Angle result) {
1280 result.setUnit(AngleUnit.RADIANS);
1281 result.setValue(getDip(position, year));
1282 }
1283
1284 /**
1285 * Returns the magnetic field dip angle from the Department of
1286 * Defense geomagnetic model and data, in radians.
1287 *
1288 * @param position a position expressed in geodetic coordinates.
1289 * @param year year expressed in decimal years.
1290 * @return the magnetic field dip angle.
1291 */
1292 public Angle getDipAsAngle(
1293 final NEDPosition position, final double year) {
1294 return new Angle(getDip(position, year), AngleUnit.RADIANS);
1295 }
1296
1297 /**
1298 * Returns the magnetic field dip angle from the Department of
1299 * Defense geomagnetic model and data, in radians.
1300 *
1301 * @param position a position expressed in geodetic coordinates.
1302 * @param calendar a calendar containing a specific instant.
1303 * @param result instance where magnetic field dip angle will be
1304 * stored.
1305 */
1306 public void getDipAsAngle(final NEDPosition position, final GregorianCalendar calendar, final Angle result) {
1307 result.setUnit(AngleUnit.RADIANS);
1308 result.setValue(getDip(position, calendar));
1309 }
1310
1311 /**
1312 * Returns the magnetic field dip angle from the Department of
1313 * Defense geomagnetic model and data, in radians.
1314 *
1315 * @param position a position expressed in geodetic coordinates.
1316 * @param calendar a calendar containing a specific instant.
1317 * @return the magnetic field dip angle.
1318 */
1319 public Angle getDipAsAngle(final NEDPosition position, final GregorianCalendar calendar) {
1320 return new Angle(getDip(position, calendar), AngleUnit.RADIANS);
1321 }
1322
1323 /**
1324 * Returns the magnetic field dip angle from the Department of
1325 * Defense geomagnetic model and data, in radians.
1326 *
1327 * @param position a position expressed in geodetic coordinates.
1328 * @param time a specific time instant.
1329 * @param result instance where magnetic field dip angle will be
1330 * stored.
1331 */
1332 public void getDipAsAngle(final NEDPosition position, final Date time, final Angle result) {
1333 result.setUnit(AngleUnit.RADIANS);
1334 result.setValue(getDip(position, time));
1335 }
1336
1337 /**
1338 * Returns the magnetic field dip angle from the Department of
1339 * Defense geomagnetic model and data, in radians.
1340 *
1341 * @param position a position expressed in geodetic coordinates.
1342 * @param time a specific time instant.
1343 * @return the magnetic field dip angle.
1344 */
1345 public Angle getDipAsAngle(final NEDPosition position, final Date time) {
1346 return new Angle(getDip(position, time), AngleUnit.RADIANS);
1347 }
1348
1349 /**
1350 * Returns the magnetic field intensity from the Department of
1351 * Defense geomagnetic model and data expressed in Teslas.
1352 * This method uses default altitude (0.0 - mean sea level) and time
1353 * (half way through the valid 5 year period of the model).
1354 *
1355 * @param latitude latitude expressed in radians.
1356 * @param longitude longitude expressed in radians.
1357 * @return magnetic field strength expressed in Teslas (T).
1358 */
1359 public double getIntensity(final double latitude, final double longitude) {
1360 final var defaultTime = model.epoch + WorldMagneticModel.LIFESPAN / 2.0;
1361 calcGeoMag(Math.toDegrees(latitude), Math.toDegrees(longitude), 0.0, defaultTime);
1362 return ti * FROM_NANO;
1363 }
1364
1365 /**
1366 * Returns the magnetic field intensity from the Department of
1367 * Defense geomagnetic model and data expressed in Teslas.
1368 * This method uses default altitude (0.0 - mean sea level) and time
1369 * (half way through the valid 5 year period of the model).
1370 *
1371 * @param latitude latitude angle.
1372 * @param longitude longitude angle.
1373 * @return magnetic field strength expressed in Teslas (T).
1374 */
1375 public double getIntensity(final Angle latitude, final Angle longitude) {
1376 return getIntensity(convertAngle(latitude), convertAngle(longitude));
1377 }
1378
1379 /**
1380 * Returns the magnetic field intensity from the Department of
1381 * Defense geomagnetic model and data expressed in Teslas.
1382 *
1383 * @param latitude latitude expressed in radians.
1384 * @param longitude longitude expressed in radians.
1385 * @param height height expressed in meters.
1386 * @param year year expressed in decimal years.
1387 * @return magnetic field strength expressed in Teslas (T).
1388 */
1389 public double getIntensity(final double latitude, final double longitude, final double height, final double year) {
1390 final var heightKm = DistanceConverter.convert(height, DistanceUnit.METER, DistanceUnit.KILOMETER);
1391 calcGeoMag(Math.toDegrees(latitude), Math.toDegrees(longitude), heightKm, year);
1392 return ti * FROM_NANO;
1393 }
1394
1395 /**
1396 * Returns the magnetic field intensity from the Department of
1397 * Defense geomagnetic model and data expressed in Teslas.
1398 *
1399 * @param latitude latitude expressed in radians.
1400 * @param longitude longitude expressed in radians.
1401 * @param height height expressed in meters.
1402 * @param calendar a calendar containing a specific instant.
1403 * @return magnetic field strength expressed in Teslas (T).
1404 */
1405 public double getIntensity(
1406 final double latitude, final double longitude, final double height, final GregorianCalendar calendar) {
1407 return getIntensity(latitude, longitude, height, convertTime(calendar));
1408 }
1409
1410 /**
1411 * Returns the magnetic field intensity from the Department of
1412 * Defense geomagnetic model and data expressed in Teslas.
1413 *
1414 * @param latitude latitude expressed in radians.
1415 * @param longitude longitude expressed in radians.
1416 * @param height height expressed in meters.
1417 * @param time a specific time instant.
1418 * @return magnetic field strength expressed in Teslas (T).
1419 */
1420 public double getIntensity(
1421 final double latitude, final double longitude, final double height, final Date time) {
1422 final var calendar = new GregorianCalendar();
1423 calendar.setTime(time);
1424 return getIntensity(latitude, longitude, height, calendar);
1425 }
1426
1427 /**
1428 * Returns the magnetic field intensity from the Department of
1429 * Defense geomagnetic model and data expressed in Teslas.
1430 *
1431 * @param latitude latitude angle.
1432 * @param longitude longitude angle.
1433 * @param height height.
1434 * @param year year expressed in decimal years.
1435 * @return magnetic field strength expressed in Teslas (T).
1436 */
1437 public double getIntensity(
1438 final Angle latitude, final Angle longitude, final Distance height, final double year) {
1439 return getIntensity(convertAngle(latitude), convertAngle(longitude), convertDistance(height), year);
1440 }
1441
1442 /**
1443 * Returns the magnetic field intensity from the Department of
1444 * Defense geomagnetic model and data expressed in Teslas.
1445 *
1446 * @param latitude latitude angle.
1447 * @param longitude longitude angle.
1448 * @param height height.
1449 * @param calendar a calendar containing a specific instant.
1450 * @return magnetic field strength expressed in Teslas (T).
1451 */
1452 public double getIntensity(
1453 final Angle latitude, final Angle longitude, final Distance height, final GregorianCalendar calendar) {
1454 return getIntensity(convertAngle(latitude), convertAngle(longitude), convertDistance(height), calendar);
1455 }
1456
1457 /**
1458 * Returns the magnetic field intensity from the Department of
1459 * Defense geomagnetic model and data expressed in Teslas.
1460 *
1461 * @param latitude latitude angle.
1462 * @param longitude longitude angle.
1463 * @param height height.
1464 * @param time a specific time instant.
1465 * @return magnetic field strength expressed in Teslas (T).
1466 */
1467 public double getIntensity(final Angle latitude, final Angle longitude, final Distance height, final Date time) {
1468 return getIntensity(convertAngle(latitude), convertAngle(longitude), convertDistance(height), time);
1469 }
1470
1471 /**
1472 * Returns the magnetic field intensity from the Department of
1473 * Defense geomagnetic model and data expressed in Teslas.
1474 *
1475 * @param position a position expressed in geodetic coordinates.
1476 * @param year year expressed in decimal years.
1477 * @return magnetic field strength expressed in Teslas (T).
1478 */
1479 public double getIntensity(final NEDPosition position, final double year) {
1480 return getIntensity(position.getLatitude(), position.getLongitude(), position.getHeight(), year);
1481 }
1482
1483 /**
1484 * Returns the magnetic field intensity from the Department of
1485 * Defense geomagnetic model and data expressed in Teslas.
1486 *
1487 * @param position a position expressed in geodetic coordinates.
1488 * @param calendar a calendar containing a specific instant.
1489 * @return magnetic field strength expressed in Teslas (T).
1490 */
1491 public double getIntensity(final NEDPosition position, final GregorianCalendar calendar) {
1492 return getIntensity(position.getLatitude(), position.getLongitude(), position.getHeight(), calendar);
1493 }
1494
1495 /**
1496 * Returns the magnetic field intensity from the Department of
1497 * Defense geomagnetic model and data expressed in Teslas.
1498 *
1499 * @param position a position expressed in geodetic coordinates.
1500 * @param time a specific time instant.
1501 * @return magnetic field strength expressed in Teslas (T).
1502 */
1503 public double getIntensity(final NEDPosition position, final Date time) {
1504 return getIntensity(position.getLatitude(), position.getLongitude(), position.getHeight(), time);
1505 }
1506
1507 /**
1508 * Returns the horizontal magnetic field intensity from the
1509 * Department of Defense geomagnetic model and data expressed
1510 * in Teslas.
1511 * This method uses default altitude (0.0 - mean sea level) and time
1512 * (half way through the valid 5 year period of the model).
1513 *
1514 * @param latitude latitude expressed in radians.
1515 * @param longitude longitude expressed in radians.
1516 * @return the horizontal magnetic field strength expressed in
1517 * Teslas (T).
1518 */
1519 public double getHorizontalIntensity(final double latitude, final double longitude) {
1520 final var defaultTime = model.epoch + WorldMagneticModel.LIFESPAN / 2.0;
1521 calcGeoMag(Math.toDegrees(latitude), Math.toDegrees(longitude), 0.0, defaultTime);
1522 return bh * FROM_NANO;
1523 }
1524
1525 /**
1526 * Returns the horizontal magnetic field intensity from the
1527 * Department of Defense geomagnetic model and data expressed
1528 * in Teslas.
1529 * This method uses default altitude (0.0 - mean sea level) and time
1530 * (half way through the valid 5 year period of the model).
1531 *
1532 * @param latitude latitude angle.
1533 * @param longitude longitude angle.
1534 * @return the horizontal magnetic field strength expressed in
1535 * Teslas (T).
1536 */
1537 public double getHorizontalIntensity(final Angle latitude, final Angle longitude) {
1538 return getHorizontalIntensity(convertAngle(latitude), convertAngle(longitude));
1539 }
1540
1541 /**
1542 * Returns the horizontal magnetic field intensity from the
1543 * Department of Defense geomagnetic model and data expressed
1544 * in Teslas.
1545 *
1546 * @param latitude latitude expressed in radians.
1547 * @param longitude longitude expressed in radians.
1548 * @param height height expressed in meters.
1549 * @param year year expressed in decimal years.
1550 * @return the horizontal magnetic field strength expressed in
1551 * Teslas (T).
1552 */
1553 public double getHorizontalIntensity(
1554 final double latitude, final double longitude, final double height, final double year) {
1555 final var heightKm = DistanceConverter.convert(height, DistanceUnit.METER, DistanceUnit.KILOMETER);
1556 calcGeoMag(Math.toDegrees(latitude), Math.toDegrees(longitude), heightKm, year);
1557 return bh * FROM_NANO;
1558 }
1559
1560 /**
1561 * Returns the horizontal magnetic field intensity from the
1562 * Department of Defense geomagnetic model and data expressed
1563 * in Teslas.
1564 *
1565 * @param latitude latitude expressed in radians.
1566 * @param longitude longitude expressed in radians.
1567 * @param height height expressed in meters.
1568 * @param calendar a calendar containing a specific instant.
1569 * @return the horizontal magnetic field strength expressed in
1570 * Teslas (T).
1571 */
1572 public double getHorizontalIntensity(
1573 final double latitude, final double longitude,
1574 final double height, final GregorianCalendar calendar) {
1575 return getHorizontalIntensity(latitude, longitude, height, convertTime(calendar));
1576 }
1577
1578 /**
1579 * Returns the horizontal magnetic field intensity from the
1580 * Department of Defense geomagnetic model and data expressed
1581 * in Teslas.
1582 *
1583 * @param latitude latitude expressed in radians.
1584 * @param longitude longitude expressed in radians.
1585 * @param height height expressed in meters.
1586 * @param time a specific time instant.
1587 * @return the horizontal magnetic field strength expressed in
1588 * Teslas (T).
1589 */
1590 public double getHorizontalIntensity(
1591 final double latitude, final double longitude, final double height, final Date time) {
1592 final var calendar = new GregorianCalendar();
1593 calendar.setTime(time);
1594 return getHorizontalIntensity(latitude, longitude, height, calendar);
1595 }
1596
1597 /**
1598 * Returns the horizontal magnetic field intensity from the
1599 * Department of Defense geomagnetic model and data expressed
1600 * in Teslas.
1601 *
1602 * @param latitude latitude angle.
1603 * @param longitude longitude angle.
1604 * @param height height.
1605 * @param year year expressed in decimal years.
1606 * @return the horizontal magnetic field strength expressed in
1607 * Teslas (T).
1608 */
1609 public double getHorizontalIntensity(
1610 final Angle latitude, final Angle longitude, final Distance height, final double year) {
1611 return getHorizontalIntensity(convertAngle(latitude), convertAngle(longitude), convertDistance(height), year);
1612 }
1613
1614 /**
1615 * Returns the horizontal magnetic field intensity from the
1616 * Department of Defense geomagnetic model and data expressed
1617 * in Teslas.
1618 *
1619 * @param latitude latitude angle.
1620 * @param longitude longitude angle.
1621 * @param height height.
1622 * @param calendar a calendar containing a specific instant.
1623 * @return the horizontal magnetic field strength expressed in
1624 * Teslas (T).
1625 */
1626 public double getHorizontalIntensity(
1627 final Angle latitude, final Angle longitude, final Distance height, final GregorianCalendar calendar) {
1628 return getHorizontalIntensity(convertAngle(latitude), convertAngle(longitude), convertDistance(height),
1629 calendar);
1630 }
1631
1632 /**
1633 * Returns the horizontal magnetic field intensity from the
1634 * Department of Defense geomagnetic model and data expressed
1635 * in Teslas.
1636 *
1637 * @param latitude latitude angle.
1638 * @param longitude longitude angle.
1639 * @param height height.
1640 * @param time a specific time instant.
1641 * @return the horizontal magnetic field strength expressed in
1642 * Teslas (T).
1643 */
1644 public double getHorizontalIntensity(
1645 final Angle latitude, final Angle longitude, final Distance height, final Date time) {
1646 return getHorizontalIntensity(convertAngle(latitude), convertAngle(longitude), convertDistance(height),
1647 time);
1648 }
1649
1650 /**
1651 * Returns the horizontal magnetic field intensity from the
1652 * Department of Defense geomagnetic model and data expressed
1653 * in Teslas.
1654 *
1655 * @param position a position expressed in geodetic coordinates.
1656 * @param year year expressed in decimal years.
1657 * @return the horizontal magnetic field strength expressed in
1658 * Teslas (T).
1659 */
1660 public double getHorizontalIntensity(final NEDPosition position, final double year) {
1661 return getHorizontalIntensity(position.getLatitude(), position.getLongitude(), position.getHeight(), year);
1662 }
1663
1664 /**
1665 * Returns the horizontal magnetic field intensity from the
1666 * Department of Defense geomagnetic model and data expressed
1667 * in Teslas.
1668 *
1669 * @param position a position expressed in geodetic coordinates.
1670 * @param calendar a calendar containing a specific instant.
1671 * @return the horizontal magnetic field strength expressed in
1672 * Teslas (T).
1673 */
1674 public double getHorizontalIntensity(final NEDPosition position, final GregorianCalendar calendar) {
1675 return getHorizontalIntensity(position.getLatitude(), position.getLongitude(), position.getHeight(), calendar);
1676 }
1677
1678 /**
1679 * Returns the horizontal magnetic field intensity from the
1680 * Department of Defense geomagnetic model and data expressed
1681 * in Teslas.
1682 *
1683 * @param position a position expressed in geodetic coordinates.
1684 * @param time a specific time instant.
1685 * @return the horizontal magnetic field strength expressed in
1686 * Teslas (T).
1687 */
1688 public double getHorizontalIntensity(final NEDPosition position, final Date time) {
1689 return getHorizontalIntensity(position.getLatitude(), position.getLongitude(), position.getHeight(), time);
1690 }
1691
1692 /**
1693 * Returns the vertical magnetic field intensity from the
1694 * Department of Defense geomagnetic model and data expressed
1695 * in Teslas.
1696 * This method uses default altitude (0.0 - mean sea level) and time
1697 * (half way through the valid 5 year period of the model).
1698 *
1699 * @param latitude latitude expressed in radians.
1700 * @param longitude longitude expressed in radians.
1701 * @return the vertical magnetic field strength expressed in
1702 * Teslas (T).
1703 */
1704 public double getVerticalIntensity(final double latitude, final double longitude) {
1705 final var defaultTime = model.epoch + WorldMagneticModel.LIFESPAN / 2.0;
1706 calcGeoMag(Math.toDegrees(latitude), Math.toDegrees(longitude), 0.0, defaultTime);
1707 return bz * FROM_NANO;
1708 }
1709
1710 /**
1711 * Returns the vertical magnetic field intensity from the
1712 * Department of Defense geomagnetic model and data expressed
1713 * in Teslas.
1714 * This method uses default altitude (0.0 - mean sea level) and time
1715 * (half way through the valid 5 year period of the model).
1716 *
1717 * @param latitude latitude angle.
1718 * @param longitude longitude angle.
1719 * @return the vertical magnetic field strength expressed in
1720 * Teslas (T).
1721 */
1722 public double getVerticalIntensity(final Angle latitude, final Angle longitude) {
1723 return getVerticalIntensity(convertAngle(latitude), convertAngle(longitude));
1724 }
1725
1726 /**
1727 * Returns the vertical magnetic field intensity from the
1728 * Department of Defense geomagnetic model and data expressed
1729 * in Teslas.
1730 *
1731 * @param latitude latitude expressed in radians.
1732 * @param longitude longitude expressed in radians.
1733 * @param height height expressed in meters.
1734 * @param year year expressed in decimal years.
1735 * @return the vertical magnetic field strength expressed in
1736 * Teslas (T).
1737 */
1738 public double getVerticalIntensity(final double latitude, final double longitude, final double height,
1739 final double year) {
1740 final var heightKm = DistanceConverter.convert(height, DistanceUnit.METER, DistanceUnit.KILOMETER);
1741 calcGeoMag(Math.toDegrees(latitude), Math.toDegrees(longitude), heightKm, year);
1742 return bz * FROM_NANO;
1743 }
1744
1745 /**
1746 * Returns the vertical magnetic field intensity from the
1747 * Department of Defense geomagnetic model and data expressed
1748 * in Teslas.
1749 *
1750 * @param latitude latitude expressed in radians.
1751 * @param longitude longitude expressed in radians.
1752 * @param height height expressed in meters.
1753 * @param calendar a calendar containing a specific instant.
1754 * @return the vertical magnetic field strength expressed in
1755 * Teslas (T).
1756 */
1757 public double getVerticalIntensity(
1758 final double latitude, final double longitude, final double height, final GregorianCalendar calendar) {
1759 return getVerticalIntensity(latitude, longitude, height, convertTime(calendar));
1760 }
1761
1762 /**
1763 * Returns the vertical magnetic field intensity from the
1764 * Department of Defense geomagnetic model and data expressed
1765 * in Teslas.
1766 *
1767 * @param latitude latitude expressed in radians.
1768 * @param longitude longitude expressed in radians.
1769 * @param height height expressed in meters.
1770 * @param time a specific time instant.
1771 * @return the vertical magnetic field strength expressed in
1772 * Teslas (T).
1773 */
1774 public double getVerticalIntensity(
1775 final double latitude, final double longitude, final double height, final Date time) {
1776 final var calendar = new GregorianCalendar();
1777 calendar.setTime(time);
1778 return getVerticalIntensity(latitude, longitude, height, calendar);
1779 }
1780
1781 /**
1782 * Returns the vertical magnetic field intensity from the
1783 * Department of Defense geomagnetic model and data expressed
1784 * in Teslas.
1785 *
1786 * @param latitude latitude angle.
1787 * @param longitude longitude angle.
1788 * @param height height.
1789 * @param year year expressed in decimal years.
1790 * @return the vertical magnetic field strength expressed in
1791 * Teslas (T).
1792 */
1793 public double getVerticalIntensity(final Angle latitude, final Angle longitude,
1794 final Distance height, final double year) {
1795 return getVerticalIntensity(convertAngle(latitude), convertAngle(longitude), convertDistance(height), year);
1796 }
1797
1798 /**
1799 * Returns the vertical magnetic field intensity from the
1800 * Department of Defense geomagnetic model and data expressed
1801 * in Teslas.
1802 *
1803 * @param latitude latitude angle.
1804 * @param longitude longitude angle.
1805 * @param height height.
1806 * @param calendar a calendar containing a specific instant.
1807 * @return the vertical magnetic field strength expressed in
1808 * Teslas (T).
1809 */
1810 public double getVerticalIntensity(
1811 final Angle latitude, final Angle longitude, final Distance height, final GregorianCalendar calendar) {
1812 return getVerticalIntensity(convertAngle(latitude), convertAngle(longitude), convertDistance(height), calendar);
1813 }
1814
1815 /**
1816 * Returns the vertical magnetic field intensity from the
1817 * Department of Defense geomagnetic model and data expressed
1818 * in Teslas.
1819 *
1820 * @param latitude latitude angle.
1821 * @param longitude longitude angle.
1822 * @param height height.
1823 * @param time a specific time instant.
1824 * @return the vertical magnetic field strength expressed in
1825 * Teslas (T).
1826 */
1827 public double getVerticalIntensity(
1828 final Angle latitude, final Angle longitude, final Distance height, final Date time) {
1829 return getVerticalIntensity(convertAngle(latitude), convertAngle(longitude), convertDistance(height), time);
1830 }
1831
1832 /**
1833 * Returns the vertical magnetic field intensity from the
1834 * Department of Defense geomagnetic model and data expressed
1835 * in Teslas.
1836 *
1837 * @param position a position expressed in geodetic coordinates.
1838 * @param year year expressed in decimal years.
1839 * @return the vertical magnetic field strength expressed in
1840 * Teslas (T).
1841 */
1842 public double getVerticalIntensity(final NEDPosition position, final double year) {
1843 return getVerticalIntensity(position.getLatitude(), position.getLongitude(), position.getHeight(), year);
1844 }
1845
1846 /**
1847 * Returns the vertical magnetic field intensity from the
1848 * Department of Defense geomagnetic model and data expressed
1849 * in Teslas.
1850 *
1851 * @param position a position expressed in geodetic coordinates.
1852 * @param calendar a calendar containing a specific instant.
1853 * @return the vertical magnetic field strength expressed in
1854 * Teslas (T).
1855 */
1856 public double getVerticalIntensity(final NEDPosition position, final GregorianCalendar calendar) {
1857 return getVerticalIntensity(position.getLatitude(), position.getLongitude(), position.getHeight(), calendar);
1858 }
1859
1860 /**
1861 * Returns the vertical magnetic field intensity from the
1862 * Department of Defense geomagnetic model and data expressed
1863 * in Teslas.
1864 *
1865 * @param position a position expressed in geodetic coordinates.
1866 * @param time a specific time instant.
1867 * @return the vertical magnetic field strength expressed in
1868 * Teslas (T).
1869 */
1870 public double getVerticalIntensity(final NEDPosition position, final Date time) {
1871 return getVerticalIntensity(position.getLatitude(), position.getLongitude(), position.getHeight(), time);
1872 }
1873
1874 /**
1875 * Returns the northerly magnetic field intensity from the
1876 * Department of Defense geomagnetic model and data expressed
1877 * in Teslas.
1878 * This method uses default altitude (0.0 - mean sea level) and time
1879 * (half way through the valid 5 year period of the model).
1880 *
1881 * @param latitude latitude expressed in radians.
1882 * @param longitude longitude expressed in radians.
1883 * @return the northerly component of the magnetic field strength
1884 * expressed in Teslas (T).
1885 */
1886 public double getNorthIntensity(final double latitude, final double longitude) {
1887 final var defaultTime = model.epoch + WorldMagneticModel.LIFESPAN / 2.0;
1888 calcGeoMag(Math.toDegrees(latitude), Math.toDegrees(longitude), 0.0, defaultTime);
1889 return bx * FROM_NANO;
1890 }
1891
1892 /**
1893 * Returns the northerly magnetic field intensity from the
1894 * Department of Defense geomagnetic model and data expressed
1895 * in Teslas.
1896 * This method uses default altitude (0.0 - mean sea level) and time
1897 * (half way through the valid 5 year period of the model).
1898 *
1899 * @param latitude latitude angle.
1900 * @param longitude longitude angle.
1901 * @return the northerly component of the magnetic field strength
1902 * expressed in Teslas (T).
1903 */
1904 public double getNorthIntensity(final Angle latitude, final Angle longitude) {
1905 return getNorthIntensity(convertAngle(latitude), convertAngle(longitude));
1906 }
1907
1908 /**
1909 * Returns the northerly magnetic field intensity from the
1910 * Department of Defense geomagnetic model and data expressed
1911 * in Teslas.
1912 *
1913 * @param latitude latitude expressed in radians.
1914 * @param longitude longitude expressed in radians.
1915 * @param height height expressed in meters.
1916 * @param year year expressed in decimal years.
1917 * @return the northerly component of the magnetic field strength
1918 * expressed in Teslas (T).
1919 */
1920 public double getNorthIntensity(
1921 final double latitude, final double longitude, final double height, final double year) {
1922 final var heightKm = DistanceConverter.convert(height, DistanceUnit.METER, DistanceUnit.KILOMETER);
1923 calcGeoMag(Math.toDegrees(latitude), Math.toDegrees(longitude), heightKm, year);
1924 return bx * FROM_NANO;
1925 }
1926
1927 /**
1928 * Returns the northerly magnetic field intensity from the
1929 * Department of Defense geomagnetic model and data expressed
1930 * in Teslas.
1931 *
1932 * @param latitude latitude expressed in radians.
1933 * @param longitude longitude expressed in radians.
1934 * @param height height expressed in meters.
1935 * @param calendar a calendar containing a specific instant.
1936 * @return the northerly component of the magnetic field strength
1937 * expressed in Teslas (T).
1938 */
1939 public double getNorthIntensity(
1940 final double latitude, final double longitude, final double height, final GregorianCalendar calendar) {
1941 return getNorthIntensity(latitude, longitude, height, convertTime(calendar));
1942 }
1943
1944 /**
1945 * Returns the northerly magnetic field intensity from the
1946 * Department of Defense geomagnetic model and data expressed
1947 * in Teslas.
1948 *
1949 * @param latitude latitude expressed in radians.
1950 * @param longitude longitude expressed in radians.
1951 * @param height height expressed in meters.
1952 * @param time a specific time instant.
1953 * @return the northerly component of the magnetic field strength
1954 * expressed in Teslas (T).
1955 */
1956 public double getNorthIntensity(
1957 final double latitude, final double longitude, final double height, final Date time) {
1958 final var calendar = new GregorianCalendar();
1959 calendar.setTime(time);
1960 return getNorthIntensity(latitude, longitude, height, calendar);
1961 }
1962
1963 /**
1964 * Returns the northerly magnetic field intensity from the
1965 * Department of Defense geomagnetic model and data expressed
1966 * in Teslas.
1967 *
1968 * @param latitude latitude angle.
1969 * @param longitude longitude angle.
1970 * @param height height.
1971 * @param year year expressed in decimal years.
1972 * @return the northerly component of the magnetic field strength
1973 * expressed in Teslas (T).
1974 */
1975 public double getNorthIntensity(
1976 final Angle latitude, final Angle longitude, final Distance height, final double year) {
1977 return getNorthIntensity(convertAngle(latitude), convertAngle(longitude), convertDistance(height), year);
1978 }
1979
1980 /**
1981 * Returns the northerly magnetic field intensity from the
1982 * Department of Defense geomagnetic model and data expressed
1983 * in Teslas.
1984 *
1985 * @param latitude latitude angle.
1986 * @param longitude longitude angle.
1987 * @param height height.
1988 * @param calendar a calendar containing a specific instant.
1989 * @return the northerly component of the magnetic field strength
1990 * expressed in Teslas (T).
1991 */
1992 public double getNorthIntensity(
1993 final Angle latitude, final Angle longitude, final Distance height, final GregorianCalendar calendar) {
1994 return getNorthIntensity(convertAngle(latitude), convertAngle(longitude), convertDistance(height), calendar);
1995 }
1996
1997 /**
1998 * Returns the northerly magnetic field intensity from the
1999 * Department of Defense geomagnetic model and data expressed
2000 * in Teslas.
2001 *
2002 * @param latitude latitude angle.
2003 * @param longitude longitude angle.
2004 * @param height height.
2005 * @param time a specific time instant.
2006 * @return the northerly component of the magnetic field strength
2007 * expressed in Teslas (T).
2008 */
2009 public double getNorthIntensity(
2010 final Angle latitude, final Angle longitude, final Distance height, final Date time) {
2011 return getNorthIntensity(convertAngle(latitude), convertAngle(longitude), convertDistance(height), time);
2012 }
2013
2014 /**
2015 * Returns the northerly magnetic field intensity from the
2016 * Department of Defense geomagnetic model and data expressed
2017 * in Teslas.
2018 *
2019 * @param position a position expressed in geodetic coordinates.
2020 * @param year year expressed in decimal years.
2021 * @return the northerly component of the magnetic field strength
2022 * expressed in Teslas (T).
2023 */
2024 public double getNorthIntensity(final NEDPosition position, final double year) {
2025 return getNorthIntensity(position.getLatitude(), position.getLongitude(), position.getHeight(), year);
2026 }
2027
2028 /**
2029 * Returns the northerly magnetic field intensity from the
2030 * Department of Defense geomagnetic model and data expressed
2031 * in Teslas.
2032 *
2033 * @param position a position expressed in geodetic coordinates.
2034 * @param calendar a calendar containing a specific instant.
2035 * @return the northerly component of the magnetic field strength
2036 * expressed in Teslas (T).
2037 */
2038 public double getNorthIntensity(final NEDPosition position, final GregorianCalendar calendar) {
2039 return getNorthIntensity(position.getLatitude(), position.getLongitude(), position.getHeight(), calendar);
2040 }
2041
2042 /**
2043 * Returns the northerly magnetic field intensity from the
2044 * Department of Defense geomagnetic model and data expressed
2045 * in Teslas.
2046 *
2047 * @param position a position expressed in geodetic coordinates.
2048 * @param time a specific time instant.
2049 * @return the northerly component of the magnetic field strength
2050 * expressed in Teslas (T).
2051 */
2052 public double getNorthIntensity(final NEDPosition position, final Date time) {
2053 return getNorthIntensity(position.getLatitude(), position.getLongitude(), position.getHeight(), time);
2054 }
2055
2056 /**
2057 * Returns the easterly magnetic field intensity from the
2058 * Department of Defense geomagnetic model and data expressed
2059 * in nano Teslas.
2060 * This method uses default altitude (0.0 - mean sea level) and time
2061 * (half way through the valid 5 year period of the model).
2062 *
2063 * @param latitude latitude expressed in radians.
2064 * @param longitude longitude expressed in radians.
2065 * @return the easterly component of the magnetic field strength
2066 * expressed in Teslas (T).
2067 */
2068 public double getEastIntensity(final double latitude, final double longitude) {
2069 final var defaultTime = model.epoch + WorldMagneticModel.LIFESPAN / 2.0;
2070 calcGeoMag(Math.toDegrees(latitude), Math.toDegrees(longitude), 0.0, defaultTime);
2071 return by * FROM_NANO;
2072 }
2073
2074 /**
2075 * Returns the easterly magnetic field intensity from the
2076 * Department of Defense geomagnetic model and data expressed
2077 * in nano Teslas.
2078 * This method uses default altitude (0.0 - mean sea level) and time
2079 * (half way through the valid 5 year period of the model).
2080 *
2081 * @param latitude latitude angle.
2082 * @param longitude longitude angle.
2083 * @return the easterly component of the magnetic field strength
2084 * expressed in Teslas (T).
2085 */
2086 public double getEastIntensity(final Angle latitude, final Angle longitude) {
2087 return getEastIntensity(convertAngle(latitude), convertAngle(longitude));
2088 }
2089
2090 /**
2091 * Returns the easterly magnetic field intensity from the
2092 * Department of Defense geomagnetic model and data expressed
2093 * in nano Teslas.
2094 *
2095 * @param latitude latitude expressed in radians.
2096 * @param longitude longitude expressed in radians.
2097 * @param height height expressed in meters.
2098 * @param year year expressed in decimal years.
2099 * @return the easterly component of the magnetic field strength
2100 * expressed in Teslas (T).
2101 */
2102 public double getEastIntensity(
2103 final double latitude, final double longitude, final double height, final double year) {
2104 final var heightKm = DistanceConverter.convert(height, DistanceUnit.METER, DistanceUnit.KILOMETER);
2105 calcGeoMag(Math.toDegrees(latitude), Math.toDegrees(longitude), heightKm, year);
2106 return by * FROM_NANO;
2107 }
2108
2109 /**
2110 * Returns the easterly magnetic field intensity from the
2111 * Department of Defense geomagnetic model and data expressed
2112 * in nano Teslas.
2113 *
2114 * @param latitude latitude expressed in radians.
2115 * @param longitude longitude expressed in radians.
2116 * @param height height expressed in meters.
2117 * @param calendar a calendar containing a specific instant.
2118 * @return the easterly component of the magnetic field strength
2119 * expressed in Teslas (T).
2120 */
2121 public double getEastIntensity(
2122 final double latitude, final double longitude, final double height, final GregorianCalendar calendar) {
2123 return getEastIntensity(latitude, longitude, height, convertTime(calendar));
2124 }
2125
2126 /**
2127 * Returns the easterly magnetic field intensity from the
2128 * Department of Defense geomagnetic model and data expressed
2129 * in nano Teslas.
2130 *
2131 * @param latitude latitude expressed in radians.
2132 * @param longitude longitude expressed in radians.
2133 * @param height height expressed in meters.
2134 * @param time a specific time instant.
2135 * @return the easterly component of the magnetic field strength
2136 * expressed in Teslas (T).
2137 */
2138 public double getEastIntensity(
2139 final double latitude, final double longitude, final double height, final Date time) {
2140 final var calendar = new GregorianCalendar();
2141 calendar.setTime(time);
2142 return getEastIntensity(latitude, longitude, height, calendar);
2143 }
2144
2145 /**
2146 * Returns the easterly magnetic field intensity from the
2147 * Department of Defense geomagnetic model and data expressed
2148 * in nano Teslas.
2149 *
2150 * @param latitude latitude angle.
2151 * @param longitude longitude angle.
2152 * @param height height.
2153 * @param year year expressed in decimal years.
2154 * @return the easterly component of the magnetic field strength
2155 * expressed in Teslas (T).
2156 */
2157 public double getEastIntensity(
2158 final Angle latitude, final Angle longitude, final Distance height, final double year) {
2159 return getEastIntensity(convertAngle(latitude), convertAngle(longitude), convertDistance(height), year);
2160 }
2161
2162 /**
2163 * Returns the easterly magnetic field intensity from the
2164 * Department of Defense geomagnetic model and data expressed
2165 * in nano Teslas.
2166 *
2167 * @param latitude latitude angle.
2168 * @param longitude longitude angle.
2169 * @param height height.
2170 * @param calendar a calendar containing a specific instant.
2171 * @return the easterly component of the magnetic field strength
2172 * expressed in Teslas (T).
2173 */
2174 public double getEastIntensity(
2175 final Angle latitude, final Angle longitude, final Distance height, final GregorianCalendar calendar) {
2176 return getEastIntensity(convertAngle(latitude), convertAngle(longitude), convertDistance(height), calendar);
2177 }
2178
2179 /**
2180 * Returns the easterly magnetic field intensity from the
2181 * Department of Defense geomagnetic model and data expressed
2182 * in nano Teslas.
2183 *
2184 * @param latitude latitude angle.
2185 * @param longitude longitude angle.
2186 * @param height height.
2187 * @param time a specific time instant.
2188 * @return the easterly component of the magnetic field strength
2189 * expressed in Teslas (T).
2190 */
2191 public double getEastIntensity(
2192 final Angle latitude, final Angle longitude, final Distance height, final Date time) {
2193 return getEastIntensity(convertAngle(latitude), convertAngle(longitude), convertDistance(height), time);
2194 }
2195
2196 /**
2197 * Returns the easterly magnetic field intensity from the
2198 * Department of Defense geomagnetic model and data expressed
2199 * in nano Teslas.
2200 *
2201 * @param position a position expressed in geodetic coordinates.
2202 * @param year year expressed in decimal years.
2203 * @return the easterly component of the magnetic field strength
2204 * expressed in Teslas (T).
2205 */
2206 public double getEastIntensity(final NEDPosition position, final double year) {
2207 return getEastIntensity(position.getLatitude(), position.getLongitude(), position.getHeight(), year);
2208 }
2209
2210 /**
2211 * Returns the easterly magnetic field intensity from the
2212 * Department of Defense geomagnetic model and data expressed
2213 * in nano Teslas.
2214 *
2215 * @param position a position expressed in geodetic coordinates.
2216 * @param calendar a calendar containing a specific instant.
2217 * @return the easterly component of the magnetic field strength
2218 * expressed in Teslas (T).
2219 */
2220 public double getEastIntensity(final NEDPosition position, final GregorianCalendar calendar) {
2221 return getEastIntensity(position.getLatitude(), position.getLongitude(), position.getHeight(), calendar);
2222 }
2223
2224 /**
2225 * Returns the easterly magnetic field intensity from the
2226 * Department of Defense geomagnetic model and data expressed
2227 * in nano Teslas.
2228 *
2229 * @param position a position expressed in geodetic coordinates.
2230 * @param time a specific time instant.
2231 * @return the easterly component of the magnetic field strength
2232 * expressed in Teslas (T).
2233 */
2234 public double getEastIntensity(final NEDPosition position, final Date time) {
2235 return getEastIntensity(position.getLatitude(), position.getLongitude(), position.getHeight(), time);
2236 }
2237
2238 /**
2239 * Estimates Earth magnetic flux density.
2240 * This method uses default altitude (0.0 - mean sea level) and time
2241 * (half way through the valid 5 year period of the model).
2242 *
2243 * @param latitude latitude expressed in radians.
2244 * @param longitude longitude expressed in radians.
2245 * @param result instance where magnetic flux will be stored resolved
2246 * around NED frame.
2247 */
2248 public void estimate(final double latitude, final double longitude, final NEDMagneticFluxDensity result) {
2249 final var defaultTime = model.epoch + WorldMagneticModel.LIFESPAN / 2.0;
2250 calcGeoMag(Math.toDegrees(latitude), Math.toDegrees(longitude), 0.0, defaultTime);
2251 final var bn = bx * FROM_NANO;
2252 final var be = by * FROM_NANO;
2253 final var bd = bz * FROM_NANO;
2254
2255 result.setCoordinates(bn, be, bd);
2256 }
2257
2258 /**
2259 * Estimates Earth magnetic flux density.
2260 * This method uses default altitude (0.0 - mean sea level) and time
2261 * (half way through the valid 5 year period of the model).
2262 *
2263 * @param latitude latitude expressed in radians.
2264 * @param longitude longitude expressed in radians.
2265 * @return Earth magnetic flux density resolved around NED frame.
2266 */
2267 public NEDMagneticFluxDensity estimate(final double latitude, final double longitude) {
2268 final var result = new NEDMagneticFluxDensity();
2269 estimate(latitude, longitude, result);
2270 return result;
2271 }
2272
2273 /**
2274 * Estimates Earth magnetic flux density.
2275 * This method uses default altitude (0.0 - mean sea level) and time
2276 * (half way through the valid 5 year period of the model).
2277 *
2278 * @param latitude latitude angle.
2279 * @param longitude longitude angle.
2280 * @param result instance where magnetic flux will be stored resolved
2281 * around NED frame.
2282 */
2283 public void estimate(final Angle latitude, final Angle longitude, final NEDMagneticFluxDensity result) {
2284 estimate(convertAngle(latitude), convertAngle(longitude), result);
2285 }
2286
2287 /**
2288 * Estimates Earth magnetic flux density.
2289 * This method uses default altitude (0.0 - mean sea level) and time
2290 * (half way through the valid 5 year period of the model).
2291 *
2292 * @param latitude latitude angle.
2293 * @param longitude longitude angle.
2294 * @return Earth magnetic flux density resolved around NED frame.
2295 */
2296 public NEDMagneticFluxDensity estimate(final Angle latitude, final Angle longitude) {
2297 return estimate(convertAngle(latitude), convertAngle(longitude));
2298 }
2299
2300 /**
2301 * Estimates Earth magnetic flux density.
2302 *
2303 * @param latitude latitude expressed in radians.
2304 * @param longitude longitude expressed in radians.
2305 * @param height height expressed in meters.
2306 * @param year year expressed in decimal years.
2307 * @param result instance where magnetic flux will be stored resolved
2308 * around NED frame.
2309 */
2310 public void estimate(final double latitude, final double longitude, final double height, final double year,
2311 final NEDMagneticFluxDensity result) {
2312
2313 final var heightKm = DistanceConverter.convert(height, DistanceUnit.METER, DistanceUnit.KILOMETER);
2314 calcGeoMag(Math.toDegrees(latitude), Math.toDegrees(longitude), heightKm, year);
2315 final var bn = bx * FROM_NANO;
2316 final var be = by * FROM_NANO;
2317 final var bd = bz * FROM_NANO;
2318
2319 result.setCoordinates(bn, be, bd);
2320 }
2321
2322 /**
2323 * Estimates Earth magnetic flux density.
2324 *
2325 * @param latitude latitude expressed in radians.
2326 * @param longitude longitude expressed in radians.
2327 * @param height height expressed in meters.
2328 * @param year year expressed in decimal years.
2329 * @return Earth magnetic flux density resolved around NED frame.
2330 */
2331 public NEDMagneticFluxDensity estimate(
2332 final double latitude, final double longitude, final double height, final double year) {
2333 final var result = new NEDMagneticFluxDensity();
2334 estimate(latitude, longitude, height, year, result);
2335 return result;
2336 }
2337
2338 /**
2339 * Estimates Earth magnetic flux density.
2340 *
2341 * @param latitude latitude expressed in radians.
2342 * @param longitude longitude expressed in radians.
2343 * @param height height expressed in meters.
2344 * @param calendar a calendar containing a specific instant.
2345 * @param result instance where magnetic flux will be stored resolved
2346 * around NED frame.
2347 */
2348 public void estimate(final double latitude, final double longitude, final double height,
2349 final GregorianCalendar calendar, final NEDMagneticFluxDensity result) {
2350 estimate(latitude, longitude, height, convertTime(calendar), result);
2351 }
2352
2353 /**
2354 * Estimates Earth magnetic flux density.
2355 *
2356 * @param latitude latitude expressed in radians.
2357 * @param longitude longitude expressed in radians.
2358 * @param height height expressed in meters.
2359 * @param calendar a calendar containing a specific instant.
2360 * @return Earth magnetic flux density resolved around NED frame.
2361 */
2362 public NEDMagneticFluxDensity estimate(
2363 final double latitude, final double longitude, final double height, final GregorianCalendar calendar) {
2364 final var result = new NEDMagneticFluxDensity();
2365 estimate(latitude, longitude, height, calendar, result);
2366 return result;
2367 }
2368
2369 /**
2370 * Estimates Earth magnetic flux density.
2371 *
2372 * @param latitude latitude expressed in radians.
2373 * @param longitude longitude expressed in radians.
2374 * @param height height expressed in meters.
2375 * @param time a specific time instant.
2376 * @param result instance where magnetic flux will be stored resolved
2377 * around NED frame.
2378 */
2379 public void estimate(final double latitude, final double longitude, final double height, final Date time,
2380 final NEDMagneticFluxDensity result) {
2381 final var calendar = new GregorianCalendar();
2382 calendar.setTime(time);
2383 estimate(latitude, longitude, height, calendar, result);
2384 }
2385
2386 /**
2387 * Estimates Earth magnetic flux density.
2388 *
2389 * @param latitude latitude expressed in radians.
2390 * @param longitude longitude expressed in radians.
2391 * @param height height expressed in meters.
2392 * @param time a specific time instant.
2393 * @return Earth magnetic flux density resolved around NED frame.
2394 */
2395 public NEDMagneticFluxDensity estimate(
2396 final double latitude, final double longitude, final double height, final Date time) {
2397 final var result = new NEDMagneticFluxDensity();
2398 estimate(latitude, longitude, height, time, result);
2399 return result;
2400 }
2401
2402 /**
2403 * Estimates Earth magnetic flux density.
2404 *
2405 * @param latitude latitude angle.
2406 * @param longitude longitude angle.
2407 * @param height height.
2408 * @param year year expressed in decimal years.
2409 * @param result instance where magnetic flux will be stored resolved
2410 * around NED frame.
2411 */
2412 public void estimate(final Angle latitude, final Angle longitude, final Distance height, final double year,
2413 final NEDMagneticFluxDensity result) {
2414 estimate(convertAngle(latitude), convertAngle(longitude), convertDistance(height), year, result);
2415 }
2416
2417 /**
2418 * Estimates Earth magnetic flux density.
2419 *
2420 * @param latitude latitude angle.
2421 * @param longitude longitude angle.
2422 * @param height height.
2423 * @param year year expressed in decimal years.
2424 * @return Earth magnetic flux density resolved around NED frame.
2425 */
2426 public NEDMagneticFluxDensity estimate(final Angle latitude, final Angle longitude, final Distance height,
2427 final double year) {
2428 final var result = new NEDMagneticFluxDensity();
2429 estimate(latitude, longitude, height, year, result);
2430 return result;
2431 }
2432
2433 /**
2434 * Estimates Earth magnetic flux density.
2435 *
2436 * @param latitude latitude angle.
2437 * @param longitude longitude angle.
2438 * @param height height.
2439 * @param calendar a calendar containing a specific instant.
2440 * @param result instance where magnetic flux will be stored resolved
2441 * around NED frame.
2442 */
2443 public void estimate(final Angle latitude, final Angle longitude, final Distance height,
2444 final GregorianCalendar calendar, final NEDMagneticFluxDensity result) {
2445 estimate(convertAngle(latitude), convertAngle(longitude), convertDistance(height), calendar, result);
2446 }
2447
2448 /**
2449 * Estimates Earth magnetic flux density.
2450 *
2451 * @param latitude latitude angle.
2452 * @param longitude longitude angle.
2453 * @param height height.
2454 * @param calendar a calendar containing a specific instant.
2455 * @return Earth magnetic flux density resolved around NED frame.
2456 */
2457 public NEDMagneticFluxDensity estimate(final Angle latitude, final Angle longitude, final Distance height,
2458 final GregorianCalendar calendar) {
2459 final var result = new NEDMagneticFluxDensity();
2460 estimate(latitude, longitude, height, calendar, result);
2461 return result;
2462 }
2463
2464 /**
2465 * Estimates Earth magnetic flux density.
2466 *
2467 * @param latitude latitude angle.
2468 * @param longitude longitude angle.
2469 * @param height height.
2470 * @param time a specific time instant.
2471 * @param result instance where magnetic flux will be stored resolved
2472 * around NED frame.
2473 */
2474 public void estimate(final Angle latitude, final Angle longitude, final Distance height, final Date time,
2475 final NEDMagneticFluxDensity result) {
2476 estimate(convertAngle(latitude), convertAngle(longitude), convertDistance(height), time, result);
2477 }
2478
2479 /**
2480 * Estimates Earth magnetic flux density.
2481 *
2482 * @param latitude latitude angle.
2483 * @param longitude longitude angle.
2484 * @param height height.
2485 * @param time a specific time instant.
2486 * @return Earth magnetic flux density resolved around NED frame.
2487 */
2488 public NEDMagneticFluxDensity estimate(final Angle latitude, final Angle longitude, final Distance height,
2489 final Date time) {
2490 final var result = new NEDMagneticFluxDensity();
2491 estimate(latitude, longitude, height, time, result);
2492 return result;
2493 }
2494
2495 /**
2496 * Estimates Earth magnetic flux density.
2497 *
2498 * @param position a position expressed in geodetic coordinates.
2499 * @param year year expressed in decimal years.
2500 * @param result instance where magnetic flux will be stored resolved
2501 * around NED frame.
2502 */
2503 public void estimate(final NEDPosition position, final double year, final NEDMagneticFluxDensity result) {
2504 estimate(position.getLatitude(), position.getLongitude(), position.getHeight(), year, result);
2505 }
2506
2507 /**
2508 * Estimates Earth magnetic flux density.
2509 *
2510 * @param position a position expressed in geodetic coordinates.
2511 * @param year year expressed in decimal years.
2512 * @return Earth magnetic flux density resolved around NED frame.
2513 */
2514 public NEDMagneticFluxDensity estimate(final NEDPosition position, final double year) {
2515 final var result = new NEDMagneticFluxDensity();
2516 estimate(position, year, result);
2517 return result;
2518 }
2519
2520 /**
2521 * Estimates Earth magnetic flux density.
2522 *
2523 * @param position a position expressed in geodetic coordinates.
2524 * @param calendar a calendar containing a specific instant.
2525 * @param result instance where magnetic flux will be stored resolved
2526 * around NED frame.
2527 */
2528 public void estimate(final NEDPosition position, final GregorianCalendar calendar,
2529 final NEDMagneticFluxDensity result) {
2530 estimate(position.getLatitude(), position.getLongitude(), position.getHeight(), calendar, result);
2531 }
2532
2533 /**
2534 * Estimates Earth magnetic flux density.
2535 *
2536 * @param position a position expressed in geodetic coordinates.
2537 * @param calendar a calendar containing a specific instant.
2538 * @return Earth magnetic flux density resolved around NED frame.
2539 */
2540 public NEDMagneticFluxDensity estimate(final NEDPosition position, final GregorianCalendar calendar) {
2541 final var result = new NEDMagneticFluxDensity();
2542 estimate(position, calendar, result);
2543 return result;
2544 }
2545
2546 /**
2547 * Estimates Earth magnetic flux density.
2548 *
2549 * @param position a position expressed in geodetic coordinates.
2550 * @param time a specific time instant.
2551 * @param result instance where magnetic flux will be stored resolved
2552 * around NED frame.
2553 */
2554 public void estimate(final NEDPosition position, final Date time, final NEDMagneticFluxDensity result) {
2555 estimate(position.getLatitude(), position.getLongitude(), position.getHeight(), time, result);
2556 }
2557
2558 /**
2559 * Estimates Earth magnetic flux density.
2560 *
2561 * @param position a position expressed in geodetic coordinates.
2562 * @param time a specific time instant.
2563 * @return Earth magnetic flux density resolved around NED frame.
2564 */
2565 public NEDMagneticFluxDensity estimate(final NEDPosition position, final Date time) {
2566 final var result = new NEDMagneticFluxDensity();
2567 estimate(position, time, result);
2568 return result;
2569 }
2570
2571 /**
2572 * Converts a time instant contained in a gregorian calendar to a
2573 * decimal year.
2574 *
2575 * @param calendar calendar containing a specific instant to be
2576 * converted.
2577 * @return converted value expressed in decimal years.
2578 */
2579 public static double convertTime(final GregorianCalendar calendar) {
2580 final var year = calendar.get(Calendar.YEAR);
2581 final double daysInYear;
2582 if (calendar.isLeapYear(year)) {
2583 daysInYear = 366.0;
2584 } else {
2585 daysInYear = 365.0;
2586 }
2587
2588 return year + Math.max(0, (calendar.get(Calendar.DAY_OF_YEAR) - 1)) / daysInYear;
2589 }
2590
2591 /**
2592 * Computes the declination (dec), inclination (dip), total intensity
2593 * (TI) and grid variation (GV - polar regions only, referenced to grid
2594 * north of polar stereographic projection) of the Earth's magnetic field
2595 * in geodetic coordinates from the coefficients of the current official
2596 * department of defense (DoD) spherical harmonic World Magnetic Model
2597 * (WMM-2010). The WMM series of models is updated every 5 years on
2598 * January 1st of those years which are divisible by 5 (i.e. 1980, 1985,
2599 * 1990, etc.) by the Naval Oceanographic Office in cooperation with the
2600 * British Geological Survey (BGS). The model is based on geomagnetic
2601 * survey measurements from aircraft, satellite and geomagnetic
2602 * observatories.
2603 * <p>
2604 * Accuracy:
2605 * In ocean areas at the Earth's surface over the entire 5 year life of
2606 * a degree and order 12 spherical harmonic model such as WMM-95, the
2607 * estimated RMS errors for the various magnetic components are:
2608 * DEC - 0.5 degrees
2609 * DIP - 0.5 degrees
2610 * TI - 200.0 nano Teslas (nT)
2611 * GV - 0.5 Degrees
2612 * <p>
2613 * Other magnetic components that can be derived from these four by
2614 * simple trigonometric relations will have the following approximate
2615 * errors over ocean areas:
2616 * X - 140 nT (North)
2617 * Y - 140 nT (East)
2618 * Z - 200 nT (Vertical) Positive is down
2619 * H - 200 nT (Horizontal)
2620 * <p>
2621 * Over land the RMS errors are expected to be somewhat higher, although
2622 * the RMS errors for DEC, DIP and GV are still estimated to be less than
2623 * 0.5 degree, for the entire 5-year life of the model at the Earth's
2624 * surface. The other component errors over land are more difficult to
2625 * estimate and so are not given.
2626 * <p>
2627 * The accuracy at any given time of all four geomagnetic parameters
2628 * depends on the geomagnetic latitude. The errors are least at the
2629 * equator and greatest at the magnetic poles.
2630 * <p>
2631 * It is very important to note that a degree and order 12 model, such
2632 * as WMM-2010 describes only the long wavelength spatial magnetic
2633 * fluctuations due to Earth's core. Not included in the WMM series
2634 * models are intermediate and short wavelength spatial fluctuations of
2635 * the geomagnetic field which originate in the Earth's mantle and crust.
2636 * Consequently, isolated angular errors at various positions on the
2637 * surface (primarily over land, in continental margins and over oceanic
2638 * seamounts, ridges and trenches) of several degrees may be expected.
2639 * Also not included in the model are nonsecular temporal fluctuations
2640 * of the geomagnetic field of magneto-spheric and ionospheric origin.
2641 * During magnetic storms, temporal fluctuations can cause substantial
2642 * deviations of the geomagnetic field from model values. In arctic and
2643 * antarctic regions, as well as in equatorial regions, deviations from
2644 * model values are both frequent and persistent.
2645 * <p>
2646 * If the required declination accuracy is more stringent than the WMM
2647 * series of models provide, then the user is advised to request special
2648 * (regional or local) surveys be performed and models prepared by the
2649 * USGS, which operates the US geomagnetic observatories.
2650 *
2651 * @param latitude the latitude in decimal degrees.
2652 * @param longitude the longitude in decimal degrees.
2653 * @param height the height (altitude) in kilometers.
2654 * @param year the date as a decimal year.
2655 */
2656 private void calcGeoMag(final double latitude, final double longitude, final double height, final double year) {
2657
2658 final var dt = year - model.epoch;
2659 final var rlon = Math.toRadians(longitude);
2660 final var rlat = Math.toRadians(latitude);
2661 final var srlon = Math.sin(rlon);
2662 final var srlat = Math.sin(rlat);
2663 final var crlon = Math.cos(rlon);
2664 final var crlat = Math.cos(rlat);
2665 final var srlat2 = srlat * srlat;
2666 final var crlat2 = crlat * crlat;
2667 sp[1] = srlon;
2668 cp[1] = crlon;
2669
2670 // Convert from geodetic coords to spherical coords.
2671 if (oldHeight == null || height != oldHeight || oldLatitude == null || latitude != oldLatitude) {
2672 final var q = Math.sqrt(a2 - c2 * srlat2);
2673 final var q1 = height * q;
2674 final var q2 = ((q1 + a2) / (q1 + b2)) * ((q1 + a2) / (q1 + b2));
2675 ct = srlat / Math.sqrt(q2 * crlat2 + srlat2);
2676 st = Math.sqrt(1.0 - (ct * ct));
2677 final var r2 = ((height * height) + 2.0 * q1 + (a4 - c4 * srlat2) / (q * q));
2678 r = Math.sqrt(r2);
2679 final var mD = Math.sqrt(a2 * crlat2 + b2 * srlat2);
2680 ca = (height + mD) / r;
2681 sa = c2 * crlat * srlat / (r * mD);
2682 }
2683 if (oldLongitude == null || longitude != oldLongitude) {
2684 for (var m = 2; m <= WorldMagneticModel.MAX_ORDER; m++) {
2685 sp[m] = sp[1] * cp[m - 1] + cp[1] * sp[m - 1];
2686 cp[m] = cp[1] * cp[m - 1] - sp[1] * sp[m - 1];
2687 }
2688 }
2689
2690 final var aor = RE_KM / r;
2691 var ar = aor * aor;
2692 var br = 0.0;
2693 var bt = 0.0;
2694 var bp = 0.0;
2695 var bpp = 0.0;
2696
2697 for (var n = 1; n <= WorldMagneticModel.MAX_ORDER; n++) {
2698 ar = ar * aor;
2699 for (int m = 0, D3 = 1, D4 = (n + m + D3) / D3; D4 > 0; D4--, m += D3) {
2700
2701 // compute unnormalized associated Legendre polynomials
2702 // and derivatives via recursion relations
2703 if (oldHeight == null || height != oldHeight || oldLatitude == null || latitude != oldLatitude) {
2704 if (n == m) {
2705 model.snorm[n + m * N] = st * model.snorm[n - 1 + (m - 1) * N];
2706 dp[m][n] = st * dp[m - 1][n - 1] + ct * model.snorm[n - 1 + (m - 1) * N];
2707 }
2708 if (n == 1 && m == 0) {
2709 model.snorm[n] = ct * model.snorm[0];
2710 dp[m][n] = ct * dp[m][n - 1] - st * model.snorm[0];
2711 }
2712 if (n > 1 && n != m) {
2713 if (m > n - 2) {
2714 model.snorm[n - 2 + m * N] = 0.0;
2715 }
2716 if (m > n - 2) {
2717 dp[m][n - 2] = 0.0;
2718 }
2719 model.snorm[n + m * N] = ct * model.snorm[n - 1 + m * N]
2720 - model.k[m][n] * model.snorm[n - 2 + m * N];
2721 dp[m][n] = ct * dp[m][n - 1] - st * model.snorm[n - 1 + m * N]
2722 - model.k[m][n] * dp[m][n - 2];
2723 }
2724 }
2725
2726 // time-adjust the Gauss coefficients
2727
2728 if (oldTime == null || year != oldTime) {
2729 tc[m][n] = model.c[m][n] + dt * model.cd[m][n];
2730
2731 if (m != 0) {
2732 tc[n][m - 1] = model.c[n][m - 1] + dt * model.cd[n][m - 1];
2733 }
2734 }
2735
2736 // accumulate terms of the spherical harmonic expansions
2737 final double temp1;
2738 final double temp2;
2739 var par = ar * model.snorm[n + m * N];
2740 if (m == 0) {
2741 temp1 = tc[m][n] * cp[m];
2742 temp2 = tc[m][n] * sp[m];
2743 } else {
2744 temp1 = tc[m][n] * cp[m] + tc[n][m - 1] * sp[m];
2745 temp2 = tc[m][n] * sp[m] - tc[n][m - 1] * cp[m];
2746 }
2747
2748 bt = bt - ar * temp1 * dp[m][n];
2749 bp += (model.fm[m] * temp2 * par);
2750 br += (model.fn[n] * temp1 * par);
2751
2752 // Special case: North/south geographic poles
2753
2754 if (st == 0.0 && m == 1) {
2755 if (n == 1) {
2756 pp[n] = pp[n - 1];
2757 } else {
2758 pp[n] = ct * pp[n - 1] - model.k[m][n] * pp[n - 2];
2759 }
2760 final double parp = ar * pp[n];
2761 bpp += (model.fm[m] * temp2 * parp);
2762 }
2763
2764 }
2765 }
2766
2767 if (st == 0.0) {
2768 bp = bpp;
2769 } else {
2770 bp /= st;
2771 }
2772
2773 // Rotate magnetic vector components from spherical to
2774 // geodetic coordinates.
2775 // by is the east-west field component
2776 // bx is the north-south field component
2777 // bz is the vertical field component.
2778 bx = -bt * ca - br * sa;
2779 by = bp;
2780 bz = bt * sa - br * ca;
2781
2782 // Compute declination (DEC), INCLINATION (DIP) and
2783 // total intensity (TI)
2784
2785 bh = Math.sqrt((bx * bx) + (by * by));
2786 ti = Math.sqrt((bh * bh) + (bz * bz));
2787 // Calculate the declination.
2788 dec = Math.toDegrees(Math.atan2(by, bx));
2789 dip = Math.toDegrees(Math.atan2(bz, bh));
2790
2791 oldTime = year;
2792 oldHeight = height;
2793 oldLatitude = latitude;
2794 oldLongitude = longitude;
2795 }
2796
2797 /**
2798 * Converts provided angle instance to radians.
2799 *
2800 * @param angle angle to be converted.
2801 * @return converted value expressed in radians.
2802 */
2803 private static double convertAngle(final Angle angle) {
2804 return AngleConverter.convert(angle.getValue().doubleValue(), angle.getUnit(), AngleUnit.RADIANS);
2805 }
2806
2807 /**
2808 * Converts provided distance to meters.
2809 *
2810 * @param distance distance to be converted.
2811 * @return converted value expressed in meters.
2812 */
2813 private static double convertDistance(final Distance distance) {
2814 return DistanceConverter.convert(distance.getValue().doubleValue(), distance.getUnit(), DistanceUnit.METER);
2815 }
2816 }