Wednesday, 17 June 2015

python - How to get appropriate CRS for a position specified in lat/lon coordinates?


I have a set of city networks whose coordinates are specified in lat/lons. I wish to run an algorithm on each network that treats them as a planar graph. I'm using python and can use pyproj to project the network appropriately, the problem is that I don't know the appropriate CRS for each city, e.g., for my test case, Harare, I define the variable: zimbabwe_reference='epsg:20935'.


What I'd like to do is call a function like crs=getCRS(lat, lon) and get the relevant string to pass on to pyproj.



Answer



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)


Alternatively, create a custom CRS on the town center:


+proj=tmerc +lat_0=(latitude) +lon_0=(longitude) +k=1 +x_0=0 +y_0=0 +ellps=WGS84 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs

No comments:

Post a Comment

arcpy - Changing output name when exporting data driven pages to JPG?

Is there a way to save the output JPG, changing the output file name to the page name, instead of page number? I mean changing the script fo...