diff --git a/crates/lox-orbits/src/analysis.rs b/crates/lox-orbits/src/analysis.rs index e6307ddd..a44be92e 100644 --- a/crates/lox-orbits/src/analysis.rs +++ b/crates/lox-orbits/src/analysis.rs @@ -1,3 +1,5 @@ +use std::f64::consts::PI; + /* * Copyright (c) 2024. Helge Eichhorn and the LOX contributors * @@ -5,15 +7,16 @@ * License, v. 2.0. If a copy of the MPL was not distributed with this * file, you can obtain one at https://mozilla.org/MPL/2.0/. */ - use lox_bodies::{RotationalElements, Spheroid}; use lox_math::roots::Brent; +use lox_math::series::{Series, SeriesError}; use lox_math::types::units::Radians; use lox_time::deltas::TimeDelta; use lox_time::julian_dates::JulianDate; use lox_time::time_scales::Tdb; use lox_time::transformations::TryToScale; use lox_time::TimeLike; +use thiserror::Error; use crate::events::{find_windows, Window}; use crate::frames::{BodyFixed, FrameTransformationProvider, Icrf, Topocentric, TryToFrame}; @@ -21,6 +24,44 @@ use crate::ground::GroundLocation; use crate::origins::{CoordinateOrigin, Origin}; use crate::trajectories::Trajectory; +#[derive(Debug, Clone, Error, PartialEq)] +pub enum ElevationMaskError { + #[error("invalid azimuth range: {}..{}", .0.to_degrees(), .1.to_degrees())] + InvalidAzimuthRange(f64, f64), + #[error("series error")] + SeriesError(#[from] SeriesError), +} + +#[derive(Debug, Clone, PartialEq)] +pub enum ElevationMask { + Fixed(f64), + Variable(Series, Vec>), +} + +impl ElevationMask { + pub fn new(azimuth: Vec, elevation: Vec) -> Result { + if !azimuth.is_empty() { + let az_min = *azimuth.iter().min_by(|a, b| a.total_cmp(b)).unwrap(); + let az_max = *azimuth.iter().max_by(|a, b| a.total_cmp(b)).unwrap(); + if az_min != -PI || az_max != PI { + return Err(ElevationMaskError::InvalidAzimuthRange(az_min, az_max)); + } + } + Ok(Self::Variable(Series::new(azimuth, elevation)?)) + } + + pub fn with_fixed_elevation(elevation: f64) -> Self { + Self::Fixed(elevation) + } + + pub fn min_elevation(&self, azimuth: f64) -> f64 { + match self { + ElevationMask::Fixed(min_elevation) => *min_elevation, + ElevationMask::Variable(series) => series.interpolate(azimuth), + } + } +} + pub fn elevation< T: TimeLike + TryToScale + Clone, O: Origin + Spheroid + RotationalElements + Clone, @@ -56,6 +97,7 @@ pub fn elevation2< >( time: T, gs: &GroundLocation, + mask: &ElevationMask, sc: &Trajectory, provider: &P, ) -> Radians { @@ -63,7 +105,7 @@ pub fn elevation2< let sc = sc.interpolate_at(time.clone()); let sc = sc.try_to_frame(body_fixed, provider).unwrap(); let obs = gs.observables(sc); - obs.elevation() + obs.elevation() - mask.min_elevation(obs.azimuth()) } pub fn visibility< @@ -72,8 +114,8 @@ pub fn visibility< P: FrameTransformationProvider, >( times: &[T], - min_elevation: Radians, gs: &GroundLocation, + mask: &ElevationMask, sc: &Trajectory, provider: &P, ) -> Vec> { @@ -92,9 +134,10 @@ pub fn visibility< elevation2( start.clone() + TimeDelta::from_decimal_seconds(t).unwrap(), gs, + mask, sc, provider, - ) - min_elevation + ) }, start.clone(), end.clone(), @@ -137,13 +180,33 @@ mod tests { } } + #[test] + fn test_elevation_mask() { + let azimuth = vec![-PI, 0.0, PI]; + let elevation = vec![-2.0, 0.0, 2.0]; + let mask = ElevationMask::new(azimuth, elevation).unwrap(); + assert_eq!(mask.min_elevation(0.0), 0.0); + } + + #[test] + fn test_elevation_mask_invalid_mask() { + let azimuth = vec![-PI, 0.0, PI / 2.0]; + let elevation = vec![-2.0, 0.0, 2.0]; + let mask = ElevationMask::new(azimuth, elevation); + assert_eq!( + mask, + Err(ElevationMaskError::InvalidAzimuthRange(-PI, PI / 2.0)) + ) + } + #[test] fn test_visibility() { let gs = location(); + let mask = ElevationMask::with_fixed_elevation(0.0); let sc = spacecraft_trajectory(); let times: Vec> = sc.states().iter().map(|s| s.time()).collect(); let expected = contacts(); - let actual = visibility(×, 0.0, &gs, &sc, &NoOpFrameTransformationProvider); + let actual = visibility(×, &gs, &mask, &sc, &NoOpFrameTransformationProvider); assert_eq!(actual.len(), expected.len()); for (actual, expected) in zip(actual, expected) { assert_close!(actual.start(), expected.start(), 0.0, 1e-4); diff --git a/crates/lox-orbits/src/python.rs b/crates/lox-orbits/src/python.rs index 7e4c5a3c..75f543a4 100644 --- a/crates/lox-orbits/src/python.rs +++ b/crates/lox-orbits/src/python.rs @@ -32,6 +32,7 @@ use lox_time::transformations::TryToScale; use lox_time::{python::time::PyTime, ut1::DeltaUt1Tai, Time}; use python::PyBody; +use crate::analysis::{ElevationMask, ElevationMaskError}; use crate::elements::{Keplerian, ToKeplerian}; use crate::events::{Event, FindEventError, Window}; use crate::frames::{BodyFixed, CoordinateSystem, Icrf, ReferenceFrame, Topocentric, TryToFrame}; @@ -966,7 +967,7 @@ impl PySgp4 { pub fn visibility( times: &Bound<'_, PyList>, gs: PyGroundLocation, - min_elevation: f64, + mask: &Bound<'_, PyElevationMask>, sc: &Bound<'_, PyTrajectory>, provider: &Bound<'_, PyUt1Provider>, ) -> PyResult> { @@ -982,15 +983,43 @@ pub fn visibility( }; let times: Vec = times.extract()?; let provider = provider.get(); + let mask = &mask.borrow().0; let sc = sc.0.with_origin_and_frame(sc_origin, Icrf); Ok( - crate::analysis::visibility(×, min_elevation, &gs.0, &sc, provider) + crate::analysis::visibility(×, &gs.0, mask, &sc, provider) .into_iter() .map(PyWindow) .collect(), ) } +impl From for PyErr { + fn from(err: ElevationMaskError) -> Self { + PyValueError::new_err(err.to_string()) + } +} + +#[pyclass(name = "ElevationMask", module = "lox_space", frozen)] +pub struct PyElevationMask(pub ElevationMask); + +#[pymethods] +impl PyElevationMask { + #[new] + fn new( + azimuth: &Bound<'_, PyArray1>, + elevation: &Bound<'_, PyArray1>, + ) -> PyResult { + let azimuth = azimuth.to_vec()?; + let elevation = elevation.to_vec()?; + Ok(PyElevationMask(ElevationMask::new(azimuth, elevation)?)) + } + + #[classmethod] + fn fixed(_cls: &Bound<'_, PyType>, min_elevation: f64) -> Self { + PyElevationMask(ElevationMask::with_fixed_elevation(min_elevation)) + } +} + #[pyclass(name = "Topocentric", module = "lox_space", frozen)] pub struct PyTopocentric(pub Topocentric); diff --git a/tools/lox-gen/Cargo.lock b/tools/lox-gen/Cargo.lock index 71e4e7ba..079b1ec4 100644 --- a/tools/lox-gen/Cargo.lock +++ b/tools/lox-gen/Cargo.lock @@ -89,12 +89,6 @@ dependencies = [ "memchr", ] -[[package]] -name = "fast-float" -version = "0.2.0" -source = "registry+https://github.com/rust-lang/crates.io-index" -checksum = "95765f67b4b18863968b4a1bd5bb576f732b29a4a28c7cd84c09fa3e2875f33c" - [[package]] name = "fast_polynomial" version = "0.1.0" @@ -194,10 +188,9 @@ dependencies = [ [[package]] name = "lox-io" -version = "0.1.0-alpha.0" +version = "0.1.0-alpha.1" dependencies = [ "csv", - "fast-float", "lox-derive", "lox-math", "nom", @@ -211,7 +204,7 @@ dependencies = [ [[package]] name = "lox-math" -version = "0.1.0-alpha.1" +version = "0.1.0-alpha.2" dependencies = [ "fast_polynomial", "float_eq",