使用3个经纬度点和3个距离的三边测量法

维基百科在这里对代数进行了一次非常深入的讨论: http://en.wikipedia.org/wiki/Trilateration

第一步是将lat/long坐标转换为笛卡尔坐标,这一点在维基百科条目中没有真正涉及:

x0 = cos( lon0 ) * cos( lat0 ) , y0 = sin( lon0 ) * cos( lat0 ) , z0 = sin( lat0 )
x1 = cos( lon1 ) * cos( lat0 ) , y1 = sin( lon1 ) * cos( lat1 ) , z1 = sin( lat1 )
x2 = cos( lon2 ) * cos( lat0 ) , y2 = sin( lon2 ) * cos( lat2 ) , z2 = sin( lat2 )

(为了使计算简单,我做了一些修改,所以我们是以“地球半径”为单位而不是以公里为单位工作的)

为了你的数据,我得到

         p0            p1           p2
X   -0.420442596  -0.420430618  -0.42040255
Y   -0.67380418   -0.673826567  -0.673825967
Z    0.607631426   0.607614889   0.607634975

维基百科文章中介绍的下一步是简化坐标,通过转换点使p0在原点,然后旋转使p1在x轴上,p2在x-y平面上。

对于翻译,只需从p1和p2中减去p0:

    p0a      p1a          p2a
X   0    1.19779E-05   4.00462E-05
Y   0   -2.23864E-05  -2.17865E-05
Z   0   -1.65372E-05   3.5486E-06

旋转不难。p1b得到(x,y)=(d,0),其中d是原点到p1a的距离(勾股定理)

对于p2b,我们需要将p2a分解成两个部分:一个平行于p1a(沿着x轴),另一个垂直于p1a(在“b”坐标系中沿着y轴)。

要做到这一点,我们需要一个朝p1a方向的单位向量,即p1a*(1/d)。以这个单位向量的点积(如果你喜欢叫它p1a_Hat)和p2a,这是p2b的x坐标。

现在Y坐标很容易。坐标变换下原点到p2的长度不能改变。所以用勾股定理计算p2a的长度,然后用勾股定理“向后”得到p2b的y坐标,以保持长度不变。这就是维基百科称之为“j”的变量。(注意,有一个含糊不清的问题,我将留给你去弄清楚j是正的还是负的)。

现在你已经得到了三个变量d,i和j,这是维基百科文章用来计算的。现在你可以将它们转换回公里,用地球半径乘以。你应该可以从这里开始计算剩下的部分

(顺便说一句,维基百科给出了不同的坐标转换计算。我喜欢尽可能避免三角函数)。

我在新形成的GIS堆栈交换上问了这个问题,也得到了一些很好的答案。

https://gis.stackexchange.com/questions/66/trilateration-using-3-latitude-and-longitude-points-and-3-distances

这里接受的答案在python中有一个(大概)工作解决方案:

https://gis.stackexchange.com/questions/66/trilateration-using-3-latitude-and-longitude-points-and-3-distances/415#415

上 Paul Bourke Geometry 页

intersection of two circles

考虑以下9个圆 A、B、C点和距离d1、d2、d3

• A的中心,半径d1
• A的中心,半径d2
• A的中心,半径d3
• B中心,半径d1
• B的中心,半径d2
• B的中心,半径d3
• C中心,半径d1
• C中心,半径d2
• C中心,半径d3

这些是你可能的圈子。现在我们可以剔除这些,因为我们知道如果d1用于a,它将不会用于b。

这使您可以输入,其中a1表示圆心为a,半径为d1的圆:

• {a1,b2,c3}
• {a1,b3,c2}
• {a2,b1,c3}
• {a2,b3,c1}
• {a3,b1,c2}
• {a3,b2,c1}

你应该能够将lat/long转换为x,y,z,知道地球的半径,以及沿着地壳弯曲距离到直线距离的距离,从那里你可以看到它们中的哪一个相交于一个公共点。记住要考虑到浮动缺陷造成的微小误差。