1

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

3

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

5

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

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

8 in I{GeographicLib version 2.2+}.

9

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

11 under 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

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

16

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

18 _1_0, _copysign_1_0, _over, _1_over

19 from pygeodesy.errors import AuxError, NN

20 from pygeodesy.fmath import hypot1 as _sc, hypot2_

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

22 from 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

25 from pygeodesy.utily import atan1

26

27 from math import asinh, copysign

28

29__all__ = ()

30__version__ = '24.03.20'

31

32

33 class 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

81 Aux = 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

93 class _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

101 def _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

125 def _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

143 def _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

167 def _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

181 def _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

188 def _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

206 def _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

214 def _sn(tx):

215 '''(INTERNAL) M{sin(x) = tan(x) / sqrt(tan(x)**2 + 1)}.

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

223 class _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.