Coverage for pygeodesy/auxilats/auxily.py: 93%

114 statements  

« prev     ^ index     » next       coverage.py v7.2.2, created at 2024-06-10 14:08 -0400

1 

2# -*- coding: utf-8 -*- 

3 

4u'''(INTERNAL) I{Auxiliary} latitudes' classes, constants and functions. 

5 

6Class L{AuxAngle} transcoded to Python from I{Karney}'s C++ class U{AuxAngle 

7<https://GeographicLib.SourceForge.io/C++/doc/classGeographicLib_1_1AuxAngle.html>} 

8in I{GeographicLib version 2.2+}. 

9 

10Copyright (C) U{Charles Karney<mailto:Karney@Alum.MIT.edu>} (2022-2023) and licensed 

11under the MIT/X11 License. For more information, see the U{GeographicLib 

12<https://GeographicLib.SourceForge.io>} documentation. 

13''' 

14# make sure int/int division yields float quotient, see .basics 

15from __future__ import division as _; del _ # PYCHOK semicolon 

16 

17from pygeodesy.constants import INF, NAN, isinf, isnan, _0_0, _0_5, \ 

18 _1_0, _copysign_1_0, _over, _1_over 

19from pygeodesy.errors import AuxError, NN 

20from pygeodesy.fmath import hypot1 as _sc, hypot2_ 

21# from pygeodesy.interns import NN # from .errors 

22from pygeodesy.karney import ADict, _ALL_DOCS, _MODS # PYCHOK used! 

23# from pygeodesy.lazily import _ALL_DOCS, _ALL_MODS as _MODS # from .karney 

24# from pygeodesy.named import ADict # from .karney 

25from pygeodesy.utily import atan1 

26 

27from math import asinh, copysign 

28 

29__all__ = () 

30__version__ = '24.03.20' 

31 

32 

33class Aux(object): 

34 '''Enum-style Aux names. 

35 ''' 

36 GEOGRAPHIC = PHI = GEODETIC = 0 

37 PARAMETRIC = BETA = REDUCED = 1 

38 GEOCENTRIC = THETA = 2 # all ... 

39 RECTIFYING = MU = 3 # use n^2 

40 CONFORMAL = CHI = 4 # use n 

41 AUTHALIC = XI = 5 # use n 

42 N = 6 

43 N2 = 36 

44 

45 def __index__(self, aux): 

46 # throws KeyError, not IndexError 

47 return _MODS.auxilats.auxLat._Aux2Greek[aux] 

48 

49 def __len__(self): 

50 return Aux.N 

51 

52 def _1d(self, auxout, auxin): 

53 '''Get the 1-d index into N^2 coeffs. 

54 ''' 

55 N = Aux.N 

56 if 0 <= auxout < N and 0 <= auxin < N: 

57 return N * auxout + auxin 

58 raise AuxError(auxout=auxout, auxin=auxin, N=N) # == AssertionError 

59 

60 def Greek(self, aux): 

61 '''Get an angle's name (C{str}). 

62 ''' 

63 return _Aux2Greek.get(aux, NN) 

64 

65 def len(self, ALorder): # PYCHOK no cover 

66 aL = ALorder # aka Lmax 

67 mu = Aux.MU * (Aux.MU + 1) 

68 ns = Aux.N2 - Aux.N 

69 return (mu * (aL * (aL + 3) - 2 * (aL // 2)) // 4 + 

70 (ns - mu) * (aL * (aL + 1)) // 2) 

71 

72 def power(self, auxout, auxin): 

73 '''Get the C{convert} exponent (C{int} or C{None}). 

74 ''' 

75 self._1d(auxout, auxin) # validate 

76 return (auxout - auxin) if max(auxin, auxout) < Aux.MU else None 

77 

78 def use_n2(self, aux): 

79 return aux not in (Aux.CHI, Aux.XI) 

80 

81Aux = Aux() # PYCHOK singleton 

82 

83_Aux2Greek = {Aux.AUTHALIC: 'Xi', 

84 Aux.CONFORMAL: 'Chi', 

85 Aux.GEOCENTRIC: 'Theta', 

86 Aux.GEODETIC: 'Phi', # == .GEOGRAPHIC 

87 Aux.PARAMETRIC: 'Beta', 

88 Aux.RECTIFYING: 'Mu'} 

89_Greek2Aux = dict(map(reversed, _Aux2Greek.items())) # PYCHOK exported 

90# _Greek2Aux.update((_g.upper(), _x) for _g, _x in _Greek2Aux.items()) 

91 

92 

93class _Coeffs(ADict): 

94 '''(INTERNAL) With C{items keys} string-ified. 

95 ''' 

96 def items(self): 

97 for n, v in ADict.items(self): 

98 yield str(n), v 

99 

100 

101def _Dasinh(x, y): 

102 d = y - x 

103 if isinf(d): # PYCHOK no cover 

104 r = _0_0 

105 elif isnan(d): # PYCHOK no cover 

106 r = NAN 

107 elif d: 

108 xy = x * y 

109 if xy > 0: 

110 hx, hy = _sc(x), _sc(y) 

111 if xy < 1: 

112 hx, hy = hy, hx 

113 else: 

114 x = _1_0 / x 

115 y = _1_0 / y 

116 r = _over(x + y, hx * x + hy * y) 

117 r = asinh(r * d) / d 

118 else: 

119 r = (asinh(y) - asinh(x)) / d 

120 else: 

121 r = _1_over(_sc(x)) 

122 return r 

123 

124 

125def _Datan(x, y): 

126 xy = x * y 

127 r = xy + _1_0 

128 if isnan(r): # PYCHOK no cover 

129 pass 

130 elif x == y: 

131 r = _1_over(r) 

132 elif x > 0 and isinf(xy): # PYCHOK no cover 

133 r = _0_0 

134 else: 

135 d = y - x 

136 if (r + xy) > 0: 

137 r = atan1(d, r) / d # atan(d / r) / d 

138 else: 

139 r = (atan1(y) - atan1(x)) / d 

140 return r 

141 

142 

143def _Dh(x, y): 

144 r = x + y 

145 if isnan(r): 

146 pass # N.B. NAN for inf-inf 

147 elif isinf(x): # PYCHOK no cover 

148 r = copysign(_0_5, x) 

149 elif isinf(y): # PYCHOK no cover 

150 r = copysign(_0_5, y) 

151 else: 

152 snx, sny = _sn(x), _sn(y) 

153 dy = sny * y 

154 dx = snx * x 

155 d = dy + dx 

156 if (d * _0_5): 

157 if (x * y) > 0: 

158 r *= hypot2_(snx / _sc(y), snx * sny, 

159 sny / _sc(x)) / (d + d) 

160 else: 

161 r = _over((dy - dx) * _0_5, y - x) 

162 else: # underflow and x == y == d == 0 

163 r *= _0_5 # PYCHOK no cover 

164 return r 

165 

166 

167def _Dlam(x, y): # Chi1.tan, Chi2.tan 

168 # I{Divided difference} of the isometric latitude 

169 # with respect to the conformal latitude 

170 if isnan(x) or isnan(y): # PYCHOK no cover 

171 r = NAN 

172 elif isinf(x) or isinf(y): # PYCHOK no cover 

173 r = INF 

174 elif x == y: 

175 r = _sc(x) 

176 else: 

177 r = _over(_Dasinh(x, y), _Datan(x, y)) 

178 return r 

179 

180 

181def _Dm(X, Y, s): # in .auxDLat, .auxDST 

182 # Return M{(X - Y) * s}, inplace X 

183 X -= Y 

184 X *= s 

185 return X # Fsum 

186 

187 

188def _Dp0Dpsi(x, y): # Chi1.tan, Chi2.tan 

189 # I{Divided difference} of the spherical rhumb area 

190 # term with respect to the isometric latitude 

191 r = x + y 

192 if isnan(r): # PYCHOK no cover 

193 pass # NAN for inf-inf 

194 elif isinf(x): # PYCHOK no cover 

195 r = _copysign_1_0(x) 

196 elif isinf(y): # PYCHOK no cover 

197 r = _copysign_1_0(y) 

198 elif x == y: 

199 r = _sn(x) 

200 else: 

201 r = _Dasinh(_h(x), _h(y)) 

202 r = _over(_Dh(x, y) * r, _Dasinh(x, y)) 

203 return r 

204 

205 

206def _h(tx): 

207 '''(INTERNAL) M{h(tan(x)) = tan(x) * sin(x) / 2} 

208 ''' 

209 if tx: 

210 tx *= _sn(tx) * _0_5 

211 return tx # preserve signed-0 

212 

213 

214def _sn(tx): 

215 '''(INTERNAL) M{tx / sqrt(tx**2 + 1)}, converting tan to sin. 

216 ''' 

217 if tx: 

218 tx = _copysign_1_0(tx) if isinf(tx) else ( 

219 NAN if isnan(tx) else (tx / _sc(tx))) 

220 return tx # preserve signed-0 

221 

222 

223class _Ufloats(dict): # in .auxilats.auxily 

224 '''(INTERNAL) "Uniquify" floats. 

225 ''' 

226 n = 0 # total number of floats 

227 

228 def __call__(self, *fs): 

229 '''Return a tuple of "uniquified" floats. 

230 ''' 

231 self.n += len(fs) 

232 _f = self.setdefault 

233 return tuple(_f(f, f) for f in map(float, fs)) # PYCHOK as attr 

234 

235 def _Coeffs(self, ALorder, coeffs): 

236 '''Return C{coeffs} (C{_Coeffs}, I{embellished}). 

237 ''' 

238# if True: 

239# n = 0 

240# for d in coeffs.values(): 

241# n += sum(d.values()) 

242# assert n == self.n 

243 Cx = _Coeffs(coeffs) 

244 Cx.set_(ALorder=ALorder, # used in .auxilats.__main__ 

245 n=self.n, # total number of floats 

246 u=len(self.keys())) # unique floats 

247 return Cx 

248 

249 

250__all__ += _ALL_DOCS(Aux.__class__) 

251 

252# **) MIT License 

253# 

254# Copyright (C) 2023-2024 -- mrJean1 at Gmail -- All Rights Reserved. 

255# 

256# Permission is hereby granted, free of charge, to any person obtaining a 

257# copy of this software and associated documentation files (the "Software"), 

258# to deal in the Software without restriction, including without limitation 

259# the rights to use, copy, modify, merge, publish, distribute, sublicense, 

260# and/or sell copies of the Software, and to permit persons to whom the 

261# Software is furnished to do so, subject to the following conditions: 

262# 

263# The above copyright notice and this permission notice shall be included 

264# in all copies or substantial portions of the Software. 

265# 

266# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS 

267# OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, 

268# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL 

269# THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR 

270# OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, 

271# ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR 

272# OTHER DEALINGS IN THE SOFTWARE.