Coverage for pygeodesy/albers.py : 95%

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 -*-
Classes L{AlbersEqualArea}, L{AlbersEqualArea2}, L{AlbersEqualArea4}, L{AlbersEqualAreaCylindrical}, L{AlbersEqualAreaNorth}, L{AlbersEqualAreaSouth} and L{AlbersError}, a transcription of I{Charles Karney}'s C++ class U{AlbersEqualArea <https://GeographicLib.SourceForge.io/html/classGeographicLib_1_1AlbersEqualArea.html>}.
See also I{Albers Equal-Area Conic Projection} in U{John P. Snyder, "Map Projections -- A Working Manual", 1987<https://Pubs.USGS.gov/pp/1395/report.pdf>}, pp 98-106 and the Albers Conical Equal-Area examples on pp 291-294. ''' # make sure int/int division yields float quotient, see .basics
_EPS__2, _EPS__4, _gamma_, _lat_, \ _lat1_, _lat2_, _lon_, _no_, _scale_, \ _x_, _y_, _0_0, _0_5, _1_0, _2_0, \ _3_0, _90_0 property_RO Meter, Scalar_
'''An L{AlbersEqualArea}, L{AlbersEqualArea2}, L{AlbersEqualArea4}, L{AlbersEqualAreaCylindrical}, L{AlbersEqualAreaNorth}, L{AlbersEqualAreaSouth} or L{Albers7Tuple} issue. '''
'''(INTERNAL) Divided differences, defined as M{Df(x, y) = (f(x) - f(y)) / (x - y)} with M{sn(x) = x / sqrt(1 + x^2)}: M{Dsn(x, y) = (x + y) / ((sn(x) + sn(y)) * (1 + x^2) * (1 + y^2))}.
@see: U{W. M. Kahan and R. J. Fateman, "Sympbolic Computation of Divided Differences"<https://People.EECS.Berkeley.EDU/~fateman/papers/divdiff.pdf>}, U{ACM SIGSAM Bulletin 33(2), 7-28 (1999)<https://DOI.org/10.1145/334714.334716>} and U{AlbersEqualArea.hpp <https://GeographicLib.sourceforge.io/html/AlbersEqualArea_8hpp_source.html>}. ''' # sx = x / hypot1(x) else: t = x - y if t: d = (sx - sy) / t
'''(INTERNAL) Scale C{B{k} >= EPS0}. '''
'''(INTERNAL) Latitude C{-90 <= B{lat} <= 90}. '''
'''(INTERNAL) Longitude C{-180 <= B{lon} <= 180}. ''' return Lon_(Error=AlbersError, **name_lon)
'''(INTERNAL) Converge tolerance. '''
'''(INTERNAL) Base class for C{AlbersEqualArea...} projections.
@see: I{Karney}'s C++ class U{AlbersEqualArea<https://GeographicLib.SourceForge.io/ html/classGeographicLib_1_1AlbersEqualArea.html>}, method C{Init}. '''
'''(INTERNAL) New C{AlbersEqualArea...} instance. '''
raise AlbersError(clat1=ca1, clat2=ca2) # determine hemisphere of tangent latitude # internally, tangent latitude positive raise AlbersError(slat1=sa1, slat2=sa2) # avoid singularities at poles
else:
# sa0m = 1 - sa0 = 1 / (sec(a0) * (tan(a0) + sec(a0))) - sa0m**2 * e2 * (2 + (_1_0 + e2) * sa0) / (e12 * sa02_) # == B + A / fsum_(s1_qZ * dg, -dD, dg * BA, g * dAB) # == u/du else: raise AlbersError(iteration=_NUMIT0, txt=_no_(Fmt.convergence(tol)))
else:
'''(INTERNAL) Compute C{sm1 / (s / qZ)} and C{C} for .__init__. '''
* (_1_0 - e2 * sa2**2)
(ca2, sa2, ta2, scb22))
'''(INTERNAL) Sum of C{sm1} terms and C{sin(xi)}s for ._s1_qZ_C2. ''' else: sa = _1_0 - sa a += (_1_0 - sxi) / sa b += scb * sa
'''Get the datum (L{Datum}). '''
'''Get the geodesic's equatorial radius, semi-axis (C{meter}). '''
'''Get the iteration number (C{int}) or C{None} if not available/applicable. '''
'''Get the geodesic's flattening (C{float}). '''
'''Convert a geodetic location to east- and northing.
@arg lat: Latitude of the location (C{degrees}). @arg lon: Longitude of the location (C{degrees}). @kwarg lon0: Optional central meridian longitude (C{degrees}). @kwarg name: Optional name for the location (C{str}).
@return: An L{Albers7Tuple}C{(x, y, lat, lon, gamma, scale, datum)}.
@note: The origin latitude is returned by C{property lat0}. No false easting or northing is added. The value of B{C{lat}} should be in the range C{[-90..90] degrees}. The returned values C{x} and C{y} will be large but finite for points projecting to infinity, i.e. one or both of the poles. '''
else: x = self._k02 * b y = _0_0
name=name or self.name)
'''Is this projection polar (C{bool})? '''
'''Get the latitude of the projection origin (C{degrees}).
This is the latitude of minimum azimuthal scale and equals the B{C{lat}} in the 1-parallel L{AlbersEqualArea} and lies between B{C{lat1}} and B{C{lat2}} for the 2-parallel L{AlbersEqualArea2} and L{AlbersEqualArea4} projections. '''
'''Get the latitude of the first parallel (C{degrees}). '''
'''Get the latitude of the second parallel (C{degrees}).
@note: The second and first parallel latitudes are the same instance for 1-parallel C{AlbersEqualArea*} projections. '''
def majoradius(self): # PYCHOK no cover '''DEPRECATED, use property C{equatoradius}.''' return self.equatoradius
'''Set the azimuthal scale for this projection.
@arg lat: Northern latitude (C{degrees}). @arg k: Azimuthal scale at latitude B{C{lat}} (C{scalar}).
@raise AlbersError: Invalid B{C{lat}} or B{C{k}}.
@note: This allows a I{latitude of conformality} to be specified. ''' k0 = _Ks(k=k) / self.forward(lat, _0_0).scale self._update(k0 != self._k0) self._k0_(k0)
'''Convert an east- and northing location to geodetic lat- and longitude.
@arg x: Easting of the location (C{meter}). @arg y: Northing of the location (C{meter}). @kwarg lon0: Optional central meridian longitude (C{degrees}). @kwarg name: Optional name for the location (C{str}). @kwarg LatLon: Class to use (C{LatLon}) or C{None}. @kwarg LatLon_kwds: Optional, additional B{C{LatLon}} keyword arguments, ignored if C{B{LatLon}=None}.
@return: The geodetic (C{LatLon}) or if B{C{LatLon}} is C{None} an L{Albers7Tuple}C{(x, y, lat, lon, gamma, scale, datum)}.
@note: The origin latitude is returned by C{property lat0}. No false easting or northing is added. The returned value of C{lon} is in the range C{[-180..180] degrees} and C{lat} is in the range C{[-90..90] degrees}. If the given B{C{x}} or B{C{y}} point is outside the valid projected space the nearest pole is returned. '''
# dsxia = scxi0 * dsxi
name=name or self.name) else: kwds = _xkwds(LatLon_kwds, datum=self.datum) r = self._xnamed(LatLon(lat, lon, **kwds), name=name)
'''Get the central scale for the projection (C{float}).
This is the azimuthal scale on the latitude of origin of the projection, see C{property lat0}. '''
'''(INTERNAL) Function M{atanhee(x)}, defined as ... atanh( E.e * x) / E.e if f > 0 # oblate atan (sqrt(-E.e2) * x) / sqrt(-E.e2) if f < 0 # prolate x if f = 0. '''
'''(INTERNAL) Function M{atanh(sqrt(x)) / sqrt(x) - 1}. ''' # x < E.e^2 = 2 * E.f use ... # x / 3 + x^2 / 5 + x^3 / 7 + ... else: s = sqrt(s) s = (atanh(s) if x > 0 else atan(s)) / s - _1_0
'''(INTERNAL) Compute the azimuthal scale C{_Ks(k0=k0)}. '''
'''(INTERNAL) Get 3-tuple C{(cos, sin, tan)} of M{xi(ta)}. '''
'''(INTERNAL) Function M{Datanhee(x, y)}, defined as M{atanhee((x - y) / (1 - E.e^2 * x * y)) / (x - y)}. ''' else: raise AlbersError(x=x, y=y, txt=_AlbersBase._Datanhee.__name__)
'''(INTERNAL) Function M{D2atanhee(x, y)}, defined as M{(Datanhee(1, y) - Datanhee(1, x)) / (y - x)}. '''
elif (_1_0 - x): s = (self._Datanhee(_1_0, y) - self._Datanhee(x, y)) / (_1_0 - x) else: raise AlbersError(x=x, y=y, txt=_AlbersBase._D2atanhee.__name__)
'''(INTERNAL) Set C{._k0}, C{._k02}, etc. '''
'''(INTERNAL) Function M{tan-phi from tan-xi}. '''
# dtxi/dta = (scxi / sca)^3 * 2 * (1 - e^2) / (qZ * (1 - e^2 * sa^2)^2) raise AlbersError(iteration=_NUMIT, txt=_no_(Fmt.convergence(tol)))
'''(INTERNAL) Function M{tan-xi from tan-phi}. '''
* (es1m1 / es2m1a + self._atanhee(es1p1)))
'''An Albers equal-area (authalic) projection with a single standard parallel.
@see: L{AlbersEqualArea2} and L{AlbersEqualArea4}. ''' '''New L{AlbersEqualArea} projection.
@arg lat: Standard parallel (C{degrees}). @kwarg k0: Azimuthal scale on the standard parallel (C{scalar}). @kwarg datum: Optional datum or ellipsoid (L{Datum}, L{Ellipsoid}, L{Ellipsoid2} or L{a_f2Tuple}). @kwarg name: Optional name for the projection (C{str}).
@raise AlbertError: Invalid B{C{lat}}, B{C{k0}} or no convergence. '''
'''An Albers equal-area (authalic) projection with two standard parallels.
@see: L{AlbersEqualArea} and L{AlbersEqualArea4}. ''' '''New L{AlbersEqualArea2} projection.
@arg lat1: First standard parallel (C{degrees}). @arg lat2: Second standard parallel (C{degrees}). @kwarg k1: Azimuthal scale on the standard parallels (C{scalar}). @kwarg datum: Optional datum or ellipsoid (L{Datum}, L{Ellipsoid}, L{Ellipsoid2} or L{a_f2Tuple}). @kwarg name: Optional name for the projection (C{str}).
@raise AlbertError: Invalid B{C{lat1}}m B{C{lat2}}, B{C{k1}} or no convergence. '''
'''An Albers equal-area (authalic) projection specified by the C{sin} and C{cos} of both standard parallels.
@see: L{AlbersEqualArea} and L{AlbersEqualArea2}. ''' '''New L{AlbersEqualArea4} projection.
@arg slat1: Sine of first standard parallel (C{scalar}). @arg clat1: Cosine of first standard parallel (non-negative C{scalar}). @arg slat2: Sine of second standard parallel (C{scalar}). @arg clat2: Cosine of second standard parallel (non-negative C{scalar}). @kwarg k1: Azimuthal scale on the standard parallels (C{scalar}). @kwarg datum: Optional datum or ellipsoid (L{Datum}, L{Ellipsoid}, L{Ellipsoid2} or L{a_f2Tuple}). @kwarg name: Optional name for the projection (C{str}).
@raise AlbertError: Negative B{C{clat1}} or B{C{clat2}}, B{C{slat1}} and B{C{slat2}} have opposite signs (hemispheres), invalid B{C{k1}} or no convergence. '''
'''An L{AlbersEqualArea} projection at C{lat=0} and C{k0=1} degenerating to the cylindrical-equal-area projection. '''
'''New L{AlbersEqualAreaCylindrical} projection.
@kwarg datum: Optional datum or ellipsoid (L{Datum}, L{Ellipsoid}, L{Ellipsoid2} or L{a_f2Tuple}). @kwarg name: Optional name for the projection (C{str}). '''
'''An azimuthal L{AlbersEqualArea} projection at C{lat=90} and C{k0=1} degenerating to the L{azimuthal} L{LambertEqualArea} projection. '''
'''New L{AlbersEqualAreaNorth} projection.
@kwarg datum: Optional datum or ellipsoid (L{Datum}, L{Ellipsoid}, L{Ellipsoid2} or L{a_f2Tuple}). @kwarg name: Optional name for the projection (C{str}). '''
'''An azimuthal L{AlbersEqualArea} projection at C{lat=-90} and C{k0=1} degenerating to the L{azimuthal} L{LambertEqualArea} projection. '''
'''New L{AlbersEqualAreaSouth} projection.
@kwarg datum: Optional datum or ellipsoid (L{Datum}, L{Ellipsoid}, L{Ellipsoid2} or L{a_f2Tuple}). @kwarg name: Optional name for the projection (C{str}). '''
'''7-Tuple C{(x, y, lat, lon, gamma, scale, datum)}, in C{meter}, C{meter}, C{degrees90}, C{degrees180}, C{degrees360}, C{scalar} and C{Datum} where C{(x, y)} is the projected, C{(lat, lon)} the geodetic location, C{gamma} the meridian convergence at point, the bearing of the y-axis measured clockwise from true North and C{scale} is the azimuthal scale of the projection at point. The radial scale is the reciprocal C{1 / scale}. '''
# **) MIT License # # Copyright (C) 2016-2021 -- 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. |