Initial Commit

backend/src/astronomy/visibility.rs

@@ -0,0 +1,217 @@
+use chrono::{DateTime, Duration, Utc};
+use serde::{Deserialize, Serialize};
+
+use crate::config::{LAT, LON, MIN_ALT_DEG};
+use super::{
+    coords::{airmass, extinction_mag, radec_to_altaz},
+    horizon::{horizon_alt, HorizonPoint},
+    lunar::{moon_altitude, moon_separation},
+    time::{julian_date, local_sidereal_time},
+};
+
+#[derive(Debug, Clone, Serialize, Deserialize)]
+pub struct CurvePoint {
+    pub utc: DateTime<Utc>,
+    pub alt_deg: f64,
+    pub az_deg: f64,
+    pub airmass: f64,
+    pub above_custom_horizon: bool,
+    pub moon_alt_deg: f64,
+}
+
+#[derive(Debug, Clone, Serialize, Deserialize)]
+pub struct VisibilitySummary {
+    pub max_alt_deg: f64,
+    pub transit_utc: Option<DateTime<Utc>>,
+    pub rise_utc: Option<DateTime<Utc>>,
+    pub set_utc: Option<DateTime<Utc>>,
+    pub best_start_utc: Option<DateTime<Utc>>,
+    pub best_end_utc: Option<DateTime<Utc>>,
+    pub usable_min: u32,
+    pub is_visible_tonight: bool,
+    pub meridian_flip_utc: Option<DateTime<Utc>>,
+    pub airmass_at_transit: f64,
+    pub extinction_at_transit: f64,
+    pub moon_sep_deg: f64,
+    pub curve: Vec<CurvePoint>,
+}
+
+pub struct TonightWindow {
+    pub dusk: DateTime<Utc>,
+    pub dawn: DateTime<Utc>,
+}
+
+pub struct MoonState {
+    pub ra_deg: f64,
+    pub dec_deg: f64,
+    pub illumination: f64,
+    pub alt_at_midnight: f64,
+}
+
+/// Compute full visibility summary for a catalog object during tonight's window.
+/// step_minutes: resolution for the altitude curve (1 for detailed view, 10 for precompute cache).
+pub fn compute_visibility(
+    ra_deg: f64,
+    dec_deg: f64,
+    window: &TonightWindow,
+    horizon: &[HorizonPoint],
+    moon: &MoonState,
+) -> VisibilitySummary {
+    compute_visibility_with_step(ra_deg, dec_deg, window, horizon, moon, 10)
+}
+
+pub fn compute_visibility_with_step(
+    ra_deg: f64,
+    dec_deg: f64,
+    window: &TonightWindow,
+    horizon: &[HorizonPoint],
+    moon: &MoonState,
+    step_minutes: i64,
+) -> VisibilitySummary {
+    let step = Duration::minutes(step_minutes);
+    let mut t = window.dusk;
+
+    let mut curve = Vec::new();
+    let mut max_alt = f64::NEG_INFINITY;
+    let mut transit_utc: Option<DateTime<Utc>> = None;
+    let mut rise_utc: Option<DateTime<Utc>> = None;
+    let mut set_utc: Option<DateTime<Utc>> = None;
+    let mut best_start: Option<DateTime<Utc>> = None;
+    let mut best_end: Option<DateTime<Utc>> = None;
+    let mut usable_min = 0u32;
+    let mut prev_alt = f64::NEG_INFINITY;
+
+    while t <= window.dawn {
+        let jd = julian_date(t);
+        let lst = local_sidereal_time(jd, LON);
+        let (alt, az) = radec_to_altaz(ra_deg, dec_deg, lst, LAT);
+        let am = airmass(alt);
+        let h_alt = horizon_alt(az, horizon);
+        let above = alt > h_alt.max(MIN_ALT_DEG);
+        let moon_alt = moon_altitude(jd, LAT, LON);
+
+        curve.push(CurvePoint {
+            utc: t,
+            alt_deg: alt,
+            az_deg: az,
+            airmass: am,
+            above_custom_horizon: above,
+            moon_alt_deg: moon_alt,
+        });
+
+        if alt > max_alt {
+            max_alt = alt;
+            transit_utc = Some(t);
+        }
+
+        // Rise: first crossing above effective horizon
+        if prev_alt <= h_alt.max(MIN_ALT_DEG) && alt > h_alt.max(MIN_ALT_DEG) && rise_utc.is_none() {
+            rise_utc = Some(t);
+        }
+        // Set: last time we were above horizon
+        if alt > h_alt.max(MIN_ALT_DEG) {
+            set_utc = Some(t);
+        }
+
+        // Best window: above 30°
+        if alt > 30.0 {
+            if best_start.is_none() {
+                best_start = Some(t);
+            }
+            best_end = Some(t);
+            usable_min += step_minutes as u32;
+        }
+
+        prev_alt = alt;
+        t += step;
+    }
+
+    let is_visible = usable_min > 0 || rise_utc.is_some();
+
+    let airmass_transit = transit_utc
+        .map(|tr| {
+            let jd = julian_date(tr);
+            let lst = local_sidereal_time(jd, LON);
+            let (alt, _) = radec_to_altaz(ra_deg, dec_deg, lst, LAT);
+            airmass(alt)
+        })
+        .unwrap_or(40.0);
+
+    let extinction_transit = extinction_mag(
+        transit_utc
+            .map(|tr| {
+                let jd = julian_date(tr);
+                let lst = local_sidereal_time(jd, LON);
+                let (alt, _) = radec_to_altaz(ra_deg, dec_deg, lst, LAT);
+                alt
+            })
+            .unwrap_or(0.0),
+    );
+
+    let moon_sep = moon_separation(moon.ra_deg, moon.dec_deg, ra_deg, dec_deg);
+
+    // Meridian flip: transit + time for HA to reach +5°
+    // 5° of HA = 5/360 * 86400 = 1200 seconds
+    let meridian_flip = transit_utc.map(|tr| tr + Duration::seconds(1200));
+
+    VisibilitySummary {
+        max_alt_deg: if max_alt == f64::NEG_INFINITY { 0.0 } else { max_alt },
+        transit_utc,
+        rise_utc,
+        set_utc,
+        best_start_utc: best_start,
+        best_end_utc: best_end,
+        usable_min,
+        is_visible_tonight: is_visible,
+        meridian_flip_utc: meridian_flip,
+        airmass_at_transit: airmass_transit,
+        extinction_at_transit: extinction_transit,
+        moon_sep_deg: moon_sep,
+        curve,
+    }
+}
+
+/// Find the longest continuous true-dark window (sun < -18° AND moon below horizon).
+pub fn true_dark_window(
+    dusk: DateTime<Utc>,
+    dawn: DateTime<Utc>,
+    lat: f64,
+    lon: f64,
+) -> Option<(DateTime<Utc>, DateTime<Utc>)> {
+    let step = Duration::minutes(5);
+    let mut t = dusk;
+    let mut best: Option<(DateTime<Utc>, DateTime<Utc>)> = None;
+    let mut current_start: Option<DateTime<Utc>> = None;
+    let mut best_duration = Duration::zero();
+
+    while t <= dawn {
+        let jd = julian_date(t);
+        let moon_alt = moon_altitude(jd, lat, lon);
+        let is_dark = moon_alt < 0.0;
+
+        if is_dark {
+            if current_start.is_none() {
+                current_start = Some(t);
+            }
+        } else if let Some(start) = current_start {
+            let dur = t - start;
+            if dur > best_duration {
+                best_duration = dur;
+                best = Some((start, t));
+            }
+            current_start = None;
+        }
+
+        t += step;
+    }
+
+    // Check if still in dark window at dawn
+    if let Some(start) = current_start {
+        let dur = dawn - start;
+        if dur > best_duration {
+            best = Some((start, dawn));
+        }
+    }
+
+    best
+}