Coverage for pygeodesy/formy.py : 98%

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 -*-
''' # make sure int/int division yields float quotient, see .basics
_mean_radius, _spherical_datum, _WGS84 # from pygeodesy.ellipsoids import Ellipsoid # from .datums LimitError, _limiterrors, _ValueError _distant_, _inside_, _near_, _null_, _opposite_, \ _outside_, _too_, _0_0, _0_125, _0_25, _0_5, \ _1_0, _2_0, _4_0, _32_0, _90_0, _180_0, _360_0 LatLon2Tuple, PhiLam2Tuple, Vector3Tuple # from pygeodesy.streprs import unstr # from .fsums Lon, Phi_, Radians, Radians_, Radius, Radius_, \ Scalar, _100km m2degrees, sincos2, sincos2_, tan_2, unroll180, \ unrollPI, wrap90, wrap180, wrapPI, wrapPI_2
'''Return the antipode, the point diametrically opposite to a given point in C{degrees}.
@arg lat: Latitude (C{degrees}). @arg lon: Longitude (C{degrees}).
@return: A L{LatLon2Tuple}C{(lat, lon)}.
@see: U{Geosphere<https://CRAN.R-Project.org/web/packages/geosphere/geosphere.pdf>}. '''
'''Return the antipode, the point diametrically opposite to a given point in C{radians}.
@arg phi: Latitude (C{radians}). @arg lam: Longitude (C{radians}).
@return: A L{PhiLam2Tuple}C{(phi, lam)}.
@see: U{Geosphere<https://CRAN.R-Project.org/web/packages/geosphere/geosphere.pdf>}. '''
'''(INTERNAL) Helper for area and spherical excess. ''' Phi_(lat1=lat1), radians(d_lon)) r *= _mean_radius(radius, lat1, lat2)**2
'''Compute the initial or final bearing (forward or reverse azimuth) between a (spherical) start and end point.
@arg lat1: Start latitude (C{degrees}). @arg lon1: Start longitude (C{degrees}). @arg lat2: End latitude (C{degrees}). @arg lon2: End longitude (C{degrees}). @kwarg options: Optional keyword arguments for function L{pygeodesy.bearing_}.
@return: Initial or final bearing (compass C{degrees360}) or zero if start and end point coincide. ''' Phi_(lat2=lat2), Lam_(lon2=lon2), **options)
'''Compute the initial or final bearing (forward or reverse azimuth) between a (spherical) start and end point.
@arg phi1: Start latitude (C{radians}). @arg lam1: Start longitude (C{radians}). @arg phi2: End latitude (C{radians}). @arg lam2: End longitude (C{radians}). @kwarg final: Return final bearing if C{True}, initial otherwise (C{bool}). @kwarg wrap: Wrap and L{pygeodesy.unrollPI} longitudes (C{bool}).
@return: Initial or final bearing (compass C{radiansPI2}) or zero if start and end point coincide.
@see: U{Bearing<https://www.Movable-Type.co.UK/scripts/latlong.html>}, U{Course between two points<https://www.EdWilliams.org/avform147.htm#Crs>} and U{Bearing Between Two Points<https://web.Archive.org/web/20020630205931/ https://MathForum.org/library/drmath/view/55417.html>}. ''' else:
'''(INTERNAL) Compute initial and final bearing. ''' except AttributeError: pass # XXX spherical version, OK for ellipsoidal ispolar? degrees(bearing_(a1, b1, a2, b2, final=True, wrap=wrap)), name=_bearingTo2.__name__)
'''Return the angle from North for the direction vector M{(lon2 - lon1, lat2 - lat1)} between two points.
Suitable only for short, not near-polar vectors up to a few hundred Km or Miles. Use function L{pygeodesy.bearing} for longer vectors.
@arg lat1: From latitude (C{degrees}). @arg lon1: From longitude (C{degrees}). @arg lat2: To latitude (C{degrees}). @arg lon2: To longitude (C{degrees}). @kwarg adjust: Adjust the longitudinal delta by the cosine of the mean latitude (C{bool}). @kwarg wrap: Wrap and L{pygeodesy.unroll180} longitudes (C{bool}).
@return: Compass angle from North (C{degrees360}).
@note: Courtesy of Martin Schultz.
@see: U{Local, flat earth approximation <https://www.EdWilliams.org/avform.htm#flat>}. '''
'''Compute the distance between two (ellipsoidal) points using the U{Andoyer-Lambert correction<https://navlib.net/wp-content/uploads/ 2013/10/admiralty-manual-of-navigation-vol-1-1964-english501c.pdf>} of the U{Law of Cosines<https://www.Movable-Type.co.UK/scripts/latlong.html#cosine-law>} fromula.
@arg lat1: Start latitude (C{degrees}). @arg lon1: Start longitude (C{degrees}). @arg lat2: End latitude (C{degrees}). @arg lon2: End longitude (C{degrees}). @kwarg datum: Datum or ellipsoid to use (L{Datum}, L{Ellipsoid}, L{Ellipsoid2} or L{a_f2Tuple}). @kwarg wrap: Wrap and L{pygeodesy.unroll180} longitudes (C{bool}).
@return: Distance (C{meter}, same units as the B{C{datum}}'s ellipsoid axes or C{radians} if B{C{datum}} is C{None}).
@raise TypeError: Invalid B{C{datum}}.
@see: Functions L{cosineAndoyerLambert_}, L{cosineForsytheAndoyerLambert}, L{cosineLaw}, L{equirectangular}, L{euclidean}, L{flatLocal}/L{hubeny}, L{flatPolar}, L{haversine}, L{thomas} and L{vincentys} and method L{Ellipsoid.distance2}. ''' *unroll180(lon1, lon2, wrap=wrap))
'''Compute the I{angular} distance between two (ellipsoidal) points using the U{Andoyer-Lambert correction<https://navlib.net/wp-content/uploads/2013/10/ admiralty-manual-of-navigation-vol-1-1964-english501c.pdf>} of the U{Law of Cosines<https://www.Movable-Type.co.UK/scripts/latlong.html#cosine-law>} fromula.
@arg phi2: End latitude (C{radians}). @arg phi1: Start latitude (C{radians}). @arg lam21: Longitudinal delta, M{end-start} (C{radians}). @kwarg datum: Datum or ellipsoid to use (L{Datum}, L{Ellipsoid}, L{Ellipsoid2} or L{a_f2Tuple}).
@return: Angular distance (C{radians}).
@raise TypeError: Invalid B{C{datum}}.
@see: Functions L{cosineAndoyerLambert}, L{cosineForsytheAndoyerLambert_}, L{cosineLaw_}, L{equirectangular_}, L{euclidean_}, L{flatLocal_}/L{hubeny_}, L{flatPolar_}, L{haversine_}, L{thomas_} and L{vincentys_} and U{Geodesy-PHP <https://GitHub.com/jtejido/geodesy-php/blob/master/src/Geodesy/ Distance/AndoyerLambert.php>}. '''
'''Compute the distance between two (ellipsoidal) points using the U{Forsythe-Andoyer-Lambert correction<https://www2.UNB.Ca/gge/Pubs/TR77.pdf>} of the U{Law of Cosines<https://www.Movable-Type.co.UK/scripts/latlong.html#cosine-law>} formula.
@arg lat1: Start latitude (C{degrees}). @arg lon1: Start longitude (C{degrees}). @arg lat2: End latitude (C{degrees}). @arg lon2: End longitude (C{degrees}). @kwarg datum: Datum or ellipsoid to use (L{Datum}, L{Ellipsoid}, L{Ellipsoid2} or L{a_f2Tuple}). @kwarg wrap: Wrap and L{pygeodesy.unroll180} longitudes (C{bool}).
@return: Distance (C{meter}, same units as the B{C{datum}}'s ellipsoid axes or C{radians} if B{C{datum}} is C{None}).
@raise TypeError: Invalid B{C{datum}}.
@see: Functions L{cosineForsytheAndoyerLambert_}, L{cosineAndoyerLambert}, L{cosineLaw}, L{equirectangular}, L{euclidean}, L{flatLocal}/L{hubeny}, L{flatPolar}, L{haversine}, L{thomas} and L{vincentys} and method L{Ellipsoid.distance2}. ''' *unroll180(lon1, lon2, wrap=wrap))
'''Compute the I{angular} distance between two (ellipsoidal) points using the U{Forsythe-Andoyer-Lambert correction<https://www2.UNB.Ca/gge/Pubs/TR77.pdf>} of the U{Law of Cosines<https://www.Movable-Type.co.UK/scripts/latlong.html#cosine-law>} formula.
@arg phi2: End latitude (C{radians}). @arg phi1: Start latitude (C{radians}). @arg lam21: Longitudinal delta, M{end-start} (C{radians}). @kwarg datum: Datum or ellipsoid to use (L{Datum}, L{Ellipsoid}, L{Ellipsoid2} or L{a_f2Tuple}).
@return: Angular distance (C{radians}).
@raise TypeError: Invalid B{C{datum}}.
@see: Functions L{cosineForsytheAndoyerLambert}, L{cosineAndoyerLambert_}, L{cosineLaw_}, L{equirectangular_}, L{euclidean_}, L{flatLocal_}/L{hubeny_}, L{flatPolar_}, L{haversine_}, L{thomas_} and L{vincentys_} and U{Geodesy-PHP <https://GitHub.com/jtejido/geodesy-php/blob/master/src/Geodesy/ Distance/ForsytheCorrection.php>}. '''
'''Compute the distance between two points using the U{spherical Law of Cosines <https://www.Movable-Type.co.UK/scripts/latlong.html#cosine-law>} formula.
@arg lat1: Start latitude (C{degrees}). @arg lon1: Start longitude (C{degrees}). @arg lat2: End latitude (C{degrees}). @arg lon2: End longitude (C{degrees}). @kwarg radius: Mean earth radius, ellipsoid or datum (C{meter}, L{Ellipsoid}, L{Ellipsoid2}, L{Datum} or L{a_f2Tuple}). @kwarg wrap: Wrap and L{pygeodesy.unroll180} longitudes (C{bool}).
@return: Distance (C{meter}, same units as B{C{radius}} or the ellipsoid or datum axes).
@raise TypeError: Invalid B{C{radius}}.
@see: Functions L{cosineLaw_}, L{cosineAndoyerLambert}, L{cosineForsytheAndoyerLambert}, L{equirectangular}, L{euclidean}, L{flatLocal}/L{hubeny}, L{flatPolar}, L{haversine}, L{thomas} and L{vincentys} and method L{Ellipsoid.distance2}.
@note: See note at function L{vincentys_}. ''' *unroll180(lon1, lon2, wrap=wrap))
'''Compute the I{angular} distance between two points using the U{spherical Law of Cosines <https://www.Movable-Type.co.UK/scripts/latlong.html#cosine-law>} formula.
@arg phi2: End latitude (C{radians}). @arg phi1: Start latitude (C{radians}). @arg lam21: Longitudinal delta, M{end-start} (C{radians}).
@return: Angular distance (C{radians}).
@see: Functions L{cosineLaw}, L{cosineAndoyerLambert_}, L{cosineForsytheAndoyerLambert_}, L{equirectangular_}, L{euclidean_}, L{flatLocal_}/L{hubeny_}, L{flatPolar_}, L{haversine_}, L{thomas_} and L{vincentys_}.
@note: See note at function L{vincentys_}. '''
'''(INTERNAL) Helper for ellipsoidal distances. ''' Phi_(lat1=lat1), radians(d_lon), datum=E)
'''(INTERNAL) Helper for spherical distances. ''' Phi_(lat1=lat1), radians(d_lon), **adjust)
'''(INTERNAL) Helper for distances. ''' earth if isinstance(earth, Datum) else _ellipsoidal_datum(earth, name=where.__name__)).ellipsoid # PYCHOK indent
'''Compute the distance between two points using the U{Equirectangular Approximation / Projection <https://www.Movable-Type.co.UK/scripts/latlong.html#equirectangular>}.
@arg lat1: Start latitude (C{degrees}). @arg lon1: Start longitude (C{degrees}). @arg lat2: End latitude (C{degrees}). @arg lon2: End longitude (C{degrees}). @kwarg radius: Mean earth radius, ellipsoid or datum (C{meter}, L{Ellipsoid}, L{Ellipsoid2}, L{Datum} or L{a_f2Tuple}). @kwarg options: Optional keyword arguments for function L{equirectangular_}.
@return: Distance (C{meter}, same units as B{C{radius}} or the ellipsoid or datum axes).
@raise TypeError: Invalid B{C{radius}}.
@see: Function L{equirectangular_} for more details, the available B{C{options}}, errors, restrictions and other, approximate or accurate distance functions. ''' Lat(lat2=lat2), Lon(lon2=lon2), **options).distance2) # PYCHOK 4 vs 2-3
adjust=True, limit=45, wrap=False): '''Compute the distance between two points using the U{Equirectangular Approximation / Projection <https://www.Movable-Type.co.UK/scripts/latlong.html#equirectangular>}.
This approximation is valid for short distance of several hundred Km or Miles, see the B{C{limit}} keyword argument and the L{LimitError}.
@arg lat1: Start latitude (C{degrees}). @arg lon1: Start longitude (C{degrees}). @arg lat2: End latitude (C{degrees}). @arg lon2: End longitude (C{degrees}). @kwarg adjust: Adjust the wrapped, unrolled longitudinal delta by the cosine of the mean latitude (C{bool}). @kwarg limit: Optional limit for lat- and longitudinal deltas (C{degrees}) or C{None} or C{0} for unlimited. @kwarg wrap: Wrap and L{pygeodesy.unroll180} longitudes (C{bool}).
@return: A L{Distance4Tuple}C{(distance2, delta_lat, delta_lon, unroll_lon2)}.
@raise LimitError: If the lat- and/or longitudinal delta exceeds the B{C{-limit..+limit}} range and L{pygeodesy.limiterrors} set to C{True}.
@see: U{Local, flat earth approximation <https://www.EdWilliams.org/avform.htm#flat>}, functions L{equirectangular}, L{cosineAndoyerLambert}, L{cosineForsytheAndoyerLambert}, L{cosineLaw}, L{euclidean}, L{flatLocal}/L{hubeny}, L{flatPolar}, L{haversine}, L{thomas} and L{vincentys} and methods L{Ellipsoid.distance2}, C{LatLon.distanceTo*} and C{LatLon.equirectangularTo}. '''
and max(abs(d_lat), abs(d_lon)) > limit > 0: lat1, lon1, lat2, lon2, limit=limit) raise LimitError('delta exceeds limit', txt=t)
'''Approximate the C{Euclidean} distance between two (spherical) points.
@arg lat1: Start latitude (C{degrees}). @arg lon1: Start longitude (C{degrees}). @arg lat2: End latitude (C{degrees}). @arg lon2: End longitude (C{degrees}). @kwarg radius: Mean earth radius, ellipsoid or datum (C{meter}, L{Ellipsoid}, L{Ellipsoid2}, L{Datum} or L{a_f2Tuple}). @kwarg adjust: Adjust the longitudinal delta by the cosine of the mean latitude (C{bool}). @kwarg wrap: Wrap and L{pygeodesy.unroll180} longitudes (C{bool}).
@return: Distance (C{meter}, same units as B{C{radius}} or the ellipsoid or datum axes).
@raise TypeError: Invalid B{C{radius}}.
@see: U{Distance between two (spherical) points <https://www.EdWilliams.org/avform.htm#Dist>}, functions L{euclid}, L{euclidean_}, L{cosineAndoyerLambert}, L{cosineForsytheAndoyerLambert}, L{cosineLaw}, L{equirectangular}, L{flatLocal}/L{hubeny}, L{flatPolar}, L{haversine}, L{thomas} and L{vincentys} and methods L{Ellipsoid.distance2}, C{LatLon.distanceTo*} and C{LatLon.equirectangularTo}. ''' *unroll180(lon1, lon2, wrap=wrap), adjust=adjust)
'''Approximate the I{angular} C{Euclidean} distance between two (spherical) points.
@arg phi2: End latitude (C{radians}). @arg phi1: Start latitude (C{radians}). @arg lam21: Longitudinal delta, M{end-start} (C{radians}). @kwarg adjust: Adjust the longitudinal delta by the cosine of the mean latitude (C{bool}).
@return: Angular distance (C{radians}).
@see: Functions L{euclid}, L{euclidean}, L{cosineAndoyerLambert_}, L{cosineForsytheAndoyerLambert_}, L{cosineLaw_}, L{equirectangular_}, L{flatLocal_}/L{hubeny_}, L{flatPolar_}, L{haversine_}, L{thomas_} and L{vincentys_}. '''
'''Compute the I{spherical excess} C{E} of a (spherical) triangle from two sides and the included angle.
@arg A: An interior triangle angle (C{radians}). @arg b: Frist adjacent triangle side (C{radians}). @arg c: Second adjacent triangle side (C{radians}).
@return: Spherical excess (C{radians}).
@raise UnitError: Invalid B{C{A}}, B{C{b}} or B{C{c}}.
@see: Function L{excessGirard}, L{excessLHuilier}, U{Spherical trigonometry<https://WikiPedia.org/wiki/Spherical_trigonometry>}. ''' Radians_(c=c) * _0_5)
'''Compute the I{spherical excess} C{E} of a (spherical) triangle using U{Girard's<https://MathWorld.Wolfram.com/GirardsSphericalExcessFormula.html>} formula.
@arg A: First interior triangle angle (C{radians}). @arg B: Second interior triangle angle (C{radians}). @arg C: Third interior triangle angle (C{radians}).
@return: Spherical excess (C{radians}).
@raise UnitError: Invalid B{C{A}}, B{C{B}} or B{C{C}}.
@see: Function L{excessLHuilier}, U{Spherical trigonometry <https://WikiPedia.org/wiki/Spherical_trigonometry>}. ''' Radians_(B=B), Radians_(C=C), -PI))
'''Compute the I{spherical excess} C{E} of a (spherical) triangle using U{L'Huilier's<https://MathWorld.Wolfram.com/LHuiliersTheorem.html>} Theorem.
@arg a: First triangle side (C{radians}). @arg b: Second triangle side (C{radians}). @arg c: Third triangle side (C{radians}).
@return: Spherical excess (C{radians}).
@raise UnitError: Invalid B{C{a}}, B{C{b}} or B{C{c}}.
@see: Function L{excessGirard}, U{Spherical trigonometry <https://WikiPedia.org/wiki/Spherical_trigonometry>}. '''
'''Compute the surface area of a (spherical) quadrilateral bounded by a segment of a great circle, two meridians and the equator using U{Karney's <https://MathOverflow.net/questions/97711/the-area-of-spherical-polygons>} method.
@arg lat1: Start latitude (C{degrees}). @arg lon1: Start longitude (C{degrees}). @arg lat2: End latitude (C{degrees}). @arg lon2: End longitude (C{degrees}). @kwarg radius: Mean earth radius, ellipsoid or datum (C{meter}, L{Ellipsoid}, L{Ellipsoid2}, L{Datum} or L{a_f2Tuple}) or C{None}. @kwarg wrap: Wrap and L{pygeodesy.unroll180} longitudes (C{bool}).
@return: Surface area, I{signed} (I{square} C{meter}, or units of B{C{radius}} I{squared}) or I{spherical excess} (C{radians}) if B{C{radius}} is C{None} or C{0}.
@raise TypeError: Invalid B{C{radius}}.
@raise UnitError: Invalid B{C{lat2}} or B{C{lat1}}.
@raise ValueError: Semi-circular longitudinal delta.
@see: Function L{excessKarney_} and L{excessQuad}. ''' *unroll180(lon1, lon2, wrap=wrap))
'''Compute the I{spherical excess} C{E} of a (spherical) quadrilateral bounded by a segment of a great circle, two meridians and the equator using U{Karney's <https://MathOverflow.net/questions/97711/the-area-of-spherical-polygons>} method.
@arg phi2: End latitude (C{radians}). @arg phi1: Start latitude (C{radians}). @arg lam21: Longitudinal delta, M{end-start} (C{radians}).
@return: Spherical excess, I{signed} (C{radians}).
@raise ValueError: Semi-circular longitudinal delta B{C{lam21}}.
@see: Function L{excessKarney}, U{Area of a spherical polygon <https://MathOverflow.net/questions/97711/the-area-of-spherical-polygons>}. ''' # from: Veness <https://www.Movable-Type.co.UK/scripts/latlong.html> Area # method due to Karney: for each edge of the polygon, # # tan(Δλ / 2) · (tan(φ1 / 2) + tan(φ2 / 2)) # tan(E / 2) = ----------------------------------------- # 1 + tan(φ1 / 2) · tan(φ2 / 2) # # where E is the spherical excess of the trapezium obtained by extending # the edge to the equator-circle vector for each edge (see also ***). _1_0 + (t1 * t2)) * _2_0)
# ***) Original post no longer available, following is a copy of the main part # <http://OSGeo-org.1560.x6.Nabble.com/Area-of-a-spherical-polygon-td3841625.html> # # The area of a polygon on a (unit) sphere is given by the spherical excess # # A = 2 * pi - sum(exterior angles) # # However this is badly conditioned if the polygon is small. In this case, use # # A = sum S12{i, i+1} over the edges of the polygon # # where S12 is the area of the quadrilateral bounded by an edge of the polygon, # two meridians and the equator, i.e. with vertices (phi1, lambda1), (phi2, # lambda2), (0, lambda1) and (0, lambda2). S12 is given by # # tan(S12 / 2) = tan(lambda21 / 2) * (tan(phi1 / 2) + tan(phi2 / 2)) / # (tan(phi1 / 2) * tan(phi2 / 2) + 1) # # = tan(lambda21 / 2) * tanh((lam(phi1) + lam(phi2)) / 2) # # where lambda21 = lambda2 - lambda1 and lam(x) is the Lambertian (or inverse # Gudermannian) function # # lam(x) = asinh(tan(x)) = atanh(sin(x)) = 2 * atanh(tan(x / 2)) # # Notes: The formula for S12 is exact, except that... # - it is indeterminate if an edge is a semi-circle # - the formula for A applies only if the polygon does not include a pole # (if it does, then add +/- 2 * pi to the result) # - in the limit of small phi and lambda, S12 reduces to the trapezoidal # formula, S12 = (lambda2 - lambda1) * (phi1 + phi2) / 2 # - I derived this result from the equation for the area of a spherical # triangle in terms of two edges and the included angle given by, e.g. # U{Todhunter, I. - Spherical Trigonometry (1871), Sec. 103, Eq. (2) # <http://Books.Google.com/books?id=3uBHAAAAIAAJ&pg=PA71>} # - I would be interested to know if this formula for S12 is already known # - Charles Karney
'''Compute the surface area of a (spherical) quadrilateral bounded by a segment of a great circle, two meridians and the equator.
@arg lat1: Start latitude (C{degrees}). @arg lon1: Start longitude (C{degrees}). @arg lat2: End latitude (C{degrees}). @arg lon2: End longitude (C{degrees}). @kwarg radius: Mean earth radius, ellipsoid or datum (C{meter}, L{Ellipsoid}, L{Ellipsoid2}, L{Datum} or L{a_f2Tuple}) or C{None}. @kwarg wrap: Wrap and L{pygeodesy.unroll180} longitudes (C{bool}).
@return: Surface area, I{signed} (I{square} C{meter}, or units of B{C{radius}} I{squared}) or I{spherical excess} (C{radians}) if B{C{radius}} is C{None} or C{0}.
@raise TypeError: Invalid B{C{radius}}.
@raise UnitError: Invalid B{C{lat2}} or B{C{lat1}}.
@see: Function L{excessQuad_} and L{excessKarney}. ''' *unroll180(lon1, lon2, wrap=wrap))
'''Compute the I{spherical excess} C{E} of a (spherical) quadrilateral bounded by a segment of a great circle, two meridians and the equator.
@arg phi2: End latitude (C{radians}). @arg phi1: Start latitude (C{radians}). @arg lam21: Longitudinal delta, M{end-start} (C{radians}).
@return: Spherical excess, I{signed} (C{radians}).
@see: Function L{excessQuad}, U{Spherical trigonometry <https://WikiPedia.org/wiki/Spherical_trigonometry>}. '''
'''Compute the distance between two (ellipsoidal) points using the U{ellipsoidal Earth to plane projection<https://WikiPedia.org/ wiki/Geographical_distance#Ellipsoidal_Earth_projected_to_a_plane>} aka U{Hubeny<https://www.OVG.AT/de/vgi/files/pdf/3781/>} formula.
@arg lat1: Start latitude (C{degrees}). @arg lon1: Start longitude (C{degrees}). @arg lat2: End latitude (C{degrees}). @arg lon2: End longitude (C{degrees}). @kwarg datum: Datum or ellipsoid to use (L{Datum}, L{Ellipsoid}, L{Ellipsoid2} or L{a_f2Tuple}). @kwarg wrap: Wrap and L{pygeodesy.unroll180} longitudes (C{bool}).
@return: Distance (C{meter}, same units as the B{C{datum}}'s ellipsoid axes).
@raise TypeError: Invalid B{C{datum}}.
@note: The meridional and prime_vertical radii of curvature are taken and scaled at the mean of both latitude.
@see: Functions L{flatLocal_} or L{hubeny_}, L{cosineLaw}, L{flatPolar}, L{cosineAndoyerLambert}, L{cosineForsytheAndoyerLambert}, L{equirectangular}, L{euclidean}, L{haversine}, L{thomas}, L{vincentys}, method L{Ellipsoid.distance2} and U{local, flat earth approximation<https://www.EdWilliams.org/avform.htm#flat>}. ''' Phi_(lat1=lat1), radians(d), datum=datum)
'''Compute the I{angular} distance between two (ellipsoidal) points using the U{ellipsoidal Earth to plane projection<https://WikiPedia.org/ wiki/Geographical_distance#Ellipsoidal_Earth_projected_to_a_plane>} aka U{Hubeny<https://www.OVG.AT/de/vgi/files/pdf/3781/>} formula.
@arg phi2: End latitude (C{radians}). @arg phi1: Start latitude (C{radians}). @arg lam21: Longitudinal delta, M{end-start} (C{radians}). @kwarg datum: Datum or ellipsoid to use (L{Datum}, L{Ellipsoid}, L{Ellipsoid2} or L{a_f2Tuple}).
@return: Angular distance (C{radians}).
@raise TypeError: Invalid B{C{datum}}.
@note: The meridional and prime_vertical radii of curvature are taken and scaled I{at the mean of both latitude}.
@see: Functions L{flatLocal} or L{hubeny}, L{cosineAndoyerLambert_}, L{cosineForsytheAndoyerLambert_}, L{cosineLaw_}, L{flatPolar_}, L{equirectangular_}, L{euclidean_}, L{haversine_}, L{thomas_} and L{vincentys_} and U{local, flat earth approximation <https://www.EdWilliams.org/avform.htm#flat>}. '''
'''Compute the distance between two (spherical) points using the U{polar coordinate flat-Earth <https://WikiPedia.org/wiki/ Geographical_distance#Polar_coordinate_flat-Earth_formula>} formula.
@arg lat1: Start latitude (C{degrees}). @arg lon1: Start longitude (C{degrees}). @arg lat2: End latitude (C{degrees}). @arg lon2: End longitude (C{degrees}). @kwarg radius: Mean earth radius, ellipsoid or datum (C{meter}, L{Ellipsoid}, L{Ellipsoid2}, L{Datum} or L{a_f2Tuple}). @kwarg wrap: Wrap and L{pygeodesy.unroll180} longitudes (C{bool}).
@return: Distance (C{meter}, same units as B{C{radius}} or the ellipsoid or datum axes).
@raise TypeError: Invalid B{C{radius}}.
@see: Functions L{flatPolar_}, L{cosineAndoyerLambert}, L{cosineForsytheAndoyerLambert},L{cosineLaw}, L{flatLocal}/L{hubeny}, L{equirectangular}, L{euclidean}, L{haversine}, L{thomas} and L{vincentys}. ''' *unroll180(lon1, lon2, wrap=wrap))
'''Compute the I{angular} distance between two (spherical) points using the U{polar coordinate flat-Earth<https://WikiPedia.org/wiki/ Geographical_distance#Polar_coordinate_flat-Earth_formula>} formula.
@arg phi2: End latitude (C{radians}). @arg phi1: Start latitude (C{radians}). @arg lam21: Longitudinal delta, M{end-start} (C{radians}).
@return: Angular distance (C{radians}).
@see: Functions L{flatPolar}, L{cosineAndoyerLambert_}, L{cosineForsytheAndoyerLambert_}, L{cosineLaw_}, L{equirectangular_}, L{euclidean_}, L{flatLocal_}/L{hubeny_}, L{haversine_}, L{thomas_} and L{vincentys_}. ''' a = _0_0
'''Compute the intersection of a Line-Of-Sight from a Point-Of-View in space with the surface of the earth.
@arg pov: Point-Of-View outside the earth (C{Cartesian}, L{Ecef9Tuple} or L{Vector3d}). @kwarg los: Line-Of-Sight, I{direction} to earth (L{Vector3d}) or C{None} to point to the earth' center. @kwarg earth: The earth model (L{Datum}, L{Ellipsoid}, L{Ellipsoid2}, L{a_f2Tuple} or C{scalar} radius in C{meter}). @kwarg LatLon_and_kwds: Optional C{LatLon} class for the intersection point plus C{LatLon} keyword arguments, include B{C{datum}} if different from B{C{earth}}.
@return: The earth intersection (L{Vector3d}, C{Cartesian type} of B{C{pov}} or B{C{LatLon}}).
@raise IntersectionError: Null B{C{pov}} or B{C{los}} vector, B{C{pov}} is inside the earth or B{C{los}} points outside the earth or points in an opposite direction.
@raise TypeError: Invalid B{C{pov}}, B{C{los}} or B{C{earth}}.
@see: U{Stephen Hartzell<https://StephenHartzell.medium.com/ satellite-line-of-sight-intersection-with-earth-d786b4a6a9b6>} and function L{pygeodesy.tyr3d} for B{C{los}}. ''' _spherical_datum(earth, name=hartzell.__name__)
raise IntersectionError(pov=pov, los=los, earth=earth, txt=_near_(_null_))
# a2 and b2 factored out, b2 == a2 and b2 / a2 == 1 c2 * u2, -u2 * z2, -w2 * x2, ux * wz * 2, -w2 * y2, -u2 * y2 * q2, -v2 * x2 * q2, ux * vy * 2 * q2) elif r < 0: # LOS pointing away from or missing the earth t = _opposite_ if max(ux, vy, wz) > 0 else _outside_ raise IntersectionError(pov=pov, los=los, earth=earth, txt=t)
t = _inside_ if min(x2 - a2, y2 - b2, z2 - c2) < EPS else _outside_ raise IntersectionError(pov=pov, los=los, earth=earth, txt=t)
raise IntersectionError(pov=pov, los=los, earth=earth, txt=_too_(_distant_))
'''Compute the distance between two (spherical) points using the U{Haversine<https://www.Movable-Type.co.UK/scripts/latlong.html>} formula.
@arg lat1: Start latitude (C{degrees}). @arg lon1: Start longitude (C{degrees}). @arg lat2: End latitude (C{degrees}). @arg lon2: End longitude (C{degrees}). @kwarg radius: Mean earth radius, ellipsoid or datum (C{meter}, L{Ellipsoid}, L{Ellipsoid2}, L{Datum} or L{a_f2Tuple}). @kwarg wrap: Wrap and L{pygeodesy.unroll180} longitudes (C{bool}).
@return: Distance (C{meter}, same units as B{C{radius}}).
@raise TypeError: Invalid B{C{radius}}.
@see: U{Distance between two (spherical) points <https://www.EdWilliams.org/avform.htm#Dist>}, functions L{cosineLaw}, L{cosineAndoyerLambert}, L{cosineForsytheAndoyerLambert}, L{equirectangular}, L{euclidean}, L{flatLocal}/L{hubeny}, L{flatPolar}, L{thomas} and L{vincentys} and methods L{Ellipsoid.distance2}, C{LatLon.distanceTo*} and C{LatLon.equirectangularTo}.
@note: See note at function L{vincentys_}. ''' *unroll180(lon1, lon2, wrap=wrap))
'''Compute the I{angular} distance between two (spherical) points using the U{Haversine<https://www.Movable-Type.co.UK/scripts/latlong.html>} formula.
@arg phi2: End latitude (C{radians}). @arg phi1: Start latitude (C{radians}). @arg lam21: Longitudinal delta, M{end-start} (C{radians}).
@return: Angular distance (C{radians}).
@see: Functions L{haversine}, L{cosineAndoyerLambert_}, L{cosineForsytheAndoyerLambert_}, L{cosineLaw_}, L{equirectangular_}, L{euclidean_}, L{flatLocal_}/L{hubeny_}, L{flatPolar_}, L{thomas_} and L{vincentys_}.
@note: See note at function L{vincentys_}. '''
'''Determine the height above the (spherical) earth' surface after traveling along a straight line at a given tilt.
@arg angle: Tilt angle above horizontal (C{degrees}). @arg distance: Distance along the line (C{meter} or same units as B{C{radius}}). @kwarg radius: Optional mean earth radius (C{meter}).
@return: Height (C{meter}, same units as B{C{distance}} and B{C{radius}}).
@raise ValueError: Invalid B{C{angle}}, B{C{distance}} or B{C{radius}}.
@see: U{MultiDop geog_lib.GeogBeamHt<https://GitHub.com/NASA/MultiDop>} (U{Shapiro et al. 2009, JTECH <https://Journals.AMetSoc.org/doi/abs/10.1175/2009JTECHA1256.1>} and U{Potvin et al. 2012, JTECH <https://Journals.AMetSoc.org/doi/abs/10.1175/JTECH-D-11-00019.1>}). ''' d, h = h, d
raise _ValueError(angle=angle, distance=distance, radius=radius)
'''Determine the distance to the horizon from a given altitude above the (spherical) earth.
@arg height: Altitude (C{meter} or same units as B{C{radius}}). @kwarg radius: Optional mean earth radius (C{meter}). @kwarg refraction: Consider atmospheric refraction (C{bool}).
@return: Distance (C{meter}, same units as B{C{height}} and B{C{radius}}).
@raise ValueError: Invalid B{C{height}} or B{C{radius}}.
@see: U{Distance to horizon<https://www.EdWilliams.org/avform.htm#Horizon>}. ''' raise _ValueError(height=height, radius=radius)
else:
lat2, lon2, radius2, datum=None, wrap=True): '''Conveniently compute the intersections of two circles each defined by a (geodetic) center point and a radius, using either ...
1) L{vector3d.intersections2} for small distances (below 100 KM or about 0.9 degrees) or if no B{C{datum}} is specified, or ...
2) L{sphericalTrigonometry.intersections2} for a spherical B{C{datum}} or if B{C{datum}} is a C{scalar} representing the earth radius, conventionally in C{meter} or ...
3) L{ellipsoidalKarney.intersections2} for an ellipsoidal B{C{datum}} and if I{Karney}'s U{geographiclib<https://PyPI.org/project/geographiclib>} is installed, otherwise ...
4) L{ellipsoidalExact.intersections2}, provided B{C{datum}} is ellipsoidal.
@arg lat1: Latitude of the first circle center (C{degrees}). @arg lon1: Longitude of the first circle center (C{degrees}). @arg radius1: Radius of the first circle (C{meter}, conventionally). @arg lat2: Latitude of the second circle center (C{degrees}). @arg lon2: Longitude of the second circle center (C{degrees}). @arg radius2: Radius of the second circle (C{meter}, same units as B{C{radius1}}). @kwarg datum: Optional ellipsoidal or spherical datum (L{Datum}, L{Ellipsoid}, L{Ellipsoid2}, L{a_f2Tuple} or C{scalar} earth radius in C{meter}, same units as B{C{radius1}} and B{C{radius2}}) or C{None}. @kwarg wrap: Wrap and unroll longitudes (C{bool}).
@return: 2-Tuple of the intersection points, each a L{LatLon2Tuple}C{(lat, lon)}. For abutting circles, the points are the same instance, aka I{radical center}.
@raise IntersectionError: Concentric, antipodal, invalid or non-intersecting circles or no convergence.
@raise TypeError: Invalid B{C{datum}}.
@raise UnitError: Invalid B{C{lat1}}, B{C{lon1}}, B{C{radius1}} B{C{lat2}}, B{C{lon2}} or B{C{radius2}}. ''' adjust=True, wrap=wrap) < _100km:
m.Vector3d(lon2, lat2, 0), r2, sphere=False, Vector=_V2T) else:
pass except ImportError: m = _MODS.ellipsoidalExact else: raise _AssertionError(datum=d)
m.LatLon(lat2, lon2, datum=d), radius2, LatLon=_LL2T, height=0, wrap=wrap)
'''Check whether two points are antipodal, on diametrically opposite sides of the earth.
@arg lat1: Latitude of one point (C{degrees}). @arg lon1: Longitude of one point (C{degrees}). @arg lat2: Latitude of the other point (C{degrees}). @arg lon2: Longitude of the other point (C{degrees}). @kwarg eps: Tolerance for near-equality (C{degrees}).
@return: C{True} if points are antipodal within the B{C{eps}} tolerance, C{False} otherwise.
@see: U{Geosphere<https://CRAN.R-Project.org/web/packages/geosphere/geosphere.pdf>}. ''' abs(abs(wrap180(lon1) - wrap180(lon2)) % _360_0 - _180_0) < eps
'''Check whether two points are antipodal, on diametrically opposite sides of the earth.
@arg phi1: Latitude of one point (C{radians}). @arg lam1: Longitude of one point (C{radians}). @arg phi2: Latitude of the other point (C{radians}). @arg lam2: Longitude of the other point (C{radians}). @kwarg eps: Tolerance for near-equality (C{radians}).
@return: C{True} if points are antipodal within the B{C{eps}} tolerance, C{False} otherwise.
@see: U{Geosphere<https://CRAN.R-Project.org/web/packages/geosphere/geosphere.pdf>}. ''' abs(wrapPI(lam1) - wrapPI(lam2)) % PI2 - PI) < eps
'''Convert lat-, longitude to C{n-vector} (normal to the earth's surface) X, Y and Z components.
@arg lat: Latitude (C{degrees}). @arg lon: Longitude (C{degrees}). @kwarg name: Optional name (C{str}).
@return: A L{Vector3Tuple}C{(x, y, z)}.
@see: Function L{philam2n_xyz}.
@note: These are C{n-vector} x, y and z components, I{NOT} geocentric ECEF x, y and z coordinates! '''
'''Convert C{n-vector} components to lat- and longitude in C{degrees}.
@arg x: X component (C{scalar}). @arg y: Y component (C{scalar}). @arg z: Z component (C{scalar}). @kwarg name: Optional name (C{str}).
@return: A L{LatLon2Tuple}C{(lat, lon)}.
@see: Function L{n_xyz2philam}. '''
'''Convert C{n-vector} components to lat- and longitude in C{radians}.
@arg x: X component (C{scalar}). @arg y: Y component (C{scalar}). @arg z: Z component (C{scalar}). @kwarg name: Optional name (C{str}).
@return: A L{PhiLam2Tuple}C{(phi, lam)}.
@see: Function L{n_xyz2latlon}. '''
'''Compare the direction of two bearings given in C{degrees}.
@arg bearing1: First bearing (compass C{degrees}). @arg bearing2: Second bearing (compass C{degrees}). @kwarg margin: Optional, interior angle bracket (C{degrees}), default C{90}.
@return: C{True} if both bearings point in opposite, C{False} if in similar or C{None} if in perpendicular directions.
@see: Function L{opposing_}. ''' True if (_180_0 - m) < d < (_180_0 + m) else None)
'''Compare the direction of two bearings given in C{radians}.
@arg radians1: First bearing (C{radians}). @arg radians2: Second bearing (C{radians}). @kwarg margin: Optional, interior angle bracket (C{radians}), default C{PI_2}.
@return: C{True} if both bearings point in opposite, C{False} if in similar or C{None} if in perpendicular directions.
@see: Function L{opposing}. ''' True if (PI - m) < r < (PI + m) else None)
'''Convert lat-, longitude to C{n-vector} (normal to the earth's surface) X, Y and Z components.
@arg phi: Latitude (C{radians}). @arg lam: Longitude (C{radians}). @kwarg name: Optional name (C{str}).
@return: A L{Vector3Tuple}C{(x, y, z)}.
@see: Function L{latlon2n_xyz}.
@note: These are C{n-vector} x, y and z components, I{NOT} geocentric ECEF x, y and z coordinates! ''' # Kenneth Gade eqn 3, but using right-handed # vector x -> 0°E,0°N, y -> 90°E,0°N, z -> 90°N
# (INTERNAL) See C{radical2} below # assert d > EPS0
'''Compute the I{radical ratio} and I{radical line} of two U{intersecting circles<https://MathWorld.Wolfram.com/ Circle-CircleIntersection.html>}.
The I{radical line} is perpendicular to the axis thru the centers of the circles at C{(0, 0)} and C{(B{distance}, 0)}.
@arg distance: Distance between the circle centers (C{scalar}). @arg radius1: Radius of the first circle (C{scalar}). @arg radius2: Radius of the second circle (C{scalar}).
@return: A L{Radical2Tuple}C{(ratio, xline)} where C{0.0 <= ratio <= 1.0} and C{xline} is along the B{C{distance}}.
@raise IntersectionError: The B{C{distance}} exceeds the sum of B{C{radius1}} and B{C{radius2}}.
@raise UnitError: Invalid B{C{distance}}, B{C{radius1}} or B{C{radius2}}.
@see: U{Circle-Circle Intersection <https://MathWorld.Wolfram.com/Circle-CircleIntersection.html>}. ''' raise IntersectionError(distance=d, radius1=r1, radius2=r2, txt=_too_(_distant_)) Radical2Tuple(_0_5, _0_0)
'''2-Tuple C{(ratio, xline)} of the I{radical} C{ratio} and I{radical} C{xline}, both C{scalar} and C{0.0 <= ratio <= 1.0} '''
# scale factor cos(mean of lats) for delta lon
# scale factor cos(mean of phis) for delta lam
'''(INTERNAL) C{sin}es, C{cos}ines and C{acos}ine. '''
'''Compute the distance between two (ellipsoidal) points using U{Thomas'<https://apps.DTIC.mil/dtic/tr/fulltext/u2/703541.pdf>} formula.
@arg lat1: Start latitude (C{degrees}). @arg lon1: Start longitude (C{degrees}). @arg lat2: End latitude (C{degrees}). @arg lon2: End longitude (C{degrees}). @kwarg datum: Datum or ellipsoid to use (L{Datum}, L{Ellipsoid}, L{Ellipsoid2} or L{a_f2Tuple}). @kwarg wrap: Wrap and L{pygeodesy.unroll180} longitudes (C{bool}).
@return: Distance (C{meter}, same units as the B{C{datum}}'s ellipsoid axes).
@raise TypeError: Invalid B{C{datum}}.
@see: Functions L{thomas_}, L{cosineAndoyerLambert}, L{cosineForsytheAndoyerLambert}, L{cosineLaw}, L{equirectangular}, L{euclidean}, L{flatLocal}/L{hubeny}, L{flatPolar}, L{haversine}, L{vincentys} and method L{Ellipsoid.distance2}. ''' *unroll180(lon1, lon2, wrap=wrap))
'''Compute the I{angular} distance between two (ellipsoidal) points using U{Thomas'<https://apps.DTIC.mil/dtic/tr/fulltext/u2/703541.pdf>} formula.
@arg phi2: End latitude (C{radians}). @arg phi1: Start latitude (C{radians}). @arg lam21: Longitudinal delta, M{end-start} (C{radians}). @kwarg datum: Datum or ellipsoid to use (L{Datum}, L{Ellipsoid}, L{Ellipsoid2} or L{a_f2Tuple}).
@return: Angular distance (C{radians}).
@raise TypeError: Invalid B{C{datum}}.
@see: Functions L{thomas}, L{cosineAndoyerLambert_}, L{cosineForsytheAndoyerLambert_}, L{cosineLaw_}, L{equirectangular_}, L{euclidean_}, L{flatLocal_}/L{hubeny_}, L{flatPolar_}, L{haversine_} and L{vincentys_} and U{Geodesy-PHP <https://GitHub.com/jtejido/geodesy-php/blob/master/src/Geodesy/ Distance/ThomasFormula.php>}. '''
'''Compute the distance between two (spherical) points using U{Vincenty's<https://WikiPedia.org/wiki/Great-circle_distance>} spherical formula.
@arg lat1: Start latitude (C{degrees}). @arg lon1: Start longitude (C{degrees}). @arg lat2: End latitude (C{degrees}). @arg lon2: End longitude (C{degrees}). @kwarg radius: Mean earth radius, ellipsoid or datum (C{meter}, L{Ellipsoid}, L{Ellipsoid2}, L{Datum} or L{a_f2Tuple}). @kwarg wrap: Wrap and L{pygeodesy.unroll180} longitudes (C{bool}).
@return: Distance (C{meter}, same units as B{C{radius}}).
@raise UnitError: Invalid B{C{radius}}.
@see: Functions L{vincentys_}, L{cosineAndoyerLambert}, L{cosineForsytheAndoyerLambert},L{cosineLaw}, L{equirectangular}, L{euclidean}, L{flatLocal}/L{hubeny}, L{flatPolar}, L{haversine} and L{thomas} and methods L{Ellipsoid.distance2}, C{LatLon.distanceTo*} and C{LatLon.equirectangularTo}.
@note: See note at function L{vincentys_}. ''' *unroll180(lon1, lon2, wrap=wrap))
'''Compute the I{angular} distance between two (spherical) points using U{Vincenty's<https://WikiPedia.org/wiki/Great-circle_distance>} spherical formula.
@arg phi2: End latitude (C{radians}). @arg phi1: Start latitude (C{radians}). @arg lam21: Longitudinal delta, M{end-start} (C{radians}).
@return: Angular distance (C{radians}).
@see: Functions L{vincentys}, L{cosineAndoyerLambert_}, L{cosineForsytheAndoyerLambert_}, L{cosineLaw_}, L{equirectangular_}, L{euclidean_}, L{flatLocal_}/L{hubeny_}, L{flatPolar_}, L{haversine_} and L{thomas_}.
@note: Functions L{vincentys_}, L{haversine_} and L{cosineLaw_} produce equivalent results, but L{vincentys_} is suitable for antipodal points and slightly more expensive (M{3 cos, 3 sin, 1 hypot, 1 atan2, 6 mul, 2 add}) than L{haversine_} (M{2 cos, 2 sin, 2 sqrt, 1 atan2, 5 mul, 1 add}) and L{cosineLaw_} (M{3 cos, 3 sin, 1 acos, 3 mul, 1 add}). '''
# **) 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. |