5

I'm looking for a smooth way to calculate the distance between two GPS Points, so I get the result like: "You have to go x meters up and y meters to the left - so I can work with a 2d-coordinate system, where I have my position as (0,0) and the other positions is showing the distance in (x, y) in meters from my position.

My idea was to calculate the distance between the points using the haversine formula. (This returns my hypotenuse)

In addition to that, I'm calculating the bearing between this two points. This is my alpha.

With this two values, I wanted to use basic trigonometry functions to resolve my problem.

So I tried to calculate:catheti_1 = sin(alpha) * hypotenuse, catheti_2 = cos(alpha) * hypotenuse.

Maybe I'm doing something wrong, but my results are useless at the moment.

So my question is: How can I calculate the distance in x and y direction between two GPS points?

I'm calculating alpha in the following procedure:

public static double bearingTo(GPSBean point1, GPSBean point2) {
    double lat1 = Math.toRadians(point1.latitude);
    double lat2 = Math.toRadians(point2.latitude);
    double lon1 = Math.toRadians(point1.longitude);
    double lon2 = Math.toRadians(point2.longitude);

    double deltaLong = lon2 - lon1;

    double y = Math.sin(deltaLong) * Math.cos(lat2);
    double x = Math.cos(lat1) * Math.sin(lat2) - Math.sin(lat1)
            * Math.cos(lat2) * Math.cos(deltaLong);
    double bearing = Math.atan2(y, x);

    return (Math.toDegrees(bearing) + 360) % 360;
}
Peter O.
  • 32,158
  • 14
  • 82
  • 96
Frame91
  • 3,670
  • 8
  • 45
  • 89
  • how are you calculating `alpha`? – Barranka Apr 08 '13 at 22:26
  • I've added the method for alpha ;) – Frame91 Apr 08 '13 at 22:28
  • @Borian No. It isn't a duplicate, because I don't need the direct distance between this two points. I know about the harversine formula, but I have to "walk" in the taxicab geometry ;) – Frame91 Apr 08 '13 at 22:30
  • maybe is something obvious, but are you sure you are passing radians to `sin()` and `cos()` (when calculating `catheti_1` and `catheti_2`)? just to discard that – Barranka Apr 08 '13 at 22:34
  • Hey, yep i'm passing radians to sin and cos :/. I'm wondering that the bearing always returns a value around 180 (+-1). Is my formula wrong? – Frame91 Apr 08 '13 at 22:43
  • 1
    You have at least one typo: your `lat2` is using `point2.longitude`. Not going to get to the right answer with that... – Floris Apr 08 '13 at 22:45
  • yep. Thanks! This was clearly wrong. I've changed it, but my results are still wrong – Frame91 Apr 08 '13 at 22:47
  • Okay well. My hyptotenuse^2 is approx. catheti_1^2 and catheti_2^2, so I guess, I'm on the correct way. However, boath values catheti_1 and catheti_2 are negative. I could simply multiply them with -1, but I think I've made a little mistake in my formula? – Frame91 Apr 08 '13 at 22:58
  • No reason why they should be positive - depends on whether point 1 is to the left/right of point 2, and above/below. See my code below for working solution. – Floris Apr 08 '13 at 23:09

3 Answers3

6

I just implemented your code, using approximate coordinates of NYC and Boston as reference points, and implementing the Haversine formula as found at http://www.movable-type.co.uk/scripts/latlong.html (which you didn't show):

long1 = -71.02; lat1 = 42.33;
long2 = -73.94; lat2 = 40.66;

lat1 *=pi/180;
lat2 *=pi/180;
long1*=pi/180;
long2*=pi/180;

dlong = (long2 - long1);
dlat  = (lat2 - lat1);

// Haversine formula:
R = 6371;
a = sin(dlat/2)*sin(dlat/2) + cos(lat1)*cos(lat2)*sin(dlong/2)*sin(dlong/2)
c = 2 * atan2( sqrt(a), sqrt(1-a) );
d = R * c;

When I run this code, I get d = 306, which agrees with the answer from the above site.

For the bearing I get 52 deg - again, close to what the site gave.

Without seeing the rest of your code it's hard to know why your answer is different.

Note: when the two points are close together, you could make all kinds of approximations, but this code should still work - the formula has good numerical stability because it's using the sin of the difference between longitudes, latitudes (rather than the difference of the sin).

Addendum:

Using your code for x, y (in your question), I get sensible values for the distance - agreeing with the "proper" answer to within 120 m (which isn't bad since one is a straight line approximation and the other follows the curvature of the earth). So I think your code is basically OK now you fixed the typo.

Richard
  • 56,349
  • 34
  • 180
  • 251
Floris
  • 45,857
  • 6
  • 70
  • 122
  • Hey thanks. I've tested it myself too, after I corrected the typo and I get some possible values. Maybe it was just the typo ;) – Frame91 Apr 08 '13 at 23:22
2

Use Haversine formula to Calculate distance (in km) between two points specified by latitude/longitude (in numeric degrees)

from: Haversine formula - R. W. Sinnott, "Virtues of the Haversine"

Sky and Telescope, vol 68, no 2, 1984

http://www.census.gov/cgi-bin/geo/gisfaq?Q5.1

Example usage from form:

result.value = LatLon.distHaversine(lat1.value.parseDeg(), long1.value.parseDeg(), * lat2.value.parseDeg(), long2.value.parseDeg());

Javascript : 

LatLon.distHaversine = function(lat1, lon1, lat2, lon2) {
   var R = 6371; // earth's mean radius in km
   var dLat = (lat2-lat1).toRad();
   var dLon = (lon2-lon1).toRad();
   lat1 = lat1.toRad(), lat2 = lat2.toRad();
   var a = Math.sin(dLat/2) * Math.sin(dLat/2) +
   Math.cos(lat1) * Math.cos(lat2) * Math.sin(dLon/2) * Math.sin(dLon/2);
   var c = 2 * Math.atan2(Math.sqrt(a), Math.sqrt(1-a));
   var d = R * c;

   return d;
}
Peter Lapisu
  • 19,915
  • 16
  • 123
  • 179
0

If anybody is interested to have a simpler formula that anyone can understand. Here is mine, it works for Sweden, but you can adapt it to work anywhere by making a more general formula for calculation of longfactor. Hope you can understand even if it is written in an odd language.

<gpsDist lat1,long1,lat2,long2> all parameters in 1/100000 degree.
Example: <getDist 5950928,1327120,5958505,1302241> => 16303
Same at https://gps-coordinates.org/distance-between-coordinates.php => 16.35 KM.

<var $latFactor,1.112>
<function getDist,
-<var $longFactor,<calc 0.638 - ($1/100000-55)*0.0171,3>>
-<var $latDist,<calc ($3-$1)*$latFactor>>
-<var $longDist,<calc ($4-$2)*$longFactor>>
-<sqrt $latDist*$latDist + $longDist*$longDist>
->

/Bertil Friman

Bertil Friman
  • 1