# Tag Info

139

The simplest way to transform coordinates in Python is pyproj, i.e. the Python interface to PROJ.4 library. In fact: from pyproj import Proj, transform inProj = Proj(init='epsg:3857') outProj = Proj(init='epsg:4326') x1,y1 = -11705274.6374,4826473.6922 x2,y2 = transform(inProj,outProj,x1,y1) print x2,y2 returns -105.150271116 39.7278572773 EDIT based on ...

63

While shapely doesn't natively understand coordinate systems, shapely.ops.transform() can do that along with pyproj. If pyproj.Proj can understand your both of your coordinate systems, then it can be made into a function that shapely can transform with. From the shapely docs: pyproj >= 2.1.0 import pyproj from shapely.geometry import Point from shapely....

37

By default the site you linked to uses the Spatial Reference System EPSG 3857 (WGS84 Web Mercator). I found this information here. You can either specify another Spatial Reference System by entering the desired EPSG into the form under Spatial Reference or you can convert the returned coordinates with Python. For instance you can use the GDAL Python ...

26

If the crs of the GeoDataFrame is known (EPSG:4326 unit=degree, here), you don't need Shapely, nor pyproj in your script because GeoPandas uses them). import geopandas as gpd test = gpd.read_file("test_wgs84.shp") print test.crs test.head(2) Now copy your GeoDataFrame and change the projection to a Cartesian system (EPSG:3857, unit= m as in the answer of ...

19

If you're using pyproj2, it's much easier to use a Transformer. Here's an example: import pyproj from shapely.ops import transform project = pyproj.Transformer.from_proj( pyproj.Proj(init='epsg:4326'), # source coordinate system pyproj.Proj(init='epsg:26913')) # destination coordinate system # g1 is a shapley Polygon g2 = transform(project....

14

While not a Shapely solution, using GeoPandas allows for relatively straightforward projection. For example, if we want to convert a shapefile to ESPG 4326: import geopandas as gpd HabModelEnviro = gpd.GeoDataFrame.from_file('data/HabModelEnviro.shp').replace({-999: None}) HabModelEnviroWGS84 = HabModelEnviro.to_crs({'proj':'longlat', 'ellps':'WGS84', '...

12

At the first question, 'epsg:32054' code has feet units. For this reason, it is necessary to use 'preserve_units=True' as parameter in inProj = Proj(init='epsg:32054') line. Now, next code works well: from pyproj import Proj, transform # wisconsing EPSG:32054 # epsg:4326 is for the entire world, wgs 84...not obvious inProj = Proj(init='epsg:32054', ...

11

from pyproj import Proj, transform print(transform(Proj(init='epsg:4326'), Proj(init='epsg:3857'), -0.1285907, 51.50809)) # longitude first, latitude second. # output (meters east of 0, meters north of 0): (-14314.651244750548, 6711665.883938471) The "trick" is to use these shortcuts for Web Mercator (EPSG 3857) and WGS 84 longitude and latitude (EPSG 4326)...

10

You can calculate the UTM zone of each town center from the longitude, starting at zone 1 from -180°E to -174°E. zone=ROUND((183+longitude)/6;0) should calculate that in one step. The EPSG code is 32600+zone for positive latitudes and 32700+zone for negatives. Together in one formula: EPSG=32700-ROUND((45+latitude)/90;0)*100+ROUND((183+longitude)/6;0) ...

10

Both coordinates are the same. As you are in the Southern Hemisphere a False Northing (of +10000000m) is usually applied to eliminate the negatives. The utm package applies the false Northing by itself. For Proj you need to specify it: >>> p = Proj(proj='utm', ellps='WGS84', zone=34, south=True) >>> p(lon,lat) >>> (261762....

10

The +init=epsg:xxxx syntax (or in form of a {'init': 'epsg:xxx'}) is being deprecated in PROJ / pyproj (see here. You are indeed not using it directly, but, the previous releases of GeoPandas (<= 0.6) were still generating this format under the hood (e.g. when reading files, as you can see from the output of df.crs). That is the reason you see this ...

10

It might be related to the proj string you specified. Use CRS instead of the proj string. from pyproj import Proj, CRS #converting lat/long - UTM lat = -22.818009494 long = -47.059235202 crs = CRS.from_epsg(32723) # 32723 -> UTM Zone 23S myProj = Proj(crs) x, y = myProj(long, lat) print (x, y) # OUT: 288652.54449622496 7475152.999493342 What is ...

9

afalciano has the right answer but wanted to include a variant usage of pyproj. It does require you know the proj4 string and is a tiny bit faster. import pyproj p = pyproj.Proj("+proj=merc +lon_0=0 +k=1 +x_0=0 +y_0=0 +a=6378137 +b=6378137 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs") lon, lat = p(x, y, inverse=True) print lat, lon

9

I believe yes. The following ought to work: gdf['geometry'].to_crs({'init': 'epsg:3395'})\ .map(lambda p: p.area / 10**6) This converts the geometry to an equal-area projection, fetches the shapely area (returned in m^2), and maps that to a km^2 (this last step is optional).

9

Simply use (Converting elevations into correct units with pyproj?, Proj4 String for NAD83(2011) / Louisiana South (ftUS), ...) preserve_units=True (as you say, pyproj assumes that your coordinates are in meters, therefore) for northing, easting, up in targets: print "Converting", easting, northing # Transform using pyproj (gives wrong answer) ...

8

The output is not a spatial/coordinate reference system, it's a pair of coordinates. You need to know what the spatial reference is to reproject the coordinates. However, that's not required in this case. Just pass an appropriate output spatial reference to the service and it will return the coordinates in Lon/Lat. Here is the page with output ...

8

For that, you must know the real Python world and the modern way to install modules. 1) The pyproj module needs the compilation of many C libraries and Windows has no compiler by default as in Linux or Mac OS X so you can't install the module with setuptools , easy_install or pip, the traditional way to install modules or unzipping the folder in C:\...

8

The proj4 library has an error identifying the difference between Israel TM Grid 2039 and IG 05/12 6991. Survey if Israel made a horrendous error and updated 2039 to include the new 7 parameter Coordinate Frame transformation, instead of the previous 3 parameter Molodansky, and left 6991 without the new and correct transformation. I've created the correct ...

8

Probably because in lon, lat, alt = 49.74761271, -113.2179781, 0 a lat of -113.2179781 is outside the valid bounds of -90 to 90.

8

If you examine the answer of afalciano in Converting projected coordinates to lat/lon using Python? 1) you define the two projections # original projection p = pyproj.Proj("+proj=stere +lat_0=90 +lat_ts=60 +lon_0=-105 +k=90 +x_0=0 +y_0=0 +a=6371200 +b=6371200 +units=m +no_defs") # resulting projection, WGS84, long, lat outProj =pyproj.Proj(init='epsg:4326')...

7

It looks like you've done everything correctly. You can evaluate the errors from each method by performing the inverse calculations to find the distance given the origin and destination coordinates, then evaluate the residuals of distances. This is a round-trip exercise. # For Vincenty's method: geopy_inv_dist = geopy.distance.vincenty(origin, destination)....

7

Tried the code suggested by Marcel Wilson and it is faster: from pyproj import Proj, transform import time import pyproj # Test 1 - 0.0006158 s start=time.time() inProj = Proj(init='epsg:3857') outProj = Proj(init='epsg:4326') x1,y1 = -11705274.6374,4826473.6922 x2,y2 = transform(inProj,outProj,x1,y1) end=time.time() print(y2,x2) print('%.7f' % (end-start))...

7

This is a great question and I will do my best to answer. To begin, the init style syntax is deprecated (https://pyproj4.github.io/pyproj/stable/gotchas.html#init-auth-auth-code-should-be-replaced-with-auth-auth-code). So, instead of CRS(init="epsg:4544"), you should use CRS("epsg:4544"). I discovered that sometimes when intialising a pyproj Proj ...

7

'<authority>:<code>' in the error means that you need just to use crs='EPSG:4326' instead of crs={'epsg:4326'}.

7

The issue is the stray comma and unnecessary latitude band K. With the proj string, +zone= must be a value between 1-60, so just use 23 without a comma: myProj = Proj("+proj=utm +zone=23 +south +ellps=WGS84 +datum=WGS84 +units=m +no_defs") x, y = myProj(longs, lats) print(x, y) # 288652.54449622496 7475152.999493342 If you are able to use pyproj ...

6

When using pyproj, note the differences from various releases in how it is used to transform data. Here are a few examples using new/old capabilities based on the question: Using pyproj >= 2.2.0 import pyproj print(pyproj.__version__) # 2.4.1 print(pyproj.proj_version_str) # 6.2.1 proj = pyproj.Transformer.from_crs(3857, 4326, always_xy=True) x1, y1 = (-...

6

Projection of this kind of files is sinusoidal. For this one: ftp://ladsweb.nascom.nasa.gov/allData/6/MOD13Q1/2016/129/MOD13Q1.A2016129.h07v06.006.2016147112419.hdf the next code can access to coordinates for 256 values for subDatasets (NDVI values). from osgeo import gdal import struct nameraster = "/home/zeito/Desktop/MOD13Q1.A2016129.h07v06....

6

Yes, simply be sure to reproject your shape in Cylindrical equal-area format with {'proj':'cea'} that preserve area measure. Then you can use .area method of your GeoDataFrame. Your also need to divide by 1000000 because .area method give area in square meters. import geopandas as gpd gdf = gpd.read_file("YOUR_SHAPE_FILE.shp") gdf = gdf['geometry']...

6

Here's some handy code from the Python GDAL/OGR cookbook that will reproject a shapefile. from osgeo import ogr, osr import os in_epsg = 5514 out_epsg = 4326 in_shp = '/path/to/input.shp' out_shp = '/path/to/reprojected.shp' driver = ogr.GetDriverByName('ESRI Shapefile') # input SpatialReference inSpatialRef = osr.SpatialReference() inSpatialRef....

6

Here is even more concise option using geopandas library: import geopandas data = geopandas.read_file(path+file_in) # change CRS to epsg 4326 data = data.to_crs(epsg=4326) # write shp file data.to_file(path+file_out) BTW Kartograph is not developed anymore according to developer git-hub page!

Only top voted, non community-wiki answers of a minimum length are eligible