Skip to main content

celestial_time/sidereal/
observatory.rs

1use super::{GAST, GMST, LAST, LMST};
2use crate::scales::{TT, UT1};
3use crate::TimeResult;
4use celestial_core::Location;
5
6#[derive(Debug, Clone, Copy)]
7pub struct ObservatoryContext<'a> {
8    ut1: &'a UT1,
9    tt: &'a TT,
10    location: &'a Location,
11}
12
13impl<'a> ObservatoryContext<'a> {
14    pub fn new(ut1: &'a UT1, tt: &'a TT, location: &'a Location) -> Self {
15        Self { ut1, tt, location }
16    }
17
18    /// Get location for a famous observatory
19    ///
20    /// Provides quick access to well-known observatory locations.
21    /// Returns an owned Location that can be used with `new()`.
22    ///
23    /// # Supported Observatories
24    /// * `"mauna_kea"` - Mauna Kea Observatory, Hawaii
25    /// * `"greenwich"` - Royal Observatory Greenwich, UK
26    /// * `"palomar"` - Palomar Observatory, California
27    /// * `"vlt"` - Very Large Telescope, Chile
28    /// * `"keck"` - W. M. Keck Observatory, Hawaii
29    ///
30    /// # Examples
31    /// ```
32    /// use celestial_time::{UT1, TT};
33    /// use celestial_time::sidereal::ObservatoryContext;
34    ///
35    /// let ut1 = UT1::j2000();
36    /// let tt = TT::j2000();
37    /// let location = ObservatoryContext::observatory_location("mauna_kea").unwrap();
38    /// let observatory = ObservatoryContext::new(&ut1, &tt, &location);
39    /// ```
40    pub fn observatory_location(observatory_name: &str) -> TimeResult<Location> {
41        match observatory_name {
42            "mauna_kea" | "keck" => Ok(Location::from_degrees(19.8283, -155.4783, 4145.0)
43                .expect("Keck Observatory coordinates are valid")),
44            "greenwich" => Ok(Location::greenwich()),
45            "palomar" => Ok(Location::from_degrees(33.3563, -116.8650, 1712.0)
46                .expect("Palomar coordinates are valid")),
47            "vlt" => Ok(Location::from_degrees(-24.6275, -70.4044, 2635.0)
48                .expect("VLT coordinates are valid")),
49            _ => Err(crate::TimeError::CalculationError(format!(
50                "Unknown observatory: {}",
51                observatory_name
52            ))),
53        }
54    }
55
56    /// Get the UT1 time
57    pub fn ut1(&self) -> &UT1 {
58        self.ut1
59    }
60
61    /// Get the TT time
62    pub fn tt(&self) -> &TT {
63        self.tt
64    }
65
66    /// Get the observer location
67    pub fn location(&self) -> &Location {
68        self.location
69    }
70
71    /// Calculate Greenwich Mean Sidereal Time (GMST)
72    ///
73    /// # Examples
74    /// ```
75    /// use celestial_time::{UT1, TT};
76    /// use celestial_time::sidereal::ObservatoryContext;
77    /// use celestial_core::Location;
78    ///
79    /// let ut1 = UT1::j2000();
80    /// let tt = TT::j2000();
81    /// let location = Location::from_degrees(19.8283, -155.4783, 4145.0).unwrap();
82    /// let observatory = ObservatoryContext::new(&ut1, &tt, &location);
83    /// let gmst = observatory.gmst().unwrap();
84    /// ```
85    pub fn gmst(&self) -> TimeResult<GMST> {
86        GMST::from_ut1_and_tt(self.ut1, self.tt)
87    }
88
89    /// Calculate Greenwich Apparent Sidereal Time (GAST)
90    ///
91    /// # Examples
92    /// ```
93    /// use celestial_time::{UT1, TT};
94    /// use celestial_time::sidereal::ObservatoryContext;
95    /// use celestial_core::Location;
96    ///
97    /// let ut1 = UT1::j2000();
98    /// let tt = TT::j2000();
99    /// let location = Location::from_degrees(19.8283, -155.4783, 4145.0).unwrap();
100    /// let observatory = ObservatoryContext::new(&ut1, &tt, &location);
101    /// let gast = observatory.gast().unwrap();
102    /// ```
103    pub fn gast(&self) -> TimeResult<GAST> {
104        GAST::from_ut1_and_tt(self.ut1, self.tt)
105    }
106
107    /// Calculate Local Mean Sidereal Time (LMST)
108    ///
109    /// # Examples
110    /// ```
111    /// use celestial_time::{UT1, TT};
112    /// use celestial_time::sidereal::ObservatoryContext;
113    /// use celestial_core::Location;
114    ///
115    /// let ut1 = UT1::j2000();
116    /// let tt = TT::j2000();
117    /// let location = Location::from_degrees(19.8283, -155.4783, 4145.0).unwrap();
118    /// let observatory = ObservatoryContext::new(&ut1, &tt, &location);
119    /// let lmst = observatory.lmst().unwrap();
120    /// ```
121    pub fn lmst(&self) -> TimeResult<LMST> {
122        LMST::from_ut1_tt_and_location(self.ut1, self.tt, self.location)
123    }
124
125    /// Calculate Local Apparent Sidereal Time (LAST)
126    ///
127    /// # Examples
128    /// ```
129    /// use celestial_time::{UT1, TT};
130    /// use celestial_time::sidereal::ObservatoryContext;
131    /// use celestial_core::Location;
132    ///
133    /// let ut1 = UT1::j2000();
134    /// let tt = TT::j2000();
135    /// let location = Location::from_degrees(19.8283, -155.4783, 4145.0).unwrap();
136    /// let observatory = ObservatoryContext::new(&ut1, &tt, &location);
137    /// let last = observatory.last().unwrap();
138    /// ```
139    pub fn last(&self) -> TimeResult<LAST> {
140        LAST::from_ut1_tt_and_location(self.ut1, self.tt, self.location)
141    }
142
143    /// Get all sidereal times at once
144    ///
145    /// Returns a tuple of (GMST, GAST, LMST, LAST) for convenience.
146    ///
147    /// # Examples
148    /// ```
149    /// use celestial_time::{UT1, TT};
150    /// use celestial_time::sidereal::ObservatoryContext;
151    /// use celestial_core::Location;
152    ///
153    /// let ut1 = UT1::j2000();
154    /// let tt = TT::j2000();
155    /// let location = Location::from_degrees(19.8283, -155.4783, 4145.0).unwrap();
156    /// let observatory = ObservatoryContext::new(&ut1, &tt, &location);
157    /// let (gmst, gast, lmst, last) = observatory.all_sidereal_times().unwrap();
158    /// ```
159    pub fn all_sidereal_times(&self) -> TimeResult<(GMST, GAST, LMST, LAST)> {
160        let gmst = self.gmst()?;
161        let gast = self.gast()?;
162        let lmst = self.lmst()?;
163        let last = self.last()?;
164        Ok((gmst, gast, lmst, last))
165    }
166
167    /// Calculate hour angle to a target right ascension
168    ///
169    /// Uses Local Apparent Sidereal Time for the most accurate hour angle calculation.
170    ///
171    /// # Arguments
172    /// * `target_ra_hours` - Target right ascension in hours
173    ///
174    /// # Returns
175    /// Hour angle in hours, normalized to [-12, 12) range
176    ///
177    /// # Examples
178    /// ```
179    /// use celestial_time::{UT1, TT};
180    /// use celestial_time::sidereal::ObservatoryContext;
181    /// use celestial_core::Location;
182    ///
183    /// let ut1 = UT1::j2000();
184    /// let tt = TT::j2000();
185    /// let location = Location::from_degrees(19.8283, -155.4783, 4145.0).unwrap();
186    /// let observatory = ObservatoryContext::new(&ut1, &tt, &location);
187    /// let hour_angle = observatory.hour_angle_to_target(6.0).unwrap(); // RA = 6h
188    /// ```
189    pub fn hour_angle_to_target(&self, target_ra_hours: f64) -> TimeResult<f64> {
190        let last = self.last()?;
191        Ok(last.hour_angle_to_target(target_ra_hours))
192    }
193
194    /// Get observatory information as a formatted string
195    ///
196    /// Useful for logging and debugging.
197    ///
198    /// # Examples
199    /// ```
200    /// use celestial_time::{UT1, TT};
201    /// use celestial_time::sidereal::ObservatoryContext;
202    /// use celestial_core::Location;
203    ///
204    /// let ut1 = UT1::j2000();
205    /// let tt = TT::j2000();
206    /// let location = Location::from_degrees(19.8283, -155.4783, 4145.0).unwrap();
207    /// let observatory = ObservatoryContext::new(&ut1, &tt, &location);
208    /// println!("{}", observatory.info());
209    /// ```
210    pub fn info(&self) -> String {
211        let lat_deg = self.location.latitude * celestial_core::constants::RAD_TO_DEG;
212        let lon_deg = self.location.longitude * celestial_core::constants::RAD_TO_DEG;
213        let height_m = self.location.height;
214
215        format!(
216            "Observatory at ({:.4}°, {:.4}°, {:.0}m) - UT1: {}, TT: {}",
217            lat_deg,
218            lon_deg,
219            height_m,
220            self.ut1.to_julian_date().jd1() + self.ut1.to_julian_date().jd2(),
221            self.tt.to_julian_date().jd1() + self.tt.to_julian_date().jd2()
222        )
223    }
224}
225
226impl<'a> std::fmt::Display for ObservatoryContext<'a> {
227    fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result {
228        write!(f, "{}", self.info())
229    }
230}
231
232#[cfg(test)]
233mod tests {
234    use super::*;
235
236    fn mauna_kea() -> Location {
237        Location::from_degrees(19.8283, -155.4783, 4145.0).unwrap()
238    }
239
240    fn greenwich() -> Location {
241        Location::greenwich()
242    }
243
244    #[test]
245    fn test_observatory_context_creation() {
246        let ut1 = UT1::j2000();
247        let tt = TT::j2000();
248        let location = mauna_kea();
249        let observatory = ObservatoryContext::new(&ut1, &tt, &location);
250
251        let obs_ut1_jd = observatory.ut1().to_julian_date();
252        let ut1_jd = ut1.to_julian_date();
253        let obs_tt_jd = observatory.tt().to_julian_date();
254        let tt_jd = tt.to_julian_date();
255        assert_eq!(
256            obs_ut1_jd.jd1() + obs_ut1_jd.jd2(),
257            ut1_jd.jd1() + ut1_jd.jd2()
258        );
259        assert_eq!(obs_tt_jd.jd1() + obs_tt_jd.jd2(), tt_jd.jd1() + tt_jd.jd2());
260        assert_eq!(observatory.location().latitude, location.latitude);
261        assert_eq!(observatory.location().longitude, location.longitude);
262        assert_eq!(observatory.location().height, location.height);
263    }
264
265    #[test]
266    fn test_all_sidereal_times() {
267        let ut1 = UT1::j2000();
268        let tt = TT::j2000();
269        let location = mauna_kea();
270        let observatory = ObservatoryContext::new(&ut1, &tt, &location);
271
272        let (gmst, gast, lmst, last) = observatory.all_sidereal_times().unwrap();
273
274        assert!(gmst.hours() >= 0.0 && gmst.hours() < 24.0);
275        assert!(gast.hours() >= 0.0 && gast.hours() < 24.0);
276        assert!(lmst.hours() >= 0.0 && lmst.hours() < 24.0);
277        assert!(last.hours() >= 0.0 && last.hours() < 24.0);
278
279        assert!((gast.hours() - gmst.hours()).abs() < 1.0);
280
281        let expected_offset = -155.4783 / 15.0;
282        let actual_offset = lmst.hours() - gmst.hours();
283        let normalized_offset = if actual_offset > 12.0 {
284            actual_offset - 24.0
285        } else if actual_offset < -12.0 {
286            actual_offset + 24.0
287        } else {
288            actual_offset
289        };
290        assert!((normalized_offset - expected_offset).abs() < 1e-10);
291
292        let last_offset = last.hours() - gast.hours();
293        let normalized_last_offset = if last_offset > 12.0 {
294            last_offset - 24.0
295        } else if last_offset < -12.0 {
296            last_offset + 24.0
297        } else {
298            last_offset
299        };
300        assert!((normalized_last_offset - expected_offset).abs() < 1e-10);
301    }
302
303    #[test]
304    fn test_hour_angle_calculation() {
305        let ut1 = UT1::j2000();
306        let tt = TT::j2000();
307        let location = mauna_kea();
308        let observatory = ObservatoryContext::new(&ut1, &tt, &location);
309
310        let target_ra = 6.0;
311        let hour_angle = observatory.hour_angle_to_target(target_ra).unwrap();
312
313        let last = observatory.last().unwrap();
314        let expected_ha = last.hour_angle_to_target(target_ra);
315        assert!((hour_angle - expected_ha).abs() < 1e-12);
316    }
317
318    #[test]
319    fn test_famous_observatories() {
320        let ut1 = UT1::j2000();
321        let tt = TT::j2000();
322
323        let observatories = ["mauna_kea", "greenwich", "palomar", "vlt", "keck"];
324
325        for name in observatories {
326            let location = ObservatoryContext::observatory_location(name).unwrap();
327            let observatory = ObservatoryContext::new(&ut1, &tt, &location);
328            let (gmst, gast, lmst, last) = observatory.all_sidereal_times().unwrap();
329
330            assert!(
331                gmst.hours() >= 0.0 && gmst.hours() < 24.0,
332                "Invalid GMST for {}",
333                name
334            );
335            assert!(
336                gast.hours() >= 0.0 && gast.hours() < 24.0,
337                "Invalid GAST for {}",
338                name
339            );
340            assert!(
341                lmst.hours() >= 0.0 && lmst.hours() < 24.0,
342                "Invalid LMST for {}",
343                name
344            );
345            assert!(
346                last.hours() >= 0.0 && last.hours() < 24.0,
347                "Invalid LAST for {}",
348                name
349            );
350        }
351    }
352
353    #[test]
354    fn test_unknown_observatory() {
355        let result = ObservatoryContext::observatory_location("unknown_observatory");
356        assert!(result.is_err());
357    }
358
359    #[test]
360    fn test_individual_sidereal_calculations() {
361        let ut1 = UT1::j2000();
362        let tt = TT::j2000();
363        let location = mauna_kea();
364        let observatory = ObservatoryContext::new(&ut1, &tt, &location);
365
366        let (gmst_batch, gast_batch, lmst_batch, last_batch) =
367            observatory.all_sidereal_times().unwrap();
368
369        let gmst_individual = observatory.gmst().unwrap();
370        let gast_individual = observatory.gast().unwrap();
371        let lmst_individual = observatory.lmst().unwrap();
372        let last_individual = observatory.last().unwrap();
373
374        assert!((gmst_individual.hours() - gmst_batch.hours()).abs() < 1e-12);
375        assert!((gast_individual.hours() - gast_batch.hours()).abs() < 1e-12);
376        assert!((lmst_individual.hours() - lmst_batch.hours()).abs() < 1e-12);
377        assert!((last_individual.hours() - last_batch.hours()).abs() < 1e-12);
378    }
379
380    #[test]
381    fn test_display_and_info() {
382        let ut1 = UT1::j2000();
383        let tt = TT::j2000();
384        let location = mauna_kea();
385        let observatory = ObservatoryContext::new(&ut1, &tt, &location);
386
387        let info = observatory.info();
388        let display = format!("{}", observatory);
389
390        assert!(info.contains("19.8283"));
391        assert!(info.contains("-155.4783"));
392        assert!(info.contains("4145"));
393
394        assert_eq!(info, display);
395    }
396
397    #[test]
398    fn test_observatory_context_hour_angle() {
399        let ut1 = UT1::j2000();
400        let tt = TT::j2000();
401        let location = greenwich();
402        let observatory = ObservatoryContext::new(&ut1, &tt, &location);
403
404        let target_ra = 12.0;
405        let hour_angle = observatory.hour_angle_to_target(target_ra).unwrap();
406
407        let last = observatory.last().unwrap();
408        let expected_ha = last.hour_angle_to_target(target_ra);
409        assert_eq!(hour_angle, expected_ha);
410
411        assert!((-12.0..12.0).contains(&hour_angle));
412    }
413}