oxiproj-transformations 0.1.2

Datum transformations and coordinate conversions for OxiProj.
Documentation
//! Epoch-aware coordinate transformation pipelines.
//!
//! Constructs end-to-end pipelines that automatically chain:
//! 1. PMM propagation from source epoch to target epoch (in the source frame,
//!    on the plate the point physically sits on)
//! 2. A frame-change Helmert **path** (possibly multi-hop, e.g.
//!    ITRF2014 → ITRF2020 → ITRF2000) evaluated at the target epoch
//!
//! A frame change never alters the coordinate epoch — only step 1 (plate
//! motion) advances the point through time. This matches the way PROJ evaluates
//! a time-dependent Helmert at the coordinate's own epoch.

use crate::epoch_transform::{plate_for_coord, plate_for_frame, propagate_epoch};
use crate::frame_chain::{apply_frame_path, find_frame_path};
use oxiproj_core::{
    epoch::{Epoch, RefFrameEpoch},
    plate_motion::{Plate, PlateMotionModel},
    ProjError, ProjResult,
};

/// Epochs closer together than this (in years) need no plate-motion propagation.
const EPOCH_EPS_YEARS: f64 = 1e-9;

/// Result of an epoch-aware pipeline transformation.
#[derive(Debug, Clone)]
pub struct EpochAwareResult {
    pub x: f64,
    pub y: f64,
    pub z: f64,
    /// Source frame name
    pub from_frame: String,
    /// Source epoch
    pub from_epoch: Epoch,
    /// Target frame name
    pub to_frame: String,
    /// Target epoch
    pub to_epoch: Epoch,
}

/// Resolve the tectonic plate for a point, preferring an explicit frame-name
/// hint and otherwise deriving it from the coordinate location.
///
/// Resolution order:
/// 1. [`plate_for_frame`] — an explicit plate token in the frame name
///    (e.g. `"NOAM_2014"`, `"EURAS_plate"`).
/// 2. [`plate_for_coord`] — a point-in-plate polygon lookup on the coordinate.
///
/// Unlike the previous implementation, there is **no** blanket fallback to the
/// Eurasian plate for bare ITRF frames: a mislabelled plate would silently
/// accumulate decimetre-level error over a decade. When neither the name nor
/// the location yields a plate, this returns `None`, and the caller reports
/// [`ProjError::NoOperation`] rather than guessing.
fn resolve_plate(frame_name: &str, x: f64, y: f64, z: f64) -> Option<Plate> {
    if let Some(p) = plate_for_frame(frame_name) {
        return Some(p);
    }
    plate_for_coord(x, y, z)
}

/// Transform geocentric ECEF coordinates from one reference frame and epoch to another.
///
/// Applies the full epoch-aware pipeline:
/// 1. **Plate-motion propagation** of the point from `from_epoch` to `to_epoch`
///    within `from_frame`, using the plate resolved for the coordinate
///    (skipped when the two epochs coincide).
/// 2. **Frame-change Helmert path** from `from_frame` to `to_frame` — a
///    breadth-first multi-hop chain through the bundled catalogue — evaluated at
///    `to_epoch`.
///
/// When `from_frame == to_frame` the frame-change path is empty and only the
/// propagation applies. When `from_epoch == to_epoch` propagation is skipped and
/// the plate is not required (so pure frame changes work anywhere, including
/// over open ocean).
///
/// # Arguments
///
/// * `x`, `y`, `z` — Geocentric Cartesian coordinates in metres at `from_epoch`.
/// * `from_frame` — Source reference frame (e.g. `"ITRF2014"`).
/// * `from_epoch` — Epoch at which the input coordinates are expressed.
/// * `to_frame` — Target reference frame (e.g. `"ITRF2020"`).
/// * `to_epoch` — Desired output epoch.
/// * `model` — Plate-motion model to use for the propagation step.
///
/// # Errors
///
/// Returns [`ProjError::NoOperation`] when:
/// - Propagation is required but the point's plate cannot be resolved (neither
///   a frame-name hint nor a point-in-plate polygon matched), or the resolved
///   plate is absent from `model`.
/// - No path of frame transforms connects `from_frame` and `to_frame`.
#[allow(clippy::too_many_arguments)]
pub fn transform_epoch_aware(
    x: f64,
    y: f64,
    z: f64,
    from_frame: &str,
    from_epoch: Epoch,
    to_frame: &str,
    to_epoch: Epoch,
    model: PlateMotionModel,
) -> ProjResult<EpochAwareResult> {
    // Step 1: plate-motion propagation in the source frame (skipped when the
    // epochs coincide, in which case the plate is irrelevant).
    let (px, py, pz) = if (to_epoch.year - from_epoch.year).abs() > EPOCH_EPS_YEARS {
        let plate = resolve_plate(from_frame, x, y, z).ok_or(ProjError::NoOperation)?;
        propagate_epoch(x, y, z, &from_epoch, &to_epoch, model, plate)?
    } else {
        (x, y, z)
    };

    // Step 2: frame-change Helmert path (multi-hop), evaluated at the target
    // (coordinate) epoch.
    let path = find_frame_path(from_frame, to_frame).ok_or(ProjError::NoOperation)?;
    let (final_x, final_y, final_z) =
        apply_frame_path(&path, px, py, pz, to_epoch.year).map_err(|_| ProjError::NoOperation)?;

    Ok(EpochAwareResult {
        x: final_x,
        y: final_y,
        z: final_z,
        from_frame: from_frame.to_string(),
        from_epoch,
        to_frame: to_frame.to_string(),
        to_epoch,
    })
}

/// Epoch-aware transform driven by two [`RefFrameEpoch`] tags (e.g. parsed from
/// `"ITRF2014@2010"` and `"ITRF2020@2020"`).
///
/// This is the ergonomic entry point for dynamic↔dynamic transforms: it reads
/// the coordinate epochs directly off the frame tags and delegates to
/// [`transform_epoch_aware`]. Both frame tags must carry an epoch (they are
/// *dynamic* frames); a static, epoch-less tag has no coordinate epoch to
/// propagate to/from.
///
/// # Errors
///
/// Returns [`ProjError::IllegalArgValue`] if either frame tag lacks an epoch,
/// and otherwise propagates the errors of [`transform_epoch_aware`].
pub fn transform_epoch_aware_frames(
    from: &RefFrameEpoch,
    to: &RefFrameEpoch,
    x: f64,
    y: f64,
    z: f64,
    model: PlateMotionModel,
) -> ProjResult<EpochAwareResult> {
    let from_epoch = from.epoch.ok_or(ProjError::IllegalArgValue)?;
    let to_epoch = to.epoch.ok_or(ProjError::IllegalArgValue)?;
    transform_epoch_aware(x, y, z, &from.frame, from_epoch, &to.frame, to_epoch, model)
}

#[cfg(test)]
mod tests {
    use super::*;
    use oxiproj_core::epoch::Epoch;

    #[test]
    fn test_itrf2014_to_itrf2020_epoch_shift() {
        let (x0, y0, z0) = (3_711_000.0, 1_023_000.0, 5_036_000.0);
        let result = transform_epoch_aware(
            x0,
            y0,
            z0,
            "ITRF2014",
            Epoch::new(2010.0).unwrap(),
            "ITRF2020",
            Epoch::new(2015.0).unwrap(),
            PlateMotionModel::Itrf2020,
        )
        .expect("should succeed");
        assert!(result.x.is_finite() && result.y.is_finite() && result.z.is_finite());
        let dx = result.x - x0;
        let dy = result.y - y0;
        let dz = result.z - z0;
        let dist = (dx * dx + dy * dy + dz * dz).sqrt();
        assert!(dist < 5.0, "displacement too large: {dist}m");
    }

    #[test]
    fn test_same_frame_epoch_propagation() {
        let result = transform_epoch_aware(
            3_711_000.0,
            1_023_000.0,
            5_036_000.0,
            "ITRF2020",
            Epoch::new(2010.0).unwrap(),
            "ITRF2020",
            Epoch::new(2020.0).unwrap(),
            PlateMotionModel::Itrf2020,
        )
        .expect("should succeed");
        assert!(result.x.is_finite());
    }

    #[test]
    fn test_unknown_frame_returns_error() {
        // The location (central Asia) resolves to a plate fine, but there is no
        // frame transform connecting "UNKNOWNFRAME" to ITRF2020, so the
        // frame-path lookup fails.
        let result = transform_epoch_aware(
            1_000_000.0,
            2_000_000.0,
            3_000_000.0,
            "UNKNOWNFRAME",
            Epoch::new(2010.0).unwrap(),
            "ITRF2020",
            Epoch::new(2015.0).unwrap(),
            PlateMotionModel::Itrf2020,
        );
        assert!(result.is_err(), "should fail for unknown frame");
    }

    #[test]
    fn test_result_metadata_preserved() {
        let result = transform_epoch_aware(
            4_627_798.0,
            119_795.0,
            4_369_668.0,
            "ITRF2020",
            Epoch::new(2015.0).unwrap(),
            "ITRF2014",
            Epoch::new(2020.0).unwrap(),
            PlateMotionModel::Itrf2014,
        )
        .expect("should succeed");
        assert_eq!(result.from_frame, "ITRF2020");
        assert_eq!(result.to_frame, "ITRF2014");
        assert!((result.from_epoch.year - 2015.0).abs() < 1e-10);
        assert!((result.to_epoch.year - 2020.0).abs() < 1e-10);
    }

    /// Same-epoch ITRF2020→ITRF2014 (no propagation) must equal the pure
    /// frame Helmert verified against `cct +init=ITRF2020:ITRF2014`:
    ///   in : 4627798.0 119795.0 4369668.0 @2015
    ///   out: 4627797.996656 119794.999050 4369667.999565
    /// This point is over open ocean (no plate polygon), proving that a pure
    /// same-epoch frame change no longer needs a plate.
    #[test]
    fn test_same_epoch_frame_change_matches_cct() {
        let result = transform_epoch_aware(
            4_627_798.0,
            119_795.0,
            4_369_668.0,
            "ITRF2020",
            Epoch::new(2015.0).unwrap(),
            "ITRF2014",
            Epoch::new(2015.0).unwrap(),
            PlateMotionModel::Itrf2020,
        )
        .expect("same-epoch frame change should succeed anywhere");
        assert!(
            (result.x - 4_627_797.996_656).abs() < 5e-6,
            "x={}",
            result.x
        );
        assert!((result.y - 119_794.999_050).abs() < 5e-6, "y={}", result.y);
        assert!(
            (result.z - 4_369_667.999_565).abs() < 5e-6,
            "z={}",
            result.z
        );
    }

    /// Multi-hop, same-epoch ITRF2014→ITRF2000 (via the ITRF2020 hub) must equal
    /// the PROJ pipeline result:
    ///   cct +proj=pipeline +step +inv +init=ITRF2020:ITRF2014
    ///                      +step +init=ITRF2020:ITRF2000
    ///   in : 4627798.0 119795.0 4369668.0 @2015 → 4627798.013556 119795.002020 4369667.976067
    #[test]
    fn test_multi_hop_itrf2014_to_itrf2000_matches_cct() {
        let result = transform_epoch_aware(
            4_627_798.0,
            119_795.0,
            4_369_668.0,
            "ITRF2014",
            Epoch::new(2015.0).unwrap(),
            "ITRF2000",
            Epoch::new(2015.0).unwrap(),
            PlateMotionModel::Itrf2020,
        )
        .expect("multi-hop should succeed");
        assert!(
            (result.x - 4_627_798.013_556).abs() < 1e-5,
            "x={}",
            result.x
        );
        assert!((result.y - 119_795.002_020).abs() < 1e-5, "y={}", result.y);
        assert!(
            (result.z - 4_369_667.976_067).abs() < 1e-5,
            "z={}",
            result.z
        );
    }

    /// A North American station must NOT be assigned Eurasian velocity. With a
    /// 10-year same-frame propagation the point moves along the North American
    /// plate; using the (wrong) Eurasian pole would give a materially different
    /// displacement vector. Denver-ish ECEF (~-105°E, ~40°N).
    #[test]
    fn test_north_american_station_uses_north_american_plate() {
        let (x0, y0, z0) = (-1_274_000.0, -4_733_000.0, 4_074_000.0);
        // Confirm the resolved plate is North American (via coordinate lookup).
        assert_eq!(
            resolve_plate("ITRF2014", x0, y0, z0),
            Some(Plate::NorthAmerican)
        );
        let via_pipeline = transform_epoch_aware(
            x0,
            y0,
            z0,
            "ITRF2014",
            Epoch::new(2010.0).unwrap(),
            "ITRF2014",
            Epoch::new(2020.0).unwrap(),
            PlateMotionModel::Itrf2014,
        )
        .expect("should succeed");
        // Reference: direct propagation with the North American plate.
        let expected = PlateMotionModel::Itrf2014
            .propagate(Plate::NorthAmerican, x0, y0, z0, 2010.0, 2020.0)
            .expect("propagate");
        assert!((via_pipeline.x - expected.0).abs() < 1e-9);
        assert!((via_pipeline.y - expected.1).abs() < 1e-9);
        assert!((via_pipeline.z - expected.2).abs() < 1e-9);
        // And the Eurasian pole would differ by clearly more than a millimetre.
        let eurasian = PlateMotionModel::Itrf2014
            .propagate(Plate::Eurasian, x0, y0, z0, 2010.0, 2020.0)
            .expect("propagate");
        let diff = ((expected.0 - eurasian.0).powi(2)
            + (expected.1 - eurasian.1).powi(2)
            + (expected.2 - eurasian.2).powi(2))
        .sqrt();
        assert!(
            diff > 1e-3,
            "N. American vs Eurasian velocity should differ by >1mm over a decade: {diff}m"
        );
    }

    #[test]
    fn test_frames_convenience_entry_point() {
        let from: RefFrameEpoch = "ITRF2020@2015".parse().unwrap();
        let to: RefFrameEpoch = "ITRF2014@2015".parse().unwrap();
        let result = transform_epoch_aware_frames(
            &from,
            &to,
            4_627_798.0,
            119_795.0,
            4_369_668.0,
            PlateMotionModel::Itrf2020,
        )
        .expect("frames entry point should succeed");
        assert!(
            (result.x - 4_627_797.996_656).abs() < 5e-6,
            "x={}",
            result.x
        );
        assert_eq!(result.from_frame, "ITRF2020");
        assert_eq!(result.to_frame, "ITRF2014");
    }

    #[test]
    fn test_frames_convenience_requires_epochs() {
        let from = RefFrameEpoch::static_frame("ITRF2020");
        let to: RefFrameEpoch = "ITRF2014@2015".parse().unwrap();
        let result = transform_epoch_aware_frames(
            &from,
            &to,
            4_627_798.0,
            119_795.0,
            4_369_668.0,
            PlateMotionModel::Itrf2020,
        );
        assert!(result.is_err(), "static source frame should be rejected");
    }
}