How to calculate distance from a point to a line segment, on a sphere?

Question

I have a line segment (great circle part) on earth. The line segment is defined by the coordinates of its ends. Obviously, two points define two line segments, so assume I am interested in the shorter one.

I am given a third point, and I am looking for the (shortest) distance between the line and the point.

All the coordinates are given in longitude\latitude (WGS 84).

How do I calculate the distance?

A solution in any reasonable programming language will do.

Answer 1

Try Distance from a Point to a Great Circle, from Ask Dr. Math. You still need to transform the longitude/latitude to spherical coordinates and scale for the earth's radius, but this seems like a good direction.

Answer 2

The shortest distance between two points on a sphere is the smaller side of the great circle passing through the two points. I'm sure you know this already. There is a similar question here http://www.physicsforums.com/archive/index.php/t-178252.html that may help you model it mathmatically.

I'm not sure how likely you are to get a coded example of this, to be honest.

Answer 3

I'm basically looking for the same thing right now, except that I strictly speaking don't care about having a segment of a great circle, but rather just want the distance to any point on the full circle.

Two links I am currently investigating:

This page mentions "Cross-track distance", which basically seems to be what you are looking for.

Also, in the following thread on the PostGIS mailing list, the attempt seems to (1) determine the closest point on the great circle with the same formula used for line-distance on a 2D-plane (with PostGIS' line_locate_point), and then (2) calculating the distance between that and the third point on a spheroid. I have no idea if mathmatically step (1) is correct, but I would be surprised.

http://postgis.refractions.net/pipermail/postgis-users/2009-July/023903.html

Finally, I just saw that the following linked under "Related":

http://stackoverflow.com/questions/1051723/distance-from-point-to-line-great-circle-functino-not-working-right-need-help

Answer 4

Here's my own solution, based on the idea in ask Dr. Math. I'd be happy to see your feedback.

Disclaimer first. This solution is correct for spheres. Earth isn't a sphere, and the coordinates system (WGS 84) doesn't assume it's a sphere. So this is just an approximation, and I can't really estimate is error. Also, for very small distances, it's probably also possible to get good approximation by assuming everything is a just a coplanar. Again I don't know how "small" the distances have to be.

Now to business. I will call the ends of the lines A,B and the third point C. Basically, the algorithm is to:

1. convert the coordinates first to Cartesian coordinates (with the origin at the center of earth) - e.g. here.

2. Calculate T, the point on the line AB that is nearest to C, using the following 3 vector products:

   G = A x B

   F = C x G

   T = G x F

3. Normalize T and multiply by the radius of earth.

4. Convert T back to longitude\latitude.
5. Calculate the distance between T and C - e.g. here.

These steps are enough if you are looking for the distance between C and the great circle defined by A and B. If like me you are interested in the distance between C and the shorter line segment, you need to take the extra step of verifying that T is indeed on this segment. If it isn't, then necessarily the nearest point is one of the ends A or B - the easiest way is to check which one.

In general terms, the idea behind the three vector products is the following. The first one (G) gives us the the plane of the great circle of A and B (so the plane containing A,B and the origin). The second (F) gives us the great circle the goes through C and is perpendicular to the G. Then T is the intersection of the great circles defined by F and G, brought to the correct length by normalization and multiplication by R.

Here's some partial Java code for doing it.

Finding the nearest point on the great circle. The inputs and output are length-2 arrays. The intermediate arrays are of length 3.

double[] nearestPointGreatCircle(double[] a, double[] b, double c[])
{
    double[] a_ = toCartsian(a);
    double[] b_ = toCartsian(b);
    double[] c_ = toCartsian(c);

    double[] G = vectorProduct(a_, b_);
    double[] F = vectorProduct(c_, G);
    double[] t = vectorProduct(G, F);
    normalize(t);
    multiplyByScalar(t, R_EARTH);
    return fromCartsian(t);
}

Finding the nearest point on the segment:

double[] nearestPointSegment (double[] a, double[] b, double[] c)
{
   double[] t= nearestPointGreatCircle(a,b,c);
   if (onSegment(a,b,t))
     return t;
   return (distance(a,c) < distance(b,c)) ? a : c;
}

This is a simple method of testing if point T, which we know is on the same great circle as A and B, is on the shorter segment of this great circle. However there are more efficient methods to do it:

   boolean onSegment (double[] a, double[] b, double[] t)
   {
     // should be   return distance(a,t)+distance(b,t)==distance(a,b), 
     // but due to rounding errors, we use: 
     return Math.abs(distance(a,b)-distance(a,t)-distance(b,t)) < PRECISION;
   }