Coverage for pygeodesy/karney.py : 97%

Hot-keys on this page
r m x p toggle line displays
j k next/prev highlighted chunk
0 (zero) top of page
1 (one) first highlighted chunk
# -*- coding: utf-8 -*-
Wrapper around Python classes C{Geodesic} and C{GeodesicLine} and several C{Math} functions from I{Karney}'s Python package U{geographiclib<https://PyPI.org/project/geographiclib>}, provided that package is installed.
The I{wrapped} class methods return a L{GDict} instance offering access to the C{dict} items either by C{key} or by C{attribute} name.
Karney-based functionality ==========================
1. The following classes and functions in C{pygeodesy}
- L{AlbersEqualArea}, L{AlbersEqualArea2}, L{AlbersEqualArea4}, L{AlbersEqualAreaCylindrical}, L{AlbersEqualAreaNorth}, L{AlbersEqualAreaSouth} -- U{AlbersEqualArea<https://GeographicLib.SourceForge.io/html/ classGeographicLib_1_1AlbersEqualArea.html>}
- L{CassiniSoldner} -- U{CassiniSoldner<https://GeographicLib.SourceForge.io/html/ classGeographicLib_1_1CassiniSoldner.html>}
- L{EcefKarney} -- U{Geocentric<https://GeographicLib.SourceForge.io/html/ classGeographicLib_1_1Geocentric.html>}
- L{Elliptic} -- U{EllipticFunction<https://GeographicLib.SourceForge.io/html/ classGeographicLib_1_1EllipticFunction.html>}
- L{EquidistantExact}, L{EquidistantGeodSolve}, L{EquidistantKarney} -- U{AzimuthalEquidistant <https://GeographicLib.SourceForge.io/html/classGeographicLib_1_1AzimuthalEquidistant.html>}
- L{Etm}, L{ExactTransverseMercator} -- U{TransverseMercatorExact <https://GeographicLib.SourceForge.io/html/classGeographicLib_1_1TransverseMercatorExact.html>}
- L{GeodesicAreaExact}, L{PolygonArea} -- U{PolygonArea<https://GeographicLib.SourceForge.io/ html/classGeographicLib_1_1PolygonAreaT.html>}
- L{GeodesicExact}, L{GeodesicLineExact} -- U{GeodesicExact<https://GeographicLib.SourceForge.io/ html/classGeographicLib_1_1GeodesicExact.html>}, U{GeodesicLineExact<https://GeographicLib.SourceForge.io/ html/classGeographicLib_1_1GeodesicLineExact.html>}
- L{GeoidKarney} -- U{Geoid<https://GeographicLib.SourceForge.io/html/geoid.html>}
- L{Georef} -- U{Georef<https://GeographicLib.SourceForge.io/html/ classGeographicLib_1_1Georef.html>}
- L{GnomonicExact}, L{GnomonicGeodSolve}, L{GnomonicKarney} -- U{Gnomonic <https://GeographicLib.SourceForge.io/html/classGeographicLib_1_1Gnomonic.html>}
- L{LocalCartesian}, L{Ltp} -- U{LocalCartesian<https://GeographicLib.SourceForge.io/html/ classGeographicLib_1_1LocalCartesian.html>}
- L{Ups} -- U{PolarStereographic<https://GeographicLib.SourceForge.io/html/ classGeographicLib_1_1PolarStereographic.html>}
- L{Utm} -- U{TransverseMercator<https://GeographicLib.SourceForge.io/html/ classGeographicLib_1_1TransverseMercator.html>}
- L{UtmUps}, L{Epsg} -- U{UTMUPS<https://GeographicLib.SourceForge.io/html/ classGeographicLib_1_1UTMUPS.html>}
- L{pygeodesy.atand}, L{pygeodesy.atan2d}, L{pygeodesy.sincos2}, L{pygeodesy.sincos2d} -- U{ Math<https://GeographicLib.SourceForge.io/html/classGeographicLib_1_1Math.html>}
are I{transcoded} from C++ classes in I{Karney}'s U{GeographicLib<https://GeographicLib.SourceForge.io/html/annotated.html>}.
2. These C{pygeodesy} modules and classes
- L{ellipsoidalGeodSolve}, L{ellipsoidalKarney}, L{geodsolve}, L{karney} - L{EquidistantKarney}, L{FrechetKarney}, L{GeodesicSolve}, L{GeodesicLineSolve}, L{GnomonicGeodSolve}, L{GnomonicKarney}, L{HeightIDWkarney}
are or use I{wrappers} around I{Karney}'s Python U{geographiclib<https://PyPI.org/project/geographiclib>} C{geodesic} or C++ utility U{GeodSolve<https://GeographicLib.SourceForge.io/html/GeodSolve.1.html>}.
3. All C{pygeodesy} functions and methods to compute I{ellipsoidal} intersections and trilaterations
- L{ellipsoidalExact.intersection3}, L{ellipsoidalExact.intersections2}, L{ellipsoidalExact.nearestOn}, L{ellipsoidalExact.LatLon.intersection3}, L{ellipsoidalExact.LatLon.intersections2}, L{ellipsoidalExact.LatLon.nearestOn}, L{ellipsoidalExact.LatLon.trilaterate5}
- L{ellipsoidalKarney.intersection3}, L{ellipsoidalKarney.intersections2}, L{ellipsoidalKarney.nearestOn}, L{ellipsoidalKarney.LatLon.intersection3}, L{ellipsoidalKarney.LatLon.intersections2}, L{ellipsoidalKarney.LatLon.nearestOn}, L{ellipsoidalKarney.LatLon.trilaterate5}
- L{ellipsoidalVincenty.intersection3}, L{ellipsoidalVincenty.intersections2}, L{ellipsoidalVincenty.nearestOn}, L{ellipsoidalVincenty.LatLon.intersection3}, L{ellipsoidalVincenty.LatLon.intersections2}, L{ellipsoidalVincenty.LatLon.nearestOn}, L{ellipsoidalVincenty.LatLon.trilaterate5}
are implementations of I{Karney}'s solution posted under U{The B{ellipsoidal} case <https://GIS.StackExchange.com/questions/48937/calculating-intersection-of-two-circles>} and in paper U{Geodesics on an ellipsoid of revolution<https://ArXiv.org/pdf/1102.1215.pdf>} (pp 20-21, section B{14. MARITIME BOUNDARIES}).
4. Spherical functions
- L{pygeodesy.excessKarney_}, L{sphericalTrigonometry.areaOf}
in C{pygeodesy} are based on I{Karney}'s post U{Area of a spherical polygon <https://MathOverflow.net/questions/97711/the-area-of-spherical-polygons>}, 3rd Answer. '''
_xgeographiclib, _xImportError, isodd # PYCHOK shared # from pygeodesy.ellipsoids import Ellipsoid2 # from .datums _0_0, _1_0, _2_0, _16_0, _180_0, _N_180_0, _360_0 _NamedTuple, _Pass, unstr # PYCHOK shared # from pygeodesy.streps import unstr # from .named Meter as _M, Meter2 as _M2, _1mm as _TOL_M # PYCHOK shared
'''(INTERNAL) Get an ellipsoid from C{(B{a_..}, B{f})} or C{B{.._ellipsoid}}. ''' _ellipsoidal_datum(a_ellipsoid, name=name, raiser=raiser).ellipsoid
'''(INTERNAL) Latitude B{C{lat}}. '''
'''(INTERNAL) Longitude B{C{lon}}. '''
def _raiseX(inst, x, *args): # PYCHOK no cover '''(INTERNAL) Throw a C{GeodesicError} for C{geographiclib} issue B{C{x}} . ''' n = _DOT_(classname(inst), callername(up=2, underOK=True)) raise GeodesicError(unstr(n, *args), txt=str(x))
'''(INTERNAL) Helper. ''' '''Convert this C{*Tuple} to a L{GDict}.
@kwarg updates: Optional items to apply (C{nam=value} pairs) '''
'''9-Tuple C{(a12, lat2, lon2, azi2, s12, m12, M12, M21, S12)} with arc length C{a12}, angles C{lat2}, C{lon2} and azimuth C{azi2} in C{degrees}, distance C{s12} and reduced length C{m12} in C{meter} and area C{S12} in C{meter} I{squared}. '''
'''Basic C{dict} with both key I{and} attribute access to the C{dict} items.
Results of all C{geodesic} methods are returned as a L{GDict} instance. ''' '''Convert this L{GDict} result to a 9-tuple, like I{Karney}'s method C{geographiclib.geodesic.Geodesic._GenDirect}.
@kwarg dflt: Default value for missing items (C{any}).
@return: L{Direct9Tuple}C{(a12, lat2, lon2, azi2, s12, m12, M12, M21, S12)} '''
'''Convert this L{GDict} result to a 12-Tuple, compatible with I{Karney}'s U{GeodSolve<https://GeographicLib.SourceForge.io/html/GeodSolve.1.html>} result.
@kwarg dflt: Default value for missing items (C{any}).
@return: L{GeodSolve12Tuple}C{(lat1, lon1, azi1, lat2, lon2, azi2, s12, a12, m12, M12, M21, S12)}. ''' return self._toTuple(GeodSolve12Tuple, dflt)
'''Convert this L{GDict} result to a 10-tuple, like I{Karney}'s method C{geographiclib.geodesic.Geodesic._GenInverse}.
@kwarg dflt: Default value for missing items (C{any}).
@return: L{Inverse10Tuple}C{(a12, s12, salp1, calp1, salp2, calp2, m12, M12, M21, S12)}. '''
'''(INTERNAL) Convert this C{GDict} to an B{C{nTuple}}. '''
'''Error raised for L{pygeodesy.geodesicx} lack of convergence or other L{pygeodesy.geodesicx} or L{pygeodesy.karney} issues. '''
'''12-Tuple C{(lat1, lon1, azi1, lat2, lon2, azi2, s12, a12, m12, M12, M21, S12)} with angles C{lat1}, C{lon1}, C{azi1}, C{lat2}, C{lon2} and C{azi2} and arc C{a12} all in C{degrees}, initial C{azi1} and final C{azi2} forward azimuths, distance C{s12} and reduced length C{m12} in C{meter}, area C{S12} in C{meter} I{squared} and geodesic scale factors C{M12} and C{M21}, both C{scalar}, see U{GeodSolve <https://GeographicLib.SourceForge.io/html/GeodSolve.1.html>}. ''' # from GeodSolve --help option -f ... lat1 lon1 azi1 lat2 lon2 azi2 s12 a12 m12 M12 M21 S12
'''10-Tuple C{(a12, s12, salp1, calp1, salp2, calp2, m12, M12, M21, S12)} with arc length C{a12} in C{degrees}, distance C{s12} and reduced length C{m12} in C{meter}, area C{S12} in C{meter} I{squared} and the sines C{salp1}, C{salp2} and cosines C{calp1}, C{calp2} of the initial C{1} and final C{2} foward azimuths. '''
'''Convert this C{Inverse10Tuple} to a L{GDict}.
@kwarg updates: Optional items to apply (C{nam=value} pairs) ''' azi2=atan2d(self.salp2, self.calp2), # PYCHOK namedTuple **updates) # PYCHOK indent
''''(INTERNAL) Wrapper for some of I{Karney}'s U{geographiclib <https://PyPI.org/project/geographiclib>} classes. '''
'''Get the I{wrapped} C{Geodesic} class, provided the U{geographiclib <https://PyPI.org/project/geographiclib>} package is installed, otherwise an C{ImportError}. '''
'''I{Karney}'s U{Geodesic<https://GeographicLib.SourceForge.io/html/ python/code.html#geographiclib.geodesic.Geodesic>} wrapper. '''
'''New C{Geodesic} instance.
@arg a_ellipsoid: An ellipsoid (L{Ellipsoid}) or datum (L{Datum}) or the equatorial radius of the ellipsoid (C{meter}). @arg f: The flattening of the ellipsoid (C{scalar}) if B{C{a_ellipsoid}) is specified as C{meter}. @kwarg name: Optional name (C{str}). ''' except (TypeError, ValueError) as x: _raiseX(self, x, *self.ellipsoid.a_f)
'''Get the C{debug} option (C{bool}). ''' return bool(self._debug)
'''Set the C{debug} option.
@arg debug: Include more details in results (C{bool}). '''
'''Return the C{Direct} result. ''' except (TypeError, ValueError) as x: _raiseX(self, x, lat1, lon1, azi1, s12, *outmask)
'''Return the destination lat, lon and reverse azimuth (final bearing) in C{degrees}.
@return: L{Destination3Tuple}C{(lat, lon, final)}. '''
'''Get this geodesic's ellipsoid (C{Ellipsoid[2]}). '''
outmask=_Geodesic.STANDARD): '''(INTERNAL) Get C{._GenDirect} result as C{GDict}. ''' except (TypeError, ValueError) as x: _raiseX(self, x, lat, lon, azi, arcmode, s12_a12, outmask)
'''(INTERNAL) Get C{._GenInverse} result as C{GDict}. ''' except (TypeError, ValueError) as x: _raiseX(self, x, lat1, lon1, lat2, lon2, outmask)
'''Return the C{Inverse} result. ''' except (TypeError, ValueError) as x: _raiseX(self, x, lat1, lon1, lat2, lon2, *outmask)
'''Return the non-negative, I{angular} distance in C{degrees}. ''' # see .FrechetKarney.distance, .HausdorffKarney._distance # and .HeightIDWkarney._distances # XXX self.DISTANCE needed for 'a12'?
'''Return the distance in C{meter} and the forward and reverse azimuths (initial and final bearing) in C{degrees}.
@return: L{Distance3Tuple}C{(distance, initial, final)}. '''
'''Set up a L{GeodesicLine} to compute several points on a single geodesic. '''
# PYCHOK .geodesicx.GeodesictExact._LineTemp and -._GDictDirect
# Geodesic.Direct.__doc__ = _Geodesic.Direct.__doc__ # Geodesic.Inverse.__doc__ = _Geodesic.Inverse.__doc__ # Geodesic.Line.__doc__ = _Geodesic.Line.__doc__
'''Get the I{wrapped} C{GeodesicLine} class, provided the U{geographiclib <https://PyPI.org/project/geographiclib>} package is installed, otherwise an C{ImportError}. '''
'''I{Karney}'s U{GeodesicLine <https://GeographicLib.SourceForge.io/html/ python/code.html#geographiclib.geodesicline.GeodesicLine>} wrapper. ''' except (TypeError, ValueError) as x: _raiseX(self, x, lat1, lon1, azi1, *caps)
except (TypeError, ValueError) as x: _raiseX(self, x, a12, *outmask)
except (TypeError, ValueError) as x: _raiseX(self, x, s12, *outmask)
# GeodesicLine.ArcPosition.__doc__ = _GeodesicLine.ArcPosition.__doc__ # GeodesicLine.Position.__doc__ = _GeodesicLine.Position.__doc__
'''Get the I{wrapped} C{Geodesic.WGS84} I{instance} provided the U{geographiclib<https://PyPI.org/project/geographiclib>} package is installed, otherwise an C{ImportError}. '''
'''Get the imported C{geographiclib}, provided the U{geographiclib <https://PyPI.org/project/geographiclib>} package is installed, otherwise an C{ImportError}. '''
'''Get the C{Math} class, provided the U{geographiclib <https://PyPI.org/project/geographiclib>} package is installed, otherwise C{None}. ''' # replace karney. with Math. functions except AttributeError: # 2.0 k._isfinite = _math_isfinite except AttributeError: pass except ImportError: M = None
'''Coarsen a scalar to zero.
@return: Coarsened value (C{float}). '''
else:
'''Compute C{deg - deg0}, reduced to C{[-180,180]} accurately.
@return: 2-Tuple C{(delta_angle, residual)} in C{degrees}. '''
d = _N_180_0
# def _Equidistant(equidistant, exact=False, geodsolve=False): # # (INTERNAL) Get the C{EquidistantExact}, C{-GeodSolve} or # # C{-Karney} class if B{C{equidistant}} in not callable. # if equidistant is None or not callable(equidistant): # if exact: # equidistant = _MODS.azimuthal.EquidistantExact # elif geodsolve: # equidistant = _MODS.azimuthal.EquidistantGeodSolve # else: # equidistant = _MODS.azimuthal.EquidistantKarney # return equidistant
'''Replace angle in C{degrees} outside [-90,90] by NAN.
@return: Angle C{degrees} or NAN. '''
'''Check finiteness of C{x}.
@return: C{True} if finite. '''
'''Reduce angle in C{degrees} to (-180,180].
@return: Reduced angle C{degrees}. '''
'''(INTERNAL) Compute the area or perimeter of a polygon, using a L{GeodesicExact}, L{GeodesicSolve} or (if the C{geographiclib} package is installed) a C{Geodesic} or C{_wrapped.Geodesic} instance. ''' raise _ValueError(wrap=wrap)
# note, lon deltas are unrolled, by default pA(p0.lat, p0.lon)
# gP.Compute returns (number_of_points, perimeter, signed area)
'''Remainder of C{x / y}.
@return: Remainder in the range M{[-y / 2, y / 2]}. ''' except AttributeError: pass
'''(INTERNAL) Remainder/modulo in the range C{[-B{y_2}, B{y_2}]}. ''' # with Python 2.7.16 and 3.7.3 on macOS 10.13.6 and # with Python 3.10.2 on macOS 12.2.1 M1 arm64 native # fmod( 0, 360) == 0.0 # fmod( 360, 360) == 0.0 # fmod(-0, 360) == 0.0 # fmod(-0.0, 360) == -0.0 # fmod(-360, 360) == -0.0 # however, using the % operator ... # 0 % 360 == 0 # 360 % 360 == 0 # 360.0 % 360 == 0.0 # -0 % 360 == 0 # -360 % 360 == 0 == (-360) % 360 # -0.0 % 360 == 0.0 == (-0.0) % 360 # -360.0 % 360 == 0.0 == (-360.0) % 360
# On Windows 32-bit with python 2.7, math.fmod(-0.0, 360) # == +0.0. This fixes this bug. See also Math::AngNormalize # in the C++ library, Math.sincosd has a similar fix. else: # PYCHOK no cover x = NAN
'''Error-free summation like C{Math::sum}.
@return: 2-Tuple C{(B{u} + B{v}, residual)}.
@note: The C{residual} can be the same as B{C{u}} or B{C{v}}.
@see: U{Algorithm 3.1<https://www.TUHH.De/ti3/paper/rump/OgRuOi05.pdf>}. ''' return M.sum(u, v)
# if Algorithm_3_1: # t = (u - t) + (v + r) # elif C_CPP: # Math::sum C/C++ # r -= u # t -= v # t += r # t = -t # else:
'''Accumulate any B{C{vs}} into a previous C{_sum2(s, t)}.
@return: 2-Tuple C{(B{s} + B{t} + B{vs}, residual)}.
@see: I{Karney's} C++ U{Accumulator<https://GeographicLib.SourceForge.io/ html/Accumulator_8hpp_source.html>} comments for more details and function C{_sum2} above.
@note: NOT "error-free", see C{pygeodesy.test/testKarney.py}. ''' # elif t: # s == 0 implies t == 0 # raise _AssertionError(t=t, txt=_not(_0_)) else: else:
'''Unroll B{C{lon2 - lon1}} like C{geodesic.Geodesic.Inverse}.
@return: 2-Tuple C{(B{lon2} - B{lon1}, B{lon2})} with B{C{lon2}} unrolled if B{C{wrap}} is C{True}, normalized otherwise. ''' else: lon2 = _norm180(lon2)
# **) MIT License # # Copyright (C) 2016-2022 -- mrJean1 at Gmail -- All Rights Reserved. # # Permission is hereby granted, free of charge, to any person obtaining a # copy of this software and associated documentation files (the "Software"), # to deal in the Software without restriction, including without limitation # the rights to use, copy, modify, merge, publish, distribute, sublicense, # and/or sell copies of the Software, and to permit persons to whom the # Software is furnished to do so, subject to the following conditions: # # The above copyright notice and this permission notice shall be included # in all copies or substantial portions of the Software. # # THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS # OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, # FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL # THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR # OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, # ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR # OTHER DEALINGS IN THE SOFTWARE. |