用纬度和经度近似距离的三元化

Question

我不太明白这件事。

我试图根据3个固定地点的3段距离来近似一个信标的位置(纬度/经度)。然而，可用的距离读数可能有高达1公里的误差。

这里(有精确的测量)、这里、这里( Java中的距离测量误差，但在lat/lon坐标中没有，也没有答案)以及其他公司也曾问过类似的问题。我还设法挖掘了处理不完美测量数据的本论文，但是它假设一个笛卡尔坐标系，而且它也是数学的，而不是接近一个可用的实现。

因此，上述任何一个链接和答案都不适用于以下问题：

• 所有可用的距离测量都是以公里为单位的近似(在这种情况下，数据最常包含的读数在-1公里到100公里之间)。
• 测量误差可达1公里。
• 3根据3个固定位置(已知纬度/经度)进行距离测量。
• 目标近似值也应该是经纬度组合。

到目前为止，我已经将这个答案调整为C#，但是我注意到由于测量的不精确性，该算法不起作用(因为算法假定三个距离圆圈彼此完美相交)：

public static class Trilateration
{
    public static GeoLocation Compute(DistanceReading point1, DistanceReading point2, DistanceReading point3)
    {
        // not my code :P
        // assuming elevation = 0
        const double earthR = 6371d;

        //using authalic sphere
        //if using an ellipsoid this step is slightly different
        //Convert geodetic Lat/Long to ECEF xyz
        //   1. Convert Lat/Long to radians
        //   2d. Convert Lat/Long(radians) to ECEF
        double xA = earthR * (Math.Cos(Radians(point1.GeoLocation.Latitude)) * Math.Cos(Radians(point1.GeoLocation.Longitude)));
        double yA = earthR * (Math.Cos(Radians(point1.GeoLocation.Latitude)) * Math.Sin(Radians(point1.GeoLocation.Longitude)));
        double zA = earthR * Math.Sin(Radians(point1.GeoLocation.Latitude));

        double xB = earthR * (Math.Cos(Radians(point2.GeoLocation.Latitude)) * Math.Cos(Radians(point2.GeoLocation.Longitude)));
        double yB = earthR * (Math.Cos(Radians(point2.GeoLocation.Latitude)) * Math.Sin(Radians(point2.GeoLocation.Longitude)));
        double zB = earthR * (Math.Sin(Radians(point2.GeoLocation.Latitude)));

        double xC = earthR * (Math.Cos(Radians(point3.GeoLocation.Latitude)) * Math.Cos(Radians(point3.GeoLocation.Longitude)));
        double yC = earthR * (Math.Cos(Radians(point3.GeoLocation.Latitude)) * Math.Sin(Radians(point3.GeoLocation.Longitude)));
        double zC = earthR * Math.Sin(Radians(point3.GeoLocation.Latitude));

        // a 64 bit Vector3 implementation :)
        Vector3_64 P1 = new(xA, yA, zA);
        Vector3_64 P2 = new(xB, yB, zB);
        Vector3_64 P3 = new(xC, yC, zC);

        //from wikipedia
        //transform to get circle 1 at origin
        //ransform to get circle 2d on x axis
        Vector3_64 ex = (P2 - P1).Normalize();
        double i = Vector3_64.Dot(ex, P3 - P1);
        Vector3_64 ey = (P3 - P1 - i * ex).Normalize();
        Vector3_64 ez = Vector3_64.Cross(ex, ey);
        double d = (P2 - P1).Length;
        double j = Vector3_64.Dot(ey, P3 - P1);

        //from wikipedia
        //plug and chug using above values
        double x = (Math.Pow(point1.DistanceKm, 2d) - Math.Pow(point2.DistanceKm, 2d) + Math.Pow(d, 2d)) / (2d * d);
        double y = ((Math.Pow(point1.DistanceKm, 2d) - Math.Pow(point3.DistanceKm, 2d) + Math.Pow(i, 2d) + Math.Pow(j, 2d)) / (2d * j)) - ((i / j) * x);

        // only one case shown here
        double z = Math.Sqrt(Math.Pow(point1.DistanceKm, 2d) - Math.Pow(x, 2d) - Math.Pow(y, 2d));

        //triPt is a vector with ECEF x,y,z of trilateration point
        Vector3_64 triPt = P1 + x * ex + y * ey + z * ez;

        //convert back to lat/long from ECEF
        //convert to degrees
        double lat = Degrees(Math.Asin(triPt.Z / earthR));
        double lon = Degrees(Math.Atan2(triPt.Y, triPt.X));

        return new GeoLocation(lat, lon);
    }

    [MethodImpl(MethodImplOptions.AggressiveInlining)]
    private static double Radians(double degrees) =>
        degrees * Math.Tau / 360d;

    [MethodImpl(MethodImplOptions.AggressiveInlining)]
    private static double Degrees(double radians) =>
        radians * 360d / Math.Tau;
}

在我的情况下，上面的代码通常不起作用，而是只返回“number”，因为它在计算最终的z值时试图取负数的平方根(由于测量不准确)。

在我的例子中，度量可能返回这样的数据(用一些随机的在线工具可视化)：

​

​

其中只有两个甚至没有一个距离圆相交：

​

​

我正在寻找的是一个算法，返回目标位置的最佳近似值，基于已知的最大误差为1公里的三种距离测量，或者我可以采取的进一步方法。

我也曾想过迭代圆上的点，然后确定到其他圆上所有点的最小平均距离，但地球的三维球面几何让我头疼。而且，也许有一种更好更简单的方法来解决这个问题，我现在还不知道。

由于这更多是一个算法问题，而不是任何特定于语言的问题，所以我很感激在任何编程语言、伪代码或自然语言方面的任何帮助。

Answer 1

如果您可以访问一个提供非线性优化实用程序的科学计算库，那么您可以尝试找到最小化以下内容的点：

(||x - p_1|| - r_1)^2 + (||x - p_2|| - r_2)^2 + (||x - p_3|| - r_3)^2 + (||x - p_earth|| - r_earth)^2

其中p_i是您测量到的i第四个位置的位置(在笛卡尔坐标中)，r_i是对应的距离读数，p_earth是地球的位置，r_earth是地球的半径，||a||是向量a的范数/长度。

表达式中的每个项都试图最小化残差半径误差。

当然，可以根据您的需要对其进行修改--例如，如果有约束优化，则可以将点位于地球表面的要求编码为约束，而不是用于优化的术语。如果球形地球模型不够精确，你可以定义一个来自地球表面的误差，或者只是将你的结果投射到地球上，如果足够准确的话。