deep_time/eop/
mod.rs

1//! Earth Orientation Parameters (EOP) data parser and interpolator.
2//!
3//! Loads UT1–UTC offsets and polar motion (x, y) from standard IERS formats
4//! (`Finals2000A`, `C04`) or custom column layouts.
5//!
6//! Provides interpolation at any MJD and integrates with [`Dt`]
7//! for UT1 time scale conversions.
8//!
9//! Tries to be planet/body agnostic, such that custom data for any world
10//! might be able to be loaded and used.
11
12#![allow(clippy::indexing_slicing)]
13#![allow(clippy::excessive_precision)]
14#![allow(clippy::approx_constant)]
15#![allow(clippy::eq_op)]
16
17use crate::{Dt, DtErr, DtErrKind, Real, Scale, an_err};
18use alloc::string::String;
19use alloc::vec::Vec;
20use core::cmp::Ordering;
21
22#[derive(Debug, Clone, Copy, Default)]
23pub enum Separator {
24    #[default]
25    Whitespace,
26    Comma,
27    Tab,
28    Pipe,
29    Semicolon,
30}
31
32/// Earth/Body Orientation Parameters Format.
33///
34/// Formats to provide to the parser, including a
35/// custom one to allow specific column indices.
36///
37/// - `Finals2000A` such as is available from
38///   https://maia.usno.navy.mil/ser7/finals2000A.all
39/// - `C04` such as is available from
40///   https://datacenter.iers.org/data/latestVersion/EOP_20u24_C04_one_file_1962-now.txt
41/// - `Custom` so you can provide your own specific column indices
42///   using [`CustomEopCols`].
43#[derive(Debug, Clone, Copy, Default)]
44pub enum EopFormat {
45    /// finals2000A.all / finals.all.iau2000.txt style files
46    #[default]
47    Finals2000A,
48    /// C04 long-term series
49    C04,
50    /// User-defined column indices (0-based)
51    Custom(CustomEopCols),
52}
53
54/// For use with [`EopFormat::Custom`].
55///
56/// Allows you to specify exactly which 0-based column contains each value
57/// when your input file does not match a standard IERS layout.
58#[derive(Debug, Clone, Copy)]
59pub struct CustomEopCols {
60    pub mjd: usize,
61    pub offset: usize,
62    pub pm_x: Option<usize>,
63    pub pm_y: Option<usize>,
64}
65
66/// A single parsed row of Earth Orientation Parameters.
67///
68/// - `mjd` — Modified Julian Date
69/// - `offset` — UT1 − UTC (or equivalent) in **seconds**
70/// - `pm_x`, `pm_y` — Polar motion in **arcseconds**
71#[derive(Debug, Clone, Copy)]
72pub struct EopDataRow {
73    pub mjd: Real,
74    /// e.g. UT1-UTC(s)
75    pub offset: Real,
76    /// polar motion x (arcsec)
77    pub pm_x: Real,
78    /// polar motion y (arcsec)
79    pub pm_y: Real,
80}
81
82/// Container for Earth/Body Orientation Parameters data.
83///
84/// - On Earth this would enable time scale conversions to and from
85///   the **UT1 time scale**.
86/// - Earth Orientation Parameters data is available from: https://maia.usno.navy.mil/ser7/finals2000A.all
87#[derive(Debug, Clone)]
88pub struct EopData {
89    rows: Vec<EopDataRow>,
90}
91
92#[cfg(feature = "std")]
93impl EopData {
94    /// Parse EOP data from any `std::io::BufRead` (file, network stream, etc.).
95    ///
96    /// Lines starting with `#` or longer than [`MAX_LINE_LEN`] are skipped.
97    /// The returned vector is always sorted by MJD.
98    pub fn data_from_reader<R: std::io::BufRead>(
99        mut reader: R,
100        format: EopFormat,
101        separator: Separator,
102    ) -> Result<Vec<EopDataRow>, DtErr> {
103        let mut line_buf = String::with_capacity(256);
104        let mut rows = Vec::new();
105
106        loop {
107            line_buf.clear();
108
109            let bytes_read = match reader.read_line(&mut line_buf) {
110                Ok(0) => break,
111                Ok(n) => n,
112                Err(e) => {
113                    return Err(an_err!(DtErrKind::IOErr, "read line: {}", e));
114                }
115            };
116
117            if bytes_read > Self::MAX_LINE_LEN {
118                continue;
119            }
120
121            let trimmed = line_buf.trim();
122            if trimmed.is_empty() || trimmed.starts_with('#') {
123                continue;
124            }
125
126            if let Some(row) = Self::try_parse_row(trimmed, format, separator) {
127                rows.push(row);
128            }
129        }
130
131        if rows.is_empty() {
132            return Err(an_err!(DtErrKind::Incomplete, "no valid rows"));
133        }
134
135        rows.sort_by(|a, b| a.mjd.partial_cmp(&b.mjd).unwrap_or(Ordering::Equal));
136        Ok(rows)
137    }
138
139    /// Returns a [`Vec`] of [`EopDataRow`] from a text file on disk.
140    ///
141    /// ## Examples
142    ///
143    /// ```rust
144    /// # #[cfg(feature = "eop-tests")]
145    /// # {
146    /// use deep_time::eop::{EopData, EopFormat, Separator};
147    ///
148    /// let path = "finals.all.iau2000.txt";
149    /// let rows = EopData::data_from_text_file(path, EopFormat::Finals2000A, Separator::Whitespace).unwrap();
150    /// # }
151    /// ```
152    ///
153    /// ## See also
154    ///
155    /// - [`EopData::from_text_file`](../eop/struct.EopData.html#method.from_text_file)
156    pub fn data_from_text_file<P: AsRef<std::path::Path>>(
157        path: P,
158        format: EopFormat,
159        separator: Separator,
160    ) -> Result<Vec<EopDataRow>, DtErr> {
161        use std::fs::File;
162        use std::io::BufReader;
163
164        let path = path.as_ref();
165        let file = File::open(path)
166            .map_err(|e| an_err!(DtErrKind::IOErr, "open file: '{}': {}", path.display(), e))?;
167
168        let reader = BufReader::new(file);
169        Self::data_from_reader(reader, format, separator)
170    }
171
172    /// Create an [`EopData`] by loading from a text file on disk.
173    ///
174    /// ## Examples
175    ///
176    /// ```rust
177    /// # #[cfg(feature = "eop-tests")]
178    /// # {
179    /// use deep_time::eop::{EopData, EopFormat, Separator};
180    ///
181    /// let path = "finals.all.iau2000.txt";
182    /// let provider = EopData::from_text_file(path, EopFormat::Finals2000A, Separator::Whitespace).unwrap();
183    /// # }
184    /// ```
185    pub fn from_text_file<P: AsRef<std::path::Path>>(
186        path: P,
187        format: EopFormat,
188        separator: Separator,
189    ) -> Result<Self, DtErr> {
190        let rows = Self::data_from_text_file(path, format, separator)?;
191        Ok(Self { rows })
192    }
193}
194
195impl EopData {
196    pub const MAX_LINE_LEN: usize = 8192;
197
198    // Small helper — parses ONE row (shared by all paths)
199    fn try_parse_row(trimmed: &str, format: EopFormat, separator: Separator) -> Option<EopDataRow> {
200        let parts: Vec<&str> = match separator {
201            Separator::Whitespace => trimmed.split_whitespace().collect(),
202            Separator::Comma => trimmed.split(',').map(|s| s.trim()).collect(),
203            Separator::Tab => trimmed.split('\t').map(|s| s.trim()).collect(),
204            Separator::Pipe => trimmed.split('|').map(|s| s.trim()).collect(),
205            Separator::Semicolon => trimmed.split(';').map(|s| s.trim()).collect(),
206        };
207
208        if parts.len() < 2 {
209            return None;
210        }
211
212        let (mjd, offset, pm_x, pm_y) = match format {
213            EopFormat::Finals2000A => {
214                let mjd_idx = parts.iter().position(|p| {
215                    p.contains('.') && p.parse::<Real>().is_ok_and(|v| v > 30000.0)
216                })?;
217
218                let mut flag_count = 0;
219                let mut offset_value: Option<Real> = None;
220                let mut pm_x_val: Real = 0.0;
221                let mut pm_y_val: Real = 0.0;
222
223                for i in (mjd_idx + 1)..parts.len() {
224                    let token = parts[i];
225
226                    let is_flag = token == "I"
227                        || token == "P"
228                        || token.starts_with("I-")
229                        || token.starts_with("P-");
230
231                    if is_flag {
232                        flag_count += 1;
233
234                        if flag_count == 1 {
235                            // After first flag: x_p and y_p
236                            if let (Some(x_str), Some(y_str)) = (parts.get(i + 1), parts.get(i + 3))
237                            {
238                                pm_x_val = x_str.parse::<Real>().unwrap_or(0.0);
239                                pm_y_val = y_str.parse::<Real>().unwrap_or(0.0);
240                            }
241                        }
242
243                        if flag_count == 2 {
244                            let value_str = if token.starts_with("I-") || token.starts_with("P-") {
245                                &token[1..]
246                            } else if i + 1 < parts.len() {
247                                parts[i + 1]
248                            } else {
249                                break;
250                            };
251                            if let Ok(val) = value_str.parse::<Real>() {
252                                offset_value = Some(val);
253                            }
254                            break;
255                        }
256                    }
257                }
258
259                let offset = offset_value?;
260                let mjd_val = parts[mjd_idx].parse::<Real>().ok()?;
261
262                (mjd_val, offset, pm_x_val, pm_y_val)
263            }
264
265            EopFormat::C04 => {
266                let mjd = parts.get(4)?.parse::<Real>().ok()?;
267                let pm_x = parts
268                    .get(5)
269                    .unwrap_or(&"0.0")
270                    .parse::<Real>()
271                    .unwrap_or(0.0);
272                let pm_y = parts
273                    .get(6)
274                    .unwrap_or(&"0.0")
275                    .parse::<Real>()
276                    .unwrap_or(0.0);
277                let offset = parts.get(7)?.parse::<Real>().ok()?;
278                (mjd, offset, pm_x, pm_y)
279            }
280
281            EopFormat::Custom(cols) => {
282                let mjd = parts.get(cols.mjd)?.parse::<Real>().ok()?;
283                let offset = parts.get(cols.offset)?.parse::<Real>().ok()?;
284                let pm_x = if let Some(pm_x_col) = cols.pm_x {
285                    parts
286                        .get(pm_x_col)
287                        .unwrap_or(&"0.0")
288                        .parse::<Real>()
289                        .ok()
290                        .unwrap_or(0.0)
291                } else {
292                    0.0
293                };
294                let pm_y = if let Some(pm_y_col) = cols.pm_y {
295                    parts
296                        .get(pm_y_col)
297                        .unwrap_or(&"0.0")
298                        .parse::<Real>()
299                        .ok()
300                        .unwrap_or(0.0)
301                } else {
302                    0.0
303                };
304                (mjd, offset, pm_x, pm_y)
305            }
306        };
307
308        Some(EopDataRow {
309            mjd,
310            offset,
311            pm_x,
312            pm_y,
313        })
314    }
315
316    fn parse_lines<'a>(
317        lines: impl Iterator<Item = &'a str>,
318        format: EopFormat,
319        separator: Separator,
320    ) -> Result<Vec<EopDataRow>, DtErr> {
321        let mut rows = Vec::new();
322
323        for line in lines {
324            let trimmed = line.trim();
325            if trimmed.is_empty() || trimmed.starts_with('#') || trimmed.len() > Self::MAX_LINE_LEN
326            {
327                continue;
328            }
329
330            if let Some(row) = Self::try_parse_row(trimmed, format, separator) {
331                rows.push(row);
332            }
333        }
334
335        if rows.is_empty() {
336            return Err(an_err!(DtErrKind::Incomplete, "no valid rows"));
337        }
338
339        rows.sort_by(|a, b| a.mjd.partial_cmp(&b.mjd).unwrap_or(Ordering::Equal));
340        Ok(rows)
341    }
342
343    /// Parse EOP data from a `&str`.
344    ///
345    /// Useful when the data is already in memory (embedded resource,
346    /// downloaded string, etc.).
347    pub fn data_from_str(
348        s: &str,
349        format: EopFormat,
350        separator: Separator,
351    ) -> Result<Vec<EopDataRow>, DtErr> {
352        Self::parse_lines(s.lines(), format, separator)
353    }
354
355    /// Parse EOP data from raw bytes.
356    ///
357    /// The bytes are interpreted as UTF-8. Invalid UTF-8 sequences
358    /// result in an empty string (and therefore an error).
359    pub fn data_from_bytes(
360        bytes: &[u8],
361        format: EopFormat,
362        separator: Separator,
363    ) -> Result<Vec<EopDataRow>, DtErr> {
364        let s = core::str::from_utf8(bytes).unwrap_or("");
365        Self::data_from_str(s, format, separator)
366    }
367
368    /// Create an [`EopData`] from a string slice.
369    pub fn from_str(s: &str, format: EopFormat, separator: Separator) -> Result<Self, DtErr> {
370        let rows = Self::data_from_str(s, format, separator)?;
371        Ok(Self { rows })
372    }
373
374    /// Create an [`EopData`] from raw bytes.
375    pub fn from_bytes(
376        bytes: &[u8],
377        format: EopFormat,
378        separator: Separator,
379    ) -> Result<Self, DtErr> {
380        let rows = Self::data_from_bytes(bytes, format, separator)?;
381        Ok(Self { rows })
382    }
383
384    /// Returns all interpolated orientation parameters (offset + polar motion)
385    /// at the given MJD.
386    ///
387    /// Returns `None` if the table is empty or the MJD is completely outside
388    /// the loaded data.
389    pub fn eop_offset(&self, mjd: Real) -> Option<EopOffset> {
390        if self.rows.is_empty() {
391            return None;
392        }
393
394        let idx = match self
395            .rows
396            .binary_search_by(|probe| probe.mjd.partial_cmp(&mjd).unwrap_or(Ordering::Equal))
397        {
398            Ok(i) => i,
399            Err(i) => {
400                if i == 0 {
401                    let row = &self.rows[0];
402                    return Some(EopOffset {
403                        offset: row.offset,
404                        pm_x: row.pm_x,
405                        pm_y: row.pm_y,
406                    });
407                }
408                if i >= self.rows.len() {
409                    let row = &self.rows[self.rows.len() - 1];
410                    return Some(EopOffset {
411                        offset: row.offset,
412                        pm_x: row.pm_x,
413                        pm_y: row.pm_y,
414                    });
415                }
416                i - 1
417            }
418        };
419
420        if idx + 1 < self.rows.len() {
421            let e0 = &self.rows[idx];
422            let e1 = &self.rows[idx + 1];
423
424            let t = (mjd - e0.mjd) / (e1.mjd - e0.mjd);
425
426            let offset = e0.offset + t * (e1.offset - e0.offset);
427            let pm_x = e0.pm_x + t * (e1.pm_x - e0.pm_x);
428            let pm_y = e0.pm_y + t * (e1.pm_y - e0.pm_y);
429
430            Some(EopOffset { offset, pm_x, pm_y })
431        } else {
432            let row = &self.rows[idx];
433            Some(EopOffset {
434                offset: row.offset,
435                pm_x: row.pm_x,
436                pm_y: row.pm_y,
437            })
438        }
439    }
440}
441
442/// Interpolated Body Orientation Parameters at a specific MJD.
443///
444/// Contains everything needed for high-precision sidereal time
445/// and polar-motion corrections.
446#[derive(Debug, Clone, Copy, Default)]
447pub struct EopOffset {
448    /// Value in **seconds** e.g. UT1 − UTC offset
449    pub offset: Real,
450    /// Polar motion x-coordinate in **arcseconds**
451    pub pm_x: Real,
452    /// Polar motion y-coordinate in **arcseconds**
453    pub pm_y: Real,
454}
455
456impl Dt {
457    /// Get an orientation parameters offset in seconds inside a struct: ([`EopOffset`])
458    /// for a particular Modified Julian Date.
459    ///
460    /// - On Earth this would be the UT1 time scale.
461    /// - Earth Orientation Parameters data is available from: https://maia.usno.navy.mil/ser7/finals2000A.all
462    pub fn mjd_to_eop_offset(mjd: Real, op_data: &EopData) -> Result<EopOffset, DtErr> {
463        let offset = op_data
464            .eop_offset(mjd)
465            .ok_or_else(|| an_err!(DtErrKind::OutOfRange, "mjd: {mjd}"))?;
466        Ok(offset)
467    }
468
469    /// Get an orientation parameters offset in seconds for a particular Modified Julian Date.
470    ///
471    /// - On Earth this would be the UT1 time scale.
472    /// - Earth Orientation Parameters data is available from: https://maia.usno.navy.mil/ser7/finals2000A.all
473    #[inline]
474    pub fn mjd_to_eop_offset_f(mjd: Real, op_data: &EopData) -> Result<Real, DtErr> {
475        Self::mjd_to_eop_offset(mjd, op_data).map(|res| res.offset)
476    }
477
478    /// Offsets a [`Dt`] using orientation parameters data.
479    ///
480    /// - On Earth this would be the UT1 time scale.
481    /// - Earth Orientation Parameters data is available from: https://maia.usno.navy.mil/ser7/finals2000A.all
482    #[inline]
483    pub fn to_eop(&self, op_data: &EopData) -> Result<Self, DtErr> {
484        Ok(self.add(Dt::from_sec_f(
485            Self::mjd_to_eop_offset_f(self.to_mjd_f(), op_data)?,
486            Scale::TAI,
487        )))
488    }
489
490    /// Convert a [`Dt`] already offset using orientation parameters data back to whatever
491    /// it was before.
492    ///
493    /// - On Earth this would be the UT1 time scale.
494    /// - Earth Orientation Parameters data is available from: https://maia.usno.navy.mil/ser7/finals2000A.all
495    pub fn from_eop(&self, op_data: &EopData) -> Result<Self, DtErr> {
496        if op_data.rows.is_empty() {
497            return Err(an_err!(DtErrKind::InternalErr, "contains no data"));
498        }
499        let mut guess = *self;
500
501        for _ in 0..8 {
502            let mjd = guess.to_mjd_f();
503            let offset = op_data
504                .eop_offset(mjd)
505                .ok_or_else(|| an_err!(DtErrKind::OutOfRange, "mjd: {mjd}"))?
506                .offset;
507
508            guess = self.sub(Dt::from_sec_f(offset, Scale::TAI)); // TODO: guess or self?
509        }
510
511        Ok(guess)
512    }
513}
deep_time/eop/mod.rs

deep_time/eop/
mod.rs
