asteroidsfactory 0.0.3

Asteroid factory — classify, build and catalogue asteroids of any type: near-Earth, main belt, trojan, centaur, binary, rubble pile, metallic, and potentially hazardous.
Documentation
use crate::config::parameters::*;

pub fn absolute_magnitude(radius: f64, albedo: f64) -> f64 {
    let diameter_km = 2.0 * radius / 1000.0;
    let p = albedo.max(0.001);
    15.618 - 5.0 * diameter_km.log10() - 2.5 * p.log10()
}

pub fn apparent_magnitude(
    h: f64,
    distance_sun_au: f64,
    distance_obs_au: f64,
    phase_angle_rad: f64,
) -> f64 {
    let phi = bowell_phase_function(phase_angle_rad);
    h + 5.0 * (distance_sun_au * distance_obs_au).log10() - 2.5 * phi.log10()
}

fn bowell_phase_function(alpha: f64) -> f64 {
    let a1 = (-3.33 * alpha.tan().powf(0.63)).exp();
    let a2 = (-1.87 * alpha.tan().powf(1.22)).exp();
    let g = 0.15;
    (1.0 - g) * a1 + g * a2
}

pub fn diameter_from_h_albedo(h: f64, albedo: f64) -> f64 {
    let p = albedo.max(0.001);
    1329.0 / p.sqrt() * 10.0_f64.powf(-0.2 * h)
}

pub fn thermal_flux(radius: f64, temperature: f64, wavelength_um: f64, distance_m: f64) -> f64 {
    let lambda = wavelength_um * 1.0e-6;
    let h_planck = 6.626e-34;
    let c = 3.0e8;
    let k_b = 1.381e-23;
    let x = h_planck * c / (lambda * k_b * temperature);
    let b_lambda = if x > 500.0 {
        0.0
    } else {
        2.0 * h_planck * c * c / (lambda.powi(5)) / (x.exp() - 1.0)
    };
    let solid_angle = PI * radius * radius / (distance_m * distance_m);
    b_lambda * solid_angle
}

pub fn geometric_albedo_from_bond(bond_albedo: f64, phase_integral: f64) -> f64 {
    if phase_integral < 1.0e-10 {
        bond_albedo
    } else {
        bond_albedo / phase_integral
    }
}