1use std::cmp::Ordering;
8use std::collections::BTreeMap;
9use std::fmt::{self, Write as _};
10
11use crate::astro::constants::time::SECONDS_PER_DAY_I64;
12use crate::astro::math::interp::lerp_ratio;
13use crate::astro::time::civil::{
14 civil_from_julian_day_number, j2000_seconds_from_split, seconds_between_splits,
15 J2000_JULIAN_DAY_NUMBER, J2000_NOON_OFFSET_S,
16};
17use crate::astro::time::model::{Instant, InstantRepr, JulianDateSplit, TimeScale};
18use crate::astro::time::scales::julian_day_number;
19use crate::constants::{
20 GPS_EPOCH_TO_J2000_S, J2000_JD, MICROSECONDS_PER_SECOND, SECONDS_PER_DAY, SECONDS_PER_HOUR,
21};
22use crate::validate::{self, FieldError};
23
24const INSTANT_SCALE_ORDER_STRIDE_S: f64 = 1.0e15;
25
26#[derive(Debug, Clone, Copy, PartialEq)]
28pub struct ClockPoint {
29 pub epoch: Instant,
31 pub bias_s: f64,
33}
34
35impl ClockPoint {
36 pub fn gps_seconds(&self) -> Option<f64> {
38 instant_to_gps_seconds(&self.epoch)
39 }
40}
41
42#[derive(Debug, Clone, Copy, PartialEq)]
44pub struct ClockEpoch {
45 pub year: i32,
47 pub month: u8,
49 pub day: u8,
51 pub hour: u8,
53 pub minute: u8,
55 pub second: f64,
57}
58
59#[derive(Debug, Clone, PartialEq)]
61pub struct RinexClock {
62 pub time_scale: TimeScale,
64 pub series: BTreeMap<String, Vec<ClockPoint>>,
66}
67
68#[derive(Debug, Clone, PartialEq, Eq)]
70pub enum RinexClockError {
71 MalformedAsRecord {
73 line: usize,
75 reason: &'static str,
77 record: String,
79 },
80 BadField {
82 line: usize,
84 field: &'static str,
86 value: String,
88 },
89 InvalidInput {
91 field: &'static str,
93 reason: &'static str,
95 },
96 UnsupportedTimeScale {
98 scale: TimeScale,
100 },
101}
102
103impl fmt::Display for RinexClockError {
104 fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
105 match self {
106 RinexClockError::MalformedAsRecord {
107 line,
108 reason,
109 record,
110 } => write!(
111 f,
112 "malformed RINEX AS clock record at line {line}: {reason}: {record}"
113 ),
114 RinexClockError::BadField { line, field, value } => write!(
115 f,
116 "bad RINEX AS clock field at line {line}: {field}={value}"
117 ),
118 RinexClockError::InvalidInput { field, reason } => {
119 write!(f, "invalid RINEX clock input {field}: {reason}")
120 }
121 RinexClockError::UnsupportedTimeScale { scale } => {
122 write!(f, "unsupported RINEX clock time scale {}", scale.abbrev())
123 }
124 }
125 }
126}
127
128impl std::error::Error for RinexClockError {}
129
130impl RinexClock {
131 pub fn parse(text: &str) -> Result<Self, RinexClockError> {
133 let time_scale = parse_time_scale(text)?;
134 let lines = data_lines(text);
135 let mut by_sat = BTreeMap::<String, Vec<(ClockPoint, usize)>>::new();
136
137 for (line_number, line) in lines {
138 if let Some((sat, point)) = parse_record(line_number, line, time_scale)? {
139 by_sat.entry(sat).or_default().push((point, line_number));
140 }
141 }
142
143 Ok(Self {
144 time_scale,
145 series: build_series(by_sat),
146 })
147 }
148
149 pub fn parse_lossy(text: &str) -> Self {
151 let time_scale = parse_time_scale(text).unwrap_or(TimeScale::Gpst);
152 let lines = data_lines(text);
153 let mut by_sat = BTreeMap::<String, Vec<(ClockPoint, usize)>>::new();
154
155 for (line_number, line) in lines {
156 if let Ok(Some((sat, point))) = parse_record(line_number, line, time_scale) {
157 by_sat.entry(sat).or_default().push((point, line_number));
158 }
159 }
160
161 Self {
162 time_scale,
163 series: build_series(by_sat),
164 }
165 }
166
167 pub fn from_series_rows(rows: Vec<(String, Vec<(f64, f64)>)>) -> Result<Self, RinexClockError> {
169 let rows = rows
170 .into_iter()
171 .map(|(sat, points)| {
172 validate::require_strictly_increasing(
173 points.iter().map(|&(gps_seconds, _)| gps_seconds),
174 "gps_seconds",
175 )
176 .map_err(map_manual_order_error)?;
177 let points = points
178 .into_iter()
179 .map(|(gps_seconds, bias_s)| {
180 validate_finite(bias_s, "bias_s")?;
181 Ok((gps_seconds_to_instant(gps_seconds), bias_s))
182 })
183 .collect::<Result<Vec<_>, RinexClockError>>()?;
184 Ok((sat, points))
185 })
186 .collect::<Result<Vec<_>, RinexClockError>>()?;
187 Self::from_instant_series_rows(TimeScale::Gpst, rows)
188 }
189
190 pub fn from_instant_series_rows(
192 time_scale: TimeScale,
193 rows: Vec<(String, Vec<(Instant, f64)>)>,
194 ) -> Result<Self, RinexClockError> {
195 let mut series = BTreeMap::new();
196 for (sat, points) in rows {
197 let mut indexed = points
198 .into_iter()
199 .enumerate()
200 .map(|(idx, (epoch, bias_s))| {
201 let point = ClockPoint { epoch, bias_s };
202 validate_clock_point(point)?;
203 Ok((point, idx))
204 })
205 .collect::<Result<Vec<_>, RinexClockError>>()?;
206 validate_instant_series_order(&indexed)?;
207 indexed.sort_by(|(a, ai), (b, bi)| {
208 compare_instants(&a.epoch, &b.epoch).then_with(|| ai.cmp(bi))
209 });
210 series.insert(sat, dedup_by_time(indexed));
211 }
212 Ok(Self { time_scale, series })
213 }
214
215 pub fn series_rows(&self) -> Vec<(String, Vec<(f64, f64)>)> {
219 self.series
220 .iter()
221 .map(|(sat, points)| {
222 (
223 sat.clone(),
224 points
225 .iter()
226 .filter_map(|point| Some((point.gps_seconds()?, point.bias_s)))
227 .collect(),
228 )
229 })
230 .collect()
231 }
232
233 pub fn instant_series_rows(&self) -> Vec<(String, Vec<(Instant, f64)>)> {
235 self.series
236 .iter()
237 .map(|(sat, points)| {
238 (
239 sat.clone(),
240 points
241 .iter()
242 .map(|point| (point.epoch, point.bias_s))
243 .collect(),
244 )
245 })
246 .collect()
247 }
248
249 pub fn clock_s(
251 &self,
252 satellite_id: &str,
253 epoch: ClockEpoch,
254 ) -> Result<Option<f64>, RinexClockError> {
255 let epoch = civil_to_clock_instant(
256 self.time_scale,
257 epoch.year,
258 epoch.month,
259 epoch.day,
260 epoch.hour,
261 epoch.minute,
262 epoch.second,
263 )
264 .ok_or_else(|| invalid_input("epoch", "invalid civil clock epoch"))?;
265 self.clock_s_at_instant(satellite_id, epoch)
266 }
267
268 pub fn clock_s_at_instant(
270 &self,
271 satellite_id: &str,
272 epoch: Instant,
273 ) -> Result<Option<f64>, RinexClockError> {
274 validate_instant(epoch, "epoch")?;
275 let Some(records) = self.series.get(satellite_id) else {
276 return Ok(None);
277 };
278 Ok(interpolate(records, epoch))
279 }
280
281 pub fn clock_s_at_gps_seconds(
283 &self,
284 satellite_id: &str,
285 gps_seconds: f64,
286 ) -> Result<Option<f64>, RinexClockError> {
287 validate_finite(gps_seconds, "gps_seconds")?;
288 self.clock_s_at_instant(satellite_id, gps_seconds_to_instant(gps_seconds))
289 }
290
291 pub fn to_rinex_string(&self) -> Result<String, RinexClockError> {
302 let mut out = String::new();
303 let label = crate::rinex_common::time_scale_rinex_label(self.time_scale).ok_or(
304 RinexClockError::UnsupportedTimeScale {
305 scale: self.time_scale,
306 },
307 )?;
308 let _ = writeln!(out, "{:<60}RINEX VERSION / TYPE", " 3.00 C");
309 let _ = writeln!(out, "{label:<60}TIME SYSTEM ID");
310 let _ = writeln!(out, "{:<60}END OF HEADER", "");
311 for (satellite, points) in &self.series {
312 for point in points {
313 validate_serializable_clock_point(self.time_scale, point)?;
314 write_as_record(&mut out, satellite, point);
315 }
316 }
317 Ok(out)
318 }
319}
320
321fn write_as_record(out: &mut String, satellite: &str, point: &ClockPoint) {
323 let (year, month, day, hour, minute, second_us) = instant_civil_microsecond(&point.epoch);
324 let second = second_us / 1_000_000;
325 let microsecond = second_us % 1_000_000;
326 let _ = writeln!(
329 out,
330 "AS {satellite:<3} {year:04} {month:02} {day:02} {hour:02} {minute:02} {second:2}.{microsecond:06} 1 {bias}",
331 bias = point.bias_s,
332 );
333}
334
335fn instant_civil_microsecond(epoch: &Instant) -> (i64, i64, i64, i64, i64, i64) {
342 let (day_number, total_us) = match epoch.repr {
343 InstantRepr::JulianDate(split) => {
344 if (-1.0 / SECONDS_PER_DAY..0.0).contains(&split.fraction) {
351 return leap_second_civil(split);
352 }
353 let day_number = (split.jd_whole + 0.5).round() as i64;
359 let total_us =
360 (split.fraction * SECONDS_PER_DAY * MICROSECONDS_PER_SECOND).round() as i64;
361 (day_number, total_us)
362 }
363 InstantRepr::Nanos(nanos) => nanos_civil_day_microsecond(nanos),
367 };
368 let (year, month, day) = civil_from_julian_day_number(day_number);
369 let hour = total_us / 3_600_000_000;
370 let rem = total_us % 3_600_000_000;
371 let minute = rem / 60_000_000;
372 let second_us = rem % 60_000_000;
373 (year, month, day, hour, minute, second_us)
374}
375
376fn leap_second_civil(split: JulianDateSplit) -> (i64, i64, i64, i64, i64, i64) {
381 let next_day_number = (split.jd_whole + 0.5).round() as i64;
382 let (year, month, day) = civil_from_julian_day_number(next_day_number - 1);
383 let remaining_s = -split.fraction * SECONDS_PER_DAY; let microsecond = ((1.0 - remaining_s) * 1_000_000.0).round() as i64;
385 (year, month, day, 23, 59, 60 * 1_000_000 + microsecond)
388}
389
390fn nanos_civil_day_microsecond(nanos: i128) -> (i64, i64) {
394 const US_PER_DAY: i128 = SECONDS_PER_DAY_I64 as i128 * 1_000_000;
395 const J2000_NOON_US: i128 = J2000_NOON_OFFSET_S as i128 * 1_000_000;
398 const J2000_DAY_NUMBER: i128 = J2000_JULIAN_DAY_NUMBER as i128;
399 let micros = (nanos + nanos.signum() * 500) / 1_000; let from_midnight = J2000_NOON_US + micros;
401 let day_offset = from_midnight.div_euclid(US_PER_DAY);
402 let us_of_day = from_midnight.rem_euclid(US_PER_DAY);
403 ((J2000_DAY_NUMBER + day_offset) as i64, us_of_day as i64)
404}
405
406pub fn civil_to_clock_instant(
408 scale: TimeScale,
409 year: i32,
410 month: u8,
411 day: u8,
412 hour: u8,
413 minute: u8,
414 second: f64,
415) -> Option<Instant> {
416 let civil = validate::civil_datetime_with_fractional_second_policy(
417 i64::from(year),
418 i64::from(month),
419 i64::from(day),
420 i64::from(hour),
421 i64::from(minute),
422 second,
423 civil_second_policy_for_time_scale(scale),
424 )
425 .ok()?;
426 civil_microsecond_to_instant(scale, civil).ok()
427}
428
429pub fn civil_to_gps_seconds(
431 year: i32,
432 month: u8,
433 day: u8,
434 hour: u8,
435 minute: u8,
436 second: f64,
437) -> Option<f64> {
438 let civil = validate::civil_datetime_with_fractional_second_policy(
439 i64::from(year),
440 i64::from(month),
441 i64::from(day),
442 i64::from(hour),
443 i64::from(minute),
444 second,
445 validate::CivilSecondPolicy::Continuous,
446 )
447 .ok()?;
448 gps_seconds_from_civil(civil)
449}
450
451fn parse_time_scale(text: &str) -> Result<TimeScale, RinexClockError> {
452 let mut time_scale = TimeScale::Gpst;
453 for (idx, line) in text.lines().enumerate() {
454 if line.contains("END OF HEADER") {
455 break;
456 }
457 if line.contains("TIME SYSTEM ID") {
458 let label = line
459 .split("TIME SYSTEM ID")
460 .next()
461 .unwrap_or(line)
462 .split_whitespace()
463 .next()
464 .unwrap_or("");
465 if label.is_empty() {
466 time_scale = TimeScale::Gpst;
467 } else {
468 time_scale = crate::rinex_common::time_scale_label(label).ok_or_else(|| {
469 RinexClockError::BadField {
470 line: idx + 1,
471 field: "time_system",
472 value: label.to_string(),
473 }
474 })?;
475 }
476 }
477 }
478 Ok(time_scale)
479}
480
481fn gps_seconds_to_instant(gps_seconds: f64) -> Instant {
482 let gps_epoch_jd = J2000_JD - GPS_EPOCH_TO_J2000_S / SECONDS_PER_DAY;
483 let days = (gps_seconds / SECONDS_PER_DAY).floor();
484 let seconds_of_day = gps_seconds - days * SECONDS_PER_DAY;
485 Instant::from_julian_date(
486 TimeScale::Gpst,
487 JulianDateSplit::new(gps_epoch_jd + days, seconds_of_day / SECONDS_PER_DAY)
488 .expect("valid split Julian date"),
489 )
490}
491
492fn validate_clock_point(point: ClockPoint) -> Result<(), RinexClockError> {
493 validate_instant(point.epoch, "epoch")?;
494 validate_finite(point.bias_s, "bias_s")
495}
496
497fn validate_serializable_clock_point(
498 product_scale: TimeScale,
499 point: &ClockPoint,
500) -> Result<(), RinexClockError> {
501 if crate::rinex_common::time_scale_rinex_label(point.epoch.scale).is_none() {
502 return Err(RinexClockError::UnsupportedTimeScale {
503 scale: point.epoch.scale,
504 });
505 }
506 if point.epoch.scale != product_scale {
507 return Err(invalid_input(
508 "epoch",
509 "epoch scale does not match clock time scale",
510 ));
511 }
512 Ok(())
513}
514
515fn validate_instant(epoch: Instant, field: &'static str) -> Result<(), RinexClockError> {
516 match epoch.repr {
517 InstantRepr::JulianDate(split) => {
518 validate_finite(split.jd_whole, field)?;
519 validate_finite(split.fraction, field)?;
520 if !(-1.0..=1.0).contains(&split.fraction) {
521 return Err(invalid_input(field, "Julian-date fraction out of range"));
522 }
523 Ok(())
524 }
525 InstantRepr::Nanos(_) => Ok(()),
526 }
527}
528
529fn validate_finite(value: f64, field: &'static str) -> Result<(), RinexClockError> {
530 if value.is_finite() {
531 Ok(())
532 } else {
533 Err(invalid_input(field, "must be finite"))
534 }
535}
536
537fn invalid_input(field: &'static str, reason: &'static str) -> RinexClockError {
538 RinexClockError::InvalidInput { field, reason }
539}
540
541fn map_manual_order_error(error: FieldError) -> RinexClockError {
542 match error {
543 FieldError::NonFinite { field } => invalid_input(field, "must be finite"),
544 FieldError::OutOfRange { field, .. } => invalid_input(field, "must be strictly increasing"),
545 _ => invalid_input(error.field(), error.reason()),
546 }
547}
548
549fn validate_instant_series_order(points: &[(ClockPoint, usize)]) -> Result<(), RinexClockError> {
550 validate::require_strictly_increasing(
551 points
552 .iter()
553 .map(|(point, _)| instant_order_key(&point.epoch)),
554 "epoch",
555 )
556 .map_err(map_manual_order_error)
557}
558
559fn instant_order_key(epoch: &Instant) -> f64 {
560 let offset_s = time_scale_rank(epoch.scale) as f64 * INSTANT_SCALE_ORDER_STRIDE_S;
561 let instant_s = match epoch.repr {
562 InstantRepr::JulianDate(split) => {
563 split.jd_whole * SECONDS_PER_DAY + split.fraction * SECONDS_PER_DAY
564 }
565 InstantRepr::Nanos(nanos) => nanos as f64 / 1.0e9,
566 };
567 offset_s + instant_s
568}
569
570fn instant_to_gps_seconds(epoch: &Instant) -> Option<f64> {
571 if epoch.scale != TimeScale::Gpst {
572 return None;
573 }
574 instant_to_j2000_seconds(epoch).map(|seconds| seconds + GPS_EPOCH_TO_J2000_S)
575}
576
577fn instant_to_j2000_seconds(epoch: &Instant) -> Option<f64> {
578 match epoch.repr {
579 InstantRepr::JulianDate(split) => {
580 Some(j2000_seconds_from_split(split.jd_whole, split.fraction))
581 }
582 InstantRepr::Nanos(_) => None,
583 }
584}
585
586fn data_lines(text: &str) -> Vec<(usize, &str)> {
587 drop_header(
588 text.lines()
589 .enumerate()
590 .map(|(idx, line)| (idx + 1, line))
591 .collect(),
592 )
593}
594
595fn drop_header(lines: Vec<(usize, &str)>) -> Vec<(usize, &str)> {
596 match lines
597 .iter()
598 .position(|(_, line)| line.contains("END OF HEADER"))
599 {
600 Some(idx) => lines.into_iter().skip(idx + 1).collect(),
601 None => lines,
602 }
603}
604
605#[derive(Debug, Clone, Copy)]
606struct ClockEpochFields<'a> {
607 year: i32,
608 month: u8,
609 day: u8,
610 hour: u8,
611 minute: u8,
612 second: &'a str,
613}
614
615fn parse_record(
616 line_number: usize,
617 line: &str,
618 time_scale: TimeScale,
619) -> Result<Option<(String, ClockPoint)>, RinexClockError> {
620 let mut fields = line.split_whitespace();
621 if fields.next() != Some("AS") {
622 return Ok(None);
623 }
624
625 let sat_field = next_as_field(&mut fields, line_number, line)?;
626 let year_field = next_as_field(&mut fields, line_number, line)?;
627 let month_field = next_as_field(&mut fields, line_number, line)?;
628 let day_field = next_as_field(&mut fields, line_number, line)?;
629 let hour_field = next_as_field(&mut fields, line_number, line)?;
630 let minute_field = next_as_field(&mut fields, line_number, line)?;
631 let second_field = next_as_field(&mut fields, line_number, line)?;
632 let _value_count_field = next_as_field(&mut fields, line_number, line)?;
633 let bias_field = next_as_field(&mut fields, line_number, line)?;
634
635 let sat = validate::strict_gnss_satellite_id(sat_field, "satellite")
636 .map_err(|error| map_field_error(line_number, error, sat_field))?
637 .to_string();
638 let year = parse_int_field::<i32>(line_number, "year", year_field)?;
639 let month = parse_int_field::<u8>(line_number, "month", month_field)?;
640 let day = parse_int_field::<u8>(line_number, "day", day_field)?;
641 let hour = parse_int_field::<u8>(line_number, "hour", hour_field)?;
642 let minute = parse_int_field::<u8>(line_number, "minute", minute_field)?;
643 let epoch = ClockEpochFields {
644 year,
645 month,
646 day,
647 hour,
648 minute,
649 second: second_field,
650 };
651 let bias_s = parse_f64_field(line_number, "bias", bias_field)?;
652 let epoch = civil_decimal_second_to_instant(time_scale, epoch)
653 .map_err(|error| map_epoch_error(line_number, error, epoch))?;
654
655 Ok(Some((sat, ClockPoint { epoch, bias_s })))
656}
657
658fn next_as_field<'a, I>(
659 fields: &mut I,
660 line_number: usize,
661 line: &str,
662) -> Result<&'a str, RinexClockError>
663where
664 I: Iterator<Item = &'a str>,
665{
666 fields
667 .next()
668 .ok_or_else(|| RinexClockError::MalformedAsRecord {
669 line: line_number,
670 reason: "expected at least 10 fields",
671 record: line.trim().to_string(),
672 })
673}
674
675fn parse_int_field<T>(
676 line_number: usize,
677 field: &'static str,
678 value: &str,
679) -> Result<T, RinexClockError>
680where
681 T: std::str::FromStr,
682{
683 validate::strict_int(value, field).map_err(|error| map_field_error(line_number, error, value))
684}
685
686fn parse_f64_field(
687 line_number: usize,
688 field: &'static str,
689 value: &str,
690) -> Result<f64, RinexClockError> {
691 validate::strict_f64(value, field).map_err(|error| map_field_error(line_number, error, value))
692}
693
694fn civil_decimal_second_to_instant(
695 scale: TimeScale,
696 epoch: ClockEpochFields<'_>,
697) -> Result<Instant, FieldError> {
698 let civil = validate::civil_datetime_with_decimal_second_policy(
699 i64::from(epoch.year),
700 i64::from(epoch.month),
701 i64::from(epoch.day),
702 i64::from(epoch.hour),
703 i64::from(epoch.minute),
704 epoch.second,
705 civil_second_policy_for_time_scale(scale),
706 )?;
707 civil_microsecond_to_instant(scale, civil)
708}
709
710fn civil_microsecond_to_instant(
711 scale: TimeScale,
712 civil: validate::ValidCivilMicrosecond,
713) -> Result<Instant, FieldError> {
714 let split = civil_microsecond_to_julian_split(scale, civil)?;
715 Ok(Instant::from_julian_date(scale, split))
716}
717
718fn civil_microsecond_to_julian_split(
719 scale: TimeScale,
720 civil: validate::ValidCivilMicrosecond,
721) -> Result<JulianDateSplit, FieldError> {
722 if civil.year < 1 {
723 return Err(FieldError::InvalidCivilDate {
724 field: "civil datetime",
725 year: civil.year,
726 month: i64::from(civil.month),
727 day: i64::from(civil.day),
728 });
729 }
730
731 let jdn = julian_day_number(civil.year as i32, civil.month as i32, civil.day as i32);
732 let jd_whole = jdn as f64 - 0.5;
733 if scale == TimeScale::Utc && civil.second == 60 {
734 let remaining_s = 1.0 - civil.microsecond as f64 / 1_000_000.0;
735 return Ok(
736 JulianDateSplit::new(jd_whole + 1.0, -remaining_s / SECONDS_PER_DAY)
737 .expect("valid leap-second split Julian date"),
738 );
739 }
740
741 let day_seconds = civil.hour as f64 * SECONDS_PER_HOUR
742 + civil.minute as f64 * 60.0
743 + civil.second as f64
744 + civil.microsecond as f64 / 1_000_000.0;
745 Ok(
746 JulianDateSplit::new(jd_whole, day_seconds / SECONDS_PER_DAY)
747 .expect("valid split Julian date"),
748 )
749}
750
751fn civil_second_policy_for_time_scale(scale: TimeScale) -> validate::CivilSecondPolicy {
752 match scale {
753 TimeScale::Utc => validate::CivilSecondPolicy::UtcLike,
754 TimeScale::Glonasst
760 | TimeScale::Tai
761 | TimeScale::Tt
762 | TimeScale::Tcg
763 | TimeScale::Tdb
764 | TimeScale::Tcb
765 | TimeScale::Gpst
766 | TimeScale::Gst
767 | TimeScale::Bdt
768 | TimeScale::Qzsst => validate::CivilSecondPolicy::Continuous,
769 }
770}
771
772fn gps_seconds_from_civil(civil: validate::ValidCivilMicrosecond) -> Option<f64> {
773 if civil.year < 1 {
774 return None;
775 }
776
777 let days = days_since_gps_epoch(civil.year as i32, civil.month as u8, civil.day as u8);
778 let whole = days as f64 * SECONDS_PER_DAY
779 + (i64::from(civil.hour) * 3_600 + i64::from(civil.minute) * 60 + i64::from(civil.second))
780 as f64;
781 Some(whole + f64::from(civil.microsecond) / 1_000_000.0)
782}
783
784fn map_field_error(line_number: usize, error: FieldError, value: &str) -> RinexClockError {
785 RinexClockError::BadField {
786 line: line_number,
787 field: error.field(),
788 value: value.to_string(),
789 }
790}
791
792fn map_epoch_error(
793 line_number: usize,
794 error: FieldError,
795 epoch: ClockEpochFields<'_>,
796) -> RinexClockError {
797 match error {
798 FieldError::FloatParse { .. }
799 | FieldError::Missing { .. }
800 | FieldError::NonFinite { .. } => RinexClockError::BadField {
801 line: line_number,
802 field: "second",
803 value: epoch.second.to_string(),
804 },
805 _ => RinexClockError::BadField {
806 line: line_number,
807 field: "epoch",
808 value: format!(
809 "{} {} {} {} {} {}",
810 epoch.year,
811 epoch.month,
812 epoch.day,
813 epoch.hour,
814 epoch.minute,
815 normalized_second_text(epoch.second)
816 ),
817 },
818 }
819}
820
821fn normalized_second_text(second: &str) -> String {
822 validate::strict_f64(second, "second")
823 .map_or_else(|_| second.to_string(), |value| value.to_string())
824}
825
826fn build_series(
827 by_sat: BTreeMap<String, Vec<(ClockPoint, usize)>>,
828) -> BTreeMap<String, Vec<ClockPoint>> {
829 by_sat
830 .into_iter()
831 .map(|(sat, mut points)| {
832 points.sort_by(|(a, ai), (b, bi)| {
833 compare_instants(&a.epoch, &b.epoch).then_with(|| ai.cmp(bi))
834 });
835 (sat, dedup_by_time(points))
836 })
837 .collect()
838}
839
840fn dedup_by_time(points: Vec<(ClockPoint, usize)>) -> Vec<ClockPoint> {
841 let mut deduped = Vec::<ClockPoint>::new();
842 for (point, _) in points {
843 match deduped.last_mut() {
844 Some(prev) if prev.epoch == point.epoch => *prev = point,
845 _ => deduped.push(point),
846 }
847 }
848 deduped
849}
850
851fn interpolate(records: &[ClockPoint], epoch: Instant) -> Option<f64> {
852 let mut prev: Option<ClockPoint> = None;
853 for point in records {
854 match compare_instants_same_scale(&point.epoch, &epoch)? {
855 Ordering::Equal => return Some(point.bias_s),
856 Ordering::Greater => {
857 let p0 = prev?;
858 let p1 = *point;
859 let span_s = seconds_between(&p1.epoch, &p0.epoch)?;
860 if span_s <= 0.0 {
861 return None;
862 }
863 let query_s = seconds_between(&epoch, &p0.epoch)?;
864 if query_s < 0.0 {
865 return None;
866 }
867 return Some(lerp_ratio(p0.bias_s, p1.bias_s, query_s, span_s));
868 }
869 Ordering::Less => prev = Some(*point),
870 }
871 }
872 None
873}
874
875fn compare_instants(a: &Instant, b: &Instant) -> Ordering {
876 time_scale_rank(a.scale)
877 .cmp(&time_scale_rank(b.scale))
878 .then_with(|| match (a.julian_date(), b.julian_date()) {
879 (Some(a), Some(b)) => compare_julian_splits(a, b),
880 _ => Ordering::Equal,
881 })
882}
883
884fn clock_timeline(scale: TimeScale) -> TimeScale {
893 match scale {
894 TimeScale::Qzsst => TimeScale::Gpst,
895 other => other,
896 }
897}
898
899fn compare_instants_same_scale(a: &Instant, b: &Instant) -> Option<Ordering> {
900 if clock_timeline(a.scale) != clock_timeline(b.scale) {
901 return None;
902 }
903 Some(compare_julian_splits(a.julian_date()?, b.julian_date()?))
904}
905
906fn compare_julian_splits(a: JulianDateSplit, b: JulianDateSplit) -> Ordering {
907 a.jd_whole
908 .partial_cmp(&b.jd_whole)
909 .unwrap_or(Ordering::Equal)
910 .then_with(|| {
911 a.fraction
912 .partial_cmp(&b.fraction)
913 .unwrap_or(Ordering::Equal)
914 })
915}
916
917fn seconds_between(later: &Instant, earlier: &Instant) -> Option<f64> {
918 if clock_timeline(later.scale) != clock_timeline(earlier.scale) {
919 return None;
920 }
921 let later = later.julian_date()?;
922 let earlier = earlier.julian_date()?;
923 let seconds = seconds_between_splits(
924 later.jd_whole,
925 later.fraction,
926 earlier.jd_whole,
927 earlier.fraction,
928 );
929 seconds.is_finite().then_some(seconds)
930}
931
932fn time_scale_rank(scale: TimeScale) -> u8 {
933 match scale {
934 TimeScale::Utc => 0,
935 TimeScale::Tai => 1,
936 TimeScale::Tt => 2,
937 TimeScale::Tcg => 3,
938 TimeScale::Tdb => 4,
939 TimeScale::Tcb => 5,
940 TimeScale::Gpst => 6,
941 TimeScale::Gst => 7,
942 TimeScale::Bdt => 8,
943 TimeScale::Glonasst => 9,
944 TimeScale::Qzsst => 10,
945 }
946}
947
948fn days_since_gps_epoch(year: i32, month: u8, day: u8) -> i64 {
949 julian_day_number(year, i32::from(month), i32::from(day)) - julian_day_number(1980, 1, 6)
950}
951
952#[cfg(test)]
953mod tests {
954 use super::*;
955
956 fn as_record(satellite: &str, bias: &str) -> String {
957 format!("AS {satellite} 2020 01 01 00 00 00.000000 1 {bias}")
958 }
959
960 #[test]
961 fn parse_rejects_non_finite_as_bias() {
962 let err = RinexClock::parse(&as_record("G01", "NaN")).unwrap_err();
963 assert_eq!(
964 err,
965 RinexClockError::BadField {
966 line: 1,
967 field: "bias",
968 value: "NaN".to_string(),
969 }
970 );
971 }
972
973 #[test]
974 fn parse_rejects_malformed_as_satellite_token() {
975 let err = RinexClock::parse(&as_record("X01", "1.0e-9")).unwrap_err();
976 assert_eq!(
977 err,
978 RinexClockError::BadField {
979 line: 1,
980 field: "satellite",
981 value: "X01".to_string(),
982 }
983 );
984 }
985
986 #[test]
987 fn explicit_utc_time_system_preserves_clock_epoch_scale() {
988 let text = " 3.00 C RINEX VERSION / TYPE\n\
989 UTC TIME SYSTEM ID\n\
990 END OF HEADER\n\
991 AS G05 2017 01 01 00 00 0.000000 1 1.0e-04\n\
992 AS G05 2017 01 01 00 00 30.000000 1 2.0e-04\n";
993 let clock = RinexClock::parse(text).expect("UTC RINEX clock");
994
995 assert_eq!(clock.time_scale, TimeScale::Utc);
996 assert_eq!(clock.series["G05"][0].epoch.scale, TimeScale::Utc);
997 let interpolated = clock
998 .clock_s(
999 "G05",
1000 ClockEpoch {
1001 year: 2017,
1002 month: 1,
1003 day: 1,
1004 hour: 0,
1005 minute: 0,
1006 second: 15.0,
1007 },
1008 )
1009 .expect("valid clock query")
1010 .expect("UTC interpolated clock");
1011 assert!((interpolated - 1.5e-4).abs() < 1.0e-18);
1012
1013 let gpst_query =
1014 civil_to_clock_instant(TimeScale::Gpst, 2017, 1, 1, 0, 0, 15.0).expect("GPST instant");
1015 assert_eq!(
1016 clock
1017 .clock_s_at_instant("G05", gpst_query)
1018 .expect("valid clock query"),
1019 None
1020 );
1021
1022 let rows = clock.instant_series_rows();
1023 assert_eq!(rows[0].1[0].0.scale, TimeScale::Utc);
1024 let rebuilt = RinexClock::from_instant_series_rows(clock.time_scale, rows)
1025 .expect("valid manual RINEX clock rows");
1026 assert_eq!(rebuilt, clock);
1027 }
1028
1029 #[test]
1030 fn manual_series_rows_reject_non_finite_inputs() {
1031 assert_eq!(
1032 RinexClock::from_series_rows(vec![("G05".to_string(), vec![(f64::NAN, 1.0e-4)])])
1033 .unwrap_err(),
1034 RinexClockError::InvalidInput {
1035 field: "gps_seconds",
1036 reason: "must be finite",
1037 }
1038 );
1039 assert_eq!(
1040 RinexClock::from_series_rows(vec![(
1041 "G05".to_string(),
1042 vec![(1_463_904_000.0, f64::INFINITY)]
1043 )])
1044 .unwrap_err(),
1045 RinexClockError::InvalidInput {
1046 field: "bias_s",
1047 reason: "must be finite",
1048 }
1049 );
1050 }
1051
1052 #[test]
1053 fn manual_series_rows_reject_unsorted_gps_seconds() {
1054 assert_eq!(
1055 RinexClock::from_series_rows(vec![(
1056 "G05".to_string(),
1057 vec![(1_463_904_030.0, 1.0e-4), (1_463_904_000.0, 2.0e-4)]
1058 )])
1059 .unwrap_err(),
1060 RinexClockError::InvalidInput {
1061 field: "gps_seconds",
1062 reason: "must be strictly increasing",
1063 }
1064 );
1065 }
1066
1067 #[test]
1068 fn manual_instant_rows_reject_non_finite_inputs() {
1069 let bad_epoch = Instant::from_julian_date(
1070 TimeScale::Gpst,
1071 JulianDateSplit {
1072 jd_whole: f64::NAN,
1073 fraction: 0.0,
1074 },
1075 );
1076 assert_eq!(
1077 RinexClock::from_instant_series_rows(
1078 TimeScale::Gpst,
1079 vec![("G05".to_string(), vec![(bad_epoch, 1.0e-4)])],
1080 )
1081 .unwrap_err(),
1082 RinexClockError::InvalidInput {
1083 field: "epoch",
1084 reason: "must be finite",
1085 }
1086 );
1087
1088 let good_epoch =
1089 civil_to_clock_instant(TimeScale::Gpst, 2026, 5, 13, 0, 0, 0.0).expect("GPST instant");
1090 assert_eq!(
1091 RinexClock::from_instant_series_rows(
1092 TimeScale::Gpst,
1093 vec![("G05".to_string(), vec![(good_epoch, f64::NAN)])],
1094 )
1095 .unwrap_err(),
1096 RinexClockError::InvalidInput {
1097 field: "bias_s",
1098 reason: "must be finite",
1099 }
1100 );
1101 }
1102
1103 #[test]
1104 fn manual_instant_rows_reject_unsorted_epochs() {
1105 let later =
1106 civil_to_clock_instant(TimeScale::Gpst, 2026, 5, 13, 0, 0, 30.0).expect("later epoch");
1107 let earlier =
1108 civil_to_clock_instant(TimeScale::Gpst, 2026, 5, 13, 0, 0, 0.0).expect("earlier epoch");
1109
1110 assert_eq!(
1111 RinexClock::from_instant_series_rows(
1112 TimeScale::Gpst,
1113 vec![("G05".to_string(), vec![(later, 1.0e-4), (earlier, 2.0e-4)])],
1114 )
1115 .unwrap_err(),
1116 RinexClockError::InvalidInput {
1117 field: "epoch",
1118 reason: "must be strictly increasing",
1119 }
1120 );
1121 }
1122
1123 #[test]
1124 fn rinex_clock_queries_reject_non_finite_inputs() {
1125 let clock = RinexClock::from_series_rows(vec![(
1126 "G05".to_string(),
1127 vec![(1_463_904_000.0, 1.0e-4)],
1128 )])
1129 .expect("valid manual RINEX clock rows");
1130 let bad_epoch = Instant::from_julian_date(
1131 TimeScale::Gpst,
1132 JulianDateSplit {
1133 jd_whole: f64::INFINITY,
1134 fraction: 0.0,
1135 },
1136 );
1137 assert_eq!(
1138 clock.clock_s_at_instant("G05", bad_epoch).unwrap_err(),
1139 RinexClockError::InvalidInput {
1140 field: "epoch",
1141 reason: "must be finite",
1142 }
1143 );
1144 assert_eq!(
1145 clock.clock_s_at_gps_seconds("G05", f64::NAN).unwrap_err(),
1146 RinexClockError::InvalidInput {
1147 field: "gps_seconds",
1148 reason: "must be finite",
1149 }
1150 );
1151 assert_eq!(
1152 clock
1153 .clock_s(
1154 "G05",
1155 ClockEpoch {
1156 year: 2026,
1157 month: 5,
1158 day: 13,
1159 hour: 0,
1160 minute: 0,
1161 second: f64::NAN,
1162 },
1163 )
1164 .unwrap_err(),
1165 RinexClockError::InvalidInput {
1166 field: "epoch",
1167 reason: "invalid civil clock epoch",
1168 }
1169 );
1170 }
1171
1172 #[test]
1173 fn interpolation_rejects_non_positive_bracket_span() {
1174 let day = 2_457_753.5;
1175 let p0 = Instant::from_julian_date(
1176 TimeScale::Utc,
1177 JulianDateSplit::new(day, 1.0).expect("valid split Julian date"),
1178 );
1179 let p1 = Instant::from_julian_date(
1180 TimeScale::Utc,
1181 JulianDateSplit::new(day + 1.0, 0.0).expect("valid split Julian date"),
1182 );
1183 let query = Instant::from_julian_date(
1184 TimeScale::Utc,
1185 JulianDateSplit::new(day + 1.0, 0.5 / SECONDS_PER_DAY)
1186 .expect("valid split Julian date"),
1187 );
1188 let records = [
1189 ClockPoint {
1190 epoch: p0,
1191 bias_s: 1.0e-4,
1192 },
1193 ClockPoint {
1194 epoch: p1,
1195 bias_s: 2.0e-4,
1196 },
1197 ];
1198
1199 assert_eq!(interpolate(&records, query), None);
1200 }
1201
1202 #[test]
1203 fn qzsst_rows_are_queryable_on_the_gpst_timeline() {
1204 let p0 = civil_to_clock_instant(TimeScale::Qzsst, 2026, 5, 13, 0, 0, 0.0)
1208 .expect("QZSST instant");
1209 let p1 = civil_to_clock_instant(TimeScale::Qzsst, 2026, 5, 13, 0, 0, 30.0)
1210 .expect("QZSST instant");
1211 let clock = RinexClock::from_instant_series_rows(
1212 TimeScale::Qzsst,
1213 vec![("J02".to_string(), vec![(p0, 1.0e-4), (p1, 3.0e-4)])],
1214 )
1215 .expect("QZSST clock builds");
1216
1217 let mid = civil_to_gps_seconds(2026, 5, 13, 0, 0, 15.0).expect("gps seconds");
1220 let bias = clock
1221 .clock_s_at_gps_seconds("J02", mid)
1222 .expect("query succeeds")
1223 .expect("QZSST row interpolates on the GPST timeline");
1224 assert!(
1225 (bias - 2.0e-4).abs() < 1.0e-12,
1226 "expected midpoint interpolation 2.0e-4, got {bias}"
1227 );
1228
1229 let start = civil_to_gps_seconds(2026, 5, 13, 0, 0, 0.0).expect("gps seconds");
1231 assert_eq!(
1232 clock
1233 .clock_s_at_gps_seconds("J02", start)
1234 .expect("query succeeds"),
1235 Some(1.0e-4)
1236 );
1237 }
1238
1239 #[test]
1240 fn to_rinex_string_round_trips_through_parse() {
1241 let text =
1245 " 3.00 C RINEX VERSION / TYPE\n\
1246 GPS TIME SYSTEM ID\n\
1247 END OF HEADER\n\
1248 AS G05 2026 05 13 00 00 0.000000 1 -2.000000000000e-04\n\
1249 AS G05 2026 05 13 00 00 30.500000 1 -2.000000600000e-04\n\
1250 AS G24 2026 05 13 00 01 0.000000 1 5.000000000000e-05\n\
1251 AS E11 2026 05 13 00 00 0.000000 1 1.234500000000e-09\n";
1252 let clock = RinexClock::parse(text).expect("parse GPST RINEX clock");
1253 let serialized = clock.to_rinex_string().expect("serialize RINEX clock");
1254 let reparsed = RinexClock::parse(&serialized).expect("re-parse serialized");
1255 assert_eq!(reparsed, clock, "serializer must round-trip through parse");
1256 assert_eq!(
1258 reparsed
1259 .to_rinex_string()
1260 .expect("serialize reparsed clock"),
1261 serialized
1262 );
1263 }
1264
1265 #[test]
1266 fn to_rinex_string_round_trips_utc_time_scale() {
1267 let text =
1269 " 3.00 C RINEX VERSION / TYPE\n\
1270 UTC TIME SYSTEM ID\n\
1271 END OF HEADER\n\
1272 AS G05 2017 01 01 00 00 0.000000 1 1.000000000000e-04\n\
1273 AS G05 2017 01 01 00 00 30.000000 1 2.000000000000e-04\n";
1274 let clock = RinexClock::parse(text).expect("parse UTC RINEX clock");
1275 assert_eq!(clock.time_scale, TimeScale::Utc);
1276 let serialized = clock.to_rinex_string().expect("serialize RINEX clock");
1277 let reparsed = RinexClock::parse(&serialized).expect("re-parse serialized");
1278 assert_eq!(reparsed.time_scale, TimeScale::Utc);
1279 assert_eq!(reparsed, clock);
1280 }
1281
1282 #[test]
1283 fn to_rinex_string_rejects_unsupported_time_scale() {
1284 let epoch =
1285 civil_to_clock_instant(TimeScale::Tcg, 2026, 5, 13, 0, 0, 0.0).expect("TCG instant");
1286 let clock = RinexClock::from_instant_series_rows(
1287 TimeScale::Tcg,
1288 vec![("G05".to_string(), vec![(epoch, 1.0e-4)])],
1289 )
1290 .expect("TCG clock builds");
1291
1292 assert_eq!(
1293 clock.to_rinex_string(),
1294 Err(RinexClockError::UnsupportedTimeScale {
1295 scale: TimeScale::Tcg
1296 })
1297 );
1298 }
1299
1300 #[test]
1301 fn to_rinex_string_rejects_unsupported_row_time_scale() {
1302 let epoch =
1303 civil_to_clock_instant(TimeScale::Tcg, 2026, 5, 13, 0, 0, 0.0).expect("TCG instant");
1304 let clock = RinexClock::from_instant_series_rows(
1305 TimeScale::Gpst,
1306 vec![("G05".to_string(), vec![(epoch, 1.0e-4)])],
1307 )
1308 .expect("mixed-scale clock builds");
1309
1310 assert_eq!(
1311 clock.to_rinex_string(),
1312 Err(RinexClockError::UnsupportedTimeScale {
1313 scale: TimeScale::Tcg
1314 })
1315 );
1316 }
1317
1318 #[test]
1319 fn nanos_repr_epoch_serializes_to_true_civil_time() {
1320 let jd_epoch =
1326 civil_to_clock_instant(TimeScale::Gpst, 2026, 5, 13, 0, 0, 30.0).expect("GPST instant");
1327 let j2000_s = instant_to_j2000_seconds(&jd_epoch).expect("J2000 seconds");
1328 let nanos = (j2000_s * 1.0e9).round() as i128;
1329 let nanos_epoch = Instant::from_nanos(TimeScale::Gpst, nanos);
1330
1331 let nanos_clock = RinexClock::from_instant_series_rows(
1332 TimeScale::Gpst,
1333 vec![("G05".to_string(), vec![(nanos_epoch, 1.0e-4)])],
1334 )
1335 .expect("nanos clock builds");
1336 let jd_clock = RinexClock::from_instant_series_rows(
1337 TimeScale::Gpst,
1338 vec![("G05".to_string(), vec![(jd_epoch, 1.0e-4)])],
1339 )
1340 .expect("jd clock builds");
1341
1342 let serialized = nanos_clock
1343 .to_rinex_string()
1344 .expect("serialize nanos RINEX clock");
1345 assert!(
1346 serialized.contains("2026 05 13 00 00 30.000000"),
1347 "Nanos epoch must serialize to its true civil time, got:\n{serialized}"
1348 );
1349 assert_eq!(
1350 serialized,
1351 jd_clock
1352 .to_rinex_string()
1353 .expect("serialize JD RINEX clock"),
1354 "Nanos- and Julian-date-repr epochs of the same instant must serialize identically"
1355 );
1356
1357 let reparsed = RinexClock::parse(&serialized).expect("re-parse serialized Nanos product");
1358 assert_eq!(reparsed, jd_clock);
1359 }
1360
1361 #[test]
1362 fn to_rinex_string_round_trips_utc_leap_second_epoch() {
1363 let text =
1368 " 3.00 C RINEX VERSION / TYPE\n\
1369 UTC TIME SYSTEM ID\n\
1370 END OF HEADER\n\
1371 AS G05 2016 12 31 23 59 60.000000 1 1.000000000000e-04\n\
1372 AS G05 2016 12 31 23 59 60.500000 1 2.000000000000e-04\n";
1373 let clock = RinexClock::parse(text).expect("parse UTC leap-second RINEX clock");
1374 let serialized = clock.to_rinex_string().expect("serialize RINEX clock");
1375 assert!(
1376 serialized.contains("23 59 60.000000"),
1377 "leap-second label must round-trip, got:\n{serialized}"
1378 );
1379 assert!(
1380 serialized.contains("23 59 60.500000"),
1381 "fractional leap second must round-trip, got:\n{serialized}"
1382 );
1383 let reparsed = RinexClock::parse(&serialized).expect("re-parse serialized leap second");
1384 assert_eq!(
1385 reparsed, clock,
1386 "leap-second epoch must round-trip bit-exact"
1387 );
1388 }
1389}