我正在尝试将状态向量(位置和速度)转换为开普勒元素,但是在尝试计算真实异常时,我遇到了负速度或位置会给我错误结果的问题。
以下是我尝试计算真实异常的不同方法:
/// <summary>
/// https://en.wikipedia.org/wiki/True_anomaly#From_state_vectors
/// </summary>
public static double TrueAnomaly(Vector4 eccentVector, Vector4 position, Vector4 velocity)
{
var dotEccPos = Vector4.Dot(eccentVector, position);
var talen = eccentVector.Length() * position.Length();
talen = dotEccPos / talen;
talen = GMath.Clamp(talen, -1, 1);
var trueAnomoly = Math.Acos(talen);
if (Vector4.Dot(position, velocity) < 0)
trueAnomoly = Math.PI * 2 - trueAnomoly;
return trueAnomoly;
}
//sgp = standard gravitational parameter
public static double TrueAnomaly(double sgp, Vector4 position, Vector4 velocity)
{
var H = Vector4.Cross(position, velocity).Length();
var R = position.Length();
var q = Vector4.Dot(position, velocity); // dot product of r*v
var TAx = H * H / (R * sgp) - 1;
var TAy = H * q / (R * sgp);
var TA = Math.Atan2(TAy, TAx);
return TA;
}
public static double TrueAnomalyFromEccentricAnomaly(double eccentricity, double eccentricAnomaly)
{
var x = Math.Sqrt(1 - Math.Pow(eccentricity, 2)) * Math.Sin(eccentricAnomaly);
var y = Math.Cos(eccentricAnomaly) - eccentricity;
return Math.Atan2(x, y);
}
public static double TrueAnomalyFromEccentricAnomaly2(double eccentricity, double eccentricAnomaly)
{
var x = Math.Cos(eccentricAnomaly) - eccentricity;
var y = 1 - eccentricity * Math.Cos(eccentricAnomaly);
return Math.Acos(x / y);
}
编辑:Spectre指出的另一种方法:
public static double TrueAnomaly(Vector4 position, double loP)
{
return Math.Atan2(position.Y, position.X) - loP;
}
位置都是相对于父体的。
如果 position.x、position.y 和 velocity.y 都是正数,这些函数都一致。如何解决这些问题,以便在位置和速度为负数时获得一致的结果?
澄清一下:我的角度似乎是正确的,只是根据位置和/或速度矢量指向错误的象限。
所以我发现了一个边缘情况,上面的大多数计算都失败了。给定位置和速度:
pos = new Vector4() { X = -0.208994076275941, Y = 0.955838328099748 };
vel = new Vector4() { X = -2.1678187689294E-07, Y = -7.93096769486992E-08 };
当我认为它应该返回 `31.1(非负数)时,我得到了一些奇怪的结果,即 ~ -31.1 度。其中一个返回〜328.8。
然而,用这个位置和速度进行测试,结果似乎还可以:
pos = new Vector4() { X = -0.25, Y = 0.25 };
vel = new Vector4() { X = Distance.KmToAU(-25), Y = Distance.KmToAU(-25) };
有关我如何测试的额外代码以及我对其他一些变量使用的数学,请参阅我的答案。
我在这个上绕圈子。这是我现有代码中的错误的结果,该错误在某些情况下出现,但在其他情况下不出现。我想现在真正的问题是为什么我会得到与我的期望或彼此不匹配的位置/速度以上的不同结果?