use sciforge::hub::prelude::constants::N_A;
use sciforge::hub::prelude::constants::elements::atomic_mass;
pub struct MolecularSpecies {
pub name: &'static str,
pub symbol: &'static str,
pub molar_mass_kg_mol: f64,
pub volume_fraction: f64,
pub refractive_index_stp: f64,
pub depolarization_factor: f64,
}
pub struct AtmosphereEndpoint {
pub planet_radius_m: f64,
pub atmosphere_height_m: f64,
pub sea_level_pressure_pa: f64,
pub sea_level_temperature_k: f64,
pub sea_level_number_density_m3: f64,
pub mean_molar_mass_kg_mol: f64,
pub rayleigh_scale_height_m: f64,
pub cloud_optical_depth: f64,
pub rayleigh_coefficients_rgb: [f64; 3],
pub species: Vec<MolecularSpecies>,
pub sun_irradiance_w_m2: f64,
}
fn co2_molar() -> f64 {
(atomic_mass(6) + 2.0 * atomic_mass(8)) * 1e-3
}
fn n2_molar() -> f64 {
2.0 * atomic_mass(7) * 1e-3
}
fn so2_molar() -> f64 {
(atomic_mass(16) + 2.0 * atomic_mass(8)) * 1e-3
}
fn ar_molar() -> f64 {
atomic_mass(18) * 1e-3
}
fn h2o_molar() -> f64 {
(2.0 * atomic_mass(1) + atomic_mass(8)) * 1e-3
}
fn co_molar() -> f64 {
(atomic_mass(6) + atomic_mass(8)) * 1e-3
}
fn he_molar() -> f64 {
atomic_mass(2) * 1e-3
}
fn mean_molar_mass() -> f64 {
crate::CO2_FRACTION * co2_molar()
+ crate::N2_FRACTION * n2_molar()
+ 0.000_150 * so2_molar()
+ crate::AR_FRACTION * ar_molar()
+ 0.000_020 * h2o_molar()
+ 0.000_017 * co_molar()
+ 0.000_012 * he_molar()
}
fn sea_level_number_density() -> f64 {
let m_air = mean_molar_mass();
let rho = crate::SURFACE_PRESSURE_PA / (crate::R_SPECIFIC * crate::MEAN_TEMPERATURE_K);
N_A * rho / m_air
}
fn rayleigh_beta_rgb(n_density: f64) -> [f64; 3] {
let n_minus_1 = 4.49e-4;
let ns2 = (2.0 * n_minus_1) * (2.0 * n_minus_1);
let coeff = 8.0 * std::f64::consts::PI.powi(3) * ns2 / (3.0 * n_density);
[
coeff / (680e-9_f64).powi(4),
coeff / (550e-9_f64).powi(4),
coeff / (440e-9_f64).powi(4),
]
}
impl AtmosphereEndpoint {
pub fn venus() -> Self {
let n_density = sea_level_number_density();
let beta = rayleigh_beta_rgb(n_density);
let species = vec![
MolecularSpecies {
name: "Carbon dioxide",
symbol: "CO2",
molar_mass_kg_mol: co2_molar(),
volume_fraction: crate::CO2_FRACTION,
refractive_index_stp: 1.000_449_0,
depolarization_factor: 0.075,
},
MolecularSpecies {
name: "Dinitrogen",
symbol: "N2",
molar_mass_kg_mol: n2_molar(),
volume_fraction: crate::N2_FRACTION,
refractive_index_stp: 1.000_298_0,
depolarization_factor: 0.030,
},
MolecularSpecies {
name: "Sulfur dioxide",
symbol: "SO2",
molar_mass_kg_mol: so2_molar(),
volume_fraction: 0.000_150,
refractive_index_stp: 1.000_686_0,
depolarization_factor: 0.05,
},
MolecularSpecies {
name: "Argon",
symbol: "Ar",
molar_mass_kg_mol: ar_molar(),
volume_fraction: crate::AR_FRACTION,
refractive_index_stp: 1.000_281_0,
depolarization_factor: 0.0,
},
MolecularSpecies {
name: "Water vapor",
symbol: "H2O",
molar_mass_kg_mol: h2o_molar(),
volume_fraction: 0.000_020,
refractive_index_stp: 1.000_256_0,
depolarization_factor: 0.17,
},
MolecularSpecies {
name: "Carbon monoxide",
symbol: "CO",
molar_mass_kg_mol: co_molar(),
volume_fraction: 0.000_017,
refractive_index_stp: 1.000_338_0,
depolarization_factor: 0.030,
},
MolecularSpecies {
name: "Helium",
symbol: "He",
molar_mass_kg_mol: he_molar(),
volume_fraction: 0.000_012,
refractive_index_stp: 1.000_035_0,
depolarization_factor: 0.0,
},
];
Self {
planet_radius_m: crate::VENUS_RADIUS,
atmosphere_height_m: 250_000.0,
sea_level_pressure_pa: crate::SURFACE_PRESSURE_PA,
sea_level_temperature_k: crate::MEAN_TEMPERATURE_K,
sea_level_number_density_m3: n_density,
mean_molar_mass_kg_mol: mean_molar_mass(),
rayleigh_scale_height_m: crate::SCALE_HEIGHT_M,
cloud_optical_depth: crate::MEAN_CLOUD_OPTICAL_DEPTH,
rayleigh_coefficients_rgb: beta,
species,
sun_irradiance_w_m2: 1_361.0 / (crate::SEMI_MAJOR_AXIS_AU * crate::SEMI_MAJOR_AXIS_AU),
}
}
pub fn rayleigh_density(&self, altitude_m: f64) -> f64 {
(-altitude_m / self.rayleigh_scale_height_m).exp()
}
pub fn sky_color(&self, sun_elevation_deg: f64) -> [f64; 3] {
let el = sun_elevation_deg.clamp(-5.0, 90.0);
let elevation_factor = ((el + 5.0) / 95.0).sqrt();
let cloud_scatter = 1.0 - (-self.cloud_optical_depth * 0.002).exp();
let r = (0.95 * cloud_scatter * elevation_factor).clamp(0.0, 1.0);
let g = (0.75 * cloud_scatter * elevation_factor).clamp(0.0, 1.0);
let b = (0.45 * cloud_scatter * elevation_factor).clamp(0.0, 1.0);
if sun_elevation_deg < 0.0 {
let fade = ((sun_elevation_deg + 5.0) / 5.0).clamp(0.0, 1.0);
[r * fade, g * fade, b * fade]
} else {
[r, g, b]
}
}
pub fn direct_transmission(&self, sun_elevation_deg: f64) -> f64 {
let el = sun_elevation_deg.clamp(0.01, 90.0);
let airmass = 1.0 / el.to_radians().sin().max(0.01);
(-self.cloud_optical_depth * airmass).exp()
}
pub fn sky_luminance(&self, observer_alt: f64, cos_zenith: f64) -> f64 {
let cos_z = cos_zenith.max(0.01);
let molecular_tau: f64 = self
.species
.iter()
.map(|s| {
let n_s = self.sea_level_number_density_m3 * s.volume_fraction;
let delta = s.depolarization_factor;
let king = (6.0 + 3.0 * delta) / (6.0 - 7.0 * delta);
let ns = s.refractive_index_stp - 1.0;
8.0 * std::f64::consts::PI.powi(3) * (2.0 * ns).powi(2) * king / (3.0 * n_s)
})
.sum::<f64>()
/ (550e-9_f64).powi(4);
let alt_factor = (-observer_alt / self.rayleigh_scale_height_m).exp();
base_luminance(self.sun_irradiance_w_m2, molecular_tau, alt_factor, cos_z)
}
pub fn barometric_density(&self, altitude_m: f64) -> f64 {
self.sea_level_number_density_m3 * (-altitude_m / self.rayleigh_scale_height_m).exp()
}
pub fn shell_volume(&self, altitude_m: f64, thickness_m: f64) -> f64 {
let r = self.planet_radius_m + altitude_m;
4.0 * std::f64::consts::PI * r * r * thickness_m
}
}
fn base_luminance(irradiance: f64, tau: f64, alt_factor: f64, cos_z: f64) -> f64 {
irradiance * tau * alt_factor * cos_z * 0.001
}