ferromagnetic/igrf/
mod.rs1use crate::{MagneticComponents, OrthogonalStrength};
2
3mod coeffs;
4mod math;
5
6pub struct IGRFresults {
7 pub result: MagneticComponents,
8 pub sv: MagneticComponents,
10}
11
12pub struct IGRF {
13 coeffs: coeffs::IGRFCoeffs,
14}
15impl Default for IGRF {
16 fn default() -> IGRF {
17 IGRF {
18 coeffs: coeffs::igrf_data(),
19 }
20 }
21}
22impl IGRF {
23 pub fn calc(&self, lat: f64, lon: f64, alt: f64, date: f64) -> IGRFresults {
24 let (start_coeffs, end_coeffs, nmax) = self.coeffs.coeffs(date);
26 let (a, b) = math::shval3(lat, lon, alt, nmax as usize, &start_coeffs, &end_coeffs);
27 let dif_a = math::Difh::from_orthognal_strength(&a);
28 let dif_b = math::Difh::from_orthognal_strength(&b);
29
30 let result = MagneticComponents {
31 declination: dif_a.declination.to_degrees(),
32 inclination: dif_a.inclination.to_degrees(),
33 horizontal_intensity: dif_a.horizontal_intensity,
34 orthogonal_strength: a,
35 total_intensity: dif_a.total_intensity,
36 };
37
38 let sv = MagneticComponents {
39 declination: {
40 let mut ddot = (dif_b.declination - dif_a.declination).to_degrees();
41 if ddot > 180.0 {
42 ddot -= 360.0
43 }
44 if ddot <= -180.0 {
45 ddot += 360.0
46 }
47 ddot * 60.0
48 },
49 inclination: (dif_b.inclination - dif_a.inclination).to_degrees() * 60.,
50 horizontal_intensity: dif_b.horizontal_intensity - dif_a.horizontal_intensity,
51 orthogonal_strength: OrthogonalStrength {
52 north: b.north - result.orthogonal_strength.north,
53 east: b.east - result.orthogonal_strength.east,
54 down: b.down - result.orthogonal_strength.down,
55 },
56 total_intensity: dif_b.total_intensity - dif_a.total_intensity,
57 };
58
59 IGRFresults { result, sv }
60 }
61}