Geopy：计算GPS方向/轴承

Question

这是我第一次发布帖子。

我正在为桥梁检查ROV八旋翼收集的GPS数据进行一些数据分析。我们使用3D扫描LIDAR、立体视觉、INS和其他一些很棒的技术，通过ROS来运行八旋翼。我目前正在使用ublox LEA-6T，类似于Doug Weibel的设置，收集原始GPS数据，如载波相位、多普勒频移和卫星星历。然后，我使用开源项目RTKLIB与本地NOAA CORS站点进行差分GPS后处理，以获得更好的姿态估计精度（用于重建桥梁的3D点云）。

无论如何，我正在使用scipy的大部分功能来对我的测试结果进行统计验证。
具体来说，对于这一部分，我只是在使用：

• python-3.3
• numpy
• geopy

我一直在使用geopy方便的距离函数来研究与测量地面真实值偏移量相关的位置协方差。稍加调整参数，可以找到与矩阵中每个标准差元素所表示的每个方向相对应的距离； 包括北、东、上和它们之间的三个方向。

然而，这些距离是绝对的，不能描述方向。
例如：正数或负数可能分别与北向或南向相关联。

我可以简单地使用纬度和经度来检测方向的极性，
但我想能够找到描述距离的点对点方位角，
因为我认为全局导航值除了当前应用外还可能有其他用途。

我发现有人提出了类似的问题，
但似乎假定了大圆近似值
而我希望至少使用WGS-84椭球模型或任何可在geopy中使用的相同模型：
跳转到计算距离

如有建议，请指教。
-ruffsl

感兴趣的来源：

• python-3.3: http://www.python.org/download/releases/3.3.0/
• numpy: http://www.numpy.org/
• geopy: https://code.google.com/p/geopy/
• scipy: http://www.scipy.org/
• ublox LEA-6T: http://www.u-blox.com/en/gps-modules/u-blox-6-timing-module/lea-6t.html
• Doug Weibel's: http://diydrones.com/profiles/blogs/proof-of-concept-test-extremely-accurate-3d-velocity-measurement
• RTKLIB: http://www.rtklib.com/
• NOAA CORS: http://geodesy.noaa.gov/CORS/
• ROS:http://www.ros.org/wiki/

- ruffsl

---

你说了很多话，但问题不太明确。你需要什么？是一个能够计算两个WG84经纬度坐标之间方位角的函数吗？ - AlexWien

@AlexWien，是的，简而言之，这正是我所追求的。我发现当人们解释他们自己的问题时，更有帮助，这为其他人提供了更相关和可搜索的关键字。 - ruffsl

但对我来说很有启发，尤其是RTKlib。 - AlexWien

Answer 1

使用Python中的geographiclib软件包。它可以在椭球体上计算距离和方位角等内容（您还可以插值路径、测量面积等）。例如，在

pip install geographiclib

你可以做

>>> from geographiclib.geodesic import Geodesic
>>> Geodesic.WGS84.Inverse(-41.32, 174.81, 40.96, -5.50)
{'lat1': -41.32, 'a12': 179.6197069334283, 's12': 19959679.26735382, 'lat2': 40.96, 'azi2': 18.825195123248392, 'azi1': 161.06766998615882, 'lon1': 174.81, 'lon2': -5.5}

这计算了从新西兰惠灵顿（南纬41.32度，东经174.81度）到西班牙萨拉曼卡（北纬40.96度，西经5.50度）的大圆航线。距离由s12给出（19959679米），初始方位角（方向）由azi1给出（相对于正北顺时针旋转161.067...度）。

- cffk

---

如果您在Geodesic.WGS84.Inverse的结果集上提供更多细节，并实际计算问题所要求的GPS航向/方位角，那将会很有帮助。 - kovac

这个链接指向文档，解释了返回的字典。https://geographiclib.sourceforge.io/1.52/python/interface.html - Levi Baguley

Answer 2

@AlexWien在Python中的回答

import math, numpy as np

def get_bearing(lat1,lon1,lat2,lon2):
    dLon = lon2 - lon1;
    y = math.sin(dLon) * math.cos(lat2);
    x = math.cos(lat1)*math.sin(lat2) - math.sin(lat1)*math.cos(lat2)*math.cos(dLon);
    brng = np.rad2deg(math.atan2(y, x));
    if brng < 0: brng+= 360
    return brng

- BBSysDyn

Answer 3

两个纬度/经度坐标之间的方位角：（lat1，lon1），（lat2，lon2）

在下面的代码中，假定lat1、lon1、lat2、lon2为弧度。
在转换之前，从度数转换为弧度。

dLon = lon2 - lon1;
y = Math.sin(dLon) * Math.cos(lat2);
x = Math.cos(lat1)*Math.sin(lat2) -
        Math.sin(lat1)*Math.cos(lat2)*Math.cos(dLon);
brng = Math.atan2(y, x).toDeg();

轴承现在在范围-180/180内。

将其标准化为指南针轴承（0-360）。

if brng < 0: brng+= 360

- AlexWien

---

我喜欢优美的数学，但我认为这是使用球坐标系，因此仍然是大圆近似。对于导航和一般方向应用（包括我自己检查极性的目的），与WG84相比，我想这已经足够了。但是对于测量和现场检查来说，这可能太粗糙了，或者不适用于本地事件估计。很棒的技术。 - ruffsl

是的，这些公式假设地球是一个椭球体。如果想要更精确，你需要使用基于椭球体的公式/实现。 - AlexWien

我感到困惑...为什么 Geopy 没有一个与计算距离的 Vincenty 函数相匹配的函数？ - Michael