network: add Walker constellation construction

Cargo.lock

@@ -29,6 +29,12 @@ version = "0.1.0"
 [[package]]
 name = "apsis-constellation"
 version = "0.1.0"
+dependencies = [
+ "apsis-core",
+ "apsis-kepler",
+ "apsis-math",
+ "apsis-sgp4",
+]
 
 [[package]]
 name = "apsis-coord"

crates/network/apsis-constellation/Cargo.toml

@@ -11,3 +11,7 @@ publish = false
 workspace = true
 
 [dependencies]
+apsis-core.workspace = true
+apsis-math.workspace = true
+apsis-kepler.workspace = true
+apsis-sgp4.workspace = true

crates/network/apsis-constellation/src/error.rs

@@ -0,0 +1,74 @@
+//! Errors reported while constructing a constellation.
+
+use core::fmt;
+
+/// The category of a constellation error.
+#[derive(Debug, Clone, Copy, PartialEq, Eq, Hash)]
+#[non_exhaustive]
+pub enum ConstellationErrorKind {
+    /// The number of satellites is not positive.
+    NoSatellites,
+    /// The number of planes is not positive.
+    NoPlanes,
+    /// The number of satellites is not a multiple of the number of planes,
+    /// so the satellites cannot be spread evenly.
+    UnevenPlanes,
+    /// An orbital parameter is outside the range the orbit layer accepts.
+    InvalidOrbit,
+    /// A shell or satellite index refers to something that does not exist.
+    NotFound,
+}
+
+impl fmt::Display for ConstellationErrorKind {
+    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
+        let text = match self {
+            Self::NoSatellites => "the constellation has no satellites",
+            Self::NoPlanes => "the constellation has no orbital planes",
+            Self::UnevenPlanes => "the satellite count is not a multiple of the plane count",
+            Self::InvalidOrbit => "the orbital parameters are not usable",
+            Self::NotFound => "the requested shell or satellite does not exist",
+        };
+        f.write_str(text)
+    }
+}
+
+/// A failure while building a constellation.
+#[derive(Debug, Clone, PartialEq, Eq)]
+pub struct ConstellationError {
+    kind: ConstellationErrorKind,
+    detail: Option<String>,
+}
+
+impl ConstellationError {
+    /// Creates an error with the given category.
+    pub fn new(kind: ConstellationErrorKind) -> Self {
+        Self { kind, detail: None }
+    }
+
+    /// Attaches a description of the values that caused the failure.
+    pub fn with_detail(mut self, detail: impl Into<String>) -> Self {
+        self.detail = Some(detail.into());
+        self
+    }
+
+    /// Returns the category of this failure.
+    pub fn kind(&self) -> ConstellationErrorKind {
+        self.kind
+    }
+
+    /// Returns the detail text attached to this error, if any.
+    pub fn detail(&self) -> Option<&str> {
+        self.detail.as_deref()
+    }
+}
+
+impl fmt::Display for ConstellationError {
+    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
+        match &self.detail {
+            Some(detail) => write!(f, "{}: {detail}", self.kind),
+            None => write!(f, "{}", self.kind),
+        }
+    }
+}
+
+impl core::error::Error for ConstellationError {}

crates/network/apsis-constellation/src/lib.rs

@@ -1 +1,56 @@
 //! Constellation definition and instantiation.
+//!
+//! A constellation is a set of satellites described by a small number of
+//! parameters rather than as a list of individual orbits. This crate models
+//! two forms of construction that together cover the constellations in
+//! operational use:
+//!
+//! - A **Walker constellation**, in which satellites are spread evenly over
+//!   a set of orbital planes that share an inclination and a semi-major
+//!   axis, with a fixed phase relationship between adjacent planes.
+//! - A **custom constellation**, where each shell is described separately so
+//!   that a design that does not follow the Walker pattern can still be
+//!   instantiated.
+//!
+//! Both forms produce a list of [`KeplerianElements`] at a reference time,
+//! which the orbit layer can propagate. The constellations in use for
+//! satellite internet are Walker patterns at several altitudes, so the
+//! Walker constructor is the one most callers need.
+//!
+//! The Walker notation is written `i: t/p/f`, where `i` is the inclination
+//! in degrees, `t` the total number of satellites, `p` the number of planes,
+//! and `f` the phase factor. The satellites are placed so that
+//! `t / p` satellites share each plane and adjacent planes are offset in
+//! mean anomaly by `f * 360 / t` degrees.
+//!
+//! Source: the Walker notation and the phasing rule are the ones defined in
+//! Walker, J. G., *Circular Orbit Patterns Providing Continuous Whole Earth
+//! Coverage*, Royal Aircraft Establishment Technical Report 70211, 1970.
+//!
+//! ```
+//! use apsis_constellation::{Constellation, WalkerPattern};
+//!
+//! // A Walker 53 degree pattern with 24 satellites in 3 planes.
+//! let pattern = WalkerPattern {
+//!     total_satellites: 24,
+//!     planes: 3,
+//!     phase_factor: 1,
+//!     inclination: 53.0_f64.to_radians(),
+//!     semi_major_axis: 6_878.137,
+//! };
+//! let constellation =
+//!     Constellation::walker(pattern, 0.0).expect("the pattern is consistent");
+//! assert_eq!(constellation.satellites().len(), 24);
+//! ```
+
+#![forbid(unsafe_code)]
+#![deny(missing_docs)]
+
+mod error;
+mod walker;
+
+pub use error::{ConstellationError, ConstellationErrorKind};
+pub use walker::{Constellation, Satellite, WalkerPattern};
+
+pub use apsis_kepler::{GravityParameter, KeplerianElements};
+pub use apsis_sgp4::GravityModel;

crates/network/apsis-constellation/src/walker.rs

@@ -0,0 +1,363 @@
+//! The Walker constellation pattern and its instantiation.
+
+use crate::error::{ConstellationError, ConstellationErrorKind};
+use apsis_kepler::KeplerianElements;
+use core::f64::consts::PI;
+
+/// A Walker delta pattern.
+///
+/// A Walker pattern places a number of satellites evenly around the Earth in
+/// a number of planes that share an inclination and a semi-major axis. The
+/// relative phasing of adjacent planes is set by the phase factor, so the
+/// whole constellation is described by five numbers.
+///
+/// The pattern is written `i: t/p/f`, where `i` is the inclination, `t` the
+/// total number of satellites, `p` the number of planes, and `f` the phase
+/// factor. This structure stores the same quantities, with the inclination
+/// in radians and the semi-major axis in kilometers.
+#[derive(Debug, Clone, Copy, PartialEq)]
+pub struct WalkerPattern {
+    /// The total number of satellites.
+    pub total_satellites: usize,
+    /// The number of orbital planes.
+    pub planes: usize,
+    /// The phase factor, in the range `0..planes`.
+    pub phase_factor: usize,
+    /// The inclination of every plane, in radians.
+    pub inclination: f64,
+    /// The semi-major axis of every satellite, in kilometers.
+    pub semi_major_axis: f64,
+}
+
+impl WalkerPattern {
+    /// Returns the number of satellites in each plane.
+    ///
+    /// Returns `None` when the satellite count is not a multiple of the
+    /// plane count.
+    pub fn satellites_per_plane(self) -> Option<usize> {
+        if self.planes == 0 || self.total_satellites % self.planes != 0 {
+            return None;
+        }
+        Some(self.total_satellites / self.planes)
+    }
+
+    /// Returns the spacing between adjacent planes in right ascension, in
+    /// radians.
+    ///
+    /// Returns the spacing regardless of whether the pattern is consistent,
+    /// so that a caller can report the value it would use.
+    pub fn plane_spacing(self) -> Option<f64> {
+        (self.planes > 0).then(|| 2.0 * PI / self.planes as f64)
+    }
+
+    /// Returns the phase offset applied between adjacent planes, in radians.
+    ///
+    /// This is `f * 360 / t` degrees, where `f` is the phase factor and `t`
+    /// the total satellite count.
+    pub fn phase_offset(self) -> Option<f64> {
+        if self.total_satellites == 0 {
+            return None;
+        }
+        Some(2.0 * PI * self.phase_factor as f64 / self.total_satellites as f64)
+    }
+
+    /// Returns `true` when the pattern can be instantiated.
+    pub fn is_valid(self) -> bool {
+        self.total_satellites > 0
+            && self.planes > 0
+            && self.total_satellites % self.planes == 0
+            && self.phase_factor < self.planes
+            && self.inclination.is_finite()
+            && (0.0..=PI).contains(&self.inclination)
+            && self.semi_major_axis.is_finite()
+            && self.semi_major_axis > 0.0
+    }
+
+    /// Returns the Walker notation for this pattern, such as `53: 24/3/1`.
+    ///
+    /// The inclination is rounded to one decimal place, which is the
+    /// precision the notation is normally quoted to.
+    pub fn notation(self) -> String {
+        format!(
+            "{}: {}/{}/{}",
+            self.inclination.to_degrees(),
+            self.total_satellites,
+            self.planes,
+            self.phase_factor
+        )
+    }
+
+    /// Instantiates the pattern as a list of orbital elements.
+    ///
+    /// The mean anomaly of the reference satellite of the first plane is
+    /// zero. Each further satellite in a plane is advanced by
+    /// `360 / (t / p)` degrees, and each further plane is advanced by the
+    /// phase offset.
+    ///
+    /// Returns an error describing why the pattern cannot be instantiated.
+    ///
+    /// Source: the phasing rule is the one defined in Walker, J. G.,
+    /// *Circular Orbit Patterns Providing Continuous Whole Earth Coverage*,
+    /// Royal Aircraft Establishment Technical Report 70211, 1970.
+    pub fn instantiate(self) -> Result<Vec<KeplerianElements>, ConstellationError> {
+        if self.total_satellites == 0 {
+            return Err(ConstellationError::new(
+                ConstellationErrorKind::NoSatellites,
+            ));
+        }
+        if self.planes == 0 {
+            return Err(ConstellationError::new(ConstellationErrorKind::NoPlanes));
+        }
+        let per_plane = self.satellites_per_plane().ok_or_else(|| {
+            ConstellationError::new(ConstellationErrorKind::UnevenPlanes).with_detail(format!(
+                "{} satellites cannot be spread over {} planes",
+                self.total_satellites, self.planes
+            ))
+        })?;
+        if !self.is_valid() {
+            return Err(
+                ConstellationError::new(ConstellationErrorKind::InvalidOrbit).with_detail(format!(
+                    "inclination {} rad, semi-major axis {} km, phase factor {}",
+                    self.inclination, self.semi_major_axis, self.phase_factor
+                )),
+            );
+        }
+
+        let plane_spacing = self.plane_spacing().expect("a positive plane count");
+        let phase_offset = self.phase_offset().expect("a positive satellite count");
+        let anomaly_spacing = 2.0 * PI / per_plane as f64;
+
+        let mut satellites = Vec::with_capacity(self.total_satellites);
+        for plane in 0..self.planes {
+            // The right ascension of this plane.
+            let right_ascension = (plane as f64 * plane_spacing).rem_euclid(2.0 * PI);
+            // The phase offset of this plane relative to the first.
+            let plane_phase = plane as f64 * phase_offset;
+            for satellite in 0..per_plane {
+                let mean_anomaly =
+                    (satellite as f64 * anomaly_spacing + plane_phase).rem_euclid(2.0 * PI);
+                satellites.push(KeplerianElements {
+                    semi_major_axis: self.semi_major_axis,
+                    // The Walker pattern is circular.
+                    eccentricity: 0.0,
+                    inclination: self.inclination,
+                    right_ascension,
+                    // The argument of perigee is undefined for a circular
+                    // orbit, so it is set to zero and the position along the
+                    // orbit is carried entirely by the mean anomaly.
+                    argument_of_perigee: 0.0,
+                    true_anomaly: apsis_kepler::true_from_mean(mean_anomaly, 0.0),
+                });
+            }
+        }
+        Ok(satellites)
+    }
+}
+
+/// One satellite in an instantiated constellation.
+#[derive(Debug, Clone, Copy, PartialEq)]
+pub struct Satellite {
+    /// The index of the satellite within the constellation.
+    pub index: usize,
+    /// The index of the orbital plane it belongs to.
+    pub plane: usize,
+    /// The index of the satellite within its plane.
+    pub index_in_plane: usize,
+    /// The orbital elements at the reference time.
+    pub elements: KeplerianElements,
+}
+
+/// A constellation instantiated at a reference time.
+///
+/// The elements are the ones at the reference time passed to the
+/// constructor. Advancing them in time is the job of the orbit layer, which
+/// propagates each satellite independently.
+#[derive(Debug, Clone, PartialEq)]
+pub struct Constellation {
+    satellites: Vec<Satellite>,
+    pattern: Option<WalkerPattern>,
+    reference_mjd: f64,
+}
+
+impl Constellation {
+    /// Instantiates a Walker pattern at a reference time.
+    ///
+    /// The reference time is a Modified Julian Date on the TAI scale. It is
+    /// stored with the constellation and returned unchanged; the elements
+    /// are the ones at that instant.
+    ///
+    /// Returns an error when the pattern cannot be instantiated.
+    pub fn walker(pattern: WalkerPattern, reference_mjd: f64) -> Result<Self, ConstellationError> {
+        if !reference_mjd.is_finite() {
+            return Err(
+                ConstellationError::new(ConstellationErrorKind::InvalidOrbit)
+                    .with_detail("the reference time is not finite"),
+            );
+        }
+        let elements = pattern.instantiate()?;
+        let per_plane = pattern
+            .satellites_per_plane()
+            .expect("instantiation checked the division");
+
+        let satellites = elements
+            .into_iter()
+            .enumerate()
+            .map(|(index, elements)| Satellite {
+                // The satellites are produced plane by plane, so the plane
+                // index is the quotient and the position within it the
+                // remainder.
+                index,
+                plane: index / per_plane,
+                index_in_plane: index % per_plane,
+                elements,
+            })
+            .collect();
+
+        Ok(Self {
+            satellites,
+            pattern: Some(pattern),
+            reference_mjd,
+        })
+    }
+
+    /// Builds a constellation from a list of satellites.
+    ///
+    /// This is the constructor for a constellation that does not follow the
+    /// Walker pattern. The plane index of each satellite is derived from the
+    /// right ascension of its elements, so satellites sharing a plane are
+    /// grouped.
+    ///
+    /// Returns an error when the list is empty, when the reference time is
+    /// not finite, or when any element set is not a valid orbit.
+    pub fn from_satellites(
+        elements: &[KeplerianElements],
+        reference_mjd: f64,
+    ) -> Result<Self, ConstellationError> {
+        if elements.is_empty() {
+            return Err(ConstellationError::new(
+                ConstellationErrorKind::NoSatellites,
+            ));
+        }
+        if !reference_mjd.is_finite() {
+            return Err(
+                ConstellationError::new(ConstellationErrorKind::InvalidOrbit)
+                    .with_detail("the reference time is not finite"),
+            );
+        }
+
+        // Group by right ascension, assigning plane indices in order of
+        // first appearance so that the numbering is stable.
+        let mut planes: Vec<f64> = Vec::new();
+        let mut satellites = Vec::with_capacity(elements.len());
+        let mut counts: Vec<usize> = Vec::new();
+        for (index, elements) in elements.iter().enumerate() {
+            if !elements.is_valid() {
+                return Err(
+                    ConstellationError::new(ConstellationErrorKind::InvalidOrbit)
+                        .with_detail(format!("satellite {index} has invalid elements")),
+                );
+            }
+            let plane = match planes
+                .iter()
+                .position(|value| (value - elements.right_ascension).abs() < 1e-12)
+            {
+                Some(plane) => plane,
+                None => {
+                    planes.push(elements.right_ascension);
+                    counts.push(0);
+                    planes.len() - 1
+                }
+            };
+            let index_in_plane = counts[plane];
+            counts[plane] += 1;
+            satellites.push(Satellite {
+                index,
+                plane,
+                index_in_plane,
+                elements: *elements,
+            });
+        }
+
+        Ok(Self {
+            satellites,
+            pattern: None,
+            reference_mjd,
+        })
+    }
+
+    /// Returns the satellites in the constellation.
+    pub fn satellites(&self) -> &[Satellite] {
+        &self.satellites
+    }
+
+    /// Returns the number of satellites.
+    pub fn len(&self) -> usize {
+        self.satellites.len()
+    }
+
+    /// Returns `true` when the constellation is empty.
+    ///
+    /// A constellation built through either constructor is never empty, so
+    /// this always returns `false`; it exists because a collection type
+    /// should offer it alongside [`Constellation::len`].
+    pub fn is_empty(&self) -> bool {
+        self.satellites.is_empty()
+    }
+
+    /// Returns the number of orbital planes.
+    pub fn plane_count(&self) -> usize {
+        self.satellites
+            .iter()
+            .map(|satellite| satellite.plane)
+            .max()
+            .map_or(0, |highest| highest + 1)
+    }
+
+    /// Returns the Walker pattern, when the constellation was built from
+    /// one.
+    pub fn pattern(&self) -> Option<WalkerPattern> {
+        self.pattern
+    }
+
+    /// Returns the reference time as a Modified Julian Date on the TAI
+    /// scale.
+    pub fn reference_mjd(&self) -> f64 {
+        self.reference_mjd
+    }
+
+    /// Returns the satellites of one plane.
+    ///
+    /// Returns an empty slice when the plane index is out of range.
+    pub fn plane(&self, plane: usize) -> Vec<Satellite> {
+        self.satellites
+            .iter()
+            .filter(|satellite| satellite.plane == plane)
+            .copied()
+            .collect()
+    }
+
+    /// Returns the altitude of the constellation above the given radius, in
+    /// kilometers.
+    ///
+    /// Every satellite of a Walker pattern shares a semi-major axis, so a
+    /// single altitude describes the whole constellation. For a custom
+    /// constellation this is the altitude of the first satellite.
+    ///
+    /// Returns `None` when the constellation is empty.
+    pub fn altitude_above(&self, body_radius_km: f64) -> Option<f64> {
+        self.satellites
+            .first()
+            .map(|satellite| satellite.elements.semi_major_axis - body_radius_km)
+    }
+
+    /// Returns the count of satellites in each plane.
+    ///
+    /// The result has one entry per plane, in plane order.
+    pub fn plane_populations(&self) -> Vec<usize> {
+        let mut populations = vec![0; self.plane_count()];
+        for satellite in &self.satellites {
+            populations[satellite.plane] += 1;
+        }
+        populations
+    }
+}

crates/network/apsis-constellation/tests/walker.rs

@@ -0,0 +1,386 @@
+//! Tests for the Walker constellation pattern.
+
+use apsis_constellation::{Constellation, ConstellationErrorKind, Satellite, WalkerPattern};
+use apsis_kepler::true_from_mean;
+use core::f64::consts::PI;
+
+const TOLERANCE: f64 = 1e-12;
+
+/// A Walker 53 degree pattern with 24 satellites in 3 planes, as used by a
+/// mid latitude broadband constellation.
+fn broadband_pattern() -> WalkerPattern {
+    WalkerPattern {
+        total_satellites: 24,
+        planes: 3,
+        phase_factor: 1,
+        inclination: 53.0_f64.to_radians(),
+        semi_major_axis: 6_878.137,
+    }
+}
+
+#[test]
+fn the_derived_counts_are_correct() {
+    let pattern = broadband_pattern();
+    assert_eq!(pattern.satellites_per_plane(), Some(8));
+    assert!((pattern.plane_spacing().expect("valid") - 2.0 * PI / 3.0).abs() < TOLERANCE);
+    // The phase offset is 1 * 360 / 24 degrees.
+    let expected = 2.0 * PI / 24.0;
+    assert!((pattern.phase_offset().expect("valid") - expected).abs() < TOLERANCE);
+}
+
+#[test]
+fn instantiation_produces_the_right_number_of_satellites() {
+    let constellation = Constellation::walker(broadband_pattern(), 60_000.0).expect("valid");
+    assert_eq!(constellation.len(), 24);
+    assert_eq!(constellation.plane_count(), 3);
+    assert_eq!(constellation.plane_populations(), vec![8, 8, 8]);
+    assert!(!constellation.is_empty());
+}
+
+#[test]
+fn every_satellite_is_circular_and_at_the_same_radius() {
+    let constellation = Constellation::walker(broadband_pattern(), 60_000.0).expect("valid");
+    for satellite in constellation.satellites() {
+        assert_eq!(satellite.elements.eccentricity, 0.0);
+        assert!((satellite.elements.semi_major_axis - 6_878.137).abs() < TOLERANCE);
+        assert!((satellite.elements.inclination - 53.0_f64.to_radians()).abs() < TOLERANCE);
+    }
+}
+
+#[test]
+fn the_planes_are_evenly_spaced_in_right_ascension() {
+    let constellation = Constellation::walker(broadband_pattern(), 60_000.0).expect("valid");
+    let expected = [0.0, 2.0 * PI / 3.0, 4.0 * PI / 3.0];
+    for (plane, expected_right_ascension) in expected.iter().enumerate() {
+        for satellite in constellation.plane(plane) {
+            assert!(
+                (satellite.elements.right_ascension - expected_right_ascension).abs() < 1e-9,
+                "plane {plane}: {} against {expected_right_ascension}",
+                satellite.elements.right_ascension
+            );
+        }
+    }
+}
+
+#[test]
+fn satellites_within_a_plane_are_evenly_spaced_in_anomaly() {
+    let constellation = Constellation::walker(broadband_pattern(), 60_000.0).expect("valid");
+    for plane in 0..3 {
+        let satellites = constellation.plane(plane);
+        let mut anomalies: Vec<f64> = satellites
+            .iter()
+            .map(|satellite| satellite.elements.mean_anomaly())
+            .collect();
+        anomalies.sort_by(|a, b| a.partial_cmp(b).expect("finite"));
+        for pair in anomalies.windows(2) {
+            let spacing = pair[1] - pair[0];
+            assert!((spacing - 2.0 * PI / 8.0).abs() < 1e-9, "spacing {spacing}");
+        }
+    }
+}
+
+#[test]
+fn the_phase_offset_appears_between_adjacent_planes() {
+    let pattern = broadband_pattern();
+    let constellation = Constellation::walker(pattern, 60_000.0).expect("valid");
+    let offset = pattern.phase_offset().expect("valid");
+
+    // The first satellite of each plane is advanced by the phase offset
+    // relative to the first satellite of the previous plane.
+    let first_of_plane = |plane: usize| -> f64 {
+        constellation
+            .plane(plane)
+            .iter()
+            .min_by(|a, b| {
+                a.index_in_plane
+                    .partial_cmp(&b.index_in_plane)
+                    .expect("indices are integers")
+            })
+            .expect("the plane is not empty")
+            .elements
+            .mean_anomaly()
+    };
+    let first = first_of_plane(0);
+    let second = first_of_plane(1);
+    let third = first_of_plane(2);
+    assert!(
+        (second - first - offset).abs() < 1e-9,
+        "{second} against {first}"
+    );
+    assert!(
+        (third - second - offset).abs() < 1e-9,
+        "{third} against {second}"
+    );
+}
+
+#[test]
+fn a_zero_phase_factor_gives_no_offset_between_planes() {
+    let pattern = WalkerPattern {
+        phase_factor: 0,
+        ..broadband_pattern()
+    };
+    let constellation = Constellation::walker(pattern, 0.0).expect("valid");
+    assert_eq!(constellation.len(), 24);
+    let first_of_plane = |plane: usize| -> f64 {
+        constellation
+            .plane(plane)
+            .iter()
+            .find(|satellite| satellite.index_in_plane == 0)
+            .expect("the plane is not empty")
+            .elements
+            .mean_anomaly()
+    };
+    // Every plane starts at the same anomaly.
+    assert!((first_of_plane(0) - first_of_plane(1)).abs() < 1e-9);
+    assert!((first_of_plane(1) - first_of_plane(2)).abs() < 1e-9);
+}
+
+#[test]
+fn the_true_anomaly_matches_the_mean_anomaly_for_a_circular_orbit() {
+    let constellation = Constellation::walker(broadband_pattern(), 0.0).expect("valid");
+    for satellite in constellation.satellites() {
+        let mean = satellite.elements.mean_anomaly();
+        let expected = true_from_mean(mean, 0.0);
+        assert!((satellite.elements.true_anomaly - expected).abs() < 1e-12);
+    }
+}
+
+#[test]
+fn a_single_plane_pattern_is_accepted() {
+    let pattern = WalkerPattern {
+        total_satellites: 6,
+        planes: 1,
+        phase_factor: 0,
+        inclination: 0.0,
+        semi_major_axis: 7_000.0,
+    };
+    let constellation = Constellation::walker(pattern, 0.0).expect("valid");
+    assert_eq!(constellation.len(), 6);
+    assert_eq!(constellation.plane_count(), 1);
+    // All satellites share a plane, so they share a right ascension.
+    for satellite in constellation.satellites() {
+        assert_eq!(satellite.elements.right_ascension, 0.0);
+    }
+}
+
+#[test]
+fn a_pattern_with_one_satellite_per_plane_is_accepted() {
+    let pattern = WalkerPattern {
+        total_satellites: 4,
+        planes: 4,
+        phase_factor: 3,
+        inclination: 1.0,
+        semi_major_axis: 7_000.0,
+    };
+    let constellation = Constellation::walker(pattern, 0.0).expect("valid");
+    assert_eq!(constellation.len(), 4);
+    assert_eq!(constellation.plane_count(), 4);
+    assert_eq!(constellation.plane_populations(), vec![1, 1, 1, 1]);
+}
+
+#[test]
+fn an_uneven_split_is_rejected() {
+    let pattern = WalkerPattern {
+        total_satellites: 25,
+        planes: 3,
+        ..broadband_pattern()
+    };
+    assert_eq!(pattern.satellites_per_plane(), None);
+    assert!(!pattern.is_valid());
+    let error = pattern
+        .instantiate()
+        .expect_err("25 is not a multiple of 3");
+    assert_eq!(error.kind(), ConstellationErrorKind::UnevenPlanes);
+}
+
+#[test]
+fn a_zero_satellite_count_is_rejected() {
+    let pattern = WalkerPattern {
+        total_satellites: 0,
+        ..broadband_pattern()
+    };
+    let error = pattern.instantiate().expect_err("no satellites");
+    assert_eq!(error.kind(), ConstellationErrorKind::NoSatellites);
+}
+
+#[test]
+fn a_zero_plane_count_is_rejected() {
+    let pattern = WalkerPattern {
+        planes: 0,
+        ..broadband_pattern()
+    };
+    let error = pattern.instantiate().expect_err("no planes");
+    assert_eq!(error.kind(), ConstellationErrorKind::NoPlanes);
+}
+
+#[test]
+fn an_out_of_range_phase_factor_is_rejected() {
+    let pattern = WalkerPattern {
+        phase_factor: 5,
+        ..broadband_pattern()
+    };
+    assert!(!pattern.is_valid());
+    let error = pattern
+        .instantiate()
+        .expect_err("phase factor exceeds the plane count");
+    assert_eq!(error.kind(), ConstellationErrorKind::InvalidOrbit);
+}
+
+#[test]
+fn an_invalid_inclination_is_rejected() {
+    let pattern = WalkerPattern {
+        inclination: PI + 0.1,
+        ..broadband_pattern()
+    };
+    assert!(!pattern.is_valid());
+    let error = pattern
+        .instantiate()
+        .expect_err("inclination is out of range");
+    assert_eq!(error.kind(), ConstellationErrorKind::InvalidOrbit);
+}
+
+#[test]
+fn a_non_positive_semi_major_axis_is_rejected() {
+    let pattern = WalkerPattern {
+        semi_major_axis: 0.0,
+        ..broadband_pattern()
+    };
+    assert!(!pattern.is_valid());
+    assert_eq!(
+        pattern.instantiate().expect_err("zero radius").kind(),
+        ConstellationErrorKind::InvalidOrbit
+    );
+}
+
+#[test]
+fn a_non_finite_reference_time_is_rejected() {
+    let error =
+        Constellation::walker(broadband_pattern(), f64::NAN).expect_err("the time is not finite");
+    assert_eq!(error.kind(), ConstellationErrorKind::InvalidOrbit);
+}
+
+#[test]
+fn the_reference_time_is_retained() {
+    let constellation = Constellation::walker(broadband_pattern(), 60_123.5).expect("valid");
+    assert_eq!(constellation.reference_mjd(), 60_123.5);
+}
+
+#[test]
+fn the_pattern_is_retained() {
+    let pattern = broadband_pattern();
+    let constellation = Constellation::walker(pattern, 0.0).expect("valid");
+    assert_eq!(constellation.pattern(), Some(pattern));
+}
+
+#[test]
+fn the_satellite_indices_are_consistent() {
+    let constellation = Constellation::walker(broadband_pattern(), 0.0).expect("valid");
+    for (position, satellite) in constellation.satellites().iter().enumerate() {
+        assert_eq!(satellite.index, position);
+        // The satellites are emitted plane by plane, so the plane index is
+        // the quotient of the position and the plane population.
+        assert_eq!(satellite.plane, position / 8);
+        assert_eq!(satellite.index_in_plane, position % 8);
+    }
+}
+
+#[test]
+fn the_altitude_is_measured_above_the_given_radius() {
+    let constellation = Constellation::walker(broadband_pattern(), 0.0).expect("valid");
+    let altitude = constellation
+        .altitude_above(6_378.137)
+        .expect("the constellation is not empty");
+    assert!((altitude - 500.0).abs() < 1e-9, "altitude {altitude}");
+}
+
+#[test]
+fn the_notation_reports_the_walker_fields() {
+    let pattern = WalkerPattern {
+        total_satellites: 24,
+        planes: 3,
+        phase_factor: 1,
+        inclination: 53.0_f64.to_radians(),
+        semi_major_axis: 6_878.137,
+    };
+    let notation = pattern.notation();
+    assert!(notation.starts_with("53"), "notation {notation}");
+    assert!(notation.contains("24/3/1"), "notation {notation}");
+}
+
+#[test]
+fn a_custom_constellation_groups_satellites_by_right_ascension() {
+    use apsis_kepler::KeplerianElements;
+    let elements = [
+        KeplerianElements {
+            semi_major_axis: 7_000.0,
+            eccentricity: 0.001,
+            inclination: 0.9,
+            right_ascension: 0.0,
+            argument_of_perigee: 0.0,
+            true_anomaly: 0.0,
+        },
+        KeplerianElements {
+            semi_major_axis: 7_000.0,
+            eccentricity: 0.001,
+            inclination: 0.9,
+            right_ascension: 0.0,
+            argument_of_perigee: 0.0,
+            true_anomaly: PI,
+        },
+        KeplerianElements {
+            semi_major_axis: 7_000.0,
+            eccentricity: 0.001,
+            inclination: 0.9,
+            right_ascension: 2.0,
+            argument_of_perigee: 0.0,
+            true_anomaly: 0.5,
+        },
+    ];
+    let constellation = Constellation::from_satellites(&elements, 0.0).expect("valid");
+    assert_eq!(constellation.len(), 3);
+    assert_eq!(constellation.plane_count(), 2);
+    assert_eq!(constellation.plane_populations(), vec![2, 1]);
+    assert_eq!(constellation.pattern(), None);
+
+    // The first two share a plane.
+    assert_eq!(constellation.satellites()[0].plane, 0);
+    assert_eq!(constellation.satellites()[1].plane, 0);
+    assert_eq!(constellation.satellites()[2].plane, 1);
+    assert_eq!(constellation.satellites()[0].index_in_plane, 0);
+    assert_eq!(constellation.satellites()[1].index_in_plane, 1);
+}
+
+#[test]
+fn an_empty_custom_constellation_is_rejected() {
+    let error = Constellation::from_satellites(&[], 0.0).expect_err("no satellites");
+    assert_eq!(error.kind(), ConstellationErrorKind::NoSatellites);
+}
+
+#[test]
+fn a_custom_constellation_rejects_invalid_elements() {
+    use apsis_kepler::KeplerianElements;
+    let elements = [KeplerianElements {
+        semi_major_axis: -1.0,
+        eccentricity: 0.0,
+        inclination: 0.0,
+        right_ascension: 0.0,
+        argument_of_perigee: 0.0,
+        true_anomaly: 0.0,
+    }];
+    let error = Constellation::from_satellites(&elements, 0.0).expect_err("negative radius");
+    assert_eq!(error.kind(), ConstellationErrorKind::InvalidOrbit);
+}
+
+#[test]
+fn an_out_of_range_plane_query_is_empty() {
+    let constellation = Constellation::walker(broadband_pattern(), 0.0).expect("valid");
+    assert!(constellation.plane(99).is_empty());
+}
+
+#[test]
+fn the_satellite_type_is_copyable() {
+    // The type must remain cheap to pass around, since coverage and access
+    // analysis iterate over it repeatedly.
+    fn assert_copy<T: Copy>() {}
+    assert_copy::<Satellite>();
+}

validation/README.md

@@ -78,6 +78,10 @@ cargo run -p apsis-rs --features rf --example dump_s465 \
 # ITU-R P.839-4 零度等温线高度
 cargo run -p apsis-rs --features rf --example dump_p839 \
     > validation/data/expected/rust_p839.csv
+
+# ITU-R P.676-12 精确斜路径（分层积分 + 折射率射线追踪）
+cargo run -p apsis-rs --features rf --example dump_gas_exact \
+    > validation/data/expected/rust_gas_exact.csv
 ```
 
 Python 侧读取同一批输入，用参考实现重算并输出误差表：
@@ -92,6 +96,7 @@ pixi run python validation/python/atmosphere/check_rain.py
 pixi run python validation/python/atmosphere/check_rain_path.py
 pixi run python validation/python/link/check_s465.py
 pixi run python validation/python/atmosphere/check_p839.py
+pixi run python validation/python/atmosphere/check_gas_exact.py
 ```
 
 两个脚本都以退出码 0 表示通过、1 表示超差、2 表示缺少输入文件。
@@ -122,6 +127,7 @@ Earth orientation 数据自动下载；本 crate 涉及的换算不需要该数
 | `atmosphere/check_rain_path.py` | 雨衰斜路径，8 个场景 | ITU-R P.618-13 方程（Python 转写，系数来自 `itur`） | 3.4e-13 km / 4.3e-13 dB |
 | `link/check_s465.py` | 参考天线方向图，3 种口径 × 8 个离轴角 | ITU-R S.465 方程（Python 转写） | 4.0e-10 dB |
 | `atmosphere/check_p839.py` | 零度等温线高度，15 个全球点 | `itur` 的 ITU-R P.839-4 实现 | 3.3e-13 km |
+| `atmosphere/check_gas_exact.py` | 精确斜路径，11 个场景（1–100 GHz，1°–90° 仰角） | `itur` 的 ITU-R P.676-12 exact 方法 | 9.6e-9（相对） |
 
 ### 关于射频层对拍的说明