我目前正在处理结合精确高度测量的 GPS 数据。 我想计算两个连续点之间的距离。有很多 关于使用 WGS84 椭球等计算两点之间距离的信息。
但是，我没有找到任何将 海拔 变化考虑在内的信息 距离计算。
有人知道描述这种方法的一些网站、论文、书籍等吗？ 谢谢
编辑：Sql Server 2008 地理扩展在计算距离时也会忽略海拔信息。
【问题讨论】：
标签： geometry gps geospatial
我使用开始和结束高度的平均值作为恒定高度实现了 WGS84 距离函数。如果您确定沿您的路径会有相对较小的高度变化，这可以很好地工作（误差与您的两个 LLA 点的高度差有关）。
这是我的代码 (C#)：
/// <summary>
/// Gets the geodesic distance between two pathpoints in the current mode's coordinate system
/// </summary>
/// <param name="point1">First point</param>
/// <param name="point2">Second point</param>
/// <param name="mode">Coordinate mode that both points are in</param>
/// <returns>Distance between the two points in the current coordinate mode</returns>
public static double GetGeodesicDistance(PathPoint point1, PathPoint point2, CoordMode mode) {
// calculate proper geodesics for LLA paths
if (mode == CoordMode.LLA) {
// meeus approximation
double f = (point1.Y + point2.Y) / 2 * LatLonAltTransformer.DEGTORAD;
double g = (point1.Y - point2.Y) / 2 * LatLonAltTransformer.DEGTORAD;
double l = (point1.X - point2.X) / 2 * LatLonAltTransformer.DEGTORAD;
double sinG = Math.Sin(g);
double sinL = Math.Sin(l);
double sinF = Math.Sin(f);
double s, c, w, r, d, h1, h2;
// not perfect but use the average altitude
double a = (LatLonAltTransformer.A + point1.Z + LatLonAltTransformer.A + point2.Z) / 2.0;
sinG *= sinG;
sinL *= sinL;
sinF *= sinF;
s = sinG * (1 - sinL) + (1 - sinF) * sinL;
c = (1 - sinG) * (1 - sinL) + sinF * sinL;
w = Math.Atan(Math.Sqrt(s / c));
r = Math.Sqrt(s * c) / w;
d = 2 * w * a;
h1 = (3 * r - 1) / 2 / c;
h2 = (3 * r + 1) / 2 / s;
return d * (1 + (1 / LatLonAltTransformer.RF) * (h1 * sinF * (1 - sinG) - h2 * (1 - sinF) * sinG));
}
PathPoint diff = new PathPoint(point2.X - point1.X, point2.Y - point1.Y, point2.Z - point1.Z, 0);
return Math.Sqrt(diff.X * diff.X + diff.Y * diff.Y + diff.Z * diff.Z);
}
在实践中，我们发现海拔差异很少会产生大的差异，我们的路径通常长 1-2 公里，海拔变化大约 100 米，与使用未修改的 WGS84 椭球相比，我们看到平均变化约 5 米.
编辑：
此外，如果您确实预计会有较大的高度变化，您可以将 WGS84 坐标转换为 ECEF（地心地球固定）并评估直线路径，如我的函数底部所示。将点转换为 ECEF 很简单：
/// <summary>
/// Converts a point in the format (Lon, Lat, Alt) to ECEF
/// </summary>
/// <param name="point">Point as (Lon, Lat, Alt)</param>
/// <returns>Point in ECEF</returns>
public static PathPoint WGS84ToECEF(PathPoint point) {
PathPoint outPoint = new PathPoint(0);
double lat = point.Y * DEGTORAD;
double lon = point.X * DEGTORAD;
double e2 = 1.0 / RF * (2.0 - 1.0 / RF);
double sinLat = Math.Sin(lat), cosLat = Math.Cos(lat);
double chi = A / Math.Sqrt(1 - e2 * sinLat * sinLat);
outPoint.X = (chi + point.Z) * cosLat * Math.Cos(lon);
outPoint.Y = (chi + point.Z) * cosLat * Math.Sin(lon);
outPoint.Z = (chi * (1 - e2) + point.Z) * sinLat;
return outPoint;
}
编辑 2：
有人问我代码中的其他一些变量：
// RF is the eccentricity of the WGS84 ellipsoid
public const double RF = 298.257223563;
// A is the radius of the earth in meters
public const double A = 6378137.0;
LatLonAltTransformer 是我用来从 LatLonAlt 坐标转换为 ECEF 坐标的类，并定义了上面的常量。
【讨论】：
您可能不关心大型 2D 距离分离的高度。因此，如果您获得的距离超过 20（或 50）公里，那么谁在乎海拔差异（取决于您的需求情况）。在 20 公里以下，输入简单的毕达哥拉斯除了高度差。顺利喂入。
【讨论】：
我建议在使用 WGS84 的任何距离上都可以显着提高准确度，因为高度差异无关紧要。在海拔差异很重要的任何距离上，您可能应该只使用直线近似。
【讨论】：
为了做到这一点，您必须解决的第一个问题是如何定义高度变化。正规方程之所以有效，是因为它们位于二维表面上，但是添加第三维意味着最短距离的简单定义不再适用，例如，现在第三维正在“发挥作用”，您的最短距离可以穿过原始椭球体。 这有点快和肮脏，但您最好的解决方案可能是假设沿椭圆体上的原始 2D 路径的高度变化率是恒定的。然后，您可以将 2D 距离计算为长度，计算出高度的变化率，然后简单地使用毕达哥拉斯计算长度的增加，三角形的一侧是 2D 距离，高度是第二个长度。
【讨论】：
首先，您需要一个模型来告诉您两点之间的直线上的高度如何变化。如果没有这样的模型，两点之间的距离就没有任何一致的定义。
如果你有一个线性模型（移动 50% 的点之间的距离也意味着你向上通过了 50% 的高度），那么你可以假设整个物体是一个直角三角形；即，为了确定高度变化如何影响距离，你表现得好像世界是平的。沿地面的距离是底边，高度变化是三角形的高度，斜边是您估计的点到点的真实行驶距离。
如果您想进一步改进，您可以注意到上面的模型非常适合无限小的距离，这意味着您可以迭代距离的各个增量，微积分风格，每次使用当前高度计算地面距离，然后使用相同的三角比来计算高度变化对行进距离的贡献。我可能会在一个包含 10 到 100 个片段的 for() 循环中执行此操作，并且可能通过反复试验找出在真实值的 epsilon 范围内所需的片段数。在这个模型下，也可以计算出线积分来计算两点之间的实际距离。
【讨论】：