diff --git a/src/orbital/utilities.py b/src/orbital/utilities.py index d6a9469..2f24710 100644 --- a/src/orbital/utilities.py +++ b/src/orbital/utilities.py @@ -284,12 +284,14 @@ def elements_from_state_vector(r, v, mu): # Argument of periapsis is the angle between # eccentricity vector and its x component. arg_pe = acos(ev.x / norm(ev)) + if ev.y < 0: + arg_pe = mod(-arg_pe, 2 * pi) else: # Right ascension of ascending node is the angle # between the node vector and its x component. raan = acos(n.x / norm(n)) if n.y < 0: - raan = 2 * pi - raan + raan = mod(-raan, 2 * pi) # Argument of periapsis is angle between # node and eccentricity vectors. @@ -301,23 +303,23 @@ def elements_from_state_vector(r, v, mu): # vector and its x component. f = acos(r.x / norm(r)) if v.x > 0: - f = 2 * pi - f + f = mod(-f, 2 * pi) else: # True anomaly is angle between node # vector and position vector. f = acos(dot(n, r) / (norm(n) * norm(r))) if dot(n, v) > 0: - f = 2 * pi - f + f = mod(-f, 2 * pi) else: if ev.z < 0: - arg_pe = 2 * pi - arg_pe + arg_pe = mod(-arg_pe, 2 * pi) # True anomaly is angle between eccentricity # vector and position vector. f = acos(dot(ev, r) / (norm(ev) * norm(r))) if dot(r, v) < 0: - f = 2 * pi - f + f = mod(-f, 2 * pi) return OrbitalElements(a=a, e=e, i=i, raan=raan, arg_pe=arg_pe, f=f) diff --git a/tests/test_orbital.py b/tests/test_orbital.py index 083d6ab..659d44a 100644 --- a/tests/test_orbital.py +++ b/tests/test_orbital.py @@ -761,15 +761,13 @@ def test_from_state_vector_elliptical_arg_pe_gt_180(self): V = Velocity(16703.9010129, 0, 0) orbit = KeplerianElements.from_state_vector(R, V, body=earth) - # XXX These do not match (they are 180° out, due to arg_pe). - # numpy.testing.assert_almost_equal(orbit.r, R) - # numpy.testing.assert_almost_equal(orbit.v, V) + numpy.testing.assert_almost_equal(orbit.r, R) + numpy.testing.assert_almost_equal(orbit.v, V) self.assertAlmostEqual(orbit.a, 10000000.0, places=3) self.assertAlmostEqual(orbit.e, 0.75) self.assertAlmostEqual(orbit.i, 0.0) self.assertAlmostEqual(orbit.raan, 0.0) - # XXX This is incorrectly calculated as 90°. - # self.assertAlmostEqual(orbit.arg_pe, radians(270.0)) + self.assertAlmostEqual(orbit.arg_pe, radians(270.0)) self.assertAlmostEqual(orbit.M0, 0.0) self.assertAlmostEqual(orbit.ref_epoch, J2000) @@ -846,7 +844,7 @@ def test_from_state_vector_roundtrips(self): # Elliptical flat (e=0.75), prograde. (Position(2500000, 0, 0), Velocity(0, 16703.901013, 0)), # Elliptical flat (e=0.75), arg_pe > 180°. - # (Position(0, -2500000, 0), Velocity(16703.901013, 0, 0)), + (Position(0, -2500000, 0), Velocity(16703.901013, 0, 0)), # Elliptical flat (e=0.75), retrograde. # (Position(2500000, 0, 0), Velocity(0, -16703.901013, 0)), # Elliptical flat (e=0.75), f > 180°. @@ -1077,19 +1075,12 @@ def test_set_v_arg_pe_gt_180(self): # arg_pe = 270°. V = Velocity(-10000, 0, 0) - def set_v(value): - orbit.v = value - - # XXX The 'r and v changed' detection logic is triggered in this case, - # causing a RuntimeError to be raised. If this was not raised, the - # following asserts would be wildly off. - self.assertRaises(RuntimeError, set_v, V) - # numpy.testing.assert_almost_equal(orbit.r, R) - # numpy.testing.assert_almost_equal(orbit.v, V) + orbit.v = V + numpy.testing.assert_almost_equal(orbit.r, R) + numpy.testing.assert_almost_equal(orbit.v, V) # arg_pe should have rotated around 180°, and M0 to match (so r is in # the same spot as it was before). - # XXX This is incorrectly calculated as 90°. - # self.assertAlmostEqual(orbit.arg_pe, radians(270.0)) + self.assertAlmostEqual(orbit.arg_pe, radians(270.0)) self.assertAlmostEqual(orbit.M0, radians(180.0)) def test_set_n(self):