pyproj

PROJ is a generic coordinate transformation software, that transforms coordinates from one coordinate reference system (CRS) to another. This includes cartographic projections as well as geodetic transformations.

PROJ.4 was started in the 1980s, and Version 4 has been around since 1995. Now there is a version 5 it is referred to as PROJ (the history of the project can be found on Wikipedia). In version 4 everything was converted by passing through WGS84. If you are interested in projections, have a read of https://proj4.org/development/migration.html for more details. PROJ4 could only guarantee metre accuracy for many transformations, so PROJ5 provides much higher accuracy. A fund raising effort was started to improve the way projections are handled at https://gdalbarn.com/ - as projections are critical to the entire open source stack.

It is used throughout GIS software, for example in ArcGIS, QGIS, and MapServer.

pyproj is a Cython wrapper to provide python interfaces to PROJ. Access to PROJ is also available through many other languages (e.g. C#, R).

Links:

• proj sourcecode at Github: https://github.com/OSGeo/proj.4/
• proj documentation at https://proj4.org/
• pyproj documentation at https://jswhit.github.io/pyproj/
• pyproj sourcode at https://github.com/jswhit/pyproj

Projection Definitions

To install pyproj you can run the following command from the command line (Windows or Linux). It is however preinstalled on OGGeoLive.

pip install pyproj

http://epsg.io/ is a very useful site to find project definitions. For example to see projections relating to Ireland you can run the following search: http://epsg.io/?q=ireland The Irish Transverse Mercator projection is found at http://epsg.io/2157 - note the structure URLs which make these easy to access programatically.

Each projection on the site hasPROJ.4 definitions

# Irish Transverse Mercator http://epsg.io/2157.proj4
+proj=tmerc +lat_0=53.5 +lon_0=-8 +k=0.99982 +x_0=600000 +y_0=750000 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
# Irish National Grid https://epsg.io/29902.proj4
+proj=tmerc +lat_0=53.5 +lon_0=-8 +k=1.000035 +x_0=200000 +y_0=250000 +ellps=mod_airy +towgs84=482.5,-130.6,564.6,-1.042,-0.214,-0.631,8.15 +units=m +no_defs

Session Code

Note - when saving the script, don’t name it pyproj.py as this will clash with the name of the pyproj library that needs to be imported and cause the script to fail!

"""
This script does the following:

+ Downloads a list of Irish towns and cities from https://data.gov.ie in a CSV format
+ Takes the ITM coordinates for each settlement and converts them to Lat/Lng

The dataset details can be found at https://data.gov.ie/dataset/settlements-ungeneralised-osi-national-statistical-boundaries/resource/c7d0e002-8fcf-4e31-820e-5d628f41eea2
Note there is no metadata indicating the coordinate projection!

As the Irish settlement names include special characters such as a fada the text needs to be encoded in UTF (Universal Text Format). 
Otherwise errors such as the following will occur: 

UnicodeEncodeError: 'ascii' codec can't encode character u'\ufeff' in position 0: ordinal not in range(128)
"""

import csv # a standard Python library
import pyproj
import requests # a third party library also installed on OSGeoLive

# show the location of the folder where all the projection definitions are stored
# On Windows: C:\VirtualEnvs\pytraining\lib\site-packages\pyproj\data
# On Linux /usr/share/proj/epsg
print("Projection definitions folder location : {}".format(pyproj.pyproj_datadir))

url = "http://data-osi.opendata.arcgis.com/datasets/059f9de9770147a2a18c5e5d710839ad_3.csv"
r = requests.get(url)

# encode the text in UTF (see above)
txt = r.text.encode(r.encoding) # r.encoding is "utf-8"

# strip any whitespace at the start and end of the text which could cause extra blank records
txt = txt.strip()

# now split the text by new line characters \n, to return a list of lines
rows = txt.split("\n") 

# get an "iterable" from the text, with a record for each row in the list
records = csv.reader(rows)

next(records, None)  # skip the headers at the start of the CSV file

# we'll create some projection objects here, outside of the loop
proj1 = pyproj.Proj(init='epsg:2157') # ITM
proj2 = pyproj.Proj(init='epsg:4326') # WGS84

# loop through the records
for r in records:
    name = r[2] # access values by list index
    x = r[3]
    y = r[4]
    # print out the ITM coords
    print("{} (X: {} Y: {})".format(name, x, y))

    # now let's reproject to Lat Lng
    x2, y2 = pyproj.transform(proj1, proj2, x, y)

    # print out the new coords
    print("{} (Lng: {} Lat: {})".format(name, x2, y2))