Coverage for magnetismi/model/cof.py: 99.07%

1import importlib.resources 

2import math 

3import pathlib 

4 

5from magnetismi import ENCODING 

6 

7CS = ('n', 'm', 'gnm', 'hnm', 'dgnm', 'dhnm') 

8STORE_PATH = pathlib.Path('magnetismi', 'model') 

9TWIN_END_TOKEN = '9' * 48 

10 

11YEARS_COVERED = tuple(y for y in range(2020, 2025 + 1)) 

12MODEL_FROM_YEAR = {y: '2020' for y in range(2020, 2025 + 1)} 

13 

14 

15class Coefficients: 

16 """The WMM coefficients relevant to the year given.""" 

17 

18 model = {} # type: ignore 

19 

20 def load(self, model_year: int) -> None: 

21 """Load the model.""" 

22 model_resource = f'wmm-{model_year}.txt' 

23 with importlib.resources.path(__package__, model_resource) as model_path: 

24 with open(model_path, 'rt', encoding=ENCODING) as handle: 

25 recs = [line.strip().split() for line in handle if line.strip() and line.strip() != TWIN_END_TOKEN] 

26 

27 epoc, model, model_date_string = recs[0] 

28 self.model['epoch'] = float(epoc) 

29 self.model['model_code'] = model 

30 self.model['model_date_string'] = model_date_string 

31 

32 self.model['coeffs'] = [ 

33 {c: int(r[s]) if c in ('n', 'm') else float(r[s]) for s, c in enumerate(CS)} 

34 for r in recs[1:] 

35 if len(r) == len(CS) 

36 ] 

37 

38 def _zeroes(self, count: int) -> list[float]: 

39 """Short and sweet.""" 

40 return [0.0 for _ in range(count)] 

41 

42 def evaluate(self) -> None: 

43 """Ja, ja, ja!""" 

44 self.maxord = self.maxdeg = 12 

45 self.tc = [self._zeroes(13) for _ in range(14)] 

46 self.sp = self._zeroes(14) 

47 self.cp = self._zeroes(14) 

48 self.cp[0] = 1.0 

49 self.pp = self._zeroes(13) 

50 self.pp[0] = 1.0 

51 self.p = [self._zeroes(14) for _ in range(14)] 

52 self.p[0][0] = 1.0 

53 self.dp = [self._zeroes(13) for _ in range(14)] 

54 self.a = 6378.137 

55 self.b = 6356.7523142 

56 self.re = 6371.2 

57 self.a2 = self.a * self.a 

58 self.b2 = self.b * self.b 

59 self.c2 = self.a2 - self.b2 

60 self.a4 = self.a2 * self.a2 

61 self.b4 = self.b2 * self.b2 

62 self.c4 = self.a4 - self.b4 

63 

64 self.c = [self._zeroes(14) for _ in range(14)] 

65 self.cd = [self._zeroes(14) for _ in range(14)] 

66 

67 for wmmnm in self.model['coeffs']: 

68 m = wmmnm['m'] 

69 n = wmmnm['n'] 

70 gnm = wmmnm['gnm'] 

71 hnm = wmmnm['hnm'] 

72 dgnm = wmmnm['dgnm'] 

73 dhnm = wmmnm['dhnm'] 

74 if m <= n: 74 ↛ 67line 74 didn't jump to line 67, because the condition on line 74 was never false

75 self.c[m][n] = gnm 

76 self.cd[m][n] = dgnm 

77 if m != 0: 

78 self.c[n][m - 1] = hnm 

79 self.cd[n][m - 1] = dhnm 

80 

81 # Convert Schmidt normalized Gauss coefficients to unnormalized 

82 self.snorm = [self._zeroes(13) for _ in range(13)] 

83 self.snorm[0][0] = 1.0 

84 self.k = [self._zeroes(13) for _ in range(13)] 

85 self.k[1][1] = 0.0 

86 self.fn = [float(n) for n in range(14) if n != 1] 

87 self.fm = [float(n) for n in range(13)] 

88 for n in range(1, self.maxord + 1): 

89 self.snorm[0][n] = self.snorm[0][n - 1] * (2.0 * n - 1) / n 

90 j = 2.0 

91 # for (m=0,D1=1,D2=(n-m+D1)/D1;D2>0;D2--,m+=D1): 

92 m = 0 

93 D1 = 1 

94 D2 = (n - m + D1) / D1 

95 while D2 > 0: 

96 self.k[m][n] = (((n - 1) * (n - 1)) - (m * m)) / ((2.0 * n - 1) * (2.0 * n - 3.0)) 

97 if m > 0: 

98 flnmj = ((n - m + 1.0) * j) / (n + m) 

99 self.snorm[m][n] = self.snorm[m - 1][n] * math.sqrt(flnmj) 

100 j = 1.0 

101 self.c[n][m - 1] = self.snorm[m][n] * self.c[n][m - 1] 

102 self.cd[n][m - 1] = self.snorm[m][n] * self.cd[n][m - 1] 

103 self.c[m][n] = self.snorm[m][n] * self.c[m][n] 

104 self.cd[m][n] = self.snorm[m][n] * self.cd[m][n] 

105 D2 = D2 - 1 

106 m = m + D1 

107 

108 def __init__(self, year: int): 

109 """Initialize the mode from file in case it covers the requested year.""" 

110 if year not in YEARS_COVERED: 

111 raise ValueError(f'requested year ({year}) not within {YEARS_COVERED}') 

112 

113 self.model = { 

114 'year_requested': year, 

115 'model_id': MODEL_FROM_YEAR[year], 

116 } 

117 self.load(self.model['model_id']) 

118 self.evaluate()