Skip to main content

Frame

Struct Frame 

pub struct Frame {
    pub ephemeris_id: i32,
    pub orientation_id: i32,
    pub mu_km3_s2: Option<f64>,
    pub shape: Option<Ellipsoid>,
}
Expand description

A Frame uniquely defined by its ephemeris center and orientation. Refer to FrameDetail for frames combined with parameters.

:type ephemeris_id: int :type orientation_id: int :type mu_km3_s2: float, optional :type shape: Ellipsoid, optional :rtype: Frame

Fields§

§ephemeris_id: i32§orientation_id: i32§mu_km3_s2: Option<f64>

Gravity parameter of this frame, only defined on celestial frames

§shape: Option<Ellipsoid>

Shape of the geoid of this frame, only defined on geodetic frames

Implementations§

§

impl Frame

pub const fn new(ephemeris_id: i32, orientation_id: i32) -> Frame

Constructs a new frame given its ephemeris and orientations IDs, without defining anything else (so this is not a valid celestial frame, although the data could be populated later).

pub const fn from_ephem_j2000(ephemeris_id: i32) -> Frame

Examples found in repository?
examples/02_jwst_covar_monte_carlo/main.rs (line 51)
51    const JWST_J2000: Frame = Frame::from_ephem_j2000(JWST_NAIF_ID);
More examples
Hide additional examples
examples/04_lro_od/main.rs (line 81)
35fn main() -> Result<(), Box<dyn Error>> {
36    pel::init();
37
38    // ====================== //
39    // === ALMANAC SET UP === //
40    // ====================== //
41
42    // Dynamics models require planetary constants and ephemerides to be defined.
43    // Let's start by grabbing those by using ANISE's MetaAlmanac.
44
45    let data_folder: PathBuf = [env!("CARGO_MANIFEST_DIR"), "examples", "04_lro_od"]
46        .iter()
47        .collect();
48
49    let meta = data_folder.join("lro-dynamics.dhall");
50
51    // Load this ephem in the general Almanac we're using for this analysis.
52    let mut almanac = MetaAlmanac::new(meta.to_string_lossy().as_ref())
53        .map_err(Box::new)?
54        .process(true)
55        .map_err(Box::new)?;
56
57    let mut moon_pc = almanac.get_planetary_data_from_id(MOON).unwrap();
58    moon_pc.mu_km3_s2 = 4902.74987;
59    almanac.set_planetary_data_from_id(MOON, moon_pc).unwrap();
60
61    let mut earth = almanac.get_planetary_data_from_id(EARTH).unwrap();
62    earth.mu_km3_s2 = 398600.436;
63    almanac.set_planetary_data_from_id(EARTH, earth).unwrap();
64
65    // Save this new kernel for reuse.
66    // In an operational context, this would be part of the "Lock" process, and should not change throughout the mission.
67    almanac
68        .planetary_data
69        .values()
70        .next()
71        .unwrap()
72        .save_as(&data_folder.join("lro-specific.pca"), true)?;
73
74    // Lock the almanac (an Arc is a read only structure).
75    let almanac = Arc::new(almanac);
76
77    // Orbit determination requires a Trajectory structure, which can be saved as parquet file.
78    // In our case, the trajectory comes from the BSP file, so we need to build a Trajectory from the almanac directly.
79    // To query the Almanac, we need to build the LRO frame in the J2000 orientation in our case.
80    // Inspecting the LRO BSP in the ANISE GUI shows us that NASA has assigned ID -85 to LRO.
81    let lro_frame = Frame::from_ephem_j2000(-85);
82
83    // To build the trajectory we need to provide a spacecraft template.
84    let sc_template = Spacecraft::builder()
85        .mass(Mass::from_dry_and_prop_masses(1018.0, 900.0)) // Launch masses
86        .srp(SRPData {
87            // SRP configuration is arbitrary, but we will be estimating it anyway.
88            area_m2: 3.9 * 2.7,
89            coeff_reflectivity: 0.96,
90        })
91        .orbit(Orbit::zero(MOON_J2000)) // Setting a zero orbit here because it's just a template
92        .build();
93    // Now we can build the trajectory from the BSP file.
94    // We'll arbitrarily set the tracking arc to 24 hours with a five second time step.
95    let traj_as_flown = Traj::from_bsp(
96        lro_frame,
97        MOON_J2000,
98        almanac.clone(),
99        sc_template,
100        5.seconds(),
101        Some(Epoch::from_str("2024-01-01 00:00:00 UTC")?),
102        Some(Epoch::from_str("2024-01-02 00:00:00 UTC")?),
103        Aberration::LT,
104        Some("LRO".to_string()),
105    )?;
106
107    println!("{traj_as_flown}");
108
109    // ====================== //
110    // === MODEL MATCHING === //
111    // ====================== //
112
113    // Set up the spacecraft dynamics.
114
115    // Specify that the orbital dynamics must account for the graviational pull of the Earth and the Sun.
116    // The gravity of the Moon will also be accounted for since the spaceraft in a lunar orbit.
117    let mut orbital_dyn = OrbitalDynamics::point_masses(vec![EARTH, SUN, JUPITER_BARYCENTER]);
118
119    // We want to include the spherical harmonics, so let's download the gravitational data from the Nyx Cloud.
120    // We're using the GRAIL JGGRX model.
121    let mut jggrx_meta = MetaFile {
122        uri: "http://public-data.nyxspace.com/nyx/models/Luna_jggrx_1500e_sha.tab.gz".to_string(),
123        crc32: Some(0x6bcacda8), // Specifying the CRC32 avoids redownloading it if it's cached.
124    };
125    // And let's download it if we don't have it yet.
126    jggrx_meta.process(true)?;
127
128    // Build the spherical harmonics.
129    // The harmonics must be computed in the body fixed frame.
130    // We're using the long term prediction of the Moon principal axes frame.
131    let moon_pa_frame = MOON_PA_FRAME.with_orient(31008);
132    let sph_harmonics = Harmonics::from_stor(
133        almanac.frame_info(moon_pa_frame)?,
134        HarmonicsMem::from_shadr(&jggrx_meta.uri, 80, 80, true)?,
135    );
136
137    // Include the spherical harmonics into the orbital dynamics.
138    orbital_dyn.accel_models.push(sph_harmonics);
139
140    // We define the solar radiation pressure, using the default solar flux and accounting only
141    // for the eclipsing caused by the Earth and Moon.
142    // Note that by default, enabling the SolarPressure model will also enable the estimation of the coefficient of reflectivity.
143    let srp_dyn = SolarPressure::new(vec![EARTH_J2000, MOON_J2000], almanac.clone())?;
144
145    // Finalize setting up the dynamics, specifying the force models (orbital_dyn) separately from the
146    // acceleration models (SRP in this case). Use `from_models` to specify multiple accel models.
147    let dynamics = SpacecraftDynamics::from_model(orbital_dyn, srp_dyn);
148
149    println!("{dynamics}");
150
151    // Now we can build the propagator.
152    let setup = Propagator::default_dp78(dynamics.clone());
153
154    // For reference, let's build the trajectory with Nyx's models from that LRO state.
155    let (sim_final, traj_as_sim) = setup
156        .with(*traj_as_flown.first(), almanac.clone())
157        .until_epoch_with_traj(traj_as_flown.last().epoch())?;
158
159    println!("SIM INIT:  {:x}", traj_as_flown.first());
160    println!("SIM FINAL: {sim_final:x}");
161    // Compute RIC difference between SIM and LRO ephem
162    let sim_lro_delta = sim_final
163        .orbit
164        .ric_difference(&traj_as_flown.last().orbit)?;
165    println!("{traj_as_sim}");
166    println!(
167        "SIM v LRO - RIC Position (m): {:.3}",
168        sim_lro_delta.radius_km * 1e3
169    );
170    println!(
171        "SIM v LRO - RIC Velocity (m/s): {:.3}",
172        sim_lro_delta.velocity_km_s * 1e3
173    );
174
175    traj_as_sim.ric_diff_to_parquet(
176        &traj_as_flown,
177        "./data/04_output/04_lro_sim_truth_error.parquet",
178        ExportCfg::default(),
179    )?;
180
181    // ==================== //
182    // === OD SIMULATOR === //
183    // ==================== //
184
185    // After quite some time trying to exactly match the model, we still end up with an oscillatory difference on the order of 150 meters between the propagated state
186    // and the truth LRO state.
187
188    // Therefore, we will actually run an estimation from a dispersed LRO state.
189    // The sc_seed is the true LRO state from the BSP.
190    let sc_seed = *traj_as_flown.first();
191
192    // Load the Deep Space Network ground stations.
193    // Nyx allows you to build these at runtime but it's pretty static so we can just load them from YAML.
194    let ground_station_file: PathBuf = [
195        env!("CARGO_MANIFEST_DIR"),
196        "examples",
197        "04_lro_od",
198        "dsn-network.yaml",
199    ]
200    .iter()
201    .collect();
202
203    let devices = GroundStation::load_named(ground_station_file)?;
204
205    let mut proc_devices = devices.clone();
206
207    // Increase the noise in the devices to accept more measurements.
208    for gs in proc_devices.values_mut() {
209        if let Some(noise) = &mut gs
210            .stochastic_noises
211            .as_mut()
212            .unwrap()
213            .get_mut(&MeasurementType::Range)
214        {
215            *noise.white_noise.as_mut().unwrap() *= 3.0;
216        }
217    }
218
219    // Typical OD software requires that you specify your own tracking schedule or you'll have overlapping measurements.
220    // Nyx can build a tracking schedule for you based on the first station with access.
221    let trkconfg_yaml: PathBuf = [
222        env!("CARGO_MANIFEST_DIR"),
223        "examples",
224        "04_lro_od",
225        "tracking-cfg.yaml",
226    ]
227    .iter()
228    .collect();
229
230    let configs: BTreeMap<String, TrkConfig> = TrkConfig::load_named(trkconfg_yaml)?;
231
232    // Build the tracking arc simulation to generate a "standard measurement".
233    let mut trk = TrackingArcSim::<Spacecraft, GroundStation>::with_seed(
234        devices.clone(),
235        traj_as_flown.clone(),
236        configs,
237        123, // Set a seed for reproducibility
238    )?;
239
240    trk.build_schedule(almanac.clone())?;
241    let arc = trk.generate_measurements(almanac.clone())?;
242    // Save the simulated tracking data
243    arc.to_parquet_simple("./data/04_output/04_lro_simulated_tracking.parquet")?;
244
245    // We'll note that in our case, we have continuous coverage of LRO when the vehicle is not behind the Moon.
246    println!("{arc}");
247
248    // Now that we have simulated measurements, we'll run the orbit determination.
249
250    // ===================== //
251    // === OD ESTIMATION === //
252    // ===================== //
253
254    let sc = SpacecraftUncertainty::builder()
255        .nominal(sc_seed)
256        .frame(LocalFrame::RIC)
257        .x_km(0.5)
258        .y_km(0.5)
259        .z_km(0.5)
260        .vx_km_s(5e-3)
261        .vy_km_s(5e-3)
262        .vz_km_s(5e-3)
263        .build();
264
265    // Build the filter initial estimate, which we will reuse in the filter.
266    let mut initial_estimate = sc.to_estimate()?;
267    initial_estimate.covar *= 3.0;
268
269    println!("== FILTER STATE ==\n{sc_seed:x}\n{initial_estimate}");
270
271    // Build the SNC in the Moon J2000 frame, specified as a velocity noise over time.
272    let process_noise = ProcessNoise3D::from_velocity_km_s(
273        &[1e-10, 1e-10, 1e-10],
274        1 * Unit::Hour,
275        10 * Unit::Minute,
276        None,
277    );
278
279    println!("{process_noise}");
280
281    // We'll set up the OD process to reject measurements whose residuals are move than 3 sigmas away from what we expect.
282    let odp = SpacecraftKalmanOD::new(
283        setup,
284        KalmanVariant::ReferenceUpdate,
285        Some(ResidRejectCrit::default()),
286        proc_devices,
287        almanac.clone(),
288    )
289    .with_process_noise(process_noise);
290
291    let od_sol = odp.process_arc(initial_estimate, &arc)?;
292
293    let final_est = od_sol.estimates.last().unwrap();
294
295    println!("{final_est}");
296
297    let ric_err = traj_as_flown
298        .at(final_est.epoch())?
299        .orbit
300        .ric_difference(&final_est.orbital_state())?;
301    println!("== RIC at end ==");
302    println!("RIC Position (m): {:.3}", ric_err.radius_km * 1e3);
303    println!("RIC Velocity (m/s): {:.3}", ric_err.velocity_km_s * 1e3);
304
305    println!(
306        "Num residuals rejected: #{}",
307        od_sol.rejected_residuals().len()
308    );
309    println!(
310        "Percentage within +/-3: {}",
311        od_sol.residual_ratio_within_threshold(3.0).unwrap()
312    );
313    println!("Ratios normal? {}", od_sol.is_normal(None).unwrap());
314
315    od_sol.to_parquet(
316        "./data/04_output/04_lro_od_results.parquet",
317        ExportCfg::default(),
318    )?;
319
320    // Create the ephemeris
321    let ephem = od_sol.to_ephemeris("LRO rebuilt".to_string());
322    let ephem_start = ephem.start_epoch().unwrap();
323    let ephem_end = ephem.end_epoch().unwrap();
324    // Check that the covariance is PSD throughout the ephemeris by interpolating it.
325    for epoch in TimeSeries::inclusive(ephem_start, ephem_end, Unit::Minute * 5) {
326        ephem
327            .covar_at(
328                epoch,
329                anise::ephemerides::ephemeris::LocalFrame::RIC,
330                &almanac,
331            )
332            .unwrap_or_else(|e| panic!("covar not PSD at {epoch}: {e}"));
333    }
334    // Export as BSP!
335    ephem
336        .write_spice_bsp(-85, "./data/04_output/04_lro_rebuilt.bsp", None)
337        .expect("could not built BSP");
338    let new_almanac = Almanac::default()
339        .load("./data/04_output/04_lro_rebuilt.bsp")
340        .unwrap();
341    new_almanac.describe(None, None, None, None, None, None, None, None);
342    let (spk_start, spk_end) = new_almanac.spk_domain(-85).unwrap();
343
344    assert!((ephem_start - spk_start).abs() < Unit::Microsecond * 1);
345    assert!((ephem_end - spk_end).abs() < Unit::Microsecond * 1);
346
347    // In our case, we have the truth trajectory from NASA.
348    // So we can compute the RIC state difference between the real LRO ephem and what we've just estimated.
349    // Export the OD trajectory first.
350    let od_trajectory = od_sol.to_traj()?;
351    // Build the RIC difference.
352    od_trajectory.ric_diff_to_parquet(
353        &traj_as_flown,
354        "./data/04_output/04_lro_od_truth_error.parquet",
355        ExportCfg::default(),
356    )?;
357
358    Ok(())
359}

pub const fn from_orient_ssb(orientation_id: i32) -> Frame

pub fn from_name(center: &str, ref_frame: &str) -> Result<Frame, AlmanacError>

Attempts to create a new frame from its center and reference frame name. This function is compatible with the CCSDS OEM names.

pub fn with_ellipsoid(self, shape: Ellipsoid) -> Frame

Define Ellipsoid shape and return a new Frame

pub fn stripped(self) -> Frame

Returns a copy of this frame with the graviational parameter and the shape information from this frame. Use this to prevent astrodynamical computations.

:rtype: None

§

impl Frame

pub const fn with_ephem(&self, new_ephem_id: i32) -> Frame

Returns a copy of this Frame whose ephemeris ID is set to the provided ID

:type new_ephem_id: int :rtype: Frame

pub const fn with_orient(&self, new_orient_id: i32) -> Frame

Returns a copy of this Frame whose orientation ID is set to the provided ID

:type new_orient_id: int :rtype: Frame

Examples found in repository?
examples/04_lro_od/main.rs (line 131)
35fn main() -> Result<(), Box<dyn Error>> {
36    pel::init();
37
38    // ====================== //
39    // === ALMANAC SET UP === //
40    // ====================== //
41
42    // Dynamics models require planetary constants and ephemerides to be defined.
43    // Let's start by grabbing those by using ANISE's MetaAlmanac.
44
45    let data_folder: PathBuf = [env!("CARGO_MANIFEST_DIR"), "examples", "04_lro_od"]
46        .iter()
47        .collect();
48
49    let meta = data_folder.join("lro-dynamics.dhall");
50
51    // Load this ephem in the general Almanac we're using for this analysis.
52    let mut almanac = MetaAlmanac::new(meta.to_string_lossy().as_ref())
53        .map_err(Box::new)?
54        .process(true)
55        .map_err(Box::new)?;
56
57    let mut moon_pc = almanac.get_planetary_data_from_id(MOON).unwrap();
58    moon_pc.mu_km3_s2 = 4902.74987;
59    almanac.set_planetary_data_from_id(MOON, moon_pc).unwrap();
60
61    let mut earth = almanac.get_planetary_data_from_id(EARTH).unwrap();
62    earth.mu_km3_s2 = 398600.436;
63    almanac.set_planetary_data_from_id(EARTH, earth).unwrap();
64
65    // Save this new kernel for reuse.
66    // In an operational context, this would be part of the "Lock" process, and should not change throughout the mission.
67    almanac
68        .planetary_data
69        .values()
70        .next()
71        .unwrap()
72        .save_as(&data_folder.join("lro-specific.pca"), true)?;
73
74    // Lock the almanac (an Arc is a read only structure).
75    let almanac = Arc::new(almanac);
76
77    // Orbit determination requires a Trajectory structure, which can be saved as parquet file.
78    // In our case, the trajectory comes from the BSP file, so we need to build a Trajectory from the almanac directly.
79    // To query the Almanac, we need to build the LRO frame in the J2000 orientation in our case.
80    // Inspecting the LRO BSP in the ANISE GUI shows us that NASA has assigned ID -85 to LRO.
81    let lro_frame = Frame::from_ephem_j2000(-85);
82
83    // To build the trajectory we need to provide a spacecraft template.
84    let sc_template = Spacecraft::builder()
85        .mass(Mass::from_dry_and_prop_masses(1018.0, 900.0)) // Launch masses
86        .srp(SRPData {
87            // SRP configuration is arbitrary, but we will be estimating it anyway.
88            area_m2: 3.9 * 2.7,
89            coeff_reflectivity: 0.96,
90        })
91        .orbit(Orbit::zero(MOON_J2000)) // Setting a zero orbit here because it's just a template
92        .build();
93    // Now we can build the trajectory from the BSP file.
94    // We'll arbitrarily set the tracking arc to 24 hours with a five second time step.
95    let traj_as_flown = Traj::from_bsp(
96        lro_frame,
97        MOON_J2000,
98        almanac.clone(),
99        sc_template,
100        5.seconds(),
101        Some(Epoch::from_str("2024-01-01 00:00:00 UTC")?),
102        Some(Epoch::from_str("2024-01-02 00:00:00 UTC")?),
103        Aberration::LT,
104        Some("LRO".to_string()),
105    )?;
106
107    println!("{traj_as_flown}");
108
109    // ====================== //
110    // === MODEL MATCHING === //
111    // ====================== //
112
113    // Set up the spacecraft dynamics.
114
115    // Specify that the orbital dynamics must account for the graviational pull of the Earth and the Sun.
116    // The gravity of the Moon will also be accounted for since the spaceraft in a lunar orbit.
117    let mut orbital_dyn = OrbitalDynamics::point_masses(vec![EARTH, SUN, JUPITER_BARYCENTER]);
118
119    // We want to include the spherical harmonics, so let's download the gravitational data from the Nyx Cloud.
120    // We're using the GRAIL JGGRX model.
121    let mut jggrx_meta = MetaFile {
122        uri: "http://public-data.nyxspace.com/nyx/models/Luna_jggrx_1500e_sha.tab.gz".to_string(),
123        crc32: Some(0x6bcacda8), // Specifying the CRC32 avoids redownloading it if it's cached.
124    };
125    // And let's download it if we don't have it yet.
126    jggrx_meta.process(true)?;
127
128    // Build the spherical harmonics.
129    // The harmonics must be computed in the body fixed frame.
130    // We're using the long term prediction of the Moon principal axes frame.
131    let moon_pa_frame = MOON_PA_FRAME.with_orient(31008);
132    let sph_harmonics = Harmonics::from_stor(
133        almanac.frame_info(moon_pa_frame)?,
134        HarmonicsMem::from_shadr(&jggrx_meta.uri, 80, 80, true)?,
135    );
136
137    // Include the spherical harmonics into the orbital dynamics.
138    orbital_dyn.accel_models.push(sph_harmonics);
139
140    // We define the solar radiation pressure, using the default solar flux and accounting only
141    // for the eclipsing caused by the Earth and Moon.
142    // Note that by default, enabling the SolarPressure model will also enable the estimation of the coefficient of reflectivity.
143    let srp_dyn = SolarPressure::new(vec![EARTH_J2000, MOON_J2000], almanac.clone())?;
144
145    // Finalize setting up the dynamics, specifying the force models (orbital_dyn) separately from the
146    // acceleration models (SRP in this case). Use `from_models` to specify multiple accel models.
147    let dynamics = SpacecraftDynamics::from_model(orbital_dyn, srp_dyn);
148
149    println!("{dynamics}");
150
151    // Now we can build the propagator.
152    let setup = Propagator::default_dp78(dynamics.clone());
153
154    // For reference, let's build the trajectory with Nyx's models from that LRO state.
155    let (sim_final, traj_as_sim) = setup
156        .with(*traj_as_flown.first(), almanac.clone())
157        .until_epoch_with_traj(traj_as_flown.last().epoch())?;
158
159    println!("SIM INIT:  {:x}", traj_as_flown.first());
160    println!("SIM FINAL: {sim_final:x}");
161    // Compute RIC difference between SIM and LRO ephem
162    let sim_lro_delta = sim_final
163        .orbit
164        .ric_difference(&traj_as_flown.last().orbit)?;
165    println!("{traj_as_sim}");
166    println!(
167        "SIM v LRO - RIC Position (m): {:.3}",
168        sim_lro_delta.radius_km * 1e3
169    );
170    println!(
171        "SIM v LRO - RIC Velocity (m/s): {:.3}",
172        sim_lro_delta.velocity_km_s * 1e3
173    );
174
175    traj_as_sim.ric_diff_to_parquet(
176        &traj_as_flown,
177        "./data/04_output/04_lro_sim_truth_error.parquet",
178        ExportCfg::default(),
179    )?;
180
181    // ==================== //
182    // === OD SIMULATOR === //
183    // ==================== //
184
185    // After quite some time trying to exactly match the model, we still end up with an oscillatory difference on the order of 150 meters between the propagated state
186    // and the truth LRO state.
187
188    // Therefore, we will actually run an estimation from a dispersed LRO state.
189    // The sc_seed is the true LRO state from the BSP.
190    let sc_seed = *traj_as_flown.first();
191
192    // Load the Deep Space Network ground stations.
193    // Nyx allows you to build these at runtime but it's pretty static so we can just load them from YAML.
194    let ground_station_file: PathBuf = [
195        env!("CARGO_MANIFEST_DIR"),
196        "examples",
197        "04_lro_od",
198        "dsn-network.yaml",
199    ]
200    .iter()
201    .collect();
202
203    let devices = GroundStation::load_named(ground_station_file)?;
204
205    let mut proc_devices = devices.clone();
206
207    // Increase the noise in the devices to accept more measurements.
208    for gs in proc_devices.values_mut() {
209        if let Some(noise) = &mut gs
210            .stochastic_noises
211            .as_mut()
212            .unwrap()
213            .get_mut(&MeasurementType::Range)
214        {
215            *noise.white_noise.as_mut().unwrap() *= 3.0;
216        }
217    }
218
219    // Typical OD software requires that you specify your own tracking schedule or you'll have overlapping measurements.
220    // Nyx can build a tracking schedule for you based on the first station with access.
221    let trkconfg_yaml: PathBuf = [
222        env!("CARGO_MANIFEST_DIR"),
223        "examples",
224        "04_lro_od",
225        "tracking-cfg.yaml",
226    ]
227    .iter()
228    .collect();
229
230    let configs: BTreeMap<String, TrkConfig> = TrkConfig::load_named(trkconfg_yaml)?;
231
232    // Build the tracking arc simulation to generate a "standard measurement".
233    let mut trk = TrackingArcSim::<Spacecraft, GroundStation>::with_seed(
234        devices.clone(),
235        traj_as_flown.clone(),
236        configs,
237        123, // Set a seed for reproducibility
238    )?;
239
240    trk.build_schedule(almanac.clone())?;
241    let arc = trk.generate_measurements(almanac.clone())?;
242    // Save the simulated tracking data
243    arc.to_parquet_simple("./data/04_output/04_lro_simulated_tracking.parquet")?;
244
245    // We'll note that in our case, we have continuous coverage of LRO when the vehicle is not behind the Moon.
246    println!("{arc}");
247
248    // Now that we have simulated measurements, we'll run the orbit determination.
249
250    // ===================== //
251    // === OD ESTIMATION === //
252    // ===================== //
253
254    let sc = SpacecraftUncertainty::builder()
255        .nominal(sc_seed)
256        .frame(LocalFrame::RIC)
257        .x_km(0.5)
258        .y_km(0.5)
259        .z_km(0.5)
260        .vx_km_s(5e-3)
261        .vy_km_s(5e-3)
262        .vz_km_s(5e-3)
263        .build();
264
265    // Build the filter initial estimate, which we will reuse in the filter.
266    let mut initial_estimate = sc.to_estimate()?;
267    initial_estimate.covar *= 3.0;
268
269    println!("== FILTER STATE ==\n{sc_seed:x}\n{initial_estimate}");
270
271    // Build the SNC in the Moon J2000 frame, specified as a velocity noise over time.
272    let process_noise = ProcessNoise3D::from_velocity_km_s(
273        &[1e-10, 1e-10, 1e-10],
274        1 * Unit::Hour,
275        10 * Unit::Minute,
276        None,
277    );
278
279    println!("{process_noise}");
280
281    // We'll set up the OD process to reject measurements whose residuals are move than 3 sigmas away from what we expect.
282    let odp = SpacecraftKalmanOD::new(
283        setup,
284        KalmanVariant::ReferenceUpdate,
285        Some(ResidRejectCrit::default()),
286        proc_devices,
287        almanac.clone(),
288    )
289    .with_process_noise(process_noise);
290
291    let od_sol = odp.process_arc(initial_estimate, &arc)?;
292
293    let final_est = od_sol.estimates.last().unwrap();
294
295    println!("{final_est}");
296
297    let ric_err = traj_as_flown
298        .at(final_est.epoch())?
299        .orbit
300        .ric_difference(&final_est.orbital_state())?;
301    println!("== RIC at end ==");
302    println!("RIC Position (m): {:.3}", ric_err.radius_km * 1e3);
303    println!("RIC Velocity (m/s): {:.3}", ric_err.velocity_km_s * 1e3);
304
305    println!(
306        "Num residuals rejected: #{}",
307        od_sol.rejected_residuals().len()
308    );
309    println!(
310        "Percentage within +/-3: {}",
311        od_sol.residual_ratio_within_threshold(3.0).unwrap()
312    );
313    println!("Ratios normal? {}", od_sol.is_normal(None).unwrap());
314
315    od_sol.to_parquet(
316        "./data/04_output/04_lro_od_results.parquet",
317        ExportCfg::default(),
318    )?;
319
320    // Create the ephemeris
321    let ephem = od_sol.to_ephemeris("LRO rebuilt".to_string());
322    let ephem_start = ephem.start_epoch().unwrap();
323    let ephem_end = ephem.end_epoch().unwrap();
324    // Check that the covariance is PSD throughout the ephemeris by interpolating it.
325    for epoch in TimeSeries::inclusive(ephem_start, ephem_end, Unit::Minute * 5) {
326        ephem
327            .covar_at(
328                epoch,
329                anise::ephemerides::ephemeris::LocalFrame::RIC,
330                &almanac,
331            )
332            .unwrap_or_else(|e| panic!("covar not PSD at {epoch}: {e}"));
333    }
334    // Export as BSP!
335    ephem
336        .write_spice_bsp(-85, "./data/04_output/04_lro_rebuilt.bsp", None)
337        .expect("could not built BSP");
338    let new_almanac = Almanac::default()
339        .load("./data/04_output/04_lro_rebuilt.bsp")
340        .unwrap();
341    new_almanac.describe(None, None, None, None, None, None, None, None);
342    let (spk_start, spk_end) = new_almanac.spk_domain(-85).unwrap();
343
344    assert!((ephem_start - spk_start).abs() < Unit::Microsecond * 1);
345    assert!((ephem_end - spk_end).abs() < Unit::Microsecond * 1);
346
347    // In our case, we have the truth trajectory from NASA.
348    // So we can compute the RIC state difference between the real LRO ephem and what we've just estimated.
349    // Export the OD trajectory first.
350    let od_trajectory = od_sol.to_traj()?;
351    // Build the RIC difference.
352    od_trajectory.ric_diff_to_parquet(
353        &traj_as_flown,
354        "./data/04_output/04_lro_od_truth_error.parquet",
355        ExportCfg::default(),
356    )?;
357
358    Ok(())
359}

pub const fn is_celestial(&self) -> bool

Returns whether this is a celestial frame

:rtype: bool

pub const fn is_geodetic(&self) -> bool

Returns whether this is a geodetic frame

:rtype: bool

pub const fn ephem_origin_id_match(&self, other_id: i32) -> bool

Returns true if the ephemeris origin is equal to the provided ID

:type other_id: int :rtype: bool

pub const fn orient_origin_id_match(&self, other_id: i32) -> bool

Returns true if the orientation origin is equal to the provided ID

:type other_id: int :rtype: bool

pub const fn ephem_origin_match(&self, other: Frame) -> bool

Returns true if the ephemeris origin is equal to the provided frame

:type other: Frame :rtype: bool

pub const fn orient_origin_match(&self, other: Frame) -> bool

Returns true if the orientation origin is equal to the provided frame

:type other: Frame :rtype: bool

pub fn strip(&mut self)

Removes the graviational parameter and the shape information from this frame. Use this to prevent astrodynamical computations.

:rtype: None

pub fn mu_km3_s2(&self) -> Result<f64, PhysicsError>

Returns the gravitational parameters of this frame, if defined

:rtype: float

pub fn with_mu_km3_s2(&self, mu_km3_s2: f64) -> Frame

Returns a copy of this frame with the graviational parameter set to the new value.

:type mu_km3_s2: float :rtype: Frame

pub fn mean_equatorial_radius_km(&self) -> Result<f64, PhysicsError>

Returns the mean equatorial radius in km, if defined

:rtype: float

pub fn semi_major_radius_km(&self) -> Result<f64, PhysicsError>

Returns the semi major radius of the tri-axial ellipoid shape of this frame, if defined

:rtype: float

pub fn flattening(&self) -> Result<f64, PhysicsError>

Returns the flattening ratio (unitless)

:rtype: float

pub fn polar_radius_km(&self) -> Result<f64, PhysicsError>

Returns the polar radius in km, if defined

:rtype: float

Trait Implementations§

§

impl Clone for Frame

§

fn clone(&self) -> Frame

Returns a duplicate of the value. Read more
1.0.0 · Source§

fn clone_from(&mut self, source: &Self)

Performs copy-assignment from source. Read more
§

impl Debug for Frame

§

fn fmt(&self, f: &mut Formatter<'_>) -> Result<(), Error>

Formats the value using the given formatter. Read more
§

impl<'de> Deserialize<'de> for Frame

§

fn deserialize<__D>( __deserializer: __D, ) -> Result<Frame, <__D as Deserializer<'de>>::Error>
where __D: Deserializer<'de>,

Deserialize this value from the given Serde deserializer. Read more
§

impl Display for Frame

§

fn fmt(&self, f: &mut Formatter<'_>) -> Result<(), Error>

Formats the value using the given formatter. Read more
§

impl From<&Frame> for FrameUid

§

fn from(frame: &Frame) -> FrameUid

Converts to this type from the input type.
§

impl From<&FrameUid> for Frame

§

fn from(uid: &FrameUid) -> Frame

Converts to this type from the input type.
§

impl From<Frame> for FrameUid

§

fn from(frame: Frame) -> FrameUid

Converts to this type from the input type.
§

impl From<FrameUid> for Frame

§

fn from(uid: FrameUid) -> Frame

Converts to this type from the input type.
§

impl LowerExp for Frame

§

fn fmt(&self, f: &mut Formatter<'_>) -> Result<(), Error>

Only prints the ephemeris name

§

impl LowerHex for Frame

§

fn fmt(&self, f: &mut Formatter<'_>) -> Result<(), Error>

Only prints the UID

§

impl Octal for Frame

§

fn fmt(&self, f: &mut Formatter<'_>) -> Result<(), Error>

Only prints the orientation name

§

impl PartialEq for Frame

§

fn eq(&self, other: &Frame) -> bool

Tests for self and other values to be equal, and is used by ==.
1.0.0 · Source§

fn ne(&self, other: &Rhs) -> bool

Tests for !=. The default implementation is almost always sufficient, and should not be overridden without very good reason.
§

impl Serialize for Frame

§

fn serialize<__S>( &self, __serializer: __S, ) -> Result<<__S as Serializer>::Ok, <__S as Serializer>::Error>
where __S: Serializer,

Serialize this value into the given Serde serializer. Read more
§

impl StaticType for Frame

§

fn static_type() -> SimpleType

Return the Dhall type that represents this type. Read more
§

impl Copy for Frame

§

impl StructuralPartialEq for Frame

Auto Trait Implementations§

§

impl Freeze for Frame

§

impl RefUnwindSafe for Frame

§

impl Send for Frame

§

impl Sync for Frame

§

impl Unpin for Frame

§

impl UnwindSafe for Frame

Blanket Implementations§

Source§

impl<T> Any for T
where T: 'static + ?Sized,

Source§

fn type_id(&self) -> TypeId

Gets the TypeId of self. Read more
Source§

impl<T> Borrow<T> for T
where T: ?Sized,

Source§

fn borrow(&self) -> &T

Immutably borrows from an owned value. Read more
Source§

impl<T> BorrowMut<T> for T
where T: ?Sized,

Source§

fn borrow_mut(&mut self) -> &mut T

Mutably borrows from an owned value. Read more
Source§

impl<T> CloneToUninit for T
where T: Clone,

Source§

unsafe fn clone_to_uninit(&self, dest: *mut u8)

🔬This is a nightly-only experimental API. (clone_to_uninit)
Performs copy-assignment from self to dest. Read more
Source§

impl<T> From<T> for T

Source§

fn from(t: T) -> T

Returns the argument unchanged.

Source§

impl<T> FromDhall for T

Source§

fn from_dhall(v: &Value) -> Result<T, Error>

§

impl<T> Instrument for T

§

fn instrument(self, span: Span) -> Instrumented<Self>

Instruments this type with the provided [Span], returning an Instrumented wrapper. Read more
§

fn in_current_span(self) -> Instrumented<Self>

Instruments this type with the current Span, returning an Instrumented wrapper. Read more
Source§

impl<T, U> Into<U> for T
where U: From<T>,

Source§

fn into(self) -> U

Calls U::from(self).

That is, this conversion is whatever the implementation of From<T> for U chooses to do.

Source§

impl<T> IntoEither for T

Source§

fn into_either(self, into_left: bool) -> Either<Self, Self>

Converts self into a Left variant of Either<Self, Self> if into_left is true. Converts self into a Right variant of Either<Self, Self> otherwise. Read more
Source§

fn into_either_with<F>(self, into_left: F) -> Either<Self, Self>
where F: FnOnce(&Self) -> bool,

Converts self into a Left variant of Either<Self, Self> if into_left(&self) returns true. Converts self into a Right variant of Either<Self, Self> otherwise. Read more
§

impl<T> Pointable for T

§

const ALIGN: usize

The alignment of pointer.
§

type Init = T

The type for initializers.
§

unsafe fn init(init: <T as Pointable>::Init) -> usize

Initializes a with the given initializer. Read more
§

unsafe fn deref<'a>(ptr: usize) -> &'a T

Dereferences the given pointer. Read more
§

unsafe fn deref_mut<'a>(ptr: usize) -> &'a mut T

Mutably dereferences the given pointer. Read more
§

unsafe fn drop(ptr: usize)

Drops the object pointed to by the given pointer. Read more
Source§

impl<T> Same for T

Source§

type Output = T

Should always be Self
§

impl<SS, SP> SupersetOf<SS> for SP
where SS: SubsetOf<SP>,

§

fn to_subset(&self) -> Option<SS>

The inverse inclusion map: attempts to construct self from the equivalent element of its superset. Read more
§

fn is_in_subset(&self) -> bool

Checks if self is actually part of its subset T (and can be converted to it).
§

fn to_subset_unchecked(&self) -> SS

Use with care! Same as self.to_subset but without any property checks. Always succeeds.
§

fn from_subset(element: &SS) -> SP

The inclusion map: converts self to the equivalent element of its superset.
Source§

impl<T> ToDhall for T
where T: Serialize,

Source§

fn to_dhall(&self, ty: Option<&SimpleType>) -> Result<Value, Error>

Source§

impl<T> ToOwned for T
where T: Clone,

Source§

type Owned = T

The resulting type after obtaining ownership.
Source§

fn to_owned(&self) -> T

Creates owned data from borrowed data, usually by cloning. Read more
Source§

fn clone_into(&self, target: &mut T)

Uses borrowed data to replace owned data, usually by cloning. Read more
Source§

impl<T> ToString for T
where T: Display + ?Sized,

Source§

fn to_string(&self) -> String

Converts the given value to a String. Read more
Source§

impl<T, U> TryFrom<U> for T
where U: Into<T>,

Source§

type Error = Infallible

The type returned in the event of a conversion error.
Source§

fn try_from(value: U) -> Result<T, <T as TryFrom<U>>::Error>

Performs the conversion.
Source§

impl<T, U> TryInto<U> for T
where U: TryFrom<T>,

Source§

type Error = <U as TryFrom<T>>::Error

The type returned in the event of a conversion error.
Source§

fn try_into(self) -> Result<U, <U as TryFrom<T>>::Error>

Performs the conversion.
§

impl<V, T> VZip<V> for T
where V: MultiLane<T>,

§

fn vzip(self) -> V

§

impl<T> WithSubscriber for T

§

fn with_subscriber<S>(self, subscriber: S) -> WithDispatch<Self>
where S: Into<Dispatch>,

Attaches the provided Subscriber to this type, returning a [WithDispatch] wrapper. Read more
§

fn with_current_subscriber(self) -> WithDispatch<Self>

Attaches the current default Subscriber to this type, returning a [WithDispatch] wrapper. Read more
§

impl<T> Allocation for T
where T: RefUnwindSafe + Send + Sync,

Source§

impl<T> DeserializeOwned for T
where T: for<'de> Deserialize<'de>,

Source§

impl<T> Scalar for T
where T: 'static + Clone + PartialEq + Debug,

Source§

impl<T> Scalar for T
where T: 'static + Clone + PartialEq + Debug,