Skip to main content

celestial_time/
julian.rs

1use crate::constants::{SECONDS_TO_DAYS, UNIX_EPOCH_JD};
2use celestial_core::constants::{J2000_JD, MJD_ZERO_POINT, SECONDS_PER_DAY_F64};
3use std::fmt;
4use std::ops::{Add, Sub};
5
6#[derive(Debug, Clone, Copy, PartialEq)]
7#[cfg_attr(feature = "serde", derive(serde::Serialize, serde::Deserialize))]
8pub struct JulianDate {
9    pub jd1: f64,
10    pub jd2: f64,
11}
12
13impl JulianDate {
14    pub fn new(jd1: f64, jd2: f64) -> Self {
15        Self { jd1, jd2 }
16    }
17
18    pub fn from_f64(jd: f64) -> Self {
19        Self::new(jd, 0.0)
20    }
21
22    pub fn j2000() -> Self {
23        Self::new(J2000_JD, 0.0)
24    }
25
26    pub fn unix_epoch() -> Self {
27        Self::new(UNIX_EPOCH_JD, 0.0)
28    }
29
30    pub fn jd1(&self) -> f64 {
31        self.jd1
32    }
33
34    pub fn jd2(&self) -> f64 {
35        self.jd2
36    }
37
38    pub fn to_f64(&self) -> f64 {
39        self.jd1 + self.jd2
40    }
41
42    pub fn add_days(&self, days: f64) -> Self {
43        Self::new(self.jd1, self.jd2 + days)
44    }
45
46    pub fn add_seconds(&self, seconds: f64) -> Self {
47        self.add_days(seconds * SECONDS_TO_DAYS)
48    }
49
50    pub fn from_calendar(year: i32, month: u8, day: u8, hour: u8, minute: u8, second: f64) -> Self {
51        // Algorithm matches ERFA's eraCal2jd + eraDtf2d convention:
52        // jd1 = full Julian Date at midnight (integer-ish)
53        // jd2 = fraction of day
54        let my = (month as i32 - 14) / 12;
55        let iypmy = year + my;
56
57        // Compute MJD for 0h of the given day (same algorithm as ERFA eraCal2jd)
58        let mjd = ((1461 * (iypmy + 4800)) / 4 + (367 * (month as i32 - 2 - 12 * my)) / 12
59            - (3 * ((iypmy + 4900) / 100)) / 4
60            + day as i32
61            - 2432076) as f64;
62
63        // Full Julian Date at midnight = MJD epoch + MJD
64        let jd1 = MJD_ZERO_POINT + mjd;
65
66        // Day fraction from time components
67        let jd2 = (60.0 * (60 * hour as i32 + minute as i32) as f64 + second) / SECONDS_PER_DAY_F64;
68
69        Self::new(jd1, jd2)
70    }
71
72    pub fn to_julian_year(&self) -> f64 {
73        const DAYS_PER_JULIAN_YEAR: f64 = 365.25;
74        2000.0 + (self.to_f64() - J2000_JD) / DAYS_PER_JULIAN_YEAR
75    }
76
77    pub fn from_julian_year(year: f64) -> Self {
78        const DAYS_PER_JULIAN_YEAR: f64 = 365.25;
79        let jd = J2000_JD + (year - 2000.0) * DAYS_PER_JULIAN_YEAR;
80        Self::from_f64(jd)
81    }
82}
83
84impl fmt::Display for JulianDate {
85    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
86        write!(f, "JD {:.9}", self.to_f64())
87    }
88}
89
90impl From<f64> for JulianDate {
91    fn from(jd: f64) -> Self {
92        Self::from_f64(jd)
93    }
94}
95
96impl Add<JulianDate> for JulianDate {
97    type Output = Self;
98
99    fn add(self, other: JulianDate) -> Self::Output {
100        Self::new(self.jd1 + other.jd1, self.jd2 + other.jd2)
101    }
102}
103
104impl Sub<JulianDate> for JulianDate {
105    type Output = Self;
106
107    fn sub(self, other: JulianDate) -> Self::Output {
108        Self::new(self.jd1 - other.jd1, self.jd2 - other.jd2)
109    }
110}
111
112#[cfg(test)]
113mod tests {
114    use super::*;
115
116    #[test]
117    fn test_julian_date_creation() {
118        let jd = JulianDate::new(J2000_JD, 0.5);
119        assert_eq!(jd.jd1(), J2000_JD);
120        assert_eq!(jd.jd2(), 0.5);
121        assert_eq!(jd.to_f64(), 2451545.5);
122    }
123
124    #[test]
125    fn test_j2000_epoch() {
126        let j2000 = JulianDate::j2000();
127        assert_eq!(j2000.to_f64(), J2000_JD);
128    }
129
130    #[test]
131    fn test_unix_epoch() {
132        let unix = JulianDate::unix_epoch();
133        assert_eq!(unix.to_f64(), crate::constants::UNIX_EPOCH_JD);
134    }
135
136    #[test]
137    fn test_arithmetic() {
138        let jd = JulianDate::new(J2000_JD, 0.0);
139        let jd_plus_day = jd.add_days(1.0);
140        assert_eq!(jd_plus_day.to_f64(), 2451546.0);
141
142        let jd_plus_hour = jd.add_seconds(3600.0);
143        assert!((jd_plus_hour.to_f64() - 2_451_545.041_666_666_5).abs() < 1e-15);
144    }
145
146    #[cfg(feature = "serde")]
147    #[test]
148    fn test_serde_round_trip() {
149        let test_cases = [
150            JulianDate::new(J2000_JD, 0.0),          // J2000.0
151            JulianDate::new(2451545.5, 0.123456789), // J2000.0 + 12h + fraction
152            JulianDate::new(2440587.5, 0.0),         // Unix epoch
153            JulianDate::new(J2000_JD, 0.999999999),  // High precision
154        ];
155
156        for original in test_cases {
157            let json = serde_json::to_string(&original).unwrap();
158            let deserialized: JulianDate = serde_json::from_str(&json).unwrap();
159
160            assert_eq!(
161                original.jd1(),
162                deserialized.jd1(),
163                "JD1 precision lost in serde round-trip"
164            );
165            assert_eq!(
166                original.jd2(),
167                deserialized.jd2(),
168                "JD2 precision lost in serde round-trip"
169            );
170            assert_eq!(
171                original, deserialized,
172                "JulianDate equality lost in serde round-trip"
173            );
174        }
175    }
176}