marss 0.0.2

Mars celestial simulation crate for the MilkyWay SolarSystem workspace
Documentation
use sciforge::hub::domain::meteorology::radiation::beer_lambert;

pub fn rayleigh_beta_rgb() -> [f64; 3] {
    [19.918e-6, 13.57e-6, 5.75e-6]
}

pub struct MarsAtmosphereParams {
    pub planet_radius: f64,
    pub atmosphere_height: f64,
    pub rayleigh_coefficients: [f64; 3],
    pub rayleigh_scale_height: f64,
    pub mie_coefficient: f64,
    pub mie_asymmetry: f64,
    pub dust_optical_depth: f64,
}

impl Default for MarsAtmosphereParams {
    fn default() -> Self {
        Self {
            planet_radius: crate::MARS_RADIUS,
            atmosphere_height: 200_000.0,
            rayleigh_coefficients: rayleigh_beta_rgb(),
            rayleigh_scale_height: crate::SCALE_HEIGHT_M,
            mie_coefficient: 4.0e-5,
            mie_asymmetry: 0.76,
            dust_optical_depth: crate::MEAN_DUST_OPTICAL_DEPTH,
        }
    }
}

impl MarsAtmosphereParams {
    pub fn rayleigh_density(&self, altitude: f64) -> f64 {
        (-altitude / self.rayleigh_scale_height).exp()
    }

    pub fn dust_density(&self, altitude: f64) -> f64 {
        (-altitude / 12_000.0).exp()
    }

    pub fn sky_color(&self, sun_elevation_deg: f64) -> [f64; 3] {
        let zenith_angle = (90.0 - sun_elevation_deg).to_radians().max(0.0);
        let airmass = if zenith_angle < 1.5 {
            1.0 / zenith_angle.cos()
        } else {
            40.0
        };

        let beta = self.rayleigh_coefficients;
        let tau_dust = self.dust_optical_depth;

        let r = (-beta[0] * airmass * self.rayleigh_scale_height).exp()
            * (-tau_dust * airmass * 0.7).exp();
        let g = (-beta[1] * airmass * self.rayleigh_scale_height).exp()
            * (-tau_dust * airmass * 0.5).exp();
        let b = (-beta[2] * airmass * self.rayleigh_scale_height).exp()
            * (-tau_dust * airmass * 0.3).exp();

        let dust_scatter = 0.3 * (1.0 - (-tau_dust * airmass).exp());
        [
            (r + dust_scatter * 0.9).min(1.0),
            (g + dust_scatter * 0.55).min(1.0),
            (b + dust_scatter * 0.25).min(1.0),
        ]
    }

    pub fn direct_transmission(&self, sun_elevation_deg: f64) -> f64 {
        let zenith = (90.0 - sun_elevation_deg).to_radians().max(0.01);
        let airmass = 1.0 / zenith.cos().max(0.025);
        let total_tau = self.dust_optical_depth + 0.01;
        beer_lambert(1.0, total_tau * airmass)
    }
}