Python-Vincenty的逆公式不收敛(查找地球上各点之间的距离) [英] Python - Vincenty's inverse formula not converging (Finding distance between points on Earth)

查看:338
本文介绍了Python-Vincenty的逆公式不收敛(查找地球上各点之间的距离)的处理方法,对大家解决问题具有一定的参考价值,需要的朋友们下面随着小编来一起学习吧!

问题描述

我正在尝试实现维基此处中所述的Vincenty逆问题.

I'm attempting to implement Vincenty's inverse problem as described on wiki HERE

问题在于lambda根本没有收敛.如果我尝试遍历公式序列,则该值保持不变,但我真的不确定为什么.也许我只是盯着自己看一个明显的问题.

The problem is that lambda is simply not converging. The value stays the same if I try to iterate over the sequence of formulas, and I'm really not sure why. Perhaps I've just stared myself blind on an obvious problem.

应该注意的是,我是Python的新手,仍然在学习该语言,所以我不确定是否滥用该语言可能会导致问题,或者我在某些计算中是否确实有错误我表演的我似乎在公式中找不到任何错误.

It should be noted that I'm new to Python and still learning the language, so I'm not sure if it's misuse of the language that might cause the problem, or if I do have some mistakes in some of the calculations that I perform. I just can't seem to find any mistakes in the formulas.

基本上,我以与Wiki文章尽可能接近的格式编写了代码,结果是这样的:

Basically, I've written in the code in as close of a format as I could to the wiki article, and the result is this:

import math

# Length of radius at equator of the ellipsoid
a = 6378137.0

# Flattening of the ellipsoid
f = 1/298.257223563

# Length of radius at the poles of the ellipsoid
b = (1 - f) * a

# Latitude points
la1, la2 = 10, 60

# Longitude points
lo1, lo2 = 5, 150

# For the inverse problem, we calculate U1, U2 and L.
# We set the initial value of lamb = L
u1 = math.atan( (1 - f) * math.tan(la1) )
u2 = math.atan( (1 - f) * math.tan(la2) )
L = (lo2 - lo1) * 0.0174532925

lamb = L

while True:
    sinArc = math.sqrt( math.pow(math.cos(u2) * math.sin(lamb),2) + math.pow(math.cos(u1) * math.sin(u2) - math.sin(u1) * math.cos(u2) * math.cos(lamb),2) )
    cosArc = math.sin(u1) * math.sin(u2) + math.cos(u1) * math.cos(u2) * math.cos(lamb)
    arc = math.atan2(sinArc, cosArc)
    sinAzimuth = ( math.cos(u1) * math.cos(u2) * math.sin(lamb) ) // ( sinArc )
    cosAzimuthSqr = 1 - math.pow(sinAzimuth, 2)
    cosProduct = cosArc - ((2 * math.sin(u1) * math.sin(u2) ) // (cosAzimuthSqr))
    C = (f//16) * cosAzimuthSqr  * (4 + f * (4 - 3 * cosAzimuthSqr))
    lamb = L + (1 - C) * f * sinAzimuth * ( arc + C * sinArc * ( cosProduct + C * cosArc * (-1 + 2 * math.pow(cosProduct, 2))))
    print(lamb)

如上所述,问题是值"lamb"(lambda)不会变小.我什至尝试将我的代码与其他实现进行比较,但是它们看起来几乎相同.

As mentioned the problem is that the value "lamb" (lambda) will not become smaller. I've even tried to compare my code to other implementations, but they looked just about the same.

我在这里做错了什么? :-)

What am I doing wrong here? :-)

谢谢大家!

推荐答案

首先,您也应该将纬度转换为弧度(已经对经度进行了转换):

First, you should convert you latitudes in radians too (you already do this for your longitudes):

u1 = math.atan( (1 - f) * math.tan(math.radians(la1)) )
u2 = math.atan( (1 - f) * math.tan(math.radians(la2)) )
L = math.radians((lo2 - lo1)) # better than * 0.0174532925

执行此操作并摆脱//(int划分)并将其替换为/(float划分),lambda停止在迭代过程中重复相同的值,并开始遵循此路径(根据您的示例坐标):

Once you do this and get rid of // (int divisions) and replace them by / (float divisions), lambda stops repeating the same value through your iterations and starts following this path (based on your example coordinates):

2.5325205864224847
2.5325167509030906
2.532516759118641
2.532516759101044
2.5325167591010813
2.5325167591010813
2.5325167591010813

您似乎希望收敛精度为10^(−12),这似乎很重要.

As you seem to expect a convergence precision of 10^(−12), it seems to make the point.

现在您可以退出循环(lambda收敛)并继续操作,直到计算出所需的测地距离s.

You can now exit the loop (lambda having converged) and keep going until you compute the desired geodesic distance s.

注意:您可以在s此处测试您的最终值.

Note: you can test your final value s here.

这篇关于Python-Vincenty的逆公式不收敛(查找地球上各点之间的距离)的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持IT屋!

查看全文
登录 关闭
扫码关注1秒登录
发送“验证码”获取 | 15天全站免登陆