[GIS] Pyproj: converting WGS84 to Robinson

coordinate systemprojpyprojpython

I have a WGS84 coordinate pair (37.7749295, -122.4194155). I have an image that is a Robinson Map Projection of the world. I'm trying to write a python script to locate my WGS84 coordinate on my Robinson Map.

As an essential step in this, I'm trying to use pyproj, the python-adaption of Proj4, to convert my WGS84 coordinate to a Robinson coordinate. Unfortunately I am new to Proj4, and somewhat confused by the documentation. I hope that someone can help me out.

Best Answer

You can do this using cartopy.

#!/usr/bin/env python

import cartopy
import cartopy.crs as ccrs
import matplotlib.pyplot as plt
from pylab import imread

im = imread('Robinson-projection.jpg')
ax = plt.axes(projection=ccrs.Robinson())
plt.imshow(im, origin='upper', extent=[-17005833.330525, 17005833.330525, -8622512.772008, 8622512.772008], interpolation='nearest')
ax.coastlines(resolution='110m', color='yellow', linewidth=1, alpha=0.7)
plt.plot(-122.4194155, 37.7749295, marker='o', color='red', transform=ccrs.Geodetic())
plt.show()

To test I used a map in the Robinson projection downloaded from Wikipedia.

http://en.wikipedia.org/wiki/List_of_cartographers#mediaviewer/File:Robinson-projection.jpg

The result is shown below.

output map

To convert coordinates between lnglat and WGS84 using pyproj:

#!/usr/bin/env python
import pyproj
crs_from = pyproj.Proj(init='EPSG:4326')
crs_to = pyproj.Proj('+proj=robin +lon_0=0 +x_0=0 +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m +no_defs')
x, y = pyproj.transform(crs_from, crs_to, -122.4194155, 37.7749295)
Related Question