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)
}
}