Coverage for pygeodesy/ltp.py: 96%
419 statements
« prev ^ index » next coverage.py v7.2.2, created at 2023-11-12 13:23 -0500
« prev ^ index » next coverage.py v7.2.2, created at 2023-11-12 13:23 -0500
2# -*- coding: utf-8 -*-
4u'''I{Local Tangent Plane} (LTP) and I{local} cartesian coordinates.
6I{Local cartesian} and I{local tangent plane} classes L{LocalCartesian}, approximations L{ChLVa}
7and L{ChLVe} and L{Ltp}, L{ChLV}, L{LocalError}, L{Attitude} and L{Frustum}.
9@see: U{Local tangent plane coordinates<https://WikiPedia.org/wiki/Local_tangent_plane_coordinates>}
10 and class L{LocalCartesian}, transcoded from I{Charles Karney}'s C++ classU{LocalCartesian
11 <https://GeographicLib.SourceForge.io/C++/doc/classGeographicLib_1_1LocalCartesian.html>}.
12'''
13# make sure int/int division yields float quotient, see .basics
14from __future__ import division as _; del _ # PYCHOK semicolon
16from pygeodesy.basics import isscalar, issubclassof, map1, _xargs_names
17from pygeodesy.constants import EPS, INT0, _umod_360, _0_0, _0_01, _0_5, _1_0, \
18 _2_0, _60_0, _90_0, _100_0, _180_0, _3600_0, \
19 _N_1_0 # PYCHOK used!
20from pygeodesy.datums import _WGS84, _xinstanceof
21from pygeodesy.ecef import _EcefBase, EcefKarney, _llhn4, _xyzn4
22from pygeodesy.errors import _NotImplementedError, _TypesError, _ValueError, \
23 _xattr, _xkwds, _xkwds_get
24from pygeodesy.fmath import fabs, fdot, Fhorner
25from pygeodesy.fsums import _floor, Fsum, fsumf_, fsum1f_
26from pygeodesy.interns import NN, _0_, _COMMASPACE_, _DOT_, _ecef_, _height_, \
27 _invalid_, _lat0_, _lon0_, _ltp_, _M_, _name_, _too_
28# from pygeodesy.lazily import _ALL_LAZY # from vector3d
29from pygeodesy.ltpTuples import Attitude4Tuple, ChLVEN2Tuple, ChLV9Tuple, \
30 ChLVYX2Tuple, Footprint5Tuple, Local9Tuple, \
31 ChLVyx2Tuple, _XyzLocals4, _XyzLocals5, Xyz4Tuple
32from pygeodesy.named import _NamedBase, notOverloaded
33from pygeodesy.namedTuples import LatLon3Tuple, LatLon4Tuple, Vector3Tuple
34from pygeodesy.props import Property, Property_RO, property_doc_, property_RO, \
35 _update_all
36from pygeodesy.streprs import Fmt, strs, unstr
37from pygeodesy.units import Bearing, Degrees, Meter
38from pygeodesy.utily import cotd, _loneg, sincos2d, sincos2d_, tand, tand_, \
39 wrap180, wrap360
40from pygeodesy.vector3d import _ALL_LAZY, Vector3d
42# from math import fabs, floor as _floor # from .fmath, .fsums
44__all__ = _ALL_LAZY.ltp
45__version__ = '23.09.22'
47_height0_ = _height_ + _0_
48_narrow_ = 'narrow'
49_wide_ = 'wide'
50_Xyz_ = 'Xyz'
53def _fov_2(**fov):
54 # Half a field-of-view angle in C{degrees}.
55 f = Degrees(Error=LocalError, **fov) * _0_5
56 if EPS < f < _90_0:
57 return f
58 t = _invalid_ if f < 0 else _too_(_wide_ if f > EPS else _narrow_)
59 raise LocalError(txt=t, **fov)
62class Attitude(_NamedBase):
63 '''The orientation of a plane or camera in space.
64 '''
65 _alt = Meter( alt =_0_0)
66 _roll = Degrees(roll=_0_0)
67 _tilt = Degrees(tilt=_0_0)
68 _yaw = Bearing(yaw =_0_0)
70 def __init__(self, alt_attitude=INT0, tilt=INT0, yaw=INT0, roll=INT0, name=NN):
71 '''New L{Attitude}.
73 @kwarg alt_attitude: An altitude (C{meter}) above earth or an attitude
74 (L{Attitude} or L{Attitude4Tuple}) with the
75 C{B{alt}itude}, B{C{tilt}}, B{C{yaw}} and B{C{roll}}.
76 @kwarg tilt: Pitch, elevation from horizontal (C{degrees180}), negative down
77 (clockwise rotation along and around the x- or East axis).
78 @kwarg yaw: Bearing, heading (compass C{degrees360}), clockwise from North
79 (counter-clockwise rotation along and around the z- or Up axis).
80 @kwarg roll: Roll, bank (C{degrees180}), positive to the right and down
81 (clockwise rotation along and around the y- or North axis).
82 @kwarg name: Optional name C{str}).
84 @raise AttitudeError: Invalid B{C{alt_attitude}}, B{C{tilt}}, B{C{yaw}} or
85 B{C{roll}}.
87 @see: U{Principal axes<https://WikiPedia.org/wiki/Aircraft_principal_axes>} and
88 U{Yaw, pitch, and roll rotations<http://MSL.CS.UIUC.edu/planning/node102.html>}.
89 '''
90 if isscalar(alt_attitude):
91 t = Attitude4Tuple(alt_attitude, tilt, yaw, roll)
92 else:
93 try:
94 t = alt_attitude.atyr
95 except AttributeError:
96 raise AttitudeError(alt=alt_attitude, tilt=tilt, yaw=yaw, rol=roll)
97 for n, v in t.items():
98 if v:
99 setattr(self, n, v)
100 n = name or t.name
101 if n:
102 self.name = n
104 @property_doc_(' altitude above earth in C{meter}.')
105 def alt(self):
106 return self._alt
108 @alt.setter # PYCHOK setter!
109 def alt(self, alt): # PYCHOK no cover
110 a = Meter(alt=alt, Error=AttitudeError)
111 if self._alt != a:
112 _update_all(self)
113 self._alt = a
115 altitude = alt
117 @Property_RO
118 def atyr(self):
119 '''Return this attitude's alt[itude], tilt, yaw and roll as an L{Attitude4Tuple}.
120 '''
121 return Attitude4Tuple(self.alt, self.tilt, self.yaw, self.roll, name=self.name)
123 @Property_RO
124 def matrix(self):
125 '''Get the 3x3 rotation matrix C{R(yaw)·R(tilt)·R(roll)}, aka I{ZYX} (C{float}, row-order).
127 @see: The matrix M of case 10 in U{Appendix A
128 <https://ntrs.NASA.gov/api/citations/19770019231/downloads/19770019231.pdf>}.
129 '''
130 def _5to3(x, y, _y, z, _z):
131 return x, fsum1f_(y, _y), fsum1f_(z, _z)
133 r0, r1, r2 = self._rows3
134 return _5to3(*r0), _5to3(*r1), r2
136 @property_doc_(' roll/bank in C{degrees180}, positive to the right and down.')
137 def roll(self):
138 return self._roll
140 @roll.setter # PYCHOK setter!
141 def roll(self, roll):
142 r = Degrees(roll=roll, wrap=wrap180, Error=AttitudeError)
143 if self._roll != r:
144 _update_all(self)
145 self._roll = r
147 bank = roll
149 @Property_RO
150 def _rows3(self):
151 # to follow the definitions of rotation angles alpha, beta and gamma:
152 # negate yaw since yaw is counter-clockwise around the z-axis, swap
153 # tilt and roll since tilt is around the x- and roll around the y-axis
154 sa, ca, sb, cb, sg, cg = sincos2d_(-self.yaw, self.roll, self.tilt)
155 return ((ca * cb, ca * sb * sg, -sa * cg, ca * sb * cg, sa * sg),
156 (sa * cb, sa * sb * sg, ca * cg, sa * sb * cg, -ca * sg),
157 ( -sb, cb * sg, cb * cg))
159 def rotate(self, x_xyz, y=None, z=None, Vector=None, **Vector_kwds):
160 '''Transform a (local) cartesian by this attitude's matrix.
162 @arg x_xyz: X component of vector (C{scalar}) or (3-D) vector
163 (C{Cartesian}, L{Vector3d} or L{Vector3Tuple}).
164 @kwarg y: Y component of vector (C{scalar}), same units as B{C{x}}.
165 @kwarg z: Z component of vector (C{scalar}), same units as B{C{x}}.
166 @kwarg Vector: Class to return transformed point (C{Cartesian},
167 L{Vector3d} or C{Vector3Tuple}) or C{None}.
168 @kwarg Vector_kwds: Optional, additional B{C{Vector}} keyword arguments,
169 ignored if C{B{Vector} is None}.
171 @return: A B{C{Vector}} instance or a L{Vector3Tuple}C{(x, y, z)} if
172 C{B{Vector}=None}.
174 @raise AttitudeError: Invalid B{C{x_xyz}}, B{C{y}} or B{C{z}}.
176 @see: U{Yaw, pitch, and roll rotations<http://MSL.CS.UIUC.edu/planning/node102.html>}.
177 '''
178 try:
179 try:
180 x, y, z = map( float, x_xyz.xyz)
181 except AttributeError:
182 x, y, z = map1(float, x_xyz, y, z)
183 except (TypeError, ValueError) as x:
184 raise AttitudeError(x_xyz=x_xyz, y=y, z=z, cause=x)
186 X, Y, Z = self._rows3
187 X = fdot(X, x, y, y, z, z)
188 Y = fdot(Y, x, y, y, z, z)
189 Z = fdot(Z, x, y, z)
190 return Vector3Tuple(X, Y, Z, name=self.name) if Vector is None else \
191 Vector(X, Y, Z, **_xkwds(Vector_kwds, name=self.name))
193 @property_doc_(' tilt/pitch/elevation from horizontal in C{degrees180}, negative down.')
194 def tilt(self):
195 return self._tilt
197 @tilt.setter # PYCHOK setter!
198 def tilt(self, tilt):
199 t = Degrees(tilt=tilt, wrap=wrap180, Error=AttitudeError)
200 if self._tilt != t:
201 _update_all(self)
202 self._tilt = t
204 elevation = pitch = tilt
206 def toStr(self, prec=6, sep=_COMMASPACE_, **unused): # PYCHOK signature
207 '''Format this attitude as string.
209 @kwarg prec: The C{float} precision, number of decimal digits (0..9).
210 Trailing zero decimals are stripped for B{C{prec}} values
211 of 1 and above, but kept for negative B{C{prec}} values.
212 @kwarg sep: Separator to join (C{str}).
214 @return: This attitude (C{str}).
215 '''
216 return self.atyr.toStr(prec=prec, sep=sep)
218 @Property_RO
219 def tyr3d(self):
220 '''Get this attitude's (3-D) directional vector (L{Vector3d}).
222 @see: U{Yaw, pitch, and roll rotations<http://MSL.CS.UIUC.edu/planning/node102.html>}.
223 '''
224 def _r2d(r):
225 return fsumf_(_N_1_0, *r)
227 return Vector3d(*map(_r2d, self._rows3), name=tyr3d.__name__)
229 @property_doc_(' yaw/bearing/heading in compass C{degrees360}, clockwise from North.')
230 def yaw(self):
231 return self._yaw
233 @yaw.setter # PYCHOK setter!
234 def yaw(self, yaw):
235 y = Bearing(yaw=yaw, Error=AttitudeError)
236 if self._yaw != y:
237 _update_all(self)
238 self._yaw = y
240 bearing = heading = yaw
243class AttitudeError(_ValueError):
244 '''An L{Attitude} or L{Attitude4Tuple} issue.
245 '''
246 pass
249class Frustum(_NamedBase):
250 '''A rectangular pyramid, typically representing a camera's I{field-of-view}
251 (fov) and the intersection with (or projection to) a I{local tangent plane}.
253 @see: U{Viewing frustum<https://WikiPedia.org/wiki/Viewing_frustum>}.
254 '''
255 _h_2 = _0_0 # half hfov in degrees
256 _ltp = None # local tangent plane
257 _tan_h_2 = _0_0 # tan(_h_2)
258 _v_2 = _0_0 # half vfov in degrees
260 def __init__(self, hfov, vfov, ltp=None):
261 '''New L{Frustum}.
263 @arg hfov: Horizontal field-of-view (C{degrees180}).
264 @arg vfov: Vertical field-of-view (C{degrees180}).
265 @kwarg ltp: Optional I{local tangent plane} (L{Ltp}).
267 @raise LocalError: Invalid B{C{hfov}} or B{C{vfov}}.
268 '''
269 self._h_2 = h = _fov_2(hfov=hfov)
270 self._v_2 = _fov_2(vfov=vfov)
272 self._tan_h_2 = tand(h, fov_2=h)
274 if ltp:
275 self._ltp = _xLtp(ltp)
277 def footprint5(self, alt_attitude, tilt=0, yaw=0, roll=0, z=_0_0, ltp=None): # MCCABE 15
278 '''Compute the center and corners of the intersection with (or projection
279 to) the I{local tangent plane} (LTP).
281 @arg alt_attitude: An altitude (C{meter}) above I{local tangent plane} or
282 an attitude (L{Attitude} or L{Attitude4Tuple}) with the
283 C{B{alt}itude}, B{C{tilt}}, B{C{yaw}} and B{C{roll}}.
284 @kwarg tilt: Pitch, elevation from horizontal (C{degrees}), negative down
285 (clockwise rotation along and around the x- or East axis).
286 @kwarg yaw: Bearing, heading (compass C{degrees}), clockwise from North
287 (counter-clockwise rotation along and around the z- or Up axis).
288 @kwarg roll: Roll, bank (C{degrees}), positive to the right and down
289 (clockwise rotation along and around the y- or North axis).
290 @kwarg z: Optional height of the footprint (C{meter}) above I{local tangent plane}.
291 @kwarg ltp: The I{local tangent plane} (L{Ltp}), overriding this
292 frustum's C{ltp}.
294 @return: A L{Footprint5Tuple}C{(center, upperleft, upperight, loweright,
295 lowerleft)} with the C{center} and 4 corners, each an L{Xyz4Tuple}.
297 @raise TypeError: Invalid B{C{ltp}}.
299 @raise UnitError: Invalid B{C{altitude}}, B{C{tilt}}, B{C{roll}} or B{C{z}}.
301 @raise ValueError: If B{C{altitude}} too low, B{C{z}} too high or B{C{tilt}}
302 or B{C{roll}} -including B{C{vfov}} respectively B{C{hfov}}-
303 over the horizon.
305 @see: U{Principal axes<https://WikiPedia.org/wiki/Aircraft_principal_axes>}.
306 '''
307 def _xy2(a, e, h_2, tan_h_2, r):
308 # left and right corners, or swapped
309 if r < EPS: # no roll
310 r = a * tan_h_2
311 l = -r # PYCHOK l is ell
312 else: # roll
313 r, l = tand_(r - h_2, r + h_2, roll_hfov=r) # PYCHOK l is ell
314 r *= -a # negate right positive
315 l *= -a # PYCHOK l is ell
316 y = a * cotd(e, tilt_vfov=e)
317 return (l, y), (r, y)
319 def _xyz5(b, xy5, z, ltp):
320 # rotate (x, y)'s by bearing, clockwise
321 s, c = sincos2d(b)
322 for x, y in xy5:
323 yield Xyz4Tuple(fsum1f_(x * c, y * s),
324 fsum1f_(y * c, -x * s), z, ltp)
326 try:
327 a, t, y, r = alt_attitude.atyr
328 except AttributeError:
329 a, t, y, r = alt_attitude, tilt, yaw, roll
331 a = Meter(altitude=a)
332 if a < EPS: # too low
333 raise _ValueError(altitude=a)
334 if z: # PYCHOK no cover
335 z = Meter(z=z)
336 a -= z
337 if a < EPS: # z above a
338 raise _ValueError(altitude_z=a)
339 else:
340 z = _0_0
342 b = Degrees(yaw=y, wrap=wrap360) # bearing
343 e = -Degrees(tilt=t, wrap=wrap180) # elevation, pitch
344 if not EPS < e < _180_0:
345 raise _ValueError(tilt=t)
346 if e > _90_0:
347 e = _loneg(e)
348 b = _umod_360(b + _180_0)
350 r = Degrees(roll=r, wrap=wrap180) # roll center
351 x = (-a * tand(r, roll=r)) if r else _0_0
352 y = a * cotd(e, tilt=t) # ground range
353 if fabs(y) < EPS:
354 y = _0_0
356 # center and corners, clockwise from upperleft, rolled
357 xy5 = ((x, y),) + _xy2(a, e - self._v_2, self._h_2, self._tan_h_2, r) \
358 + _xy2(a, e + self._v_2, -self._h_2, -self._tan_h_2, r) # swapped
359 # turn center and corners by yaw, clockwise
360 p = self.ltp if ltp is None else ltp # None OK
361 return Footprint5Tuple(_xyz5(b, xy5, z, p)) # *_xyz5
363 @Property_RO
364 def hfov(self):
365 '''Get the horizontal C{fov} (C{degrees}).
366 '''
367 return Degrees(hfov=self._h_2 * _2_0)
369 @Property_RO
370 def ltp(self):
371 '''Get the I{local tangent plane} (L{Ltp}) or C{None}.
372 '''
373 return self._ltp
375 def toStr(self, prec=3, fmt=Fmt.F, sep=_COMMASPACE_): # PYCHOK signature
376 '''Convert this frustum to a "hfov, vfov, ltp" string.
378 @kwarg prec: Number of (decimal) digits, unstripped (0..8 or C{None}).
379 @kwarg fmt: Optional, C{float} format (C{str}).
380 @kwarg sep: Separator to join (C{str}).
382 @return: Frustum in the specified form (C{str}).
383 '''
384 t = self.hfov, self.vfov
385 if self.ltp:
386 t += self.ltp,
387 t = strs(t, prec=prec, fmt=fmt)
388 return sep.join(t) if sep else t
390 @Property_RO
391 def vfov(self):
392 '''Get the vertical C{fov} (C{degrees}).
393 '''
394 return Degrees(vfov=self._v_2 * _2_0)
397class LocalError(_ValueError):
398 '''A L{LocalCartesian} or L{Ltp} related issue.
399 '''
400 pass
403class LocalCartesian(_NamedBase):
404 '''Conversion between geodetic C{(lat, lon, height)} and I{local
405 cartesian} C{(x, y, z)} coordinates with I{geodetic} origin
406 C{(lat0, lon0, height0)}, transcoded from I{Karney}'s C++ class
407 U{LocalCartesian<https://GeographicLib.SourceForge.io/C++/doc/
408 classGeographicLib_1_1LocalCartesian.html>}.
410 The C{z} axis is normal to the ellipsoid, the C{y} axis points due
411 North. The plane C{z = -height0} is tangent to the ellipsoid.
413 The conversions all take place via geocentric coordinates using a
414 geocentric L{EcefKarney}, by default the WGS84 datum/ellipsoid.
416 @see: Class L{Ltp}.
417 '''
418 _ecef = EcefKarney(_WGS84)
419 _EcefClass = EcefKarney
420 _lon00 = INT0 # self.lon0
421 _t0 = None # origin (..., lat0, lon0, height0, ...) L{Ecef9Tuple}
422 _9Tuple = Local9Tuple
424 def __init__(self, latlonh0=INT0, lon0=INT0, height0=INT0, ecef=None, name=NN, **lon00):
425 '''New L{LocalCartesian} converter.
427 @kwarg latlonh0: The (geodetic) origin (C{LatLon}, L{LatLon4Tuple}, L{Ltp}
428 L{LocalCartesian} or L{Ecef9Tuple}) or the C{scalar}
429 latitude of the (goedetic) origin (C{degrees}).
430 @kwarg lon0: Longitude of the (goedetic) origin (C{degrees}) for C{scalar}
431 B{C{latlonh0}}, ignored otherwise.
432 @kwarg height0: Optional height (C{meter}, conventionally) at the (goedetic)
433 origin perpendicular to and above (or below) the ellipsoid's
434 surface and for C{scalar} B{C{latlonh0}}, ignored otherwise.
435 @kwarg ecef: An ECEF converter (L{EcefKarney} I{only}) for C{scalar}
436 B{C{latlonh0}}, ignored otherwise.
437 @kwarg name: Optional name (C{str}).
438 @kwarg lon00: An arbitrary, I{polar} longitude (C{degrees}), overriding
439 the default C{B{lon00}=B{lon0}}, see method C{reverse}.
441 @raise LocalError: If B{C{latlonh0}} not C{LatLon}, L{LatLon4Tuple}, L{Ltp},
442 L{LocalCartesian} or L{Ecef9Tuple} or B{C{latlonh0}},
443 B{C{lon0}}, B{C{height0}} or B{C{lon00}} invalid.
445 @raise TypeError: Invalid B{C{ecef}} or not L{EcefKarney}.
447 @note: If BC{latlonh0} is an L{Ltp} or L{LocalCartesian}, only C{lat0}, C{lon0},
448 C{height0} and I{polar} C{lon00} are copied, I{not} the ECEF converter.
449 '''
450 self.reset(latlonh0, lon0=lon0, height0=height0, ecef=ecef, name=name, **lon00)
452 def __eq__(self, other):
453 '''Compare this and an other instance.
455 @arg other: The other ellipsoid (L{LocalCartesian} or L{Ltp}).
457 @return: C{True} if equal, C{False} otherwise.
458 '''
459 return other is self or (isinstance(other, self.__class__) and
460 other.ecef == self.ecef and
461 other._t0 == self._t0)
463 @Property_RO
464 def datum(self):
465 '''Get the ECEF converter's datum (L{Datum}).
466 '''
467 return self.ecef.datum
469 @Property_RO
470 def ecef(self):
471 '''Get the ECEF converter (L{EcefKarney}).
472 '''
473 return self._ecef
475 def _ecef2local(self, ecef, Xyz, Xyz_kwds):
476 '''(INTERNAL) Convert geocentric/geodetic to local, like I{forward}.
478 @arg ecef: Geocentric (and geodetic) (L{Ecef9Tuple}).
479 @arg Xyz: An L{XyzLocal}, L{Enu} or L{Ned} I{class} or C{None}.
480 @arg Xyz_kwds: B{C{Xyz}} keyword arguments, ignored if C{B{Xyz} is None}.
482 @return: An C{B{Xyz}(x, y, z, ltp, **B{Xyz_kwds}} instance or if
483 C{B{Xyz} is None}, a L{Local9Tuple}C{(x, y, z, lat, lon,
484 height, ltp, ecef, M)} with this C{ltp}, B{C{ecef}}
485 (L{Ecef9Tuple}) converted to this C{datum} and C{M=None},
486 always.
487 '''
488 ltp = self
489 if ecef.datum != ltp.datum:
490 ecef = ecef.toDatum(ltp.datum)
491 x, y, z = self.M.rotate(ecef.xyz, *ltp._t0_xyz)
492 r = Local9Tuple(x, y, z, ecef.lat, ecef.lon, ecef.height,
493 ltp, ecef, None, name=ecef.name)
494 if Xyz:
495 if not issubclassof(Xyz, *_XyzLocals4): # Vector3d
496 raise _TypesError(_Xyz_, Xyz, *_XyzLocals4)
497 r = r.toXyz(Xyz=Xyz, **Xyz_kwds)
498 return r
500 @Property_RO
501 def ellipsoid(self):
502 '''Get the ECEF converter's ellipsoid (L{Ellipsoid}).
503 '''
504 return self.ecef.datum.ellipsoid
506 def forward(self, latlonh, lon=None, height=0, M=False, name=NN):
507 '''Convert I{geodetic} C{(lat, lon, height)} to I{local} cartesian
508 C{(x, y, z)}.
510 @arg latlonh: Either a C{LatLon}, L{Ltp}, L{Ecef9Tuple} or C{scalar}
511 (geodetic) latitude (C{degrees}).
512 @kwarg lon: Optional C{scalar} (geodetic) longitude for C{scalar}
513 B{C{latlonh}} (C{degrees}).
514 @kwarg height: Optional height (C{meter}, conventionally) perpendicular
515 to and above (or below) the ellipsoid's surface.
516 @kwarg M: Optionally, return the I{concatenated} rotation L{EcefMatrix},
517 iff available (C{bool}).
518 @kwarg name: Optional name (C{str}).
520 @return: A L{Local9Tuple}C{(x, y, z, lat, lon, height, ltp, ecef, M)}
521 with I{local} C{x}, C{y}, C{z}, I{geodetic} C{(lat}, C{lon},
522 C{height}, this C{ltp}, C{ecef} (L{Ecef9Tuple}) with
523 I{geocentric} C{x}, C{y}, C{z} (and I{geodetic} C{lat},
524 C{lon}, C{height}) and the I{concatenated} rotation matrix
525 C{M} (L{EcefMatrix}) if requested.
527 @raise LocalError: If B{C{latlonh}} not C{scalar}, C{LatLon}, L{Ltp},
528 L{Ecef9Tuple} or invalid or if B{C{lon}} not
529 C{scalar} for C{scalar} B{C{latlonh}} or invalid
530 or if B{C{height}} invalid.
531 '''
532 lat, lon, h, n = _llhn4(latlonh, lon, height, Error=LocalError, name=name)
533 t = self.ecef._forward(lat, lon, h, n, M=M)
534 x, y, z = self.M.rotate(t.xyz, *self._t0_xyz)
535 m = self.M.multiply(t.M) if M else None
536 return self._9Tuple(x, y, z, lat, lon, h, self, t, m, name=n or self.name)
538 @Property_RO
539 def height0(self):
540 '''Get the origin's height (C{meter}).
541 '''
542 return self._t0.height
544 @Property_RO
545 def lat0(self):
546 '''Get the origin's latitude (C{degrees}).
547 '''
548 return self._t0.lat
550 @Property_RO
551 def latlonheight0(self):
552 '''Get the origin's lat-, longitude and height (L{LatLon3Tuple}C{(lat, lon, height)}).
553 '''
554 return LatLon3Tuple(self.lat0, self.lon0, self.height0, name=self.name)
556 def _local2ecef(self, local, nine=False, M=False):
557 '''(INTERNAL) Convert I{local} to geocentric/geodetic, like I{.reverse}.
559 @arg local: Local (L{XyzLocal}, L{Enu}, L{Ned}, L{Aer} or L{Local9Tuple}).
560 @kwarg nine: Return 3- or 9-tuple (C{bool}).
561 @kwarg M: Include the rotation matrix (C{bool}).
563 @return: A I{geocentric} 3-tuple C{(x, y, z)} or if C{B{nine}=True},
564 an L{Ecef9Tuple}C{(x, y, z, lat, lon, height, C, M, datum)},
565 optionally including rotation matrix C{M} or C{None}.
566 '''
567 t = self.M.unrotate(local.xyz, *self._t0_xyz)
568 if nine:
569 t = self.ecef.reverse(*t, M=M)
570 return t
572 @Property_RO
573 def lon0(self):
574 '''Get the origin's longitude (C{degrees}).
575 '''
576 return self._t0.lon
578 @Property
579 def lon00(self):
580 '''Get the arbitrary, I{polar} longitude (C{degrees}).
581 '''
582 return self._lon00
584 @lon00.setter # PYCHOK setter!
585 def lon00(self, lon00):
586 '''Set the arbitrary, I{polar} longitude (C{degrees}).
587 '''
588 # lon00 <https://GitHub.com/mrJean1/PyGeodesy/issues/77>
589 self._lon00 = Degrees(lon00=lon00)
591 @Property_RO
592 def M(self):
593 '''Get the rotation matrix (C{EcefMatrix}).
594 '''
595 return self._t0.M
597 def reset(self, latlonh0=INT0, lon0=INT0, height0=INT0, ecef=None, name=NN, **lon00):
598 '''Reset this converter, see L{LocalCartesian.__init__} and L{Ltp.__init__} for more details.
599 '''
600 if isinstance(latlonh0, LocalCartesian):
601 if self._t0:
602 _update_all(self)
603 self._ecef = latlonh0.ecef
604 self._lon00 = latlonh0.lon00
605 self._t0 = latlonh0._t0
606 n = name or latlonh0.name
607 else:
608 lat0, lon0, height0, n = _llhn4(latlonh0, lon0, height0, suffix=_0_,
609 Error=LocalError, name=name or self.name)
610 if ecef: # PYCHOK no cover
611 _xinstanceof(self._EcefClass, ecef=ecef)
612 _update_all(self)
613 self._ecef = ecef
614 elif self._t0:
615 _update_all(self)
616 self._t0 = self.ecef._forward(lat0, lon0, height0, n, M=True)
617 self.lon00 = _xattr(latlonh0, lon00=_xkwds_get(lon00, lon00=lon0))
618 if n:
619 self.rename(n)
621 def reverse(self, xyz, y=None, z=None, M=False, name=NN, **lon00):
622 '''Convert I{local} C{(x, y, z)} to I{geodetic} C{(lat, lon, height)}.
624 @arg xyz: A I{local} (L{XyzLocal}, L{Enu}, L{Ned}, L{Aer}, L{Local9Tuple}) or
625 local C{x} coordinate (C{scalar}).
626 @kwarg y: Local C{y} coordinate for C{scalar} B{C{xyz}} and B{C{z}} (C{meter}).
627 @kwarg z: Local C{z} coordinate for C{scalar} B{C{xyz}} and B{C{y}} (C{meter}).
628 @kwarg M: Optionally, return the I{concatenated} rotation L{EcefMatrix}, iff
629 available (C{bool}).
630 @kwarg name: Optional name (C{str}).
631 @kwarg lon00: An arbitrary, I{polar} longitude (C{degrees}), returned for local
632 C{B{x}=0} and C{B{y}=0} at I{polar} latitudes C{abs(B{lat0}) == 90},
633 overriding property C{lon00} and default C{B{lon00}=B{lon0}}.
635 @return: An L{Local9Tuple}C{(x, y, z, lat, lon, height, ltp, ecef, M)} with
636 I{local} C{x}, C{y}, C{z}, I{geodetic} C{lat}, C{lon}, C{height},
637 this C{ltp}, an C{ecef} (L{Ecef9Tuple}) with the I{geocentric} C{x},
638 C{y}, C{z} (and I{geodetic} C{lat}, C{lon}, C{height}) and the
639 I{concatenated} rotation matrix C{M} (L{EcefMatrix}) if requested.
641 @raise LocalError: Invalid B{C{xyz}} or C{scalar} C{x} or B{C{y}} and/or B{C{z}}
642 not C{scalar} for C{scalar} B{C{xyz}}.
643 '''
644 x, y, z, n = _xyzn4(xyz, y, z, _XyzLocals5, Error=LocalError, name=name)
645 c = self.M.unrotate((x, y, z), *self._t0_xyz)
646 t = self.ecef.reverse(*c, M=M, lon00=_xkwds_get(lon00, lon00=self.lon00))
647 m = self.M.multiply(t.M) if M else None
648 return self._9Tuple(x, y, z, t.lat, t.lon, t.height, self, t, m, name=n or self.name)
650 @Property_RO
651 def _t0_xyz(self):
652 '''(INTERNAL) Get C{(x0, y0, z0)} as L{Vector3Tuple}.
653 '''
654 return self._t0.xyz
656 def toStr(self, prec=9, **unused): # PYCHOK signature
657 '''Return this L{LocalCartesian} as a string.
659 @kwarg prec: Precision, number of (decimal) digits (0..9).
661 @return: This L{LocalCartesian} representation (C{str}).
662 '''
663 return self.attrs(_lat0_, _lon0_, _height0_, _M_, _ecef_, _name_, prec=prec)
666class Ltp(LocalCartesian):
667 '''A I{local tangent plan} (LTP), a sub-class of C{LocalCartesian} with
668 (re-)configurable ECEF converter.
669 '''
670 _EcefClass = _EcefBase
672 def __init__(self, latlonh0=INT0, lon0=INT0, height0=INT0, ecef=None, name=NN, **lon00):
673 '''New C{Ltp}, see L{LocalCartesian.__init__} for more details.
675 @kwarg ecef: Optional ECEF converter (L{EcefKarney}, L{EcefFarrell21},
676 L{EcefFarrell22}, L{EcefSudano}, L{EcefVeness} or
677 L{EcefYou} I{instance}), overriding the default
678 L{EcefKarney}C{(datum=Datums.WGS84)} for C{scalar}.
680 @raise TypeError: Invalid B{C{ecef}}.
681 '''
682 LocalCartesian.reset(self, latlonh0, lon0=lon0, height0=height0,
683 ecef=ecef, name=name, **lon00)
685 @Property
686 def ecef(self):
687 '''Get this LTP's ECEF converter (C{Ecef...} I{instance}).
688 '''
689 return self._ecef
691 @ecef.setter # PYCHOK setter!
692 def ecef(self, ecef):
693 '''Set this LTP's ECEF converter (C{Ecef...} I{instance}).
695 @raise TypeError: Invalid B{C{ecef}}.
696 '''
697 _xinstanceof(_EcefBase, ecef=ecef)
698 if ecef != self._ecef: # PYCHOK no cover
699 self.reset(self._t0)
700 self._ecef = ecef
703class _ChLV(object):
704 '''(INTERNAL) Base class for C{ChLV*} classes.
705 '''
706 _03_falsing = ChLVyx2Tuple(0.6e6, 0.2e6)
707# _92_falsing = ChLVYX2Tuple(2.0e6, 1.0e6) # _95_ - _03_
708 _95_falsing = ChLVEN2Tuple(2.6e6, 1.2e6)
710 def _ChLV9Tuple(self, fw, M, name, *Y_X_h_lat_lon_h):
711 '''(INTERNAL) Helper for C{ChLVa/e.forward} and C{.reverse}.
712 '''
713 if bool(M): # PYCHOK no cover
714 m = self.forward if fw else self.reverse # PYCHOK attr
715 n = _DOT_(self.__class__.__name__, m.__name__)
716 raise _NotImplementedError(unstr(n, M=M), txt=None)
717 t = Y_X_h_lat_lon_h + (self, self._t0, None) # PYCHOK _t0
718 return ChLV9Tuple(t, name=name)
720 @property_RO
721 def _enh_n_h(self):
722 '''(INTERNAL) Get C{ChLV*.reverse} args[1:4] names, I{once}.
723 '''
724 t = _xargs_names(_ChLV.reverse)[1:4]
725 _ChLV._enh_n_h = t # overwrite this property_RO
726 # assert _xargs_names( ChLV.reverse)[1:4] == t
727 # assert _xargs_names(ChLVa.reverse)[1:4] == t
728 # assert _xargs_names(ChLVe.reverse)[1:4] == t
729 return t
731 def forward(self, latlonh, lon=None, height=0, M=None, name=NN): # PYCHOK no cover
732 '''Convert WGS84 geodetic to I{Swiss} projection coordinates. I{Must be overloaded}.
734 @arg latlonh: Either a C{LatLon}, L{Ltp} or C{scalar} (geodetic) latitude (C{degrees}).
735 @kwarg lon: Optional, C{scalar} (geodetic) longitude for C{scalar} B{C{latlonh}} (C{degrees}).
736 @kwarg height: Optional, height, vertically above (or below) the surface of the ellipsoid
737 (C{meter}) for C{scalar} B{C{latlonh}} and B{C{lon}}.
738 @kwarg M: If C{True}, return the I{concatenated} rotation L{EcefMatrix} iff available
739 for C{ChLV} only, C{None} otherwise (C{bool}).
740 @kwarg name: Optional name (C{str}).
742 @return: A L{ChLV9Tuple}C{(Y, X, h_, lat, lon, height, ltp, ecef, M)} with the unfalsed
743 I{Swiss Y, X} coordinates, I{Swiss h_} height, the given I{geodetic} C{lat},
744 C{lon} and C{height}, this C{ChLV*} instance and C{ecef} (L{Ecef9Tuple}) at
745 I{Bern, Ch} and rotation matrix C{M}. The returned C{ltp} is this C{ChLV},
746 C{ChLVa} or C{ChLVe} instance.
748 @raise LocalError: Invalid or non-C{scalar} B{C{latlonh}}, B{C{lon}} or B{C{height}}.
749 '''
750 notOverloaded(self, latlonh, lon=lon, height=height, M=M, name=name)
752 def reverse(self, enh_, n=None, h_=0, M=None, **name): # PYCHOK no cover
753 '''Convert I{Swiss} projection to WGS84 geodetic coordinates. I{Must be overloaded}.
755 @arg enh_: A Swiss projection (L{ChLV9Tuple}) or the C{scalar}, falsed I{Swiss E_LV95}
756 or I{y_LV03} easting (C{meter}).
757 @kwarg n: Falsed I{Swiss N_LV85} or I{x_LV03} northing for C{scalar} B{C{enh_}} and
758 B{C{h_}} (C{meter}).
759 @kwarg h_: I{Swiss h'} height for C{scalar} B{C{enh_}} and B{C{n}} (C{meter}).
760 @kwarg M: If C{True}, return the I{concatenated} rotation L{EcefMatrix} iff available
761 for C{ChLV} only, C{None} otherwise (C{bool}).
762 @kwarg name: Optional name (C{str}).
764 @return: A L{ChLV9Tuple}C{(Y, X, h_, lat, lon, height, ltp, ecef, M)} with the unfalsed
765 I{Swiss Y, X} coordinates, I{Swiss h_} height, the given I{geodetic} C{lat},
766 C{lon} and C{height}, this C{ChLV*} instance and C{ecef} (L{Ecef9Tuple}) at
767 I{Bern, Ch} and rotation matrix C{M}. The returned C{ltp} is this C{ChLV},
768 C{ChLVa} or C{ChLVe} instance.
770 @raise LocalError: Invalid or non-C{scalar} B{C{enh_}}, B{C{n}} or B{C{h_}}.
771 '''
772 notOverloaded(self, enh_, n=n, h_=h_, M=M, **name)
774 @staticmethod
775 def _falsing2(LV95):
776 '''(INTERNAL) Get the C{LV95} or C{LV03} falsing.
777 '''
778 return _ChLV._95_falsing if LV95 in (True, 95) else (
779 _ChLV._03_falsing if LV95 in (False, 3) else ChLVYX2Tuple(0, 0))
781 @staticmethod
782 def _llh2abh_3(lat, lon, h):
783 '''(INTERNAL) Helper for C{ChLVa/e.forward}.
784 '''
785 def _deg2ab(deg, sLL):
786 # convert degrees to arc-seconds
787 def _dms(ds, p, q, swap):
788 d = _floor(ds)
789 t = (ds - d) * p
790 m = _floor(t)
791 s = (t - m) * p
792 if swap:
793 d, s = s, d
794 return d + (m + s * q) * q
796 s = _dms(deg, _60_0, _0_01, False) # deg2sexag
797 s = _dms( s, _100_0, _60_0, True) # sexag2asec
798 return (s - sLL) / ChLV._s_ab
800 a = _deg2ab(lat, ChLV._sLat) # phi', lat_aux
801 b = _deg2ab(lon, ChLV._sLon) # lam', lng_aux
802 h_ = fsumf_(h, -ChLV.Bern.height, 2.73 * b, 6.94 * a)
803 return a, b, h_
805 @staticmethod
806 def _YXh_2abh3(Y, X, h_):
807 '''(INTERNAL) Helper for C{ChLVa/e.reverse}.
808 '''
809 def _YX2ab(YX):
810 return YX * ChLV._ab_m
812 a, b = map1(_YX2ab, Y, X)
813 h = fsumf_(h_, ChLV.Bern.height, -12.6 * a, -22.64 * b)
814 return a, b, h
816 def _YXh_n4(self, enh_, n, h_, **name):
817 '''(INTERNAL) Helper for C{ChLV*.reverse}.
818 '''
819 Y, X, h_, name = _xyzn4(enh_, n, h_, ChLV9Tuple,
820 _xyz_y_z_names=self._enh_n_h, **name)
821 if isinstance(enh_, ChLV9Tuple):
822 Y, X = enh_.Y, enh_.X
823 else: # isscalar(enh_)
824 Y, X = ChLV.unfalse2(Y, X) # PYCHOK ChLVYX2Tuple
825 return Y, X, h_, name
828class ChLV(_ChLV, Ltp):
829 '''Conversion between I{WGS84 geodetic} and I{Swiss} projection coordinates using
830 L{pygeodesy.EcefKarney}'s Earth-Centered, Earth-Fixed (ECEF) methods.
832 @see: U{Swiss projection formulas<https://www.SwissTopo.admin.CH/en/maps-data-online/
833 calculation-services.html>}, page 7ff, U{NAVREF<https://www.SwissTopo.admin.CH/en/
834 maps-data-online/calculation-services/navref.html>}, U{REFRAME<https://www.SwissTopo.admin.CH/
835 en/maps-data-online/calculation-services/reframe.html>} and U{SwissTopo Scripts GPS WGS84
836 <-> LV03<https://GitHub.com/ValentinMinder/Swisstopo-WGS84-LV03>}.
837 '''
838 _9Tuple = ChLV9Tuple
840 _ab_d = 0.36 # a, b units per degree, ...
841 _ab_m = 1.0e-6 # ... per meter and ...
842 _ab_M = _1_0 # ... per 1,000 kilometer
843 _s_d = _3600_0 # arc-seconds per degree ...
844 _s_ab = _s_d / _ab_d # ... and per a, b unit
845 _sLat = 169028.66 # Bern, Ch in ...
846 _sLon = 26782.5 # ... arc-seconds ...
847 # lat, lon, height == 46°57'08.66", 7°26'22.50", 49.55m ("new" 46°57'07.89", 7°26'22.335")
848 Bern = LatLon4Tuple(_sLat / _s_d, _sLon / _s_d, 49.55, _WGS84, name='Bern')
850 def __init__(self, latlonh0=Bern, **other_Ltp_kwds):
851 '''New ECEF-based I{WGS84-Swiss} L{ChLV} converter, centered at I{Bern, Ch}.
853 @kwarg latlonh0: The I{geodetic} origin and height, overriding C{Bern, Ch}.
854 @kwarg other_Ltp_kwds: Optional, other L{Ltp.__init__} keyword arguments.
856 @see: L{Ltp.__init__} for more information.
857 '''
858 Ltp.__init__(self, latlonh0, **_xkwds(other_Ltp_kwds, ecef=None, name=ChLV.Bern.name))
860 def forward(self, latlonh, lon=None, height=0, M=None, name=NN): # PYCHOK unused M
861 # overloaded for the _ChLV.forward.__doc__
862 return Ltp.forward(self, latlonh, lon=lon, height=height, M=M, name=name)
864 def reverse(self, enh_, n=None, h_=0, M=None, **name): # PYCHOK signature
865 # overloaded for the _ChLV.reverse.__doc__
866 Y, X, h_, name = self._YXh_n4(enh_, n, h_, **name)
867 return Ltp.reverse(self, Y, X, h_, M=M, name=name)
869 @staticmethod
870 def false2(Y, X, LV95=True, name=NN):
871 '''Add the I{Swiss LV95} or I{LV03} falsing.
873 @arg Y: Unfalsed I{Swiss Y} easting (C{meter}).
874 @arg X: Unfalsed I{Swiss X} northing (C{meter}).
875 @kwarg LV95: If C{True} add C{LV95} falsing, if C{False} add
876 C{LV03} falsing, otherwise leave unfalsed.
877 @kwarg name: Optional name (C{str}).
879 @return: A L{ChLVEN2Tuple}C{(E_LV95, N_LV95)} or a
880 L{ChLVyx2Tuple}C{(y_LV03, x_LV03)} with falsed B{C{Y}}
881 and B{C{X}}, otherwise a L{ChLVYX2Tuple}C{(Y, X)}
882 with B{C{Y}} and B{C{X}} as-is.
883 '''
884 e, n = t = _ChLV._falsing2(LV95)
885 return t.classof(e + Y, n + X, name=name)
887 @staticmethod
888 def isLV03(e, n):
889 '''Is C{(B{e}, B{n})} a valid I{Swiss LV03} projection?
891 @arg e: Falsed (or unfalsed) I{Swiss} easting (C{meter}).
892 @arg n: Falsed (or unfalsed) I{Swiss} northing (C{meter}).
894 @return: C{True} if C{(B{e}, B{n})} is a valid, falsed I{Swiss
895 LV03}, projection C{False} otherwise.
896 '''
897 # @see: U{Map<https://www.SwissTopo.admin.CH/en/knowledge-facts/
898 # surveying-geodesy/reference-frames/local/lv95.html>}
899 return 400.0e3 < e < 900.0e3 and 40.0e3 < n < 400.0e3
901 @staticmethod
902 def isLV95(e, n, raiser=True):
903 '''Is C{(B{e}, B{n})} a valid I{Swiss LV95} or I{LV03} projection?
905 @arg e: Falsed (or unfalsed) I{Swiss} easting (C{meter}).
906 @arg n: Falsed (or unfalsed) I{Swiss} northing (C{meter}).
907 @kwarg raiser: If C{True}, throw a L{LocalError} if B{C{e}} and
908 B{C{n}} are invalid I{Swiss LV95} nor I{LV03}.
910 @return: C{True} or C{False} if C{(B{e}, B{n})} is a valid I{Swiss
911 LV95} respectively I{LV03} projection, C{None} otherwise.
912 '''
913 if ChLV.isLV03(e, n):
914 return False
915 elif ChLV.isLV03(e - 2.0e6, n - 1.0e6): # _92_falsing = _95_ - _03_
916 return True
917 elif raiser: # PYCHOK no cover
918 raise LocalError(unstr(ChLV.isLV95, e=e, n=n))
919 return None
921 @staticmethod
922 def unfalse2(e, n, LV95=None, name=NN):
923 '''Remove the I{Swiss LV95} or I{LV03} falsing.
925 @arg e: Falsed I{Swiss E_LV95} or I{y_LV03} easting (C{meter}).
926 @arg n: Falsed I{Swiss N_LV95} or I{x_LV03} northing (C{meter}).
927 @kwarg LV95: If C{True} remove I{LV95} falsing, if C{False} remove
928 I{LV03} falsing, otherwise use method C{isLV95(B{e},
929 B{n})}.
930 @kwarg name: Optional name (C{str}).
932 @return: A L{ChLVYX2Tuple}C{(Y, X)} with the unfalsed B{C{e}}
933 respectively B{C{n}}.
934 '''
935 Y, X = _ChLV._falsing2(ChLV.isLV95(e, n) if LV95 is None else LV95)
936 return ChLVYX2Tuple(e - Y, n - X, name=name)
939class ChLVa(_ChLV, LocalCartesian):
940 '''Conversion between I{WGS84 geodetic} and I{Swiss} projection coordinates
941 using the U{Approximate<https://www.SwissTopo.admin.CH/en/maps-data-online/
942 calculation-services.html>} formulas, page 13.
944 @see: Older U{references<https://GitHub.com/alphasldiallo/Swisstopo-WGS84-LV03>}.
945 '''
946 def __init__(self, name=ChLV.Bern.name):
947 '''New I{Approximate WGS84-Swiss} L{ChLVa} converter, centered at I{Bern, Ch}.
949 @kwarg name: Optional name (C{str}), overriding C{Bern.name}.
950 '''
951 LocalCartesian.__init__(self, latlonh0=ChLV.Bern, name=name)
953 def forward(self, latlonh, lon=None, height=0, M=None, name=NN):
954 # overloaded for the _ChLV.forward.__doc__
955 lat, lon, h, name = _llhn4(latlonh, lon, height, name=name)
956 a, b, h_ = _ChLV._llh2abh_3(lat, lon, h)
957 a2, b2 = a**2, b**2
959 Y = fsumf_( 72.37, 211455.93 * b,
960 -10938.51 * b * a,
961 -0.36 * b * a2,
962 -44.54 * b * b2) # + 600_000
963 X = fsumf_(147.07, 308807.95 * a,
964 3745.25 * b2,
965 76.63 * a2,
966 -194.56 * b2 * a,
967 119.79 * a2 * a) # + 200_000
968 return self._ChLV9Tuple(True, M, name, Y, X, h_, lat, lon, h)
970 def reverse(self, enh_, n=None, h_=0, M=None, **name): # PYCHOK signature
971 # overloaded for the _ChLV.reverse.__doc__
972 Y, X, h_, name = self._YXh_n4(enh_, n, h_, **name)
973 a, b, h = _ChLV._YXh_2abh3(Y, X, h_)
974 ab_d, a2, b2 = ChLV._ab_d, a**2, b**2
976 lat = Fsum(16.9023892, 3.238272 * b,
977 -0.270978 * a2,
978 -0.002528 * b2,
979 -0.0447 * a2 * b,
980 -0.014 * b2 * b).fover(ab_d)
981 lon = Fsum( 2.6779094, 4.728982 * a,
982 0.791484 * a * b,
983 0.1306 * a * b2,
984 -0.0436 * a * a2).fover(ab_d)
985 return self._ChLV9Tuple(False, M, name, Y, X, h_, lat, lon, h)
988class ChLVe(_ChLV, LocalCartesian):
989 '''Conversion between I{WGS84 geodetic} and I{Swiss} projection coordinates
990 using the U{Ellipsoidal approximate<https://www.SwissTopo.admin.CH/en/
991 maps-data-online/calculation-services.html>} formulas, pp 10-11 and U{Bolliger,
992 J.<https://eMuseum.GGGS.CH/literatur-lv/liste-Dateien/1967_Bolliger_a.pdf>}
993 pp 148-151 (also U{GGGS<https://eMuseum.GGGS.CH/literatur-lv/liste.htm>}).
995 @note: Methods L{ChLVe.forward} and L{ChLVe.reverse} have an additional keyword
996 argument C{B{gamma}=False} to approximate the I{meridian convergence}.
997 If C{B{gamma}=True} a 2-tuple C{(t, gamma)} is returned with C{t} the
998 usual result (C{ChLV9Tuple}) and C{gamma}, the I{meridian convergence}
999 (decimal C{degrees}). To convert C{gamma} to C{grades} or C{gons},
1000 use function L{pygeodesy.degrees2grades}.
1002 @see: Older U{references<https://GitHub.com/alphasldiallo/Swisstopo-WGS84-LV03>}.
1003 '''
1004 def __init__(self, name=ChLV.Bern.name):
1005 '''New I{Approximate WGS84-Swiss} L{ChLVe} converter, centered at I{Bern, Ch}.
1007 @kwarg name: Optional name (C{str}), overriding C{Bern.name}.
1008 '''
1009 LocalCartesian.__init__(self, latlonh0=ChLV.Bern, name=name)
1011 def forward(self, latlonh, lon=None, height=0, M=None, name=NN, gamma=False): # PYCHOK gamma
1012 # overloaded for the _ChLV.forward.__doc__
1013 lat, lon, h, name = _llhn4(latlonh, lon, height, name=name)
1014 a, b, h_ = _ChLV._llh2abh_3(lat, lon, h)
1015 ab_M, z, _F = ChLV._ab_M, 0, Fhorner
1017 B1 = _F(a, 211428.533991, -10939.608605, -2.658213, -8.539078, -0.00345, -0.007992)
1018 B3 = _F(a, -44.232717, 4.291740, -0.309883, 0.013924)
1019 B5 = _F(a, 0.019784, -0.004277)
1020 Y = _F(b, z, B1, z, B3, z, B5).fover(ab_M) # 1,000 Km!
1022 B0 = _F(a, z, 308770.746371, 75.028131, 120.435227, 0.009488, 0.070332, -0.00001)
1023 B2 = _F(a, 3745.408911, -193.792705, 4.340858, -0.376174, 0.004053)
1024 B4 = _F(a, -0.734684, 0.144466, -0.011842)
1025 B6 = 0.000488
1026 X = _F(b, B0, z, B2, z, B4, z, B6).fover(ab_M) # 1,000 Km!
1028 t = self._ChLV9Tuple(True, M, name, Y, X, h_, lat, lon, h)
1029 if gamma:
1030 U1 = _F(a, 2255515.207166, 2642.456961, 1.284180, 2.577486, 0.001165)
1031 U3 = _F(a, -412.991934, 64.106344, -2.679566, 0.123833)
1032 U5 = _F(a, 0.204129, -0.037725)
1033 g = _F(b, z, U1, z, U3, z, U5).fover(ChLV._ab_m) # * ChLV._ab_d degrees?
1034 t = t, g
1035 return t
1037 def reverse(self, enh_, n=None, h_=0, M=None, name=NN, gamma=False): # PYCHOK gamma
1038 # overloaded for the _ChLV.reverse.__doc__
1039 Y, X, h_, name = self._YXh_n4(enh_, n, h_, name=name)
1040 a, b, h = _ChLV._YXh_2abh3(Y, X, h_)
1041 s_d, _F, z = ChLV._s_d, Fhorner, 0
1043 A0 = _F(b, ChLV._sLat, 32386.4877666, -25.486822, -132.457771, 0.48747, 0.81305, -0.0069)
1044 A2 = _F(b, -2713.537919, -450.442705, -75.53194, -14.63049, -2.7604)
1045 A4 = _F(b, 24.42786, 13.20703, 4.7476)
1046 A6 = -0.4249
1047 lat = _F(a, A0, z, A2, z, A4, z, A6).fover(s_d)
1049 A1 = _F(b, 47297.3056722, 7925.714783, 1328.129667, 255.02202, 48.17474, 9.0243)
1050 A3 = _F(b, -442.709889, -255.02202, -96.34947, -30.0808)
1051 A5 = _F(b, 9.63495, 9.0243)
1052 lon = _F(a, ChLV._sLon, A1, z, A3, z, A5).fover(s_d)
1053 # == (ChLV._sLon + a * (A1 + a**2 * (A3 + a**2 * A5))) / s_d
1055 t = self._ChLV9Tuple(False, M, name, Y, X, h_, lat, lon, h)
1056 if gamma:
1057 U1 = _F(b, 106679.792202, 17876.57022, 4306.5241, 794.87772, 148.1545, 27.8725)
1058 U3 = _F(b, -1435.508, -794.8777, -296.309, -92.908)
1059 U5 = _F(b, 29.631, 27.873)
1060 g = _F(a, z, U1, z, U3, z, U5).fover(ChLV._s_ab) # degrees
1061 t = t, g
1062 return t
1065def tyr3d(tilt=INT0, yaw=INT0, roll=INT0, Vector=Vector3d, **Vector_kwds):
1066 '''Convert an attitude oriention into a (3-D) direction vector.
1068 @kwarg tilt: Pitch, elevation from horizontal (C{degrees}), negative down
1069 (clockwise rotation along and around the x-axis).
1070 @kwarg yaw: Bearing, heading (compass C{degrees360}), clockwise from North
1071 (counter-clockwise rotation along and around the z-axis).
1072 @kwarg roll: Roll, bank (C{degrees}), positive to the right and down
1073 (clockwise rotation along and around the y-axis).
1075 @return: A named B{C{Vector}} instance or if B{C{Vector}} is C{None},
1076 a named L{Vector3Tuple}C{(x, y, z)}.
1078 @see: U{Yaw, pitch, and roll rotations<http://MSL.CS.UIUC.edu/planning/node102.html>}
1079 and function L{pygeodesy.hartzell} argument C{los}.
1080 '''
1081 d = Attitude4Tuple(_0_0, tilt, yaw, roll).tyr3d
1082 return d if Vector is type(d) else (
1083 Vector3Tuple(d.x, d.y, d.z, name=d.name) if Vector is None else
1084 Vector(d.x, d.y, d.z, **_xkwds(Vector_kwds, name=d.name))) # PYCHOK indent
1087def _xLtp(ltp, *dflt):
1088 '''(INTERNAL) Validate B{C{ltp}}.
1089 '''
1090 if dflt and ltp is None:
1091 ltp = dflt[0]
1092 if isinstance(ltp, (LocalCartesian, Ltp)):
1093 return ltp
1094 raise _TypesError(_ltp_, ltp, Ltp, LocalCartesian)
1096# **) MIT License
1097#
1098# Copyright (C) 2016-2023 -- mrJean1 at Gmail -- All Rights Reserved.
1099#
1100# Permission is hereby granted, free of charge, to any person obtaining a
1101# copy of this software and associated documentation files (the "Software"),
1102# to deal in the Software without restriction, including without limitation
1103# the rights to use, copy, modify, merge, publish, distribute, sublicense,
1104# and/or sell copies of the Software, and to permit persons to whom the
1105# Software is furnished to do so, subject to the following conditions:
1106#
1107# The above copyright notice and this permission notice shall be included
1108# in all copies or substantial portions of the Software.
1109#
1110# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
1111# OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
1112# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
1113# THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR
1114# OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE,
1115# ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
1116# OTHER DEALINGS IN THE SOFTWARE.