Add Geosynchronous Orbit

src/tatc/analysis/dop.py

@@ -5,6 +5,7 @@
 @author: Michael P. Jones <mpj@mit.edu>
 @author: Paul T. Grogan <paul.grogan@asu.edu>
 """
+
 import warnings
 from datetime import datetime
 from enum import Enum

src/tatc/constants.py

@@ -10,7 +10,6 @@
 import numpy as np
 from skyfield.api import load, Loader
 
-
 # load ephemeris file
 resources_dir = os.path.join(os.path.dirname(__file__), "resources")
 de421_loader = Loader(resources_dir)

src/tatc/schemas/__init__.py

@@ -11,6 +11,7 @@
     KeplerianOrbit,
     MolniyaOrbit,
     TundraOrbit,
+    GeosynchronousOrbit,
 )
 from .point import Point, GroundStation
 from .satellite import (

src/tatc/schemas/orbit.py

@@ -1181,3 +1181,83 @@ def get_derived_orbit(
             eccentricity=self.eccentricity,
             perigee_argument=self.perigee_argument,
         )
+
+
+class GeosynchronousOrbit(CircularOrbit):
+    """
+    Geosynchronous orbit defined by longitude to be observed.
+    Geostationary orbit can be specified with inclination equal to 0.
+    """
+
+    type: Literal["geosynchronous"] = Field(
+        "geosynchronous", description="Orbit type discriminator."
+    )
+    altitude: float = Field(35786000, description="Mean altitude (meters).")
+    longitude: float = Field(0, description="Longitude (degrees).", ge=0, lt=360)
+
+    def get_true_anomaly(self) -> float:
+        """
+        Gets the true anomaly for the satellite over a given longitude at a given epoch.
+
+        Returns:
+            float: True anomaly in degrees.
+        """
+        t = constants.timescale.from_datetime(self.epoch)
+        gst = t.gmst * 15
+        return (self.longitude + gst) % 360
+
+    def get_derived_orbit(
+        self, delta_mean_anomaly: float, delta_raan: float
+    ) -> GeosynchronousOrbit:
+        """
+        Gets a derived geosynchronous orbit with pertubations to the
+        mean anomaly (interpreted as a shift in observed longitude) and
+        right ascension of ascending node.
+
+        Args:
+            delta_mean_anomaly (float): Delta mean anomaly (degrees).
+            delta_raan (float): Delta right ascension of ascending node (degrees).
+
+        Returns:
+            GeosynchronousOrbit: the derived orbit
+        """
+        true_anomaly = utils.mean_anomaly_to_true_anomaly(
+            np.mod(self.get_mean_anomaly() + delta_mean_anomaly, 360)
+        )
+        longitude = np.mod(self.longitude + delta_mean_anomaly, 360)
+        raan = np.mod(self.right_ascension_ascending_node + delta_raan, 360)
+        return GeosynchronousOrbit(
+            altitude=self.altitude,
+            true_anomaly=true_anomaly,
+            epoch=self.epoch,
+            inclination=self.inclination,
+            right_ascension_ascending_node=raan,
+            longitude=longitude,
+        )
+
+    def to_tle(self, lazy_load: bool = None) -> TwoLineElements:
+        """
+        Converts this orbit to a two line elements representation.
+
+        Args:
+            lazy_load (bool): True, if this tle should be lazy-loaded.
+
+        Returns:
+            TwoLineElements: the two line elements orbit
+        """
+        if lazy_load is None:
+            lazy_load = config.rc.orbit_tle_lazy_load
+        if lazy_load:
+            tle = self.__dict__.get("tle")
+        else:
+            tle = None
+        if tle is None:
+            tle = CircularOrbit(
+                altitude=self.altitude,
+                inclination=self.inclination,
+                right_ascension_ascending_node=self.right_ascension_ascending_node,
+                true_anomaly=self.get_true_anomaly(),
+                epoch=self.epoch,
+            ).to_tle()
+            self.__dict__["tle"] = tle
+        return tle

src/tatc/schemas/satellite.py

@@ -4,6 +4,7 @@
 
 @author: Paul T. Grogan <paul.grogan@asu.edu>
 """
+
 from __future__ import annotations
 
 import copy

src/tatc/utils.py

@@ -4,6 +4,7 @@
 
 @author: Paul T. Grogan <paul.grogan@asu.edu>
 """
+
 import re
 from typing import List, Union
 
@@ -40,7 +41,7 @@ def mean_anomaly_to_true_anomaly(mean_anomaly: float, eccentricity: float = 0) -
 
     Args:
         mean_anomaly (float): The mean anomaly (degrees).
-        true_anomaly (float): The orbit eccentricity.
+        eccentricity (float): The orbit eccentricity.
 
     Returns:
         float: The true anomaly (degrees).

tests/schemas/test_geosynchronous_orbit.py

@@ -0,0 +1,122 @@
+import unittest
+
+from datetime import datetime, timezone
+
+from tatc import constants
+from tatc.schemas import GeosynchronousOrbit
+
+
+class TestGeosynchronousOrbit(unittest.TestCase):
+    def setUp(self):
+        self.test_data = {
+            "altitude": 35786000,
+            "epoch": datetime(2022, 1, 1, 12, tzinfo=timezone.utc),
+            "inclination": 0.0,
+            "right_ascension_ascending_node": 0.0,
+            "longitude": 10.0,
+        }
+        self.test_orbit = GeosynchronousOrbit(**self.test_data)
+
+    def test_good_data(self):
+        self.assertEqual(self.test_orbit.altitude, self.test_data.get("altitude"))
+        self.assertEqual(self.test_orbit.epoch, self.test_data.get("epoch"))
+        self.assertEqual(self.test_orbit.inclination, self.test_data.get("inclination"))
+        self.assertEqual(
+            self.test_orbit.right_ascension_ascending_node,
+            self.test_data.get("right_ascension_ascending_node"),
+        )
+        self.assertEqual(self.test_orbit.longitude, self.test_data.get("longitude"))
+
+    def test_default_altitude(self):
+        o = GeosynchronousOrbit(longitude=0.0)
+        self.assertEqual(o.altitude, 35786000)
+
+    def test_good_data_iso8601_datetime(self):
+        good_data = {
+            "altitude": 35786000,
+            "epoch": "2022-01-01T12:00:00Z",
+            "inclination": 0.0,
+            "right_ascension_ascending_node": 0.0,
+            "longitude": 10.0,
+        }
+        o = GeosynchronousOrbit(**good_data)
+        self.assertEqual(o.altitude, good_data.get("altitude"))
+        self.assertEqual(o.epoch, datetime(2022, 1, 1, 12, tzinfo=timezone.utc))
+        self.assertEqual(o.inclination, good_data.get("inclination"))
+        self.assertEqual(
+            o.right_ascension_ascending_node,
+            good_data.get("right_ascension_ascending_node"),
+        )
+        self.assertEqual(o.longitude, good_data.get("longitude"))
+
+    def test_get_true_anomaly_matches_sidereal_time(self):
+        t = constants.timescale.from_datetime(self.test_data.get("epoch"))
+        expected = (self.test_data.get("longitude") + t.gmst * 15) % 360
+        self.assertAlmostEqual(self.test_orbit.get_true_anomaly(), expected, delta=1e-6)
+
+    def test_get_derived_orbit(self):
+        derived_orbit = self.test_orbit.get_derived_orbit(20, 10)
+        self.assertAlmostEqual(
+            derived_orbit.longitude,
+            (self.test_orbit.longitude + 20) % 360,
+            delta=0.001,
+        )
+        self.assertAlmostEqual(
+            derived_orbit.right_ascension_ascending_node,
+            self.test_orbit.right_ascension_ascending_node + 10,
+            delta=0.001,
+        )
+
+    def test_to_tle(self):
+        tle = self.test_orbit.to_tle()
+        self.assertAlmostEqual(
+            tle.get_altitude(), self.test_data.get("altitude"), delta=1.0
+        )
+        self.assertAlmostEqual(
+            tle.get_epoch().timestamp(),
+            self.test_data.get("epoch").timestamp(),
+            delta=1,
+        )
+        self.assertEqual(
+            tle.get_inclination(),
+            self.test_data.get("inclination"),
+        )
+        self.assertAlmostEqual(
+            tle.get_right_ascension_ascending_node(),
+            self.test_data.get("right_ascension_ascending_node"),
+            delta=0.001,
+        )
+
+    def test_to_tle_subpoint_matches_longitude(self):
+        tle = self.test_orbit.to_tle()
+        sat = tle.as_skyfield()
+        t = constants.timescale.from_datetime(self.test_data.get("epoch"))
+        subpoint = sat.at(t).subpoint()
+        self.assertAlmostEqual(
+            subpoint.longitude.degrees, self.test_data.get("longitude"), delta=0.1
+        )
+        self.assertAlmostEqual(subpoint.latitude.degrees, 0.0, delta=0.1)
+
+    def test_to_tle_subpoint_multiple_longitudes(self):
+        for longitude in (0.0, 45.0, 137.5, 270.0):
+            with self.subTest(longitude=longitude):
+                o = GeosynchronousOrbit(
+                    epoch=self.test_data.get("epoch"),
+                    longitude=longitude,
+                    inclination=0.0,
+                )
+                tle = o.to_tle()
+                sat = tle.as_skyfield()
+                t = constants.timescale.from_datetime(self.test_data.get("epoch"))
+                subpoint = sat.at(t).subpoint()
+                sub_lon = subpoint.longitude.degrees % 360
+                self.assertAlmostEqual(sub_lon, longitude, delta=0.1)
+                self.assertAlmostEqual(subpoint.latitude.degrees, 0.0, delta=0.1)
+
+    def test_bad_longitude_negative(self):
+        with self.assertRaises(Exception):
+            GeosynchronousOrbit(longitude=-1.0)
+
+    def test_bad_longitude_too_large(self):
+        with self.assertRaises(Exception):
+            GeosynchronousOrbit(longitude=360.0)