In my last project I used (standalone) Python for geoprocessing. Since ArcGIS or something like that was not there to help with projecting geodata from one coordinate system to another, I wrote a function which converts well-known WGS1984 to Swiss national grid coordinates (Swissgrid, CH1903 LV03).

Swisstopo has the formulas of an approximate transformation and also funtions in C#, Javascript and PHP, but not in Python. So here is to your disposal:

[sourcecode language=”python” gutter=”false”]

def projWGS1984ToCH1903(x, y):

“””Converts coordinates from WGS1984 to CH1903_LV03 using the

method from http://www.swisstopo.admin.ch/internet/swisstopo/

de/home/products/software/products/skripts.html, returns a

list [x,y]

:param x: x coordinates in degrees in WGS1984 :

:param y: y coordinates in degrees in WGS1984 “””

# Transformation into sexagesimal seconds

x = x * 3600

y = y * 3600

# Latitude and longitude difference to Bern

x_fact = (x – 26782.5) / 10000 # LAMBDA

y_fact = (y – 169028.66) / 10000 # PHI

x = 600072.37 + 211455.93 * x_fact – 10938.51 * x_fact * y_fact – 0.36 * x_fact * y_fact**2 – 44.54 * x_fact**3

y = 200147.07 + 308807.95 * y_fact + 3745.25 * x_fact**2 + 76.63 * y_fact**2 – 194.56 * x_fact**2 * y_fact + 119.79 * y_fact**3

return [x, y][/sourcecode]