如何计算一个球体上一个点到一个线段的距离？

我在地球上有一个线段（大圆圈部分）。 线段由其端点的坐标定义。 显然，两点限定了两条线段，所以假设我对较短的线段感兴趣。

我得到第三点，我正在寻找线和点之间的（最短的）距离。

所有坐标都以经度/纬度（WGS 84）给出。

我如何计算距离？

任何合理的编程语言的解决scheme将做。

这是我自己的解决scheme，基于问Math博士的想法。 我很乐意看到您的反馈。

免责声明。 这个解决scheme对于领域是正确的。 地球不是一个球体，坐标系统（WGS 84）并不认为它是一个球体。 所以这只是一个近似值，我不能真正估计是错误的。 另外，对于非常小的距离，通过假定所有东西都是共面的，也可能得到好的近似值。 再次，我不知道如何“小”的距离必须是。

现在去做生意。 我将调用行A，B和第三点C的结束。基本上，algorithm是：

  1. 首先将坐标转换为笛卡尔坐标（原点位于地球中心） – 例如这里 。

  2. 使用以下3个vector产品计算T，AB上距离C最近的点：

    G = A×B

    F = C×G

    T = G×F

  3. 标准化T并乘以地球半径。

  4. 将T转换回经度\纬度。
  5. 计算T和C之间的距离 – 例如这里 。

如果您正在查找C和由A和B定义的大圆之间的距离，则这些步骤就足够了。如果像我一样，您对C和较短线段之间的距离感兴趣，则需要采取额外步骤来validationT确实在这个部分。 如果不是，那么最近的点就是A或B的最后一个点 – 最简单的方法就是检查哪一个。

一般而言，三个向量产品背后的想法如下。 第一个（G）给我们A和B的大圆的平面（所以包含A，B和原点的平面）。 第二个（F）给出我们通过C的大圆，并且垂直于G.然后T是由F和G定义的大圆的交点，通过归一化和乘以R得到正确的长度。

这里有一些不完整的Java代码。

find大圆上的最近点。 input和输出是长度为2的数组。 中间arrays的长度为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); }

查找细分受众群最近的点： 

 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; }

这是一个简单的testing方法，如果我们知道与A和B在同一个大圆上的T点在这个大圆的较短部分上。 但是有更有效的方法来做到这一点： 

  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; }

尝试从math博士算起，尝试从一个点到一个大圆的距离。 您仍然需要将经度/纬度转换为球形坐标，并将地球半径缩放，但这似乎是一个好方向。

这是接受答案的完整代码作为ideone小提琴（在这里find）： 

 import java.util.*; import java.lang.*; import java.io.*; /* Name of the class has to be "Main" only if the class is public. */ class Ideone { private static final double _eQuatorialEarthRadius = 6378.1370D; private static final double _d2r = (Math.PI / 180D); private static double PRECISION = 0.1; // Haversine Algorithm // source: http://stackoverflow.com/questions/365826/calculate-distance-between-2-gps-coordinates private static double HaversineInM(double lat1, double long1, double lat2, double long2) { return (1000D * HaversineInKM(lat1, long1, lat2, long2)); } private static double HaversineInKM(double lat1, double long1, double lat2, double long2) { double dlong = (long2 - long1) * _d2r; double dlat = (lat2 - lat1) * _d2r; double a = Math.pow(Math.sin(dlat / 2D), 2D) + Math.cos(lat1 * _d2r) * Math.cos(lat2 * _d2r) * Math.pow(Math.sin(dlong / 2D), 2D); double c = 2D * Math.atan2(Math.sqrt(a), Math.sqrt(1D - a)); double d = _eQuatorialEarthRadius * c; return d; } // Distance between a point and a line public static void pointLineDistanceTest() { //line //double [] a = {50.174315,19.054743}; //double [] b = {50.176019,19.065042}; double [] a = {52.00118, 17.53933}; double [] b = {52.00278, 17.54008}; //point //double [] c = {50.184373,19.054657}; double [] c = {52.008308, 17.542927}; double[] nearestNode = nearestPointGreatCircle(a, b, c); System.out.println("nearest node: " + Double.toString(nearestNode[0]) + "," + Double.toString(nearestNode[1])); double result = HaversineInM(c[0], c[1], nearestNode[0], nearestNode[1]); System.out.println("result: " + Double.toString(result)); } // source: http://stackoverflow.com/questions/1299567/how-to-calculate-distance-from-a-point-to-a-line-segment-on-a-sphere private static 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); return fromCartsian(multiplyByScalar(normalize(t), _eQuatorialEarthRadius)); } @SuppressWarnings("unused") private static double[] nearestPointSegment (double[] a, double[] b, double[] c) { double[] t= nearestPointGreatCircle(a,b,c); if (onSegment(a,b,t)) return t; return (HaversineInKM(a[0], a[1], c[0], c[1]) < HaversineInKM(b[0], b[1], c[0], c[1])) ? a : b; } private static 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(HaversineInKM(a[0], a[1], b[0], b[1])-HaversineInKM(a[0], a[1], t[0], t[1])-HaversineInKM(b[0], b[1], t[0], t[1])) < PRECISION; } // source: http://stackoverflow.com/questions/1185408/converting-from-longitude-latitude-to-cartesian-coordinates private static double[] toCartsian(double[] coord) { double[] result = new double[3]; result[0] = _eQuatorialEarthRadius * Math.cos(Math.toRadians(coord[0])) * Math.cos(Math.toRadians(coord[1])); result[1] = _eQuatorialEarthRadius * Math.cos(Math.toRadians(coord[0])) * Math.sin(Math.toRadians(coord[1])); result[2] = _eQuatorialEarthRadius * Math.sin(Math.toRadians(coord[0])); return result; } private static double[] fromCartsian(double[] coord){ double[] result = new double[2]; result[0] = Math.toDegrees(Math.asin(coord[2] / _eQuatorialEarthRadius)); result[1] = Math.toDegrees(Math.atan2(coord[1], coord[0])); return result; } // Basic functions private static double[] vectorProduct (double[] a, double[] b){ double[] result = new double[3]; result[0] = a[1] * b[2] - a[2] * b[1]; result[1] = a[2] * b[0] - a[0] * b[2]; result[2] = a[0] * b[1] - a[1] * b[0]; return result; } private static double[] normalize(double[] t) { double length = Math.sqrt((t[0] * t[0]) + (t[1] * t[1]) + (t[2] * t[2])); double[] result = new double[3]; result[0] = t[0]/length; result[1] = t[1]/length; result[2] = t[2]/length; return result; } private static double[] multiplyByScalar(double[] normalize, double k) { double[] result = new double[3]; result[0] = normalize[0]*k; result[1] = normalize[1]*k; result[2] = normalize[2]*k; return result; } public static void main(String []args){ System.out.println("Hello World"); Ideone.pointLineDistanceTest(); } }

它适用于评论数据： 

 //line double [] a = {50.174315,19.054743}; double [] b = {50.176019,19.065042}; //point double [] c = {50.184373,19.054657};

最近的节点是：50.17493121381319,19.05846668493702

但是我有这个数据的问题： 

 double [] a = {52.00118, 17.53933}; double [] b = {52.00278, 17.54008}; //point double [] c = {52.008308, 17.542927};

最近的节点是：52.00834987257176,17.542691313436357这是错误的。

我认为由两点指定的线不是封闭的线段。

如果有人需要它这是loleksy答案移植到C＃ 

  private static double _eQuatorialEarthRadius = 6378.1370D; private static double _d2r = (Math.PI / 180D); private static double PRECISION = 0.1; // Haversine Algorithm // source: http://stackoverflow.com/questions/365826/calculate-distance-between-2-gps-coordinates private static double HaversineInM(double lat1, double long1, double lat2, double long2) { return (1000D * HaversineInKM(lat1, long1, lat2, long2)); } private static double HaversineInKM(double lat1, double long1, double lat2, double long2) { double dlong = (long2 - long1) * _d2r; double dlat = (lat2 - lat1) * _d2r; double a = Math.Pow(Math.Sin(dlat / 2D), 2D) + Math.Cos(lat1 * _d2r) * Math.Cos(lat2 * _d2r) * Math.Pow(Math.Sin(dlong / 2D), 2D); double c = 2D * Math.Atan2(Math.Sqrt(a), Math.Sqrt(1D - a)); double d = _eQuatorialEarthRadius * c; return d; } // Distance between a point and a line static double pointLineDistanceGEO(double[] a, double[] b, double[] c) { double[] nearestNode = nearestPointGreatCircle(a, b, c); double result = HaversineInKM(c[0], c[1], nearestNode[0], nearestNode[1]); return result; } // source: http://stackoverflow.com/questions/1299567/how-to-calculate-distance-from-a-point-to-a-line-segment-on-a-sphere private static 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); return fromCartsian(multiplyByScalar(normalize(t), _eQuatorialEarthRadius)); } private static double[] nearestPointSegment (double[] a, double[] b, double[] c) { double[] t= nearestPointGreatCircle(a,b,c); if (onSegment(a,b,t)) return t; return (HaversineInKM(a[0], a[1], c[0], c[1]) < HaversineInKM(b[0], b[1], c[0], c[1])) ? a : b; } private static bool 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(HaversineInKM(a[0], a[1], b[0], b[1])-HaversineInKM(a[0], a[1], t[0], t[1])-HaversineInKM(b[0], b[1], t[0], t[1])) < PRECISION; } // source: http://stackoverflow.com/questions/1185408/converting-from-longitude-latitude-to-cartesian-coordinates private static double[] toCartsian(double[] coord) { double[] result = new double[3]; result[0] = _eQuatorialEarthRadius * Math.Cos(deg2rad(coord[0])) * Math.Cos(deg2rad(coord[1])); result[1] = _eQuatorialEarthRadius * Math.Cos(deg2rad(coord[0])) * Math.Sin(deg2rad(coord[1])); result[2] = _eQuatorialEarthRadius * Math.Sin(deg2rad(coord[0])); return result; } private static double[] fromCartsian(double[] coord){ double[] result = new double[2]; result[0] = rad2deg(Math.Asin(coord[2] / _eQuatorialEarthRadius)); result[1] = rad2deg(Math.Atan2(coord[1], coord[0])); return result; } // Basic functions private static double[] vectorProduct (double[] a, double[] b){ double[] result = new double[3]; result[0] = a[1] * b[2] - a[2] * b[1]; result[1] = a[2] * b[0] - a[0] * b[2]; result[2] = a[0] * b[1] - a[1] * b[0]; return result; } private static double[] normalize(double[] t) { double length = Math.Sqrt((t[0] * t[0]) + (t[1] * t[1]) + (t[2] * t[2])); double[] result = new double[3]; result[0] = t[0]/length; result[1] = t[1]/length; result[2] = t[2]/length; return result; } private static double[] multiplyByScalar(double[] normalize, double k) { double[] result = new double[3]; result[0] = normalize[0]*k; result[1] = normalize[1]*k; result[2] = normalize[2]*k; return result; }

对于长达几千米的距离，我会简化从飞机到飞机的问题。 然后，这个问题很简单，因为可以使用一个简单的三angular形计算：

我们有A点和B点，并在AB线寻找距离X. 然后： 

 Location a; Location b; Location x; double ax = a.distanceTo(x); double alfa = (Math.abs(a.bearingTo(b) - a.bearingTo(x))) / 180 * Math.PI; double distance = Math.sin(alfa) * ax;

球体上两点之间的最短距离是通过两点的大圆的较小侧。 我相信你已经知道了 在这里有一个类似的问题http://www.physicsforums.com/archive/index.php/t-178252.html ，可以帮助你用math模型。

我不确定你有多大可能得到一个这样的代码示例，说实话。

我现在基本上都在寻找同样的东西，只不过我严格地说并不关心是否有一个大圆的一部分，而只是想把距离定在整个圆上的任何一点。

我正在调查的两个链接：

这个页面提到“跨轨道距离”，这基本上是你正在寻找的东西。

另外，在PostGIS邮件列表中的以下线程中，尝试似乎（1）使用与2D平面上的线距相同的公式（使用PostGIS'line_locate_point）确定大圆上的最近点，然后（2）计算球体与第三点之间的距离。 我不知道math步骤（1）是否正确，但我会感到惊讶。

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

最后，我刚才看到下面链接的“相关”：

距离点到线大圆圈function不能正常工作。

Interesting Posts