计算测地距离时考虑高度 [英] Taking altitude into account when calculating geodesic distance
问题描述
我目前正在处理GPS数据和精确的高度测量。
我想计算两个连续点之间的距离。关于使用WGS84椭球等计算两点之间的距离，有很多
的信息。然而，我没有找到任何信息
没有人知道一些网站，论文，书籍等，它们描述了这个
距离计算中的高度更改。这样的方法？
感谢
编辑：计算距离时，Sql Server 2008的地理分机也忽略了海拔高度信息。
我使用开始和结束高度的平均值作为恒定高度，实现了WGS84距离函数。如果你确定沿着你的路线会有相对较小的高度变化，这个效果可以接受的很好（错误是相对于你的两个LLA点的高度差异而言）。
这里是我的代码（C＃）：
///< summary>
///获取当前模式坐标系中两个路径点之间的测地距离
///< / summary>
///< param name =point1>第一点< / param>
///< param name =point2>第二点< / param>
///< param name =mode>两点都在&的坐标模式< / param>
///<返回>当前坐标模式中的两点之间的距离< /返回>
public static double GetGeodesicDistance（PathPoint point1，PathPoint point2，CoordMode mode）{
//计算LLA路径的正确测地线
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;
//不完美，但使用平均高度
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米左右，我们发现平均约5米的变化与未修改的WGS84椭圆体相比。
编辑：
为了补充说明，如果您确实希望大幅度改变高度，您可以将WGS84坐标转换为ECEF（以地球为中心的固定地球）并评估直线路径我的功能的底部。将点转换为ECEF很简单：
///< summary>
///将格式（Lon，Lat，Alt）中的点转换为ECEF
///< / summary>
///< param name =point>指向（Lon，Lat，Alt）< / param>
///<返回> ECEF中的点< /返回>
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;
返回outPoint;
$ / code> 编辑2：
<
// RF是WGS84椭圆体的偏心率
public const double RF = 298.257223563;
// A是以米为单位的地球半径
public const double A = 6378137.0;
LatLonAltTransformer 是我用来转换的类从LatLonAlt坐标到ECEF坐标，并定义上述常量。
i´m currently dealing with gps data combined with precise altitude measurement. I want to calculate the distance between two consecuting points. There is a lot of information out there about calculating distance between two points using the WGS84 ellipsoid and so on.
however, i did not find any information that takes Altitude changes into account for this distance calculation.
does anyone know about some websites, papers, books etc. that describes such a method? thanks
edit: Sql Server 2008 geographic extensions also neglect altitude information when calculating distance.
I implemented a WGS84 distance function using the average of the start and end altitude as the constant altitude. If you are certain that there will be relatively little altitude variation along your path this works acceptably well (error is relative to the altitude difference of your two LLA points).
Here's my code (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);
}
In practice we've found that the altitude difference rarely makes a large difference, our paths are typically 1-2km long with altitude varying on the order of 100m and we see about ~5m change on average versus using the WGS84 ellipsoid unmodified.
Edit:
To add to this, if you do expect large altitude changes, you can convert your WGS84 coordinates to ECEF (earth centered earth fixed) and evaluate straight-line paths as shown at the bottom of my function. Converting a point to ECEF is simple to do:
/// <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;
}
Edit 2:
I was asked about some of the other variables in my code:
// 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 is a class I used to convert from LatLonAlt coordinates to ECEF coordinates and defines the constants above.
这篇关于计算测地距离时考虑高度的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！
