iron_shapes/
edge.rs

1// Copyright (c) 2018-2020 Thomas Kramer.
2// SPDX-FileCopyrightText: 2018-2022 Thomas Kramer
3//
4// SPDX-License-Identifier: AGPL-3.0-or-later
5
6//! An edge is a line segment from a start-point to a end-point.
7
8use crate::point::Point;
9use crate::rect::Rect;
10use crate::vector::Vector;
11
12use crate::CoordinateType;
13
14use num_traits::cast::NumCast;
15use num_traits::Float;
16
17use crate::traits::MapPointwise;
18pub use crate::traits::{BoundingBox, RotateOrtho, TryBoundingBox, TryCastCoord};
19pub use crate::types::Angle;
20
21pub use crate::types::{ContainsResult, Side};
22
23/// Get the endpoints of an edge.
24pub trait EdgeEndpoints<T> {
25    /// Get the start point of the edge.
26    fn start(&self) -> Point<T>;
27    /// Get the end point of the edge.
28    fn end(&self) -> Point<T>;
29}
30
31/// Define the intersection between two edges (i.e. line segments).
32pub trait EdgeIntersect
33where
34    Self: Sized,
35{
36    /// Numeric type used for expressing the end-point coordinates of the edge.
37    type Coord;
38    /// Numeric type used for expressing an intersection-point of two edges.
39    /// Often this might be the same as `Coord`.
40    type IntersectionCoord;
41
42    /// Compute intersection of two edges.
43    fn edge_intersection(
44        &self,
45        other: &Self,
46    ) -> EdgeIntersection<Self::Coord, Self::IntersectionCoord, Self>;
47}
48
49/// Iterate over edges.
50/// For an n-gon this would produce n edges.
51pub trait IntoEdges<T> {
52    /// Type of edge which will be returned.
53    type Edge: EdgeEndpoints<T>;
54    /// Iterator type.
55    type EdgeIter: Iterator<Item = Self::Edge>;
56
57    /// Get an iterator over edges.
58    fn into_edges(self) -> Self::EdgeIter;
59}
60
61/// Return type for the edge-edge intersection functions.
62/// Stores all possible results of an edge to edge intersection.
63///
64/// # Note on coordinate types:
65/// There are two coordinate types (which may be the same concrete type):
66/// * `OriginalCoord` is the coordinate type used to define the edge end-points. An intersection at the end-points
67/// can be expressed with this coordinate type.
68/// * `IntersectionCoord` is the coordinate type used to express intersection points somewhere in the middle of the edge.
69/// This may differ from the coordinate type of the end-points. For example if the end-points are stored in integer coordinates
70/// the intersection may require rational coordinates. But in special cases such as axis-aligned edges, the intersection point
71/// can indeed be expressed in integer coordinates.
72#[derive(Clone, Copy, PartialEq, Eq, Debug)]
73pub enum EdgeIntersection<IntersectionCoord, OriginalCoord, Edge> {
74    /// No intersection.
75    None,
76    /// Intersection in a single point but not on an endpoint of an edge.
77    Point(Point<IntersectionCoord>),
78    /// Intersection in an endpoint of an edge.
79    EndPoint(Point<OriginalCoord>),
80    /// Full or partial overlap.
81    Overlap(Edge),
82}
83
84/// Return type for the line-line intersection functions.
85/// Stores all possible results of a line to line intersection.
86#[derive(Clone, Copy, PartialEq, Debug)]
87pub enum LineIntersection<TP: CoordinateType, TE: CoordinateType> {
88    /// No intersection at all.
89    None,
90    /// Intersection in a single point.
91    /// Besides the intersection point also an other expression for the intersection point is given.
92    /// The three values `(a, b, c)` describe the intersection point in terms of a starting point (the starting point
93    /// of the edge which defines the line) and the direction of the edge multiplied by a fraction.
94    ///
95    /// `edge.start + edge.vector()*a/c == p` and
96    /// `other_edge.start + other_edge.vector()*b/c == p`.
97    Point(Point<TP>, (TE, TE, TE)),
98    /// Lines are collinear.
99    Collinear,
100}
101
102impl<T: Copy> From<Edge<T>> for (Point<T>, Point<T>) {
103    fn from(e: Edge<T>) -> Self {
104        (e.start, e.end)
105    }
106}
107
108impl<T: Copy> From<&Edge<T>> for (Point<T>, Point<T>) {
109    fn from(e: &Edge<T>) -> Self {
110        (e.start, e.end)
111    }
112}
113
114impl<T: CoordinateType> From<(Point<T>, Point<T>)> for Edge<T> {
115    fn from(points: (Point<T>, Point<T>)) -> Self {
116        Edge::new(points.0, points.1)
117    }
118}
119
120impl<T: CoordinateType> From<[Point<T>; 2]> for Edge<T> {
121    fn from(points: [Point<T>; 2]) -> Self {
122        Edge::new(points[0], points[1])
123    }
124}
125
126/// An edge (line segment) is represented by its starting point and end point.
127#[derive(Clone, Copy, PartialEq, Eq, Hash, Debug)]
128#[cfg_attr(feature = "serde", derive(Serialize, Deserialize))]
129pub struct Edge<T> {
130    /// Start-point of the edge.
131    pub start: Point<T>,
132    /// End-point of the edge.
133    pub end: Point<T>,
134}
135
136impl<T: Copy> EdgeEndpoints<T> for Edge<T> {
137    fn start(&self) -> Point<T> {
138        self.start
139    }
140
141    fn end(&self) -> Point<T> {
142        self.end
143    }
144}
145
146impl<T: Copy> Edge<T> {
147    /// Create a new `Edge` from two arguments that implement `Into<Point>`.
148    pub fn new<C>(start: C, end: C) -> Self
149    where
150        C: Into<Point<T>>,
151    {
152        Edge {
153            start: start.into(),
154            end: end.into(),
155        }
156    }
157
158    /// Return the same edge but with the two points swapped.
159    pub fn reversed(&self) -> Self {
160        Edge {
161            start: self.end,
162            end: self.start,
163        }
164    }
165}
166
167impl<T: PartialEq> Edge<T> {
168    /// Check if edge is degenerate.
169    /// An edge is degenerate if start point and end point are equal.
170    pub fn is_degenerate(&self) -> bool {
171        self.start == self.end
172    }
173
174    /// Test if this edge is either horizontal or vertical.
175    pub fn is_rectilinear(&self) -> bool {
176        !self.is_degenerate() && (self.start.x == self.end.x || self.start.y == self.end.y)
177    }
178
179    /// Test if this edge is horizontal.
180    pub fn is_horizontal(&self) -> bool {
181        !self.is_degenerate() && self.start.y == self.end.y
182    }
183
184    /// Test if this edge is vertical.
185    pub fn is_vertical(&self) -> bool {
186        !self.is_degenerate() && self.start.x == self.end.x
187    }
188}
189
190impl<T: CoordinateType> Edge<T> {
191    /// Returns the vector from `self.start` to `self.end`.
192    pub fn vector(&self) -> Vector<T> {
193        self.end - self.start
194    }
195
196    /// Tells on which side of the edge a point is.
197    ///
198    /// # Panics
199    /// Panics if the edge is degenerate.
200    ///
201    /// Returns `Side::Left` if the point is on the left side,
202    /// `Side::Right` if the point is on the right side
203    /// or `Side::Center` if the point lies exactly on the line.
204    pub fn side_of(&self, point: Point<T>) -> Side {
205        assert!(!self.is_degenerate());
206
207        let a = self.vector();
208        let b = point - self.start;
209
210        // Handle rectilinear cases.
211        if a.x.is_zero() {
212            let side = if point.x < self.start.x {
213                Side::Left
214            } else if point.x > self.start.x {
215                Side::Right
216            } else {
217                Side::Center
218            };
219
220            // Maybe flip the side depending on the orientation of the edge.
221            let side = if self.start.y < self.end.y {
222                side
223            } else {
224                side.other()
225            };
226
227            return side;
228        } else if a.y.is_zero() {
229            let side = if point.y < self.start.y {
230                Side::Right
231            } else if point.y > self.start.y {
232                Side::Left
233            } else {
234                Side::Center
235            };
236
237            // Maybe flip the side depending on the orientation of the edge.
238            let side = if self.start.x < self.end.x {
239                side
240            } else {
241                side.other()
242            };
243
244            return side;
245        }
246
247        let area = b.cross_prod(a);
248
249        if area.is_zero() {
250            Side::Center
251        } else if area < T::zero() {
252            Side::Left
253        } else {
254            Side::Right
255        }
256    }
257
258    //    /// Find minimal distance between two edges.
259    //    pub fn distance_to_edge(self, other: Edge<T>) -> f64 {
260    //        let d1 = self.distance(other.start);
261    //        let d2 = self.distance(other.end);
262    //        let d3 = other.distance(self.start);
263    //        let d4 = other.distance(self.end);
264    //
265    //        let min1 = d1.min(d2);
266    //        let min2 = d3.min(d4);
267    //
268    //        // TODO: if intersects
269    //
270    //        min1.min(min2)
271    //    }
272
273    /// Test if point lies on the edge.
274    /// Includes start and end points of edge.
275    pub fn contains_point(&self, point: Point<T>) -> ContainsResult {
276        if self.start == point || self.end == point {
277            ContainsResult::OnBounds
278        } else if self.is_degenerate() {
279            ContainsResult::No
280        } else {
281            let a = self.start - point;
282            let b = self.end - point;
283            // Check if the triangle point-start-end has zero area and that a and b have opposite directions.
284            // If a and b have opposite directions then point is between start and end.
285            if a.cross_prod(b).is_zero() && a.dot(b) <= T::zero() {
286                ContainsResult::WithinBounds
287            } else {
288                ContainsResult::No
289            }
290        }
291    }
292
293    /// Test if point lies on the line defined by the edge.
294    pub fn line_contains_point(&self, point: Point<T>) -> bool {
295        if self.is_degenerate() {
296            self.start == point
297        } else {
298            let l = self.vector();
299            let b = point - self.start;
300
301            l.cross_prod(b).is_zero()
302        }
303    }
304
305    /// Test if two edges are parallel.
306    pub fn is_parallel(&self, other: &Edge<T>) -> bool {
307        if self.is_degenerate() || other.is_degenerate() {
308            false
309        } else {
310            let a = self.vector();
311            let b = other.vector();
312
313            if a.x.is_zero() {
314                // Both vertical?
315                b.x.is_zero()
316            } else if a.y.is_zero() {
317                // Both horizontal?
318                b.y.is_zero()
319            } else if b.x.is_zero() || b.y.is_zero() {
320                false
321            } else {
322                a.cross_prod(b).is_zero()
323            }
324        }
325    }
326
327    /// Test if two edges are collinear, i.e. are on the same line.
328    pub fn is_collinear(&self, other: &Edge<T>) -> bool
329    where
330        T: CoordinateType,
331    {
332        if self.is_degenerate() || other.is_degenerate() {
333            false
334        } else {
335            let v = self.vector();
336            let a = other.start - self.start;
337            let b = other.end - self.start;
338
339            if self.start.x == self.end.x {
340                // Self is vertical.
341                // The edges are collinear iff `other` has the same x coordinates as `self`.
342                other.start.x == self.start.x && other.end.x == self.start.x
343            } else if self.start.y == self.end.y {
344                // `self` is horizontal
345                other.start.y == self.start.y && other.end.y == self.start.y
346            } else if other.start.x == other.end.x {
347                self.start.x == other.start.x && self.end.x == other.start.x
348            } else if other.start.y == other.end.y {
349                // `other` is horizontal
350                self.start.y == other.start.y && self.end.y == other.start.y
351            } else {
352                v.cross_prod(a).is_zero() && v.cross_prod(b).is_zero()
353            }
354        }
355    }
356
357    /// Test edges for coincidence.
358    /// Two edges are coincident if they are oriented the same way
359    /// and share more than one point (implies that they must be parallel).
360    pub fn is_coincident(&self, other: &Edge<T>) -> bool {
361        // TODO: use is_collinear
362        // TODO: approximate version for floating point types
363
364        let self_start_in_other = other.contains_point(self.start).inclusive_bounds();
365        let self_end_in_other = other.contains_point(self.end).inclusive_bounds();
366        let other_start_in_self = self.contains_point(other.start).inclusive_bounds();
367        let other_end_in_self = self.contains_point(other.end).inclusive_bounds();
368
369        let share_more_than_one_point = self.end != other.start
370            && self.start != other.end
371            && (self_start_in_other || other_start_in_self)
372            && (other_end_in_self || self_end_in_other);
373
374        // Sharing more than one point should imply that the edges are parallel.
375        debug_assert!(if share_more_than_one_point {
376            self.is_parallel(other)
377        } else {
378            true
379        });
380
381        let oriented_the_same_way = self.vector().dot(other.vector()) > T::zero();
382
383        share_more_than_one_point && oriented_the_same_way
384    }
385
386    /// Test if two edges are approximately parallel.
387    /// To be used for float coordinates.
388    /// Inspired by algorithm on page 241 of "Geometric Tools for Computer Graphics".
389    pub fn is_parallel_approx(&self, other: &Edge<T>, epsilon_squared: T) -> bool {
390        if self.is_degenerate() || other.is_degenerate() {
391            false
392        } else {
393            let d1 = self.vector();
394            let d2 = other.vector();
395
396            let len1_sqr = d1.norm2_squared();
397            let len2_sqr = d2.norm2_squared();
398
399            let cross = d1.cross_prod(d2);
400            let cross_sqr = cross * cross;
401
402            // TODO: require square tolerance form caller?
403            cross_sqr <= len1_sqr * len2_sqr * epsilon_squared
404        }
405    }
406
407    /// Test if two edges are approximately collinear, i.e. are on the same line.
408    /// Inspired by algorithm on page 241 of "Geometric Tools for Computer Graphics".
409    pub fn is_collinear_approx(&self, other: &Edge<T>, epsilon_squared: T) -> bool {
410        if self.is_degenerate() || other.is_degenerate() {
411            false
412        } else {
413            let d1 = self.vector();
414            let d2 = other.vector();
415
416            let len1_sqr = d1.norm2_squared();
417            let len2_sqr = d2.norm2_squared();
418
419            let cross = d1.cross_prod(d2);
420            let cross_sqr = cross * cross;
421
422            let approx_parallel = cross_sqr <= len1_sqr * len2_sqr * epsilon_squared;
423
424            if approx_parallel {
425                let e = other.start - self.start;
426                let len_e_sqrt = e.norm2_squared();
427                let cross = e.cross_prod(d1);
428                let cross_sqr = cross * cross;
429                cross_sqr <= len1_sqr * len_e_sqrt * epsilon_squared
430            } else {
431                false
432            }
433        }
434    }
435
436    /// Test if lines defined by the edges intersect.
437    /// If the lines are collinear they are also considered intersecting.
438    pub fn lines_intersect_approx(&self, other: &Edge<T>, epsilon_squared: T) -> bool {
439        !self.is_parallel_approx(other, epsilon_squared)
440            || self.is_collinear_approx(other, epsilon_squared)
441    }
442
443    /// Test if this edge is crossed by the line defined by the other edge.
444    ///
445    /// Returns `WithinBounds` if start and end point of this edge lie on different sides
446    /// of the line defined by the `other` edge or `OnBounds` if at least one of the points
447    /// lies on the line.
448    pub fn crossed_by_line(&self, other: &Edge<T>) -> ContainsResult {
449        // TODO: Handle degenerate cases.
450        let side1 = other.side_of(self.start);
451
452        if side1 == Side::Center {
453            ContainsResult::OnBounds
454        } else {
455            let side2 = other.side_of(self.end);
456
457            if side2 == Side::Center {
458                ContainsResult::OnBounds
459            } else if side1 == side2 {
460                ContainsResult::No
461            } else {
462                ContainsResult::WithinBounds
463            }
464        }
465    }
466
467    /// Test if lines defined by the edges intersect.
468    /// If the lines are collinear they are also considered intersecting.
469    pub fn lines_intersect(&self, other: &Edge<T>) -> bool {
470        !self.is_parallel(other) || self.is_collinear(other)
471    }
472
473    /// Test if two edges intersect.
474    /// If the edges coincide, they also intersect.
475    pub fn edges_intersect(&self, other: &Edge<T>) -> ContainsResult {
476        // Two edges intersect if the start and end point of one edge
477        // lie on opposite sides of the other edge.
478
479        if self.is_degenerate() {
480            other.contains_point(self.start)
481        } else if other.is_degenerate() {
482            self.contains_point(other.start)
483        } else if !self.bounding_box().touches(&other.bounding_box()) {
484            ContainsResult::No
485        } else {
486            // TODO:
487            //        else if self.is_ortho() && other.is_ortho() {
488            //            // We know now that the bounding boxes touch each other.
489            //            // For rectilinear edges this implies that they touch somewhere or overlap.
490            //            true
491            //        } else {
492            match self.crossed_by_line(other) {
493                ContainsResult::No => ContainsResult::No,
494                r => r.min(other.crossed_by_line(self)),
495            }
496        }
497    }
498}
499
500impl<T: CoordinateType + NumCast> Edge<T> {
501    /// Test if point lies on the line defined by the edge.
502    pub fn line_contains_point_approx<F: Float + NumCast>(
503        &self,
504        point: Point<T>,
505        tolerance: F,
506    ) -> bool {
507        if self.is_degenerate() {
508            self.start == point
509        } else {
510            self.distance_to_line_abs_approx::<F>(point) <= tolerance
511        }
512    }
513
514    /// Compute the intersection point of the lines defined by the two edges.
515    ///
516    /// Degenerate lines don't intersect by definition.
517    ///
518    /// Returns `LineIntersection::None` iff the two lines don't intersect.
519    /// Returns `LineIntersection::Collinear` iff both lines are equal.
520    /// Returns `LineIntersection::Point(p,(a,b,c))` iff the lines intersect in exactly one point `p`.
521    /// `f` is a value such that `self.start + self.vector()*a/c == p` and
522    /// `other.start + other.vector()*b/c == p`.
523    ///
524    /// # Examples
525    ///
526    /// ```
527    /// use iron_shapes::point::Point;
528    /// use iron_shapes::edge::*;
529    ///
530    /// let e1 = Edge::new((0, 0), (2, 2));
531    /// let e2 = Edge::new((0, 2), (2, 0));
532    ///
533    /// assert_eq!(e1.line_intersection_approx(&e2, 1e-6),
534    ///     LineIntersection::Point(Point::new(1., 1.), (4, 4, 8)));
535    ///
536    /// assert_eq!(Point::zero() + e1.vector().cast() * 0.5, Point::new(1., 1.));
537    /// ```
538    ///
539    pub fn line_intersection_approx<F: Float>(
540        &self,
541        other: &Edge<T>,
542        tolerance: F,
543    ) -> LineIntersection<F, T> {
544        // TODO: implement algorithm on page 241 of Geometric Tools for Computer Graphics.
545
546        debug_assert!(tolerance >= F::zero(), "Tolerance cannot be negative.");
547
548        if self.is_degenerate() || other.is_degenerate() {
549            LineIntersection::None
550        } else {
551            // TODO: faster implementation if both lines are orthogonal
552
553            let ab = self.vector();
554            let cd = other.vector();
555
556            // Assert that the vectors have a non-zero length. This should already be the case
557            // because the degenerate cases are handled before.
558            debug_assert!(ab.norm2_squared() > T::zero());
559            debug_assert!(cd.norm2_squared() > T::zero());
560
561            let s = ab.cross_prod(cd);
562            let s_float = F::from(s).unwrap();
563
564            // TODO: What if approximate zero due to rounding error?
565            if s_float.abs() <= tolerance * ab.length() * cd.length() {
566                // Lines are parallel
567                // TODO use assertion
568                //                debug_assert!(self.is_parallel_approx(other, tolerance));
569
570                // TODO: check more efficiently for collinear lines.
571                if self.line_contains_point_approx(other.start, tolerance) {
572                    // If the line defined by `self` contains at least one point of `other` then they are equal.
573                    // TODO use assertion
574                    //                    debug_assert!(self.is_collinear_approx(other, tolerance));
575                    LineIntersection::Collinear
576                } else {
577                    LineIntersection::None
578                }
579            } else {
580                let ac = other.start - self.start;
581                let ac_cross_cd = ac.cross_prod(cd);
582                let i = F::from(ac_cross_cd).unwrap() / s_float;
583
584                let p: Point<F> = self.start.cast() + ab.cast() * i;
585
586                let ca_cross_ab = ac.cross_prod(ab);
587
588                // Check that the intersection point lies on the lines indeed.
589                debug_assert!(self
590                    .cast()
591                    .line_contains_point_approx(p, tolerance + tolerance));
592                debug_assert!(other
593                    .cast()
594                    .line_contains_point_approx(p, tolerance + tolerance));
595
596                debug_assert!({
597                    let j = F::from(ca_cross_ab).unwrap() / s_float;
598                    let p2: Point<F> = other.start.cast() + cd.cast() * j;
599                    (p - p2).norm2_squared() < tolerance + tolerance
600                });
601
602                let positions = if s < T::zero() {
603                    (
604                        T::zero() - ac_cross_cd,
605                        T::zero() - ca_cross_ab,
606                        T::zero() - s,
607                    )
608                } else {
609                    (ac_cross_cd, ca_cross_ab, s)
610                };
611
612                LineIntersection::Point(p, positions)
613            }
614        }
615    }
616
617    /// Compute the intersection with another edge.
618    pub fn edge_intersection_approx<F: Float>(
619        &self,
620        other: &Edge<T>,
621        tolerance: F,
622    ) -> EdgeIntersection<F, T, Edge<T>> {
623        debug_assert!(tolerance >= F::zero(), "Tolerance cannot be negative.");
624
625        // Swap direction of other edge such that both have the same direction.
626        let other = if (self.start < self.end) != (other.start < other.end) {
627            other.reversed()
628        } else {
629            *other
630        };
631
632        // Check endpoints for coincidence.
633        // This must be handled separately because equality of the intersection point and endpoints
634        // will not necessarily be detected due to rounding errors.
635        let same_start_start = self.start == other.start;
636        let same_start_end = self.start == other.end;
637
638        let same_end_start = self.end == other.start;
639        let same_end_end = self.end == other.end;
640
641        // Are the edges equal but not degenerate?
642        let fully_coincident =
643            (same_start_start & same_end_end) ^ (same_start_end & same_end_start);
644
645        let result = if self.is_degenerate() {
646            // First degenerate case
647            if other.contains_point_approx(self.start, tolerance) {
648                EdgeIntersection::EndPoint(self.start)
649            } else {
650                EdgeIntersection::None
651            }
652        } else if other.is_degenerate() {
653            // Second degenerate case
654            if self.contains_point_approx(other.start, tolerance) {
655                EdgeIntersection::EndPoint(other.start)
656            } else {
657                EdgeIntersection::None
658            }
659        } else if fully_coincident {
660            EdgeIntersection::Overlap(*self)
661        } else if !self.bounding_box().touches(&other.bounding_box()) {
662            // TODO: Add tolerance here.
663            // If bounding boxes do not touch, then intersection is impossible.
664            EdgeIntersection::None
665        } else {
666            // Compute the intersection of the lines defined by the two edges.
667            let line_intersection = self.line_intersection_approx(&other, tolerance);
668
669            // Then check if the intersection point is on both edges
670            // or find the intersection if the edges overlap.
671
672            match line_intersection {
673                LineIntersection::None => EdgeIntersection::None,
674
675                // Coincident at an endpoint:
676                LineIntersection::Point(_, _) if same_start_start || same_start_end => {
677                    EdgeIntersection::EndPoint(self.start)
678                }
679
680                // Coincident at an endpoint:
681                LineIntersection::Point(_, _) if same_end_start || same_end_end => {
682                    EdgeIntersection::EndPoint(self.end)
683                }
684
685                // Intersection in one point:
686                LineIntersection::Point(p, (pos1, pos2, len)) => {
687                    let len_f = F::from(len).unwrap();
688                    let zero_tol = -tolerance;
689                    let len_tol = len_f + tolerance;
690
691                    let pos1 = F::from(pos1).unwrap();
692                    let pos2 = F::from(pos2).unwrap();
693
694                    // Check if the intersection point `p` lies on the edge.
695                    if pos1 >= zero_tol && pos1 <= len_tol && pos2 >= zero_tol && pos2 <= len_tol {
696                        // Intersection
697
698                        // TODO: rework
699
700                        let zero_tol = tolerance;
701                        let len_tol = len_f - tolerance;
702
703                        if pos1 <= zero_tol {
704                            EdgeIntersection::EndPoint(self.start)
705                        } else if pos1 >= len_tol {
706                            EdgeIntersection::EndPoint(self.end)
707                        } else if pos2 <= zero_tol {
708                            EdgeIntersection::EndPoint(other.start)
709                        } else if pos2 >= len_tol {
710                            EdgeIntersection::EndPoint(other.end)
711                        } else {
712                            // Intersection in-between the endpoints.
713                            debug_assert!(
714                                {
715                                    let e1 = self.cast();
716                                    let e2 = other.cast();
717                                    p != e1.start && p != e1.end && p != e2.start && p != e2.end
718                                },
719                                "Intersection should not be an endpoint."
720                            );
721                            EdgeIntersection::Point(p)
722                        }
723                    } else {
724                        // No intersection.
725                        EdgeIntersection::None
726                    }
727                }
728                LineIntersection::Collinear => {
729                    // TODO use assertion
730                    //                    debug_assert!(self.is_collinear_approx(other, tolerance));
731
732                    // Project all points of the two edges on the line defined by the first edge
733                    // (scaled by the length of the first edge).
734                    // This allows to calculate the interval of overlap in one dimension.
735
736                    let (pa, pb) = self.into();
737                    let (pc, pd) = other.into();
738
739                    let b = pb - pa;
740                    let c = pc - pa;
741                    let d = pd - pa;
742
743                    let dist_a = T::zero();
744                    let dist_b = b.dot(b);
745
746                    let dist_c = b.dot(c);
747                    let dist_d = b.dot(d);
748
749                    let start1 = (dist_a, pa);
750                    let end1 = (dist_b, pb);
751
752                    // Sort end points of other edge.
753                    let (start2, end2) = if dist_c < dist_d {
754                        ((dist_c, pc), (dist_d, pd))
755                    } else {
756                        ((dist_d, pd), (dist_c, pc))
757                    };
758
759                    // Find maximum by distance.
760                    let start = if start1.0 < start2.0 { start2 } else { start1 };
761
762                    // Find minimum by distance.
763                    let end = if end1.0 < end2.0 { end1 } else { end2 };
764
765                    // Check if the edges overlap in more than one point, in exactly one point or
766                    // in zero points.
767                    if start.0 < end.0 {
768                        EdgeIntersection::Overlap(Edge::new(start.1, end.1))
769                    } else if start.0 == end.0 {
770                        EdgeIntersection::EndPoint(start.1)
771                    } else {
772                        EdgeIntersection::None
773                    }
774                }
775            }
776        };
777
778        // Check that the result is consistent with the edge intersection test.
779        //        debug_assert_eq!(
780        //            result != EdgeIntersection::None,
781        //            self.edges_intersect(other)
782        //        );
783
784        debug_assert!(if self.edges_intersect(&other).inclusive_bounds() {
785            result != EdgeIntersection::None
786        } else {
787            true
788        });
789
790        result
791    }
792}
793
794impl<T: CoordinateType + NumCast> Edge<T> {
795    /// Try to cast into other data type.
796    /// When the conversion fails `None` is returned.
797    pub fn try_cast<Target>(&self) -> Option<Edge<Target>>
798    where
799        Target: CoordinateType + NumCast,
800    {
801        match (self.start.try_cast(), self.end.try_cast()) {
802            (Some(start), Some(end)) => Some(Edge { start, end }),
803            _ => None,
804        }
805    }
806
807    /// Cast to other data type.
808    ///
809    /// # Panics
810    /// Panics when the conversion fails.
811    pub fn cast<Target>(&self) -> Edge<Target>
812    where
813        Target: CoordinateType + NumCast,
814    {
815        Edge {
816            start: self.start.cast(),
817            end: self.end.cast(),
818        }
819    }
820
821    /// Cast to float.
822    ///
823    /// # Panics
824    /// Panics when the conversion fails.
825    pub fn cast_to_float<Target>(&self) -> Edge<Target>
826    where
827        Target: CoordinateType + NumCast + Float,
828    {
829        Edge {
830            start: self.start.cast_to_float(),
831            end: self.end.cast_to_float(),
832        }
833    }
834
835    /// Calculate the distance from the point to the line given by the edge.
836    ///
837    /// Distance will be positive if the point lies on the right side of the edge and negative
838    /// if the point is on the left side.
839    pub fn distance_to_line<F: Float>(&self, point: Point<T>) -> F {
840        assert!(!self.is_degenerate());
841
842        let a = self.vector();
843        let b = point - self.start;
844
845        let area = b.cross_prod(a);
846
847        F::from(area).unwrap() / a.length() // distance
848    }
849
850    /// Calculate distance from point to the edge.
851    pub fn distance<F: Float>(&self, point: Point<T>) -> F {
852        let dist_ortho: F = self.distance_to_line_abs_approx(point);
853        let dist_a: F = (point - self.start).length();
854        let dist_b = (point - self.end).length();
855
856        dist_ortho.max(dist_a.min(dist_b))
857    }
858
859    /// Find the perpendicular projection of a point onto the line of the edge.
860    pub fn projection_approx<F: Float>(&self, point: Point<T>) -> Point<F> {
861        assert!(!self.is_degenerate());
862
863        let p = (point - self.start).cast();
864
865        let d = (self.vector()).cast();
866
867        // Orthogonal projection of point onto the line.
868
869        self.start.cast() + d * d.dot(p) / d.norm2_squared()
870    }
871
872    /// Find the mirror image of `point`.
873    pub fn reflection_approx<F: Float>(&self, point: Point<T>) -> Point<F> {
874        let proj = self.projection_approx(point);
875        proj + proj - point.cast().v()
876    }
877
878    /// Calculate the absolute distance from the point onto the unbounded line coincident with this edge.
879    pub fn distance_to_line_abs_approx<F: Float>(&self, point: Point<T>) -> F {
880        let dist: F = self.distance_to_line(point);
881        dist.abs()
882    }
883
884    /// Test if point lies approximately on the edge.
885    /// Returns true if `point` is up to `tolerance` away from the edge
886    /// and lies between start and end points (inclusive).
887    pub fn contains_point_approx<F: Float>(&self, point: Point<T>, tolerance: F) -> bool {
888        debug_assert!(tolerance >= F::zero());
889        if self.is_degenerate() {
890            let l = (self.start - point).norm2_squared();
891            F::from(l).unwrap() <= tolerance
892        } else {
893            let p = point - self.start;
894
895            let v = self.vector();
896            // Rotate by -pi/4 to get the normal.
897            let n: Vector<F> = v.rotate_ortho(Angle::R270).cast() / v.length();
898            // Project on to normal
899            let o = n.dot(p.cast());
900
901            if o.abs() >= tolerance {
902                // Point is not close enough to line.
903                false
904            } else {
905                // Is point between start and end?
906
907                // Project on to line
908                let l = F::from(v.dot(p)).unwrap();
909
910                l >= -tolerance && l <= F::from(v.norm2_squared()).unwrap() + tolerance
911            }
912        }
913    }
914}
915
916impl<T: Copy> MapPointwise<T> for Edge<T> {
917    fn transform<F: Fn(Point<T>) -> Point<T>>(&self, tf: F) -> Self {
918        Edge {
919            start: tf(self.start),
920            end: tf(self.end),
921        }
922    }
923}
924
925impl<T: Copy + PartialOrd> BoundingBox<T> for Edge<T> {
926    fn bounding_box(&self) -> Rect<T> {
927        Rect::new(self.start, self.end)
928    }
929}
930
931impl<T: Copy + PartialOrd> TryBoundingBox<T> for Edge<T> {
932    /// Get bounding box of edge (always exists).
933    fn try_bounding_box(&self) -> Option<Rect<T>> {
934        Some(self.bounding_box())
935    }
936}
937
938impl<T: Copy + NumCast, Dst: Copy + NumCast> TryCastCoord<T, Dst> for Edge<T> {
939    type Output = Edge<Dst>;
940
941    fn try_cast(&self) -> Option<Self::Output> {
942        match (self.start.try_cast(), self.end.try_cast()) {
943            (Some(s), Some(e)) => Some(Edge::new(s, e)),
944            _ => None,
945        }
946    }
947}
948
949#[cfg(test)]
950mod tests {
951    extern crate rand;
952
953    use super::*;
954    use crate::edge::Edge;
955    use crate::point::Point;
956    use crate::types::*;
957    use std::f64;
958
959    use self::rand::distributions::{Bernoulli, Distribution, Uniform};
960    use self::rand::rngs::StdRng;
961    use self::rand::SeedableRng;
962
963    #[test]
964    fn test_is_parallel() {
965        let e1 = Edge::new((1, 2), (3, 4));
966        let e2 = Edge::new((1 - 10, 2 - 20), (5 - 10, 6 - 20));
967        let e3 = Edge::new((1 - 10, 2 - 20), (5 - 10, 6 - 20 + 1));
968
969        assert!(e1.is_parallel(&e2));
970        assert!(!e1.is_parallel(&e3));
971
972        assert!(e1.is_parallel_approx(&e2, 0));
973        assert!(!e1.is_parallel_approx(&e3, 0));
974
975        assert!(e1.is_parallel_approx(&e3, 1));
976    }
977
978    #[test]
979    fn test_is_collinear() {
980        let e1 = Edge::new((0, 0), (1, 2));
981        let e2 = Edge::new((10, 20), (100, 200));
982        assert!(e1.is_collinear(&e2));
983        assert!(e2.is_collinear(&e1));
984        assert!(e1.is_collinear_approx(&e2, 0));
985        assert!(e2.is_collinear_approx(&e1, 0));
986
987        // Not collinear.
988        let e1 = Edge::new((0i64, 0), (1, 2));
989        let e2 = Edge::new((10, 20), (1000, 2001));
990        assert!(!e1.is_collinear(&e2));
991        assert!(!e2.is_collinear(&e1));
992        assert!(!e1.is_collinear_approx(&e2, 0));
993        assert!(!e2.is_collinear_approx(&e1, 0));
994
995        assert!(e1.is_collinear_approx(&e2, 1));
996        assert!(e2.is_collinear_approx(&e1, 1));
997    }
998
999    #[test]
1000    fn test_distance_to_line() {
1001        let e1 = Edge::new((1, 0), (2, 1));
1002        let p0 = Point::new(2, 0);
1003
1004        let d: f64 = e1.distance_to_line(p0);
1005        let diff = (d - f64::sqrt(2.) / 2.).abs();
1006
1007        assert!(diff < 1e-9);
1008    }
1009
1010    #[test]
1011    fn test_line_contains_point() {
1012        let e1 = Edge::new((1, 2), (5, 6));
1013        let p0 = Point::new(1, 2);
1014        let p1 = Point::new(3, 4);
1015        let p2 = Point::new(5, 6);
1016        let p3 = Point::new(6, 7);
1017        let p4 = Point::new(0, 1);
1018        let p5 = Point::new(0, 0);
1019
1020        assert!(e1.line_contains_point(p0));
1021        assert!(e1.line_contains_point(p1));
1022        assert!(e1.line_contains_point(p2));
1023        assert!(e1.line_contains_point(p3));
1024        assert!(e1.line_contains_point(p4));
1025        assert!(!e1.line_contains_point(p5));
1026
1027        let e1 = Edge::new((0, 0), (1, 1));
1028        assert!(e1.line_contains_point(Point::new(2, 2)));
1029        assert!(e1.line_contains_point(Point::new(-1, -1)));
1030    }
1031
1032    #[test]
1033    fn test_contains() {
1034        let e1 = Edge::new((1, 2), (5, 6));
1035        let p0 = Point::new(1, 2);
1036        let p1 = Point::new(3, 4);
1037        let p2 = Point::new(5, 6);
1038        let p3 = Point::new(0, 0);
1039        let p4 = Point::new(6, 7);
1040
1041        assert!(e1.contains_point(p0).inclusive_bounds());
1042        assert!(e1.contains_point(p1).inclusive_bounds());
1043        assert!(e1.contains_point(p2).inclusive_bounds());
1044        assert!(!e1.contains_point(p3).inclusive_bounds());
1045        assert!(!e1.contains_point(p4).inclusive_bounds());
1046
1047        let tol = 1e-6;
1048        assert!(e1.contains_point_approx(p0, tol));
1049        assert!(e1.contains_point_approx(p1, tol));
1050        assert!(e1.contains_point_approx(p2, tol));
1051        assert!(!e1.contains_point_approx(p3, tol));
1052        assert!(!e1.contains_point_approx(p4, tol));
1053
1054        let e1 = Edge::new((0, 0), (1, 1));
1055        let p0 = Point::new(2, 2);
1056        assert!(!e1.contains_point(p0).inclusive_bounds());
1057    }
1058
1059    #[test]
1060    fn test_projection() {
1061        let e1 = Edge::new((-6., -5.), (4., 7.));
1062        let p1 = Point::new(1., 2.);
1063        let p2 = Point::new(-10., 10.);
1064
1065        let proj1 = e1.projection_approx(p1);
1066        let proj2 = e1.projection_approx(p2);
1067
1068        assert!(e1.contains_point_approx(proj1, PREC_DISTANCE));
1069        assert!(e1.contains_point_approx(proj2, PREC_DISTANCE));
1070    }
1071
1072    #[test]
1073    fn test_side_of() {
1074        let e1 = Edge::new((1, 0), (4, 4));
1075        let p1 = Point::new(-10, 0);
1076        let p2 = Point::new(10, -10);
1077        let p3 = Point::new(1, 0);
1078
1079        assert_eq!(e1.side_of(p1), Side::Left);
1080        assert_eq!(e1.side_of(p2), Side::Right);
1081        assert_eq!(e1.side_of(p3), Side::Center);
1082    }
1083
1084    #[test]
1085    fn test_crossed_by() {
1086        let e1 = Edge::new((1, 0), (4, 4));
1087        let e2 = Edge::new((1 + 1, 0), (4 + 1, 4));
1088        let e3 = Edge::new((1, 0), (4, 5));
1089        let e4 = Edge::new((2, -2), (0, 0));
1090
1091        // Coincident lines
1092        assert!(e1.crossed_by_line(&e1).inclusive_bounds());
1093
1094        // Parallel but not coincident.
1095        assert!(e1.is_parallel(&e2));
1096        assert!(!e1.crossed_by_line(&e2).inclusive_bounds());
1097
1098        // Crossing lines.
1099        assert!(e1.crossed_by_line(&e3).inclusive_bounds());
1100
1101        // crossed_by is not commutative
1102        assert!(!e1.crossed_by_line(&e4).inclusive_bounds());
1103        assert!(e4.crossed_by_line(&e1).inclusive_bounds());
1104    }
1105
1106    #[test]
1107    fn test_intersect() {
1108        let e1 = Edge::new((0, 0), (4, 4));
1109        let e2 = Edge::new((1, 0), (0, 1));
1110        let e3 = Edge::new((0, -1), (-1, 1));
1111
1112        assert!(e1.edges_intersect(&e1).inclusive_bounds());
1113        assert!(e1.edges_intersect(&e2).inclusive_bounds());
1114        assert!(!e1.crossed_by_line(&e3).inclusive_bounds());
1115        assert!(e3.crossed_by_line(&e1).inclusive_bounds());
1116        assert!(!e1.edges_intersect(&e3).inclusive_bounds());
1117
1118        // Bounding boxes overlap but edges are not touching.
1119        let e1 = Edge::new((0, 0), (8, 8));
1120        let e2 = Edge::new((4, 0), (3, 1));
1121        assert!(!e1.is_rectilinear());
1122        assert!(!e2.is_rectilinear());
1123        assert!(e1.bounding_box().touches(&e2.bounding_box()));
1124        assert!(!e1.edges_intersect(&e2).inclusive_bounds());
1125        assert!(e1.lines_intersect(&e2));
1126
1127        // Intersection at an endpoint.
1128        let e1 = Edge::new((0, 0), (4, 4));
1129        let e2 = Edge::new((1, 1), (2, 0));
1130        assert_eq!(e1.edges_intersect(&e2), ContainsResult::OnBounds);
1131        assert_eq!(e2.edges_intersect(&e1), ContainsResult::OnBounds);
1132    }
1133
1134    #[test]
1135    fn test_line_intersection() {
1136        let tol = 1e-12;
1137        let e1 = Edge::new((0, 0), (2, 2));
1138        let e2 = Edge::new((1, 0), (0, 1));
1139        let e3 = Edge::new((1, 0), (3, 2));
1140
1141        assert_eq!(
1142            e1.line_intersection_approx(&e2, tol),
1143            LineIntersection::Point(Point::new(0.5, 0.5), (1, 2, 4))
1144        );
1145        // Parallel lines should not intersect
1146        assert_eq!(
1147            e1.line_intersection_approx(&e3, tol),
1148            LineIntersection::None
1149        );
1150
1151        let e4 = Edge::new((-320., 2394.), (94., -448382.));
1152        let e5 = Edge::new((71., 133.), (-13733., 1384.));
1153
1154        if let LineIntersection::Point(intersection, _) = e4.line_intersection_approx(&e5, tol) {
1155            assert!(e4.distance_to_line::<f64>(intersection).abs() < tol);
1156            assert!(e5.distance_to_line::<f64>(intersection).abs() < tol);
1157        } else {
1158            assert!(false);
1159        }
1160
1161        // Collinear lines.
1162        let e1 = Edge::new((0., 0.), (2., 2.));
1163        let e2 = Edge::new((4., 4.), (8., 8.));
1164        assert!(!e1.is_coincident(&e2));
1165        assert!(e1.is_parallel_approx(&e2, tol));
1166        assert_eq!(
1167            e1.line_intersection_approx(&e2, tol),
1168            LineIntersection::Collinear
1169        );
1170    }
1171
1172    #[test]
1173    fn test_edge_intersection() {
1174        let tol = 1e-20;
1175        // Point intersection inside both edges.
1176        let e1 = Edge::new((0, 0), (2, 2));
1177        let e2 = Edge::new((2, 0), (0, 2));
1178        assert!(e1.edges_intersect(&e2).inclusive_bounds());
1179        assert_eq!(
1180            e1.edge_intersection_approx(&e2, tol),
1181            EdgeIntersection::Point(Point::new(1f64, 1f64))
1182        );
1183        assert_eq!(
1184            e2.edge_intersection_approx(&e1, tol),
1185            EdgeIntersection::Point(Point::new(1f64, 1f64))
1186        );
1187
1188        // Point intersection on the end of one edge.
1189        let e1 = Edge::new((0, 0), (2, 2));
1190        let e2 = Edge::new((2, 0), (1, 1));
1191        assert!(e1.edges_intersect(&e2).inclusive_bounds());
1192        assert_eq!(
1193            e1.edge_intersection_approx(&e2, tol),
1194            EdgeIntersection::EndPoint(Point::new(1, 1))
1195        );
1196        assert_eq!(
1197            e2.edge_intersection_approx(&e1, tol),
1198            EdgeIntersection::EndPoint(Point::new(1, 1))
1199        );
1200
1201        // No intersection
1202        let e1 = Edge::new((0, 0), (4, 4));
1203        let e2 = Edge::new((3, 0), (2, 1));
1204        assert!(!e1.edges_intersect(&e2).inclusive_bounds());
1205        assert_eq!(
1206            e1.edge_intersection_approx(&e2, tol),
1207            EdgeIntersection::None
1208        );
1209        assert_eq!(
1210            e2.edge_intersection_approx(&e1, tol),
1211            EdgeIntersection::None
1212        );
1213    }
1214
1215    #[test]
1216    fn test_edge_intersection_overlap() {
1217        let tol = 1e-12;
1218        // Overlapping edges. Same orientations.
1219        let e1 = Edge::new((0, 0), (2, 0));
1220        let e2 = Edge::new((1, 0), (3, 0));
1221        assert!(e1.edges_intersect(&e2).inclusive_bounds());
1222        assert!(e1.is_collinear(&e2));
1223
1224        assert_eq!(
1225            e1.edge_intersection_approx(&e2, tol),
1226            EdgeIntersection::Overlap(Edge::new((1, 0), (2, 0)))
1227        );
1228        assert_eq!(
1229            e2.edge_intersection_approx(&e1, tol),
1230            EdgeIntersection::Overlap(Edge::new((1, 0), (2, 0)))
1231        );
1232
1233        // Overlapping edges. Opposing orientations.
1234        let e1 = Edge::new((0, 0), (2, 2));
1235        let e2 = Edge::new((3, 3), (1, 1));
1236        assert!(e1.edges_intersect(&e2).inclusive_bounds());
1237        assert!(e1.is_collinear(&e2));
1238
1239        assert_eq!(
1240            e1.edge_intersection_approx(&e2, tol),
1241            EdgeIntersection::Overlap(Edge::new((1, 1), (2, 2)))
1242        );
1243        assert_eq!(
1244            e2.edge_intersection_approx(&e1, tol),
1245            EdgeIntersection::Overlap(Edge::new((2, 2), (1, 1)))
1246        );
1247
1248        // Overlapping edges. One is fully contained in the other. Same orientations.
1249        let e1 = Edge::new((0, 0), (4, 4));
1250        let e2 = Edge::new((1, 1), (2, 2));
1251        assert!(e1.edges_intersect(&e2).inclusive_bounds());
1252        assert!(e1.is_collinear(&e2));
1253
1254        assert_eq!(
1255            e1.edge_intersection_approx(&e2, tol),
1256            EdgeIntersection::Overlap(e2)
1257        );
1258        assert_eq!(
1259            e2.edge_intersection_approx(&e1, tol),
1260            EdgeIntersection::Overlap(e2)
1261        );
1262
1263        // Overlapping edges. One is fully contained in the other. Opposing orientations.
1264        let e1 = Edge::new((0, 0), (4, 4));
1265        let e2 = Edge::new((2, 2), (1, 1));
1266        assert!(e1.edges_intersect(&e2).inclusive_bounds());
1267        assert!(e1.is_collinear(&e2));
1268
1269        assert_eq!(
1270            e1.edge_intersection_approx(&e2, tol),
1271            EdgeIntersection::Overlap(e2.reversed())
1272        );
1273        assert_eq!(
1274            e2.edge_intersection_approx(&e1, tol),
1275            EdgeIntersection::Overlap(e2)
1276        );
1277
1278        // Collinear edge, touch in exactly one point.
1279        let e1 = Edge::new((0, 0), (1, 1));
1280        let e2 = Edge::new((1, 1), (2, 2));
1281        assert!(e1.edges_intersect(&e2).inclusive_bounds());
1282        assert!(e1.is_collinear(&e2));
1283
1284        assert_eq!(
1285            e1.edge_intersection_approx(&e2, tol),
1286            EdgeIntersection::EndPoint(Point::new(1, 1))
1287        );
1288        assert_eq!(
1289            e2.edge_intersection_approx(&e1, tol),
1290            EdgeIntersection::EndPoint(Point::new(1, 1))
1291        );
1292
1293        // Edges touch in exactly one point. Floating point.
1294        let e1 = Edge::new((1.1, 1.2), (0., 0.));
1295        let e2 = Edge::new((1.1, 1.2), (2., 2.));
1296        assert!(e1.edges_intersect(&e2).inclusive_bounds());
1297        assert!(e2.edges_intersect(&e1).inclusive_bounds());
1298
1299        assert_eq!(
1300            e1.edge_intersection_approx(&e2, tol),
1301            EdgeIntersection::EndPoint(Point::new(1.1, 1.2))
1302        );
1303        assert_eq!(
1304            e2.edge_intersection_approx(&e1, tol),
1305            EdgeIntersection::EndPoint(Point::new(1.1, 1.2))
1306        );
1307
1308        // Intersection on endpoint with floats.
1309        let e1 = Edge {
1310            start: Point::new(20.90725737794763, 65.33301386746126),
1311            end: Point::new(32.799556599584776, 63.13131890182373),
1312        };
1313        let e2 = Edge {
1314            start: Point::new(22.217533978705163, 70.84660296990562),
1315            end: Point::new(32.799556599584776, 63.13131890182373),
1316        };
1317        assert!(e1.edges_intersect(&e2).inclusive_bounds());
1318        assert!(e2.edges_intersect(&e1).inclusive_bounds());
1319
1320        assert_eq!(
1321            e1.edge_intersection_approx(&e2, tol),
1322            EdgeIntersection::EndPoint(Point::new(32.799556599584776, 63.13131890182373))
1323        );
1324        assert_eq!(
1325            e2.edge_intersection_approx(&e1, tol),
1326            EdgeIntersection::EndPoint(Point::new(32.799556599584776, 63.13131890182373))
1327        );
1328    }
1329
1330    #[test]
1331    fn test_coincident() {
1332        let e1 = Edge::new((0, 0), (2, 2));
1333        let e2 = Edge::new((1, 1), (3, 3));
1334        assert!(e1.is_coincident(&e2));
1335        assert!(e2.is_coincident(&e1));
1336
1337        let e1 = Edge::new((0, 0), (3, 3));
1338        let e2 = Edge::new((1, 1), (2, 2));
1339        assert!(e1.is_coincident(&e2));
1340        assert!(e2.is_coincident(&e1));
1341
1342        let e1 = Edge::new((0, 0), (1, 1));
1343        let e2 = Edge::new((1, 1), (2, 2));
1344        assert!(!e1.is_coincident(&e2));
1345        assert!(!e2.is_coincident(&e1));
1346    }
1347
1348    #[test]
1349    fn test_intersection_of_random_edges() {
1350        let tol = 1e-12;
1351        let seed1 = [1u8; 32];
1352        let seed2 = [2u8; 32];
1353
1354        let between = Uniform::from(-1.0..1.0);
1355        let mut rng = StdRng::from_seed(seed1);
1356
1357        let mut rand_edge = || -> Edge<f64> {
1358            let points: Vec<(f64, f64)> = (0..2)
1359                .into_iter()
1360                .map(|_| (between.sample(&mut rng), between.sample(&mut rng)))
1361                .collect();
1362
1363            Edge::new(points[0], points[1])
1364        };
1365
1366        let bernoulli_rare = Bernoulli::new(0.05).unwrap();
1367        let bernoulli_50 = Bernoulli::new(0.5).unwrap();
1368        let mut rng = StdRng::from_seed(seed2);
1369
1370        for _i in 0..1000 {
1371            // Create a random pair of edges. 5% of the edge pairs share an endpoint.
1372            let (a, b) = {
1373                let a = rand_edge();
1374
1375                let b = {
1376                    let b = rand_edge();
1377
1378                    if bernoulli_rare.sample(&mut rng) {
1379                        // Share a point with the other edge.
1380                        let shared = if bernoulli_50.sample(&mut rng) {
1381                            a.start
1382                        } else {
1383                            a.end
1384                        };
1385
1386                        let result = Edge::new(b.start, shared);
1387                        if bernoulli_50.sample(&mut rng) {
1388                            result
1389                        } else {
1390                            result.reversed()
1391                        }
1392                    } else {
1393                        b
1394                    }
1395                };
1396                (a, b)
1397            };
1398
1399            let intersection_ab = a.edge_intersection_approx(&b, tol);
1400            assert_eq!(
1401                intersection_ab != EdgeIntersection::None,
1402                a.edges_intersect(&b).inclusive_bounds()
1403            );
1404
1405            let intersection_ba = b.edge_intersection_approx(&a, tol);
1406            assert_eq!(
1407                intersection_ba != EdgeIntersection::None,
1408                b.edges_intersect(&a).inclusive_bounds()
1409            );
1410        }
1411    }
1412}