axiolid_linear_intersection/
line_line.rs

1//! Line/line classification in the plane.
2//!
3//! The parallel-versus-coincident distinction is decided by certified
4//! predicates, never by comparing a determinant against a global epsilon. Two
5//! lines that are parallel and two that are the same line are different
6//! topological facts, and a rounding artefact must not be able to swap them.
7
8use axiolid_core::{Point2, Scalar, Tolerance, Vec2};
9use axiolid_guarantees::Sign;
10use axiolid_linear::Line2;
11use axiolid_predicates::orient2d;
12
13use crate::error::{InputSide, LinearIntersectionError};
14use crate::validate::{finite_point, finite_vector, nonzero_direction};
15
16/// How two unbounded planar lines relate.
17#[non_exhaustive]
18#[derive(Debug, Clone, Copy, PartialEq)]
19pub enum LineLineIntersection2 {
20    /// The lines cross at exactly one point.
21    Point {
22        /// The intersection point.
23        point: Point2,
24        /// Parameter on the left line.
25        left_parameter: Scalar,
26        /// Parameter on the right line.
27        right_parameter: Scalar,
28    },
29    /// The lines have the same direction and never meet.
30    Parallel,
31    /// The lines are the same line.
32    Coincident,
33}
34
35/// Classify two unbounded planar lines.
36///
37/// `tolerance` governs the residual acceptance of the computed point, not the
38/// topological parallel/coincident decision: that decision uses certified
39/// orientation predicates so it cannot be flipped by a scaled epsilon.
40pub fn line_line2(
41    left: Line2,
42    right: Line2,
43    tolerance: Tolerance,
44) -> Result<LineLineIntersection2, LinearIntersectionError> {
45    validate(left, InputSide::Left)?;
46    validate(right, InputSide::Right)?;
47
48    // Parallelism is a property of the DIRECTIONS alone. Sampling a point from
49    // each line and comparing sides is wrong for unbounded lines: two crossing
50    // lines can easily place both sample points on the same side.
51    let direction_origin = left.origin;
52    let left_tip = translate(direction_origin, left.direction)?;
53    let right_tip = translate(direction_origin, right.direction)?;
54    let directions_parallel =
55        sign(orient2d(direction_origin, left_tip, right_tip).sign())? == Sign::Zero;
56
57    if directions_parallel {
58        // Same direction: either the same line, or a disjoint translate of it.
59        let left_second = translate(left.origin, left.direction)?;
60        let right_on_left = sign(orient2d(left.origin, left_second, right.origin).sign())?;
61        return Ok(if right_on_left == Sign::Zero {
62            LineLineIntersection2::Coincident
63        } else {
64            LineLineIntersection2::Parallel
65        });
66    }
67
68    let determinant = left.direction.x * right.direction.y - left.direction.y * right.direction.x;
69    if !determinant.is_finite() {
70        return Err(LinearIntersectionError::ArithmeticOverflow);
71    }
72    if determinant == 0.0 {
73        // The predicate proved the directions are not collinear, so a zero
74        // float determinant means the parameters are not computable in f64.
75        return Err(LinearIntersectionError::ArithmeticOverflow);
76    }
77
78    let delta_x = right.origin.x - left.origin.x;
79    let delta_y = right.origin.y - left.origin.y;
80    let left_parameter = (delta_x * right.direction.y - delta_y * right.direction.x) / determinant;
81    let right_parameter = (delta_x * left.direction.y - delta_y * left.direction.x) / determinant;
82    if !left_parameter.is_finite() || !right_parameter.is_finite() {
83        return Err(LinearIntersectionError::ArithmeticOverflow);
84    }
85
86    let point = Point2 {
87        x: left.origin.x + left_parameter * left.direction.x,
88        y: left.origin.y + left_parameter * left.direction.y,
89    };
90    if !point.x.is_finite() || !point.y.is_finite() {
91        return Err(LinearIntersectionError::ArithmeticOverflow);
92    }
93
94    // Residual check against the caller's tolerance: the classification is
95    // certified, but the returned coordinate is still a rounded value and the
96    // caller asked for a specific accuracy.
97    let residual_x = point.x - (right.origin.x + right_parameter * right.direction.x);
98    let residual_y = point.y - (right.origin.y + right_parameter * right.direction.y);
99    let residual = residual_x.hypot(residual_y);
100    if !residual.is_finite() {
101        return Err(LinearIntersectionError::ArithmeticOverflow);
102    }
103    let scale = point.x.abs().max(point.y.abs()).max(1.0);
104    if residual > tolerance.linear() * scale && residual > f64::EPSILON * scale * 16.0 {
105        return Err(LinearIntersectionError::ArithmeticOverflow);
106    }
107
108    Ok(LineLineIntersection2::Point {
109        point,
110        left_parameter,
111        right_parameter,
112    })
113}
114
115fn validate(line: Line2, side: InputSide) -> Result<(), LinearIntersectionError> {
116    finite_point(line.origin, side)?;
117    finite_vector(line.direction, side)?;
118    nonzero_direction(line.direction, side)
119}
120
121/// Offset a point by a vector, refusing a non-finite result rather than
122/// letting an infinity reach a predicate.
123fn translate(point: Point2, offset: Vec2) -> Result<Point2, LinearIntersectionError> {
124    let moved = Point2 {
125        x: point.x + offset.x,
126        y: point.y + offset.y,
127    };
128    if moved.x.is_finite() && moved.y.is_finite() {
129        Ok(moved)
130    } else {
131        Err(LinearIntersectionError::ArithmeticOverflow)
132    }
133}
134
135/// A certified predicate must produce a sign; uncertainty is a refusal.
136fn sign(sign: Option<Sign>) -> Result<Sign, LinearIntersectionError> {
137    sign.ok_or(LinearIntersectionError::ArithmeticOverflow)
138}
axiolid_linear_intersection/line_line.rs

axiolid_linear_intersection/
line_line.rs
