from_state_vector inverts arg_pe when Z coordinate is extremely small and negative

[mgiuca]

Reproduction case:

from orbital import Position, Velocity, KeplerianElements, earth
R = Position(0, 2500000, -0.01)
V = Velocity(-16703.901013, 0, 0)
print(KeplerianElements.from_state_vector(R, V, body=earth))

Expected output:

KeplerianElements:
    Semimajor axis (a)                           =  10000.000 km
    Eccentricity (e)                             =      0.750000
    Inclination (i)                              =      0.0 deg
    Right ascension of the ascending node (raan) =      0.0 deg
    Argument of perigee (arg_pe)                 =     90.0 deg
    Mean anomaly at reference epoch (M0)         =      0.0 deg
    Period (T)                                   = 2:45:52.014055
    Reference epoch (ref_epoch)                  = J2000.000
        Mean anomaly (M)                         =      0.0 deg
        Time (t)                                 = 0:00:00
        Epoch (epoch)                            = J2000.000

Actual output:

KeplerianElements:
    Semimajor axis (a)                           =  10000.000 km
    Eccentricity (e)                             =      0.750000
    Inclination (i)                              =      0.0 deg
    Right ascension of the ascending node (raan) =      0.0 deg
    Argument of perigee (arg_pe)                 =    270.0 deg
    Mean anomaly at reference epoch (M0)         =      0.0 deg
    Period (T)                                   = 2:45:52.014055
    Reference epoch (ref_epoch)                  = J2000.000
        Mean anomaly (M)                         =      0.0 deg
        Time (t)                                 = 0:00:00
        Epoch (epoch)                            = J2000.000

(Note: arg_pe is 270°, it should be 90°)

Analysis:

The incorrectness of the output can be seen by changing the R.z of -0.01 to 0.0, which then changes arg_pe (correctly) to 90°. Such a small change in the position should not rotate the body 180° around its orbit. Alternatively, changing R.z to 1 (correctly) causes raan to become 180° and arg_pe 270°. Either of these would be fine - it doesn't matter what raan and arg_pe are when the inclination is very small, but they have to be 90° apart.

The problem is this code in utilities.elements_from_state_vector:

        if ev.z < 0:
            arg_pe = 2 * pi - arg_pe

This is executed whenever e is not close to 0, inverting arg_pe to compensate for the fact that raan also flips around when the body goes below the horizontal. That is the correct behaviour when i is also not close to 0, but when i is close to 0, raan is hard-coded to 0, so it is not correct to flip arg_pe. So the above code should only run in the "e is not small and i is also not small" case.

Note that if i is exactly 0, this isn't an issue, since ev.z will also be 0. The bug is when the orbit is very very slightly inclined, but i counts as close enough to 0.

I am working on a fix for this alongside several other bugs in from_state_vector and will send a PR soon.