I am trying to convert lat,lon to x,y coordinates in France, using the pyproj library and what I am doing is obviously wrong.

Here is some of my code:

import pyproj
import math

PROJ='+proj=utm +zone=31, +north +ellps=WGS84 +datum=WGS84 +units=m +no_defs'

def LatLon_To_XY(Lat,Lon):

def XY_To_LatLon(x,y):

def distance(x1,y1,x2,y2):

Problem is when I am validating this computing the distance between Paris and Toulouse, I get 816km where the distance is supposed to be around 580km.

I guess I have wrong parameters for the projection definition.

Here are the values displayed:

  • For Paris lat,lon (48.856614, 2.3522219) and x,y: 6253348.594037977 374425.11288885673

  • For Toulouse lat,lon (43.604652, 1.444209) and x,y: 5453729.038474649 210744.18988320036

  • thanks. so I can trust my computed x and y ? – Stéphane Oct 15 '18 at 12:45
  • made it into an answer...I do get some slighty different x & y when transforming in PostGIS, but the relative distance is the same – geozelot Oct 15 '18 at 12:50
  • I think the Latitude and Longitude are swapped in your code, try with (lon,lat). Normally distances in UTM within a specific zone aren't distorted much more than 1 in 1000 – FSimardGIS Oct 15 '18 at 16:23

Nope, there's nothing wrong with your calculations; you are looking at projection distortions!

I´d suggest to use the built in Geod class in pyproj for geodetic computations. e.g.

import pyproj
import math

P = pyproj.Proj(proj='utm', zone=31, ellps='WGS84', preserve_units=True)
G = pyproj.Geod(ellps='WGS84')

def LatLon_To_XY(Lat,Lon):

def XY_To_LatLon(x,y):
    return P(x,y,inverse=True)    

def distance(Lat1, Lon1, Lat2, Lon2):
    return G.inv(Lon1, Lat1, Lon2, Lat2)[2]

or a framework that implements geodetic algebra, e.g. vincenty, GeoPy...

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.