Coverage for pygeodesy/auxilats/auxily.py: 93%
110 statements
« prev ^ index » next coverage.py v7.2.2, created at 2023-08-06 15:27 -0400
« prev ^ index » next coverage.py v7.2.2, created at 2023-08-06 15:27 -0400
2# -*- coding: utf-8 -*-
4u'''(INTERNAL) Private I{Auxiliary} base classes, constants and functions.
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+}.
10Copyright (C) U{Charles Karney<mailto:Charles@Karney.com>} (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
17from pygeodesy.constants import INF, NAN, _0_0, _0_5, _1_0, _2_0, \
18 _copysign_1_0, isinf, isnan, _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 _ALL_DOCS, _Dict, _MODS # PYCHOK used!
23# from pygeodesy.lazily import _ALL_DOCS, _ALL_MODS as _MODS # from .karney
24# from pygeodesy.named import _Dict # from .karney
26from math import asinh, atan, copysign
28__all__ = ()
29__version__ = '23.08.06'
32class Aux(object):
33 '''Enum-style Aux names.
34 '''
35 GEOGRAPHIC = PHI = GEODETIC = 0
36 PARAMETRIC = BETA = REDUCED = 1
37 GEOCENTRIC = THETA = 2 # all ...
38 RECTIFYING = MU = 3 # use n^2
39 CONFORMAL = CHI = 4 # use n
40 AUTHALIC = XI = 5 # use n
41 N = 6
42 N2 = 36
44 def __index__(self, aux):
45 # throws KeyError, not IndexError
46 return _MODS.auxilats.auxLat._Aux2Greek[aux]
48 def __len__(self):
49 return Aux.N
51 def _1d(self, auxout, auxin):
52 '''Get the 1-d index into N^2 coeffs.
53 '''
54 N = Aux.N
55 if 0 <= auxout < N and 0 <= auxin < N:
56 return N * auxout + auxin
57 raise AuxError(auxout=auxout, auxin=auxin, N=N) # == AssertionError
59 def Greek(self, aux):
60 '''Get an angle's name (C{str}).
61 '''
62 return _Aux2Greek.get(aux, NN)
64 def len(self, ALorder): # PYCHOK no cover
65 aL = ALorder # aka Lmax
66 mu = Aux.MU * (Aux.MU + 1)
67 ns = Aux.N2 - Aux.N
68 return (mu * (aL * (aL + 3) - 2 * (aL // 2)) // 4 +
69 (ns - mu) * (aL * (aL + 1)) // 2)
71 def power(self, auxout, auxin):
72 '''Get the C{convert} exponent (C{int} or C{None}).
73 '''
74 self._1d(auxout, auxin) # validate
75 return (auxout - auxin) if max(auxin, auxout) < Aux.MU else None
77 def use_n2(self, aux):
78 return aux not in (Aux.CHI, Aux.XI)
80Aux = Aux() # PYCHOK singleton
82_Aux2Greek = {Aux.AUTHALIC: 'Xi',
83 Aux.CONFORMAL: 'Chi',
84 Aux.GEOCENTRIC: 'Theta',
85 Aux.GEODETIC: 'Phi', # == .GEOGRAPHIC
86 Aux.PARAMETRIC: 'Beta',
87 Aux.RECTIFYING: 'Mu'}
88_Greek2Aux = dict(map(reversed, _Aux2Greek.items())) # PYCHOK exported
91class _Coeffs(_Dict):
92 '''(INTERNAL) With C{items keys} string-ified.
93 '''
94 def items(self):
95 for n, v in _Dict.items(self):
96 yield str(n), v
99class _Ufloats(dict): # in .auxilats.auxily
100 '''(INTERNAL) "Uniquify" floats.
101 '''
102 n = 0 # total number of floats
104 def __call__(self, *fs):
105 '''Return a tuple of "uniquified" floats.
106 '''
107 self.n += len(fs)
108 _f = self.setdefault
109 return tuple(_f(f, f) for f in map(float, fs)) # PYCHOK as attr
111 def _Coeffs(self, ALorder, coeffs):
112 '''Return C{coeffs} (C{_Coeffs}, I{embellished}).
113 '''
114# if True:
115# n = 0
116# for d in coeffs.values():
117# for t in d.values():
118# n += len(t)
119# assert n == self.n
120 Cx = _Coeffs(coeffs)
121 Cx.set_(ALorder=ALorder, # used in .auxilats.__main__
122 n=self.n, # total number of floats
123 u=len(self.keys())) # unique floats
124 return Cx
127def _2cos2x(cosx, sinx):
128 # M{2 * cos(2 * x)} from cos(x) and sin(x)
129 r = cosx - sinx
130 if r:
131 r *= (cosx + sinx) * _2_0
132 return r
135def _Dasinh(x, y):
136 d = y - x
137 if isinf(d): # PYCHOK no cover
138 r = _0_0
139 elif isnan(d): # PYCHOK no cover
140 r = NAN
141 elif d:
142 xy = x * y
143 if xy > 0:
144 hx, hy = _sc(x), _sc(y)
145 if xy < 1:
146 hx, hy = hy, hx
147 else:
148 x = _1_0 / x
149 y = _1_0 / y
150 r = _over(x + y, hx * x + hy * y)
151 r = asinh(r * d) / d
152 else:
153 r = (asinh(y) - asinh(x)) / d
154 else:
155 r = _1_over(_sc(x))
156 return r
159def _Datan(x, y):
160 xy = x * y
161 r = xy + _1_0
162 if isnan(r): # PYCHOK no cover
163 pass
164 elif x == y:
165 r = _1_over(r)
166 elif x > 0 and isinf(xy): # PYCHOK no cover
167 r = _0_0
168 else:
169 d = y - x
170 if (r + xy) > 0:
171 r = atan(_over(d, r)) / d # atan2(d, r) / d
172 else:
173 r = (atan(y) - atan(x)) / d
174 return r
177def _Dh(x, y):
178 r = x + y
179 if isnan(r):
180 pass # N.B. NAN for inf-inf
181 elif isinf(x): # PYCHOK no cover
182 r = copysign(_0_5, x)
183 elif isinf(y): # PYCHOK no cover
184 r = copysign(_0_5, y)
185 else:
186 snx, sny = _sn(x), _sn(y)
187 d = sny * y + snx * x
188 if (d * _0_5):
189 if (x * y) > 0:
190 r *= hypot2_(snx / _sc(y), snx * sny,
191 sny / _sc(x)) / (d + d)
192 else:
193 r = _over(_h(y) - _h(x), y - x)
194 else: # underflow and x == y == d == 0
195 r *= _0_5 # PYCHOK no cover
196 return r
199def _Dlam(x, y): # Chi1.tan, Chi2.tan
200 # I{Divided difference} of the isometric latitude
201 # with respect to the conformal latitude
202 if isnan(x) or isnan(y): # PYCHOK no cover
203 r = NAN
204 elif isinf(x) or isinf(y): # PYCHOK no cover
205 r = INF
206 elif x == y:
207 r = _sc(x)
208 else:
209 r = _over(_Dasinh(x, y), _Datan(x, y))
210 return r
213def _Dp0Dpsi(x, y): # Chi1.tan, Chi2.tan
214 # I{Divided difference} of the spherical rhumb area
215 # term with respect to the isometric latitude
216 r = x + y
217 if isnan(r): # PYCHOK no cover
218 pass # NAN for inf-inf
219 elif isinf(x): # PYCHOK no cover
220 r = _copysign_1_0(x)
221 elif isinf(y): # PYCHOK no cover
222 r = _copysign_1_0(y)
223 elif x == y:
224 r = _sn(x)
225 else:
226 r = _Dasinh(_h(x), _h(y))
227 r = _over(_Dh(x, y) * r, _Dasinh(x, y))
228 return r
231def _h(tx):
232 '''(INTERNAL) M{h(tan(x)) = tan(x) * sin(x) / 2}
233 '''
234 return (_sn(tx) * tx * _0_5) if tx else _0_0
237def _sn(tx):
238 '''(INTERNAL) M{tx / sqrt(tx**2 + 1)}, converting tan to sin.
239 '''
240 if tx:
241 tx = _copysign_1_0(tx) if isinf(tx) else (
242 NAN if isnan(tx) else (tx / _sc(tx)))
243 return tx # preserve signed-0
246__all__ += _ALL_DOCS(Aux.__class__)
248# **) MIT License
249#
250# Copyright (C) 2023-2023 -- mrJean1 at Gmail -- All Rights Reserved.
251#
252# Permission is hereby granted, free of charge, to any person obtaining a
253# copy of this software and associated documentation files (the "Software"),
254# to deal in the Software without restriction, including without limitation
255# the rights to use, copy, modify, merge, publish, distribute, sublicense,
256# and/or sell copies of the Software, and to permit persons to whom the
257# Software is furnished to do so, subject to the following conditions:
258#
259# The above copyright notice and this permission notice shall be included
260# in all copies or substantial portions of the Software.
261#
262# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
263# OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
264# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
265# THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR
266# OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE,
267# ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
268# OTHER DEALINGS IN THE SOFTWARE.