gistools/space/util/
time.rs

1use crate::util::Date;
2use core::f64::consts::TAU;
3use libm::floor;
4
5/// Time object to track month-day-hour-minute-second
6#[derive(Debug, Default, Clone, Copy, PartialEq)]
7pub struct TimeStamp {
8    /// month
9    pub mon: f64,
10    /// day
11    pub day: f64,
12    /// hour
13    pub hr: f64,
14    /// minute
15    pub min: f64,
16    /// second
17    pub sec: f64,
18}
19
20/// procedure days2mdhms
21///
22/// this procedure converts the day of the year, days, to the equivalent month
23/// day, hour, minute and second.
24///
25/// algorithm     set up array for the number of days per month
26/// find leap year - use 1900 because 2000 is a leap year
27/// loop through a temp value while the value is < the days
28/// perform int conversions to the correct day and month
29/// convert remainder into h m s using type conversions
30///
31/// author        david vallado                  719-573-2600    1 mar 2001
32///
33/// ## Parameters
34/// - `year`: year to date
35/// - `days`: day of year
36///
37/// ## Returns
38/// Decomposed information into year, month, day, hour, minute and second
39pub fn days2mdhms(year: u16, days: f64) -> TimeStamp {
40    let lmonth = [
41        31.,
42        if year.is_multiple_of(4) { 29. } else { 28. },
43        31.,
44        30.,
45        31.,
46        30.,
47        31.,
48        31.,
49        30.,
50        31.,
51        30.,
52        31.,
53    ];
54    let dayofyr = floor(days);
55
56    //  ----------------- find month and day of month ----------------
57    let mut i = 1;
58    let mut inttemp = 0.;
59    while dayofyr > (inttemp + lmonth[i - 1]) && i < 12 {
60        inttemp += lmonth[i - 1];
61        i += 1;
62    }
63
64    let mon = i;
65    let day = dayofyr - inttemp;
66
67    //  ----------------- find hours minutes and seconds -------------
68    let mut temp = (days - dayofyr) * 24.0;
69    let hr = floor(temp);
70    temp = (temp - hr) * 60.0;
71    let min = floor(temp);
72    let sec = (temp - min) * 60.0;
73
74    TimeStamp { mon: mon as f64, day, hr, min, sec }
75}
76
77/// procedure jday
78///
79/// this procedure finds the julian date given the year, month, day, and time.
80/// the julian date is defined by each elapsed day since noon, jan 1, 4713 bc.
81///
82/// algorithm     calculate the answer in one step for efficiency
83///
84/// author         david vallado                  719-573-2600    1 mar 2001
85///
86/// ## Parameters
87/// - `year`: year
88/// - `mon`: month
89/// - `day`: day
90/// - `hr`: hour
91/// - `min`: minute
92/// - `sec`: seconds
93/// - `msec`: milliseconds
94///
95/// ## Returns
96/// Julian date
97pub fn jday_internal(
98    year: f64,
99    mon: f64,
100    day: f64,
101    hr: f64,
102    min: f64,
103    sec: f64,
104    msec: Option<f64>,
105) -> f64 {
106    let msec = msec.unwrap_or(0.);
107    367.0 * year - floor(7. * (year + floor((mon + 9.) / 12.0)) * 0.25)
108        + floor((275. * mon) / 9.0)
109        + day
110        + 1721013.5
111        + ((msec / 60000. + sec / 60.0 + min) / 60.0 + hr) / 24.0
112}
113
114/// Builds a julian date from the given parameters
115///
116/// ## Parameters
117/// - `date`: date
118///
119/// ## Returns
120/// Julian date
121pub fn jday(date: &Date) -> f64 {
122    jday_internal(
123        date.year as f64,
124        date.month as f64,
125        date.day as f64,
126        date.hour as f64,
127        date.minute as f64,
128        date.second as f64,
129        None,
130    )
131}
132
133/// function gstime
134///
135/// this function finds the greenwich sidereal time.
136///
137/// author         david vallado                  719-573-2600    1 mar 2001
138///
139/// ## Parameters
140/// - `jdut1`: julian date
141///
142/// ## Returns
143/// Greenwich sidereal time
144fn gstime_internal(jdut1: f64) -> f64 {
145    let tut1 = (jdut1 - 2451545.0) / 36525.0;
146
147    let mut temp = -6.2e-6 * tut1 * tut1 * tut1
148        + 0.093104 * tut1 * tut1
149        + (876600.0 * 3600.0 + 8640184.812866) * tut1
150        + 67310.54841; // # sec
151    temp = ((temp.to_radians()) / 240.0) % TAU; // 360/86400 = 1/240, to deg, to rad
152
153    //  ------------------------ check quadrants ---------------------
154    if temp < 0.0 {
155        temp += TAU;
156    }
157
158    temp
159}
160
161/// find greenwich sidereal time
162///
163/// ## Parameters
164/// - `time`: julian date
165///
166/// ## Returns
167/// Greenwich sidereal time
168pub fn gstime(time: &Date) -> f64 {
169    gstime_internal(jday(time))
170}