Commit e86a797942586cd229a79ca7920214626a7d1b22 - python-pyproj

New upstream version 2.1.0+ds Bas Couwenberg 5 years ago

16 changed file(s) with 1066 addition(s) and 813 deletion(s). Raw diff Collapse all Expand all

-1

.travis.yml less more

79	79
80	80	script:
81	81	- python -c "import pyproj; pyproj.Proj(init='epsg:4269')"
82		- nose2 -v --coverage pyproj --coverage-report xml
	82	- nose2 -v
	83	- nose2 --coverage pyproj --coverage-report xml # coverage report on separate line as does not show failures
83	84
84	85	after_success:
85	86	- coveralls

-4

README.md less more

0	0	pyproj
1	1	======
2	2
3		[![Build Status](https://travis-ci.org/jswhit/pyproj.svg)](https://travis-ci.org/jswhit/pyproj)
	3	[![Build Status](https://travis-ci.org/pyproj4/pyproj.svg)](https://travis-ci.org/pyproj4/pyproj)
4	4
5	5	[![Build status](https://ci.appveyor.com/api/projects/status/8xkka4s97uwhkc64/branch/master?svg=true
6	6	)](https://ci.appveyor.com/project/jswhit/pyproj)
	7
	8	[![Coverage Status](https://coveralls.io/repos/github/pyproj4/pyproj/badge.svg?branch=master)](https://coveralls.io/github/pyproj4/pyproj?branch=master)
7	9
8	10	[![PyPI version](https://badge.fury.io/py/pyproj.svg)](https://badge.fury.io/py/pyproj)
9	11

14	16	Python interface to [PROJ.4](https://github.com/OSGeo/proj.4).
15	17
16	18
17		[Installation](https://jswhit.github.io/pyproj/html/installation.html)
	19	[Installation](https://pyproj4.github.io/pyproj/html/installation.html)
18	20	------------
19	21
20		[Documentation](http://jswhit.github.io/pyproj)
	22	[Documentation](http://pyproj4.github.io/pyproj)
21	23	-------------
22	24
23	25	Bugs/Questions
24	26	--------------
25		Report bugs/ask questions at https://github.com/jswhit/pyproj/issues.
	27	Report bugs/ask questions at https://github.com/pyproj4/pyproj/issues.

+18

-551

pyproj/__init__.py less more

25	25
26	26	Download: http://python.org/pypi/pyproj
27	27
28		Requirements: Python 2.6, 2.7, 3.2 or higher version.
	28	Requirements: Python 2.7 or 3.5+.
29	29
30	30	Example scripts are in 'test' subdirectory of source distribution.
31	31	The 'test()' function will run the examples in the docstrings.

46	46	LOSS OF USE, DATA OR PROFITS, WHETHER IN AN ACTION OF CONTRACT,
47	47	NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT OF OR IN
48	48	CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE. """
49		__version__ = "2.0.2"
50		__all__ = ["Proj", "Geod", "CRS", "transform", "itransform", "pj_ellps", "pj_list"]
	49	__version__ = "2.1.0"
	50	__all__ = [
	51	"Proj",
	52	"Geod",
	53	"CRS",
	54	"Transformer",
	55	"transform",
	56	"itransform",
	57	"pj_ellps",
	58	"pj_list",
	59	]
51	60
52		import re
53		import sys
54		import warnings
55		from array import array
56		from itertools import chain, islice
57
58		from pyproj import _proj
59		from pyproj.compat import cstrencode, pystrdecode
60	61	from pyproj.crs import CRS
61	62	from pyproj.exceptions import ProjError
62	63	from pyproj.geod import Geod, geodesic_version_str, pj_ellps
63		from pyproj.utils import _convertback, _copytobuffer
64
65		try:
66		from future_builtins import zip # python 2.6+
67		except ImportError:
68		pass # python 3.x
69
70
71		# import numpy as np
72		proj_version_str = _proj.proj_version_str
73
74		pj_list = {
75		"aea": "Albers Equal Area",
76		"aeqd": "Azimuthal Equidistant",
77		"affine": "Affine transformation",
78		"airy": "Airy",
79		"aitoff": "Aitoff",
80		"alsk": "Mod. Stererographics of Alaska",
81		"apian": "Apian Globular I",
82		"august": "August Epicycloidal",
83		"bacon": "Bacon Globular",
84		"bertin1953": "Bertin 1953",
85		"bipc": "Bipolar conic of western hemisphere",
86		"boggs": "Boggs Eumorphic",
87		"bonne": "Bonne (Werner lat_1=90)",
88		"calcofi": "Cal Coop Ocean Fish Invest Lines/Stations",
89		"cart": "Geodetic/cartesian conversions",
90		"cass": "Cassini",
91		"cc": "Central Cylindrical",
92		"ccon": "Central Conic",
93		"cea": "Equal Area Cylindrical",
94		"chamb": "Chamberlin Trimetric",
95		"collg": "Collignon",
96		"comill": "Compact Miller",
97		"crast": "Craster Parabolic (Putnins P4)",
98		"deformation": "Kinematic grid shift",
99		"denoy": "Denoyer Semi-Elliptical",
100		"eck1": "Eckert I",
101		"eck2": "Eckert II",
102		"eck3": "Eckert III",
103		"eck4": "Eckert IV",
104		"eck5": "Eckert V",
105		"eck6": "Eckert VI",
106		"eqc": "Equidistant Cylindrical (Plate Caree)",
107		"eqdc": "Equidistant Conic",
108		"euler": "Euler",
109		"etmerc": "Extended Transverse Mercator",
110		"fahey": "Fahey",
111		"fouc": "Foucaut",
112		"fouc_s": "Foucaut Sinusoidal",
113		"gall": "Gall (Gall Stereographic)",
114		"geoc": "Geocentric Latitude",
115		"geocent": "Geocentric",
116		"geogoffset": "Geographic Offset",
117		"geos": "Geostationary Satellite View",
118		"gins8": "Ginsburg VIII (TsNIIGAiK)",
119		"gn_sinu": "General Sinusoidal Series",
120		"gnom": "Gnomonic",
121		"goode": "Goode Homolosine",
122		"gs48": "Mod. Stererographics of 48 U.S.",
123		"gs50": "Mod. Stererographics of 50 U.S.",
124		"hammer": "Hammer & Eckert-Greifendorff",
125		"hatano": "Hatano Asymmetrical Equal Area",
126		"healpix": "HEALPix",
127		"rhealpix": "rHEALPix",
128		"helmert": "3- and 7-parameter Helmert shift",
129		"hgridshift": "Horizontal grid shift",
130		"horner": "Horner polynomial evaluation",
131		"igh": "Interrupted Goode Homolosine",
132		"imw_p": "International Map of the World Polyconic",
133		"isea": "Icosahedral Snyder Equal Area",
134		"kav5": "Kavraisky V",
135		"kav7": "Kavraisky VII",
136		"krovak": "Krovak",
137		"labrd": "Laborde",
138		"laea": "Lambert Azimuthal Equal Area",
139		"lagrng": "Lagrange",
140		"larr": "Larrivee",
141		"lask": "Laskowski",
142		"lonlat": "Lat/long (Geodetic)",
143		"latlon": "Lat/long (Geodetic alias)",
144		"latlong": "Lat/long (Geodetic alias)",
145		"longlat": "Lat/long (Geodetic alias)",
146		"lcc": "Lambert Conformal Conic",
147		"lcca": "Lambert Conformal Conic Alternative",
148		"leac": "Lambert Equal Area Conic",
149		"lee_os": "Lee Oblated Stereographic",
150		"loxim": "Loximuthal",
151		"lsat": "Space oblique for LANDSAT",
152		"mbt_s": "McBryde-Thomas Flat-Polar Sine",
153		"mbt_fps": "McBryde-Thomas Flat-Pole Sine (No. 2)",
154		"mbtfpp": "McBride-Thomas Flat-Polar Parabolic",
155		"mbtfpq": "McBryde-Thomas Flat-Polar Quartic",
156		"mbtfps": "McBryde-Thomas Flat-Polar Sinusoidal",
157		"merc": "Mercator",
158		"mil_os": "Miller Oblated Stereographic",
159		"mill": "Miller Cylindrical",
160		"misrsom": "Space oblique for MISR",
161		"moll": "Mollweide",
162		"molobadekas": "Molodensky-Badekas transform",
163		"molodensky": "Molodensky transform",
164		"murd1": "Murdoch I",
165		"murd2": "Murdoch II",
166		"murd3": "Murdoch III",
167		"natearth": "Natural Earth",
168		"natearth2": "Natural Earth 2",
169		"nell": "Nell",
170		"nell_h": "Nell-Hammer",
171		"nicol": "Nicolosi Globular",
172		"nsper": "Near-sided perspective",
173		"nzmg": "New Zealand Map Grid",
174		"ob_tran": "General Oblique Transformation",
175		"ocea": "Oblique Cylindrical Equal Area",
176		"oea": "Oblated Equal Area",
177		"omerc": "Oblique Mercator",
178		"ortel": "Ortelius Oval",
179		"ortho": "Orthographic",
180		"patterson": "Patterson Cylindrical",
181		"pconic": "Perspective Conic",
182		"pipeline": "Transformation pipeline manager",
183		"poly": "Polyconic (American)",
184		"pop": "Retrieve coordinate value from pipeline stack",
185		"putp1": "Putnins P1",
186		"putp2": "Putnins P2",
187		"putp3": "Putnins P3",
188		"putp3p": "Putnins P3'",
189		"putp4p": "Putnins P4'",
190		"putp5": "Putnins P5",
191		"putp5p": "Putnins P5'",
192		"putp6": "Putnins P6",
193		"putp6p": "Putnins P6'",
194		"qua_aut": "Quartic Authalic",
195		"robin": "Robinson",
196		"rouss": "Roussilhe Stereographic",
197		"rpoly": "Rectangular Polyconic",
198		"sch": "Spherical Cross-track Height",
199		"sinu": "Sinusoidal (Sanson-Flamsteed)",
200		"somerc": "Swiss. Obl. Mercator",
201		"stere": "Stereographic",
202		"sterea": "Oblique Stereographic Alternative",
203		"gstmerc": "Gauss-Schreiber Transverse Mercator (aka Gauss-Laborde Reunion)",
204		"tcc": "Transverse Central Cylindrical",
205		"tcea": "Transverse Cylindrical Equal Area",
206		"times": "Times",
207		"tissot": "Tissot Conic",
208		"tmerc": "Transverse Mercator",
209		"tpeqd": "Two Point Equidistant",
210		"tpers": "Tilted perspective",
211		"ups": "Universal Polar Stereographic",
212		"urm5": "Urmaev V",
213		"urmfps": "Urmaev Flat-Polar Sinusoidal",
214		"utm": "Universal Transverse Mercator (UTM)",
215		"vandg": "van der Grinten (I)",
216		"vandg2": "van der Grinten II",
217		"vandg3": "van der Grinten III",
218		"vandg4": "van der Grinten IV",
219		"vitk1": "Vitkovsky I",
220		"wag1": "Wagner I (Kavraisky VI)",
221		"wag2": "Wagner II",
222		"wag3": "Wagner III",
223		"wag4": "Wagner IV",
224		"wag5": "Wagner V",
225		"wag6": "Wagner VI",
226		"wag7": "Wagner VII",
227		"webmerc": "Web Mercator / Pseudo Mercator",
228		"weren": "Werenskiold I",
229		"wink1": "Winkel I",
230		"wink2": "Winkel II",
231		"wintri": "Winkel Tripel",
232		}
233
234
235		class Proj(_proj.Proj):
236		"""
237		performs cartographic transformations (converts from
238		longitude,latitude to native map projection x,y coordinates and
239		vice versa) using proj (https://github.com/OSGeo/proj.4/wiki).
240
241		A Proj class instance is initialized with proj map projection
242		control parameter key/value pairs. The key/value pairs can
243		either be passed in a dictionary, or as keyword arguments,
244		or as a proj4 string (compatible with the proj command). See
245		http://www.remotesensing.org/geotiff/proj_list for examples of
246		key/value pairs defining different map projections.
247
248		Calling a Proj class instance with the arguments lon, lat will
249		convert lon/lat (in degrees) to x/y native map projection
250		coordinates (in meters). If optional keyword 'inverse' is True
251		(default is False), the inverse transformation from x/y to
252		lon/lat is performed. If optional keyword 'errcheck' is True (default is
253		False) an exception is raised if the transformation is invalid.
254		If errcheck=False and the transformation is invalid, no
255		exception is raised and 1.e30 is returned. If the optional keyword
256		'preserve_units' is True, the units in map projection coordinates
257		are not forced to be meters.
258
259		Works with numpy and regular python array objects, python
260		sequences and scalars.
261		"""
262
263		def __init__(self, projparams=None, preserve_units=True, **kwargs):
264		"""
265		initialize a Proj class instance.
266
267		See the proj documentation (https://github.com/OSGeo/proj.4/wiki)
268		for more information about projection parameters.
269
270		Parameters
271		----------
272		projparams: int, str, dict, pyproj.CRS
273		A proj.4 or WKT string, proj.4 dict, EPSG integer, or a pyproj.CRS instnace.
274		preserve_units: bool
275		If false, will ensure +units=m.
276		**kwargs:
277		proj.4 projection parameters.
278
279
280		Example usage:
281
282		>>> from pyproj import Proj
283		>>> p = Proj(proj='utm',zone=10,ellps='WGS84', preserve_units=False) # use kwargs
284		>>> x,y = p(-120.108, 34.36116666)
285		>>> 'x=%9.3f y=%11.3f' % (x,y)
286		'x=765975.641 y=3805993.134'
287		>>> 'lon=%8.3f lat=%5.3f' % p(x,y,inverse=True)
288		'lon=-120.108 lat=34.361'
289		>>> # do 3 cities at a time in a tuple (Fresno, LA, SF)
290		>>> lons = (-119.72,-118.40,-122.38)
291		>>> lats = (36.77, 33.93, 37.62 )
292		>>> x,y = p(lons, lats)
293		>>> 'x: %9.3f %9.3f %9.3f' % x
294		'x: 792763.863 925321.537 554714.301'
295		>>> 'y: %9.3f %9.3f %9.3f' % y
296		'y: 4074377.617 3763936.941 4163835.303'
297		>>> lons, lats = p(x, y, inverse=True) # inverse transform
298		>>> 'lons: %8.3f %8.3f %8.3f' % lons
299		'lons: -119.720 -118.400 -122.380'
300		>>> 'lats: %8.3f %8.3f %8.3f' % lats
301		'lats: 36.770 33.930 37.620'
302		>>> p2 = Proj('+proj=utm +zone=10 +ellps=WGS84', preserve_units=False) # use proj4 string
303		>>> x,y = p2(-120.108, 34.36116666)
304		>>> 'x=%9.3f y=%11.3f' % (x,y)
305		'x=765975.641 y=3805993.134'
306		>>> p = Proj(init="epsg:32667", preserve_units=False)
307		>>> 'x=%12.3f y=%12.3f (meters)' % p(-114.057222, 51.045)
308		'x=-1783506.250 y= 6193827.033 (meters)'
309		>>> p = Proj("+init=epsg:32667")
310		>>> 'x=%12.3f y=%12.3f (feet)' % p(-114.057222, 51.045)
311		'x=-5851386.754 y=20320914.191 (feet)'
312		>>> # test data with radian inputs
313		>>> p1 = Proj(init="epsg:4214")
314		>>> x1, y1 = p1(116.366, 39.867)
315		>>> '{:.3f} {:.3f}'.format(x1, y1)
316		'2.031 0.696'
317		>>> x2, y2 = p1(x1, y1, inverse=True)
318		>>> '{:.3f} {:.3f}'.format(x2, y2)
319		'116.366 39.867'
320		"""
321		self.crs = CRS.from_user_input(projparams if projparams is not None else kwargs)
322		# make sure units are meters if preserve_units is False.
323		if not preserve_units and self.crs.is_projected:
324		projstring = self.crs.to_proj4(4)
325		projstring = re.sub(r"\s\+units=[\w-]+", "", projstring)
326		projstring += " +units=m"
327		self.crs = CRS(projstring)
328		super(Proj, self).__init__(
329		cstrencode(self.crs.to_proj4().replace("+type=crs", "").strip())
330		)
331
332		def __call__(self, args, *kw):
333		# ,lon,lat,inverse=False,errcheck=False):
334		"""
335		Calling a Proj class instance with the arguments lon, lat will
336		convert lon/lat (in degrees) to x/y native map projection
337		coordinates (in meters). If optional keyword 'inverse' is True
338		(default is False), the inverse transformation from x/y to
339		lon/lat is performed. If optional keyword 'errcheck' is True (default is
340		False) an exception is raised if the transformation is invalid.
341		If errcheck=False and the transformation is invalid, no
342		exception is raised and 1.e30 is returned.
343
344		Inputs should be doubles (they will be cast to doubles if they
345		are not, causing a slight performance hit).
346
347		Works with numpy and regular python array objects, python
348		sequences and scalars, but is fastest for array objects.
349		"""
350		inverse = kw.get("inverse", False)
351		errcheck = kw.get("errcheck", False)
352		# if len(args) == 1:
353		# latlon = np.array(args[0], copy=True,
354		# order='C', dtype=float, ndmin=2)
355		# if inverse:
356		# _proj.Proj._invn(self, latlon, radians=radians, errcheck=errcheck)
357		# else:
358		# _proj.Proj._fwdn(self, latlon, radians=radians, errcheck=errcheck)
359		# return latlon
360		lon, lat = args
361		# process inputs, making copies that support buffer API.
362		inx, xisfloat, xislist, xistuple = _copytobuffer(lon)
363		iny, yisfloat, yislist, yistuple = _copytobuffer(lat)
364		# call proj4 functions. inx and iny modified in place.
365		if inverse:
366		self._inv(inx, iny, errcheck=errcheck)
367		else:
368		self._fwd(inx, iny, errcheck=errcheck)
369		# if inputs were lists, tuples or floats, convert back.
370		outx = _convertback(xisfloat, xislist, xistuple, inx)
371		outy = _convertback(yisfloat, yislist, xistuple, iny)
372		return outx, outy
373
374		def is_latlong(self):
375		"""
376		Returns
377		-------
378		bool: True if projection in geographic (lon/lat) coordinates.
379		"""
380		warnings.warn("'is_latlong()' is deprecated. Please use 'crs.is_geographic'.")
381		return self.crs.is_geographic
382
383		def is_geocent(self):
384		"""
385		Returns
386		-------
387		bool: True if projection in geocentric (x/y) coordinates
388		"""
389		warnings.warn("'is_geocent()' is deprecated. Please use 'crs.is_geocent'.")
390		return self.is_geocent
391
392		def definition_string(self):
393		"""Returns formal definition string for projection
394
395		>>> Proj('+init=epsg:4326').definition_string()
396		'proj=longlat datum=WGS84 no_defs ellps=WGS84 towgs84=0,0,0'
397		>>>
398		"""
399		return pystrdecode(self.definition)
400
401		def to_latlong_def(self):
402		"""return the definition string of the geographic (lat/lon)
403		coordinate version of the current projection"""
404		# This is a little hacky way of getting a latlong proj object
405		# Maybe instead of this function the __cinit__ function can take a
406		# Proj object and a type (where type = "geographic") as the libproj
407		# java wrapper
408		return self.crs.to_geodetic().to_proj4(4)
409
410		# deprecated : using in transform raised a TypeError in release 1.9.5.1
411		# reported in issue #53, resolved in #73.
412		def to_latlong(self):
413		"""return a new Proj instance which is the geographic (lat/lon)
414		coordinate version of the current projection"""
415		return Proj(self.crs.to_geodetic())
416
417
418		def transform(p1, p2, x, y, z=None, radians=False):
419		"""
420		x2, y2, z2 = transform(p1, p2, x1, y1, z1)
421
422		Transform points between two coordinate systems defined by the
423		Proj instances p1 and p2.
424
425		The points x1,y1,z1 in the coordinate system defined by p1 are
426		transformed to x2,y2,z2 in the coordinate system defined by p2.
427
428		z1 is optional, if it is not set it is assumed to be zero (and
429		only x2 and y2 are returned). If the optional keyword
430		'radians' is True (default is False), then all input and
431		output coordinates will be in radians instead of the default
432		of degrees for geographic input/output projections.
433
434		In addition to converting between cartographic and geographic
435		projection coordinates, this function can take care of datum
436		shifts (which cannot be done using the __call__ method of the
437		Proj instances). It also allows for one of the coordinate
438		systems to be geographic (proj = 'latlong').
439
440		x,y and z can be numpy or regular python arrays, python
441		lists/tuples or scalars. Arrays are fastest. For projections in
442		geocentric coordinates, values of x and y are given in meters.
443		z is always meters.
444
445		Example usage:
446
447		>>> # projection 1: UTM zone 15, grs80 ellipse, NAD83 datum
448		>>> # (defined by epsg code 26915)
449		>>> p1 = Proj(init='epsg:26915', preserve_units=False)
450		>>> # projection 2: UTM zone 15, clrk66 ellipse, NAD27 datum
451		>>> p2 = Proj(init='epsg:26715', preserve_units=False)
452		>>> # find x,y of Jefferson City, MO.
453		>>> x1, y1 = p1(-92.199881,38.56694)
454		>>> # transform this point to projection 2 coordinates.
455		>>> x2, y2 = transform(p1,p2,x1,y1)
456		>>> '%9.3f %11.3f' % (x1,y1)
457		'569704.566 4269024.671'
458		>>> '%9.3f %11.3f' % (x2,y2)
459		'569722.342 4268814.028'
460		>>> '%8.3f %5.3f' % p2(x2,y2,inverse=True)
461		' -92.200 38.567'
462		>>> # process 3 points at a time in a tuple
463		>>> lats = (38.83,39.32,38.75) # Columbia, KC and StL Missouri
464		>>> lons = (-92.22,-94.72,-90.37)
465		>>> x1, y1 = p1(lons,lats)
466		>>> x2, y2 = transform(p1,p2,x1,y1)
467		>>> xy = x1+y1
468		>>> '%9.3f %9.3f %9.3f %11.3f %11.3f %11.3f' % xy
469		'567703.344 351730.944 728553.093 4298200.739 4353698.725 4292319.005'
470		>>> xy = x2+y2
471		>>> '%9.3f %9.3f %9.3f %11.3f %11.3f %11.3f' % xy
472		'567721.149 351747.558 728569.133 4297989.112 4353489.645 4292106.305'
473		>>> lons, lats = p2(x2,y2,inverse=True)
474		>>> xy = lons+lats
475		>>> '%8.3f %8.3f %8.3f %5.3f %5.3f %5.3f' % xy
476		' -92.220 -94.720 -90.370 38.830 39.320 38.750'
477		>>> # test datum shifting, installation of extra datum grid files.
478		>>> p1 = Proj(proj='latlong',datum='WGS84')
479		>>> x1 = -111.5; y1 = 45.25919444444
480		>>> p2 = Proj(proj="utm",zone=10,datum='NAD27', preserve_units=False)
481		>>> x2, y2 = transform(p1, p2, x1, y1)
482		>>> "%s %s" % (str(x2)[:9],str(y2)[:9])
483		'1402285.9 5076292.4'
484		>>> from pyproj import CRS
485		>>> c1 = CRS(proj='latlong',datum='WGS84')
486		>>> x1 = -111.5; y1 = 45.25919444444
487		>>> c2 = CRS(proj="utm",zone=10,datum='NAD27')
488		>>> x2, y2 = transform(c1, c2, x1, y1)
489		>>> "%s %s" % (str(x2)[:9],str(y2)[:9])
490		'1402291.0 5076289.5'
491		>>> x3, y3 = transform("epsg:4326", "epsg:3857", 33, 98)
492		>>> "%.3f %.3f" % (x3, y3)
493		'10909310.098 3895303.963'
494		>>> pj = Proj(init="epsg:4214")
495		>>> pjx, pjy = pj(116.366, 39.867)
496		>>> xr, yr = transform(pj, Proj(4326), pjx, pjy, radians=True)
497		>>> "%.3f %.3f" % (xr, yr)
498		'2.031 0.696'
499		"""
500		# check that p1 and p2 are valid
501		if not isinstance(p1, Proj):
502		p1 = CRS.from_user_input(p1)
503		if not isinstance(p2, Proj):
504		p2 = CRS.from_user_input(p2)
505
506		# process inputs, making copies that support buffer API.
507		inx, xisfloat, xislist, xistuple = _copytobuffer(x)
508		iny, yisfloat, yislist, yistuple = _copytobuffer(y)
509		if z is not None:
510		inz, zisfloat, zislist, zistuple = _copytobuffer(z)
511		else:
512		inz = None
513		# call pj_transform. inx,iny,inz buffers modified in place.
514		_proj._transform(p1, p2, inx, iny, inz, radians)
515		# if inputs were lists, tuples or floats, convert back.
516		outx = _convertback(xisfloat, xislist, xistuple, inx)
517		outy = _convertback(yisfloat, yislist, xistuple, iny)
518		if inz is not None:
519		outz = _convertback(zisfloat, zislist, zistuple, inz)
520		return outx, outy, outz
521		else:
522		return outx, outy
523
524
525		def itransform(p1, p2, points, switch=False, radians=False):
526		"""
527		points2 = transform(p1, p2, points1)
528		Iterator/generator version of the function pyproj.transform.
529		Transform points between two coordinate systems defined by the
530		Proj instances p1 and p2. This function can be used as an alternative
531		to pyproj.transform when there is a need to transform a big number of
532		coordinates lazily, for example when reading and processing from a file.
533		Points1 is an iterable/generator of coordinates x1,y1(,z1) or lon1,lat1(,z1)
534		in the coordinate system defined by p1. Points2 is an iterator that returns tuples
535		of x2,y2(,z2) or lon2,lat2(,z2) coordinates in the coordinate system defined by p2.
536		z are provided optionally.
537
538		Points1 can be:
539		- a tuple/list of tuples/lists i.e. for 2d points: [(xi,yi),(xj,yj),....(xn,yn)]
540		- a Nx3 or Nx2 2d numpy array where N is the point number
541		- a generator of coordinates (xi,yi) for 2d points or (xi,yi,zi) for 3d
542
543		If optional keyword 'switch' is True (default is False) then x, y or lon,lat coordinates
544		of points are switched to y, x or lat, lon. If the optional keyword 'radians' is True
545		(default is False), then all input and output coordinates will be in radians instead
546		of the default of degrees for geographic input/output projections.
547
548
549		Example usage:
550
551		>>> # projection 1: WGS84
552		>>> # (defined by epsg code 4326)
553		>>> p1 = Proj(init='epsg:4326', preserve_units=False)
554		>>> # projection 2: GGRS87 / Greek Grid
555		>>> p2 = Proj(init='epsg:2100', preserve_units=False)
556		>>> # Three points with coordinates lon, lat in p1
557		>>> points = [(22.95, 40.63), (22.81, 40.53), (23.51, 40.86)]
558		>>> # transform this point to projection 2 coordinates.
559		>>> for pt in itransform(p1,p2,points): '%6.3f %7.3f' % pt
560		'411200.657 4498214.742'
561		'399210.500 4487264.963'
562		'458703.102 4523331.451'
563		>>> for pt in itransform(4326, 2100, points): '{:.3f} {:.3f}'.format(*pt)
564		'2221638.801 2637034.372'
565		'2212924.125 2619851.898'
566		'2238294.779 2703763.736'
567		>>> pj = Proj(init="epsg:4214")
568		>>> pjx, pjy = pj(116.366, 39.867)
569		>>> for pt in itransform(pj, Proj(4326), [(pjx, pjy)]): '{:.3f} {:.3f}'.format(*pt)
570		'2.031 0.696'
571		"""
572		if not isinstance(p1, Proj):
573		p1 = CRS.from_user_input(p1)
574		if not isinstance(p2, Proj):
575		p2 = CRS.from_user_input(p2)
576
577		it = iter(points) # point iterator
578		# get first point to check stride
579		try:
580		fst_pt = next(it)
581		except StopIteration:
582		raise ValueError("iterable must contain at least one point")
583
584		stride = len(fst_pt)
585		if stride not in (2, 3):
586		raise ValueError("points can contain up to 3 coordinates")
587
588		# create a coordinate sequence generator etc. x1,y1,z1,x2,y2,z2,....
589		# chain so the generator returns the first point that was already acquired
590		coord_gen = chain(fst_pt, (coords[c] for coords in it for c in range(stride)))
591
592		while True:
593		# create a temporary buffer storage for the next 64 points (64stride8 bytes)
594		buff = array("d", islice(coord_gen, 0, 64 * stride))
595		if len(buff) == 0:
596		break
597
598		_proj._transform_sequence(p1, p2, stride, buff, switch, radians)
599
600		for pt in zip(([iter(buff)] stride)):
601		yield pt
	64	from pyproj.proj import Proj, pj_list, proj_version_str
	65	from pyproj.transformer import Transformer, itransform, transform
602	66
603	67
604	68	def test(**kwargs):

607	71	import pyproj
608	72
609	73	verbose = kwargs.get("verbose")
610		failure_count, test_count = doctest.testmod(pyproj, verbose=verbose)
	74	failure_count, test_count = doctest.testmod(pyproj.proj, verbose=verbose)
611	75	failure_count_crs, test_count_crs = doctest.testmod(pyproj.crs, verbose=verbose)
612	76	failure_count_geod, test_count_geod = doctest.testmod(pyproj.geod, verbose=verbose)
	77	failure_count_transform, test_count_transform = doctest.testmod(
	78	pyproj.transformer, verbose=verbose
	79	)
613	80
614	81	return failure_count + failure_count_crs + failure_count_geod
615	82

-241

pyproj/_proj.pyx less more

0		# cimport c_numpy
1		# c_numpy.import_array()
2	0	include "base.pxi"
3	1
4		from pyproj.crs import CRS
5	2	from pyproj.compat import cstrencode, pystrdecode
6	3	from pyproj.datadir import get_data_dir
7	4	from pyproj.exceptions import ProjError

174	171	xdatab[i] = projlonlatout.uv.u
175	172	ydatab[i] = projlonlatout.uv.v
176	173
177		# def _fwdn(self, c_numpy.ndarray lonlat, radians=False, errcheck=False):
178		# """
179		# forward transformation - lons,lats to x,y (done in place).
180		# Uses ndarray of shape ...,2.
181		# if radians=True, lons/lats are radians instead of degrees.
182		# if errcheck=True, an exception is raised if the forward
183		# transformation is invalid.
184		# if errcheck=False and the forward transformation is
185		# invalid, no exception is
186		# raised and 1.e30 is returned.
187		# """
188		# cdef PJ_UV projxyout, projlonlatin
189		# cdef PJ_UV *llptr
190		# cdef int err
191		# cdef Py_ssize_t npts, i
192		# npts = c_numpy.PyArray_SIZE(lonlat)//2
193		# llptr = <PJ_UV *>lonlat.data
194		# for i from 0 <= i < npts:
195		# if radians:
196		# projlonlatin = llptr[i]
197		# else:
198		# projlonlatin.u = _DG2RAD*llptr[i].u
199		# projlonlatin.v = _DG2RAD*llptr[i].v
200		# projxyout = pj_fwd(projlonlatin,self.projpj)
201		#
202		# if errcheck:
203		# err = proj_context_errno(self.projctx)
204		# if err != 0:
205		# raise ProjError(proj_errno_string(err))
206		# # since HUGE_VAL can be 'inf',
207		# # change it to a real (but very large) number.
208		# if projxyout.u == HUGE_VAL:
209		# llptr[i].u = 1.e30
210		# else:
211		# llptr[i].u = projxyout.u
212		# if projxyout.v == HUGE_VAL:
213		# llptr[i].u = 1.e30
214		# else:
215		# llptr[i].v = projxyout.v
216		#
217		# def _invn(self, c_numpy.ndarray xy, radians=False, errcheck=False):
218		# """
219		# inverse transformation - x,y to lons,lats (done in place).
220		# Uses ndarray of shape ...,2.
221		# if radians=True, lons/lats are radians instead of degrees.
222		# if errcheck=True, an exception is raised if the inverse transformation is invalid.
223		# if errcheck=False and the inverse transformation is invalid, no exception is
224		# raised and 1.e30 is returned.
225		# """
226		# cdef PJ_UV projxyin, projlonlatout
227		# cdef PJ_UV *llptr
228		# cdef Py_ssize_t npts, i
229		# npts = c_numpy.PyArray_SIZE(xy)//2
230		# llptr = <PJ_UV *>xy.data
231		#
232		# for i from 0 <= i < npts:
233		# projxyin = llptr[i]
234		# projlonlatout = pj_inv(projxyin, self.projpj)
235		# if errcheck:
236		# err = proj_context_errno(self.projctx)
237		# if err != 0:
238		# raise ProjError(proj_errno_string(err))
239		# # since HUGE_VAL can be 'inf',
240		# # change it to a real (but very large) number.
241		# if projlonlatout.u == HUGE_VAL:
242		# llptr[i].u = 1.e30
243		# elif radians:
244		# llptr[i].u = projlonlatout.u
245		# else:
246		# llptr[i].u = _RAD2DG*projlonlatout.u
247		# if projlonlatout.v == HUGE_VAL:
248		# llptr[i].v = 1.e30
249		# elif radians:
250		# llptr[i].v = projlonlatout.v
251		# else:
252		# llptr[i].v = _RAD2DG*projlonlatout.v
253
254	174	def __repr__(self):
255	175	return "Proj('{srs}', preserve_units=True)".format(srs=self.srs)
256
257		cdef class TransProj:
258		def __cinit__(self):
259		self.projpj = NULL
260		self.projctx = NULL
261
262		def __init__(self, p1, p2):
263		# set up the context
264		self.projctx = proj_context_create()
265		py_data_dir = cstrencode(get_data_dir())
266		cdef const char* data_dir = py_data_dir
267		proj_context_set_search_paths(self.projctx, 1, &data_dir)
268
269		self.projpj = proj_create_crs_to_crs(
270		self.projctx,
271		TransProj.definition_from_object(p1),
272		TransProj.definition_from_object(p2),
273		NULL)
274		if self.projpj is NULL:
275		raise ProjError("Error creating CRS to CRS.")
276
277		def __dealloc__(self):
278		"""destroy projection definition"""
279		if self.projpj is not NULL:
280		proj_destroy(self.projpj)
281		if self.projctx is not NULL:
282		proj_context_destroy(self.projctx)
283
284		@staticmethod
285		def definition_from_object(in_proj):
286		"""
287		Parameters
288		----------
289		in_proj: :obj:`pyproj.Proj` or :obj:`pyproj.CRS`
290
291		Returns
292		-------
293		char*: Definition string for `proj_create_crs_to_crs`.
294
295		"""
296		if isinstance(in_proj, Proj):
297		return in_proj.definition
298		return cstrencode(CRS.from_user_input(in_proj).to_wkt())
299
300
301		def is_geographic(proj):
302		if hasattr(proj, "crs"):
303		proj = proj.crs
304		return proj.is_geographic
305
306
307		def _transform(p1, p2, inx, iny, inz, radians):
308		pj_trans = TransProj(p1, p2)
309		# private function to call pj_transform
310		cdef void *xdata
311		cdef void *ydata
312		cdef void *zdata
313		cdef double *xx
314		cdef double *yy
315		cdef double *zz
316		cdef Py_ssize_t buflenx, bufleny, buflenz, npts, i
317		cdef int err
318		if PyObject_AsWriteBuffer(inx, &xdata, &buflenx) <> 0:
319		raise ProjError
320		if PyObject_AsWriteBuffer(iny, &ydata, &bufleny) <> 0:
321		raise ProjError
322		if inz is not None:
323		if PyObject_AsWriteBuffer(inz, &zdata, &buflenz) <> 0:
324		raise ProjError
325		else:
326		buflenz = bufleny
327		if not (buflenx == bufleny == buflenz):
328		raise ProjError('x,y and z must be same size')
329		xx = <double *>xdata
330		yy = <double *>ydata
331		if inz is not None:
332		zz = <double *>zdata
333		else:
334		zz = NULL
335		npts = buflenx//8
336
337		if radians and is_geographic(p1):
338		for i from 0 <= i < npts:
339		xx[i] = xx[i]*_RAD2DG
340		yy[i] = yy[i]*_RAD2DG
341
342		proj_trans_generic(
343		pj_trans.projpj,
344		PJ_FWD,
345		xx, _DOUBLESIZE, npts,
346		yy, _DOUBLESIZE, npts,
347		zz, _DOUBLESIZE, npts,
348		NULL, 0, 0,
349		)
350		cdef int errno = proj_errno(pj_trans.projpj)
351		if errno:
352		raise ProjError("proj_trans_generic error: {}".format(
353		pystrdecode(proj_errno_string(errno))))
354
355		if radians and is_geographic(p2):
356		for i from 0 <= i < npts:
357		xx[i] = xx[i]*_DG2RAD
358		yy[i] = yy[i]*_DG2RAD
359
360		def _transform_sequence(p1, p2, Py_ssize_t stride, inseq, bint switch, radians):
361		pj_trans = TransProj(p1, p2)
362		# private function to itransform function
363		cdef:
364		void *buffer
365		double *coords
366		double *x
367		double *y
368		double *z
369		Py_ssize_t buflen, npts, i, j
370		int err
371
372		if stride < 2:
373		raise ProjError("coordinates must contain at least 2 values")
374		if PyObject_AsWriteBuffer(inseq, &buffer, &buflen) <> 0:
375		raise ProjError("object does not provide the python buffer writeable interface")
376
377		coords = <double*>buffer
378		npts = buflen // (stride * _DOUBLESIZE)
379
380		if radians and is_geographic(p1):
381		for i from 0 <= i < npts:
382		j = stride*i
383		coords[j] *= _RAD2DG
384		coords[j+1] *= _RAD2DG
385
386		if not switch:
387		x = coords
388		y = coords + 1
389		else:
390		x = coords + 1
391		y = coords
392
393		if stride == 2:
394		z = NULL
395		else:
396		z = coords + 2
397
398		proj_trans_generic (
399		pj_trans.projpj,
400		PJ_FWD,
401		x, stride*_DOUBLESIZE, npts,
402		y, stride*_DOUBLESIZE, npts,
403		z, stride*_DOUBLESIZE, npts,
404		NULL, 0, 0,
405		)
406
407		cdef int errno = proj_errno(pj_trans.projpj)
408		if errno:
409		raise ProjError("proj_trans_generic error: {}".format(
410		proj_errno_string(errno)))
411
412		if radians and is_geographic(p2):
413		for i from 0 <= i < npts:
414		j = stride*i
415		coords[j] *= _DG2RAD
416		coords[j+1] *= _DG2RAD

+10

-0

pyproj/_transformer.pxd less more

	0	include "proj.pxi"
	1
	2	cdef class _Transformer:
	3	cdef PJ * projpj
	4	cdef PJ_CONTEXT * projctx
	5	cdef public object input_geographic
	6	cdef public object output_geographic
	7	cdef public object input_radians
	8	cdef public object output_radians
	9	cdef public object is_pipeline

+233

-0

pyproj/_transformer.pyx less more

	0	include "base.pxi"
	1
	2	from pyproj.crs import CRS
	3	from pyproj.proj import Proj
	4	from pyproj.compat import cstrencode, pystrdecode
	5	from pyproj.datadir import get_data_dir
	6	from pyproj.exceptions import ProjError
	7
	8	cdef class _Transformer:
	9	def __cinit__(self):
	10	self.projpj = NULL
	11	self.projctx = NULL
	12	self.input_geographic = False
	13	self.output_geographic = False
	14	self.input_radians = False
	15	self.output_radians = False
	16	self.is_pipeline = False
	17
	18	def __init__(self):
	19	# set up the context
	20	self.projctx = proj_context_create()
	21	py_data_dir = cstrencode(get_data_dir())
	22	cdef const char* data_dir = py_data_dir
	23	proj_context_set_search_paths(self.projctx, 1, &data_dir)
	24
	25	def __dealloc__(self):
	26	"""destroy projection definition"""
	27	if self.projpj is not NULL:
	28	proj_destroy(self.projpj)
	29	if self.projctx is not NULL:
	30	proj_context_destroy(self.projctx)
	31
	32	@staticmethod
	33	def _init_crs_to_crs(proj_from, proj_to):
	34	cdef _Transformer transformer = _Transformer()
	35	transformer.projpj = proj_create_crs_to_crs(
	36	transformer.projctx,
	37	_Transformer._definition_from_object(proj_from),
	38	_Transformer._definition_from_object(proj_to),
	39	NULL)
	40	if transformer.projpj is NULL:
	41	raise ProjError("Error creating CRS to CRS.")
	42	transformer.input_radians = proj_angular_input(transformer.projpj, PJ_FWD)
	43	transformer.output_radians = proj_angular_output(transformer.projpj, PJ_FWD)
	44	transformer.is_pipeline = False
	45	return transformer
	46
	47	@staticmethod
	48	def from_proj(proj_from, proj_to):
	49	if not isinstance(proj_from, Proj):
	50	proj_from = Proj(proj_from)
	51	if not isinstance(proj_to, Proj):
	52	proj_to = Proj(proj_to)
	53	transformer = _Transformer._init_crs_to_crs(proj_from, proj_to)
	54	transformer.input_geographic = proj_from.crs.is_geographic
	55	transformer.output_geographic = proj_to.crs.is_geographic
	56	return transformer
	57
	58	@staticmethod
	59	def from_crs(crs_from, crs_to):
	60	if not isinstance(crs_from, CRS):
	61	crs_from = CRS.from_user_input(crs_from)
	62	if not isinstance(crs_to, CRS):
	63	crs_to = CRS.from_user_input(crs_to)
	64	transformer = _Transformer._init_crs_to_crs(crs_from, crs_to)
	65	transformer.input_geographic = crs_from.is_geographic
	66	transformer.output_geographic = crs_to.is_geographic
	67	return transformer
	68
	69	@staticmethod
	70	def from_pipeline(const char *proj_pipeline):
	71	cdef _Transformer transformer = _Transformer()
	72
	73	# initialize projection
	74	transformer.projpj = proj_create(transformer.projctx, proj_pipeline)
	75	if transformer.projpj is NULL:
	76	raise ProjError("Invalid projection {}.".format(proj_pipeline))
	77	transformer.input_radians = proj_angular_input(transformer.projpj, PJ_FWD)
	78	transformer.output_radians = proj_angular_output(transformer.projpj, PJ_FWD)
	79	transformer.is_pipeline = True
	80	return transformer
	81
	82	@staticmethod
	83	def _definition_from_object(in_proj):
	84	"""
	85	Parameters
	86	----------
	87	in_proj: :obj:`pyproj.Proj` or :obj:`pyproj.CRS`
	88
	89	Returns
	90	-------
	91	char*: Definition string for `proj_create_crs_to_crs`.
	92
	93	"""
	94	if isinstance(in_proj, Proj):
	95	return cstrencode(in_proj.srs)
	96	return cstrencode(in_proj.to_wkt())
	97
	98	def _transform(self, inx, iny, inz, radians):
	99	# private function to call pj_transform
	100	cdef void *xdata
	101	cdef void *ydata
	102	cdef void *zdata
	103	cdef double *xx
	104	cdef double *yy
	105	cdef double *zz
	106	cdef Py_ssize_t buflenx, bufleny, buflenz, npts, i
	107	cdef int err
	108	if PyObject_AsWriteBuffer(inx, &xdata, &buflenx) <> 0:
	109	raise ProjError
	110	if PyObject_AsWriteBuffer(iny, &ydata, &bufleny) <> 0:
	111	raise ProjError
	112	if inz is not None:
	113	if PyObject_AsWriteBuffer(inz, &zdata, &buflenz) <> 0:
	114	raise ProjError
	115	else:
	116	buflenz = bufleny
	117	if not (buflenx == bufleny == buflenz):
	118	raise ProjError('x,y and z must be same size')
	119	xx = <double *>xdata
	120	yy = <double *>ydata
	121	if inz is not None:
	122	zz = <double *>zdata
	123	else:
	124	zz = NULL
	125	npts = buflenx//8
	126
	127	# degrees to radians
	128	if not self.is_pipeline and not radians and self.input_radians:
	129	for i from 0 <= i < npts:
	130	xx[i] = xx[i]*_DG2RAD
	131	yy[i] = yy[i]*_DG2RAD
	132	# radians to degrees
	133	elif not self.is_pipeline and radians and not self.input_radians and self.input_geographic:
	134	for i from 0 <= i < npts:
	135	xx[i] = xx[i]*_RAD2DG
	136	yy[i] = yy[i]*_RAD2DG
	137
	138	proj_trans_generic(
	139	self.projpj,
	140	PJ_FWD,
	141	xx, _DOUBLESIZE, npts,
	142	yy, _DOUBLESIZE, npts,
	143	zz, _DOUBLESIZE, npts,
	144	NULL, 0, 0,
	145	)
	146	cdef int errno = proj_errno(self.projpj)
	147	if errno:
	148	raise ProjError("proj_trans_generic error: {}".format(
	149	pystrdecode(proj_errno_string(errno))))
	150
	151	# radians to degrees
	152	if not self.is_pipeline and not radians and self.output_radians:
	153	for i from 0 <= i < npts:
	154	xx[i] = xx[i]*_RAD2DG
	155	yy[i] = yy[i]*_RAD2DG
	156	# degrees to radians
	157	elif not self.is_pipeline and radians and not self.output_radians and self.output_geographic:
	158	for i from 0 <= i < npts:
	159	xx[i] = xx[i]*_DG2RAD
	160	yy[i] = yy[i]*_DG2RAD
	161
	162
	163	def _transform_sequence(self, Py_ssize_t stride, inseq, bint switch, radians):
	164	# private function to itransform function
	165	cdef:
	166	void *buffer
	167	double *coords
	168	double *x
	169	double *y
	170	double *z
	171	Py_ssize_t buflen, npts, i, j
	172	int err
	173
	174	if stride < 2:
	175	raise ProjError("coordinates must contain at least 2 values")
	176	if PyObject_AsWriteBuffer(inseq, &buffer, &buflen) <> 0:
	177	raise ProjError("object does not provide the python buffer writeable interface")
	178
	179	coords = <double*>buffer
	180	npts = buflen // (stride * _DOUBLESIZE)
	181
	182	# degrees to radians
	183	if not self.is_pipeline and not radians and self.input_radians:
	184	for i from 0 <= i < npts:
	185	j = stride*i
	186	coords[j] *= _DG2RAD
	187	coords[j+1] *= _DG2RAD
	188	# radians to degrees
	189	elif not self.is_pipeline and radians and not self.input_radians and self.input_geographic:
	190	for i from 0 <= i < npts:
	191	j = stride*i
	192	coords[j] *= _RAD2DG
	193	coords[j+1] *= _RAD2DG
	194
	195	if not switch:
	196	x = coords
	197	y = coords + 1
	198	else:
	199	x = coords + 1
	200	y = coords
	201
	202	if stride == 2:
	203	z = NULL
	204	else:
	205	z = coords + 2
	206
	207	proj_trans_generic (
	208	self.projpj,
	209	PJ_FWD,
	210	x, stride*_DOUBLESIZE, npts,
	211	y, stride*_DOUBLESIZE, npts,
	212	z, stride*_DOUBLESIZE, npts,
	213	NULL, 0, 0,
	214	)
	215
	216	cdef int errno = proj_errno(self.projpj)
	217	if errno:
	218	raise ProjError("proj_trans_generic error: {}".format(
	219	proj_errno_string(errno)))
	220
	221	# radians to degrees
	222	if not self.is_pipeline and not radians and self.output_radians:
	223	for i from 0 <= i < npts:
	224	j = stride*i
	225	coords[j] *= _RAD2DG
	226	coords[j+1] *= _RAD2DG
	227	# degrees to radians
	228	elif not self.is_pipeline and radians and not self.output_radians and self.output_geographic:
	229	for i from 0 <= i < npts:
	230	j = stride*i
	231	coords[j] *= _DG2RAD
	232	coords[j+1] *= _DG2RAD

+389

-0

pyproj/proj.py less more

	0	# -- coding: utf-8 --
	1	"""
	2	Cython wrapper to provide python interfaces to
	3	PROJ.4 (https://github.com/OSGeo/proj.4/wiki) functions.
	4
	5	Performs cartographic transformations and geodetic computations.
	6
	7	The Proj class can convert from geographic (longitude,latitude)
	8	to native map projection (x,y) coordinates and vice versa, or
	9	from one map projection coordinate system directly to another.
	10	The module variable pj_list is a dictionary containing all the
	11	available projections and their descriptions.
	12
	13	Input coordinates can be given as python arrays, lists/tuples,
	14	scalars or numpy/Numeric/numarray arrays. Optimized for objects
	15	that support the Python buffer protocol (regular python and
	16	numpy array objects).
	17
	18	Download: http://python.org/pypi/pyproj
	19
	20	Contact: Jeffrey Whitaker <jeffrey.s.whitaker@noaa.gov
	21
	22	copyright (c) 2006 by Jeffrey Whitaker.
	23
	24	Permission to use, copy, modify, and distribute this software
	25	and its documentation for any purpose and without fee is hereby
	26	granted, provided that the above copyright notice appear in all
	27	copies and that both the copyright notice and this permission
	28	notice appear in supporting documentation. THE AUTHOR DISCLAIMS
	29	ALL WARRANTIES WITH REGARD TO THIS SOFTWARE, INCLUDING ALL
	30	IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT
	31	SHALL THE AUTHOR BE LIABLE FOR ANY SPECIAL, INDIRECT OR
	32	CONSEQUENTIAL DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM
	33	LOSS OF USE, DATA OR PROFITS, WHETHER IN AN ACTION OF CONTRACT,
	34	NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT OF OR IN
	35	CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE. """
	36	import re
	37	import warnings
	38
	39	from pyproj import _proj
	40	from pyproj.compat import cstrencode, pystrdecode
	41	from pyproj.crs import CRS
	42	from pyproj.utils import _convertback, _copytobuffer
	43
	44	# import numpy as np
	45	proj_version_str = _proj.proj_version_str
	46
	47	pj_list = {
	48	"aea": "Albers Equal Area",
	49	"aeqd": "Azimuthal Equidistant",
	50	"affine": "Affine transformation",
	51	"airy": "Airy",
	52	"aitoff": "Aitoff",
	53	"alsk": "Mod. Stererographics of Alaska",
	54	"apian": "Apian Globular I",
	55	"august": "August Epicycloidal",
	56	"bacon": "Bacon Globular",
	57	"bertin1953": "Bertin 1953",
	58	"bipc": "Bipolar conic of western hemisphere",
	59	"boggs": "Boggs Eumorphic",
	60	"bonne": "Bonne (Werner lat_1=90)",
	61	"calcofi": "Cal Coop Ocean Fish Invest Lines/Stations",
	62	"cart": "Geodetic/cartesian conversions",
	63	"cass": "Cassini",
	64	"cc": "Central Cylindrical",
	65	"ccon": "Central Conic",
	66	"cea": "Equal Area Cylindrical",
	67	"chamb": "Chamberlin Trimetric",
	68	"collg": "Collignon",
	69	"comill": "Compact Miller",
	70	"crast": "Craster Parabolic (Putnins P4)",
	71	"deformation": "Kinematic grid shift",
	72	"denoy": "Denoyer Semi-Elliptical",
	73	"eck1": "Eckert I",
	74	"eck2": "Eckert II",
	75	"eck3": "Eckert III",
	76	"eck4": "Eckert IV",
	77	"eck5": "Eckert V",
	78	"eck6": "Eckert VI",
	79	"eqc": "Equidistant Cylindrical (Plate Caree)",
	80	"eqdc": "Equidistant Conic",
	81	"euler": "Euler",
	82	"etmerc": "Extended Transverse Mercator",
	83	"fahey": "Fahey",
	84	"fouc": "Foucaut",
	85	"fouc_s": "Foucaut Sinusoidal",
	86	"gall": "Gall (Gall Stereographic)",
	87	"geoc": "Geocentric Latitude",
	88	"geocent": "Geocentric",
	89	"geogoffset": "Geographic Offset",
	90	"geos": "Geostationary Satellite View",
	91	"gins8": "Ginsburg VIII (TsNIIGAiK)",
	92	"gn_sinu": "General Sinusoidal Series",
	93	"gnom": "Gnomonic",
	94	"goode": "Goode Homolosine",
	95	"gs48": "Mod. Stererographics of 48 U.S.",
	96	"gs50": "Mod. Stererographics of 50 U.S.",
	97	"hammer": "Hammer & Eckert-Greifendorff",
	98	"hatano": "Hatano Asymmetrical Equal Area",
	99	"healpix": "HEALPix",
	100	"rhealpix": "rHEALPix",
	101	"helmert": "3- and 7-parameter Helmert shift",
	102	"hgridshift": "Horizontal grid shift",
	103	"horner": "Horner polynomial evaluation",
	104	"igh": "Interrupted Goode Homolosine",
	105	"imw_p": "International Map of the World Polyconic",
	106	"isea": "Icosahedral Snyder Equal Area",
	107	"kav5": "Kavraisky V",
	108	"kav7": "Kavraisky VII",
	109	"krovak": "Krovak",
	110	"labrd": "Laborde",
	111	"laea": "Lambert Azimuthal Equal Area",
	112	"lagrng": "Lagrange",
	113	"larr": "Larrivee",
	114	"lask": "Laskowski",
	115	"lonlat": "Lat/long (Geodetic)",
	116	"latlon": "Lat/long (Geodetic alias)",
	117	"latlong": "Lat/long (Geodetic alias)",
	118	"longlat": "Lat/long (Geodetic alias)",
	119	"lcc": "Lambert Conformal Conic",
	120	"lcca": "Lambert Conformal Conic Alternative",
	121	"leac": "Lambert Equal Area Conic",
	122	"lee_os": "Lee Oblated Stereographic",
	123	"loxim": "Loximuthal",
	124	"lsat": "Space oblique for LANDSAT",
	125	"mbt_s": "McBryde-Thomas Flat-Polar Sine",
	126	"mbt_fps": "McBryde-Thomas Flat-Pole Sine (No. 2)",
	127	"mbtfpp": "McBride-Thomas Flat-Polar Parabolic",
	128	"mbtfpq": "McBryde-Thomas Flat-Polar Quartic",
	129	"mbtfps": "McBryde-Thomas Flat-Polar Sinusoidal",
	130	"merc": "Mercator",
	131	"mil_os": "Miller Oblated Stereographic",
	132	"mill": "Miller Cylindrical",
	133	"misrsom": "Space oblique for MISR",
	134	"moll": "Mollweide",
	135	"molobadekas": "Molodensky-Badekas transform",
	136	"molodensky": "Molodensky transform",
	137	"murd1": "Murdoch I",
	138	"murd2": "Murdoch II",
	139	"murd3": "Murdoch III",
	140	"natearth": "Natural Earth",
	141	"natearth2": "Natural Earth 2",
	142	"nell": "Nell",
	143	"nell_h": "Nell-Hammer",
	144	"nicol": "Nicolosi Globular",
	145	"nsper": "Near-sided perspective",
	146	"nzmg": "New Zealand Map Grid",
	147	"ob_tran": "General Oblique Transformation",
	148	"ocea": "Oblique Cylindrical Equal Area",
	149	"oea": "Oblated Equal Area",
	150	"omerc": "Oblique Mercator",
	151	"ortel": "Ortelius Oval",
	152	"ortho": "Orthographic",
	153	"patterson": "Patterson Cylindrical",
	154	"pconic": "Perspective Conic",
	155	"pipeline": "Transformation pipeline manager",
	156	"poly": "Polyconic (American)",
	157	"pop": "Retrieve coordinate value from pipeline stack",
	158	"putp1": "Putnins P1",
	159	"putp2": "Putnins P2",
	160	"putp3": "Putnins P3",
	161	"putp3p": "Putnins P3'",
	162	"putp4p": "Putnins P4'",
	163	"putp5": "Putnins P5",
	164	"putp5p": "Putnins P5'",
	165	"putp6": "Putnins P6",
	166	"putp6p": "Putnins P6'",
	167	"qua_aut": "Quartic Authalic",
	168	"robin": "Robinson",
	169	"rouss": "Roussilhe Stereographic",
	170	"rpoly": "Rectangular Polyconic",
	171	"sch": "Spherical Cross-track Height",
	172	"sinu": "Sinusoidal (Sanson-Flamsteed)",
	173	"somerc": "Swiss. Obl. Mercator",
	174	"stere": "Stereographic",
	175	"sterea": "Oblique Stereographic Alternative",
	176	"gstmerc": "Gauss-Schreiber Transverse Mercator (aka Gauss-Laborde Reunion)",
	177	"tcc": "Transverse Central Cylindrical",
	178	"tcea": "Transverse Cylindrical Equal Area",
	179	"times": "Times",
	180	"tissot": "Tissot Conic",
	181	"tmerc": "Transverse Mercator",
	182	"tpeqd": "Two Point Equidistant",
	183	"tpers": "Tilted perspective",
	184	"ups": "Universal Polar Stereographic",
	185	"urm5": "Urmaev V",
	186	"urmfps": "Urmaev Flat-Polar Sinusoidal",
	187	"utm": "Universal Transverse Mercator (UTM)",
	188	"vandg": "van der Grinten (I)",
	189	"vandg2": "van der Grinten II",
	190	"vandg3": "van der Grinten III",
	191	"vandg4": "van der Grinten IV",
	192	"vitk1": "Vitkovsky I",
	193	"wag1": "Wagner I (Kavraisky VI)",
	194	"wag2": "Wagner II",
	195	"wag3": "Wagner III",
	196	"wag4": "Wagner IV",
	197	"wag5": "Wagner V",
	198	"wag6": "Wagner VI",
	199	"wag7": "Wagner VII",
	200	"webmerc": "Web Mercator / Pseudo Mercator",
	201	"weren": "Werenskiold I",
	202	"wink1": "Winkel I",
	203	"wink2": "Winkel II",
	204	"wintri": "Winkel Tripel",
	205	}
	206
	207
	208	class Proj(_proj.Proj):
	209	"""
	210	performs cartographic transformations (converts from
	211	longitude,latitude to native map projection x,y coordinates and
	212	vice versa) using proj (https://github.com/OSGeo/proj.4/wiki).
	213
	214	A Proj class instance is initialized with proj map projection
	215	control parameter key/value pairs. The key/value pairs can
	216	either be passed in a dictionary, or as keyword arguments,
	217	or as a proj4 string (compatible with the proj command). See
	218	http://www.remotesensing.org/geotiff/proj_list for examples of
	219	key/value pairs defining different map projections.
	220
	221	Calling a Proj class instance with the arguments lon, lat will
	222	convert lon/lat (in degrees) to x/y native map projection
	223	coordinates (in meters). If optional keyword 'inverse' is True
	224	(default is False), the inverse transformation from x/y to
	225	lon/lat is performed. If optional keyword 'errcheck' is True (default is
	226	False) an exception is raised if the transformation is invalid.
	227	If errcheck=False and the transformation is invalid, no
	228	exception is raised and 1.e30 is returned. If the optional keyword
	229	'preserve_units' is True, the units in map projection coordinates
	230	are not forced to be meters.
	231
	232	Works with numpy and regular python array objects, python
	233	sequences and scalars.
	234	"""
	235
	236	def __init__(self, projparams=None, preserve_units=True, **kwargs):
	237	"""
	238	initialize a Proj class instance.
	239
	240	See the proj documentation (https://github.com/OSGeo/proj.4/wiki)
	241	for more information about projection parameters.
	242
	243	Parameters
	244	----------
	245	projparams: int, str, dict, pyproj.CRS
	246	A proj.4 or WKT string, proj.4 dict, EPSG integer, or a pyproj.CRS instnace.
	247	preserve_units: bool
	248	If false, will ensure +units=m.
	249	**kwargs:
	250	proj.4 projection parameters.
	251
	252
	253	Example usage:
	254
	255	>>> from pyproj import Proj
	256	>>> p = Proj(proj='utm',zone=10,ellps='WGS84', preserve_units=False) # use kwargs
	257	>>> x,y = p(-120.108, 34.36116666)
	258	>>> 'x=%9.3f y=%11.3f' % (x,y)
	259	'x=765975.641 y=3805993.134'
	260	>>> 'lon=%8.3f lat=%5.3f' % p(x,y,inverse=True)
	261	'lon=-120.108 lat=34.361'
	262	>>> # do 3 cities at a time in a tuple (Fresno, LA, SF)
	263	>>> lons = (-119.72,-118.40,-122.38)
	264	>>> lats = (36.77, 33.93, 37.62 )
	265	>>> x,y = p(lons, lats)
	266	>>> 'x: %9.3f %9.3f %9.3f' % x
	267	'x: 792763.863 925321.537 554714.301'
	268	>>> 'y: %9.3f %9.3f %9.3f' % y
	269	'y: 4074377.617 3763936.941 4163835.303'
	270	>>> lons, lats = p(x, y, inverse=True) # inverse transform
	271	>>> 'lons: %8.3f %8.3f %8.3f' % lons
	272	'lons: -119.720 -118.400 -122.380'
	273	>>> 'lats: %8.3f %8.3f %8.3f' % lats
	274	'lats: 36.770 33.930 37.620'
	275	>>> p2 = Proj('+proj=utm +zone=10 +ellps=WGS84', preserve_units=False) # use proj4 string
	276	>>> x,y = p2(-120.108, 34.36116666)
	277	>>> 'x=%9.3f y=%11.3f' % (x,y)
	278	'x=765975.641 y=3805993.134'
	279	>>> p = Proj(init="epsg:32667", preserve_units=False)
	280	>>> 'x=%12.3f y=%12.3f (meters)' % p(-114.057222, 51.045)
	281	'x=-1783506.250 y= 6193827.033 (meters)'
	282	>>> p = Proj("+init=epsg:32667")
	283	>>> 'x=%12.3f y=%12.3f (feet)' % p(-114.057222, 51.045)
	284	'x=-5851386.754 y=20320914.191 (feet)'
	285	>>> # test data with radian inputs
	286	>>> p1 = Proj(init="epsg:4214")
	287	>>> x1, y1 = p1(116.366, 39.867)
	288	>>> '{:.3f} {:.3f}'.format(x1, y1)
	289	'2.031 0.696'
	290	>>> x2, y2 = p1(x1, y1, inverse=True)
	291	>>> '{:.3f} {:.3f}'.format(x2, y2)
	292	'116.366 39.867'
	293	"""
	294	self.crs = CRS.from_user_input(projparams if projparams is not None else kwargs)
	295	# make sure units are meters if preserve_units is False.
	296	if not preserve_units and self.crs.is_projected:
	297	projstring = self.crs.to_proj4(4)
	298	projstring = re.sub(r"\s\+units=[\w-]+", "", projstring)
	299	projstring += " +units=m"
	300	self.crs = CRS(projstring)
	301	super(Proj, self).__init__(
	302	cstrencode(self.crs.to_proj4().replace("+type=crs", "").strip())
	303	)
	304
	305	def __call__(self, args, *kw):
	306	# ,lon,lat,inverse=False,errcheck=False):
	307	"""
	308	Calling a Proj class instance with the arguments lon, lat will
	309	convert lon/lat (in degrees) to x/y native map projection
	310	coordinates (in meters). If optional keyword 'inverse' is True
	311	(default is False), the inverse transformation from x/y to
	312	lon/lat is performed. If optional keyword 'errcheck' is True (default is
	313	False) an exception is raised if the transformation is invalid.
	314	If errcheck=False and the transformation is invalid, no
	315	exception is raised and 1.e30 is returned.
	316
	317	Inputs should be doubles (they will be cast to doubles if they
	318	are not, causing a slight performance hit).
	319
	320	Works with numpy and regular python array objects, python
	321	sequences and scalars, but is fastest for array objects.
	322	"""
	323	inverse = kw.get("inverse", False)
	324	errcheck = kw.get("errcheck", False)
	325	# if len(args) == 1:
	326	# latlon = np.array(args[0], copy=True,
	327	# order='C', dtype=float, ndmin=2)
	328	# if inverse:
	329	# _proj.Proj._invn(self, latlon, radians=radians, errcheck=errcheck)
	330	# else:
	331	# _proj.Proj._fwdn(self, latlon, radians=radians, errcheck=errcheck)
	332	# return latlon
	333	lon, lat = args
	334	# process inputs, making copies that support buffer API.
	335	inx, xisfloat, xislist, xistuple = _copytobuffer(lon)
	336	iny, yisfloat, yislist, yistuple = _copytobuffer(lat)
	337	# call proj4 functions. inx and iny modified in place.
	338	if inverse:
	339	self._inv(inx, iny, errcheck=errcheck)
	340	else:
	341	self._fwd(inx, iny, errcheck=errcheck)
	342	# if inputs were lists, tuples or floats, convert back.
	343	outx = _convertback(xisfloat, xislist, xistuple, inx)
	344	outy = _convertback(yisfloat, yislist, xistuple, iny)
	345	return outx, outy
	346
	347	def is_latlong(self):
	348	"""
	349	Returns
	350	-------
	351	bool: True if projection in geographic (lon/lat) coordinates.
	352	"""
	353	warnings.warn("'is_latlong()' is deprecated. Please use 'crs.is_geographic'.")
	354	return self.crs.is_geographic
	355
	356	def is_geocent(self):
	357	"""
	358	Returns
	359	-------
	360	bool: True if projection in geocentric (x/y) coordinates
	361	"""
	362	warnings.warn("'is_geocent()' is deprecated. Please use 'crs.is_geocent'.")
	363	return self.is_geocent
	364
	365	def definition_string(self):
	366	"""Returns formal definition string for projection
	367
	368	>>> Proj('+init=epsg:4326').definition_string()
	369	'proj=longlat datum=WGS84 no_defs ellps=WGS84 towgs84=0,0,0'
	370	>>>
	371	"""
	372	return pystrdecode(self.definition)
	373
	374	def to_latlong_def(self):
	375	"""return the definition string of the geographic (lat/lon)
	376	coordinate version of the current projection"""
	377	# This is a little hacky way of getting a latlong proj object
	378	# Maybe instead of this function the __cinit__ function can take a
	379	# Proj object and a type (where type = "geographic") as the libproj
	380	# java wrapper
	381	return self.crs.to_geodetic().to_proj4(4)
	382
	383	# deprecated : using in transform raised a TypeError in release 1.9.5.1
	384	# reported in issue #53, resolved in #73.
	385	def to_latlong(self):
	386	"""return a new Proj instance which is the geographic (lat/lon)
	387	coordinate version of the current projection"""
	388	return Proj(self.crs.to_geodetic())

+354

-0

pyproj/transformer.py less more

	0	# -- coding: utf-8 --
	1	"""
	2	The transformer module is for performing cartographic transformations.
	3
	4	Copyright (c) 2019 pyproj Contributors.
	5
	6	Permission to use, copy, modify, and distribute this software
	7	and its documentation for any purpose and without fee is hereby
	8	granted, provided that the above copyright notice appear in all
	9	copies and that both the copyright notice and this permission
	10	notice appear in supporting documentation. THE AUTHOR DISCLAIMS
	11	ALL WARRANTIES WITH REGARD TO THIS SOFTWARE, INCLUDING ALL
	12	IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT
	13	SHALL THE AUTHOR BE LIABLE FOR ANY SPECIAL, INDIRECT OR
	14	CONSEQUENTIAL DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM
	15	LOSS OF USE, DATA OR PROFITS, WHETHER IN AN ACTION OF CONTRACT,
	16	NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT OF OR IN
	17	CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE."""
	18
	19	from array import array
	20	from itertools import chain, islice
	21
	22	from pyproj._transformer import _Transformer
	23	from pyproj.compat import cstrencode
	24	from pyproj.utils import _convertback, _copytobuffer
	25
	26	try:
	27	from future_builtins import zip # python 2.6+
	28	except ImportError:
	29	pass # python 3.x
	30
	31
	32	class Transformer(object):
	33	"""
	34	The Transformer class is for facilitating re-using
	35	transforms without needing to re-create them. The goal
	36	is to make repeated transforms faster.
	37
	38	Additionally, it provides multiple methods for initialization.
	39	"""
	40
	41	@staticmethod
	42	def from_proj(proj_from, proj_to):
	43	"""Make a Transformer from a :obj:`pyproj.Proj` or input used to create one.
	44
	45	Parameters
	46	----------
	47	proj_from: :obj:`pyproj.Proj` or input used to create one
	48	Projection of input data.
	49	proj_from: :obj:`pyproj.Proj` or input used to create one
	50	Projection of output data.
	51
	52	Returns
	53	-------
	54	:obj:`pyproj.Transformer`
	55
	56	"""
	57
	58	transformer = Transformer()
	59	transformer._transformer = _Transformer.from_proj(proj_from, proj_to)
	60	return transformer
	61
	62	@staticmethod
	63	def from_crs(crs_from, crs_to):
	64	"""Make a Transformer from a :obj:`pyproj.CRS` or input used to create one.
	65
	66	Parameters
	67	----------
	68	proj_from: :obj:`pyproj.CRS` or input used to create one
	69	Projection of input data.
	70	proj_from: :obj:`pyproj.CRS` or input used to create one
	71	Projection of output data.
	72
	73	Returns
	74	-------
	75	:obj:`pyproj.Transformer`
	76
	77	"""
	78	transformer = Transformer()
	79	transformer._transformer = _Transformer.from_crs(crs_from, crs_to)
	80	return transformer
	81
	82	@staticmethod
	83	def from_pipeline(proj_pipeline):
	84	"""Make a Transformer from a PROJ pipeline string.
	85
	86	https://proj4.org/operations/pipeline.html
	87
	88	Parameters
	89	----------
	90	proj_pipeline: str
	91	Projection pipeline string.
	92
	93	Returns
	94	-------
	95	:obj:`pyproj.Transformer`
	96
	97	"""
	98	transformer = Transformer()
	99	transformer._transformer = _Transformer.from_pipeline(cstrencode(proj_pipeline))
	100	return transformer
	101
	102	def transform(self, xx, yy, zz=None, radians=False):
	103	"""
	104	Transform points between two coordinate systems.
	105
	106	Parameters
	107	----------
	108	xx: scalar or array (numpy or python)
	109	Input x coordinate(s).
	110	yy: scalar or array (numpy or python)
	111	Input y coordinate(s).
	112	zz: scalar or array (numpy or python), optional
	113	Input z coordinate(s).
	114	radians: boolean, optional
	115	If True, will expect input data to be in radians and will return radians
	116	if the projection is geographic. Default is False (degrees). Ignored for
	117	pipeline transformations.
	118
	119
	120	Example:
	121
	122	>>> from pyproj import Transformer
	123	>>> transformer = Transformer.from_crs("epsg:4326", "epsg:3857")
	124	>>> x3, y3 = transformer.transform(33, 98)
	125	>>> "%.3f %.3f" % (x3, y3)
	126	'10909310.098 3895303.963'
	127	>>> pipeline_str = "+proj=pipeline +step +proj=longlat +ellps=WGS84 +step +proj=unitconvert +xy_in=rad +xy_out=deg"
	128	>>> pipe_trans = Transformer.from_pipeline(pipeline_str)
	129	>>> xt, yt = pipe_trans.transform(2.1, 0.001)
	130	>>> "%.3f %.3f" % (xt, yt)
	131	'120.321 0.057'
	132	>>> transproj = Transformer.from_proj({"proj":'geocent', "ellps":'WGS84', "datum":'WGS84'}, '+init=EPSG:4326')
	133	>>> xpj, ypj, zpj = transproj.transform(-2704026.010, -4253051.810, 3895878.820, radians=True)
	134	>>> "%.3f %.3f %.3f" % (xpj, ypj, zpj)
	135	'-2.137 0.661 -20.531'
	136	>>> transprojr = Transformer.from_proj('+init=EPSG:4326', {"proj":'geocent', "ellps":'WGS84', "datum":'WGS84'})
	137	>>> xpjr, ypjr, zpjr = transprojr.transform(xpj, ypj, zpj, radians=True)
	138	>>> "%.3f %.3f %.3f" % (xpjr, ypjr, zpjr)
	139	'-2704026.010 -4253051.810 3895878.820'
	140
	141	"""
	142	# process inputs, making copies that support buffer API.
	143	inx, xisfloat, xislist, xistuple = _copytobuffer(xx)
	144	iny, yisfloat, yislist, yistuple = _copytobuffer(yy)
	145	if zz is not None:
	146	inz, zisfloat, zislist, zistuple = _copytobuffer(zz)
	147	else:
	148	inz = None
	149	# call pj_transform. inx,iny,inz buffers modified in place.
	150	self._transformer._transform(inx, iny, inz, radians)
	151	# if inputs were lists, tuples or floats, convert back.
	152	outx = _convertback(xisfloat, xislist, xistuple, inx)
	153	outy = _convertback(yisfloat, yislist, xistuple, iny)
	154	if inz is not None:
	155	outz = _convertback(zisfloat, zislist, zistuple, inz)
	156	return outx, outy, outz
	157	else:
	158	return outx, outy
	159
	160	def itransform(self, points, switch=False, radians=False):
	161	"""
	162	Iterator/generator version of the function pyproj.Transformer.transform.
	163
	164
	165	Parameters
	166	----------
	167	points: list
	168	List of point tuples.
	169	switch: boolean, optional
	170	If True x, y or lon,lat coordinates of points are switched to y, x
	171	or lat, lon. Default is False.
	172	radians: boolean, optional
	173	If True, will expect input data to be in radians and will return radians
	174	if the projection is geographic. Default is False (degrees). Ignored for
	175	pipeline transformations.
	176
	177
	178	Example:
	179
	180	>>> from pyproj import Transformer
	181	>>> transformer = Transformer.from_crs(4326, 2100)
	182	>>> points = [(22.95, 40.63), (22.81, 40.53), (23.51, 40.86)]
	183	>>> for pt in transformer.itransform(points): '{:.3f} {:.3f}'.format(*pt)
	184	'2221638.801 2637034.372'
	185	'2212924.125 2619851.898'
	186	'2238294.779 2703763.736'
	187	>>> pipeline_str = "+proj=pipeline +step +proj=longlat +ellps=WGS84 +step +proj=unitconvert +xy_in=rad +xy_out=deg"
	188	>>> pipe_trans = Transformer.from_pipeline(pipeline_str)
	189	>>> for pt in pipe_trans.itransform([(2.1, 0.001)]): '{:.3f} {:.3f}'.format(*pt)
	190	'120.321 0.057'
	191	>>> transproj = Transformer.from_proj({"proj":'geocent', "ellps":'WGS84', "datum":'WGS84'}, '+init=EPSG:4326')
	192	>>> for pt in transproj.itransform([(-2704026.010, -4253051.810, 3895878.820)], radians=True): '{:.3f} {:.3f} {:.3f}'.format(*pt)
	193	'-2.137 0.661 -20.531'
	194	>>> transprojr = Transformer.from_proj('+init=EPSG:4326', {"proj":'geocent', "ellps":'WGS84', "datum":'WGS84'})
	195	>>> for pt in transprojr.itransform([(-2.137, 0.661, -20.531)], radians=True): '{:.3f} {:.3f} {:.3f}'.format(*pt)
	196	'-2704214.394 -4254414.478 3894270.731'
	197
	198	"""
	199	it = iter(points) # point iterator
	200	# get first point to check stride
	201	try:
	202	fst_pt = next(it)
	203	except StopIteration:
	204	raise ValueError("iterable must contain at least one point")
	205
	206	stride = len(fst_pt)
	207	if stride not in (2, 3):
	208	raise ValueError("points can contain up to 3 coordinates")
	209
	210	# create a coordinate sequence generator etc. x1,y1,z1,x2,y2,z2,....
	211	# chain so the generator returns the first point that was already acquired
	212	coord_gen = chain(fst_pt, (coords[c] for coords in it for c in range(stride)))
	213
	214	while True:
	215	# create a temporary buffer storage for the next 64 points (64stride8 bytes)
	216	buff = array("d", islice(coord_gen, 0, 64 * stride))
	217	if len(buff) == 0:
	218	break
	219
	220	self._transformer._transform_sequence(stride, buff, switch, radians)
	221
	222	for pt in zip(([iter(buff)] stride)):
	223	yield pt
	224
	225
	226	def transform(p1, p2, x, y, z=None, radians=False):
	227	"""
	228	x2, y2, z2 = transform(p1, p2, x1, y1, z1)
	229
	230	Transform points between two coordinate systems defined by the
	231	Proj instances p1 and p2.
	232
	233	The points x1,y1,z1 in the coordinate system defined by p1 are
	234	transformed to x2,y2,z2 in the coordinate system defined by p2.
	235
	236	z1 is optional, if it is not set it is assumed to be zero (and
	237	only x2 and y2 are returned). If the optional keyword
	238	'radians' is True (default is False), then all input and
	239	output coordinates will be in radians instead of the default
	240	of degrees for geographic input/output projections.
	241
	242	In addition to converting between cartographic and geographic
	243	projection coordinates, this function can take care of datum
	244	shifts (which cannot be done using the __call__ method of the
	245	Proj instances). It also allows for one of the coordinate
	246	systems to be geographic (proj = 'latlong').
	247
	248	x,y and z can be numpy or regular python arrays, python
	249	lists/tuples or scalars. Arrays are fastest. For projections in
	250	geocentric coordinates, values of x and y are given in meters.
	251	z is always meters.
	252
	253	Example usage:
	254
	255	>>> from pyproj import Proj, transform
	256	>>> # projection 1: UTM zone 15, grs80 ellipse, NAD83 datum
	257	>>> # (defined by epsg code 26915)
	258	>>> p1 = Proj(init='epsg:26915', preserve_units=False)
	259	>>> # projection 2: UTM zone 15, clrk66 ellipse, NAD27 datum
	260	>>> p2 = Proj(init='epsg:26715', preserve_units=False)
	261	>>> # find x,y of Jefferson City, MO.
	262	>>> x1, y1 = p1(-92.199881,38.56694)
	263	>>> # transform this point to projection 2 coordinates.
	264	>>> x2, y2 = transform(p1,p2,x1,y1)
	265	>>> '%9.3f %11.3f' % (x1,y1)
	266	'569704.566 4269024.671'
	267	>>> '%9.3f %11.3f' % (x2,y2)
	268	'569722.342 4268814.028'
	269	>>> '%8.3f %5.3f' % p2(x2,y2,inverse=True)
	270	' -92.200 38.567'
	271	>>> # process 3 points at a time in a tuple
	272	>>> lats = (38.83,39.32,38.75) # Columbia, KC and StL Missouri
	273	>>> lons = (-92.22,-94.72,-90.37)
	274	>>> x1, y1 = p1(lons,lats)
	275	>>> x2, y2 = transform(p1,p2,x1,y1)
	276	>>> xy = x1+y1
	277	>>> '%9.3f %9.3f %9.3f %11.3f %11.3f %11.3f' % xy
	278	'567703.344 351730.944 728553.093 4298200.739 4353698.725 4292319.005'
	279	>>> xy = x2+y2
	280	>>> '%9.3f %9.3f %9.3f %11.3f %11.3f %11.3f' % xy
	281	'567721.149 351747.558 728569.133 4297989.112 4353489.645 4292106.305'
	282	>>> lons, lats = p2(x2,y2,inverse=True)
	283	>>> xy = lons+lats
	284	>>> '%8.3f %8.3f %8.3f %5.3f %5.3f %5.3f' % xy
	285	' -92.220 -94.720 -90.370 38.830 39.320 38.750'
	286	>>> # test datum shifting, installation of extra datum grid files.
	287	>>> p1 = Proj(proj='latlong',datum='WGS84')
	288	>>> x1 = -111.5; y1 = 45.25919444444
	289	>>> p2 = Proj(proj="utm",zone=10,datum='NAD27', preserve_units=False)
	290	>>> x2, y2 = transform(p1, p2, x1, y1)
	291	>>> "%s %s" % (str(x2)[:9],str(y2)[:9])
	292	'1402291.0 5076289.5'
	293	>>> from pyproj import CRS
	294	>>> c1 = CRS(proj='latlong',datum='WGS84')
	295	>>> x1 = -111.5; y1 = 45.25919444444
	296	>>> c2 = CRS(proj="utm",zone=10,datum='NAD27')
	297	>>> x2, y2 = transform(c1, c2, x1, y1)
	298	>>> "%s %s" % (str(x2)[:9],str(y2)[:9])
	299	'1402291.0 5076289.5'
	300	>>> pj = Proj(init="epsg:4214")
	301	>>> pjx, pjy = pj(116.366, 39.867)
	302	>>> xr, yr = transform(pj, Proj(4326), pjx, pjy, radians=True)
	303	>>> "%.3f %.3f" % (xr, yr)
	304	'2.031 0.696'
	305	"""
	306	return Transformer.from_proj(p1, p2).transform(x, y, z, radians)
	307
	308
	309	def itransform(p1, p2, points, switch=False, radians=False):
	310	"""
	311	points2 = itransform(p1, p2, points1)
	312	Iterator/generator version of the function pyproj.transform.
	313	Transform points between two coordinate systems defined by the
	314	Proj instances p1 and p2. This function can be used as an alternative
	315	to pyproj.transform when there is a need to transform a big number of
	316	coordinates lazily, for example when reading and processing from a file.
	317	Points1 is an iterable/generator of coordinates x1,y1(,z1) or lon1,lat1(,z1)
	318	in the coordinate system defined by p1. Points2 is an iterator that returns tuples
	319	of x2,y2(,z2) or lon2,lat2(,z2) coordinates in the coordinate system defined by p2.
	320	z are provided optionally.
	321
	322	Points1 can be:
	323	- a tuple/list of tuples/lists i.e. for 2d points: [(xi,yi),(xj,yj),....(xn,yn)]
	324	- a Nx3 or Nx2 2d numpy array where N is the point number
	325	- a generator of coordinates (xi,yi) for 2d points or (xi,yi,zi) for 3d
	326
	327	If optional keyword 'switch' is True (default is False) then x, y or lon,lat coordinates
	328	of points are switched to y, x or lat, lon. If the optional keyword 'radians' is True
	329	(default is False), then all input and output coordinates will be in radians instead
	330	of the default of degrees for geographic input/output projections.
	331
	332
	333	Example usage:
	334
	335	>>> from pyproj import Proj, itransform
	336	>>> # projection 1: WGS84
	337	>>> # (defined by epsg code 4326)
	338	>>> p1 = Proj(init='epsg:4326', preserve_units=False)
	339	>>> # projection 2: GGRS87 / Greek Grid
	340	>>> p2 = Proj(init='epsg:2100', preserve_units=False)
	341	>>> # Three points with coordinates lon, lat in p1
	342	>>> points = [(22.95, 40.63), (22.81, 40.53), (23.51, 40.86)]
	343	>>> # transform this point to projection 2 coordinates.
	344	>>> for pt in itransform(p1,p2,points): '%6.3f %7.3f' % pt
	345	'411200.657 4498214.742'
	346	'399210.500 4487264.963'
	347	'458703.102 4523331.451'
	348	>>> pj = Proj(init="epsg:4214")
	349	>>> pjx, pjy = pj(116.366, 39.867)
	350	>>> for pt in itransform(pj, Proj(4326), [(pjx, pjy)], radians=True): '{:.3f} {:.3f}'.format(*pt)
	351	'2.031 0.696'
	352	"""
	353	return Transformer.from_proj(p1, p2).itransform(points, switch, radians)

-2

setup.py less more

80	80	cythonize_options["compiler_directives"] = {"linetrace": True}
81	81	cythonize_options["annotate"] = True
82	82
83
84	83	proj_libdir = os.environ.get("PROJ_LIBDIR")
85	84	libdirs = []
86	85	if proj_libdir is None:

100	99	os.path.join(BASE_INTERNAL_PROJ_DIR, "lib")
101	100	):
102	101	package_data["pyproj"].append(os.path.join(BASE_INTERNAL_PROJ_DIR, "lib", "*"))
103
104	102
105	103	proj_incdir = os.environ.get("PROJ_INCDIR")
106	104	incdirs = []

134	132	Extension("pyproj._proj", ["pyproj/_proj.pyx"], **ext_options),
135	133	Extension("pyproj._geod", ["pyproj/_geod.pyx"], **ext_options),
136	134	Extension("pyproj._crs", ["pyproj/_crs.pyx"], **ext_options),
	135	Extension(
	136	"pyproj._transformer", ["pyproj/_transformer.pyx"], **ext_options
	137	),
137	138	],
138	139	quiet=True,
139	140	**cythonize_options

-0

sphinx/api/index.rst less more

7	7	proj
8	8	crs
9	9	geod
	10	transformer
10	11	datadir

-12

sphinx/api/proj.rst less more

4	4	pyproj.Proj
5	5	-----------
6	6
7		.. autoclass:: pyproj.Proj
	7	.. autoclass:: pyproj.proj.Proj
8	8	:members:
9	9	:inherited-members:
10	10	:special-members: __init__, __call__
11	11	:show-inheritance:
12
13		pyproj.transform
14		----------------
15
16		.. autofunction:: pyproj.transform
17
18
19		pyproj.itransform
20		-----------------
21
22		.. autofunction:: pyproj.itransform

+19

-0

sphinx/api/transformer.rst less more

	0	Transformer
	1	===========
	2
	3	pyproj.Transformer
	4	------------------
	5
	6	.. autoclass:: pyproj.transformer.Transformer
	7	:members:
	8
	9	pyproj.transform
	10	----------------
	11
	12	.. autofunction:: pyproj.transformer.transform
	13
	14
	15	pyproj.itransform
	16	-----------------
	17
	18	.. autofunction:: pyproj.transformer.itransform⏎

-0

sphinx/history.rst less more

0	0	Change Log
1	1	==========
	2
	3	2.1.0
	4	~~~~~
	5	* Added :class:`pyproj.Transformer` to make repetitive transformations more efficient (issue #187)
	6	* Added fix for using local datumgrids with transform (issue #191)
	7	* Added :class:`pyproj.Transformer.from_pipeline` to support pipeline transformations.
	8	* Added fix for conversion between radians/degrees for transformations (issues #192 & #195)
	9
2	10	2.0.2
3	11	~~~~~
4	12	* add filter for boolean values in dict2string so "no_rot=True" works (issue #183).

-1

sphinx/index.rst less more

2	2
3	3	Python interface to `PROJ.4 <https://proj4.org/>`_.
4	4
5		GitHub Repository: https://github.com/jswhit/pyproj
	5	GitHub Repository: https://github.com/pyproj4/pyproj
6	6
7	7
8	8	.. toctree::

-1

sphinx/installation.rst less more

91	91
92	92	.. code-block:: bash
93	93
94		pip install git+https://github.com/jswhit/pyproj.git
	94	pip install git+https://github.com/pyproj4/pyproj.git
95	95
96	96	From cloned GitHub repo for development:
97	97	^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^

+20

-0

unittest/test_transformer.py less more

	0	import pyproj
	1
	2
	3	def test_tranform_wgs84_to_custom():
	4	custom_proj = pyproj.Proj(
	5	"+proj=geos +lon_0=0.000000 +lat_0=0 +h=35807.414063"
	6	" +a=6378.169000 +b=6356.583984"
	7	)
	8	wgs84 = pyproj.Proj("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
	9	lat, lon = 51.04715, 3.23406
	10	xx, yy = pyproj.transform(wgs84, custom_proj, lon, lat)
	11	assert "{:.3f} {:.3f}".format(xx, yy) == "212.623 4604.975"
	12
	13
	14	def test_transform_wgs84_to_alaska():
	15	lat_lon_proj = pyproj.Proj(init="epsg:4326", preserve_units=False)
	16	alaska_aea_proj = pyproj.Proj(init="epsg:2964", preserve_units=False)
	17	test = (-179.72638, 49.752533)
	18	xx, yy = pyproj.transform(lat_lon_proj, alaska_aea_proj, *test)
	19	assert "{:.3f} {:.3f}".format(xx, yy) == "-1824924.495 330822.800"