1#![allow(dead_code)]
9
10use core::str::FromStr;
11
12use crate::id::GnssSatelliteId;
13
14#[derive(Debug, Clone, Copy, PartialEq, Eq, thiserror::Error)]
15#[error("{field} length is {actual}, expected {expected}")]
16pub(crate) struct LengthError {
17 pub(crate) field: &'static str,
18 pub(crate) expected: usize,
19 pub(crate) actual: usize,
20}
21
22#[derive(Debug, Clone, Copy, PartialEq, Eq, thiserror::Error)]
23#[error("{field} integer arithmetic overflow")]
24pub(crate) struct ArithmeticError {
25 pub(crate) field: &'static str,
26}
27
28#[derive(Debug, Clone, PartialEq, thiserror::Error)]
29pub enum FieldError {
30 #[error("{field} is missing")]
31 Missing { field: &'static str },
32 #[error("{field} is not a finite number")]
33 NonFinite { field: &'static str },
34 #[error("{field} must be positive")]
35 NotPositive { field: &'static str },
36 #[error("{field} must be non-negative")]
37 Negative { field: &'static str },
38 #[error("{field} is out of range")]
39 OutOfRange {
40 field: &'static str,
41 min: f64,
42 max: f64,
43 upper_inclusive: bool,
44 },
45 #[error("{field} is not a valid float: {value:?}")]
46 FloatParse { field: &'static str, value: String },
47 #[error("{field} is not a valid integer: {value:?}")]
48 IntParse { field: &'static str, value: String },
49 #[error("{field} is not a valid civil date: {year:04}-{month:02}-{day:02}")]
50 InvalidCivilDate {
51 field: &'static str,
52 year: i64,
53 month: i64,
54 day: i64,
55 },
56 #[error("{field} is not a valid civil time: {hour:02}:{minute:02}:{second}")]
57 InvalidCivilTime {
58 field: &'static str,
59 hour: i64,
60 minute: i64,
61 second: f64,
62 },
63}
64
65impl FieldError {
66 pub const fn field(&self) -> &'static str {
67 match self {
68 Self::Missing { field }
69 | Self::NonFinite { field }
70 | Self::NotPositive { field }
71 | Self::Negative { field }
72 | Self::OutOfRange { field, .. }
73 | Self::FloatParse { field, .. }
74 | Self::IntParse { field, .. }
75 | Self::InvalidCivilDate { field, .. }
76 | Self::InvalidCivilTime { field, .. } => field,
77 }
78 }
79
80 pub const fn reason(&self) -> &'static str {
81 match self {
82 Self::Missing { .. } => "missing",
83 Self::NonFinite { .. } => "not finite",
84 Self::NotPositive { .. } => "not positive",
85 Self::Negative { .. } => "negative",
86 Self::OutOfRange { .. } => "out of range",
87 Self::FloatParse { .. } => "invalid float",
88 Self::IntParse { .. } => "invalid integer",
89 Self::InvalidCivilDate { .. } => "invalid civil date",
90 Self::InvalidCivilTime { .. } => "invalid civil time",
91 }
92 }
93}
94
95#[derive(Debug, Clone, Copy, PartialEq)]
96pub(crate) struct ValidCivil {
97 pub(crate) year: i64,
98 pub(crate) month: u32,
99 pub(crate) day: u32,
100 pub(crate) hour: u32,
101 pub(crate) minute: u32,
102 pub(crate) second: f64,
103}
104
105#[derive(Debug, Clone, Copy, PartialEq, Eq)]
106pub(crate) struct ValidCivilMicrosecond {
107 pub(crate) year: i64,
108 pub(crate) month: u32,
109 pub(crate) day: u32,
110 pub(crate) hour: u32,
111 pub(crate) minute: u32,
112 pub(crate) second: u32,
113 pub(crate) microsecond: u32,
114}
115
116#[derive(Debug, Clone, Copy, PartialEq, Eq)]
117pub(crate) struct ValidCivilFemtosecond {
118 pub(crate) year: i64,
119 pub(crate) month: u32,
120 pub(crate) day: u32,
121 pub(crate) hour: u32,
122 pub(crate) minute: u32,
123 pub(crate) second: u32,
124 pub(crate) microsecond: u32,
126 pub(crate) femtosecond: u32,
129}
130
131#[derive(Debug, Clone, Copy, PartialEq, Eq)]
132pub(crate) enum CivilSecondPolicy {
133 UtcLike,
135 Continuous,
138}
139
140#[derive(Debug, Clone, Copy)]
141struct CivilMinute {
142 year: i64,
143 month: i64,
144 day: i64,
145 hour: i64,
146 minute: i64,
147}
148
149const CIVIL_YEAR_MIN: i64 = 0;
150const CIVIL_YEAR_MAX: i64 = 9999;
151
152const MICROSECOND_FRACTION_DIGITS: usize = 6;
153const FEMTOSECOND_FRACTION_DIGITS: usize = 15;
154const FEMTOSECONDS_PER_SECOND: i128 = 1_000_000_000_000_000;
155const FEMTOSECONDS_PER_MICROSECOND: i128 = 1_000_000_000;
156
157impl CivilSecondPolicy {
158 const fn allows_leap_second_label(self) -> bool {
159 match self {
160 Self::UtcLike => true,
161 Self::Continuous => false,
162 }
163 }
164}
165
166pub(crate) fn finite(x: f64, field: &'static str) -> Result<f64, FieldError> {
167 if x.is_finite() {
168 Ok(x)
169 } else {
170 Err(FieldError::NonFinite { field })
171 }
172}
173
174pub(crate) fn finite_positive(x: f64, field: &'static str) -> Result<f64, FieldError> {
175 finite(x, field).and_then(|x| {
176 if x > 0.0 {
177 Ok(x)
178 } else {
179 Err(FieldError::NotPositive { field })
180 }
181 })
182}
183
184pub(crate) fn finite_nonneg(x: f64, field: &'static str) -> Result<f64, FieldError> {
185 finite(x, field).and_then(|x| {
186 if x >= 0.0 {
187 Ok(x)
188 } else {
189 Err(FieldError::Negative { field })
190 }
191 })
192}
193
194pub(crate) fn finite_in_range(
195 x: f64,
196 min: f64,
197 max: f64,
198 field: &'static str,
199) -> Result<f64, FieldError> {
200 finite_in_range_impl(x, min, max, true, field)
201}
202
203pub(crate) fn finite_in_range_exclusive_upper(
204 x: f64,
205 min: f64,
206 max: f64,
207 field: &'static str,
208) -> Result<f64, FieldError> {
209 finite_in_range_impl(x, min, max, false, field)
210}
211
212pub(crate) fn fraction(x: f64, field: &'static str) -> Result<f64, FieldError> {
213 finite_in_range(x, 0.0, 1.0, field)
214}
215
216pub(crate) fn second_of_day(x: f64, field: &'static str) -> Result<f64, FieldError> {
217 finite_in_range_exclusive_upper(x, 0.0, crate::constants::SECONDS_PER_DAY, field)
218}
219
220pub(crate) fn positive_step(x: f64, field: &'static str) -> Result<f64, FieldError> {
221 finite_positive(x, field)
222}
223
224pub(crate) fn range_order(lo: f64, hi: f64, field: &'static str) -> Result<(), FieldError> {
225 if lo <= hi {
226 Ok(())
227 } else {
228 Err(FieldError::OutOfRange {
229 field,
230 min: lo,
231 max: hi,
232 upper_inclusive: true,
233 })
234 }
235}
236
237pub(crate) fn clamp_magnitude(x: f64, max_magnitude: f64) -> f64 {
238 debug_assert!(max_magnitude.is_finite());
239 debug_assert!(max_magnitude > 0.0);
240
241 x.clamp(-max_magnitude, max_magnitude)
242}
243
244fn finite_in_range_impl(
245 x: f64,
246 min: f64,
247 max: f64,
248 upper_inclusive: bool,
249 field: &'static str,
250) -> Result<f64, FieldError> {
251 debug_assert!(min.is_finite());
252 debug_assert!(max.is_finite());
253 debug_assert!(min <= max);
254
255 let x = finite(x, field)?;
256 let upper_ok = if upper_inclusive { x <= max } else { x < max };
257 if x >= min && upper_ok {
258 Ok(x)
259 } else {
260 Err(FieldError::OutOfRange {
261 field,
262 min,
263 max,
264 upper_inclusive,
265 })
266 }
267}
268
269pub(crate) fn finite_vec3(v: [f64; 3], field: &'static str) -> Result<[f64; 3], FieldError> {
270 finite_slice(&v, field).map(|()| v)
271}
272
273pub(crate) fn finite_slice(xs: &[f64], field: &'static str) -> Result<(), FieldError> {
274 if xs.iter().all(|x| x.is_finite()) {
275 Ok(())
276 } else {
277 Err(FieldError::NonFinite { field })
278 }
279}
280
281pub(crate) fn require_strictly_increasing<I>(
282 values: I,
283 field: &'static str,
284) -> Result<(), FieldError>
285where
286 I: IntoIterator<Item = f64>,
287{
288 let mut previous = None;
289 for value in values {
290 finite(value, field)?;
291 if let Some(previous) = previous {
292 if value <= previous {
293 return Err(FieldError::OutOfRange {
294 field,
295 min: previous,
296 max: value,
297 upper_inclusive: false,
298 });
299 }
300 }
301 previous = Some(value);
302 }
303 Ok(())
304}
305
306#[allow(clippy::needless_range_loop)]
307pub(crate) fn validate_covariance_psd<const N: usize>(
308 m: &[[f64; N]; N],
309 field: &'static str,
310) -> Result<(), FieldError> {
311 for row in m {
312 finite_slice(row, field)?;
313 }
314
315 let scale = matrix_scale(m);
316 let tol = covariance_matrix_tolerance(N, scale);
317 for i in 0..N {
318 for j in (i + 1)..N {
319 if (m[i][j] - m[j][i]).abs() > tol {
320 return Err(FieldError::NotPositive { field });
321 }
322 }
323 }
324
325 let mut symmetric_part = [[0.0_f64; N]; N];
326 for i in 0..N {
327 symmetric_part[i][i] = m[i][i];
328 for j in (i + 1)..N {
329 let value = 0.5 * (m[i][j] + m[j][i]);
330 symmetric_part[i][j] = value;
331 symmetric_part[j][i] = value;
332 }
333 }
334
335 if symmetric_min_eigenvalue(&mut symmetric_part, tol) < -tol {
336 return Err(FieldError::NotPositive { field });
337 }
338
339 Ok(())
340}
341
342#[allow(clippy::needless_range_loop)]
343pub(crate) fn validate_covariance_psd_rows(
344 rows: &[&[f64]],
345 field: &'static str,
346) -> Result<(), FieldError> {
347 let n = rows.len();
348 for row in rows {
349 finite_slice(row, field)?;
350 debug_assert_eq!(row.len(), n);
351 }
352
353 let scale = matrix_rows_scale(rows);
354 let tol = covariance_matrix_tolerance(n, scale);
355 for i in 0..n {
356 for j in (i + 1)..n {
357 if (rows[i][j] - rows[j][i]).abs() > tol {
358 return Err(FieldError::NotPositive { field });
359 }
360 }
361 }
362
363 let mut symmetric_part = vec![vec![0.0_f64; n]; n];
364 for i in 0..n {
365 symmetric_part[i][i] = rows[i][i];
366 for j in (i + 1)..n {
367 let value = 0.5 * (rows[i][j] + rows[j][i]);
368 symmetric_part[i][j] = value;
369 symmetric_part[j][i] = value;
370 }
371 }
372
373 if symmetric_rows_min_eigenvalue(&mut symmetric_part, tol) < -tol {
374 return Err(FieldError::NotPositive { field });
375 }
376
377 Ok(())
378}
379
380fn matrix_scale<const N: usize>(m: &[[f64; N]; N]) -> f64 {
381 let mut scale = 1.0_f64;
382 for row in m {
383 for value in row {
384 scale = scale.max(value.abs());
385 }
386 }
387 scale
388}
389
390fn matrix_rows_scale(rows: &[&[f64]]) -> f64 {
391 let mut scale = 1.0_f64;
392 for row in rows {
393 for value in *row {
394 scale = scale.max(value.abs());
395 }
396 }
397 scale
398}
399
400fn covariance_matrix_tolerance(n: usize, scale: f64) -> f64 {
401 let scale = scale.max(1.0);
402 (128.0 * f64::EPSILON * (n.max(1) as f64) * scale).max(1.0e-9 * scale)
403}
404
405#[allow(clippy::needless_range_loop)]
406fn symmetric_min_eigenvalue<const N: usize>(a: &mut [[f64; N]; N], tol: f64) -> f64 {
407 let max_sweeps = (16 * N * N).max(32);
408 for _ in 0..max_sweeps {
409 let mut p = 0usize;
410 let mut q = 0usize;
411 let mut max_offdiag = 0.0_f64;
412 for i in 0..N {
413 for j in (i + 1)..N {
414 let offdiag = a[i][j].abs();
415 if offdiag > max_offdiag {
416 max_offdiag = offdiag;
417 p = i;
418 q = j;
419 }
420 }
421 }
422
423 if max_offdiag <= tol {
424 break;
425 }
426
427 let app = a[p][p];
428 let aqq = a[q][q];
429 let apq = a[p][q];
430 if apq == 0.0 {
431 break;
432 }
433
434 let tau = (aqq - app) / (2.0 * apq);
435 let t = if tau >= 0.0 {
436 1.0 / (tau + (1.0 + tau * tau).sqrt())
437 } else {
438 -1.0 / (-tau + (1.0 + tau * tau).sqrt())
439 };
440 let c = 1.0 / (1.0 + t * t).sqrt();
441 let s = t * c;
442
443 for k in 0..N {
444 if k != p && k != q {
445 let akp = a[k][p];
446 let akq = a[k][q];
447 let new_kp = c * akp - s * akq;
448 let new_kq = s * akp + c * akq;
449 a[k][p] = new_kp;
450 a[p][k] = new_kp;
451 a[k][q] = new_kq;
452 a[q][k] = new_kq;
453 }
454 }
455
456 a[p][p] = c * c * app - 2.0 * s * c * apq + s * s * aqq;
457 a[q][q] = s * s * app + 2.0 * s * c * apq + c * c * aqq;
458 a[p][q] = 0.0;
459 a[q][p] = 0.0;
460 }
461
462 let mut min = f64::INFINITY;
463 for (i, row) in a.iter().enumerate() {
464 min = min.min(row[i]);
465 }
466 min
467}
468
469#[allow(clippy::needless_range_loop)]
470fn symmetric_rows_min_eigenvalue(a: &mut [Vec<f64>], tol: f64) -> f64 {
471 let n = a.len();
472 let max_sweeps = (16 * n * n).max(32);
473 for _ in 0..max_sweeps {
474 let mut p = 0usize;
475 let mut q = 0usize;
476 let mut max_offdiag = 0.0_f64;
477 for (i, row) in a.iter().enumerate() {
478 for (j, value) in row.iter().enumerate().skip(i + 1) {
479 let offdiag = value.abs();
480 if offdiag > max_offdiag {
481 max_offdiag = offdiag;
482 p = i;
483 q = j;
484 }
485 }
486 }
487
488 if max_offdiag <= tol {
489 break;
490 }
491
492 let app = a[p][p];
493 let aqq = a[q][q];
494 let apq = a[p][q];
495 if apq == 0.0 {
496 break;
497 }
498
499 let tau = (aqq - app) / (2.0 * apq);
500 let t = if tau >= 0.0 {
501 1.0 / (tau + (1.0 + tau * tau).sqrt())
502 } else {
503 -1.0 / (-tau + (1.0 + tau * tau).sqrt())
504 };
505 let c = 1.0 / (1.0 + t * t).sqrt();
506 let s = t * c;
507
508 for k in 0..n {
509 if k != p && k != q {
510 let akp = a[k][p];
511 let akq = a[k][q];
512 let new_kp = c * akp - s * akq;
513 let new_kq = s * akp + c * akq;
514 a[k][p] = new_kp;
515 a[p][k] = new_kp;
516 a[k][q] = new_kq;
517 a[q][k] = new_kq;
518 }
519 }
520
521 a[p][p] = c * c * app - 2.0 * s * c * apq + s * s * aqq;
522 a[q][q] = s * s * app + 2.0 * s * c * apq + c * c * aqq;
523 a[p][q] = 0.0;
524 a[q][p] = 0.0;
525 }
526
527 let mut min = f64::INFINITY;
528 for (i, row) in a.iter().enumerate() {
529 min = min.min(row[i]);
530 }
531 min
532}
533
534pub(crate) fn present<T>(value: Option<T>, field: &'static str) -> Result<T, FieldError> {
535 value.ok_or(FieldError::Missing { field })
536}
537
538pub(crate) fn exact_len<T>(
539 xs: &[T],
540 expected: usize,
541 field: &'static str,
542) -> Result<(), LengthError> {
543 if xs.len() == expected {
544 Ok(())
545 } else {
546 Err(LengthError {
547 field,
548 expected,
549 actual: xs.len(),
550 })
551 }
552}
553
554pub(crate) fn checked_i64_add(
555 lhs: i64,
556 rhs: i64,
557 field: &'static str,
558) -> Result<i64, ArithmeticError> {
559 lhs.checked_add(rhs).ok_or(ArithmeticError { field })
560}
561
562pub(crate) fn checked_i64_sub(
563 lhs: i64,
564 rhs: i64,
565 field: &'static str,
566) -> Result<i64, ArithmeticError> {
567 lhs.checked_sub(rhs).ok_or(ArithmeticError { field })
568}
569
570pub(crate) fn checked_i64_mul(
571 lhs: i64,
572 rhs: i64,
573 field: &'static str,
574) -> Result<i64, ArithmeticError> {
575 lhs.checked_mul(rhs).ok_or(ArithmeticError { field })
576}
577
578pub(crate) fn strict_f64(s: &str, field: &'static str) -> Result<f64, FieldError> {
579 let value = s.trim();
580 if value.is_empty() {
581 return Err(FieldError::Missing { field });
582 }
583 let normalized = value.replace(['D', 'd'], "e");
584 let parsed = normalized
585 .parse::<f64>()
586 .map_err(|_| FieldError::FloatParse {
587 field,
588 value: value.to_string(),
589 })?;
590 finite(parsed, field)
591}
592
593pub(crate) fn strict_int<T>(s: &str, field: &'static str) -> Result<T, FieldError>
594where
595 T: FromStr,
596{
597 let value = s.trim();
598 if value.is_empty() {
599 return Err(FieldError::Missing { field });
600 }
601 value.parse::<T>().map_err(|_| FieldError::IntParse {
602 field,
603 value: value.to_string(),
604 })
605}
606
607pub(crate) fn strict_gnss_satellite_id(
608 s: &str,
609 field: &'static str,
610) -> Result<GnssSatelliteId, FieldError> {
611 let value = s.trim();
612 if value.is_empty() {
613 return Err(FieldError::Missing { field });
614 }
615 value
616 .parse::<GnssSatelliteId>()
617 .map_err(|_| FieldError::IntParse {
618 field,
619 value: value.to_string(),
620 })
621}
622
623pub(crate) fn civil_datetime_with_second_policy(
624 year: i64,
625 month: i64,
626 day: i64,
627 hour: i64,
628 minute: i64,
629 second: f64,
630 second_policy: CivilSecondPolicy,
631) -> Result<ValidCivil, FieldError> {
632 const FIELD: &str = "civil datetime";
633 finite(second, FIELD)?;
634
635 validate_civil_year(year, month, day)?;
636 if !(1..=12).contains(&month) {
637 return Err(FieldError::InvalidCivilDate {
638 field: FIELD,
639 year,
640 month,
641 day,
642 });
643 }
644 let last_day = crate::astro::time::civil::days_in_month(year, month);
645 if !(1..=last_day).contains(&day) {
646 return Err(FieldError::InvalidCivilDate {
647 field: FIELD,
648 year,
649 month,
650 day,
651 });
652 }
653 let time_fields_valid = (0..=23).contains(&hour) && (0..=59).contains(&minute);
654 let seconds_valid = (0.0..60.0).contains(&second)
655 || (second_policy.allows_leap_second_label()
656 && (60.0..61.0).contains(&second)
657 && time_fields_valid
658 && is_positive_utc_leap_second_label(year, month, day, hour, minute));
659 if !time_fields_valid || !seconds_valid {
660 return Err(FieldError::InvalidCivilTime {
661 field: FIELD,
662 hour,
663 minute,
664 second,
665 });
666 }
667
668 Ok(ValidCivil {
669 year,
670 month: month as u32,
671 day: day as u32,
672 hour: hour as u32,
673 minute: minute as u32,
674 second,
675 })
676}
677
678fn validate_civil_year(year: i64, month: i64, day: i64) -> Result<(), FieldError> {
679 if (CIVIL_YEAR_MIN..=CIVIL_YEAR_MAX).contains(&year) {
680 Ok(())
681 } else {
682 Err(FieldError::InvalidCivilDate {
683 field: "civil datetime",
684 year,
685 month,
686 day,
687 })
688 }
689}
690
691pub(crate) fn civil_datetime_with_decimal_second_policy(
692 year: i64,
693 month: i64,
694 day: i64,
695 hour: i64,
696 minute: i64,
697 second: &str,
698 second_policy: CivilSecondPolicy,
699) -> Result<ValidCivilMicrosecond, FieldError> {
700 const FIELD: &str = "civil datetime";
701 let second = second.trim();
702 if second.is_empty() {
703 return Err(FieldError::Missing { field: FIELD });
704 }
705
706 let (whole_second, microsecond) =
707 decimal_second_to_subseconds(second, FIELD, MICROSECOND_FRACTION_DIGITS)?;
708 civil_datetime_with_whole_microsecond_policy(
709 CivilMinute {
710 year,
711 month,
712 day,
713 hour,
714 minute,
715 },
716 whole_second,
717 microsecond as i64,
718 second_policy,
719 )
720}
721
722pub(crate) fn civil_datetime_with_femtosecond_policy(
728 year: i64,
729 month: i64,
730 day: i64,
731 hour: i64,
732 minute: i64,
733 second: &str,
734 second_policy: CivilSecondPolicy,
735) -> Result<ValidCivilFemtosecond, FieldError> {
736 const FIELD: &str = "civil datetime";
737 let second = second.trim();
738 if second.is_empty() {
739 return Err(FieldError::Missing { field: FIELD });
740 }
741
742 let (whole_second, femtosecond) =
743 decimal_second_to_subseconds(second, FIELD, FEMTOSECOND_FRACTION_DIGITS)?;
744 civil_datetime_with_whole_femtosecond_policy(
745 CivilMinute {
746 year,
747 month,
748 day,
749 hour,
750 minute,
751 },
752 whole_second,
753 femtosecond,
754 second_policy,
755 )
756}
757
758pub(crate) fn civil_datetime_with_fractional_second_policy(
759 year: i64,
760 month: i64,
761 day: i64,
762 hour: i64,
763 minute: i64,
764 second: f64,
765 second_policy: CivilSecondPolicy,
766) -> Result<ValidCivilMicrosecond, FieldError> {
767 let civil =
768 civil_datetime_with_second_policy(year, month, day, hour, minute, second, second_policy)?;
769
770 let whole_second = civil.second.trunc() as i64;
771 let microsecond = ((civil.second - whole_second as f64) * 1_000_000.0).round() as i64;
772 civil_datetime_with_whole_microsecond_policy(
773 CivilMinute {
774 year: civil.year,
775 month: i64::from(civil.month),
776 day: i64::from(civil.day),
777 hour: i64::from(civil.hour),
778 minute: i64::from(civil.minute),
779 },
780 whole_second,
781 microsecond,
782 second_policy,
783 )
784}
785
786fn civil_datetime_with_whole_microsecond_policy(
787 minute_parts: CivilMinute,
788 whole_second: i64,
789 microsecond: i64,
790 second_policy: CivilSecondPolicy,
791) -> Result<ValidCivilMicrosecond, FieldError> {
792 let civil = civil_datetime_with_whole_femtosecond_policy(
793 minute_parts,
794 whole_second,
795 i128::from(microsecond) * FEMTOSECONDS_PER_MICROSECOND,
796 second_policy,
797 )?;
798 Ok(ValidCivilMicrosecond {
799 year: civil.year,
800 month: civil.month,
801 day: civil.day,
802 hour: civil.hour,
803 minute: civil.minute,
804 second: civil.second,
805 microsecond: civil.microsecond,
806 })
807}
808
809fn civil_datetime_with_whole_femtosecond_policy(
810 minute_parts: CivilMinute,
811 whole_second: i64,
812 mut femtosecond: i128,
813 second_policy: CivilSecondPolicy,
814) -> Result<ValidCivilFemtosecond, FieldError> {
815 const FIELD: &str = "civil datetime";
816 if !(0..=FEMTOSECONDS_PER_SECOND).contains(&femtosecond) {
817 return Err(FieldError::InvalidCivilTime {
818 field: FIELD,
819 hour: minute_parts.hour,
820 minute: minute_parts.minute,
821 second: whole_second as f64 + femtosecond as f64 / FEMTOSECONDS_PER_SECOND as f64,
822 });
823 }
824
825 let civil = civil_datetime_with_second_policy(
826 minute_parts.year,
827 minute_parts.month,
828 minute_parts.day,
829 minute_parts.hour,
830 minute_parts.minute,
831 whole_second as f64,
832 second_policy,
833 )?;
834
835 let mut rounded_second = whole_second;
836 let rounded_from_subsecond = femtosecond == FEMTOSECONDS_PER_SECOND;
837 if rounded_from_subsecond {
838 rounded_second += 1;
839 femtosecond = 0;
840 }
841 let microsecond = (femtosecond / FEMTOSECONDS_PER_MICROSECOND) as u32;
842 let femtosecond = (femtosecond % FEMTOSECONDS_PER_MICROSECOND) as u32;
843
844 let leap_second_label_allowed = second_policy.allows_leap_second_label()
845 && is_positive_utc_leap_second_label(
846 minute_parts.year,
847 minute_parts.month,
848 minute_parts.day,
849 minute_parts.hour,
850 minute_parts.minute,
851 );
852 if rounded_second <= 59 || (rounded_second == 60 && leap_second_label_allowed) {
853 return Ok(ValidCivilFemtosecond {
854 year: civil.year,
855 month: civil.month,
856 day: civil.day,
857 hour: civil.hour,
858 minute: civil.minute,
859 second: rounded_second as u32,
860 microsecond,
861 femtosecond,
862 });
863 }
864
865 if rounded_second == 60
866 || (rounded_second == 61 && rounded_from_subsecond && leap_second_label_allowed)
867 {
868 let (year, month, day, hour, minute) =
869 carry_to_next_minute(civil.year, civil.month, civil.day, civil.hour, civil.minute)?;
870 return Ok(ValidCivilFemtosecond {
871 year,
872 month,
873 day,
874 hour,
875 minute,
876 second: 0,
877 microsecond,
878 femtosecond,
879 });
880 }
881
882 Err(FieldError::InvalidCivilTime {
883 field: FIELD,
884 hour: minute_parts.hour,
885 minute: minute_parts.minute,
886 second: rounded_second as f64
887 + (i128::from(microsecond) * FEMTOSECONDS_PER_MICROSECOND + i128::from(femtosecond))
888 as f64
889 / FEMTOSECONDS_PER_SECOND as f64,
890 })
891}
892
893fn is_positive_utc_leap_second_label(
894 year: i64,
895 month: i64,
896 day: i64,
897 hour: i64,
898 minute: i64,
899) -> bool {
900 let (Ok(year), Ok(month), Ok(day), Ok(hour), Ok(minute)) = (
901 i32::try_from(year),
902 i32::try_from(month),
903 i32::try_from(day),
904 i32::try_from(hour),
905 i32::try_from(minute),
906 ) else {
907 return false;
908 };
909 crate::astro::time::scales::is_positive_leap_second_label(year, month, day, hour, minute)
910}
911
912fn decimal_second_to_subseconds(
917 second: &str,
918 field: &'static str,
919 fraction_digits: usize,
920) -> Result<(i64, i128), FieldError> {
921 let (whole, fraction) = second
922 .split_once('.')
923 .map_or((second, None), |(whole, fraction)| (whole, Some(fraction)));
924 if whole.is_empty() {
925 return Err(FieldError::FloatParse {
926 field,
927 value: second.to_string(),
928 });
929 }
930 let negative_zero_whole = whole.starts_with('-');
931 let whole = whole.parse::<i64>().map_err(|_| FieldError::FloatParse {
932 field,
933 value: second.to_string(),
934 })?;
935 let negative_zero_whole = negative_zero_whole && whole == 0;
936
937 let Some(fraction) = fraction else {
938 return Ok((whole, 0));
939 };
940 if fraction.is_empty() || !fraction.bytes().all(|b| b.is_ascii_digit()) {
941 return Err(FieldError::FloatParse {
942 field,
943 value: second.to_string(),
944 });
945 }
946
947 let mut subsecond = 0i128;
948 for i in 0..fraction_digits {
949 subsecond *= 10;
950 subsecond += fraction
951 .as_bytes()
952 .get(i)
953 .map_or(0, |b| i128::from(b - b'0'));
954 }
955 if fraction
956 .as_bytes()
957 .get(fraction_digits)
958 .is_some_and(|b| b - b'0' >= 5)
959 {
960 subsecond += 1;
961 }
962 if negative_zero_whole {
963 subsecond = -subsecond.max(1);
964 }
965 Ok((whole, subsecond))
966}
967
968fn carry_to_next_minute(
969 mut year: i64,
970 mut month: u32,
971 mut day: u32,
972 mut hour: u32,
973 mut minute: u32,
974) -> Result<(i64, u32, u32, u32, u32), FieldError> {
975 minute += 1;
976 if minute < 60 {
977 return Ok((year, month, day, hour, minute));
978 }
979 minute = 0;
980 hour += 1;
981 if hour < 24 {
982 return Ok((year, month, day, hour, minute));
983 }
984 hour = 0;
985 day += 1;
986 if i64::from(day) <= crate::astro::time::civil::days_in_month(year, i64::from(month)) {
987 return Ok((year, month, day, hour, minute));
988 }
989 day = 1;
990 month += 1;
991 if month <= 12 {
992 return Ok((year, month, day, hour, minute));
993 }
994 month = 1;
995 year = year
996 .checked_add(1)
997 .ok_or_else(|| FieldError::InvalidCivilDate {
998 field: "civil datetime",
999 year,
1000 month: i64::from(month),
1001 day: i64::from(day),
1002 })?;
1003 validate_civil_year(year, i64::from(month), i64::from(day))?;
1004 Ok((year, month, day, hour, minute))
1005}
1006
1007#[cfg(test)]
1008mod tests {
1009 use super::*;
1010
1011 #[test]
1012 fn bounded_quantity_helpers_accept_valid_inputs() {
1013 assert_eq!(finite_in_range(1.0, -1.0, 1.0, "bounded"), Ok(1.0));
1014 assert_eq!(fraction(0.0, "fraction"), Ok(0.0));
1015 assert_eq!(fraction(1.0, "fraction"), Ok(1.0));
1016 assert_eq!(second_of_day(86_399.999, "second_of_day"), Ok(86_399.999));
1017 assert_eq!(positive_step(0.25, "step"), Ok(0.25));
1018 assert_eq!(range_order(-1.0, 1.0, "range"), Ok(()));
1019 assert_eq!(range_order(1.0, 1.0, "range"), Ok(()));
1020 assert_eq!(clamp_magnitude(3.0, 2.0), 2.0);
1021 assert_eq!(clamp_magnitude(-3.0, 2.0), -2.0);
1022 assert_eq!(clamp_magnitude(1.5, 2.0), 1.5);
1023 }
1024
1025 #[test]
1026 fn bounded_quantity_helpers_reject_bad_inputs() {
1027 assert!(matches!(
1028 fraction(50.0, "fraction"),
1029 Err(FieldError::OutOfRange {
1030 field: "fraction",
1031 min: 0.0,
1032 max: 1.0,
1033 upper_inclusive: true
1034 })
1035 ));
1036 assert!(matches!(
1037 second_of_day(crate::constants::SECONDS_PER_DAY, "second_of_day"),
1038 Err(FieldError::OutOfRange {
1039 field: "second_of_day",
1040 min: 0.0,
1041 max: crate::constants::SECONDS_PER_DAY,
1042 upper_inclusive: false
1043 })
1044 ));
1045 assert!(matches!(
1046 positive_step(0.0, "step"),
1047 Err(FieldError::NotPositive { field: "step" })
1048 ));
1049 assert!(matches!(
1050 range_order(2.0, 1.0, "range"),
1051 Err(FieldError::OutOfRange {
1052 field: "range",
1053 min: 2.0,
1054 max: 1.0,
1055 upper_inclusive: true
1056 })
1057 ));
1058 assert!(matches!(
1059 finite_in_range(f64::NAN, -1.0, 1.0, "bounded"),
1060 Err(FieldError::NonFinite { field: "bounded" })
1061 ));
1062 }
1063
1064 #[test]
1065 fn covariance_psd_helper_accepts_positive_semidefinite_matrices() {
1066 let semidefinite = [[1.0, 1.0], [1.0, 1.0]];
1067 assert_eq!(validate_covariance_psd(&semidefinite, "covariance"), Ok(()));
1068
1069 let scaled = [
1070 [4.0e12, 2.0e12, 0.0],
1071 [2.0e12, 1.0e12, 0.0],
1072 [0.0, 0.0, 3.0],
1073 ];
1074 assert_eq!(validate_covariance_psd(&scaled, "covariance"), Ok(()));
1075 }
1076
1077 #[test]
1078 fn covariance_psd_helper_rejects_malformed_matrices() {
1079 let non_finite = [[1.0, f64::NAN], [f64::NAN, 1.0]];
1080 assert!(matches!(
1081 validate_covariance_psd(&non_finite, "covariance"),
1082 Err(FieldError::NonFinite {
1083 field: "covariance"
1084 })
1085 ));
1086
1087 let asymmetric = [[1.0, 1.0e-5], [0.0, 1.0]];
1088 assert!(matches!(
1089 validate_covariance_psd(&asymmetric, "covariance"),
1090 Err(FieldError::NotPositive {
1091 field: "covariance"
1092 })
1093 ));
1094
1095 let indefinite = [[1.0, 2.0], [2.0, 1.0]];
1096 assert!(matches!(
1097 validate_covariance_psd(&indefinite, "covariance"),
1098 Err(FieldError::NotPositive {
1099 field: "covariance"
1100 })
1101 ));
1102 }
1103
1104 #[test]
1105 fn covariance_psd_rows_helper_matches_array_helper() {
1106 let valid = [vec![2.0, 0.25], vec![0.25, 1.0]];
1107 let valid_rows: Vec<&[f64]> = valid.iter().map(Vec::as_slice).collect();
1108 assert_eq!(
1109 validate_covariance_psd_rows(&valid_rows, "covariance"),
1110 Ok(())
1111 );
1112
1113 let indefinite = [vec![1.0, 2.0], vec![2.0, 1.0]];
1114 let indefinite_rows: Vec<&[f64]> = indefinite.iter().map(Vec::as_slice).collect();
1115 assert!(matches!(
1116 validate_covariance_psd_rows(&indefinite_rows, "covariance"),
1117 Err(FieldError::NotPositive {
1118 field: "covariance"
1119 })
1120 ));
1121 }
1122
1123 #[test]
1124 fn strictly_increasing_helper_rejects_nonfinite_and_nonmonotonic_values() {
1125 assert_eq!(require_strictly_increasing([1.0, 2.0, 3.0], "time"), Ok(()));
1126 assert!(matches!(
1127 require_strictly_increasing([1.0, 1.0], "time"),
1128 Err(FieldError::OutOfRange {
1129 field: "time",
1130 min: 1.0,
1131 max: 1.0,
1132 upper_inclusive: false
1133 })
1134 ));
1135 assert!(matches!(
1136 require_strictly_increasing([1.0, f64::NAN], "time"),
1137 Err(FieldError::NonFinite { field: "time" })
1138 ));
1139 }
1140
1141 #[test]
1142 fn civil_datetime_uses_time_system_leap_second_policy() {
1143 let utc = civil_datetime_with_second_policy(
1144 2016,
1145 12,
1146 31,
1147 23,
1148 59,
1149 60.0,
1150 CivilSecondPolicy::UtcLike,
1151 )
1152 .expect("UTC-like civil time accepts a leap second label");
1153 assert_eq!(utc.second, 60.0);
1154
1155 let utc_ordinary = civil_datetime_with_second_policy(
1156 2016,
1157 12,
1158 30,
1159 23,
1160 59,
1161 59.0,
1162 CivilSecondPolicy::UtcLike,
1163 )
1164 .expect("UTC-like civil time accepts an ordinary final minute second");
1165 assert_eq!(utc_ordinary.second, 59.0);
1166
1167 assert!(matches!(
1168 civil_datetime_with_second_policy(
1169 2016,
1170 12,
1171 30,
1172 23,
1173 59,
1174 60.0,
1175 CivilSecondPolicy::UtcLike
1176 ),
1177 Err(FieldError::InvalidCivilTime {
1178 field: "civil datetime",
1179 hour: 23,
1180 minute: 59,
1181 second: 60.0
1182 })
1183 ));
1184
1185 let gps = civil_datetime_with_second_policy(
1186 2016,
1187 12,
1188 31,
1189 23,
1190 59,
1191 59.0,
1192 CivilSecondPolicy::Continuous,
1193 )
1194 .expect("GPS-like civil time accepts an ordinary final minute second");
1195 assert_eq!(gps.second, 59.0);
1196
1197 assert!(matches!(
1198 civil_datetime_with_second_policy(
1199 2016,
1200 12,
1201 31,
1202 23,
1203 59,
1204 60.0,
1205 CivilSecondPolicy::Continuous
1206 ),
1207 Err(FieldError::InvalidCivilTime {
1208 field: "civil datetime",
1209 hour: 23,
1210 minute: 59,
1211 second: 60.0
1212 })
1213 ));
1214 assert!(matches!(
1215 civil_datetime_with_second_policy(
1216 2016,
1217 12,
1218 31,
1219 23,
1220 59,
1221 61.0,
1222 CivilSecondPolicy::UtcLike
1223 ),
1224 Err(FieldError::InvalidCivilTime {
1225 field: "civil datetime",
1226 hour: 23,
1227 minute: 59,
1228 second: 61.0
1229 })
1230 ));
1231 }
1232
1233 #[test]
1234 fn decimal_second_parser_carries_rounded_microseconds() {
1235 let ordinary = civil_datetime_with_decimal_second_policy(
1236 2026,
1237 6,
1238 17,
1239 4,
1240 32,
1241 "52.9999995",
1242 CivilSecondPolicy::Continuous,
1243 )
1244 .expect("rounded fractional second carries into second");
1245 assert_eq!(
1246 ordinary,
1247 ValidCivilMicrosecond {
1248 year: 2026,
1249 month: 6,
1250 day: 17,
1251 hour: 4,
1252 minute: 32,
1253 second: 53,
1254 microsecond: 0,
1255 }
1256 );
1257
1258 let day_boundary = civil_datetime_with_decimal_second_policy(
1259 2026,
1260 6,
1261 17,
1262 23,
1263 59,
1264 "59.9999995",
1265 CivilSecondPolicy::Continuous,
1266 )
1267 .expect("rounded fractional second carries across day boundary");
1268 assert_eq!(
1269 day_boundary,
1270 ValidCivilMicrosecond {
1271 year: 2026,
1272 month: 6,
1273 day: 18,
1274 hour: 0,
1275 minute: 0,
1276 second: 0,
1277 microsecond: 0,
1278 }
1279 );
1280
1281 let utc_ordinary_boundary = civil_datetime_with_decimal_second_policy(
1282 2026,
1283 6,
1284 17,
1285 23,
1286 59,
1287 "59.9999995",
1288 CivilSecondPolicy::UtcLike,
1289 )
1290 .expect("rounded UTC ordinary second carries across day boundary");
1291 assert_eq!(
1292 utc_ordinary_boundary,
1293 ValidCivilMicrosecond {
1294 year: 2026,
1295 month: 6,
1296 day: 18,
1297 hour: 0,
1298 minute: 0,
1299 second: 0,
1300 microsecond: 0,
1301 }
1302 );
1303
1304 let utc_leap_boundary = civil_datetime_with_decimal_second_policy(
1305 2016,
1306 12,
1307 31,
1308 23,
1309 59,
1310 "59.9999995",
1311 CivilSecondPolicy::UtcLike,
1312 )
1313 .expect("rounded UTC leap-second label stays in the leap minute");
1314 assert_eq!(
1315 utc_leap_boundary,
1316 ValidCivilMicrosecond {
1317 year: 2016,
1318 month: 12,
1319 day: 31,
1320 hour: 23,
1321 minute: 59,
1322 second: 60,
1323 microsecond: 0,
1324 }
1325 );
1326 }
1327
1328 #[test]
1329 fn decimal_second_parser_rejects_bad_fraction() {
1330 assert!(matches!(
1331 civil_datetime_with_decimal_second_policy(
1332 2026,
1333 6,
1334 17,
1335 4,
1336 32,
1337 "52.x",
1338 CivilSecondPolicy::Continuous
1339 ),
1340 Err(FieldError::FloatParse {
1341 field: "civil datetime",
1342 ..
1343 })
1344 ));
1345 }
1346
1347 #[test]
1348 fn decimal_second_parser_rejects_negative_fractional_zero_second() {
1349 assert!(matches!(
1350 civil_datetime_with_decimal_second_policy(
1351 2026,
1352 6,
1353 17,
1354 4,
1355 32,
1356 "-0.1",
1357 CivilSecondPolicy::Continuous
1358 ),
1359 Err(FieldError::InvalidCivilTime {
1360 field: "civil datetime",
1361 hour: 4,
1362 minute: 32,
1363 second
1364 }) if (second + 0.1).abs() < f64::EPSILON
1365 ));
1366
1367 let fractional_zero_second = civil_datetime_with_decimal_second_policy(
1368 2026,
1369 6,
1370 17,
1371 4,
1372 32,
1373 "0.1",
1374 CivilSecondPolicy::Continuous,
1375 )
1376 .expect("positive fractional zero second parses");
1377 assert_eq!(
1378 fractional_zero_second,
1379 ValidCivilMicrosecond {
1380 year: 2026,
1381 month: 6,
1382 day: 17,
1383 hour: 4,
1384 minute: 32,
1385 second: 0,
1386 microsecond: 100_000,
1387 }
1388 );
1389
1390 let positive_second = civil_datetime_with_decimal_second_policy(
1391 2026,
1392 6,
1393 17,
1394 4,
1395 32,
1396 "52.1",
1397 CivilSecondPolicy::Continuous,
1398 )
1399 .expect("positive decimal second parses");
1400 assert_eq!(
1401 positive_second,
1402 ValidCivilMicrosecond {
1403 year: 2026,
1404 month: 6,
1405 day: 17,
1406 hour: 4,
1407 minute: 32,
1408 second: 52,
1409 microsecond: 100_000,
1410 }
1411 );
1412 }
1413}