Coverage for pygeodesy/utily.py: 90%
323 statements
« prev ^ index » next coverage.py v7.2.2, created at 2023-09-20 13:43 -0400
« prev ^ index » next coverage.py v7.2.2, created at 2023-09-20 13:43 -0400
2# -*- coding: utf-8 -*-
4u'''Various utility functions.
6After I{(C) Chris Veness 2011-2015} published under the same MIT Licence**, see
7U{Latitude/Longitude<https://www.Movable-Type.co.UK/scripts/latlong.html>} and
8U{Vector-based geodesy<https://www.Movable-Type.co.UK/scripts/latlong-vectors.html>}.
9'''
10# make sure int/int division yields float quotient, see .basics
11from __future__ import division as _; del _ # PYCHOK semicolon
13from pygeodesy.basics import copysign0, isint, isstr
14from pygeodesy.constants import EPS, EPS0, INF, NAN, NEG0, NINF, PI, PI2, PI_2, R_M, \
15 _float as _F, _isfinite, isnan, isnear0, isneg0, _M_KM, \
16 _M_NM, _M_SM, _0_0, _1__90, _0_5, _1_0, _N_1_0, _2__PI, \
17 _10_0, _90_0, _N_90_0, _180_0, _N_180_0, _360_0, _400_0
18from pygeodesy.errors import _ValueError, _xkwds, _xkwds_get
19from pygeodesy.interns import _edge_, _radians_, _semi_circular_, _SPACE_
20from pygeodesy.lazily import _ALL_LAZY, _ALL_MODS as _MODS
21from pygeodesy.units import Degrees, Feet, Float, Lam, Lam_, Meter, Meter2, Radians
23from math import acos, asin, atan2, cos, degrees, fabs, radians, sin, tan # pow
25__all__ = _ALL_LAZY.utily
26__version__ = '23.09.19'
28# read constant name "_M_Unit" as "meter per Unit"
29_M_CHAIN = _F( 20.1168) # yard2m(1) * 22
30_M_FATHOM = _F( 1.8288) # yard2m(1) * 2 or _M_NM * 1e-3
31_M_FOOT = _F( 0.3048) # Int'l (1 / 3.2808398950131 = 10_000 / (254 * 12))
32_M_FOOT_GE = _F( 0.31608) # German Fuss (1 / 3.1637560111364)
33_M_FOOT_FR = _F( 0.3248406) # French Pied-du-Roi or pied (1 / 3.0784329298739)
34_M_FOOT_USVY = _F( 0.3048006096012192) # US Survey (1200 / 3937)
35_M_FURLONG = _F( 201.168) # 220 * yard2m(1) == 10 * m2chain(1)
36# _M_KM = _F(1000.0) # kilo meter
37# _M_NM = _F(1852.0) # nautical mile
38# _M_SM = _F(1609.344) # statute mile
39_M_TOISE = _F( 1.9490436) # French toise, 6 pieds (6 / 3.0784329298739)
40_M_YARD_UK = _F( 0.9144) # 254 * 12 * 3 / 10_000 == 3 * ft2m(1) Int'l
43def acos1(x):
44 '''Return C{math.acos(max(-1, min(1, B{x})))}.
45 '''
46 return acos(x) if fabs(x) < _1_0 else (PI if x < 0 else _0_0)
49def acre2ha(acres):
50 '''Convert acres to hectare.
52 @arg acres: Value in acres (C{scalar}).
54 @return: Value in C{hectare} (C{float}).
56 @raise ValueError: Invalid B{C{acres}}.
57 '''
58 # 0.40468564224 == acre2m2(1) / 10_000
59 return Float(ha=Float(acres) * 0.40468564224)
62def acre2m2(acres):
63 '''Convert acres to I{square} meter.
65 @arg acres: Value in acres (C{scalar}).
67 @return: Value in C{meter^2} (C{float}).
69 @raise ValueError: Invalid B{C{acres}}.
70 '''
71 # 4046.8564224 == chain2m(1) * furlong2m(1)
72 return Meter2(Float(acres) * 4046.8564224)
75def asin1(x):
76 '''Return C{math.asin(max(-1, min(1, B{x})))}.
77 '''
78 return asin(x) if fabs(x) < _1_0 else (PI_2 if x > 0 else -PI_2) # not PI3_2!
81def atand(y_x):
82 '''Return C{atan(B{y_x})} angle in C{degrees}.
84 @see: Function L{pygeodesy.atan2d}.
85 '''
86 return atan2d(y_x, _1_0)
89def atan2b(y, x):
90 '''Return C{atan2(B{y}, B{x})} in degrees M{[0..+360]}.
92 @see: Function L{pygeodesy.atan2d}.
93 '''
94 d = atan2d(y, x)
95 if d < 0:
96 d += _360_0
97 return d
100def atan2d(y, x, reverse=False):
101 '''Return C{atan2(B{y}, B{x})} in degrees M{[-180..+180]},
102 optionally I{reversed} (by 180 degrees for C{azi2}).
104 @see: I{Karney}'s C++ function U{Math.atan2d
105 <https://GeographicLib.SourceForge.io/C++/doc/classGeographicLib_1_1Math.html>}.
106 '''
107 if fabs(y) > fabs(x) > 0:
108 if y < 0: # q = 3
109 d = degrees(atan2(x, -y)) - _90_0
110 else: # q = 2
111 d = _90_0 - degrees(atan2(x, y))
112 elif x < 0: # q = 1
113 d = copysign0(_180_0, y) - degrees(atan2(y, -x))
114 elif x > 0: # q = 0
115 d = degrees(atan2(y, x)) if y else _0_0
116 elif isnan(x) or isnan(y):
117 return NAN
118 else: # x == 0
119 d = _N_90_0 if y < 0 else (_90_0 if y > 0 else _0_0)
120 return _azireversed(d) if reverse else d
123def _azireversed(azimuth): # in .rhumbBase
124 '''(INTERNAL) Return the I{reverse} B{C{azimuth}}.
125 '''
126 return azimuth + (_N_180_0 if azimuth > 0 else _180_0)
129def chain2m(chains):
130 '''Convert I{UK} chains to meter.
132 @arg chains: Value in chains (C{scalar}).
134 @return: Value in C{meter} (C{float}).
136 @raise ValueError: Invalid B{C{chains}}.
137 '''
138 return Meter(Float(chains=chains) * _M_CHAIN)
141def circle4(earth, lat):
142 '''Get the equatorial or a parallel I{circle of latitude}.
144 @arg earth: The earth radius, ellipsoid or datum
145 (C{meter}, L{Ellipsoid}, L{Ellipsoid2},
146 L{Datum} or L{a_f2Tuple}).
147 @arg lat: Geodetic latitude (C{degrees90}, C{str}).
149 @return: A L{Circle4Tuple}C{(radius, height, lat, beta)}
150 instance.
152 @raise RangeError: Latitude B{C{lat}} outside valid range and
153 L{pygeodesy.rangerrors} set to C{True}.
155 @raise TypeError: Invalid B{C{earth}}.
157 @raise ValueError: B{C{earth}} or B{C{lat}}.
158 '''
159 E = _MODS.datums._spherical_datum(earth).ellipsoid
160 return E.circle4(lat)
163def cot(rad, **error_kwds):
164 '''Return the C{cotangent} of an angle in C{radians}.
166 @arg rad: Angle (C{radians}).
167 @kwarg error_kwds: Error to raise (C{ValueError}).
169 @return: C{cot(B{rad})}.
171 @raise ValueError: L{pygeodesy.isnear0}C{(sin(B{rad})}.
172 '''
173 s, c = sincos2(rad)
174 if isnear0(s):
175 raise _valueError(cot, rad, **error_kwds)
176 return c / s
179def cot_(*rads, **error_kwds):
180 '''Return the C{cotangent} of angle(s) in C{radiansresection}.
182 @arg rads: One or more angles (C{radians}).
183 @kwarg error_kwds: Error to raise (C{ValueError}).
185 @return: Yield the C{cot(B{rad})} for each angle.
187 @raise ValueError: See L{pygeodesy.cot}.
188 '''
189 try:
190 for r in rads:
191 yield cot(r)
192 except ValueError:
193 raise _valueError(cot_, r, **error_kwds)
196def cotd(deg, **error_kwds):
197 '''Return the C{cotangent} of an angle in C{degrees}.
199 @arg deg: Angle (C{degrees}).
200 @kwarg error_kwds: Error to raise (C{ValueError}).
202 @return: C{cot(B{deg})}.
204 @raise ValueError: L{pygeodesy.isnear0}C{(sin(B{deg})}.
205 '''
206 s, c = sincos2d(deg)
207 if isnear0(s):
208 raise _valueError(cotd, deg, **error_kwds)
209 return c / s
212def cotd_(*degs, **error_kwds):
213 '''Return the C{cotangent} of angle(s) in C{degrees}.
215 @arg degs: One or more angles (C{degrees}).
216 @kwarg error_kwds: Error to raise (C{ValueError}).
218 @return: Yield the C{cot(B{deg})} for each angle.
220 @raise ValueError: See L{pygeodesy.cotd}.
221 '''
222 try:
223 for d in degs:
224 yield cotd(d)
225 except ValueError:
226 raise _valueError(cotd_, d, **error_kwds)
229def degrees90(rad):
230 '''Convert radians to degrees and wrap M{[-270..+90]}.
232 @arg rad: Angle (C{radians}).
234 @return: Angle, wrapped (C{degrees90}).
235 '''
236 return _wrap(degrees(rad), _90_0, _360_0)
239def degrees180(rad):
240 '''Convert radians to degrees and wrap M{[-180..+180]}.
242 @arg rad: Angle (C{radians}).
244 @return: Angle, wrapped (C{degrees180}).
245 '''
246 return _wrap(degrees(rad), _180_0, _360_0)
249def degrees360(rad):
250 '''Convert radians to degrees and wrap M{[0..+360)}.
252 @arg rad: Angle (C{radians}).
254 @return: Angle, wrapped (C{degrees360}).
255 '''
256 return _wrap(degrees(rad), _360_0, _360_0)
259def degrees2grades(deg):
260 '''Convert degrees to I{grades} (aka I{gons} or I{gradians}).
262 @arg deg: Angle (C{degrees}).
264 @return: Angle (C{grades}).
265 '''
266 return Degrees(deg) * _400_0 / _360_0
269def degrees2m(deg, radius=R_M, lat=0):
270 '''Convert an angle to a distance along the equator or
271 along the parallel at an other (geodetic) latitude.
273 @arg deg: The angle (C{degrees}).
274 @kwarg radius: Mean earth radius, ellipsoid or datum
275 (C{meter}, L{Ellipsoid}, L{Ellipsoid2},
276 L{Datum} or L{a_f2Tuple}).
277 @kwarg lat: Parallel latitude (C{degrees90}, C{str}).
279 @return: Distance (C{meter}, same units as B{C{radius}}
280 or equatorial and polar radii) or C{0.0} for
281 near-polar B{C{lat}}.
283 @raise RangeError: Latitude B{C{lat}} outside valid range and
284 L{pygeodesy.rangerrors} set to C{True}.
286 @raise TypeError: Invalid B{C{radius}}.
288 @raise ValueError: Invalid B{C{deg}}, B{C{radius}} or
289 B{C{lat}}.
291 @see: Function L{radians2m} and L{m2degrees}.
292 '''
293 return _radians2m(Lam_(deg=deg, clip=0), radius, lat)
296def fathom2m(fathoms):
297 '''Convert I{Imperial} fathom to meter.
299 @arg fathoms: Value in fathoms (C{scalar}).
301 @return: Value in C{meter} (C{float}).
303 @raise ValueError: Invalid B{C{fathoms}}.
305 @see: Function L{toise2m}, U{Fathom<https://WikiPedia.org/wiki/Fathom>}
306 and U{Klafter<https://WikiPedia.org/wiki/Klafter>}.
307 '''
308 return Meter(Float(fathoms=fathoms) * _M_FATHOM)
311def ft2m(feet, usurvey=False, pied=False, fuss=False):
312 '''Convert I{International}, I{US Survey}, I{French} or I{German}
313 B{C{feet}} to C{meter}.
315 @arg feet: Value in feet (C{scalar}).
316 @kwarg usurvey: If C{True}, convert I{US Survey} foot else ...
317 @kwarg pied: If C{True}, convert French I{pied-du-Roi} else ...
318 @kwarg fuss: If C{True}, convert German I{Fuss}, otherwise
319 I{International} foot to C{meter}.
321 @return: Value in C{meter} (C{float}).
323 @raise ValueError: Invalid B{C{feet}}.
324 '''
325 return Meter(Feet(feet) * (_M_FOOT_USVY if usurvey else
326 (_M_FOOT_FR if pied else
327 (_M_FOOT_GE if fuss else _M_FOOT))))
330def furlong2m(furlongs):
331 '''Convert a furlong to meter.
333 @arg furlongs: Value in furlongs (C{scalar}).
335 @return: Value in C{meter} (C{float}).
337 @raise ValueError: Invalid B{C{furlongs}}.
338 '''
339 return Meter(Float(furlongs=furlongs) * _M_FURLONG)
342def grades(rad):
343 '''Convert radians to I{grades} (aka I{gons} or I{gradians}).
345 @arg rad: Angle (C{radians}).
347 @return: Angle (C{grades}).
348 '''
349 return Float(grades=Float(rad=rad) * _400_0 / PI2)
352def grades400(rad):
353 '''Convert radians to I{grades} (aka I{gons} or I{gradians}) and wrap M{[0..+400)}.
355 @arg rad: Angle (C{radians}).
357 @return: Angle, wrapped (C{grades}).
358 '''
359 return Float(grades400=_wrap(grades(rad), _400_0, _400_0))
362def grades2degrees(gon):
363 '''Convert I{grades} (aka I{gons} or I{gradians}) to C{degrees}.
365 @arg gon: Angle (C{grades}).
367 @return: Angle (C{degrees}).
368 '''
369 return Degrees(Float(gon=gon) * _360_0 / _400_0)
372def grades2radians(gon):
373 '''Convert I{grades} (aka I{gons} or I{gradians}) to C{radians}.
375 @arg gon: Angle (C{grades}).
377 @return: Angle (C{radians}).
378 '''
379 return Radians(Float(gon=gon) * PI2 / _400_0)
382def km2m(km):
383 '''Convert kilo meter to meter (m).
385 @arg km: Value in kilo meter (C{scalar}).
387 @return: Value in meter (C{float}).
389 @raise ValueError: Invalid B{C{km}}.
390 '''
391 return Meter(Float(km=km) * _M_KM)
394def _loneg(lon):
395 '''(INTERNAL) "Complement" of C{lon}.
396 '''
397 return _180_0 - lon
400def m2chain(meter):
401 '''Convert meter to I{UK} chains.
403 @arg meter: Value in meter (C{scalar}).
405 @return: Value in C{chains} (C{float}).
407 @raise ValueError: Invalid B{C{meter}}.
408 '''
409 return Float(chain=Meter(meter) / _M_CHAIN) # * 0.049709695378986715
412def m2degrees(distance, radius=R_M, lat=0):
413 '''Convert a distance to an angle along the equator or
414 along the parallel at an other (geodetic) latitude.
416 @arg distance: Distance (C{meter}, same units as B{C{radius}}).
417 @kwarg radius: Mean earth radius, ellipsoid or datum (C{meter},
418 an L{Ellipsoid}, L{Ellipsoid2}, L{Datum} or
419 L{a_f2Tuple}).
420 @kwarg lat: Parallel latitude (C{degrees90}, C{str}).
422 @return: Angle (C{degrees}) or C{INF} for near-polar B{C{lat}}.
424 @raise RangeError: Latitude B{C{lat}} outside valid range and
425 L{pygeodesy.rangerrors} set to C{True}.
427 @raise TypeError: Invalid B{C{radius}}.
429 @raise ValueError: Invalid B{C{distance}}, B{C{radius}}
430 or B{C{lat}}.
432 @see: Function L{m2radians} and L{degrees2m}.
433 '''
434 return degrees(m2radians(distance, radius=radius, lat=lat))
437def m2fathom(meter):
438 '''Convert meter to I{Imperial} fathoms.
440 @arg meter: Value in meter (C{scalar}).
442 @return: Value in C{fathoms} (C{float}).
444 @raise ValueError: Invalid B{C{meter}}.
446 @see: Function L{m2toise}, U{Fathom<https://WikiPedia.org/wiki/Fathom>}
447 and U{Klafter<https://WikiPedia.org/wiki/Klafter>}.
448 '''
449 return Float(fathom=Meter(meter) / _M_FATHOM) # * 0.546806649
452def m2ft(meter, usurvey=False, pied=False, fuss=False):
453 '''Convert meter to I{International}, I{US Survey}, I{French} or
454 or I{German} feet (C{ft}).
456 @arg meter: Value in meter (C{scalar}).
457 @kwarg usurvey: If C{True}, convert to I{US Survey} foot else ...
458 @kwarg pied: If C{True}, convert to French I{pied-du-Roi} else ...
459 @kwarg fuss: If C{True}, convert to German I{Fuss}, otherwise to
460 I{International} foot.
462 @return: Value in C{feet} (C{float}).
464 @raise ValueError: Invalid B{C{meter}}.
465 '''
466 # * 3.2808333333333333, US Survey 3937 / 1200
467 # * 3.2808398950131235, Int'l 10_000 / (254 * 12)
468 return Float(feet=Meter(meter) / (_M_FOOT_USVY if usurvey else
469 (_M_FOOT_FR if pied else
470 (_M_FOOT_GE if fuss else _M_FOOT))))
473def m2furlong(meter):
474 '''Convert meter to furlongs.
476 @arg meter: Value in meter (C{scalar}).
478 @return: Value in C{furlongs} (C{float}).
480 @raise ValueError: Invalid B{C{meter}}.
481 '''
482 return Float(furlong=Meter(meter) / _M_FURLONG) # * 0.00497096954
485def m2km(meter):
486 '''Convert meter to kilo meter (Km).
488 @arg meter: Value in meter (C{scalar}).
490 @return: Value in Km (C{float}).
492 @raise ValueError: Invalid B{C{meter}}.
493 '''
494 return Float(km=Meter(meter) / _M_KM)
497def m2NM(meter):
498 '''Convert meter to nautical miles (NM).
500 @arg meter: Value in meter (C{scalar}).
502 @return: Value in C{NM} (C{float}).
504 @raise ValueError: Invalid B{C{meter}}.
505 '''
506 return Float(NM=Meter(meter) / _M_NM) # * 5.39956804e-4
509def m2radians(distance, radius=R_M, lat=0):
510 '''Convert a distance to an angle along the equator or
511 along the parallel at an other (geodetic) latitude.
513 @arg distance: Distance (C{meter}, same units as B{C{radius}}).
514 @kwarg radius: Mean earth radius, ellipsoid or datum (C{meter},
515 an L{Ellipsoid}, L{Ellipsoid2}, L{Datum} or
516 L{a_f2Tuple}).
517 @kwarg lat: Parallel latitude (C{degrees90}, C{str}).
519 @return: Angle (C{radians}) or C{INF} for near-polar B{C{lat}}.
521 @raise RangeError: Latitude B{C{lat}} outside valid range and
522 L{pygeodesy.rangerrors} set to C{True}.
524 @raise TypeError: Invalid B{C{radius}}.
526 @raise ValueError: Invalid B{C{distance}}, B{C{radius}}
527 or B{C{lat}}.
529 @see: Function L{m2degrees} and L{radians2m}.
530 '''
531 m = circle4(radius, lat).radius
532 return INF if m < EPS0 else Radians(Float(distance=distance) / m)
535def m2SM(meter):
536 '''Convert meter to statute miles (SM).
538 @arg meter: Value in meter (C{scalar}).
540 @return: Value in C{SM} (C{float}).
542 @raise ValueError: Invalid B{C{meter}}.
543 '''
544 return Float(SM=Meter(meter) / _M_SM) # * 6.21369949e-4 == 1 / 1609.344
547def m2toise(meter):
548 '''Convert meter to French U{toises<https://WikiPedia.org/wiki/Toise>}.
550 @arg meter: Value in meter (C{scalar}).
552 @return: Value in C{toises} (C{float}).
554 @raise ValueError: Invalid B{C{meter}}.
556 @see: Function L{m2fathom}.
557 '''
558 return Float(toise=Meter(meter) / _M_TOISE) # * 0.513083632632119
561def m2yard(meter):
562 '''Convert meter to I{UK} yards.
564 @arg meter: Value in meter (C{scalar}).
566 @return: Value in C{yards} (C{float}).
568 @raise ValueError: Invalid B{C{meter}}.
569 '''
570 return Float(yard=Meter(meter) / _M_YARD_UK) # * 1.0936132983377078
573def NM2m(nm):
574 '''Convert nautical miles to meter (m).
576 @arg nm: Value in nautical miles (C{scalar}).
578 @return: Value in meter (C{float}).
580 @raise ValueError: Invalid B{C{nm}}.
581 '''
582 return Meter(Float(nm=nm) * _M_NM)
585def _passarg(arg): # in .auxilats.auxLat, .formy
586 '''(INTERNAL) Helper, no-op.
587 '''
588 return arg
591def _passargs(*args): # in .formy
592 '''(INTERNAL) Helper, no-op.
593 '''
594 return args
597def radians2m(rad, radius=R_M, lat=0):
598 '''Convert an angle to a distance along the equator or
599 along the parallel at an other (geodetic) latitude.
601 @arg rad: The angle (C{radians}).
602 @kwarg radius: Mean earth radius, ellipsoid or datum
603 (C{meter}, L{Ellipsoid}, L{Ellipsoid2},
604 L{Datum} or L{a_f2Tuple}).
605 @kwarg lat: Parallel latitude (C{degrees90}, C{str}).
607 @return: Distance (C{meter}, same units as B{C{radius}}
608 or equatorial and polar radii) or C{0.0} for
609 near-polar B{C{lat}}.
611 @raise RangeError: Latitude B{C{lat}} outside valid range and
612 L{pygeodesy.rangerrors} set to C{True}.
614 @raise TypeError: Invalid B{C{radius}}.
616 @raise ValueError: Invalid B{C{rad}}, B{C{radius}} or
617 B{C{lat}}.
619 @see: Function L{degrees2m} and L{m2radians}.
620 '''
621 return _radians2m(Lam(rad=rad, clip=0), radius, lat)
624def _radians2m(rad, radius, lat):
625 '''(INTERNAL) Helper for C{degrees2m} and C{radians2m}.
626 '''
627 m = circle4(radius, lat).radius
628 return _0_0 if m < EPS0 else (rad * m)
631def radiansPI(deg):
632 '''Convert and wrap degrees to radians M{[-PI..+PI]}.
634 @arg deg: Angle (C{degrees}).
636 @return: Radians, wrapped (C{radiansPI})
637 '''
638 return _wrap(radians(deg), PI, PI2)
641def radiansPI2(deg):
642 '''Convert and wrap degrees to radians M{[0..+2PI)}.
644 @arg deg: Angle (C{degrees}).
646 @return: Radians, wrapped (C{radiansPI2})
647 '''
648 return _wrap(radians(deg), PI2, PI2)
651def radiansPI_2(deg):
652 '''Convert and wrap degrees to radians M{[-3PI/2..+PI/2]}.
654 @arg deg: Angle (C{degrees}).
656 @return: Radians, wrapped (C{radiansPI_2})
657 '''
658 return _wrap(radians(deg), PI_2, PI2)
661def _sin0cos2(q, r, sign):
662 '''(INTERNAL) 2-tuple (C{sin(r), cos(r)}) in quadrant C{0 <= B{q} <= 3}
663 and C{sin} zero I{signed} with B{C{sign}}.
664 '''
665 if r < PI_2:
666 s, c = sin(r), cos(r)
667 t = s, c, -s, -c, s
668 else: # r == PI_2
669 t = _1_0, _0_0, _N_1_0, _0_0, _1_0
670# else: # r == 0, testUtility failures
671# t = _0_0, _1_0, _0_0, _N_1_0, _0_0
672# q &= 3
673 s = t[q] or (NEG0 if sign < 0 else _0_0)
674 c = t[q + 1] or _0_0
675 return s, c
678def sincos2(rad):
679 '''Return the C{sine} and C{cosine} of an angle in C{radians}.
681 @arg rad: Angle (C{radians}).
683 @return: 2-Tuple (C{sin(B{rad})}, C{cos(B{rad})}).
685 @see: U{GeographicLib<https://GeographicLib.SourceForge.io/C++/doc/
686 classGeographicLib_1_1Math.html#sincosd>} function U{sincosd
687 <https://SourceForge.net/p/geographiclib/code/ci/release/tree/
688 python/geographiclib/geomath.py#l155>} and C++ U{sincosd
689 <https://SourceForge.net/p/geographiclib/code/ci/release/tree/
690 include/GeographicLib/Math.hpp#l558>}.
691 '''
692 if _isfinite(rad):
693 q = int(rad * _2__PI) # int(math.floor)
694 if q < 0:
695 q -= 1
696 t = _sin0cos2(q & 3, rad - q * PI_2, rad)
697 else:
698 t = NAN, NAN
699 return t
702def sincos2_(*rads):
703 '''Return the C{sine} and C{cosine} of angle(s) in C{radians}.
705 @arg rads: One or more angles (C{radians}).
707 @return: Yield the C{sin(B{rad})} and C{cos(B{rad})} for each angle.
709 @see: function L{sincos2}.
710 '''
711 for r in rads:
712 s, c = sincos2(r)
713 yield s
714 yield c
717def sincos2d(deg, **adeg):
718 '''Return the C{sine} and C{cosine} of an angle in C{degrees}.
720 @arg deg: Angle (C{degrees}).
721 @kwarg adeg: Optional correction (C{degrees}).
723 @return: 2-Tuple (C{sin(B{deg_})}, C{cos(B{deg_})}, C{B{deg_} =
724 B{deg} + B{adeg}}).
726 @see: U{GeographicLib<https://GeographicLib.SourceForge.io/C++/doc/
727 classGeographicLib_1_1Math.html#sincosd>} function U{sincosd
728 <https://SourceForge.net/p/geographiclib/code/ci/release/tree/
729 python/geographiclib/geomath.py#l155>} and C++ U{sincosd
730 <https://SourceForge.net/p/geographiclib/code/ci/release/tree/
731 include/GeographicLib/Math.hpp#l558>}.
732 '''
733 if _isfinite(deg):
734 q = int(deg * _1__90) # int(math.floor)
735 if q < 0:
736 q -= 1
737 d = deg - q * _90_0
738 if adeg:
739 t = _xkwds_get(adeg, adeg=_0_0)
740 d = _MODS.karney._around(d + t)
741 t = _sin0cos2(q & 3, radians(d), deg)
742 else:
743 t = NAN, NAN
744 return t
747def SinCos2(x):
748 '''Get C{sin} and C{cos} of I{typed} angle.
750 @arg x: Angle (L{Degrees}, L{Radians} or C{radians}).
752 @return: 2-Tuple (C{sin(B{x})}, C{cos(B{x})}).
753 '''
754 return sincos2d(x) if isinstance(x, Degrees) else (
755 sincos2(x) if isinstance(x, Radians) else
756 sincos2(float(x))) # assume C{radians}
759def sincos2d_(*degs):
760 '''Return the C{sine} and C{cosine} of angle(s) in C{degrees}.
762 @arg degs: One or more angles (C{degrees}).
764 @return: Yield the C{sin(B{deg})} and C{cos(B{deg})} for each angle.
766 @see: Function L{sincos2d}.
767 '''
768 for d in degs:
769 s, c = sincos2d(d)
770 yield s
771 yield c
774def sincostan3(rad):
775 '''Return the C{sine}, C{cosine} and C{tangent} of an angle in C{radians}.
777 @arg rad: Angle (C{radians}).
779 @return: 3-Tuple (C{sin(B{rad})}, C{cos(B{rad})}, C{tan(B{rad})}).
781 @see: Function L{sincos2}.
782 '''
783 rad %= PI2 # see ._wrap comments
784 if rad:
785 s, c = sincos2(rad)
786 t = (s / c) if c else (NINF if s < 0
787 else (INF if s > 0 else s))
788 else: # sin(-0.0) == tan(-0.0) = -0.0
789 c = _1_0
790 s = t = NEG0 if isneg0(rad) else _0_0
791 return s, c, t
794def SM2m(sm):
795 '''Convert statute miles to meter (m).
797 @arg sm: Value in statute miles (C{scalar}).
799 @return: Value in meter (C{float}).
801 @raise ValueError: Invalid B{C{sm}}.
802 '''
803 return Meter(Float(sm=sm) * _M_SM)
806def tan_2(rad, **semi): # edge=1
807 '''Compute the tangent of half angle.
809 @arg rad: Angle (C{radians}).
810 @kwarg semi: Angle or edge name and index
811 for semi-circular error.
813 @return: M{tan(rad / 2)} (C{float}).
815 @raise ValueError: If B{C{rad}} is semi-circular
816 and B{C{semi}} is given.
817 '''
818 # .formy.excessKarney_, .sphericalTrigonometry.areaOf
819 if semi and isnear0(fabs(rad) - PI):
820 for n, v in semi.items():
821 break
822 n = _SPACE_(n, _radians_) if not isint(v) else \
823 _SPACE_(_MODS.streprs.Fmt.SQUARE(**semi), _edge_)
824 raise _ValueError(n, rad, txt=_semi_circular_)
826 return tan(rad * _0_5)
829def tand(deg, **error_kwds):
830 '''Return the C{tangent} of an angle in C{degrees}.
832 @arg deg: Angle (C{degrees}).
833 @kwarg error_kwds: Error to raise (C{ValueError}).
835 @return: C{tan(B{deg})}.
837 @raise ValueError: If L{pygeodesy.isnear0}C{(cos(B{deg})}.
838 '''
839 s, c = sincos2d(deg)
840 if isnear0(c):
841 raise _valueError(tand, deg, **error_kwds)
842 return s / c
845def tand_(*degs, **error_kwds):
846 '''Return the C{tangent} of angle(s) in C{degrees}.
848 @arg degs: One or more angles (C{degrees}).
849 @kwarg error_kwds: Error to raise (C{ValueError}).
851 @return: Yield the C{tan(B{deg})} for each angle.
853 @raise ValueError: See L{pygeodesy.tand}.
854 '''
855 for d in degs:
856 yield tand(d, **error_kwds)
859def tanPI_2_2(rad):
860 '''Compute the tangent of half angle, 90 degrees rotated.
862 @arg rad: Angle (C{radians}).
864 @return: M{tan((rad + PI/2) / 2)} (C{float}).
865 '''
866 return tan((rad + PI_2) * _0_5)
869def toise2m(toises):
870 '''Convert French U{toises<https://WikiPedia.org/wiki/Toise>} to meter.
872 @arg toises: Value in toises (C{scalar}).
874 @return: Value in C{meter} (C{float}).
876 @raise ValueError: Invalid B{C{toises}}.
878 @see: Function L{fathom2m}.
879 '''
880 return Meter(Float(toises=toises) * _M_TOISE)
883def truncate(x, ndigits=None):
884 '''Truncate to the given number of digits.
886 @arg x: Value to truncate (C{scalar}).
887 @kwarg ndigits: Number of digits (C{int}),
888 aka I{precision}.
890 @return: Truncated B{C{x}} (C{float}).
892 @see: Python function C{round}.
893 '''
894 if isint(ndigits):
895 p = _10_0**ndigits
896 x = int(x * p) / p
897 return x
900def unroll180(lon1, lon2, wrap=True):
901 '''Unroll longitudinal delta and wrap longitude in degrees.
903 @arg lon1: Start longitude (C{degrees}).
904 @arg lon2: End longitude (C{degrees}).
905 @kwarg wrap: If C{True}, wrap and unroll to the M{(-180..+180]}
906 range (C{bool}).
908 @return: 2-Tuple C{(B{lon2}-B{lon1}, B{lon2})} unrolled (C{degrees},
909 C{degrees}).
911 @see: Capability C{LONG_UNROLL} in U{GeographicLib
912 <https://GeographicLib.SourceForge.io/html/python/interface.html#outmask>}.
913 '''
914 d = lon2 - lon1
915 if wrap and fabs(d) > _180_0:
916 u = _wrap(d, _180_0, _360_0)
917 if u != d:
918 return u, (lon1 + u)
919 return d, lon2
922def _unrollon(p1, p2, wrap=False): # unroll180 == .karney._unroll2
923 '''(INTERNAL) Wrap/normalize, unroll and replace longitude.
924 '''
925 lat, lon = p2.lat, p2.lon
926 if wrap and _Wrap.normal:
927 lat, lon = _Wrap.latlon(lat, lon)
928 _, lon = unroll180(p1.lon, lon, wrap=True)
929 if lat != p2.lat or fabs(lon - p2.lon) > EPS:
930 p2 = p2.dup(lat=lat, lon=wrap180(lon))
931 # p2 = p2.copy(); p2.latlon = lat, wrap180(lon)
932 return p2
935def _unrollon3(p1, p2, p3, wrap=False):
936 '''(INTERNAL) Wrap/normalize, unroll 2 points.
937 '''
938 w = wrap
939 if w:
940 w = _Wrap.normal
941 p2 = _unrollon(p1, p2, wrap=w)
942 p3 = _unrollon(p1, p3, wrap=w)
943 p2 = _unrollon(p2, p3)
944 return p2, p3, w # was wrapped?
947def unrollPI(rad1, rad2, wrap=True):
948 '''Unroll longitudinal delta and wrap longitude in radians.
950 @arg rad1: Start longitude (C{radians}).
951 @arg rad2: End longitude (C{radians}).
952 @kwarg wrap: If C{True}, wrap and unroll to the M{(-PI..+PI]}
953 range (C{bool}).
955 @return: 2-Tuple C{(B{rad2}-B{rad1}, B{rad2})} unrolled
956 (C{radians}, C{radians}).
958 @see: Capability C{LONG_UNROLL} in U{GeographicLib
959 <https://GeographicLib.SourceForge.io/html/python/interface.html#outmask>}.
960 '''
961 r = rad2 - rad1
962 if wrap and fabs(r) > PI:
963 u = _wrap(r, PI, PI2)
964 if u != r:
965 return u, (rad1 + u)
966 return r, rad2
969def _valueError(where, x, **kwds):
970 '''(INTERNAL) Return a C{_ValueError}.
971 '''
972 x = _MODS.streprs.Fmt.PAREN(where.__name__, x)
973 return _ValueError(x, **kwds)
976class _Wrap(object):
978 _normal = False # default
980 @property
981 def normal(self):
982 '''Get the current L{normal} setting (C{True},
983 C{False} or C{None}).
984 '''
985 return self._normal
987 @normal.setter # PYCHOK setter!
988 def normal(self, setting):
989 '''Set L{normal} to C{True}, C{False} or C{None}.
990 '''
991 t = {True: (_MODS.formy.normal, _MODS.formy.normal_),
992 False: (self.wraplatlon, self.wraphilam),
993 None: (_passargs, _passargs)}.get(setting, ())
994 if t:
995 self.latlon, self.philam = t
996 self._normal = setting
998 def latlonDMS2(self, lat, lon, **DMS2_kwds):
999 if isstr(lat) or isstr(lon):
1000 kwds = _xkwds(DMS2_kwds, clipLon=0, clipLat=0)
1001 lat, lon = _MODS.dms.parseDMS2(lat, lon, **kwds)
1002 return self.latlon(lat, lon)
1004# def normalatlon(self, *latlon):
1005# return _MODS.formy.normal(*latlon)
1007# def normalamphi(self, *philam):
1008# return _MODS.formy.normal_(*philam)
1010 def wraplatlon(self, lat, lon):
1011 return wrap90(lat), wrap180(lon)
1013 latlon = wraplatlon # default
1015 def latlon3(self, lon1, lat2, lon2, wrap):
1016 if wrap:
1017 lat2, lon2 = self.latlon(lat2, lon2)
1018 lon21, lon2 = unroll180(lon1, lon2)
1019 else:
1020 lon21 = lon2 - lon1
1021 return lon21, lat2, lon2
1023 def _latlonop(self, wrap):
1024 if wrap and self._normal is not None:
1025 return self.latlon
1026 else:
1027 return _passargs
1029 def wraphilam(self, phi, lam):
1030 return wrapPI_2(phi), wrapPI(lam)
1032 philam = wraphilam # default
1034 def philam3(self, lam1, phi2, lam2, wrap):
1035 if wrap:
1036 phi2, lam2 = self.philam(phi2, lam2)
1037 lam21, lam2 = unrollPI(lam1, lam2)
1038 else:
1039 lam21 = lam2 - lam1
1040 return lam21, phi2, lam2
1042 def _philamop(self, wrap):
1043 if wrap and self._normal is not None:
1044 return self.philam
1045 else:
1046 return _passargs
1048 def wraphilam(self, phi, lam):
1049 return wrapPI_2(phi), wrapPI(lam)
1051 def point(self, ll, wrap=True): # in .points._fractional, -.PointsIter.iterate, ...
1052 '''Return C{ll} or a copy, I{normalized} or I{wrap}'d.
1053 '''
1054 if wrap and self._normal is not None:
1055 lat, lon = ll.latlon
1056 if fabs(lon) > 180 or fabs(lat) > 90:
1057 _n = self.latlon
1058 ll = ll.copy(name=_n.__name__)
1059 ll.latlon = _n(lat, lon)
1060 return ll
1062_Wrap = _Wrap() # PYCHOK singleton
1065def _wrap(angle, wrap, modulo):
1066 '''(INTERNAL) Angle wrapper M{((wrap-modulo)..+wrap]}.
1068 @arg angle: Angle (C{degrees}, C{radians} or C{grades}).
1069 @arg wrap: Range (C{degrees}, C{radians} or C{grades}).
1070 @arg modulo: Upper limit (360 C{degrees}, PI2 C{radians} or 400 C{grades}).
1072 @return: The B{C{angle}}, wrapped (C{degrees}, C{radians} or C{grades}).
1073 '''
1074 a = float(angle)
1075 if not (wrap - modulo) <= a < wrap:
1076 # math.fmod(-1.5, 3.14) == -1.5, but -1.5 % 3.14 == 1.64
1077 # math.fmod(-1.5, 360) == -1.5, but -1.5 % 360 == 358.5
1078 a %= modulo
1079 if a > wrap:
1080 a -= modulo
1081 return a
1084def wrap90(deg):
1085 '''Wrap degrees to M{[-270..+90]}.
1087 @arg deg: Angle (C{degrees}).
1089 @return: Degrees, wrapped (C{degrees90}).
1090 '''
1091 return _wrap(deg, _90_0, _360_0)
1094def wrap180(deg):
1095 '''Wrap degrees to M{[-180..+180]}.
1097 @arg deg: Angle (C{degrees}).
1099 @return: Degrees, wrapped (C{degrees180}).
1100 '''
1101 return _wrap(deg, _180_0, _360_0)
1104def wrap360(deg): # see .streprs._umod_360
1105 '''Wrap degrees to M{[0..+360)}.
1107 @arg deg: Angle (C{degrees}).
1109 @return: Degrees, wrapped (C{degrees360}).
1110 '''
1111 return _wrap(deg, _360_0, _360_0)
1114def wrapPI(rad):
1115 '''Wrap radians to M{[-PI..+PI]}.
1117 @arg rad: Angle (C{radians}).
1119 @return: Radians, wrapped (C{radiansPI}).
1120 '''
1121 return _wrap(rad, PI, PI2)
1124def wrapPI2(rad):
1125 '''Wrap radians to M{[0..+2PI)}.
1127 @arg rad: Angle (C{radians}).
1129 @return: Radians, wrapped (C{radiansPI2}).
1130 '''
1131 return _wrap(rad, PI2, PI2)
1134def wrapPI_2(rad):
1135 '''Wrap radians to M{[-3PI/2..+PI/2]}.
1137 @arg rad: Angle (C{radians}).
1139 @return: Radians, wrapped (C{radiansPI_2}).
1140 '''
1141 return _wrap(rad, PI_2, PI2)
1144def wrap_normal(*normal):
1145 '''Define the operation for the keyword argument C{B{wrap}=True},
1146 across L{pygeodesy}: I{wrap}, I{normalize} or I{no-op}. For
1147 backward compatibility, the default is I{wrap}.
1149 @arg normal: If C{True}, I{normalize} lat- and longitude using
1150 L{normal} or L{normal_}, if C{False}, I{wrap} the
1151 lat- and longitude individually by L{wrap90} or
1152 L{wrapPI_2} respectively L{wrap180}, L{wrapPI} or
1153 if C{None}, leave lat- and longitude I{unchanged}.
1154 Do not supply any value to get the current setting.
1156 @return: The previous L{wrap_normal} setting (C{bool} or C{None}).
1157 '''
1158 t = _Wrap.normal
1159 if normal:
1160 _Wrap.normal = normal[0]
1161 return t
1164def yard2m(yards):
1165 '''Convert I{UK} yards to meter.
1167 @arg yards: Value in yards (C{scalar}).
1169 @return: Value in C{meter} (C{float}).
1171 @raise ValueError: Invalid B{C{yards}}.
1172 '''
1173 return Float(yards=yards) * _M_YARD_UK
1175# **) MIT License
1176#
1177# Copyright (C) 2016-2023 -- mrJean1 at Gmail -- All Rights Reserved.
1178#
1179# Permission is hereby granted, free of charge, to any person obtaining a
1180# copy of this software and associated documentation files (the "Software"),
1181# to deal in the Software without restriction, including without limitation
1182# the rights to use, copy, modify, merge, publish, distribute, sublicense,
1183# and/or sell copies of the Software, and to permit persons to whom the
1184# Software is furnished to do so, subject to the following conditions:
1185#
1186# The above copyright notice and this permission notice shall be included
1187# in all copies or substantial portions of the Software.
1188#
1189# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
1190# OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
1191# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
1192# THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR
1193# OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE,
1194# ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
1195# OTHER DEALINGS IN THE SOFTWARE.