In [3]:

import gdal
import osr

gdal.UseExceptions()
gdal.VersionInfo()

Out[3]:

'3000000'

Define projections¶

In [4]:

srs_4326 = osr.SpatialReference()
srs_4326.ImportFromEPSG(4326)

Out[4]:

In [5]:

srs_3857 = osr.SpatialReference()
srs_3857.ImportFromEPSG(3857)

Out[5]:

Define transformations¶

In [6]:

ct_3857_to_4326 = osr.CoordinateTransformation(srs_3857, srs_4326)
ct_4326_to_3857 = osr.CoordinateTransformation(srs_4326, srs_3857)

Forward transform (x, y)¶

In [7]:

lat = 50.0
lon = -120.0

# Incorrect output, and different compared to GDAL 2.4
mapx, mapy, z = ct_4326_to_3857.TransformPoint(lon, lat)
mapx, mapy

Out[7]:

(inf, inf)

Forward transform (y, x)¶

In [11]:

lat = 50.0
lon = -120.0

# this output is equivalent to the cmd (gdaltransform.exe) (but (x,y) -> (y,x))
mapx, mapy, z = ct_4326_to_3857.TransformPoint(lat, lon)
mapx, mapy

Out[11]:

(-13358338.895192828, 6446275.841017158)

Backward transform (x, y)¶

In [14]:

mapx = -13358338.895192828
mapy = 6446275.841017158

# Wrong output order (?), but correct values # (x,y) -> (y,x)
ct_3857_to_4326.TransformPoint(mapx, mapy)

Out[14]:

(49.99999999999999, -119.99999999999999, 0.0)

Backward transform (y, x)¶

In [15]:

mapx = -13358338.895192828
mapy = 6446275.841017158

# Incorrect results
ct_3857_to_4326.TransformPoint(mapy, mapx)

Out[15]:

(-75.95934455387825, 57.90788113636135, 0.0)

CMD output reference¶

In [17]:

# (py37) C:\Users\Rutger>gdaltransform --version
# GDAL 3.0.0, released 2019/05/05

# (py37) C:\Users\Rutger>gdaltransform -s_srs EPSG:4326 -t_srs EPSG:3857
# -120 50
# -13358338.8951928 6446275.84101716 0

CoordinateTransform docstring¶

In [18]:

# docstring suggests (x,y,z)
ct_4326_to_3857.TransformPoint?

Signature: ct_4326_to_3857.TransformPoint(*args)
Docstring:
TransformPoint(CoordinateTransformation self, double [3] inout)
TransformPoint(CoordinateTransformation self, double [4] inout)
TransformPoint(CoordinateTransformation self, double x, double y, double z=0.0)
TransformPoint(CoordinateTransformation self, double x, double y, double z, double t)
File:      c:\miniconda3\envs\py37\lib\site-packages\osgeo\osr.py
Type:      method

In [ ]: