zwo_mount_control 0.2.1

//! Coordinate conversion utilities for telescope mount control.
//!
//! This module provides functions for converting between various coordinate systems
//! used in astronomy:
//!
//! - **Equatorial coordinates**: Right Ascension (RA) and Declination (Dec)
//! - **Horizontal coordinates**: Azimuth (Az) and Altitude (Alt)
//! - **HMS/DMS**: Hours-Minutes-Seconds and Degrees-Minutes-Seconds formats
//!
//! # Examples
//!
//! ```
//! use zwo_mount_control::coordinates::Coordinates;
//!
//! // Convert HMS to decimal hours
//! let ra = Coordinates::hms_to_decimal(18, 36, 56.0);
//! assert!((ra - 18.6155556).abs() < 0.0001);
//!
//! // Convert DMS to decimal degrees
//! let dec = Coordinates::dms_to_decimal(38, 47, 1.0);
//! assert!((dec - 38.7836111).abs() < 0.0001);
//! ```

use std::f64::consts::PI;

/// Coordinate conversion utilities.
pub struct Coordinates;

impl Coordinates {
    // ============================================
    // HMS (Hours-Minutes-Seconds) Conversions
    // ============================================

    /// Convert hours, minutes, seconds to decimal hours.
    ///
    /// # Arguments
    /// * `hours` - Hours (0-23)
    /// * `minutes` - Minutes (0-59)
    /// * `seconds` - Seconds (0-59.999...)
    ///
    /// # Returns
    /// Decimal hours (0.0 to 24.0)
    pub fn hms_to_decimal(hours: u8, minutes: u8, seconds: f64) -> f64 {
        hours as f64 + minutes as f64 / 60.0 + seconds / 3600.0
    }

    /// Convert decimal hours to hours, minutes, seconds.
    ///
    /// # Arguments
    /// * `decimal_hours` - Decimal hours (will be normalized to 0-24 range)
    ///
    /// # Returns
    /// Tuple of (hours, minutes, seconds)
    pub fn decimal_to_hms(decimal_hours: f64) -> (u8, u8, f64) {
        let normalized = decimal_hours.rem_euclid(24.0);
        let hours = normalized.floor() as u8;
        let remaining = (normalized - hours as f64) * 60.0;
        let minutes = remaining.floor() as u8;
        let seconds = (remaining - minutes as f64) * 60.0;
        (hours, minutes, seconds)
    }

    /// Format decimal hours as HH:MM:SS string.
    pub fn decimal_hours_to_string(decimal_hours: f64) -> String {
        let (h, m, s) = Self::decimal_to_hms(decimal_hours);
        format!("{:02}:{:02}:{:05.2}", h, m, s)
    }

    // ============================================
    // DMS (Degrees-Minutes-Seconds) Conversions
    // ============================================

    /// Convert degrees, minutes, seconds to decimal degrees.
    ///
    /// # Arguments
    /// * `degrees` - Degrees (can be negative for southern/western coordinates)
    /// * `minutes` - Minutes (0-59)
    /// * `seconds` - Seconds (0-59.999...)
    ///
    /// # Returns
    /// Decimal degrees
    pub fn dms_to_decimal(degrees: i16, minutes: u8, seconds: f64) -> f64 {
        let sign = if degrees >= 0 { 1.0 } else { -1.0 };
        sign * (degrees.abs() as f64 + minutes as f64 / 60.0 + seconds / 3600.0)
    }

    /// Convert decimal degrees to degrees, minutes, seconds.
    ///
    /// # Arguments
    /// * `decimal_degrees` - Decimal degrees (can be negative)
    ///
    /// # Returns
    /// Tuple of (degrees, minutes, seconds) where degrees is signed
    pub fn decimal_to_dms(decimal_degrees: f64) -> (i16, u8, f64) {
        let sign = if decimal_degrees >= 0.0 { 1 } else { -1 };
        let abs_value = decimal_degrees.abs();
        let degrees = abs_value.floor() as i16 * sign;
        let remaining = (abs_value - abs_value.floor()) * 60.0;
        let minutes = remaining.floor() as u8;
        let seconds = (remaining - remaining.floor()) * 60.0;
        (degrees, minutes, seconds)
    }

    /// Format decimal degrees as sDD*MM:SS string.
    pub fn decimal_degrees_to_string(decimal_degrees: f64) -> String {
        let (d, m, s) = Self::decimal_to_dms(decimal_degrees);
        let sign = if decimal_degrees >= 0.0 { '+' } else { '-' };
        format!("{}{}°{:02}'{:05.2}\"", sign, d.abs(), m, s)
    }

    // ============================================
    // RA/Dec Conversions
    // ============================================

    /// Convert Right Ascension from decimal hours to degrees.
    pub fn ra_hours_to_degrees(ra_hours: f64) -> f64 {
        ra_hours * 15.0
    }

    /// Convert Right Ascension from degrees to decimal hours.
    pub fn ra_degrees_to_hours(ra_degrees: f64) -> f64 {
        ra_degrees / 15.0
    }

    /// Normalize Right Ascension to 0-24 hours range.
    pub fn normalize_ra(ra_hours: f64) -> f64 {
        ra_hours.rem_euclid(24.0)
    }

    /// Normalize Declination to -90 to +90 degrees range.
    pub fn normalize_dec(dec_degrees: f64) -> f64 {
        dec_degrees.clamp(-90.0, 90.0)
    }

    // ============================================
    // Az/Alt Conversions
    // ============================================

    /// Normalize Azimuth to 0-360 degrees range.
    pub fn normalize_azimuth(azimuth: f64) -> f64 {
        azimuth.rem_euclid(360.0)
    }

    /// Normalize Altitude to -90 to +90 degrees range.
    pub fn normalize_altitude(altitude: f64) -> f64 {
        altitude.clamp(-90.0, 90.0)
    }

    // ============================================
    // Equatorial <-> Horizontal Conversions
    // ============================================

    /// Convert equatorial coordinates (RA/Dec) to horizontal coordinates (Az/Alt).
    ///
    /// # Arguments
    /// * `ra_hours` - Right Ascension in decimal hours
    /// * `dec_degrees` - Declination in decimal degrees
    /// * `latitude` - Observer's latitude in decimal degrees
    /// * `lst_hours` - Local Sidereal Time in decimal hours
    ///
    /// # Returns
    /// Tuple of (azimuth, altitude) in decimal degrees
    pub fn equatorial_to_horizontal(
        ra_hours: f64,
        dec_degrees: f64,
        latitude: f64,
        lst_hours: f64,
    ) -> (f64, f64) {
        // Convert to radians
        let dec = dec_degrees.to_radians();
        let lat = latitude.to_radians();

        // Calculate Hour Angle
        let ha_hours = lst_hours - ra_hours;
        let ha = (ha_hours * 15.0).to_radians();

        // Calculate altitude
        let sin_alt = dec.sin() * lat.sin() + dec.cos() * lat.cos() * ha.cos();
        let alt = sin_alt.asin();

        // Calculate azimuth
        let cos_az = (dec.sin() - alt.sin() * lat.sin()) / (alt.cos() * lat.cos());
        let cos_az_clamped = cos_az.clamp(-1.0, 1.0);
        let mut az = cos_az_clamped.acos();

        // Adjust azimuth based on hour angle
        if ha.sin() > 0.0 {
            az = 2.0 * PI - az;
        }

        (az.to_degrees(), alt.to_degrees())
    }

    /// Convert horizontal coordinates (Az/Alt) to equatorial coordinates (RA/Dec).
    ///
    /// # Arguments
    /// * `azimuth` - Azimuth in decimal degrees (0 = North, 90 = East)
    /// * `altitude` - Altitude in decimal degrees
    /// * `latitude` - Observer's latitude in decimal degrees
    /// * `lst_hours` - Local Sidereal Time in decimal hours
    ///
    /// # Returns
    /// Tuple of (ra_hours, dec_degrees)
    pub fn horizontal_to_equatorial(
        azimuth: f64,
        altitude: f64,
        latitude: f64,
        lst_hours: f64,
    ) -> (f64, f64) {
        // Convert to radians
        let az = azimuth.to_radians();
        let alt = altitude.to_radians();
        let lat = latitude.to_radians();

        // Calculate declination
        let sin_dec = alt.sin() * lat.sin() + alt.cos() * lat.cos() * az.cos();
        let dec = sin_dec.asin();

        // Calculate hour angle
        let cos_ha = (alt.sin() - dec.sin() * lat.sin()) / (dec.cos() * lat.cos());
        let cos_ha_clamped = cos_ha.clamp(-1.0, 1.0);
        let mut ha = cos_ha_clamped.acos();

        // Adjust hour angle based on azimuth
        if az.sin() > 0.0 {
            ha = 2.0 * PI - ha;
        }

        // Calculate RA from hour angle and LST
        let ha_hours = ha.to_degrees() / 15.0;
        let ra_hours = Self::normalize_ra(lst_hours - ha_hours);

        (ra_hours, dec.to_degrees())
    }

    // ============================================
    // Local Sidereal Time
    // ============================================

    /// Calculate Local Sidereal Time from Julian Date and longitude.
    ///
    /// # Arguments
    /// * `jd` - Julian Date
    /// * `longitude` - Observer's longitude in decimal degrees (positive = East)
    ///
    /// # Returns
    /// Local Sidereal Time in decimal hours
    pub fn julian_to_lst(jd: f64, longitude: f64) -> f64 {
        // Calculate Julian centuries from J2000.0
        let t = (jd - 2451545.0) / 36525.0;

        // Greenwich Mean Sidereal Time at 0h UT
        let gmst = 280.46061837 + 360.98564736629 * (jd - 2451545.0) + 0.000387933 * t * t
            - t * t * t / 38710000.0;

        // Convert to hours and add longitude
        let lst = (gmst + longitude) / 15.0;

        Self::normalize_ra(lst)
    }

    /// Calculate Julian Date from year, month, day, and fractional day.
    ///
    /// # Arguments
    /// * `year` - Year (e.g., 2024)
    /// * `month` - Month (1-12)
    /// * `day` - Day with fractional part for time
    ///
    /// # Returns
    /// Julian Date
    pub fn to_julian_date(year: i32, month: u8, day: f64) -> f64 {
        let (y, m) = if month <= 2 {
            (year - 1, month + 12)
        } else {
            (year, month)
        };

        let a = (y as f64 / 100.0).floor();
        let b = 2.0 - a + (a / 4.0).floor();

        (365.25 * (y as f64 + 4716.0)).floor() + (30.6001 * (m as f64 + 1.0)).floor() + day + b
            - 1524.5
    }

    // ============================================
    // Angular Separation
    // ============================================

    /// Calculate angular separation between two points on the celestial sphere.
    ///
    /// # Arguments
    /// * `ra1`, `dec1` - First position (RA in hours, Dec in degrees)
    /// * `ra2`, `dec2` - Second position (RA in hours, Dec in degrees)
    ///
    /// # Returns
    /// Angular separation in degrees
    pub fn angular_separation(ra1: f64, dec1: f64, ra2: f64, dec2: f64) -> f64 {
        let ra1_rad = Self::ra_hours_to_degrees(ra1).to_radians();
        let dec1_rad = dec1.to_radians();
        let ra2_rad = Self::ra_hours_to_degrees(ra2).to_radians();
        let dec2_rad = dec2.to_radians();

        let cos_sep = dec1_rad.sin() * dec2_rad.sin()
            + dec1_rad.cos() * dec2_rad.cos() * (ra1_rad - ra2_rad).cos();

        cos_sep.clamp(-1.0, 1.0).acos().to_degrees()
    }

    // ============================================
    // Rate Conversions (for satellite tracking)
    // ============================================

    /// Convert angular rate from degrees per second to arcseconds per second.
    pub fn deg_per_sec_to_arcsec_per_sec(deg_per_sec: f64) -> f64 {
        deg_per_sec * 3600.0
    }

    /// Convert angular rate from arcseconds per second to degrees per second.
    pub fn arcsec_per_sec_to_deg_per_sec(arcsec_per_sec: f64) -> f64 {
        arcsec_per_sec / 3600.0
    }

    /// Sidereal rate in arcseconds per second (approximately 15.041 "/s).
    pub const SIDEREAL_RATE_ARCSEC_PER_SEC: f64 = 15.041067;

    /// Sidereal rate in degrees per hour.
    pub const SIDEREAL_RATE_DEG_PER_HOUR: f64 = 15.0;
}

/// Represents a position in equatorial coordinates.
#[derive(Debug, Clone, Copy, PartialEq)]
pub struct EquatorialPosition {
    /// Right Ascension in decimal hours (0-24)
    pub ra: f64,
    /// Declination in decimal degrees (-90 to +90)
    pub dec: f64,
}

impl EquatorialPosition {
    /// Create a new equatorial position.
    pub fn new(ra: f64, dec: f64) -> Self {
        Self {
            ra: Coordinates::normalize_ra(ra),
            dec: Coordinates::normalize_dec(dec),
        }
    }

    /// Create from HMS and DMS components.
    pub fn from_hms_dms(ra_h: u8, ra_m: u8, ra_s: f64, dec_d: i16, dec_m: u8, dec_s: f64) -> Self {
        Self::new(
            Coordinates::hms_to_decimal(ra_h, ra_m, ra_s),
            Coordinates::dms_to_decimal(dec_d, dec_m, dec_s),
        )
    }

    /// Convert to horizontal coordinates.
    pub fn to_horizontal(&self, latitude: f64, lst_hours: f64) -> HorizontalPosition {
        let (az, alt) =
            Coordinates::equatorial_to_horizontal(self.ra, self.dec, latitude, lst_hours);
        HorizontalPosition::new(az, alt)
    }

    /// Get RA as HMS components.
    pub fn ra_hms(&self) -> (u8, u8, f64) {
        Coordinates::decimal_to_hms(self.ra)
    }

    /// Get Dec as DMS components.
    pub fn dec_dms(&self) -> (i16, u8, f64) {
        Coordinates::decimal_to_dms(self.dec)
    }

    /// Calculate angular separation to another position.
    pub fn angular_separation(&self, other: &EquatorialPosition) -> f64 {
        Coordinates::angular_separation(self.ra, self.dec, other.ra, other.dec)
    }
}

impl std::fmt::Display for EquatorialPosition {
    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
        let (ra_h, ra_m, ra_s) = self.ra_hms();
        let (dec_d, dec_m, dec_s) = self.dec_dms();
        let dec_sign = if self.dec >= 0.0 { '+' } else { '-' };
        write!(
            f,
            "RA {:02}h{:02}m{:05.2}s, Dec {}{}°{:02}'{:05.2}\"",
            ra_h,
            ra_m,
            ra_s,
            dec_sign,
            dec_d.abs(),
            dec_m,
            dec_s
        )
    }
}

/// Represents a position in horizontal (Alt-Az) coordinates.
#[derive(Debug, Clone, Copy, PartialEq)]
pub struct HorizontalPosition {
    /// Azimuth in decimal degrees (0-360, 0 = North, 90 = East)
    pub azimuth: f64,
    /// Altitude in decimal degrees (-90 to +90, 0 = horizon, 90 = zenith)
    pub altitude: f64,
}

impl HorizontalPosition {
    /// Create a new horizontal position.
    pub fn new(azimuth: f64, altitude: f64) -> Self {
        Self {
            azimuth: Coordinates::normalize_azimuth(azimuth),
            altitude: Coordinates::normalize_altitude(altitude),
        }
    }

    /// Create from DMS components.
    pub fn from_dms(az_d: u16, az_m: u8, az_s: f64, alt_d: i8, alt_m: u8, alt_s: f64) -> Self {
        Self::new(
            az_d as f64 + az_m as f64 / 60.0 + az_s / 3600.0,
            Coordinates::dms_to_decimal(alt_d as i16, alt_m, alt_s),
        )
    }

    /// Convert to equatorial coordinates.
    pub fn to_equatorial(&self, latitude: f64, lst_hours: f64) -> EquatorialPosition {
        let (ra, dec) =
            Coordinates::horizontal_to_equatorial(self.azimuth, self.altitude, latitude, lst_hours);
        EquatorialPosition::new(ra, dec)
    }

    /// Get azimuth as DMS components.
    pub fn azimuth_dms(&self) -> (u16, u8, f64) {
        let degrees = self.azimuth.floor() as u16;
        let remaining = (self.azimuth - degrees as f64) * 60.0;
        let minutes = remaining.floor() as u8;
        let seconds = (remaining - minutes as f64) * 60.0;
        (degrees, minutes, seconds)
    }

    /// Get altitude as DMS components.
    pub fn altitude_dms(&self) -> (i8, u8, f64) {
        let (d, m, s) = Coordinates::decimal_to_dms(self.altitude);
        (d as i8, m, s)
    }

    /// Check if the position is above the horizon.
    pub fn is_above_horizon(&self) -> bool {
        self.altitude > 0.0
    }

    /// Check if the position is above a given elevation limit.
    pub fn is_above_limit(&self, limit_degrees: f64) -> bool {
        self.altitude > limit_degrees
    }
}

impl std::fmt::Display for HorizontalPosition {
    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
        let (az_d, az_m, az_s) = self.azimuth_dms();
        let (alt_d, alt_m, alt_s) = self.altitude_dms();
        let alt_sign = if self.altitude >= 0.0 { '+' } else { '-' };
        write!(
            f,
            "Az {}°{:02}'{:05.2}\", Alt {}{}°{:02}'{:05.2}\"",
            az_d,
            az_m,
            az_s,
            alt_sign,
            alt_d.abs(),
            alt_m,
            alt_s
        )
    }
}

/// Represents tracking rates in RA and Dec.
#[derive(Debug, Clone, Copy, PartialEq)]
pub struct TrackingRates {
    /// RA rate in arcseconds per second (positive = westward/tracking direction)
    pub ra_rate: f64,
    /// Dec rate in arcseconds per second (positive = northward)
    pub dec_rate: f64,
}

impl TrackingRates {
    /// Create new tracking rates.
    pub fn new(ra_rate: f64, dec_rate: f64) -> Self {
        Self { ra_rate, dec_rate }
    }

    /// Sidereal tracking rate (for stars).
    pub fn sidereal() -> Self {
        Self {
            ra_rate: Coordinates::SIDEREAL_RATE_ARCSEC_PER_SEC,
            dec_rate: 0.0,
        }
    }

    /// No tracking (stationary).
    pub fn zero() -> Self {
        Self {
            ra_rate: 0.0,
            dec_rate: 0.0,
        }
    }

    /// Check if rates are within typical satellite tracking range.
    pub fn is_valid_for_satellite(&self) -> bool {
        // Satellites typically move at rates up to ~1-2 degrees per second
        let max_rate = 7200.0; // 2 deg/s in arcsec/s
        self.ra_rate.abs() <= max_rate && self.dec_rate.abs() <= max_rate
    }
}

impl Default for TrackingRates {
    fn default() -> Self {
        Self::sidereal()
    }
}

#[cfg(test)]
mod tests {
    use super::*;

    #[test]
    fn test_hms_to_decimal() {
        let ra = Coordinates::hms_to_decimal(12, 30, 0.0);
        assert!((ra - 12.5).abs() < 0.0001);

        let ra = Coordinates::hms_to_decimal(18, 36, 56.0);
        assert!((ra - 18.6155556).abs() < 0.0001);
    }

    #[test]
    fn test_decimal_to_hms() {
        let (h, m, s) = Coordinates::decimal_to_hms(12.5);
        assert_eq!(h, 12);
        assert_eq!(m, 30);
        assert!(s.abs() < 0.01);
    }

    #[test]
    fn test_dms_to_decimal() {
        let dec = Coordinates::dms_to_decimal(45, 30, 0.0);
        assert!((dec - 45.5).abs() < 0.0001);

        let dec = Coordinates::dms_to_decimal(-23, 26, 21.0);
        assert!((dec - (-23.439167)).abs() < 0.0001);
    }

    #[test]
    fn test_decimal_to_dms() {
        let (d, m, s) = Coordinates::decimal_to_dms(45.5);
        assert_eq!(d, 45);
        assert_eq!(m, 30);
        assert!(s.abs() < 0.01);

        let (d, m, s) = Coordinates::decimal_to_dms(-23.5);
        assert_eq!(d, -23);
        assert_eq!(m, 30);
        assert!(s.abs() < 0.01);
    }

    #[test]
    fn test_equatorial_position() {
        let pos = EquatorialPosition::from_hms_dms(18, 36, 56.0, 38, 47, 1.0);
        assert!((pos.ra - 18.6155556).abs() < 0.0001);
        assert!((pos.dec - 38.7836111).abs() < 0.0001);
    }

    #[test]
    fn test_horizontal_position() {
        let pos = HorizontalPosition::new(180.0, 45.0);
        assert!(pos.is_above_horizon());
        assert!(pos.is_above_limit(30.0));
        assert!(!pos.is_above_limit(60.0));
    }

    #[test]
    fn test_angular_separation() {
        // Same position should have 0 separation
        let sep = Coordinates::angular_separation(12.0, 45.0, 12.0, 45.0);
        assert!(sep.abs() < 0.0001);

        // Test known separation
        let sep = Coordinates::angular_separation(0.0, 0.0, 6.0, 0.0);
        assert!((sep - 90.0).abs() < 0.001);
    }

    #[test]
    fn test_normalize_ra() {
        assert!((Coordinates::normalize_ra(25.0) - 1.0).abs() < 0.0001);
        assert!((Coordinates::normalize_ra(-1.0) - 23.0).abs() < 0.0001);
    }

    #[test]
    fn test_tracking_rates() {
        let sidereal = TrackingRates::sidereal();
        assert!((sidereal.ra_rate - Coordinates::SIDEREAL_RATE_ARCSEC_PER_SEC).abs() < 0.001);
        assert!(sidereal.dec_rate.abs() < 0.001);
    }
}