egml-core 0.0.2-alpha.8

use crate::error::Error;
use crate::model::base::{AsAbstractGml, Id};
use crate::model::common::{ApplyTransform, ComputeEnvelope, IterGeometries};
use crate::model::geometry::primitives::{AbstractRing, AsAbstractRing, AsAbstractRingMut};
use crate::model::geometry::refs::AbstractGeometryKindRef;
use crate::model::geometry::{DirectPosition, Envelope};
use crate::util::distance::point_segment_distance_2d;
use crate::{impl_abstract_ring_mut_traits, impl_abstract_ring_traits, impl_has_geometry_type};
use nalgebra::{Isometry3, Rotation3, Scale3, Transform3, Vector2, Vector3};

const MINIMUM_NUMBER_OF_POINTS: usize = 3;

/// An implicitly closed ring of at least 3 distinct, non-adjacent positions.
///
/// Corresponds to `gml:LinearRing` in [OGC 07-036 §10.5.8](https://docs.ogc.org/is/07-036/07-036.pdf).  The ring is
/// implicitly closed: the last position is not repeated.
///
/// # Invariants
///
/// - At least 3 positions.
/// - No two adjacent positions are equal.
/// - First and last positions are not equal.
#[derive(Debug, Clone, PartialEq, Default)]
pub struct LinearRing {
    pub abstract_ring: AbstractRing,
    /// Boxed slice rather than `Vec`: points are written once and then only
    /// read or transformed in place, so the capacity field would be dead weight.
    points: Box<[DirectPosition]>,
}

impl LinearRing {
    /// Creates a new `LinearRing` from an ordered list of positions.
    ///
    /// # Errors
    ///
    /// Returns [`Error::TooFewElements`] if `points` has fewer than 3 entries.
    /// Returns [`Error::AdjacentDuplicatePositions`] if adjacent positions are equal.
    /// Returns [`Error::RepeatedClosingVertex`] if the first and last positions are equal.
    pub fn new(points: impl IntoIterator<Item = DirectPosition>) -> Result<Self, Error> {
        let points: Vec<DirectPosition> = points.into_iter().collect();
        Self::validate_points(&points, None)?;

        Ok(Self {
            abstract_ring: AbstractRing::default(),
            points: points.into_boxed_slice(),
        })
    }

    pub fn from_abstract_ring(
        abstract_ring: AbstractRing,
        points: impl IntoIterator<Item = DirectPosition>,
    ) -> Result<Self, Error> {
        let points: Vec<DirectPosition> = points.into_iter().collect();
        Self::validate_points(&points, None)?;
        Ok(Self {
            abstract_ring,
            points: points.into_boxed_slice(),
        })
    }

    fn validate_points(points: &[DirectPosition], id: Option<&Id>) -> Result<(), Error> {
        if let Some((index, window)) = points.windows(2).enumerate().find(|(_, w)| w[0] == w[1]) {
            return Err(Error::AdjacentDuplicatePositions {
                index,
                position: window[0],
            });
        }

        if points.len() < MINIMUM_NUMBER_OF_POINTS {
            let detail = if id.is_none() {
                Some(format!(
                    "points: {}",
                    points
                        .iter()
                        .map(|p| p.to_string())
                        .collect::<Vec<String>>()
                        .join(", ")
                ))
            } else {
                None
            };

            return Err(Error::TooFewElements {
                geometry: "gml:LinearRing",
                minimum: 3,
                spec: Some("OGC 07-036 §10.5.8"),
                id: id.cloned(),
                detail,
            });
        }

        let first = *points.first().expect("non-empty validated above");
        if first == *points.last().expect("non-empty validated above") {
            return Err(Error::RepeatedClosingVertex { position: first });
        }

        Ok(())
    }

    /// Returns the positions of this ring.
    pub fn points(&self) -> &[DirectPosition] {
        &self.points
    }

    /// Replaces the positions of this ring.
    ///
    /// # Errors
    ///
    /// Returns the same errors as [`new`](Self::new).
    pub fn set_points(&mut self, val: Vec<DirectPosition>) -> Result<(), Error> {
        Self::validate_points(&val, self.id())?;

        self.points = val.into_boxed_slice();
        Ok(())
    }
}

impl AsAbstractRing for LinearRing {
    fn abstract_ring(&self) -> &AbstractRing {
        &self.abstract_ring
    }
}

impl AsAbstractRingMut for LinearRing {
    fn abstract_ring_mut(&mut self) -> &mut AbstractRing {
        &mut self.abstract_ring
    }
}

impl_abstract_ring_traits!(LinearRing);
impl_abstract_ring_mut_traits!(LinearRing);
impl_has_geometry_type!(LinearRing, LinearRing);

impl LinearRing {
    /// Returns the 3D perimeter of this ring.
    ///
    /// Sums the distances between consecutive vertices, including the implicit
    /// closing segment from the last vertex back to the first.
    pub fn length_3d(&self) -> f64 {
        let n = self.points.len();
        (0..n)
            .map(|i| {
                let vi: Vector3<f64> = self.points[i].into();
                let vj: Vector3<f64> = self.points[(i + 1) % n].into();
                (vj - vi).norm()
            })
            .sum()
    }

    /// Returns the 3D area_3d of this ring using the cross-product summation formula.
    ///
    /// Computes `|Σ (vᵢ × vᵢ₊₁)| / 2` over all consecutive vertex pairs (with wrap-around),
    /// which gives the correct planar area_3d regardless of orientation in 3D space.
    pub fn area_3d(&self) -> f64 {
        let n = self.points.len();
        let mut cross_sum = Vector3::zeros();
        for i in 0..n {
            let vi: Vector3<f64> = self.points[i].into();
            let vj: Vector3<f64> = self.points[(i + 1) % n].into();
            cross_sum += vi.cross(&vj);
        }
        cross_sum.norm() * 0.5
    }

    /// Returns the minimum Euclidean distance between `point` and the boundary
    /// of this ring, measured in the XY plane (z coordinates are ignored).
    ///
    /// Includes the implicit closing segment from the last vertex back to the
    /// first. The ring is treated as a curve, not an area: a point inside the
    /// ring has a positive distance to its nearest edge.
    pub fn distance_2d(&self, point: &DirectPosition) -> f64 {
        let p = Vector2::new(point.x(), point.y());
        let n = self.points.len();
        (0..n)
            .map(|i| {
                let a = &self.points[i];
                let b = &self.points[(i + 1) % n];
                point_segment_distance_2d(p, Vector2::new(a.x(), a.y()), Vector2::new(b.x(), b.y()))
            })
            .fold(f64::INFINITY, f64::min)
    }
}

impl ApplyTransform for LinearRing {
    fn apply_transform(&mut self, transform: Transform3<f64>) {
        self.points.iter_mut().for_each(|p| {
            p.apply_transform(transform);
        });
    }

    fn apply_isometry(&mut self, isometry: Isometry3<f64>) {
        self.points.iter_mut().for_each(|p| {
            p.apply_isometry(isometry);
        });
    }

    fn apply_translation(&mut self, vector: Vector3<f64>) {
        self.points.iter_mut().for_each(|p| {
            p.apply_translation(vector);
        });
    }

    fn apply_rotation(&mut self, rotation: Rotation3<f64>) {
        self.points.iter_mut().for_each(|p| {
            p.apply_rotation(rotation);
        });
    }

    fn apply_scale(&mut self, scale: Scale3<f64>) {
        self.points.iter_mut().for_each(|p| {
            p.apply_scale(scale);
        });
    }
}

impl ComputeEnvelope for LinearRing {
    /// Returns the axis-aligned bounding box of all positions in this ring.
    fn compute_envelope(&self) -> Option<Envelope> {
        Some(Envelope::from_points(&self.points).expect("linear ring must have valid points"))
    }
}

impl IterGeometries for LinearRing {
    fn iter_geometries(&self) -> Box<dyn Iterator<Item = AbstractGeometryKindRef<'_>> + '_> {
        Box::new(std::iter::once(self.into()))
    }
}

#[cfg(test)]
mod tests {
    use super::*;
    use nalgebra::{Isometry3, Vector3};

    #[test]
    fn length_3d_unit_square_includes_closing_segment() {
        let ring = LinearRing::new([
            DirectPosition::new(0.0, 0.0, 0.0).unwrap(),
            DirectPosition::new(1.0, 0.0, 0.0).unwrap(),
            DirectPosition::new(1.0, 1.0, 0.0).unwrap(),
            DirectPosition::new(0.0, 1.0, 0.0).unwrap(),
        ])
        .unwrap();
        assert!((ring.length_3d() - 4.0).abs() < 1e-10);
    }

    #[test]
    fn length_3d_right_triangle_in_3d() {
        // 3-4-5 right triangle in the XZ plane — perimeter 12.
        let ring = LinearRing::new([
            DirectPosition::new(0.0, 0.0, 0.0).unwrap(),
            DirectPosition::new(3.0, 0.0, 0.0).unwrap(),
            DirectPosition::new(3.0, 0.0, 4.0).unwrap(),
        ])
        .unwrap();
        assert!((ring.length_3d() - 12.0).abs() < 1e-10);
    }

    #[test]
    fn area_3d_unit_square_xy() {
        let ring = LinearRing::new([
            DirectPosition::new(0.0, 0.0, 0.0).unwrap(),
            DirectPosition::new(1.0, 0.0, 0.0).unwrap(),
            DirectPosition::new(1.0, 1.0, 0.0).unwrap(),
            DirectPosition::new(0.0, 1.0, 0.0).unwrap(),
        ])
        .unwrap();
        assert!((ring.area_3d() - 1.0).abs() < 1e-10);
    }

    #[test]
    fn area_3d_tilted_rectangle() {
        // Rectangle with sides 1 and sqrt(5) tilted in 3D — area_3d should be sqrt(5).
        let ring = LinearRing::new([
            DirectPosition::new(0.0, 0.0, 0.0).unwrap(),
            DirectPosition::new(1.0, 0.0, 0.0).unwrap(),
            DirectPosition::new(1.0, 1.0, 2.0).unwrap(),
            DirectPosition::new(0.0, 1.0, 2.0).unwrap(),
        ])
        .unwrap();
        let expected = 5.0_f64.sqrt();
        assert!((ring.area_3d() - expected).abs() < 1e-10);
    }

    #[test]
    fn area_3d_triangle() {
        // Right triangle with legs 3 and 4 — area_3d should be 6.
        let ring = LinearRing::new([
            DirectPosition::new(0.0, 0.0, 0.0).unwrap(),
            DirectPosition::new(3.0, 0.0, 0.0).unwrap(),
            DirectPosition::new(0.0, 4.0, 0.0).unwrap(),
        ])
        .unwrap();
        assert!((ring.area_3d() - 6.0).abs() < 1e-10);
    }

    #[test]
    fn linear_ring_construction_test() {
        let points = vec![
            DirectPosition::new(
                601.927_914_447_452_5,
                1130.4631113024607,
                9.013_090_391_538_235,
            )
            .unwrap(),
            DirectPosition::new(
                601.927_918_328_473_4,
                1130.4631032795705,
                9.013_090_723_310_274,
            )
            .unwrap(),
            DirectPosition::new(
                601.927_918_328_473_4,
                1130.4631032795705,
                9.013_090_723_310_274,
            )
            .unwrap(),
        ];
        let result = LinearRing::new(points);

        assert!(matches!(
            result,
            Err(Error::AdjacentDuplicatePositions { .. })
        ));
    }

    #[test]
    fn offset_linear_ring_test() {
        let mut linear_ring = LinearRing::new([
            DirectPosition::new(1.0, 2.0, 3.0).unwrap(),
            DirectPosition::new(2.0, 4.0, 6.0).unwrap(),
            DirectPosition::new(4.0, 7.0, 9.0).unwrap(),
        ])
        .unwrap();
        //let offset = nalgebra::Vector3::<f64>::new(1.0, -1.0, 3.0);
        let isometry: Isometry3<f64> =
            Isometry3::new(Vector3::new(1.0, -1.0, 3.0), Default::default());
        let expected_linear_ring = LinearRing::new([
            DirectPosition::new(2.0, 1.0, 6.0).unwrap(),
            DirectPosition::new(3.0, 3.0, 9.0).unwrap(),
            DirectPosition::new(5.0, 6.0, 12.0).unwrap(),
        ])
        .unwrap();

        linear_ring.apply_isometry(isometry);

        assert_eq!(linear_ring, expected_linear_ring);
    }

    fn pos(x: f64, y: f64, z: f64) -> DirectPosition {
        DirectPosition::new(x, y, z).unwrap()
    }

    fn square_ring() -> LinearRing {
        LinearRing::new([
            pos(0.0, 0.0, 0.0),
            pos(10.0, 0.0, 0.0),
            pos(10.0, 10.0, 0.0),
            pos(0.0, 10.0, 0.0),
        ])
        .unwrap()
    }

    #[test]
    fn distance_2d_outside() {
        assert!((square_ring().distance_2d(&pos(5.0, -3.0, 0.0)) - 3.0).abs() < 1e-10);
    }

    #[test]
    fn distance_2d_uses_implicit_closing_segment() {
        // Nearest edge is (0,10)→(0,0), which is only present via wrap-around.
        assert!((square_ring().distance_2d(&pos(-2.0, 5.0, 0.0)) - 2.0).abs() < 1e-10);
    }

    #[test]
    fn distance_2d_inside_measures_to_boundary() {
        // The ring is a curve, not an area: interior points have positive distance.
        assert!((square_ring().distance_2d(&pos(5.0, 4.0, 0.0)) - 4.0).abs() < 1e-10);
    }

    #[test]
    fn distance_2d_clamps_to_corner() {
        assert!((square_ring().distance_2d(&pos(13.0, 14.0, 0.0)) - 5.0).abs() < 1e-10);
    }

    #[test]
    fn distance_2d_on_boundary_is_zero() {
        assert!(square_ring().distance_2d(&pos(10.0, 7.0, 0.0)).abs() < 1e-10);
    }

    #[test]
    fn distance_2d_ignores_z_and_vertical_wall_ring() {
        // A vertical wall ring collapses to a line segment from (0,0) to (10,0) in 2D.
        let wall = LinearRing::new([
            pos(0.0, 0.0, 0.0),
            pos(10.0, 0.0, 0.0),
            pos(10.0, 0.0, 5.0),
            pos(0.0, 0.0, 5.0),
        ])
        .unwrap();
        let d = wall.distance_2d(&pos(5.0, 3.0, 100.0));
        assert!(!d.is_nan());
        assert!((d - 3.0).abs() < 1e-10);
    }
}