0
$\begingroup$

I'm trying to calculate future spacecraft's position and velocity at time t. I need to calculate the true anomaly at time t for this. I've already made it for the elliptic orbits and it is working good. But I also need it for the hyperbolic orbits. I write some code to do it, but looks like it is working wrong (especially if future time t is far that 1s).

How to calculate the future true anomaly in case of hyperbolic orbit? Did I make errors in my hyperbolic calculations?

There is my code of future true anomaly calculation:

void calculate_keplerian_params(float mu, godot::Vector2 r_vec, godot::Vector2 v_vec, float t)
{
    float h = r_vec.cross(v_vec);

    float r = r_vec.length();
    float v = v_vec.length();

    float E = 0.5f * (v * v) - mu / r;
    float a = -mu / (2.f * E);
    float e = std::sqrtf(1.f - (h * h) / (a * mu));

    godot::Vector2 e_vec = (v * v / mu - 1.f / r) * r_vec - (r_vec.dot(v_vec) / mu) * v_vec;

    float i = 0.f;
    float omega_LAN = e_vec.angle() - M_PI;
    float omega_RAN = e_vec.angle();

    float p = a * (1.f - e * e);
    float nu = 0.f;
    if (e < 1.f)
    {
        nu = std::atan2f(std::sqrtf(p / mu) * r_vec.dot(v_vec), p - r);
    }
    else
    {
        nu = std::acosf((p - r) / (e * r));
    }

    float EA = 0.f;
    float M0 = 0.f;

    if (e < 1.f)
    {
        EA = 2.f * std::atanf(sqrt((1.f - e) / (1.f + e)) * std::tanf(nu / 2.f));
        M0 = EA - e * sin(EA);
    }
    else
    {
        EA = 2.f * std::atanhf(std::sqrtf((e - 1.f) / (e + 1.f)) * std::tanf(nu / 2.f));
        M0 = e * std::sinhf(EA) - EA;
    }

    float n = std::sqrtf(mu / std::fabs(a * a * a));
    float Mt = M0 + n * t;
    float EAt = newtons_method(Mt, e, Mt);

    float nut = 0.f;
    if (e < 1.f)
    {
        nut = 2.f * std::atanf(sqrt((1.f + e) / (1.f - e)) * std::tanf(EAt / 2.f));
    }
    else
    {
        nut = 2.f * std::atanhf(std::sqrtf((e + 1.f) / (e - 1.f)) * std::tanf(EAt / 2.f));
    }

    // write params to store struct

}

Where:

float godot::Orbit::fn(float E, float e, float M) const
{
    if (e < 1.f) {
        return E - e * std::sinf(E) - M;
    } else {
        //return e * std::sinhf(E) - E - M;
        return M - e * std::sinhf(E) + E;
    }
}

float godot::Orbit::de(float E, float e) const
{
    if (e < 1.f) {
        return 1.f - e * std::cosf(E);
    } else {
        return e * std::coshf(E) - 1.f;
    }
}

float godot::Orbit::newtons_method(float E0, float e, float M, float max_error, uint32_t max_iter)
{
    float value = E0, i = 0, root;

    float sign = e < 1.f ? 1.f : -1.f;

    do {
        i++;
        value = value - sign * (fn(value, e, M) / de(value, e));
    } while (std::fabs(fn(value, e, M)) > max_error && i < max_iter);

    return value;
}

I used this formulas to calculate position and velocity from the true anomaly:

godot::Vector2 position = e_vec.normalized().rotated(predictedTA) * radiust;
godot::Vector2 velocity = godot::Vector2(
            position.x * h * e / (radiust * p) * std::sinf(predictedTA) - h / radiust * (std::cosf(omega_RAN) * std::sinf(predictedTA) + std::sinf(omega_RAN) * std::cosf(predictedTA)),
            position.y * h * e / (radiust * p) * std::sinf(predictedTA) - h / radiust * (std::sinf(omega_RAN) * std::sinf(predictedTA) - std::cosf(omega_RAN) * std::cosf(predictedTA))
        );

Can I use the same formulas to calculate velocity and position in case of hyperbolic movement? Should I change it to something different?

$\endgroup$
5
  • 1
    $\begingroup$ In float godot::Orbit::fn(float E, float e, float M) const , the line return e * std::sinhf(E) - E - M should not be commented out; it is the correct function relating Hyperbolic anomaly to Mean anomaly for Newton's Method. $\endgroup$
    – notovny
    Nov 16 at 17:27
  • 1
    $\begingroup$ "Can I use the same formulas...?" You haven't shown any formula at all, just big long chunks of code. They are not the same thing and they don't always behave the same way because code can have all kinds of bugs or unexpected behavior. It is also hard to read other people's code when looking for mathematics rather than programming. If you want to ask a question about formulas, remove all the code and just post the actual formula in question using MathJax. You can learn how to use it with this tutorial or just reviewing other posts. $\endgroup$
    – uhoh
    Nov 16 at 20:51
  • $\begingroup$ I think you may turn out to be quite an active user here in Space SE. While the first time we use MathJax it takes a little time, it's actually quite easy to get used to and you may find it quite worthwhile. $\endgroup$
    – uhoh
    Nov 16 at 21:01
  • 1
    $\begingroup$ @uhoh Sorry, I will try to use MathJax $\endgroup$
    – Robotex
    Nov 16 at 21:20
  • 1
    $\begingroup$ As I've posted in another thread, there are validated algorithms in C++ and Rust to compute all of the orbital elements might want to compute. You may also check out this reference which details, in maths, how to compute these: nyxspace.com/MathSpec/celestial/orbital_elements . $\endgroup$
    – ChrisR
    Nov 17 at 3:47

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.