计算球面上圆弧交点时纬度值偏差过大的原因排查

核心错误原因 #

1. 笛卡尔坐标转换公式混淆了地理纬度与极角

你使用的坐标转换公式把地理纬度当成了极角（从z轴正方向，即北极点的角度），与标准地理坐标系定义完全相反：

// 你的错误公式
x = sin(toRadians(p.lat)) * cos(toRadians(p.lon));
y = sin(toRadians(p.lat)) * sin(toRadians(p.lon));
z = cos(toRadians(p.lat));

标准地理坐标（lat为北纬/南纬，lon为东经/西经）转单位球面笛卡尔坐标的正确公式应为：

// 正确公式
x = cos(toRadians(p.lat)) * cos(toRadians(p.lon));
y = cos(toRadians(p.lat)) * sin(toRadians(p.lon));
z = sin(toRadians(p.lat));

两者的差异在于：你的公式将输入的lat视为极角（从北极点到点的夹角），而地理纬度是从赤道到点的夹角，两者关系为极角 = 90° - 地理纬度，这直接导致所有输入点的笛卡尔坐标完全错位。

2. 球面坐标转回地理纬度的公式错误

由于坐标转换的前置错误，你后续转回球面坐标的公式也不符合地理纬度定义：

// 你的错误公式
lat = toDegrees(atan2(sqrt(x * x + y * y), z));

该公式计算的是极角，直接赋值给纬度会出现超过±90°的无效值（比如你得到的120.77°纬度）。正确的转换公式（基于标准笛卡尔坐标）应为：

// 正确公式：地理纬度范围-90°~90°
lat = toDegrees(atan2(z, sqrt(x * x + y * y)));
// 经度范围-180°~180°
lon = toDegrees(atan2(y, x));

3. 算法适用场景错误（小圆误用大圆逻辑）

你当前的算法仅能求解两个大圆的交点（大圆是经过球心的平面与球面的交线），但测试用例中的第一条圆弧是小圆（北纬60°的纬线，非赤道的纬线所在平面不经过球心）。用大圆的平面法向量计算逻辑处理小圆，会得到完全错误的平面参数，自然无法得到正确交点。

测试用例修正验证 #

你的测试用例：

• 圆弧1：(60°N,30°E) → (60°N,60°E)（北纬60°纬线，小圆）
• 圆弧2：(40°N,50°E) → (60°N,50°E)（东经50°经线，大圆）

正确交点应为(60°N,50°E)，修正后可通过以下方式计算：

1. 经线（大圆）的平面方程：法向量为(sin(50°), -cos(50°), 0)，平面方程为y*cos(50°) - x*sin(50°) = 0
2. 北纬60°纬线（小圆）的平面方程：所有点的z坐标等于sin(60°)（标准转换公式），即z ≈ 0.8660
3. 联立球面方程x²+y²+z²=1，解出的点对应地理坐标(60°N,50°E)，完全符合预期。

通用解决方案 #

1. 修正坐标转换逻辑：严格遵循标准地理坐标系与笛卡尔坐标的转换公式，明确区分地理纬度和极角的定义。
2. 区分大圆与小圆处理：
   • 大圆：使用两点叉乘获取平面法向量（大圆平面经过球心，两点与球心共面），再通过两个法向量叉乘得到交线方向，单位化后得到交点。
   • 小圆：先确定小圆的平面方程（可通过圆心点、半径或三个共面点计算），再联立球面方程和两个平面方程求解交点。
3. 修正球面坐标转换公式：确保转回的纬度、经度符合地理坐标范围（纬度-90°_{90°，经度-180°}180°）。

内容的提问来源于stack exchange，提问作者Remerd