2

I have a set of lat/lon coordinates for different nodes.
Since I have to represent these positions in a simulator, I have to convert them to cartesian coordinates.
For this I use the python package pyproj.

lat1 = 50.0
lon1 = 10.0
lat2 = 50.5
lon2 = 10.5

x1, y1, z1 = wgs84_to_utm32(lat1, lon1, 0.0)
x2, y2, z2 = wgs84_to_utm32(lat2, lon2, 0.0)

print(distance_haversine(lat1, lon1, lat2, lon2))
print(distance_cartesian(x1, y1, x2, y2))

the ouput is:

66012.5130481
102485.874713

that is a difference of over 36 km.

So my question is, how can I convert lat/lon coordinates so that the distances are preserved. I don't care for minimal error.

EDIT:

#
# Convert WGS84 (Geographic) coordinates to UTM32 (Germany)
#
def wgs84_to_utm(lat, lon, alt):
   wgs84 = Proj(init='epsg:4326')
   utm_de = Proj(init='epsg:32632')  # utm germany = 32

   return transform(wgs84, utm_de, lat, lon, alt)

EDIT2:
Okay to be clear, I know that I am trying to compare a distance on a sphere to a distance on a flat surface.

But since I am simulating WLAN nodes, the distance between those nodes is crucial. But I have no other information than their lat/lon coordinates.

How would I go about representing those lat/lon coordinates on a flat surface, so that the distances are preserved?

EDIT3:

def distance_haversine(lat1, lon1, lat2, lon2):
    # approximate radius of earth in km
    R = 6373.0 * 1000

    lat1 = radians(lat1)
    lon1 = radians(lon1)
    lat2 = radians(lat2)
    lon2 = radians(lon2)

    dlon = lon2 - lon1
    dlat = lat2 - lat1

    a = sin(dlat / 2) ** 2 + cos(lat1) * cos(lat2) * sin(dlon / 2) ** 2
    c = 2 * atan2(sqrt(a), sqrt(1 - a))

    distance = R * c

    return distance


def distance_cartesian(x1, y1, x2, y2):
    dx = x1 - x2
    dy = y1 - y2

    return sqrt(dx * dx + dy * dy)
stulleman
  • 163
  • 2
  • 14
  • What does wgs84_to_utm32 do? Also your two points are the same - surely the distance between them is zero? – chasmani Oct 20 '17 at 15:35
  • FWIW, using [geographiclib](https://pypi.python.org/pypi/geographiclib/1.49), the Geodesic.WGS84 distance between those locations is 66067.75789575798 metres. – PM 2Ring Oct 20 '17 at 15:36
  • @chasmani I added the code for the function. What do you mean with, the points are the same? They're not. – stulleman Oct 20 '17 at 15:41
  • Maybe this helps: https://stackoverflow.com/questions/31365848/how-to-convert-longitude-latitude-elevation-to-cartesian-coordinates – mrCarnivore Oct 20 '17 at 15:41
  • @mrCarnivore no, because the transformation is not the problem. I just need a transformation where distances are preserved. – stulleman Oct 20 '17 at 15:44
  • @stulleman: don't you think missing the real elevation could be the explanation? – mrCarnivore Oct 20 '17 at 15:49
  • Haversine and cartesian calculations will give different results, because Haversine distance is across the surface of the sphere, whereas cartesian is straight-line. But that would only be a minimal difference over those points – chasmani Oct 20 '17 at 15:51
  • I added another edit – stulleman Oct 20 '17 at 15:55
  • I checked the result for Haversine and that looks right to me. What is the distance_cartesian method? Can you show that? – chasmani Oct 20 '17 at 16:01
  • I don't know if you get notifications when I edit my post. So I edited my post – stulleman Oct 20 '17 at 16:10
  • Your utm coordinates were coming out wrong, not sure why. You could fix that within pyproj, or use the utm module instead, as I did in my answer – chasmani Oct 20 '17 at 16:21

1 Answers1

2

There is a mistake somewhere in the conversion to utm. Try the utm module instead

import utm

lat1 = 50.0
lon1 = 10.0
lat2 = 50.5
lon2 = 10.5

x1, y1, z1, u = utm.from_latlon(lat1, lon1)
x2, y2, z2, u = utm.from_latlon(lat2, lon2)

print(distance_haversine(lat1, lon1, lat2, lon2))
print(distance_cartesian(x1, y1, x2, y2))

Output:

66012.51304805529
66047.84250743943
chasmani
  • 2,362
  • 2
  • 23
  • 35
  • This still leads to a bit of error, but much better than my error. Do you have any idea why this works? – stulleman Oct 20 '17 at 16:26
  • Yeah. The UTM co-ordinates from your original wgs84_to_utm method were coming out wrong. Everything else was fine. I don't know anything about the pyproj module so I just used the utm module instead – chasmani Oct 20 '17 at 16:28
  • Okay so most probably I'm using pyproj wrong. Thank you, I'll look into what I'm doing wrong and post it here if I find anything. – stulleman Oct 20 '17 at 16:30