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

deep_time/eop/
mod.rs
