diff --git a/docs/nyxspace/tutorials/index.md b/docs/nyxspace/tutorials/index.md new file mode 100644 index 0000000..d68b919 --- /dev/null +++ b/docs/nyxspace/tutorials/index.md @@ -0,0 +1,9 @@ +# Nyx Tutorials + +Welcome to the Nyx tutorials section! These guides provide practical, step-by-step instructions on how to use Nyx Space for various astrodynamics and orbit determination tasks in Python. + +## Tutorials + +* [Mission Design](mission_design.md): Learn how to propagate orbits with perturbations and execute simple Monte Carlo simulations. +* [Orbit Determination](orbit_determination.md): Understand how to execute an extended Kalman filter (EKF) and smoother with tracking data. +* [Spacecraft Definition](spacecraft.md): Explore the basic definitions of spacecraft, including properties like mass, drag, solar radiation pressure (SRP), and defining dispersed spacecraft for Monte Carlo simulations. diff --git a/docs/nyxspace/tutorials/mission_design.md b/docs/nyxspace/tutorials/mission_design.md new file mode 100644 index 0000000..e4c6a72 --- /dev/null +++ b/docs/nyxspace/tutorials/mission_design.md @@ -0,0 +1,203 @@ +# Mission Design + +This tutorial demonstrates how to configure propagators, handle perturbations like atmospheric drag and solar radiation pressure (SRP), use spherical harmonics, and perform simple Monte Carlo propagations in Nyx. + +## Propagate with Perturbations + +Here is a guide on how to propagate a spacecraft's orbit over time while accounting for Solar Radiation Pressure (SRP) and atmospheric drag. Note that this demonstrates simple propagation segments, not the more complex sequencing available in Nyx Premium. + +```python +import logging +from nyx_space import DragData, ExportCfg, Mass, Spacecraft, SRPData +from nyx_space.anise import MetaAlmanac +from nyx_space.anise.analysis import Condition, Event, OrbitalElement, ScalarExpr +from nyx_space.anise.astro import Orbit +from nyx_space.anise.constants import CelestialObjects, Frames +from nyx_space.mission_design import ( + AccelModels, + AtmDensity, + Drag, + Dynamics, + ForceModels, + GravityFieldConfig, + IntegratorMethod, + IntegratorOptions, + PointMasses, + Propagator, + SolarPressure, +) +from nyx_space.monte_carlo import MvnSpacecraft, StateDispersion, StateParameter +from nyx_space.time import Duration, Epoch, Unit + +# Step 1: Define Integrator Options +# We default to a variable step integrator. This ensures numerical accuracy +# without wasting CPU cycles on regions of the orbit where dynamics change slowly. +opts = IntegratorOptions() + +# Step 2: Define the Environment +# Define which celestial bodies exert a point-mass gravitational pull. +accel_models = AccelModels( + point_masses=PointMasses( + celestial_objects=[ + CelestialObjects.EARTH, + CelestialObjects.MOON, + CelestialObjects.SUN, + ] + ) +) + +# Step 3: Configure Non-Gravitational Force Models +almanac = MetaAlmanac.latest() + +# Define Atmospheric Drag. +drag = Drag( + AtmDensity.earth_exponential(), + almanac.frame_info(Frames.IAU_EARTH_FRAME), + estimate=False, +) + +# Define Solar Radiation Pressure (SRP). +# We configure this to account for eclipses cast by the Earth and the Moon. +srp = SolarPressure( + [Frames.EARTH_J2000, Frames.MOON_J2000], almanac, estimate=False +) +force_models = ForceModels(drag=drag, srp=srp) + +# Step 4: Assemble the Dynamics +# The dynamics aggregate the gravitational and non-gravitational models. +dynamics = Dynamics(accel_models, force_models) + +# Initialize the propagator with the defined dynamics and the Runge-Kutta 8(9) integrator. +# The Almanac provides the necessary ephemerides for the celestial bodies and frame transformations. +propagator = Propagator(dynamics, almanac, IntegratorMethod.RungeKutta89, opts) + +# Step 5: Define the Initial Spacecraft State +eme2k = almanac.frame_info(Frames.EME2000) +orbit = Orbit.from_keplerian( + sma_km=6800.0, + ecc=1e-4, + inc_deg=45.0, + raan_deg=60.0, + aop_deg=75.0, + ta_deg=90.0, + epoch=Epoch("2030-12-29 01:02:03 TDB"), + frame=eme2k, +) + +# Define the physical characteristics of the spacecraft necessary to compute drag and SRP. +my_sc = Spacecraft( + orbit, + Mass(dry_mass_kg=123.0), + SRPData(area_m2=10.0, coeff_reflectivity=1.2), + DragData(area_m2=10.0, coeff_drag=2.0), +) + +# Step 6: Execute Propagations +# Scenario A: Propagate until a specific temporal epoch. +target_epoch = Epoch("2030-12-30 01:02:03 UTC") +result_epoch = propagator.until_epoch(my_sc, target_epoch, trajectory=True) + +# The trajectory may be converted to an ANISE Ephemeris, e.g. for serialization to a BSP/OEM file. +ephem = result_epoch.trajectory.to_ephemeris("MySpacecraftTraj", ExportCfg()) + +# Scenario B: Propagate until a specific orbital event (apoapsis). +# This demonstrates the root-finding capabilities of Nyx to stop exactly at an event. +result_event = propagator.until_event( + my_sc, Event.apoapsis(), max_duration=Duration("1 day") +) +``` + +## Configure Spherical Harmonics + +This example configures a high-fidelity spherical harmonic gravity field model for precise near-Earth propagation. + +```python +# Assuming 'almanac' is loaded +eme2k = almanac.frame_info(Frames.EME2000) + +# Define a low-Earth orbit using Cartesian state vectors +orbit = Orbit( + 5442.1625926801835, + -4068.9498468206248, + -13.456851447751518, + 2.8581975428173836, + 3.8097859312745794, + 6.002126693122689, + Epoch("2025-08-25 11:55:44 UTC"), + eme2k, +) +spacecraft = Spacecraft(orbit) + +# Step 1: Define the Gravity Field Configuration +# We apply the EGM2008 gravity model (degree 10, order 10). +# The frame must be Earth-fixed (e.g., ITRF93) for the field to rotate with the planet. +accel_models = AccelModels( + point_masses=PointMasses( + celestial_objects=[CelestialObjects.MOON, CelestialObjects.SUN] + ), + gravity_field=GravityFieldConfig( + degree=10, + order=10, + filepath="../data/01_planetary/EGM2008_to2190_TideFree.gz", + frame=Frames.EARTH_ITRF93.to_frameuid(), + ), +) + +# Step 2: Propagate +# For high-precision gravity fields, Dormand-Prince 7-8 is a robust choice. +dynamics = Dynamics(accel_models) +propagator = Propagator(dynamics, almanac, IntegratorMethod.DormandPrince78) + +final_state = propagator.for_duration(spacecraft, Unit.Day * 1) +``` + +## Execute Simple Monte Carlo + +Generate a multi-variate dispersion of initial orbital states and propagate them in parallel to analyze trajectory divergence. + +```python +import polars as pl + +# Define the multivariate normal dispersion +disp = [ + StateDispersion.zero_mean( + StateParameter.Element(OrbitalElement.SemiMajorAxis), 15.0 + ), + StateDispersion.zero_mean( + StateParameter.Element(OrbitalElement.Eccentricity), 1e-6 + ), +] + +# Generate dispersed variants based on the distributions +mvn = MvnSpacecraft(spacecraft, disp) +dispersed_spacecraft = mvn.sample(100, seed=123) + +# Configure the propagator +accel_models = AccelModels( + point_masses=PointMasses( + celestial_objects=[CelestialObjects.EARTH, CelestialObjects.MOON] + ), +) +srp = SolarPressure([Frames.EARTH_J2000], almanac, estimate=False) +dynamics = Dynamics(accel_models, force_models=ForceModels(srp)) +propagator = Propagator(dynamics, almanac, IntegratorMethod.RungeKutta89) + +# Execute the single segment Monte Carlo +# Propagate each spacecraft until it hits its 25th occurrence of Local Solar Time equaling 6.0 hours. +results = propagator.many_until_event( + dispersed_spacecraft, + Event(ScalarExpr.LocalSolarTime(), Condition.Equals(6.0), Duration("1 ms")), + trajectory=False, + trigger=25, + max_duration=Duration("1 day"), +) + +# Analyze Results +df_data = {"Epoch (UTC)": [], "LTAN (deg)": []} +for rslt in results: + df_data["Epoch (UTC)"] += [rslt.state.orbit.epoch.to_datetime()] + df_data["LTAN (deg)"] += [rslt.state.orbit.ltan_deg()] + +df = pl.DataFrame(df_data) +print(df.describe()) +``` diff --git a/docs/nyxspace/tutorials/orbit_determination.md b/docs/nyxspace/tutorials/orbit_determination.md new file mode 100644 index 0000000..82179ca --- /dev/null +++ b/docs/nyxspace/tutorials/orbit_determination.md @@ -0,0 +1,158 @@ +# Orbit Determination + +This tutorial explains how to run tracking data through a scalar orbit determination filter to estimate a spacecraft's position, velocity, and coefficient of reflectivity. + +## Executing an Orbit Determination Filter + +```python +from nyx_space import ExportCfg, Spacecraft +from nyx_space.anise import MetaAlmanac +from nyx_space.anise.analysis import OrbitalElement +from nyx_space.anise.astro import Orbit +from nyx_space.anise.constants import CelestialObjects, Frames +from nyx_space.mission_design import ( + AccelModels, + Dynamics, + GravityFieldConfig, + PointMasses, + Propagator, +) +from nyx_space.monte_carlo import StateDispersion, StateParameter +from nyx_space.orbit_determination import ( + GroundStation, + GroundTrackingArcSim, + Handoff, + KalmanVariant, + LocalFrame, + Location, + MeasurementType, + ProcessNoise, + Scheduler, + SigmaRejection, + SpacecraftEstimate, + SpacecraftODProcess, + StochasticNoise, + TrkConfig, +) +from nyx_space.time import Epoch, Unit + +almanac = MetaAlmanac.latest() +eme2k = almanac.frame_info(Frames.EME2000) + +# Step 1: Build a reference trajectory +orbit = Orbit( + 5442.1625926801835, + -4068.9498468206248, + -13.456851447751518, + 2.8581975428173836, + 3.8097859312745794, + 6.002126693122689, + Epoch("2025-08-25 11:55:44 UTC"), + eme2k, +) +spacecraft = Spacecraft(orbit) + +# Step 2: Disperse this spacecraft state to create an initial filter error +disp = [ + StateDispersion.zero_mean( + StateParameter.Element(OrbitalElement.SemiMajorAxis), 1.0 + ), + StateDispersion.zero_mean( + StateParameter.Element(OrbitalElement.Eccentricity), 1e-6 + ), +] + +# Builds a zero mean estimate from these dispersions +estimate = SpacecraftEstimate.from_dispersions( + nominal_state=spacecraft, dispersions=disp, seed=123 +) + +# Propagate to generate the ground-truth trajectory. +accel_models = AccelModels( + point_masses=PointMasses( + celestial_objects=[CelestialObjects.MOON, CelestialObjects.SUN] + ), + gravity_field=GravityFieldConfig( + degree=10, + order=10, + filepath="../data/01_planetary/EGM2008_to2190_TideFree.gz", + frame=Frames.IAU_EARTH_FRAME.to_frameuid(), + ), +) +dynamics = Dynamics(accel_models) +propagator = Propagator(dynamics, almanac) +traj = propagator.for_duration(spacecraft, Unit.Day * 1, trajectory=True).trajectory + +# Step 3: Configure Stochastic Measurement Noise, generate tracking data. +stochastic_noises = { + MeasurementType.Range: StochasticNoise(name="Range"), + MeasurementType.Doppler: StochasticNoise(name="Doppler"), +} + +# Step 4: Build the ground station network +gs0 = GroundStation( + name="Paris, FR", + location=Location( + latitude_deg=48.8566, + longitude_deg=2.3522, + height_km=0.4, + frame=Frames.IAU_EARTH_FRAME.to_frameuid(), + terrain_mask=[], + terrain_mask_ignored=True, + ), + stochastic_noises=stochastic_noises, +) + +gs1 = GroundStation( + "Denver, CO", + Location( + 39.7420, -104.9915, 1.8, Frames.IAU_EARTH_FRAME.to_frameuid(), [], True + ), + stochastic_noises, +) + +# Step 5: Define Tracking Schedules +configs = { + "Paris, FR": TrkConfig(Scheduler(Handoff.Greedy)), + "Denver, CO": TrkConfig(Scheduler(Handoff.Greedy)), +} + +network = {"Paris, FR": gs0, "Denver, CO": gs1} +trk_sim = GroundTrackingArcSim(network, traj, configs) +trk_sim.build_schedule(almanac) + +# Step 6: Generate synthetic measurements +trk_arc = trk_sim.generate_measurements(almanac) + +# Step 7: Run the orbit determination filter +od_proc_deviation = SpacecraftODProcess( + propagator, KalmanVariant.DeviationTracking, network +) +od_dev_sol = od_proc_deviation.process_arc(estimate, trk_arc) + +# Export the whole orbit determination solution into a single Parquet file. +od_dev_sol.to_parquet("od_dev.pq", ExportCfg(False)) + +# Step 8: Running the smoother +smoothed = od_dev_sol.smooth(almanac) + +# Reference update algorithm setup +process_noise = ProcessNoise.from_velocity_m_s( + vx_m_s=1e-6, + vy_m_s=1e-6, + vz_m_s=1e-6, + noise_duration=Unit.Second * 1, + disable_time=Unit.Hour * 2, + local_frame=LocalFrame.Inertial, +) +od_proc = SpacecraftODProcess( + propagator, KalmanVariant.ReferenceUpdate, network, process_noise=process_noise +) +od_sol = od_proc.process_arc(estimate, trk_arc) +od_sol.to_parquet("od_ref_update.pq", ExportCfg(False)) + +# Export definitive ephemeris to OEM +definitive_ephem = od_sol.to_ephemeris("Test OD Spacecraft") +oem_filepath = "definitive_ephem.oem" +definitive_ephem.write_ccsds_oem(oem_filepath) +``` diff --git a/docs/nyxspace/tutorials/spacecraft.md b/docs/nyxspace/tutorials/spacecraft.md new file mode 100644 index 0000000..57306a1 --- /dev/null +++ b/docs/nyxspace/tutorials/spacecraft.md @@ -0,0 +1,58 @@ +# Spacecraft Definition + +This tutorial explores the basic definitions of spacecraft, including properties like mass, drag, solar radiation pressure (SRP), and defining dispersed spacecraft for Monte Carlo simulations. + +## Defining a Spacecraft + +```python +from nyx_space import DragData, Mass, Spacecraft, SRPData +from nyx_space.anise import MetaAlmanac +from nyx_space.anise.analysis import OrbitalElement +from nyx_space.anise.astro import Orbit +from nyx_space.anise.constants import Frames +from nyx_space.anise.time import Epoch +from nyx_space.monte_carlo import MvnSpacecraft, StateDispersion, StateParameter + +almanac = MetaAlmanac.latest() +eme2k = almanac.frame_info(Frames.EME2000) + +# Define an orbit +orbit = Orbit.from_keplerian( + 6800.0, 1e-4, 45.0, 60.0, 75.0, 90.0, Epoch("2020-02-29 01:02:03 TDB"), eme2k +) + +# A simple spacecraft with just an orbit +sc = Spacecraft(orbit) + +# A spacecraft with defined mass (in kg) +sc_mass = Spacecraft(orbit, Mass(123.0)) + +# A spacecraft with mass, SRP (area, coefficient of reflectivity), and Drag (area, drag coefficient) +sc_mass_srp_drag = Spacecraft( + orbit, Mass(123.0), SRPData(10.0, 1.2), DragData(10.0, 2.0) +) + +# Test serialization and deserialization to ASN1 format +sc_mass2 = Spacecraft.from_asn1(sc_mass.to_asn1()) +sc_mass_srp_drag2 = Spacecraft.from_asn1(sc_mass_srp_drag.to_asn1()) +``` + +## Multivariate Normal Spacecraft + +Nyx supports defining dispersions for Monte Carlo simulations using a multivariate normal distribution. + +```python +# Define state dispersions +disp = [ + StateDispersion.zero_mean( + StateParameter.Element(OrbitalElement.SemiMajorAxis), 15.0 + ), + StateDispersion.zero_mean(StateParameter.Element(OrbitalElement.RAAN), 5.0), +] + +# Create a Multivariate Normal Spacecraft +mvn = MvnSpacecraft(sc, disp) + +# Sample 1000 dispersed spacecraft based on the defined distributions +dispersed = mvn.sample(1000, 123) +``` diff --git a/mkdocs.yml b/mkdocs.yml index 65f3304..e189ad1 100644 --- a/mkdocs.yml +++ b/mkdocs.yml @@ -46,6 +46,11 @@ nav: # - Orbit definition and manipulation [WIP]: nyxspace/user_guide/orbits.md # - Propagators and Monte Carlo runs [WIP]: nyxspace/user_guide/prop-mc.md # - High fidelity models [WIP]: nyxspace/user_guide/hifimodels.md + - Tutorials: + - Tutorials index: nyxspace/tutorials/index.md + - Mission design: nyxspace/tutorials/mission_design.md + - Orbit determination: nyxspace/tutorials/orbit_determination.md + - Spacecraft: nyxspace/tutorials/spacecraft.md - MathSpec: - Introduction: nyxspace/MathSpec/index.md - Celestial: