python - How to call gdaltransform with its input?

I want to convert between different coordinates. If I'm in the command line, it's as simple as:

gdaltransform -s_srs EPSG:4326 -t_srs EPSG:32616
-122 46

And I receive the output

-2193988.77788563 5724599.39928024 0

I'm trying to use this in a python program, and it seems like I need the input in the same command. Something like this (although it doesn't work)

gdaltransform -s_srs EPSG:4326 -122 46 -t_srs EPSG:32616

In my python program, I have tried to make separate calls:

command = "gdaltransform -s_srs EPSG:4326 -t_srs EPSG:" + data['crs']
os.system(command)
command = data['lon'] + ' ' + data['lat']

output = subprocess.check_output(command, shell=True)

but get an error:

CalledProcessError: Command '-122 46' returned non-zero exit status 2

I'm not sure of the correct way of doing this.

Answer

You can do this in Python without a call to an external process using the GDAL python bindings.

Here's an example:

from osgeo import osr


src = osr.SpatialReference()
tgt = osr.SpatialReference()
src.ImportFromEPSG(4326)
tgt.ImportFromEPSG(int(data['crs']))

transform = osr.CoordinateTransformation(src, tgt)
coords = transform.TransformPoint(-122, 46)
x,y = coords[0:2]