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

deep_time/eop/
mod.rs
