# How to use a specific transformation for a coordinate projection with Python

I try to use Python to project single coordinates from Gauß-Krueger-coordinate system (EPSG:31467) to ETRS89-UTM (EPSG:25832). I'm using the way Antonio Falciano showed here.

So my script looks like this:

``````from pyproj import Proj, transform

x_in = 3468511.6749226563
y_in = 5428925.356241324

input_proj = Proj(init='epsg:31467')
output_proj = Proj(init='epsg:25832')

x_out, y_out = transform(input_proj, output_proj, x_in, y_in)
print x_out, y_out

returns 468448.880312 5427193.93989
``````

This works for me, but now I want to add a specific geographic (datum) transfomation (DHDN_To_WGS_1984_4_NTv2). How can I manage this. I would prefere to do this without using a GIS but I cannot find a way so far.

After the tip of mkenedy I found out, that I have to use cs2cs and a .gsb-file. I think the following line has to put in my code somehow, but I really have no idea how. If I try, it'll cause a SyntaxError: invalid syntax:

``````cs2cs +proj=tmerc +lat_0=0 +lon_0=9 +k=1 +x_0=3500000 +y_0=0 +ellps=bessel + units=m + nadgrids=BETA2007.gsb +to +proj=utm +ellps=GRS80 +zone=32 +nadgrids=@null x_in x_out
``````

Although i have been reading internet pages about cs2cs for hours now, I don't understand how to use cs2cs correctly. Can someone help me with this?

Now I tried to do an reprojection first like it is done on this site. It looks like this:

``````from pyproj import Proj, transform
import pyproj

# GK3 = EPSG: 31467
input_EPSG = '31467'
# UTM32 = EPSG: 25832
output_EPSG = '25832'
gsb_filename = 'BWTA2017.gsb'
x_in = 3468511.6749226563
y_in = 5428925.356241324

gk = pyproj.Proj("+proj=tmerc +lat_0=0 +lon_0=9 +k=1 +x_0=3500000 +y_0=0 +ellps=bessel +datum=potsdam +nadgrids=" + gsb_filename + "+units=m + no_defs")
x, y = gk(x_in, y_in)
print x, y

returns: 1e+30 1e+30
RuntimeError: latitude or longitude exceeded limits
``````
• Explain please what this specific transfomation means? I found this is a Geographic(datum) transformation, do you have proj4 definition of this or anything else? – Jane Jul 3 at 12:16
• You would like to transform from DHDN (EPSG:4314) to WGS84 (EPSG:4326) -? – Jane Jul 3 at 12:22
• No, sorry, I apologize for the confusion. The input system is the DHDN (EPSG: 31467) and it has to be in ETRS89-UTM (EPSG: 25832) – AnPy Jul 3 at 14:13
• I don't know for sure, but I think you would have to change your Proj definitions to use the proj4 strings and include +nadgrids=dhdnfilename.gsb parameter which means you have to have that grid file downloaded and available. Partial info at proj.maptools.org/gen_parms.html – mkennedy Jul 3 at 16:57
• @AnPy, EPSG: 31467 and EPSG: 25832 are in your code, and you mentioned "This works for you", So what specific transformation you would like to add? – Jane Jul 4 at 8:41

!!!Solved!!! The way that worked in the end:

``````import pyproj

x_in = 3468072.147662042
y_in = 5430007.049974753

input_EPSG = "31467"
output_EPSG = "25832"

def proj_coord(x_in, y_in):
"""
projection process - In this method the projection of a single pair of coordinates takes place.
:param x_in: X - coordinate-value
:param y_in: Y - coordinate-value
:return: x_out, y_out (projected and transformed coordinate)
"""
input_proj = pyproj.Proj(init="EPSG:"+input_EPSG, nadgrids='..\\' + gsb_filename)
output_proj = pyproj.Proj(init="EPSG:"+output_EPSG)
print "Coordinates are currently in EPSG ", input_EPSG

x_out, y_out = pyproj.transform(input_proj, output_proj, x_in, y_in)

print "Now the coordinates are projected to EPSG", output_EPSG
print "\t The new coordinate - X: ", x_out
print "\t The new coordinate - Y: ", y_out
return x_out, y_out

returns:
Coordinates are currently in EPSG  31467
Now the coordinates are projected to EPSG 25832
The new coordinate - X:  468009.548038
The new coordinate - Y:  5428274.49364
``````

Thanks to all of you. Your questions and tips helped bring to find the right ideas and this procedural method.