Coverage for pygeodesy / triaxials / bases.py: 91%

516 statements  

« prev     ^ index     » next       coverage.py v7.14.0, created at 2026-08-27 13:47 -0400

1 

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

3 

4u'''(INTERNAL) Base classes for I{ordered} triaxial ellipsoid classes L{Conformal}, L{Conformal3}, 

5L{Triaxial}, L{Triaxial3} and I{unordered} L{Triaxial_}. 

6 

7Transcoded to pure Python from I{Karney}'s GeographicLib 2.7 C++ classes U{Ellipsoid3<https:// 

8GeographicLib.SourceForge.io/C++/doc/classGeographicLib_1_1Triaxial_1_1Ellipsoid3.html>}, 

9U{Cartesian3<https://GeographicLib.SourceForge.io/C++/doc/classGeographicLib_1_1Triaxial_1_1Cartesian3.html>} and 

10U{Conformal3<https://GeographicLib.SourceForge.io/C++/doc/classGeographicLib_1_1Triaxial_1_1Conformal3.html>}. 

11 

12GeographicLib 2.5.2 C++ class U{JacobiConformal<https://GeographicLib.SourceForge.io/C++/doc/ 

13classGeographicLib_1_1JacobiConformal.html#details>}. 

14 

15Copyright (C) U{Charles Karney<mailto:Karney@Alum.MIT.edu>} (2008-2024, 2025) and licensed under the MIT/X11 License. 

16For more information, see the U{GeographicLib 2.5.2 and 2.7<https://GeographicLib.SourceForge.io/>} documentation. 

17 

18Enum-like C{Lat-/Longitude Kinds (LLK)}, see I{Karney}'s U{coord<https://GeographicLib.SourceForge.io/ 

19C++/doc/classGeographicLib_1_1Triaxial_1_1Cartesian3.html>}: 

20 

21@var LLK.CONFORMAL: Jacobi conformal X and Y projection 

22@var LLK.ELLIPSOIDAL: Ellipsoidal lat-, longitude and heading C{bet}, C{omg}, C{alp} (L{Ang}) 

23@var LLK.GEOCENTRIC: Geocentric lat-, longitude and heading C{phi}", C{lam}" and C{zet} (L{Ang}) 

24@var LLK.GEOCENTRIC_X: Geocentric with pole along major X axis 

25@var LLK.GEODETIC: Geodetic lat-, longitude and heading C{phi}, C{lam} and C{zet} (L{Ang}) 

26@var LLK.GEODETIC_X: Geodetic with pole along major X axis 

27@var LLK.GEODETIC_LON0: Geodetic lat-, longitude I{- lon0} and heading C{phi}, C{lam} and C{zet} (L{Ang}) 

28@var LLK.GEOGRAPHIC = LLK.GEODETIC 

29@var LLK.PARAMETRIC: Parametric lat-, longitude and heading C{phi}', C{lam}' and C{zet} (L{Ang}) 

30@var LLK.PARAMETRIC_X: Parametric with pole along major X axis 

31@var LLK.PLANETODETIC = LLK.GEODETIC 

32@var LLK.PLANETOCENTRIC = LLK.GEOCENTRIC 

33''' 

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

35from __future__ import division as _; del _ # noqa: E702 ; 

36 

37# from pygeodesy.angles import Ang # _MODS 

38# from pygeodesy.basics import map1 # from .namedTuples 

39from pygeodesy.constants import EPS, EPS0, EPS02, EPS4, INT0, NAN, PI_3, PI2, PI4, \ 

40 _EPS2e4, _isfinite, float0_, _1_over, _0_0, _1_0, \ 

41 _N_1_0, _3_0, _4_0 # PYCHOK used! 

42# from pygeodesy.ellipses import Ellipse, _isFlat # _MODS 

43# from pygeodesy.ellipsoids import Ellipsoid, _EWGS84 # _MODS 

44# from pygeodesy.elliptic import Elliptic # _MODS 

45# from pygeodesy.errors import _ValueError, _xkwds # from .utily 

46from pygeodesy.fmath import cbrt, fmean_, hypot, norm2, sqrt0, fabs, sqrt 

47from pygeodesy.fsums import _Fsumf_, fsumf_ 

48# from pygeodesy.internals import typename # _MODS 

49from pygeodesy.interns import _a_, _b_, _c_, _h_, _inside_, _not_, _NOTEQUAL_, _null_, \ 

50 _outside_, _scale_, _SPACE_, _spherical_, _x_, _y_, _z_ 

51from pygeodesy.lazily import _ALL_DOCS, _ALL_LAZY, _ALL_MODS as _MODS, _FOR_DOCS 

52from pygeodesy.named import _NamedEnum, _NamedEnumItem, _NamedTuple, _Pass 

53# from pygeodesy.named import _lazyNamedEnumItem as _lazy # _MODS 

54from pygeodesy.namedTuples import Ellipse5Tuple, Vector4Tuple, map1 

55from pygeodesy.props import Property_RO, property_doc_, property_RO, \ 

56 deprecated_method, deprecated_property_RO 

57# from pygeodesy.streprs import Fmt # _MODS 

58from pygeodesy.units import Degrees, Easting, Float, Height, Height_, _Lat0, \ 

59 Meter, Meter2, Meter3, Northing, Radius_, Scalar 

60from pygeodesy.utily import asin1, km2m, m2km, _ValueError, _xkwds 

61from pygeodesy.vector3d import _otherV3d, Vector3d 

62 

63# from math import fabs, sqrt # from .fmath 

64 

65__all__ = _ALL_LAZY.triaxials_bases 

66__version__ = '26.08.08' 

67 

68_bet_ = 'bet' # PYCHOK shared 

69_llk_ = 'llk' # PYCHOK shared 

70_KTpFlat = 1.5849625007 

71_MAXIT = 33 # 20 # PYCHOK shared 

72_not_ordered_ = _not_('ordered') 

73_omg_ = 'omg' # PYCHOK shared 

74 

75 

76class Conformal5Tuple(_NamedTuple): # see .Forward4Tuple 

77 '''5-Tuple C{(x, y, z, scale, llk)} with the easting C{x} and 

78 northing C{y} projection, C{scale} or C{NAN} I{but with} 

79 C{z=INT0} I{and kind} C{llk=LLK.CONFORMAL} I{always}. 

80 ''' 

81 _Names_ = (_x_, _y_, _z_, _scale_, _llk_) 

82 _Units_ = ( Easting, Northing, _Pass, Scalar, _Pass) 

83 

84 def __new__(cls, x, y, z=INT0, scale=NAN, llk=None, **kwds): # **iteration_name 

85 args = x, y, (z or INT0), scale, (llk or LLK.CONFORMAL) 

86 return _NamedTuple.__new__(cls, args, **kwds) 

87 

88 

89class _LLK(str): 

90 '''(INTERNAL) Lat-/Longitude Kind. 

91 ''' 

92 def __init__(self, llk): # aka C++ alt 

93 self._X = bool(llk.endswith('_X')) 

94 str.__init__(llk) 

95 

96 

97class LLK(object): 

98 '''Enum-like C{Lat-/Longitude Kinds (LLK)}, see U{coord<https://GeographicLib. 

99 SourceForge.io/C++/doc/classGeographicLib_1_1Triaxial_1_1Cartesian3.html>}. 

100 ''' 

101 CONFORMAL = _LLK('CONFORMAL') 

102 

103 ELLIPSOIDAL = _LLK('ELLIPSOIDAL') # bet, omg, alp 

104 GEOCENTRIC = _LLK('GEOCENTRIC') # phi2p, lam2p, zet 

105 GEOCENTRIC_X = _LLK('GEOCENTRIC_X') 

106 GEODETIC = _LLK('GEODETIC') # phi, lam, zet 

107 GEODETIC_LON0 = _LLK('GEODETIC_LON0') 

108 GEODETIC_X = _LLK('GEODETIC_X') 

109 GEOGRAPHIC = GEODETIC 

110 PARAMETRIC = _LLK('PARAMETRIC') # phi1p, lam1p, zet 

111 PARAMETRIC_X = _LLK('PARAMETRIC_X') 

112 PLANETODETIC = GEODETIC 

113 PLANETOCENTRIC = GEOCENTRIC 

114 

115 _CENTRICS = (GEOCENTRIC, GEOCENTRIC_X, PLANETOCENTRIC) 

116 _DETICS = (GEODETIC, GEODETIC_X, GEODETIC_LON0, GEOGRAPHIC, PLANETODETIC) 

117 _METRICS = (PARAMETRIC, PARAMETRIC_X) 

118 _NOIDAL = (None, ELLIPSOIDAL) 

119# _XCLUDE = (CONFORMAL, GEOGRAPHIC, PLANETOCENTRIC, PLANETODETIC) 

120 

121 def __getitem__(self, name): 

122 llk = self.get(name, None) 

123 if llk is None: 

124 t = _MODS.internals.typename(self) 

125 t = _MODS.streprs.Fmt.SQUARE(t, name) 

126 raise _ValueError(t, name) 

127 return llk 

128 

129 def get(self, name, dflt=None): 

130 '''Get an C{LLK} by C{name}. 

131 ''' 

132 llk = getattr(self, name, None) 

133 return llk if isinstance(llk, _LLK) else dflt 

134 

135 def items(self): 

136 '''Yield all C{LLK (name, value)} pairs. 

137 ''' 

138 for n, llk in LLK.__class__.__dict__.items(): 

139 if isinstance(llk, _LLK): 

140 yield n, llk 

141 

142 def keys(self): 

143 '''Yield all C{LLK} names. 

144 ''' 

145 for n, _ in self.items(): 

146 yield n 

147 

148 def values(self): 

149 '''Yield all C{LLK} values. 

150 ''' 

151 for _, llk in self.items(): 

152 yield llk 

153 

154if not _FOR_DOCS: # PYCHOK force epydoc 

155 LLK = LLK() # singleton 

156del _FOR_DOCS 

157 

158 

159def _HeightINT0(h, name=_h_, **kwds): # Error=... 

160 '''(INTERNAL) Return C{INT0} or C{Height(h=h, **kwds)}. 

161 ''' 

162 return h if h is INT0 else Height(h, name=name, **kwds) 

163 

164 

165class TriaxialError(_ValueError): 

166 '''Raised for any C{triaxial} issue. 

167 ''' 

168 pass # ... 

169 

170 

171class _UnOrderedTriaxialBase(_NamedEnumItem): 

172 '''(INTERNAL) Base class for all I{unordered} triaxial classes. 

173 ''' 

174 _ijk = _kji = None 

175 _unordered = True 

176 

177 def __init__(self, a_triaxial, b=None, c=None, **name): 

178 '''New I{unordered} C{Triaxial_}. 

179 

180 @arg a_triaxial: Large, C{X} semi-axis (C{scalar}, conventionally in 

181 C{meter}) or an other L{Triaxial}, L{Triaxial_} or 

182 L{TriaxialB} instance. 

183 @kwarg b: Middle, C{Y} semi-axis (C{meter}, same units as B{C{a}}), 

184 required if C{B{a_triaxial} is scalar}, ignored otherwise. 

185 @kwarg c: Small, C{Z} semi-axis (C{meter}, like B{C{b}}). 

186 @kwarg name: Optional C{B{name}=NN} (C{str}). 

187 

188 @raise TriaxialError: Invalid semi-axis or -axes. 

189 ''' 

190 try: 

191 try: 

192 a = a_triaxial 

193 t = a._abc3 

194 name = _xkwds(name, name=a.name) 

195 except AttributeError: 

196 t = Radius_(a=a), Radius_(b=b), Radius_(c=c) 

197 except (TypeError, ValueError) as x: 

198 raise TriaxialError(a=a, b=b, c=c, cause=x) 

199 if name: 

200 self.name = name 

201 

202 a, b, c = self._abc3 = t 

203 if self._unordered: # == not isinstance(self, Triaxial) 

204 s, _, t = sorted(t) 

205 if not (_isfinite(t) and _isfinite(s) and s > 0): 

206 raise TriaxialError(a=a, b=b, c=c) # txt=_invalid_ 

207 elif not (_isfinite(a) and a >= b >= c > 0): # see TriaxialB 

208 raise TriaxialError(a=a, b=b, c=c, txt=_not_ordered_) 

209 elif not (a > c and self._a2c2 > 0 and self.e2ac > 0): 

210 raise TriaxialError(a=a, c=c, e2ac=self.e2ac, txt=_spherical_) 

211 

212 def __repr__(self): 

213 '''Default C{repr(self)}. 

214 ''' 

215 return self.toRepr(terse=0) 

216 

217# def __str__(self): # in _NamedEnumItem 

218# return self.toStr() 

219 

220 @Property_RO 

221 def a(self): 

222 '''Get the C{largest, x} semi-axis (C{meter}, conventionally). 

223 ''' 

224 a, _, _ = self._abc3 

225 return a 

226 

227 @Property_RO 

228 def a2(self): 

229 '''Get C{a**2}. 

230 ''' 

231 return self.a**2 

232 

233 @Property_RO 

234 def _a2b2(self): 

235 '''(INTERNAL) Get C{a**2 - b**2} == E_sub_e**2. 

236 ''' 

237 a, b, _ = self._abc3 

238 d = a - b 

239 return (d * (a + b)) if d else _0_0 

240 

241 @Property_RO 

242 def _a2_b2(self): 

243 '''(INTERNAL) Get C{(a / b)**2}. 

244 ''' 

245 a, b, _ = self._abc3 

246 return (a / b)**2 if a != b else _1_0 

247 

248 @Property_RO 

249 def abc3(self): # in geed3solve._a12d 

250 '''Get the semi-axes as 3-tuple C{(a, b, c)}. 

251 ''' 

252 return self._abc3 

253 

254 @Property_RO 

255 def _a2b2c23(self): # in .triaxials.triaxial3 

256 '''(INTERNAL) Get 3-tuple C{(a**2, b**2, c**2)}. 

257 ''' 

258 return self.a2, self.b2, self.c2 

259 

260 @Property_RO 

261 def _a2c2(self): 

262 '''(INTERNAL) Get C{a**2 - c**2} == E_sub_x**2. 

263 ''' 

264 a, _, c = self._abc3 

265 d = a - c 

266 return (d * (a + c)) if d else _0_0 

267 

268 @Property_RO 

269 def area(self): 

270 '''Get the surface area (C{meter} I{squared}). 

271 ''' 

272 return self.areaKT(_KTpFlat) if self.isFlat else self.areaRG 

273 

274 def areaKT(self, *p): 

275 '''I{Approximate} the surface area using U{Knud Thomson's 

276 <https://WikiPedia.org/wiki/Ellipsoid#Approximate_formula>} 

277 formula (C{meter} I{squared}). 

278 

279 @arg p: Exponent (C{scalar} > 0), 1.6075 for near-spherical 

280 or 1.5849625007 for "near-flat" triaxials. 

281 ''' 

282 a, b, c = self._abc3 

283 _p = pow 

284 p = p[0] if p else (_KTpFlat if self.isFlat else 1.6075) 

285 a = _p(fmean_(_p(a * b, p), _p(a * c, p), _p(b * c, p)), _1_over(p)) 

286 return Meter2(areaKT=a * PI4) 

287 

288 @deprecated_method 

289 def area_p(self, p=1.6075): 

290 '''DEPRECATED on 2026-02-15, use method L{areaKT<Triaxial_.areaKT>}.''' 

291 return Meter2(area_p=self.areaKT(p)) 

292 

293 @Property_RO 

294 def areaRG(self): 

295 '''Get the surface area using Carlson's U{symmetric RG 

296 <https://WikiPedia.org/wiki/Ellipsoid#Surface_Area>} 

297 form (C{meter} I{squared}), see also C{Elliptic.fRG} 

298 ''' 

299 t = sorted(self._a2b2c23) # all non-zero 

300 r = _MODS.elliptic._rG3(*map(_1_over, t)) 

301 return Meter2(areaRG=self.volume * r * _3_0) 

302 

303 @Property_RO 

304 def b(self): 

305 '''Get the C{middle, y} semi-axis (C{meter}, same units as B{C{a}}). 

306 ''' 

307 _, b, _ = self._abc3 

308 return b 

309 

310 @Property_RO 

311 def b2(self): 

312 '''Get C{b**2}. 

313 ''' 

314 return self.b**2 

315 

316 @Property_RO 

317 def _b2_a2(self): 

318 '''(INTERNAL) Get C{(b / a)**2}. 

319 ''' 

320 a, b, _ = self._abc3 

321 return (b / a)**2 if a != b else _1_0 

322 

323 @Property_RO 

324 def _b2c2(self): 

325 '''(INTERNAL) Get C{b**2 - c**2} == E_sub_y**2. 

326 ''' 

327 _, b, c = self._abc3 

328 d = b - c 

329 return (d * (b + c)) if d else _0_0 

330 

331 @Property_RO 

332 def c(self): 

333 '''Get the C{smallest, z} semi-axis (C{meter}, same units as B{C{a}}). 

334 ''' 

335 _, _, c = self._abc3 

336 return c 

337 

338 @Property_RO 

339 def c2(self): 

340 '''Get C{c**2}. 

341 ''' 

342 return self.c**2 

343 

344 @Property_RO 

345 def _c2_a2(self): 

346 '''(INTERNAL) Get C{(c / a)**2}. 

347 ''' 

348 a, _, c = self._abc3 

349 return (c / a)**2 if a != c else _1_0 

350 

351 @Property_RO 

352 def _c2_b2(self): 

353 '''(INTERNAL) Get C{(c / b)**2}. 

354 ''' 

355 _, b, c = self._abc3 

356 return (c / b)**2 if b != c else _1_0 

357 

358 @Property_RO 

359 def e2ab(self): 

360 '''Get the C{ab} ellipse' I{(1st) eccentricity squared} (C{scalar}), M{1 - (b/a)**2}. 

361 ''' 

362 return Float(e2ab=(_1_0 - self._b2_a2) or _0_0) 

363 

364# _1e2ab = _b2_a2 # == C{1 - e2ab} == C{(b/a)**2} 

365 

366 @Property_RO 

367 def e2ac(self): 

368 '''Get the C{ac} ellipse' I{(1st) eccentricity squared} (C{scalar}), M{1 - (c/a)**2}. 

369 ''' 

370 return Float(e2ac=(_1_0 - self._c2_a2) or _0_0) 

371 

372# _1e2ac = _c2_a2 # == C{1 - e2ac} == C{(c/a)**2} 

373 

374 @Property_RO 

375 def e2bc(self): 

376 '''Get the C{bc} ellipse' I{(1st) eccentricity squared} (C{scalar}), M{1 - (c/b)**2}. 

377 ''' 

378 return Float(e2bc=(_1_0 - self._c2_b2) or _0_0) 

379 

380# _1e2bc = _c2_b2 # == C{1 - e2bc} == C{(c/b)**2} 

381 

382 def ellipse5(self, lat): 

383 '''Get the equatorial or a parallel I{ellipse of lattitude}. 

384 

385 @arg lat: Geodetic latitude (C{degrees90}, C{str} or C{Ang}). 

386 

387 @return: An L{Ellipse5Tuple}C{(a, b, height, lat, beta)} with C{a}, 

388 C{b} and C{height} measured along this triaxial's semi-axis 

389 C{a}, C{b} and C{c}, respectively. 

390 

391 @see: Method L{Ellipsoid.circle4<pygeodesy.Ellipsoid.circle4>} for 

392 further details. 

393 ''' 

394 a, b, c = self._abc3 

395 lat = _Lat0(lat) 

396 if lat and c > 0: 

397 E = _MODS.ellipsoids.Ellipsoid 

398 if a > b: 

399 r, z, lat, B = E(a, b=c).circle4(lat) 

400 b *= r / a 

401 a = r 

402 elif b > a: 

403 r, z, lat, B = E(b, b=c).circle4(lat) 

404 a *= r / b 

405 b = r 

406 else: # a == b 

407 r, z, lat, B = E(a, b=c).circle4(lat) 

408 a = b = r 

409 else: # equatorial or "flat" 

410 z = lat = B = _0_0 

411 return Ellipse5Tuple(a, b, z, lat, B) 

412 

413 def hartzell4(self, pov, los=False, **name): 

414 '''Compute the intersection of this triaxial's surface with a Line-Of-Sight 

415 from a Point-Of-View in space. 

416 

417 @see: Function L{hartzell4<triaxials.triaxial5.hartzell4>} for further details. 

418 ''' 

419 return _MODS.triaxials.hartzell4(pov, los=los, tri_biax=self, **name) 

420 

421 def height4(self, x_xyz, y=None, z=None, normal=True, eps=EPS, **name): 

422 '''Compute the projection on and the height above or below this triaxial's surface. 

423 

424 @see: Function L{height4<triaxials.triaxial5.height4>} for further details. 

425 ''' 

426 return _MODS.triaxials.height4(x_xyz, y=y, z=z, tri_biax=self, normal=normal, eps=eps, **name) 

427 

428 @Property_RO 

429 def isFlat(self): 

430 '''Is this triaxial "flat", too pro-/oblate (C{bool})? 

431 ''' 

432 _f = _MODS.ellipses._isFlat 

433 c, b, a = sorted(self._abc3) 

434 return _f(a, c) or _f(b, c) or _f(a, b) 

435 

436 @Property_RO 

437 def isOblate(self): 

438 '''Is this triaxial oblate (C{bool})? 

439 ''' 

440 return not (self.isProlate or self.isSpherical) 

441 

442 @Property_RO 

443 def isOrdered(self): 

444 '''Is this triaxial I{ordered} and I{not spherical} (C{bool})? 

445 ''' 

446 a, b, c = self._abc3 

447 return bool(a >= b > c) # b > c! 

448 

449 @Property_RO 

450 def isProlate(self): 

451 '''Is this triaxial prolate (C{bool})? 

452 ''' 

453 a, b, c = self._abc3 

454 return a < b or b < c or a < c 

455 

456 @Property_RO 

457 def isSpherical(self): 

458 '''Is this triaxial I{spherical} (C{Radius} or INT0)? 

459 ''' 

460 a, b, c = self._abc3 

461 return a if a == b == c else INT0 

462 

463 def _norm2(self, s, c, *a): 

464 '''(INTERNAL) Normalize C{s} and C{c} iff not already. 

465 ''' 

466 if fabs(_hypot2_1(s, c)) > EPS02: 

467 s, c = norm2(s, c) 

468 if a: 

469 s, c = norm2(s * self.b, c * a[0]) 

470 return float0_(s, c) 

471 

472 def normal3d(self, x_xyz, y=None, z=None, length=_1_0): 

473 '''Get a 3-D vector I{on and perpendicular to} this triaxial's surface. 

474 

475 @arg x_xyz: X component (C{scalar}) or a cartesian (C{Cartesian}, 

476 L{Ecef9Tuple}, L{Vector3d}, L{Vector3Tuple} or L{Vector4Tuple}). 

477 @kwarg y: Y component (C{scalar}), required if B{C{x_xyz}} if C{scalar}, ignored 

478 otherwise. 

479 @kwarg z: Z component (C{scalar}), like B{C{y}}. 

480 @kwarg length: Optional, signed length in out-/inward direction (C{scalar}). 

481 

482 @return: A C{Vector3d(x_, y_, z_)} normalized to B{C{length}}, pointing out- 

483 or inward for postive respectively negative B{C{length}}. 

484 

485 @raise TriaxialError: Zero length cartesian or vector. 

486 

487 @note: Cartesian C{(B{x}, B{y}, B{z})} I{must be on} this triaxial's surface, 

488 use method L{Triaxial.sideOf} to validate. 

489 

490 @see: Methods L{Triaxial.height4} and L{Triaxial.sideOf}. 

491 ''' 

492 # n = 2 * (x / a2, y / b2, z / c2) 

493 # == 2 * (x, y * a2 / b2, z * a2 / c2) / a2 # iff ordered 

494 # == 2 * (x, y / _b2_a2, z / _c2_a2) / a2 

495 # == unit(x, y / _b2_a2, z / _c2_a2).times(length) 

496 x, y, z = _otherV3d_(x_xyz, y, z).xyz3 

497 n = Vector3d(x, y / self._b2_a2, 

498 z / self._c2_a2, name__=self.normal3d) 

499 u = n.length 

500 if u < EPS0: 

501 raise TriaxialError(x=x_xyz, y=y, z=z, txt=_null_) 

502 return n.times(length / u) 

503 

504 def normal4(self, x_xyz, y=None, z=None, height=0, normal=True): 

505 '''Compute a cartesian at a B{C{height}} above or below this triaxial's surface. 

506 

507 @arg x_xyz: X component (C{scalar}) or a cartesian (C{Cartesian}, L{Ecef9Tuple}, 

508 L{Vector3d}, L{Vector3Tuple} or L{Vector4Tuple}). 

509 @kwarg y: Y component (C{scalar}), required if B{C{x_xyz}} if C{scalar}, ignored 

510 otherwise. 

511 @kwarg z: Z component (C{scalar}), like B{C{y}}. 

512 @kwarg normal: If C{True}, the B{C{height}} is I{perpendicular, plumb} to the 

513 triaxial's surface, otherwise C{radially} to the center of this 

514 triaxial (C{bool}). 

515 

516 @return: L{Vector4Tuple}C{(x, y, z, h)} with the cartesian coordinates C{x}, 

517 C{y} and C{z} and C{h} the I{signed, normal distance} to the triaxial's 

518 surface in C{meter}, conventionally. Positive C{h} indicates, the 

519 cartesian is outside the triaxial, negative C{h} means inside. 

520 

521 @raise TriaxialError: Zero length cartesian or vector. 

522 

523 @note: Cartesian C{(B{x}, B{y}, B{z})} I{must be on} this triaxial's surface, 

524 use method L{Triaxial.sideOf} to validate. 

525 

526 @see: Methods L{Triaxial.normal3d} and L{Triaxial.height4}. 

527 ''' 

528 v, h = _otherV3d_(x_xyz, y, z), Height_(height, low=None) 

529 if h: 

530 if v.length < EPS0: 

531 raise TriaxialError(x=x_xyz, y=y, z=z, txt=_null_) 

532 if normal: 

533 n = self.normal3d(v, length=h) 

534 h = n.length 

535 n += v 

536 else: 

537 h = h / v.length # /= chokes PyChecker 

538 n = v.times(h + _1_0) 

539 else: 

540 n = v 

541 return Vector4Tuple(n.x, n.y, n.z, h, name__=self.normal4) 

542 

543 def _order3(self, *abc, **reverse): # reverse=False 

544 '''(INTERNAL) Un-/Order C{a}, C{b} and C{c}. 

545 

546 @return: 3-Tuple C{(a, b, c)} ordered by or un-ordered 

547 (reverse-ordered) C{ijk} if C{B{reverse}=True}. 

548 ''' 

549 ijk = self._order_ijk(**reverse) 

550 return _getitems(abc, *ijk) if ijk else abc 

551 

552 def _order3d(self, v, **reverse): # reverse=False 

553 '''(INTERNAL) Un-/Order a C{Vector3d}. 

554 

555 @return: Vector3d(x, y, z) un-/ordered. 

556 ''' 

557 ijk = self._order_ijk(**reverse) 

558 return v.classof(*_getitems(v.xyz3, *ijk)) if ijk else v 

559 

560 @Property_RO 

561 def _ordered4(self): 

562 '''(INTERNAL) Helper for C{_hartzell3} and C{_plumbTo5}. 

563 ''' 

564 def _order2(reverse, a, b, c): 

565 '''(INTERNAL) Un-Order C{a}, C{b} and C{c}. 

566 

567 @return: 2-Tuple C{((a, b, c), ijk)} with C{a} >= C{b} >= C{c} 

568 and C{ijk} a 3-tuple with the initial indices. 

569 ''' 

570 i, j, k = range(3) 

571 if a < b: 

572 a, b, i, j = b, a, j, i 

573 if a < c: 

574 a, c, i, k = c, a, k, i 

575 if b < c: 

576 b, c, j, k = c, b, k, j 

577 # reverse (k, j, i) since (a, b, c) is reversed-sorted 

578 ijk = (k, j, i) if reverse else (None if i < j < k else (i, j, k)) 

579 return (a, b, c), ijk 

580 

581 abc, T = self._abc3, self 

582 if not self.isOrdered: 

583 abc, ijk = _order2(False, *abc) 

584 if ijk: 

585 _, kji = _order2(True, *ijk) 

586 T = _UnOrderedTriaxialBase(*abc) 

587 T._ijk, T._kji = ijk, kji 

588 return abc + (T,) 

589 

590 def _order_ijk(self, reverse=False): 

591 '''(INTERNAL) Get the un-/order indices. 

592 ''' 

593 return self._kji if reverse else self._ijk 

594 

595 @deprecated_property_RO 

596 def perimeter4ab(self): 

597 '''DEPRECATED on 2026.02.09, use property L{Ellipse<pygeodesy.Ellipse>}C{(a, b).perimeter2k_}.''' 

598 a, b, _ = self._abc3 

599 return Meter(perimeter4ab=_MODS.ellipses.Ellipse(a, b).perimeter2k_) 

600 

601 @deprecated_property_RO 

602 def perimeter4ac(self): 

603 '''DEPRECATED on 2026.02.09, use property L{Ellipse<pygeodesy.Ellipse>}C{(a, c).perimeter2k_}.''' 

604 a, _, c = self._abc3 

605 return Meter(perimeter4ac=_MODS.ellipses.Ellipse(a, c).perimeter2k_) 

606 

607 @deprecated_property_RO 

608 def perimeter4bc(self): 

609 '''DEPRECATED on 2026.02.09, use property L{Ellipse<pygeodesy.Ellipse>}C{(b, c).perimeter2k_}.''' 

610 _, b, c = self._abc3 

611 return Meter(perimeter4bc=_MODS.ellipses.Ellipse(b, c).perimeter2k_) 

612 

613 @Property_RO 

614 def R2(self): 

615 '''Get the I{authalic} earth radius (C{meter}), M{sqrt(area / PI4)}. 

616 ''' 

617 r = self.isSpherical 

618 return Meter(R2=r if r else sqrt(self.area / PI4)) # Radius 

619 

620 Rauthalic = R2 

621 

622 @Property_RO 

623 def R3(self): 

624 '''Get the I{volumetric} earth radius (C{meter}), M{(a * b * c)**(1/3)}. 

625 ''' 

626 a, b, c = self._abc3 

627 return Meter(R3=a if a == b == c else cbrt(a * b * c)) # Radius 

628 

629 Rvolumetric = R3 

630 

631 def _radialTo3(self, sbeta, cbeta, somega, comega): 

632 '''(INTERNAL) I{Unordered} helper for C{.height4}. 

633 ''' 

634 def _rphi(a, b, sphi, cphi): 

635 # <https://WikiPedia.org/wiki/Ellipse#Polar_form_relative_to_focus> 

636 # polar form: radius(phi) = a * b / hypot(a * sphi, b * cphi) 

637 return (b / hypot(sphi, b / a * cphi)) if a > b else ( 

638 (a / hypot(cphi, a / b * sphi)) if a < b else a) 

639 

640 sa, ca = self._norm2(sbeta, cbeta) 

641 sb, cb = self._norm2(somega, comega) 

642 

643 a, b, c = self._abc3 

644 if a != b: 

645 a = _rphi(a, b, sb, cb) 

646 if a != c: 

647 c = _rphi(a, c, sa, ca) 

648 t = c * ca 

649 return (t * cb), (t * sb), (c * sa) 

650 

651 def sideOf(self, x_xyz, y=None, z=None, eps=EPS4): 

652 '''Is a cartesian on, above or below the surface of this triaxial? 

653 

654 @arg x_xyz: X component (C{scalar}) or a cartesian (C{Cartesian}, 

655 L{Ecef9Tuple}, L{Vector3d}, L{Vector3Tuple} or L{Vector4Tuple}). 

656 @kwarg y: Y component (C{scalar}), required if B{C{x_xyz}} is C{scalar}, 

657 ignored otherwise. 

658 @kwarg z: Z component (C{scalar}), like B{C{y}}. 

659 @kwarg eps: On-surface tolerance (C{scalar}, distance I{squared}). 

660 

661 @return: Signed, radial distance I{squared} to this triangle's surface 

662 (C{scalar}), C{INT0} if within tolerance B{C{eps}}, positive 

663 if outside or negative if inside this triaxial. 

664 

665 @see: Methods L{Triaxial.height4} and L{Triaxial.normal3d}. 

666 ''' 

667 v = _otherV3d_(x_xyz, y, z) 

668 s2 = fsumf_(_N_1_0, *map(_over02, v.xyz3, self._abc3)) 

669 return INT0 if fabs(s2) < eps else Scalar(sideOf=s2) 

670 

671 def _sideOn(self, v, eps=_EPS2e4, Error=TriaxialError): # in pyaxqg 

672 s = self.sideOf(v.xyz, eps=eps) 

673 if s and Error: # PYCHOK no cover 

674 t = _SPACE_((_inside_ if s < 0 else _outside_), repr(self)) 

675 raise Error(eps=eps, sideOf=s, x=v.x, y=v.y, z=v.z, txt=t) 

676 return s 

677 

678 def toEllipsoid(self, **name): 

679 '''Convert this triaxial to a I{biaxial} L{Ellipsoid}, provided 2 axes match. 

680 

681 @kwarg name: Optional C{B{name}=NN} (C{str}). 

682 

683 @return: An L{Ellipsoid} with north along this C{Z} axis if C{a == b}, 

684 this C{Y} axis if C{a == c} or this C{X} axis if C{b == c}. 

685 

686 @raise TriaxialError: This C{a != b}, C{b != c} and C{c != a}. 

687 

688 @see: Method L{Ellipsoid.toTriaxial}. 

689 ''' 

690 a, b, c = self._abc3 

691 if a == b: 

692 b = c # N = c-Z 

693 elif b == c: # N = a-X 

694 a, b = b, a 

695 elif a != c: # N = b-Y 

696 t = _SPACE_(_a_, _NOTEQUAL_, _b_, _NOTEQUAL_, _c_) 

697 raise TriaxialError(a=a, b=b, c=c, txt=t) 

698 return _MODS.ellipsoids.Ellipsoid(a, b=b, name=self._name__(name)) 

699 

700 toBiaxial = toEllipsoid 

701 

702 def toStr(self, prec=9, terse=-3, **name): # PYCHOK signature 

703 '''Return this C{Triaxial} as a string. 

704 

705 @kwarg prec: Precision, number of decimal digits (0..9). 

706 @kwarg terse: Limit the number of items (C{int}, 3..11), 

707 use C{B{terse}=0} or C{=None} for all. 

708 @kwarg name: Optional name (C{str}), to override or C{None} 

709 to exclude this triaxial's name. 

710 

711 @return: This C{Triaxial}'s attributes (C{str}). 

712 ''' 

713 T = _UnOrderedTriaxialBase 

714 m = _MODS.triaxials 

715 C = m.Triaxial3B 

716 k = dict(**name) 

717 if isinstance(self, C): 

718 t = T.b, C.e2, C.k2, C.kp2 

719 else: 

720 t = T.a, # props 

721 C = m.ConformalSphere 

722 t += (C.ab, C.bc) if isinstance(self, C) else (T.b, T.c) 

723 C = _Triaxial3Base 

724 if isinstance(self, C): 

725 t += C.k2, C.kp2 

726 pm = self.Lon0 # PYCHOK attr 

727 if pm: 

728 k.update(Lon0=pm.degrees) 

729 else: 

730 t += T.e2ab, T.e2bc, T.e2ac 

731 for C in (m.Conformal, m.Conformal3): 

732 if isinstance(self, C): 

733 t += C.xyQ2, 

734 break 

735 t += T.volume, T.area, T.R2 

736 if terse: 

737 t = t[:terse] 

738 return self._instr(prec=prec, props=t, **k) 

739 

740 @Property_RO 

741 def unOrdered(self): 

742 '''Is this triaxial I{un-ordered} and I{not spherical} (C{bool})? 

743 ''' 

744 return not (self.isOrdered or bool(self.isSpherical)) 

745 

746 @Property_RO 

747 def volume(self): 

748 '''Get the volume (C{meter**3}), M{4 / 3 * PI * a * b * c}. 

749 ''' 

750 a, b, c = self._abc3 

751 return Meter3(volume=a * b * c * PI_3 * _4_0) 

752 

753 

754class _OrderedTriaxialBase(_UnOrderedTriaxialBase): 

755 '''(INTERNAL) Base class for all I{ordered} triaxial classes. 

756 ''' 

757 _unordered = False 

758 

759 def __init__(self, a_triaxial, b=None, c=None, **name): 

760 '''New I{ordered} L{Triaxial}, L{Triaxial3}, L{Conformal} or L{Conformal3}. 

761 

762 @arg a_triaxial: Largest semi-axis (C{scalar}, conventionally in C{meter}) 

763 or an other L{Triaxial} or L{Triaxial_} instance. 

764 @kwarg b: Middle semi-axis (C{meter}, same units as B{C{a}}), required 

765 if C{B{a_triaxial} is scalar}, ignored otherwise. 

766 @kwarg c: Smallest semi-axis (C{meter}, like B{C{b}}). 

767 @kwarg name: Optional C{B{name}=NN} (C{str}). 

768 

769 @note: The semi-axes must be ordered as C{B{a} >= B{b} >= B{c} > 0} and 

770 must be ellipsoidal, C{B{a} > B{c}}. 

771 

772 @raise TriaxialError: Semi-axes unordered, spherical or invalid. 

773 ''' 

774 _UnOrderedTriaxialBase.__init__(self, a_triaxial, b=b, c=c, **name) 

775 

776 @Property_RO 

777 def _a2b2_a2c2(self): 

778 '''@see: Methods C{.forwardBetaOmega} and property C{._k2_kp2E}. 

779 ''' 

780 s = self._a2c2 

781 if s: 

782 s = self._a2b2 / s 

783 return s or _0_0 

784 

785 @Property_RO 

786 def area(self): 

787 '''Get the surface area (C{meter} I{squared}). 

788 

789 @see: U{Surface area<https://WikiPedia.org/wiki/Ellipsoid#Surface_area>}. 

790 ''' 

791 a = self._areax 

792 if a is None: 

793 a = _UnOrderedTriaxialBase(self).area # or self.area21k 

794 return a 

795 

796 @Property_RO 

797 def area21k(self): 

798 '''Get the surface area using incomplete elliptic integrals of the 

799 2nd and 1st kind (C{meter} I{squared}), see also C{Elliptic.fE} 

800 respectively C{Elliptic.fF}. 

801 ''' 

802 a = self._areax 

803 if a is None: 

804 k2, kp2 = t = self._k2_kp2E 

805 if self.e2ac < EPS or min(t) < EPS or max(t) > _1_0: 

806 a = self.areaKT() # "flat" or near-spherical 

807 else: 

808 aE = _MODS.elliptic.Elliptic(k2=kp2, kp2=k2) # swapped! 

809 s = sqrt(self.e2ac) # == sin(phi) 

810 t = self._c2_a2 / s # == cos(phi)**2 / sin(phi) 

811 r = asin1(s) # phi 

812 a, b, c = self._abc3 

813 t = (aE.fE(r) * s + aE.fF(r) * t) * a * b 

814 a = Meter2(area21k=(c**2 + t) * PI2) 

815 return a 

816 

817 @Property_RO 

818 def _areax(self): 

819 '''(INTERNAL) Get the area as ellipsoidal or C{None}. 

820 ''' 

821 a, b, c = self._abc3 

822 return None if a != b else \ 

823 _MODS.ellipsoids.Ellipsoid(a, b=c).areax 

824 

825 @Property_RO 

826 def _k2_kp2E(self): 

827 '''(INTERNAL) Get elliptic C{k2} and C{kp2} for C{._xE}, C{._yE} and C{.areaE}. 

828 ''' 

829 # k2 = a2b2 / a2c2 * c2_b2 

830 # kp2 = b2c2 / a2c2 * a2_b2 

831 # b2 = b**2 

832 # xE = Elliptic(k2, -a2b2 / b2, kp2, a2_b2) 

833 # yE = Elliptic(kp2, b2c2 / b2, k2, c2_b2) 

834 # aE = Elliptic(kp2, 0, k2, 1) 

835 k2 = (self._c2_b2 * self._a2b2_a2c2) or _0_0 

836 kp2 = (self._a2_b2 * self._b2c2 / self._a2c2) if k2 else _1_0 

837 return k2, kp2 

838 

839 def _radialTo3(self, sbeta, cbeta, somega, comega): 

840 '''(INTERNAL) Convert I{ellipsoidal} lat- C{beta} and longitude 

841 C{omega} to a cartesian I{on this triaxial's surface}, also 

842 I{ordered} helper for C{.height4 with normal=False}. 

843 ''' 

844 sa, ca = self._norm2(sbeta, cbeta) 

845 sb, cb = self._norm2(somega, comega) 

846 

847 b2_a2 = self._b2_a2 # == (b/a)**2 

848 c2_a2 = -self._c2_a2 # == -(c/a)**2 

849 a2c2_a2 = self. e2ac # (a**2 - c**2) / a**2 == 1 - (c/a)**2 

850 

851 x2 = _Fsumf_(_1_0, -b2_a2 * sa**2, c2_a2 * ca**2).fover(a2c2_a2) 

852 z2 = _Fsumf_(c2_a2, sb**2, b2_a2 * cb**2).fover(a2c2_a2) 

853 

854 x, y, z = self._abc3 

855 x *= cb * _sqrt0(x2) 

856 y *= ca * sb 

857 z *= sa * _sqrt0(z2) 

858 return x, y, z 

859 

860 

861class _Triaxial3Base(_OrderedTriaxialBase): 

862 '''(INTERNAL) Base class for I{unordered} triaxialC{3} classes. 

863 ''' 

864 _e2_k2_kp2 = None 

865 _Lon0 = None 

866 _lon0WGS84_3 = -(1493 / 100) # in pyaxqg 

867 

868 @Property_RO 

869 def e2(self): 

870 '''Get the I{squared eccentricity} (C{scalar}), M{(a**2 - c**2) / b**2}. 

871 ''' 

872 if self._e2_k2_kp2: 

873 e2, _, _ = self._e2_k2_kp2 

874 else: 

875 e2 = self._a2c2 / self.b2 

876 return Float(e2=e2) 

877 

878 def _init_abc3_e2_k2_kp2(self, b, e2, k2, kp2, **name): 

879 '''(INTERNAL) C{Triaxial3B.__init__}. 

880 ''' 

881 if name: 

882 self.name = name 

883 s = k2 + kp2 

884 if s > 0 and s != _1_0: 

885 k2 = k2 / s # /= chokes PyChecker 

886 kp2 = kp2 / s 

887 if min(e2, k2, kp2) < 0 or not s > 0: 

888 raise TriaxialError(e2=e2, k2=k2, kp2=kp2) 

889 if e2: 

890 a = Radius_(a=_sqrt0(_1_0 + e2 * kp2) * b) if kp2 else b 

891 c = Radius_(c=_sqrt0(_1_0 - e2 * k2) * b) if k2 else b 

892 else: # spherical 

893 a = c = b 

894 if not (_isfinite(b) and a >= b >= c > 0): 

895 raise TriaxialError(b=b, a=a, c=c, e2=e2, 

896 k2=k2, kp2=kp2, txt=_not_ordered_) 

897 self._abc3 = a, b, c 

898 self._e2_k2_kp2 = e2, k2, kp2 

899 

900 @property_RO 

901 def isBiaxial(self): 

902 '''Is this triaxial I{biaxial} (C{bool}), C{a} == C{b} or C{b} == C{c} or C{a} == C{c}? 

903 ''' 

904 return self.isOblate or self.isProlate or self.a == self.c 

905 

906 @property_RO 

907 def isOblate(self): 

908 '''Is this triaxial I{oblate} (C{bool}), C{a} == C{b}? 

909 ''' 

910 return bool(self.kp2 == 0) 

911 

912 @property_RO 

913 def isProlate(self): 

914 '''Is this triaxial I{prolate} (C{bool}), C{b} == C{c}? 

915 ''' 

916 return bool(self.k2 == 0) 

917 

918 @property_RO 

919 def isTriaxial(self): 

920 '''Is this triaxial I{triaxial} (C{bool}), C{a} != C{b} and C{b} != C{c} and c{a} != C{c}? 

921 ''' 

922 return not self.isBiaxial 

923 

924 @Property_RO 

925 def _k_kp(self): 

926 '''(INTERNAL) Get the oblate C{k} and prolate C{kp} parameters. 

927 ''' 

928 return map1(_sqrt0, *self._k2_kp2) 

929 

930 @Property_RO 

931 def k2(self): 

932 '''(INTERNAL) Get the oblate C{k2} parameter I{squared}. 

933 ''' 

934 k2, _ = self._k2_kp2 

935 return k2 

936 

937 @Property_RO 

938 def _k2_kp2(self): 

939 '''(INTERNAL) Get the oblate C{k2} and prolate C{kp2} parameters I{squared}. 

940 ''' 

941 if self._e2_k2_kp2: 

942 _, k2, kp2 = self._e2_k2_kp2 

943 else: 

944 s = self._a2c2 

945 k2 = (self._b2c2 / s) if s else _1_0 

946 kp2 = (self._a2b2 / s) if s else _0_0 

947 return k2, kp2 

948 

949 @Property_RO 

950 def kp2(self): 

951 '''(INTERNAL) Get the prolate C{kp2} parameter I{squared}. 

952 ''' 

953 _, kp2 = self._k2_kp2 

954 return kp2 

955 

956 @Property_RO 

957 def _lcc23(self): 

958 return self._a2c2, self._b2c2, _0_0 

959 

960 @property_doc_(" prime-meridian rotation, longitude of the I{earth}'s major semi-axis C{a}, (L{Ang}), Karney's C{Triaxial_Earth_lon0}.") 

961 def Lon0(self): 

962 if self._Lon0 is None: 

963 WGS84_3 = self.name.startswith('WGS84_3') 

964 self.Lon0 = self._lon0WGS84_3 if WGS84_3 else 0 

965 return self._Lon0 

966 

967 @Lon0.setter # PYCHOK setter! 

968 def Lon0(self, lon0): 

969 A = _MODS.angles.Ang 

970 n = _Triaxial3Base.Lon0.name 

971 self._Lon0 = A(lon0, unit=Degrees, name=n) 

972 

973 @Property_RO 

974 def _xE(self): 

975 '''(INTERNAL) Get the x-elliptic function. 

976 ''' 

977 return self._xyE(self.e2, self.k2, self.kp2) 

978 

979 def _xyE(self, e2, k2, kp2): 

980 '''(INTERNAL) Helper for C{._xE} and C{._yE}. 

981 ''' 

982 if e2: 

983 a2 = -kp2 * e2 

984 ap2 = _1_0 - a2 

985 kp2 *= _1_0 - k2 * e2 

986 k2 *= ap2 

987 else: 

988 a2, ap2 = _0_0, _1_0 

989 return _MODS.elliptic.Elliptic(kp2, a2, k2, ap2) 

990 

991 @Property_RO 

992 def _yE(self): 

993 '''(INTERNAL) Get the y-elliptic function. 

994 ''' 

995 return self._xyE(-self.e2, self.kp2, self.k2) 

996 

997 

998class _TriaxialsBase(_NamedEnum): 

999 '''(INTERNAL) C{Triaxial*} registry, I{must} be a sub-class 

1000 to accommodate the L{_LazyNamedEnumItem} properties. 

1001 ''' 

1002 _assert_kwds = {} # like propertyROnce 

1003 _Triaxial = None # must be overloaded 

1004 

1005 def _Lazy(self, *abc, **name): 

1006 '''(INTERNAL) Instantiate the C{self._Triaxial}. 

1007 ''' 

1008 return self._Triaxial(*abc, **name) 

1009 

1010 def _assert(self): # PYCHOK signature 

1011 kwds = _TriaxialsBase._assert_kwds 

1012 if not kwds: 

1013 _lazy = _MODS.named._lazyNamedEnumItem 

1014 EWGS84 = _MODS.ellipsoids._EWGS84 

1015 abc84_35 = map1(m2km, EWGS84.a + 35, EWGS84.a - 35, EWGS84.b) 

1016 # <https://ArxIV.org/pdf/1909.06452.pdf> Table 1 Semi-axes in Km 

1017 # <https://www.JPS.NASA.gov/education/images/pdf/ss-moons.pdf> 

1018 # <https://link.Springer.com/article/10.1007/s00190-022-01650-9> 

1019 # <https://GeographicLib.SourceForge.io/C++/doc/classGeographicLib_1_1Constants.html> 

1020 # <https://www.ResearchGate.net/publication/344992491_Fitting_a_triaxial_ellipsoid_to_a_geoid_model> 

1021 for n, abc in dict( # a (Km) b (Km) c (Km) planet 

1022 Amalthea= (125.0, 73.0, 64.0), # Jupiter 

1023 Ariel= (581.1, 577.9, 577.7), # Uranus 

1024 Earth= (6378.173435, 6378.1039, 6356.7544), 

1025 Enceladus=(256.6, 251.4, 248.3), # Saturn 

1026 Europa= (1564.13, 1561.23, 1560.93), # Jupiter 

1027 Io= (1829.4, 1819.3, 1815.7), # Jupiter 

1028 Mars= (3394.6, 3393.3, 3376.3), 

1029 Mimas= (207.4, 196.8, 190.6), # Saturn 

1030 Miranda= (240.4, 234.2, 232.9), # Uranus 

1031 Moon= (1735.55, 1735.324, 1734.898), # Earth 

1032 Tethys= (535.6, 528.2, 525.8), # Saturn 

1033 WGS84_3= (6378.17136, 6378.10161, 6356.75184), # C++ 

1034 WGS84_3r=(6378.172, 6378.102, 6356.752), # C++, rounded 

1035# Panou= (6378.17188, 6378.10203, 6356.75224), # et.al. Fitting ... 

1036 WGS84_35=abc84_35).items(): 

1037 kwds[n] = _lazy(n, *map(km2m, abc)) 

1038 _NamedEnum._assert(self, **kwds) 

1039 

1040 

1041def _getitems(items, *indices): 

1042 '''(INTERNAL) Get the C{items} at the given I{indices}. 

1043 

1044 @return: C{Type(items[i] for i in indices)} with 

1045 C{Type = type(items)}, any C{type} having 

1046 the special method C{__getitem__}. 

1047 ''' 

1048 return type(items)(map(items.__getitem__, indices)) 

1049 

1050 

1051def _hypot2_1(x, y, z=0): 

1052 '''(INTERNAL) Compute M{x**2 + y**2 + z**2 - 1} with C{max(fabs(x), fabs(y), 

1053 fabs(z))} rarely greater than 1.0. 

1054 ''' 

1055 return fsumf_(_N_1_0, x*x, y*y, z*z) 

1056 

1057 

1058def _otherV3d_(x_xyz, y, z, **name): 

1059 '''(INTERNAL) Get a Vector3d from C{x_xyz}, C{y} and C{z}. 

1060 ''' 

1061 return _otherV3d(x_xyz=x_xyz, **name) if y is z is None else \ 

1062 Vector3d(x_xyz, y, z, **name) 

1063 

1064 

1065def _over0(p, q): 

1066 '''(INTERNAL) Return C{p / q} or C{0}. 

1067 ''' 

1068 return (p / q) if q > fabs(p) else _0_0 

1069 

1070 

1071def _over02(p, q): 

1072 '''(INTERNAL) Return C{(p / q)**2} or C{0}. 

1073 ''' 

1074 return (p / q)**2 if p and q else _0_0 

1075 

1076 

1077def _sqrt0(x): 

1078 '''(INTERNAL) C{sqrt0} with C{TriaxialError}. 

1079 ''' 

1080 return sqrt0(x, Error=TriaxialError) 

1081 

1082 

1083__all__ += _ALL_DOCS(_OrderedTriaxialBase, _Triaxial3Base, _UnOrderedTriaxialBase) 

1084 

1085# **) MIT License 

1086# 

1087# Copyright (C) 2025-2026 -- mrJean1 at Gmail -- All Rights Reserved. 

1088# 

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

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

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

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

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

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

1095# 

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

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

1098# 

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

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

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

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

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

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

1105# OTHER DEALINGS IN THE SOFTWARE.