Coverage for pygeodesy/utily.py: 91%
335 statements
« prev ^ index » next coverage.py v7.2.2, created at 2024-01-26 16:28 -0500
« prev ^ index » next coverage.py v7.2.2, created at 2024-01-26 16:28 -0500
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 _copysign, isinstanceof, isint, isstr, neg, _passargs
14from pygeodesy.constants import EPS, EPS0, INF, NAN, NEG0, PI, PI2, PI_2, R_M, \
15 _float as _F, _isfinite, isnan, isnear0, _over, \
16 _umod_360, _umod_PI2, _M_KM, _M_NM, _M_SM, _0_0, \
17 _1__90, _0_5, _1_0, _N_1_0, _2__PI, _10_0, _90_0, \
18 _N_90_0, _180_0, _N_180_0, _360_0, _400_0
19from pygeodesy.errors import _ValueError, _xkwds, _xkwds_get, _ALL_LAZY, _MODS
20from pygeodesy.interns import _edge_, _radians_, _semi_circular_, _SPACE_
21# from pygeodesy.lazily import _ALL_LAZY, _ALL_MODS as _MODS # from .errors
22from pygeodesy.units import Degrees, Degrees_, Feet, Float, Lam, Lam_, \
23 Meter, Meter2, Radians, Radians_
25from math import acos, asin, atan2, cos, degrees, fabs, radians, sin, tan # pow
27__all__ = _ALL_LAZY.utily
28__version__ = '24.01.12'
30# read constant name "_M_Unit" as "meter per Unit"
31_M_CHAIN = _F( 20.1168) # yard2m(1) * 22
32_M_FATHOM = _F( 1.8288) # yard2m(1) * 2 or _M_NM * 1e-3
33_M_FOOT = _F( 0.3048) # Int'l (1 / 3.2808398950131 = 10_000 / (254 * 12))
34_M_FOOT_GE = _F( 0.31608) # German Fuss (1 / 3.1637560111364)
35_M_FOOT_FR = _F( 0.3248406) # French Pied-du-Roi or pied (1 / 3.0784329298739)
36_M_FOOT_USVY = _F( 0.3048006096012192) # US Survey (1200 / 3937)
37_M_FURLONG = _F( 201.168) # 220 * yard2m(1) == 10 * m2chain(1)
38# _M_KM = _F(1000.0) # kilo meter
39# _M_NM = _F(1852.0) # nautical mile
40# _M_SM = _F(1609.344) # statute mile
41_M_TOISE = _F( 1.9490436) # French toise, 6 pieds (6 / 3.0784329298739)
42_M_YARD_UK = _F( 0.9144) # 254 * 12 * 3 / 10_000 == 3 * ft2m(1) Int'l
45def _abs1nan(x):
46 '''(INTERNAL) Bracket C{x}.
47 '''
48 return _N_1_0 < x < _1_0 or isnan(x)
51def acos1(x):
52 '''Return C{math.acos(max(-1, min(1, B{x})))}.
53 '''
54 return acos(x) if _abs1nan(x) else (PI if x < 0 else _0_0)
57def acre2ha(acres):
58 '''Convert acres to hectare.
60 @arg acres: Value in acres (C{scalar}).
62 @return: Value in C{hectare} (C{float}).
64 @raise ValueError: Invalid B{C{acres}}.
65 '''
66 # 0.40468564224 == acre2m2(1) / 10_000
67 return Float(ha=Float(acres) * 0.40468564224)
70def acre2m2(acres):
71 '''Convert acres to I{square} meter.
73 @arg acres: Value in acres (C{scalar}).
75 @return: Value in C{meter^2} (C{float}).
77 @raise ValueError: Invalid B{C{acres}}.
78 '''
79 # 4046.8564224 == chain2m(1) * furlong2m(1)
80 return Meter2(Float(acres) * 4046.8564224)
83def asin1(x):
84 '''Return C{math.asin(max(-1, min(1, B{x})))}.
85 '''
86 return asin(x) if _abs1nan(x) else _copysign(PI_2, x)
89def atan1(y, x=_1_0):
90 '''Return C{atan(B{y} / B{x})} angle in C{radians} M{[-PI/2..+PI/2]}
91 using C{atan2} for consistency and to avoid C{ZeroDivisionError}.
92 '''
93 return atan2(-y, -x) if x < 0 else atan2(y, x or _0_0) # -0. to 0.
96def atan1d(y, x=_1_0):
97 '''Return C{atan(B{y} / B{x})} angle in C{degrees} M{[-90..+90]}
98 using C{atan2d} for consistency and to avoid C{ZeroDivisionError}.
100 @see: Function L{pygeodesy.atan2d}.
101 '''
102 return atan2d(-y, -x) if x < 0 else atan2d(y, x or _0_0) # -0. to 0.
105def atan2b(y, x):
106 '''Return C{atan2(B{y}, B{x})} in degrees M{[0..+360]}.
108 @see: Function L{pygeodesy.atan2d}.
109 '''
110 b = atan2d(y, x)
111 if b < 0:
112 b += _360_0
113 return b
116def atan2d(y, x, reverse=False):
117 '''Return C{atan2(B{y}, B{x})} in degrees M{[-180..+180]},
118 optionally I{reversed} (by 180 degrees for C{azimuth}s).
120 @see: I{Karney}'s C++ function U{Math.atan2d
121 <https://GeographicLib.SourceForge.io/C++/doc/classGeographicLib_1_1Math.html>}.
122 '''
123 if fabs(y) > fabs(x) > 0:
124 if y < 0: # q = 3
125 d = degrees(atan2(x, -y)) - _90_0
126 else: # q = 2
127 d = _90_0 - degrees(atan2(x, y))
128 elif isnan(x) or isnan(y):
129 return NAN
130 elif y:
131 if x > 0: # q = 0
132 d = degrees(atan2(y, x))
133 elif x < 0: # q = 1
134 d = _copysign(_180_0, y) - degrees(atan2(y, -x))
135 else: # x == 0
136 d = _copysign(_90_0, y)
137 else:
138 d = _180_0 if x < 0 else _0_0
139 return _azireversed(d) if reverse else d
142def _azireversed(azimuth): # in .rhumbBase
143 '''(INTERNAL) Return the I{reverse} B{C{azimuth}} in degrees M{[-180..+180]}.
144 '''
145 return azimuth + (_N_180_0 if azimuth > 0 else _180_0)
148def chain2m(chains):
149 '''Convert I{UK} chains to meter.
151 @arg chains: Value in chains (C{scalar}).
153 @return: Value in C{meter} (C{float}).
155 @raise ValueError: Invalid B{C{chains}}.
156 '''
157 return Meter(Float(chains=chains) * _M_CHAIN)
160def circle4(earth, lat):
161 '''Get the equatorial or a parallel I{circle of latitude}.
163 @arg earth: The earth radius, ellipsoid or datum
164 (C{meter}, L{Ellipsoid}, L{Ellipsoid2},
165 L{Datum} or L{a_f2Tuple}).
166 @arg lat: Geodetic latitude (C{degrees90}, C{str}).
168 @return: A L{Circle4Tuple}C{(radius, height, lat, beta)}
169 instance.
171 @raise RangeError: Latitude B{C{lat}} outside valid range and
172 L{pygeodesy.rangerrors} set to C{True}.
174 @raise TypeError: Invalid B{C{earth}}.
176 @raise ValueError: B{C{earth}} or B{C{lat}}.
177 '''
178 E = _MODS.datums._earth_ellipsoid(earth)
179 return E.circle4(lat)
182def cot(rad, **error_kwds):
183 '''Return the C{cotangent} of an angle in C{radians}.
185 @arg rad: Angle (C{radians}).
186 @kwarg error_kwds: Error to raise (C{ValueError}).
188 @return: C{cot(B{rad})}.
190 @raise ValueError: L{pygeodesy.isnear0}C{(sin(B{rad})}.
191 '''
192 s, c = sincos2(rad)
193 if c:
194 if isnear0(s):
195 raise _valueError(cot, rad, **error_kwds)
196 c = c / s # /= chokes PyChecker
197 return c
200def cot_(*rads, **error_kwds):
201 '''Return the C{cotangent} of angle(s) in C{radiansresection}.
203 @arg rads: One or more angles (C{radians}).
204 @kwarg error_kwds: Error to raise (C{ValueError}).
206 @return: Yield the C{cot(B{rad})} for each angle.
208 @raise ValueError: See L{pygeodesy.cot}.
209 '''
210 try:
211 for r in rads:
212 yield cot(r)
213 except ValueError:
214 raise _valueError(cot_, r, **error_kwds)
217def cotd(deg, **error_kwds):
218 '''Return the C{cotangent} of an angle in C{degrees}.
220 @arg deg: Angle (C{degrees}).
221 @kwarg error_kwds: Error to raise (C{ValueError}).
223 @return: C{cot(B{deg})}.
225 @raise ValueError: L{pygeodesy.isnear0}C{(sin(B{deg})}.
226 '''
227 s, c = sincos2d(deg)
228 if c:
229 if isnear0(s):
230 raise _valueError(cotd, deg, **error_kwds)
231 c = c / s # /= chokes PyChecker
232 elif s < 0:
233 c = -c # negate-0
234 return c
237def cotd_(*degs, **error_kwds):
238 '''Return the C{cotangent} of angle(s) in C{degrees}.
240 @arg degs: One or more angles (C{degrees}).
241 @kwarg error_kwds: Error to raise (C{ValueError}).
243 @return: Yield the C{cot(B{deg})} for each angle.
245 @raise ValueError: See L{pygeodesy.cotd}.
246 '''
247 try:
248 for d in degs:
249 yield cotd(d)
250 except ValueError:
251 raise _valueError(cotd_, d, **error_kwds)
254def degrees90(rad):
255 '''Convert radians to degrees and wrap M{[-90..+90)}.
257 @arg rad: Angle (C{radians}).
259 @return: Angle, wrapped (C{degrees90}).
260 '''
261 return wrap90(degrees(rad))
264def degrees180(rad):
265 '''Convert radians to degrees and wrap M{[-180..+180)}.
267 @arg rad: Angle (C{radians}).
269 @return: Angle, wrapped (C{degrees180}).
270 '''
271 return wrap180(degrees(rad))
274def degrees360(rad):
275 '''Convert radians to degrees and wrap M{[0..+360)}.
277 @arg rad: Angle (C{radians}).
279 @return: Angle, wrapped (C{degrees360}).
280 '''
281 return _umod_360(degrees(rad))
284def degrees2grades(deg):
285 '''Convert degrees to I{grades} (aka I{gons} or I{gradians}).
287 @arg deg: Angle (C{degrees}).
289 @return: Angle (C{grades}).
290 '''
291 return Degrees(deg) * _400_0 / _360_0
294def degrees2m(deg, radius=R_M, lat=0):
295 '''Convert an angle to a distance along the equator or
296 along the parallel at an other (geodetic) latitude.
298 @arg deg: The angle (C{degrees}).
299 @kwarg radius: Mean earth radius, ellipsoid or datum
300 (C{meter}, L{Ellipsoid}, L{Ellipsoid2},
301 L{Datum} or L{a_f2Tuple}).
302 @kwarg lat: Parallel latitude (C{degrees90}, C{str}).
304 @return: Distance (C{meter}, same units as B{C{radius}}
305 or equatorial and polar radii) or C{0.0} for
306 near-polar B{C{lat}}.
308 @raise RangeError: Latitude B{C{lat}} outside valid range and
309 L{pygeodesy.rangerrors} set to C{True}.
311 @raise TypeError: Invalid B{C{radius}}.
313 @raise ValueError: Invalid B{C{deg}}, B{C{radius}} or
314 B{C{lat}}.
316 @see: Function L{radians2m} and L{m2degrees}.
317 '''
318 return _Radians2m(Lam_(deg=deg, clip=0), radius, lat)
321def fathom2m(fathoms):
322 '''Convert I{Imperial} fathom to meter.
324 @arg fathoms: Value in fathoms (C{scalar}).
326 @return: Value in C{meter} (C{float}).
328 @raise ValueError: Invalid B{C{fathoms}}.
330 @see: Function L{toise2m}, U{Fathom<https://WikiPedia.org/wiki/Fathom>}
331 and U{Klafter<https://WikiPedia.org/wiki/Klafter>}.
332 '''
333 return Meter(Float(fathoms=fathoms) * _M_FATHOM)
336def ft2m(feet, usurvey=False, pied=False, fuss=False):
337 '''Convert I{International}, I{US Survey}, I{French} or I{German}
338 B{C{feet}} to C{meter}.
340 @arg feet: Value in feet (C{scalar}).
341 @kwarg usurvey: If C{True}, convert I{US Survey} foot else ...
342 @kwarg pied: If C{True}, convert French I{pied-du-Roi} else ...
343 @kwarg fuss: If C{True}, convert German I{Fuss}, otherwise
344 I{International} foot to C{meter}.
346 @return: Value in C{meter} (C{float}).
348 @raise ValueError: Invalid B{C{feet}}.
349 '''
350 return Meter(Feet(feet) * (_M_FOOT_USVY if usurvey else
351 (_M_FOOT_FR if pied else
352 (_M_FOOT_GE if fuss else _M_FOOT))))
355def furlong2m(furlongs):
356 '''Convert a furlong to meter.
358 @arg furlongs: Value in furlongs (C{scalar}).
360 @return: Value in C{meter} (C{float}).
362 @raise ValueError: Invalid B{C{furlongs}}.
363 '''
364 return Meter(Float(furlongs=furlongs) * _M_FURLONG)
367def grades(rad):
368 '''Convert radians to I{grades} (aka I{gons} or I{gradians}).
370 @arg rad: Angle (C{radians}).
372 @return: Angle (C{grades}).
373 '''
374 return Float(grades=Float(rad=rad) * _400_0 / PI2)
377def grades400(rad):
378 '''Convert radians to I{grades} (aka I{gons} or I{gradians}) and wrap M{[0..+400)}.
380 @arg rad: Angle (C{radians}).
382 @return: Angle, wrapped (C{grades}).
383 '''
384 return Float(grades400=(grades(rad) % _400_0) or _0_0) # _umod_400
387def grades2degrees(gon):
388 '''Convert I{grades} (aka I{gons} or I{gradians}) to C{degrees}.
390 @arg gon: Angle (C{grades}).
392 @return: Angle (C{degrees}).
393 '''
394 return Degrees(Float(gon=gon) * _360_0 / _400_0)
397def grades2radians(gon):
398 '''Convert I{grades} (aka I{gons} or I{gradians}) to C{radians}.
400 @arg gon: Angle (C{grades}).
402 @return: Angle (C{radians}).
403 '''
404 return Radians(Float(gon=gon) * PI2 / _400_0)
407def km2m(km):
408 '''Convert kilo meter to meter (m).
410 @arg km: Value in kilo meter (C{scalar}).
412 @return: Value in meter (C{float}).
414 @raise ValueError: Invalid B{C{km}}.
415 '''
416 return Meter(Float(km=km) * _M_KM)
419def _loneg(lon):
420 '''(INTERNAL) "Complement" of C{lon}.
421 '''
422 return _180_0 - lon
425def m2chain(meter):
426 '''Convert meter to I{UK} chains.
428 @arg meter: Value in meter (C{scalar}).
430 @return: Value in C{chains} (C{float}).
432 @raise ValueError: Invalid B{C{meter}}.
433 '''
434 return Float(chain=Meter(meter) / _M_CHAIN) # * 0.049709695378986715
437def m2degrees(distance, radius=R_M, lat=0):
438 '''Convert a distance to an angle along the equator or
439 along the parallel at an other (geodetic) latitude.
441 @arg distance: Distance (C{meter}, same units as B{C{radius}}).
442 @kwarg radius: Mean earth radius, ellipsoid or datum (C{meter},
443 an L{Ellipsoid}, L{Ellipsoid2}, L{Datum} or
444 L{a_f2Tuple}).
445 @kwarg lat: Parallel latitude (C{degrees90}, C{str}).
447 @return: Angle (C{degrees}) or C{INF} for near-polar B{C{lat}}.
449 @raise RangeError: Latitude B{C{lat}} outside valid range and
450 L{pygeodesy.rangerrors} set to C{True}.
452 @raise TypeError: Invalid B{C{radius}}.
454 @raise ValueError: Invalid B{C{distance}}, B{C{radius}}
455 or B{C{lat}}.
457 @see: Function L{m2radians} and L{degrees2m}.
458 '''
459 return degrees(m2radians(distance, radius=radius, lat=lat))
462def m2fathom(meter):
463 '''Convert meter to I{Imperial} fathoms.
465 @arg meter: Value in meter (C{scalar}).
467 @return: Value in C{fathoms} (C{float}).
469 @raise ValueError: Invalid B{C{meter}}.
471 @see: Function L{m2toise}, U{Fathom<https://WikiPedia.org/wiki/Fathom>}
472 and U{Klafter<https://WikiPedia.org/wiki/Klafter>}.
473 '''
474 return Float(fathom=Meter(meter) / _M_FATHOM) # * 0.546806649
477def m2ft(meter, usurvey=False, pied=False, fuss=False):
478 '''Convert meter to I{International}, I{US Survey}, I{French} or
479 or I{German} feet (C{ft}).
481 @arg meter: Value in meter (C{scalar}).
482 @kwarg usurvey: If C{True}, convert to I{US Survey} foot else ...
483 @kwarg pied: If C{True}, convert to French I{pied-du-Roi} else ...
484 @kwarg fuss: If C{True}, convert to German I{Fuss}, otherwise to
485 I{International} foot.
487 @return: Value in C{feet} (C{float}).
489 @raise ValueError: Invalid B{C{meter}}.
490 '''
491 # * 3.2808333333333333, US Survey 3937 / 1200
492 # * 3.2808398950131235, Int'l 10_000 / (254 * 12)
493 return Float(feet=Meter(meter) / (_M_FOOT_USVY if usurvey else
494 (_M_FOOT_FR if pied else
495 (_M_FOOT_GE if fuss else _M_FOOT))))
498def m2furlong(meter):
499 '''Convert meter to furlongs.
501 @arg meter: Value in meter (C{scalar}).
503 @return: Value in C{furlongs} (C{float}).
505 @raise ValueError: Invalid B{C{meter}}.
506 '''
507 return Float(furlong=Meter(meter) / _M_FURLONG) # * 0.00497096954
510def m2km(meter):
511 '''Convert meter to kilo meter (Km).
513 @arg meter: Value in meter (C{scalar}).
515 @return: Value in Km (C{float}).
517 @raise ValueError: Invalid B{C{meter}}.
518 '''
519 return Float(km=Meter(meter) / _M_KM)
522def m2NM(meter):
523 '''Convert meter to nautical miles (NM).
525 @arg meter: Value in meter (C{scalar}).
527 @return: Value in C{NM} (C{float}).
529 @raise ValueError: Invalid B{C{meter}}.
530 '''
531 return Float(NM=Meter(meter) / _M_NM) # * 5.39956804e-4
534def m2radians(distance, radius=R_M, lat=0):
535 '''Convert a distance to an angle along the equator or along the
536 parallel at an other (geodetic) latitude.
538 @arg distance: Distance (C{meter}, same units as B{C{radius}}).
539 @kwarg radius: Mean earth radius, ellipsoid or datum (C{meter},
540 an L{Ellipsoid}, L{Ellipsoid2}, L{Datum} or
541 L{a_f2Tuple}).
542 @kwarg lat: Parallel latitude (C{degrees90}, C{str}).
544 @return: Angle (C{radians}) or C{INF} for near-polar B{C{lat}}.
546 @raise RangeError: Latitude B{C{lat}} outside valid range and
547 L{pygeodesy.rangerrors} set to C{True}.
549 @raise TypeError: Invalid B{C{radius}}.
551 @raise ValueError: Invalid B{C{distance}}, B{C{radius}}
552 or B{C{lat}}.
554 @see: Function L{m2degrees} and L{radians2m}.
555 '''
556 m = circle4(radius, lat).radius
557 return INF if m < EPS0 else Radians(Float(distance=distance) / m)
560def m2SM(meter):
561 '''Convert meter to statute miles (SM).
563 @arg meter: Value in meter (C{scalar}).
565 @return: Value in C{SM} (C{float}).
567 @raise ValueError: Invalid B{C{meter}}.
568 '''
569 return Float(SM=Meter(meter) / _M_SM) # * 6.21369949e-4 == 1 / 1609.344
572def m2toise(meter):
573 '''Convert meter to French U{toises<https://WikiPedia.org/wiki/Toise>}.
575 @arg meter: Value in meter (C{scalar}).
577 @return: Value in C{toises} (C{float}).
579 @raise ValueError: Invalid B{C{meter}}.
581 @see: Function L{m2fathom}.
582 '''
583 return Float(toise=Meter(meter) / _M_TOISE) # * 0.513083632632119
586def m2yard(meter):
587 '''Convert meter to I{UK} yards.
589 @arg meter: Value in meter (C{scalar}).
591 @return: Value in C{yards} (C{float}).
593 @raise ValueError: Invalid B{C{meter}}.
594 '''
595 return Float(yard=Meter(meter) / _M_YARD_UK) # * 1.0936132983377078
598def NM2m(nm):
599 '''Convert nautical miles to meter (m).
601 @arg nm: Value in nautical miles (C{scalar}).
603 @return: Value in meter (C{float}).
605 @raise ValueError: Invalid B{C{nm}}.
606 '''
607 return Meter(Float(nm=nm) * _M_NM)
610def radians2m(rad, radius=R_M, lat=0):
611 '''Convert an angle to a distance along the equator or
612 along the parallel at an other (geodetic) latitude.
614 @arg rad: The angle (C{radians}).
615 @kwarg radius: Mean earth radius, ellipsoid or datum
616 (C{meter}, L{Ellipsoid}, L{Ellipsoid2},
617 L{Datum} or L{a_f2Tuple}).
618 @kwarg lat: Parallel latitude (C{degrees90}, C{str}).
620 @return: Distance (C{meter}, same units as B{C{radius}}
621 or equatorial and polar radii) or C{0.0} for
622 near-polar B{C{lat}}.
624 @raise RangeError: Latitude B{C{lat}} outside valid range and
625 L{pygeodesy.rangerrors} set to C{True}.
627 @raise TypeError: Invalid B{C{radius}}.
629 @raise ValueError: Invalid B{C{rad}}, B{C{radius}} or
630 B{C{lat}}.
632 @see: Function L{degrees2m} and L{m2radians}.
633 '''
634 return _Radians2m(Lam(rad=rad, clip=0), radius, lat)
637def _Radians2m(rad, radius, lat):
638 '''(INTERNAL) Helper for C{degrees2m} and C{radians2m}.
639 '''
640 m = circle4(radius, lat).radius
641 return _0_0 if m < EPS0 else (rad * m)
644def radiansPI(deg):
645 '''Convert and wrap degrees to radians M{[-PI..+PI]}.
647 @arg deg: Angle (C{degrees}).
649 @return: Radians, wrapped (C{radiansPI})
650 '''
651 return wrapPI(radians(deg))
654def radiansPI2(deg):
655 '''Convert and wrap degrees to radians M{[0..+2PI)}.
657 @arg deg: Angle (C{degrees}).
659 @return: Radians, wrapped (C{radiansPI2})
660 '''
661 return _umod_PI2(radians(deg))
664def radiansPI_2(deg):
665 '''Convert and wrap degrees to radians M{[-3PI/2..+PI/2]}.
667 @arg deg: Angle (C{degrees}).
669 @return: Radians, wrapped (C{radiansPI_2})
670 '''
671 return wrapPI_2(radians(deg))
674def _sin0cos2(q, r, sign):
675 '''(INTERNAL) 2-tuple (C{sin(r), cos(r)}) in quadrant C{0 <= B{q} <= 3}
676 and C{sin} zero I{signed} with B{C{sign}}.
677 '''
678 if r < PI_2:
679 s, c = sin(r), cos(r)
680 t = s, c, -s, -c, s
681 else: # r == PI_2
682 t = _1_0, _0_0, _N_1_0, _0_0, _1_0
683# else: # r == 0, testUtility failures
684# t = _0_0, _1_0, _0_0, _N_1_0, _0_0
685# q &= 3
686 s = t[q] or (NEG0 if sign < 0 else _0_0)
687 c = t[q + 1] or _0_0
688 return s, c
691def SinCos2(x):
692 '''Get C{sin} and C{cos} of I{typed} angle.
694 @arg x: Angle (L{Degrees}, L{Degrees_}, L{Radians}, L{Radians_}
695 or scalar C{radians}).
697 @return: 2-Tuple (C{sin(B{x})}, C{cos(B{x})}).
698 '''
699 return sincos2d(x) if isinstanceof(x, Degrees, Degrees_) else (
700 sincos2(x) if isinstanceof(x, Radians, Radians_) else
701 sincos2(float(x))) # assume C{radians}
704def sincos2(rad):
705 '''Return the C{sine} and C{cosine} of an angle in C{radians}.
707 @arg rad: Angle (C{radians}).
709 @return: 2-Tuple (C{sin(B{rad})}, C{cos(B{rad})}).
711 @see: U{GeographicLib<https://GeographicLib.SourceForge.io/C++/doc/
712 classGeographicLib_1_1Math.html#sincosd>} function U{sincosd
713 <https://SourceForge.net/p/geographiclib/code/ci/release/tree/
714 python/geographiclib/geomath.py#l155>} and C++ U{sincosd
715 <https://SourceForge.net/p/geographiclib/code/ci/release/tree/
716 include/GeographicLib/Math.hpp#l558>}.
717 '''
718 if _isfinite(rad):
719 q = int(rad * _2__PI) # int(math.floor)
720 if q < 0:
721 q -= 1
722 t = _sin0cos2(q & 3, rad - q * PI_2, rad)
723 else:
724 t = NAN, NAN
725 return t
728def sincos2_(*rads):
729 '''Return the C{sine} and C{cosine} of angle(s) in C{radians}.
731 @arg rads: One or more angles (C{radians}).
733 @return: Yield the C{sin(B{rad})} and C{cos(B{rad})} for each angle.
735 @see: function L{sincos2}.
736 '''
737 for r in rads:
738 s, c = sincos2(r)
739 yield s
740 yield c
743def sincos2d(deg, **adeg):
744 '''Return the C{sine} and C{cosine} of an angle in C{degrees}.
746 @arg deg: Angle (C{degrees}).
747 @kwarg adeg: Optional correction (C{degrees}).
749 @return: 2-Tuple (C{sin(B{deg_})}, C{cos(B{deg_})}, C{B{deg_} =
750 B{deg} + B{adeg}}).
752 @see: U{GeographicLib<https://GeographicLib.SourceForge.io/C++/doc/
753 classGeographicLib_1_1Math.html#sincosd>} function U{sincosd
754 <https://SourceForge.net/p/geographiclib/code/ci/release/tree/
755 python/geographiclib/geomath.py#l155>} and C++ U{sincosd
756 <https://SourceForge.net/p/geographiclib/code/ci/release/tree/
757 include/GeographicLib/Math.hpp#l558>}.
758 '''
759 if _isfinite(deg):
760 q = int(deg * _1__90) # int(math.floor)
761 if q < 0:
762 q -= 1
763 d = deg - q * _90_0
764 if adeg:
765 t = _xkwds_get(adeg, adeg=_0_0)
766 d = _MODS.karney._around(d + t)
767 t = _sin0cos2(q & 3, radians(d), deg)
768 else:
769 t = NAN, NAN
770 return t
773def sincos2d_(*degs):
774 '''Return the C{sine} and C{cosine} of angle(s) in C{degrees}.
776 @arg degs: One or more angles (C{degrees}).
778 @return: Yield the C{sin(B{deg})} and C{cos(B{deg})} for each angle.
780 @see: Function L{sincos2d}.
781 '''
782 for d in degs:
783 s, c = sincos2d(d)
784 yield s
785 yield c
788def sincostan3(rad):
789 '''Return the C{sine}, C{cosine} and C{tangent} of an angle in C{radians}.
791 @arg rad: Angle (C{radians}).
793 @return: 3-Tuple (C{sin(B{rad})}, C{cos(B{rad})}, C{tan(B{rad})}).
795 @see: Function L{sincos2}.
796 '''
797 s, c = sincos2(float(rad))
798 t = NAN if s is NAN else (_over(s, c) if s else neg(s, neg0=c < 0))
799 return s, c, t
802def SM2m(sm):
803 '''Convert statute miles to meter (m).
805 @arg sm: Value in statute miles (C{scalar}).
807 @return: Value in meter (C{float}).
809 @raise ValueError: Invalid B{C{sm}}.
810 '''
811 return Meter(Float(sm=sm) * _M_SM)
814def tan_2(rad, **semi): # edge=1
815 '''Compute the tangent of half angle.
817 @arg rad: Angle (C{radians}).
818 @kwarg semi: Angle or edge name and index
819 for semi-circular error.
821 @return: M{tan(rad / 2)} (C{float}).
823 @raise ValueError: If B{C{rad}} is semi-circular
824 and B{C{semi}} is given.
825 '''
826 # .formy.excessKarney_, .sphericalTrigonometry.areaOf
827 if semi and isnear0(fabs(rad) - PI):
828 for n, v in semi.items():
829 break
830 n = _SPACE_(n, _radians_) if not isint(v) else \
831 _SPACE_(_MODS.streprs.Fmt.SQUARE(**semi), _edge_)
832 raise _ValueError(n, rad, txt=_semi_circular_)
834 return tan(rad * _0_5) if _isfinite(rad) else NAN
837def tand(deg, **error_kwds):
838 '''Return the C{tangent} of an angle in C{degrees}.
840 @arg deg: Angle (C{degrees}).
841 @kwarg error_kwds: Error to raise (C{ValueError}).
843 @return: C{tan(B{deg})}.
845 @raise ValueError: If L{pygeodesy.isnear0}C{(cos(B{deg})}.
846 '''
847 s, c = sincos2d(deg)
848 if s:
849 if isnear0(c):
850 raise _valueError(tand, deg, **error_kwds)
851 s = s / c # /= chokes PyChecker
852 elif c < 0:
853 s = -s # negate-0
854 return s
857def tand_(*degs, **error_kwds):
858 '''Return the C{tangent} of angle(s) in C{degrees}.
860 @arg degs: One or more angles (C{degrees}).
861 @kwarg error_kwds: Error to raise (C{ValueError}).
863 @return: Yield the C{tan(B{deg})} for each angle.
865 @raise ValueError: See L{pygeodesy.tand}.
866 '''
867 for d in degs:
868 yield tand(d, **error_kwds)
871def tanPI_2_2(rad):
872 '''Compute the tangent of half angle, 90 degrees rotated.
874 @arg rad: Angle (C{radians}).
876 @return: M{tan((rad + PI/2) / 2)} (C{float}).
877 '''
878 return tan((rad + PI_2) * _0_5) if _isfinite(rad) else (
879 NAN if isnan(rad) else (_N_90_0 if rad < 0 else _90_0))
882def toise2m(toises):
883 '''Convert French U{toises<https://WikiPedia.org/wiki/Toise>} to meter.
885 @arg toises: Value in toises (C{scalar}).
887 @return: Value in C{meter} (C{float}).
889 @raise ValueError: Invalid B{C{toises}}.
891 @see: Function L{fathom2m}.
892 '''
893 return Meter(Float(toises=toises) * _M_TOISE)
896def truncate(x, ndigits=None):
897 '''Truncate to the given number of digits.
899 @arg x: Value to truncate (C{scalar}).
900 @kwarg ndigits: Number of digits (C{int}),
901 aka I{precision}.
903 @return: Truncated B{C{x}} (C{float}).
905 @see: Python function C{round}.
906 '''
907 if isint(ndigits):
908 p = _10_0**ndigits
909 x = int(x * p) / p
910 return x
913def unroll180(lon1, lon2, wrap=True):
914 '''Unroll longitudinal delta and wrap longitude in degrees.
916 @arg lon1: Start longitude (C{degrees}).
917 @arg lon2: End longitude (C{degrees}).
918 @kwarg wrap: If C{True}, wrap and unroll to the M{(-180..+180]}
919 range (C{bool}).
921 @return: 2-Tuple C{(B{lon2}-B{lon1}, B{lon2})} unrolled (C{degrees},
922 C{degrees}).
924 @see: Capability C{LONG_UNROLL} in U{GeographicLib
925 <https://GeographicLib.SourceForge.io/html/python/interface.html#outmask>}.
926 '''
927 d = lon2 - lon1
928 if wrap:
929 u = wrap180(d)
930 if u != d:
931 return u, (lon1 + u)
932 return d, lon2
935def _unrollon(p1, p2, wrap=False): # unroll180 == .karney._unroll2
936 '''(INTERNAL) Wrap/normalize, unroll and replace longitude.
937 '''
938 lat, lon = p2.lat, p2.lon
939 if wrap and _Wrap.normal:
940 lat, lon = _Wrap.latlon(lat, lon)
941 _, lon = unroll180(p1.lon, lon, wrap=True)
942 if lat != p2.lat or fabs(lon - p2.lon) > EPS:
943 p2 = p2.dup(lat=lat, lon=wrap180(lon))
944 # p2 = p2.copy(); p2.latlon = lat, wrap180(lon)
945 return p2
948def _unrollon3(p1, p2, p3, wrap=False):
949 '''(INTERNAL) Wrap/normalize, unroll 2 points.
950 '''
951 w = wrap
952 if w:
953 w = _Wrap.normal
954 p2 = _unrollon(p1, p2, wrap=w)
955 p3 = _unrollon(p1, p3, wrap=w)
956 p2 = _unrollon(p2, p3)
957 return p2, p3, w # was wrapped?
960def unrollPI(rad1, rad2, wrap=True):
961 '''Unroll longitudinal delta and wrap longitude in radians.
963 @arg rad1: Start longitude (C{radians}).
964 @arg rad2: End longitude (C{radians}).
965 @kwarg wrap: If C{True}, wrap and unroll to the M{(-PI..+PI]}
966 range (C{bool}).
968 @return: 2-Tuple C{(B{rad2}-B{rad1}, B{rad2})} unrolled
969 (C{radians}, C{radians}).
971 @see: Capability C{LONG_UNROLL} in U{GeographicLib
972 <https://GeographicLib.SourceForge.io/html/python/interface.html#outmask>}.
973 '''
974 r = rad2 - rad1
975 if wrap:
976 u = wrapPI(r)
977 if u != r:
978 return u, (rad1 + u)
979 return r, rad2
982def _valueError(where, x, **kwds):
983 '''(INTERNAL) Return a C{_ValueError}.
984 '''
985 t = _MODS.streprs.Fmt.PAREN(where.__name__, x)
986 return _ValueError(t, **kwds)
989class _Wrap(object):
991 _normal = False # default
993 @property
994 def normal(self):
995 '''Get the current L{normal} setting (C{True},
996 C{False} or C{None}).
997 '''
998 return self._normal
1000 @normal.setter # PYCHOK setter!
1001 def normal(self, setting):
1002 '''Set L{normal} to C{True}, C{False} or C{None}.
1003 '''
1004 t = {True: (_MODS.formy.normal, _MODS.formy.normal_),
1005 False: (self.wraplatlon, self.wraphilam),
1006 None: (_passargs, _passargs)}.get(setting, ())
1007 if t:
1008 self.latlon, self.philam = t
1009 self._normal = setting
1011 def latlonDMS2(self, lat, lon, **DMS2_kwds):
1012 if isstr(lat) or isstr(lon):
1013 kwds = _xkwds(DMS2_kwds, clipLon=0, clipLat=0)
1014 lat, lon = _MODS.dms.parseDMS2(lat, lon, **kwds)
1015 return self.latlon(lat, lon)
1017# def normalatlon(self, *latlon):
1018# return _MODS.formy.normal(*latlon)
1020# def normalamphi(self, *philam):
1021# return _MODS.formy.normal_(*philam)
1023 def wraplatlon(self, lat, lon):
1024 return wrap90(lat), wrap180(lon)
1026 latlon = wraplatlon # default
1028 def latlon3(self, lon1, lat2, lon2, wrap):
1029 if wrap:
1030 lat2, lon2 = self.latlon(lat2, lon2)
1031 lon21, lon2 = unroll180(lon1, lon2)
1032 else:
1033 lon21 = lon2 - lon1
1034 return lon21, lat2, lon2
1036 def _latlonop(self, wrap):
1037 if wrap and self._normal is not None:
1038 return self.latlon
1039 else:
1040 return _passargs
1042 def wraphilam(self, phi, lam):
1043 return wrapPI_2(phi), wrapPI(lam)
1045 philam = wraphilam # default
1047 def philam3(self, lam1, phi2, lam2, wrap):
1048 if wrap:
1049 phi2, lam2 = self.philam(phi2, lam2)
1050 lam21, lam2 = unrollPI(lam1, lam2)
1051 else:
1052 lam21 = lam2 - lam1
1053 return lam21, phi2, lam2
1055 def _philamop(self, wrap):
1056 if wrap and self._normal is not None:
1057 return self.philam
1058 else:
1059 return _passargs
1061 def point(self, ll, wrap=True): # in .points._fractional, -.PointsIter.iterate, ...
1062 '''Return C{ll} or a copy, I{normalized} or I{wrap}'d.
1063 '''
1064 if wrap and self._normal is not None:
1065 lat, lon = ll.latlon
1066 if fabs(lon) > 180 or fabs(lat) > 90:
1067 _n = self.latlon
1068 ll = ll.copy(name=_n.__name__)
1069 ll.latlon = _n(lat, lon)
1070 return ll
1072_Wrap = _Wrap() # PYCHOK singleton
1075# def _wrap(angle, wrap, modulo):
1076# '''(INTERNAL) Angle wrapper M{((wrap-modulo)..+wrap]}.
1077#
1078# @arg angle: Angle (C{degrees}, C{radians} or C{grades}).
1079# @arg wrap: Range (C{degrees}, C{radians} or C{grades}).
1080# @arg modulo: Upper limit (360 C{degrees}, PI2 C{radians} or 400 C{grades}).
1081#
1082# @return: The B{C{angle}}, wrapped (C{degrees}, C{radians} or C{grades}).
1083# '''
1084# a = float(angle)
1085# if not (wrap - modulo) <= a < wrap:
1086# # math.fmod(-1.5, 3.14) == -1.5, but -1.5 % 3.14 == 1.64
1087# # math.fmod(-1.5, 360) == -1.5, but -1.5 % 360 == 358.5
1088# a %= modulo
1089# if a > wrap:
1090# a -= modulo
1091# return a
1094def wrap90(deg):
1095 '''Wrap degrees to M{[-90..+90]}.
1097 @arg deg: Angle (C{degrees}).
1099 @return: Degrees, wrapped (C{degrees90}).
1100 '''
1101 w = wrap180(deg)
1102 return (w - _180_0) if w > 90 else ((w + _180_0) if w < -90 else w)
1105def wrap180(deg):
1106 '''Wrap degrees to M{[-180..+180]}.
1108 @arg deg: Angle (C{degrees}).
1110 @return: Degrees, wrapped (C{degrees180}).
1111 '''
1112 d = float(deg)
1113 w = _umod_360(d)
1114 if w > 180:
1115 w -= _360_0
1116 elif d < 0 and w == 180:
1117 w = -w
1118 return w
1121def wrap360(deg): # see .streprs._umod_360
1122 '''Wrap degrees to M{[0..+360)}.
1124 @arg deg: Angle (C{degrees}).
1126 @return: Degrees, wrapped (C{degrees360}).
1127 '''
1128 return _umod_360(float(deg))
1131def wrapPI(rad):
1132 '''Wrap radians to M{[-PI..+PI]}.
1134 @arg rad: Angle (C{radians}).
1136 @return: Radians, wrapped (C{radiansPI}).
1137 '''
1138 r = float(rad)
1139 w = _umod_PI2(r)
1140 if w > PI:
1141 w -= PI2
1142 elif r < 0 and w == PI:
1143 w = -PI
1144 return w
1147def wrapPI2(rad):
1148 '''Wrap radians to M{[0..+2PI)}.
1150 @arg rad: Angle (C{radians}).
1152 @return: Radians, wrapped (C{radiansPI2}).
1153 '''
1154 return _umod_PI2(float(rad))
1157def wrapPI_2(rad):
1158 '''Wrap radians to M{[-PI/2..+PI/2]}.
1160 @arg rad: Angle (C{radians}).
1162 @return: Radians, wrapped (C{radiansPI_2}).
1163 '''
1164 w = wrapPI(rad)
1165 return (w - PI) if w > PI_2 else ((w + PI) if w < (-PI_2) else w)
1168# def wraplatlon(lat, lon):
1169# '''Both C{wrap90(B{lat})} and C{wrap180(B{lon})}.
1170# '''
1171# return wrap90(lat), wrap180(lon)
1174def wrap_normal(*normal):
1175 '''Define the operation for the keyword argument C{B{wrap}=True},
1176 across L{pygeodesy}: I{wrap}, I{normalize} or I{no-op}. For
1177 backward compatibility, the default is I{wrap}.
1179 @arg normal: If C{True}, I{normalize} lat- and longitude using
1180 L{normal} or L{normal_}, if C{False}, I{wrap} the
1181 lat- and longitude individually by L{wrap90} or
1182 L{wrapPI_2} respectively L{wrap180}, L{wrapPI} or
1183 if C{None}, leave lat- and longitude I{unchanged}.
1184 Do not supply any value to get the current setting.
1186 @return: The previous L{wrap_normal} setting (C{bool} or C{None}).
1187 '''
1188 t = _Wrap.normal
1189 if normal:
1190 _Wrap.normal = normal[0]
1191 return t
1194# def wraphilam(phi, lam,):
1195# '''Both C{wrapPI_2(B{phi})} and C{wrapPI(B{lam})}.
1196# '''
1197# return wrapPI_2(phi), wrapPI(lam)
1200def yard2m(yards):
1201 '''Convert I{UK} yards to meter.
1203 @arg yards: Value in yards (C{scalar}).
1205 @return: Value in C{meter} (C{float}).
1207 @raise ValueError: Invalid B{C{yards}}.
1208 '''
1209 return Float(yards=yards) * _M_YARD_UK
1211# **) MIT License
1212#
1213# Copyright (C) 2016-2024 -- mrJean1 at Gmail -- All Rights Reserved.
1214#
1215# Permission is hereby granted, free of charge, to any person obtaining a
1216# copy of this software and associated documentation files (the "Software"),
1217# to deal in the Software without restriction, including without limitation
1218# the rights to use, copy, modify, merge, publish, distribute, sublicense,
1219# and/or sell copies of the Software, and to permit persons to whom the
1220# Software is furnished to do so, subject to the following conditions:
1221#
1222# The above copyright notice and this permission notice shall be included
1223# in all copies or substantial portions of the Software.
1224#
1225# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
1226# OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
1227# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
1228# THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR
1229# OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE,
1230# ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
1231# OTHER DEALINGS IN THE SOFTWARE.