jupiters 0.0.3

Jupiter celestial simulation crate for the MilkyWay SolarSystem workspace
Documentation
use crate::temporal::epoch::J2000EPOCH;
use sciforge::hub::prelude::constants::AU;

pub const JUPITERAXIALTILTDEG: f64 = 3.13;

pub struct SolarPosition {
    pub azimuthdeg: f64,
    pub elevationdeg: f64,
    pub distanceau: f64,
    pub direction: [f64; 3],
}

impl SolarPosition {
    pub fn compute(juliandate: f64, latdeg: f64, londeg: f64) -> Self {
        let d = juliandate - J2000EPOCH;
        let jupiterorbitalperioddays = 4332.59;
        let meananomaly = (20.020 + 360.0 * d / jupiterorbitalperioddays) % 360.0;
        let mrad = meananomaly.to_radians();
        let equationofcenter =
            5.5549 * mrad.sin() + 0.1726 * (2.0 * mrad).sin() + 0.004 * (3.0 * mrad).sin();
        let eclipticlon = (meananomaly + equationofcenter + 180.0 + 14.331) % 360.0;
        let eclrad = eclipticlon.to_radians();
        let obliquityrad = JUPITERAXIALTILTDEG.to_radians();
        let declination = (obliquityrad.sin() * eclrad.sin()).asin();
        let rightascension = (obliquityrad.cos() * eclrad.sin()).atan2(eclrad.cos());
        let jupitersiderealrotperiodh = 9.925;
        let siderealtime =
            (280.16 + 360.0 / (jupitersiderealrotperiodh / 24.0) * d + londeg) % 360.0;
        let hourangle = siderealtime.to_radians() - rightascension;
        let latrad = latdeg.to_radians();
        let altitude = (latrad.sin() * declination.sin()
            + latrad.cos() * declination.cos() * hourangle.cos())
        .asin();
        let azimuth = (declination.sin() - altitude.sin() * latrad.sin())
            .atan2(altitude.cos() * latrad.cos());
        let dirx = altitude.cos() * azimuth.sin();
        let diry = altitude.sin();
        let dirz = altitude.cos() * azimuth.cos();
        let distanceau = 5.2044 - 0.2543 * mrad.cos();
        Self {
            azimuthdeg: azimuth.to_degrees(),
            elevationdeg: altitude.to_degrees(),
            distanceau,
            direction: [dirx, diry, dirz],
        }
    }

    pub fn distancem(&self) -> f64 {
        self.distanceau * AU
    }

    pub fn isabovehorizon(&self) -> bool {
        self.elevationdeg > 0.0
    }

    pub fn solarfluxwm2(&self) -> f64 {
        1361.0 / (self.distanceau * self.distanceau)
    }
}