我有一个航班数据集。飞机与地面站通信，通过 GPS，我可以获得当前飞机位置和地面站位置之间的距离，以及各自的经度、纬度和高度。 我现在需要做的是将 long、lat、alt 位置数据转换为 x、y、z 坐标系。
由于地面站处于静态位置，我希望它位于坐标系的原点（0,0,0）。 我正在寻找一个可以处理这种转换的 python 模块。
作为示例，我提供了来自地面站的实际数据以及需要转换的飞机的 1 个位置：
地面站 [应该是 (0, 0, 0)]
纬度：41.947.694
经度：3.209.083
海拔高度（米）：379.41
飞机位置： 纬度：419.763.015.750 长：34.890.851.772 海拔[米]: 971.32
将地面站移动到原点，如果我理解正确的话，我必须从飞机位置中减去它的位置，这意味着：
地面站： 纬度: 0 经度: 0 海拔高度（米）：0
飞机位置：
纬度：419.721.068.056
长：34.887.642.689
海拔[米]: 591.91
由于我是处理此类数据的新手，并尝试了一些简单的解决方案，但没有找到任何转换此数据的内容，我希望您能帮助我找到解决方案。
我需要转换后的数据来模拟 3D 椭球体以找到可能的交点。 这也可以使用与 python 不同的方法来完成。所以如果我从来没有使用过的话，matlab 也很好。
提前致谢！
转换纬度、经度和高度 GPS 坐标， 这是大地坐标（垂直方向是 垂直于椭球体而不是地球中心的方向） 基于 WGS84 参考椭球体，以“xyz”坐标，地心， 在所谓的地心地球固定框架中是一个函数 在多个库中以多种语言实现。自从Matlab 提到过，它的实现为
geodetic2ecef一些文本和网站也对此进行了描述，甚至在 维基百科， 这样我们就可以安全地从头开始实施它。
如果将站点坐标和平面坐标都变换为这个 ECEF x,y,z 坐标，你可以计算它们的差值，然后你 已经有一个 XYZ 坐标系统，即 (0, 0, 0) 车站。
请注意，确实，正如评论中提到的，您不能执行 纬度或经度本身与任何位置的差异 有意义的解释。找出两个之间的距离 点我们只能减去大小，如向量或笛卡尔坐标； 纬度和经度与这些不线性相关。
现在，按照前面描述计算的 xyz 差异并不是很大 本身很重要，因为它们的参考方向 是全局的 - z 轴平行于旋转轴 地球和 xz 平面平行于本初子午线。
为此，我建议通过添加（双）来“本地化”这些差异 旋转这些全局方向以使它们指向 与参考站相关的方向：新的 z 轴 为大地垂直线（朝向天顶）， 新的 x 轴从车站向北，y 轴向东。
同样，这些计算在文献中有详细记录，并且 在线，在这里从头开始实现它们是安全的。
这是一个例子，可能是用python实现的，效率不是很高 因为它重复计算三角函数 基础，但我的目的是尽可能清楚地表达，而不是 高效。
import math
a = 6378137 # semi major axis, WGS84
# f = 1 / 298.257223563 # WGS84
e2 = 6.69437999014e-3 # eccentricity, WGS84 e2 = 1 - (1-f)*(1-f)
def gps_to_geocentric_ecef(lat, lon, alt):
# converts geodetic latitude, longitude, altitude based on WGS84 (GPS coordinates)
# to geocentric x, y, z in meters (ECEF), sqrt(x^2+y^2+z^2) = geocentric radius
# |z| distance to equatorial plane, |x| distance to the plane of the prime meridian
# https://en.wikipedia.org/wiki/Geodetic_coordinates#Conversion
# https://github.com/proj4js/proj4js/blob/master/lib/datumUtils.js#L68C5-L70C44
# verify: https://tool-online.com/en/coordinate-converter.php (right select WORLD and XYZ (geocentric))
# or: https://www.oc.nps.edu/oc2902w/coord/llhxyz.htm - paste lat, lon, height, then press "LLH to ECEF"
lat = math.radians(lat)
lon = math.radians(lon)
sin_lat = math.sin(lat)
cos_lat = math.cos(lat)
cos_lon = math.cos(lon)
sin_lon = math.sin(lon)
rn = a / (math.sqrt(1 - e2 * sin_lat * sin_lat))
x = (rn + alt) * cos_lat * cos_lon
y = (rn + alt) * cos_lat * sin_lon
z = ((rn * (1 - e2)) + alt) * sin_lat
return [x, y, z]
def rotate_ecef_to_local_xyz(dx, dy, dz, lat, lon):
# rotates the coordinate differences dx, dy, dz fromm the global frame (ECEF) to
# a local frame that has the geodetic normal to (lat, lon) as z-axis, while
# the new x-axis is towards North and the y-axis towards East
lat = math.radians(lat)
lon = math.radians(lon)
cos_lon = math.cos(lon)
sin_lon = math.sin(lon)
cos_lat = math.cos(lat)
sin_lat = math.sin(lat)
dx, dy, dz = (- dx * cos_lon * sin_lat - dy * sin_lon * sin_lat + dz * cos_lat,
- dx * sin_lon + dy * cos_lon,
dx * cos_lon * cos_lat + dy * sin_lon * cos_lat + dz * sin_lat)
return [dx, dy, dz]
def main():
bgr_lat, bgr_lon, bgr_alt = [41.947694, 3.209083, 379.41]
airplane_gps = [
[41.9763015750, 3.4890851772, 971.32]
]
# ecef coordinates for the station:
bgr_xyz = gps_to_geocentric_ecef(bgr_lat, bgr_lon, bgr_alt)
# ecef coordinates for the (first) airplane:
plane_xyz = gps_to_geocentric_ecef(*airplane_gps[0])
# xyz position of the airplane relative to the station, ecef orientation:
diff_xyz = [plane_xyz[i] - bgr_xyz[i] for i in range(3)]
# xyz position of the airplane relative to the station, station-local orientation
local_xyz = rotate_ecef_to_local_xyz(*diff_xyz, bgr_lat, bgr_lon)
print(local_xyz) # final result, in meters (x-towards North, y-towards East, z-towards Zenith)
if __name__ == '__main__':
main()
由于海拔高度以米为单位，因此我使用了 WGS84 半长轴也以米为单位，结果以米为单位： 该平面向北 3215m，向东 23210m， 高 549m（从地球切线平面测量） 经过参考站 BGR）
虽然我通过多种方式彻底验证了这段代码，但我仍然希望 没有人会基于它来驾驶飞机:)