use super::vertical_profile::{ValidatedZenithDomain, VerticalEmissionProfileError};
use crate::error::{NsbError, Result};
use crate::units::ScaleFactors;
use siderust::atmosphere::{van_rhijn_factor, ScatteringFactor};
use siderust::qtty::{unit::Radian, Degrees, Kilometers, Quantity};
pub(crate) const VAN_RHIJN_IMPLEMENTATION_VERSION: &str = "siderust-0.11.0-mean-earth-radius";
pub(crate) const DEFAULT_VAN_RHIJN_EMISSION_HEIGHT_KM: Kilometers = Kilometers::new(88.0);
#[derive(Debug, Clone, Copy, PartialEq)]
pub struct VanRhijnConfig {
emission_height_km: Kilometers,
}
impl VanRhijnConfig {
pub fn new(
emission_height_km: Kilometers,
) -> std::result::Result<Self, VerticalEmissionProfileError> {
if !emission_height_km.is_finite() || emission_height_km <= Kilometers::new(0.0) {
return Err(VerticalEmissionProfileError::Invalid(
"Van Rhijn emission height must be finite and positive".into(),
));
}
Ok(Self { emission_height_km })
}
pub(crate) const fn from_continuum_height(emission_height_km: Kilometers) -> Self {
Self { emission_height_km }
}
pub const fn emission_height_km(self) -> Kilometers {
self.emission_height_km
}
pub const fn validated_zenith_domain(self) -> ValidatedZenithDomain {
ValidatedZenithDomain {
min: Degrees::new(0.0),
max: Degrees::new(90.0),
}
}
pub(crate) fn geometry_factor(self, zenith: Degrees) -> Result<ScaleFactors> {
if !zenith.is_finite() || zenith < Degrees::new(0.0) || zenith > Degrees::new(90.0) {
return Err(NsbError::OutOfRange(format!(
"Van Rhijn zenith angle must be in [0, 90] deg, got {}",
zenith.value()
)));
}
let factor = van_rhijn_factor(zenith.to::<Radian>(), self.emission_height_km());
if !factor.is_finite() || factor <= Quantity::<ScatteringFactor>::new(0.0) {
return Err(NsbError::Unsupported(
"Van Rhijn configuration produced a non-finite geometry factor".into(),
));
}
Ok(factor.to())
}
}
impl Default for VanRhijnConfig {
fn default() -> Self {
Self {
emission_height_km: DEFAULT_VAN_RHIJN_EMISSION_HEIGHT_KM,
}
}
}
#[cfg(test)]
mod tests {
use super::super::AirglowGeometryModel;
use super::*;
use siderust::coordinates::centers::Geodetic;
use siderust::coordinates::frames::ECEF;
use siderust::qtty::Meters;
fn observer(height_m: f64) -> Geodetic<ECEF> {
Geodetic::new_raw(
Degrees::new(12.345),
Degrees::new(-43.21),
Meters::new(height_m),
)
}
#[test]
fn van_rhijn_default_is_the_palace_representative_configuration() {
let default = AirglowGeometryModel::default();
let explicit =
AirglowGeometryModel::VanRhijn(VanRhijnConfig::new(Kilometers::new(88.0)).unwrap());
for zenith in [0.0, 30.0, 60.0, 80.0, 90.0] {
let default_factor = default
.geometry_factor(observer(2_635.0), Degrees::new(zenith))
.unwrap();
let explicit_factor = explicit
.geometry_factor(observer(2_635.0), Degrees::new(zenith))
.unwrap();
assert_eq!(
default_factor.value().to_bits(),
explicit_factor.value().to_bits()
);
}
}
#[test]
fn invalid_emission_height_is_rejected() {
assert!(VanRhijnConfig::new(Kilometers::new(0.0)).is_err());
assert!(VanRhijnConfig::new(Kilometers::new(-1.0)).is_err());
assert!(VanRhijnConfig::new(Kilometers::new(f64::NAN)).is_err());
}
#[test]
fn zenith_outside_domain_is_rejected() {
let config = VanRhijnConfig::default();
assert!(config.geometry_factor(Degrees::new(-0.1)).is_err());
assert!(config.geometry_factor(Degrees::new(90.1)).is_err());
}
}