Oregami
Repositories/oxedyne/fe2o3

oxedyne/fe2o3/fe2o3_geom/src/planar.rs

47.2 KiB, 9 runs

created by r1870400018:17290, which is this file's identity for as long as the history lasts, whatever it is later renamed to

download · who wrote it · its history

1//! Floating-point planar geometry: primitives and predicates for a dynamic-geometry editor.
2//!
3//! This module provides the continuous, real-valued counterpart to the integer UI-layout
4//! types elsewhere in the crate. Where `dim`, `rect` and `shape` serve pixel-aligned widget
5//! layout, this module serves geometric construction and constraint solving: the kind of
6//! work a compass-and-straightedge editor performs when it places a point on a line, hangs a
7//! circle off two others, or asks where two loci meet.
8//!
9//! The primitives are `Pt` (a point, also exported as `Point`), `Vec2` (a free vector),
10//! `Line` (an infinite line), `Ray` (a half-line), `Segment` (a bounded line), `Circle` and
11//! `Arc`. `Angle` and `Vec2` are first-class *types* rather than incidental pairs of floats,
12//! so that a constraint and a back-solver can pass them around without re-deriving their
13//! invariants each time.
14//!
15//! The predicates are the ones a constraint check and a back-solve both call: incidence of a
16//! point on a line, intersection of line with line, line with circle, and circle with circle,
17//! the foot of a perpendicular, the projection of a point onto a line or a (clamped) segment,
18//! and the angle between two rays.
19//!
20//! # Tolerance
21//!
22//! Floating-point coordinates never compare exactly, so every predicate that answers a
23//! yes/no or counts intersections takes an explicit `eps` tolerance, measured in the same
24//! units as the coordinates (a *distance*, not a raw coordinate difference). A sensible
25//! default is exposed as [`EPSILON`]. Because line and ray directions are normalised on
26//! construction, the tolerance stays scale-independent: it is compared against genuine
27//! perpendicular distances and against sines of angles, not against raw cross products.
28
29use oxedyne_fe2o3_core::prelude::*;
30
31use std::{
32 f64::consts::PI,
33 ops::{
34 Add,
35 Div,
36 Mul,
37 Neg,
38 Sub,
39 },
40};
41
42/// The default distance tolerance for planar predicates.
43///
44/// Two points closer than this are treated as coincident, a point this near a line is treated
45/// as lying on it, and an intersection this close to a tangent is treated as a single point.
46pub const EPSILON: f64 = 1.0e-9;
47
48/// Returns `true` when `a` and `b` are within `eps` of one another.
49///
50/// A small free helper so callers need not repeat the absolute-difference idiom.
51pub fn approx_eq(a: f64, b: f64, eps: f64) -> bool {
52 (a - b).abs() <= eps
53}
54
55// ---------------------------------------------------------------------------------------------
56// Vec2
57// ---------------------------------------------------------------------------------------------
58
59/// A free vector in the plane, with `f64` components.
60///
61/// Distinct from [`Pt`]: a `Vec2` is a displacement or direction, not a location. The type
62/// distinction lets the arithmetic express intent -- subtracting two points yields a `Vec2`,
63/// and adding a `Vec2` to a point yields a point.
64#[derive(Clone, Copy, Debug, PartialEq)]
65pub struct Vec2 {
66 /// Horizontal component.
67 pub x: f64,
68 /// Vertical component.
69 pub y: f64,
70}
71
72impl Vec2 {
73 /// Creates a new vector from its components.
74 pub fn new(x: f64, y: f64) -> Self {
75 Self { x, y }
76 }
77
78 /// Returns the dot (scalar) product with `other`.
79 pub fn dot(&self, other: &Vec2) -> f64 {
80 self.x * other.x + self.y * other.y
81 }
82
83 /// Returns the `z` component of the 3D cross product, i.e. the signed area of the
84 /// parallelogram spanned by the two vectors.
85 ///
86 /// Positive when `other` lies anticlockwise of `self`. Zero (within tolerance) means the
87 /// two vectors are parallel.
88 pub fn cross(&self, other: &Vec2) -> f64 {
89 self.x * other.y - self.y * other.x
90 }
91
92 /// Returns the squared length, avoiding the square root where only comparison is needed.
93 pub fn length_sq(&self) -> f64 {
94 self.x * self.x + self.y * self.y
95 }
96
97 /// Returns the Euclidean length (magnitude).
98 pub fn length(&self) -> f64 {
99 self.length_sq().sqrt()
100 }
101
102 /// Returns the unit vector in the same direction.
103 ///
104 /// # Errors
105 /// Fails when the vector is shorter than [`EPSILON`], as a zero vector has no direction.
106 pub fn normalise(&self) -> Outcome<Vec2> {
107 let len = self.length();
108 if len < EPSILON {
109 return Err(err!(
110 "Cannot normalise a zero-length vector.";
111 Invalid, Numeric, Range));
112 }
113 Ok(Vec2::new(self.x / len, self.y / len))
114 }
115
116 /// Returns the vector rotated a quarter turn anticlockwise (the left normal).
117 pub fn perp(&self) -> Vec2 {
118 Vec2::new(-self.y, self.x)
119 }
120
121 /// Returns the direction of the vector as an [`Angle`] measured from the positive `x` axis.
122 pub fn angle(&self) -> Angle {
123 Angle::from_radians(self.y.atan2(self.x))
124 }
125
126 /// Returns `true` when both components are within `eps` of `other`'s.
127 pub fn approx_eq(&self, other: &Vec2, eps: f64) -> bool {
128 (*self - *other).length() <= eps
129 }
130}
131
132impl Add for Vec2 {
133 type Output = Vec2;
134
135 fn add(self, other: Vec2) -> Vec2 {
136 Vec2::new(self.x + other.x, self.y + other.y)
137 }
138}
139
140impl Sub for Vec2 {
141 type Output = Vec2;
142
143 fn sub(self, other: Vec2) -> Vec2 {
144 Vec2::new(self.x - other.x, self.y - other.y)
145 }
146}
147
148impl Neg for Vec2 {
149 type Output = Vec2;
150
151 fn neg(self) -> Vec2 {
152 Vec2::new(-self.x, -self.y)
153 }
154}
155
156impl Mul<f64> for Vec2 {
157 type Output = Vec2;
158
159 fn mul(self, scalar: f64) -> Vec2 {
160 Vec2::new(self.x * scalar, self.y * scalar)
161 }
162}
163
164impl Div<f64> for Vec2 {
165 type Output = Vec2;
166
167 fn div(self, scalar: f64) -> Vec2 {
168 Vec2::new(self.x / scalar, self.y / scalar)
169 }
170}
171
172// ---------------------------------------------------------------------------------------------
173// Pt
174// ---------------------------------------------------------------------------------------------
175
176/// A point (location) in the plane, with `f64` coordinates.
177///
178/// Also exported as [`Point`] for callers who prefer the fuller name.
179#[derive(Clone, Copy, Debug, PartialEq)]
180pub struct Pt {
181 /// Horizontal coordinate.
182 pub x: f64,
183 /// Vertical coordinate.
184 pub y: f64,
185}
186
187/// A fuller alias for [`Pt`].
188pub type Point = Pt;
189
190impl Pt {
191 /// Creates a new point from its coordinates.
192 pub fn new(x: f64, y: f64) -> Self {
193 Self { x, y }
194 }
195
196 /// Returns the point as a position vector from the origin.
197 pub fn to_vec(&self) -> Vec2 {
198 Vec2::new(self.x, self.y)
199 }
200
201 /// Returns the squared distance to `other`, avoiding the square root where only
202 /// comparison is needed.
203 pub fn distance_sq(&self, other: &Pt) -> f64 {
204 (*self - *other).length_sq()
205 }
206
207 /// Returns the Euclidean distance to `other`.
208 pub fn distance(&self, other: &Pt) -> f64 {
209 (*self - *other).length()
210 }
211
212 /// Returns the midpoint between this point and `other`.
213 pub fn midpoint(&self, other: &Pt) -> Pt {
214 Pt::new((self.x + other.x) / 2.0, (self.y + other.y) / 2.0)
215 }
216
217 /// Returns `true` when the two points are within `eps` of one another.
218 pub fn approx_eq(&self, other: &Pt, eps: f64) -> bool {
219 self.distance(other) <= eps
220 }
221}
222
223impl Sub for Pt {
224 type Output = Vec2;
225
226 /// Point minus point is the displacement between them.
227 fn sub(self, other: Pt) -> Vec2 {
228 Vec2::new(self.x - other.x, self.y - other.y)
229 }
230}
231
232impl Add<Vec2> for Pt {
233 type Output = Pt;
234
235 /// Point plus vector is the translated point.
236 fn add(self, v: Vec2) -> Pt {
237 Pt::new(self.x + v.x, self.y + v.y)
238 }
239}
240
241impl Sub<Vec2> for Pt {
242 type Output = Pt;
243
244 /// Point minus vector is the point translated in the opposite direction.
245 fn sub(self, v: Vec2) -> Pt {
246 Pt::new(self.x - v.x, self.y - v.y)
247 }
248}
249
250// ---------------------------------------------------------------------------------------------
251// Angle
252// ---------------------------------------------------------------------------------------------
253
254/// An angle, stored internally in radians.
255///
256/// A first-class type so that a bare `f64` cannot be mistaken for degrees, and so that
257/// normalisation and trigonometry live in one place.
258#[derive(Clone, Copy, Debug, PartialEq)]
259pub struct Angle {
260 /// The angle in radians.
261 rad: f64,
262}
263
264impl Angle {
265 /// Creates an angle from a value in radians.
266 pub fn from_radians(rad: f64) -> Self {
267 Self { rad }
268 }
269
270 /// Creates an angle from a value in degrees.
271 pub fn from_degrees(deg: f64) -> Self {
272 Self { rad: deg * PI / 180.0 }
273 }
274
275 /// Returns the angle in radians.
276 pub fn radians(&self) -> f64 {
277 self.rad
278 }
279
280 /// Returns the angle in degrees.
281 pub fn degrees(&self) -> f64 {
282 self.rad * 180.0 / PI
283 }
284
285 /// Returns the sine of the angle.
286 pub fn sin(&self) -> f64 {
287 self.rad.sin()
288 }
289
290 /// Returns the cosine of the angle.
291 pub fn cos(&self) -> f64 {
292 self.rad.cos()
293 }
294
295 /// Returns the angle normalised into the half-open interval `[0, 2π)`.
296 pub fn normalised(&self) -> Angle {
297 let two_pi = 2.0 * PI;
298 let mut r = self.rad % two_pi;
299 if r < 0.0 {
300 r += two_pi;
301 }
302 Angle::from_radians(r)
303 }
304
305 /// Returns `true` when the two angles are within `eps` radians of one another, comparing
306 /// on the circle so that values astride the `2π` wrap are still recognised as equal.
307 pub fn approx_eq(&self, other: &Angle, eps: f64) -> bool {
308 let two_pi = 2.0 * PI;
309 let mut d = (self.rad - other.rad).abs() % two_pi;
310 if d > PI {
311 d = two_pi - d;
312 }
313 d <= eps
314 }
315}
316
317// ---------------------------------------------------------------------------------------------
318// Line
319// ---------------------------------------------------------------------------------------------
320
321/// An infinite line, stored as a point on the line and a unit direction.
322///
323/// The direction is normalised on construction, which keeps perpendicular-distance and
324/// angle tolerances scale-independent.
325#[derive(Clone, Copy, Debug)]
326pub struct Line {
327 /// A point through which the line passes.
328 pub origin: Pt,
329 /// The unit direction of the line.
330 pub dir: Vec2,
331}
332
333impl Line {
334 /// Creates a line through `origin` in direction `dir`.
335 ///
336 /// # Errors
337 /// Fails when `dir` has zero length, as a line needs a direction.
338 pub fn new(origin: Pt, dir: Vec2) -> Outcome<Self> {
339 let unit = res!(dir.normalise());
340 Ok(Self { origin, dir: unit })
341 }
342
343 /// Creates a line through the two given points.
344 ///
345 /// # Errors
346 /// Fails when the two points are within [`EPSILON`] of each other, as coincident points
347 /// do not determine a line.
348 pub fn through(a: Pt, b: Pt) -> Outcome<Self> {
349 if a.approx_eq(&b, EPSILON) {
350 return Err(err!(
351 "Cannot build a line through two coincident points {:?} and {:?}.", a, b;
352 Invalid, Input, Range));
353 }
354 Line::new(a, b - a)
355 }
356
357 /// Returns the signed perpendicular distance from `pt` to the line.
358 ///
359 /// The sign follows the left normal: positive when `pt` lies anticlockwise of the
360 /// direction.
361 pub fn signed_distance(&self, pt: &Pt) -> f64 {
362 let w = *pt - self.origin;
363 w.cross(&self.dir).neg()
364 }
365
366 /// Returns the (unsigned) perpendicular distance from `pt` to the line.
367 pub fn distance_to(&self, pt: &Pt) -> f64 {
368 let w = *pt - self.origin;
369 w.cross(&self.dir).abs()
370 }
371
372 /// Returns `true` when `pt` lies on the line, within perpendicular distance `eps`.
373 pub fn contains(&self, pt: &Pt, eps: f64) -> bool {
374 self.distance_to(pt) <= eps
375 }
376
377 /// Returns the foot of the perpendicular dropped from `pt` onto the line: the closest
378 /// point of the line to `pt`.
379 pub fn foot_of_perpendicular(&self, pt: &Pt) -> Pt {
380 let w = *pt - self.origin;
381 let t = w.dot(&self.dir); // Projection parameter along the unit direction.
382 self.origin + self.dir * t
383 }
384
385 /// Returns the orthogonal projection of `pt` onto the line.
386 ///
387 /// For an infinite line this is exactly the foot of the perpendicular; the alias exists
388 /// so callers can name the operation as a projection.
389 pub fn project_point(&self, pt: &Pt) -> Pt {
390 self.foot_of_perpendicular(pt)
391 }
392
393 /// Intersects this line with `other`.
394 ///
395 /// Returns a single point, or reports the lines parallel (never meeting) or coincident
396 /// (meeting everywhere). Parallelism is judged by the sine of the angle between the unit
397 /// directions against `eps`; coincidence additionally requires `other`'s origin to lie
398 /// within `eps` of this line.
399 pub fn intersect_line(&self, other: &Line, eps: f64) -> LineIntersect {
400 let denom = self.dir.cross(&other.dir); // sin of the angle between unit directions.
401 if denom.abs() <= eps {
402 // Directions parallel; decide coincident versus strictly parallel.
403 if self.distance_to(&other.origin) <= eps {
404 return LineIntersect::Coincident;
405 }
406 return LineIntersect::Parallel;
407 }
408 let w = other.origin - self.origin;
409 let t = w.cross(&other.dir) / denom; // Parameter along this line.
410 LineIntersect::Point(self.origin + self.dir * t)
411 }
412
413 /// Intersects this line with `circle`, returning zero, one (tangent) or two points.
414 ///
415 /// The point pair is ordered along the line's direction (the smaller parameter first).
416 pub fn intersect_circle(&self, circle: &Circle, eps: f64) -> LineCircleIntersect {
417 // Parameter of the foot of the perpendicular from the centre onto the line.
418 let w = circle.centre - self.origin;
419 let t0 = w.dot(&self.dir);
420 let foot = self.origin + self.dir * t0;
421 let d = foot.distance(&circle.centre); // Distance from centre to the line.
422 if d > circle.radius + eps {
423 return LineCircleIntersect::None;
424 }
425 if approx_eq(d, circle.radius, eps) {
426 return LineCircleIntersect::Tangent(foot);
427 }
428 // Half-chord length; clamp the radicand to guard against a tiny negative from noise.
429 let h = (circle.radius * circle.radius - d * d).max(0.0).sqrt();
430 let p0 = self.origin + self.dir * (t0 - h);
431 let p1 = self.origin + self.dir * (t0 + h);
432 LineCircleIntersect::Secant(p0, p1)
433 }
434}
435
436/// The outcome of intersecting two [`Line`]s.
437#[derive(Clone, Copy, Debug, PartialEq)]
438pub enum LineIntersect {
439 /// The lines cross at exactly one point.
440 Point(Pt),
441 /// The lines are parallel and distinct, so they never meet.
442 Parallel,
443 /// The lines are the same line, so they meet everywhere.
444 Coincident,
445}
446
447/// The outcome of intersecting a [`Line`] with a [`Circle`].
448#[derive(Clone, Copy, Debug, PartialEq)]
449pub enum LineCircleIntersect {
450 /// The line misses the circle.
451 None,
452 /// The line touches the circle at a single point.
453 Tangent(Pt),
454 /// The line cuts the circle at two points, ordered along the line's direction.
455 Secant(Pt, Pt),
456}
457
458// ---------------------------------------------------------------------------------------------
459// Ray
460// ---------------------------------------------------------------------------------------------
461
462/// A half-line: an origin and a unit direction, extending infinitely one way.
463#[derive(Clone, Copy, Debug)]
464pub struct Ray {
465 /// The endpoint from which the ray extends.
466 pub origin: Pt,
467 /// The unit direction of the ray.
468 pub dir: Vec2,
469}
470
471impl Ray {
472 /// Creates a ray from `origin` in direction `dir`.
473 ///
474 /// # Errors
475 /// Fails when `dir` has zero length, as a ray needs a direction.
476 pub fn new(origin: Pt, dir: Vec2) -> Outcome<Self> {
477 let unit = res!(dir.normalise());
478 Ok(Self { origin, dir: unit })
479 }
480
481 /// Creates a ray from `origin` pointing towards `through`.
482 ///
483 /// # Errors
484 /// Fails when the two points are within [`EPSILON`] of each other.
485 pub fn towards(origin: Pt, through: Pt) -> Outcome<Self> {
486 if origin.approx_eq(&through, EPSILON) {
487 return Err(err!(
488 "Cannot build a ray from {:?} towards a coincident point {:?}.", origin, through;
489 Invalid, Input, Range));
490 }
491 Ray::new(origin, through - origin)
492 }
493
494 /// Returns the angle between this ray and `other`, in the range `[0, π]`.
495 ///
496 /// Both rays carry unit directions, so this is the arccosine of their dot product,
497 /// clamped to guard against floating-point drift past the ends of the domain.
498 pub fn angle_between(&self, other: &Ray) -> Angle {
499 let d = self.dir.dot(&other.dir).clamp(-1.0, 1.0);
500 Angle::from_radians(d.acos())
501 }
502}
503
504// ---------------------------------------------------------------------------------------------
505// Segment
506// ---------------------------------------------------------------------------------------------
507
508/// A bounded line between two endpoints.
509#[derive(Clone, Copy, Debug, PartialEq)]
510pub struct Segment {
511 /// The start endpoint.
512 pub a: Pt,
513 /// The end endpoint.
514 pub b: Pt,
515}
516
517impl Segment {
518 /// Creates a segment between the two endpoints.
519 ///
520 /// # Errors
521 /// Fails when the endpoints are within [`EPSILON`] of each other, as a degenerate segment
522 /// has no direction and cannot be projected onto meaningfully.
523 pub fn new(a: Pt, b: Pt) -> Outcome<Self> {
524 if a.approx_eq(&b, EPSILON) {
525 return Err(err!(
526 "Cannot build a segment between two coincident points {:?} and {:?}.", a, b;
527 Invalid, Input, Range));
528 }
529 Ok(Self { a, b })
530 }
531
532 /// Returns the length of the segment.
533 pub fn length(&self) -> f64 {
534 self.a.distance(&self.b)
535 }
536
537 /// Returns the projection of `pt` onto the segment, clamped to lie between the endpoints.
538 ///
539 /// Unlike the projection onto an infinite line, the parameter is clamped to `[0, 1]`, so
540 /// a point beyond an end projects to that end.
541 pub fn project_point(&self, pt: &Pt) -> Pt {
542 let ab = self.b - self.a;
543 let l2 = ab.length_sq(); // Non-zero: the constructor forbids a degenerate segment.
544 let t = ((*pt - self.a).dot(&ab) / l2).clamp(0.0, 1.0);
545 self.a + ab * t
546 }
547
548 /// The point a fraction `t` of the way from `a` to `b`, unclamped.
549 pub fn point_at(&self, t: f64) -> Pt {
550 self.a + (self.b - self.a) * t
551 }
552
553 /// The infinite line this segment lies on.
554 ///
555 /// # Errors
556 /// Fails when the segment is degenerate, which its constructor already forbids.
557 pub fn to_line(&self) -> Outcome<Line> {
558 Line::through(self.a, self.b)
559 }
560
561 /// Crosses this segment with an infinite line, returning where along the *segment* they
562 /// meet as a fraction from `a` to `b`, or `None` when they do not meet within it.
563 ///
564 /// This is the primitive a scanline needs: sweeping a family of parallel lines across a
565 /// shape and asking each of its edges where it is cut. A crossing exactly at an endpoint
566 /// counts, since dropping it would open a gap in the sweep; an edge lying *along* the line
567 /// does not, since it has no single crossing to report.
568 pub fn cross_line(&self, line: &Line, eps: f64) -> Option<f64> {
569 let ab = self.b - self.a;
570 let denom = line.dir.cross(&ab);
571 if denom.abs() <= eps {
572 return None; // Parallel to the line, coincident or not.
573 }
574 let t = line.dir.cross(&(line.origin - self.a)) / denom;
575 if t < -eps || t > 1.0 + eps {
576 return None;
577 }
578 Some(t.clamp(0.0, 1.0))
579 }
580
581 /// Crosses this segment with another, returning the meeting point when both are cut
582 /// within their own extents.
583 ///
584 /// Parallel and collinear pairs report nothing: a collinear overlap has no single point
585 /// to name, and a caller wanting the overlap wants a different question answered.
586 pub fn intersect_segment(&self, other: &Segment, eps: f64) -> Option<Pt> {
587 let (r, s) = (self.b - self.a, other.b - other.a);
588 let denom = r.cross(&s);
589 if denom.abs() <= eps {
590 return None;
591 }
592 let w = other.a - self.a;
593 let t = w.cross(&s) / denom;
594 let u = w.cross(&r) / denom;
595 if t < -eps || t > 1.0 + eps || u < -eps || u > 1.0 + eps {
596 return None;
597 }
598 Some(self.point_at(t.clamp(0.0, 1.0)))
599 }
600}
601
602// ---------------------------------------------------------------------------------------------
603// Triangle
604// ---------------------------------------------------------------------------------------------
605
606/// A triangle, stored as its three corners in the order they were given.
607///
608/// Its reason for existing here is [`Triangle::barycentric`] and its inverse
609/// [`Triangle::point_at`], which together are the standard way to carry a point from one
610/// triangle into another: read the point's weights in the first, then rebuild it from the
611/// second's corners. That is how a triangulated region warps whatever lies inside it, and it
612/// is also how a value sampled at three points is interpolated across the space between them.
613#[derive(Clone, Copy, Debug, PartialEq)]
614pub struct Triangle {
615 /// The first corner.
616 pub a: Pt,
617 /// The second corner.
618 pub b: Pt,
619 /// The third corner.
620 pub c: Pt,
621}
622
623impl Triangle {
624 /// Creates a triangle from three corners.
625 ///
626 /// # Errors
627 /// Fails when the three corners are collinear to within [`EPSILON`], since a degenerate
628 /// triangle has no interior and no barycentric coordinates.
629 pub fn new(a: Pt, b: Pt, c: Pt) -> Outcome<Self> {
630 let t = Self { a, b, c };
631 if t.signed_area().abs() <= EPSILON {
632 return Err(err!(
633 "Cannot build a triangle from the collinear points {:?}, {:?} and {:?}.", a, b, c;
634 Invalid, Input, Range));
635 }
636 Ok(t)
637 }
638
639 /// Returns twice the signed area, positive when the corners wind anticlockwise in a
640 /// y-up frame. This is the cross product of two edges, so it is also the quantity a
641 /// degeneracy test and an orientation test both want.
642 pub fn signed_area2(&self) -> f64 {
643 (self.b - self.a).cross(&(self.c - self.a))
644 }
645
646 /// Returns the signed area, positive when the corners wind anticlockwise in a y-up frame.
647 pub fn signed_area(&self) -> f64 {
648 self.signed_area2() / 2.0
649 }
650
651 /// Returns the area.
652 pub fn area(&self) -> f64 {
653 self.signed_area().abs()
654 }
655
656 /// Returns the centroid, the average of the three corners.
657 pub fn centroid(&self) -> Pt {
658 Pt::new(
659 (self.a.x + self.b.x + self.c.x) / 3.0,
660 (self.a.y + self.b.y + self.c.y) / 3.0,
661 )
662 }
663
664 /// Returns the barycentric coordinates of `pt`: the three weights that sum to one and
665 /// rebuild the point from the corners.
666 ///
667 /// All three are non-negative exactly when the point lies inside the triangle or on its
668 /// boundary, which is what [`Triangle::contains`] tests. A negative weight says which
669 /// side the point fell out of, so the value is useful beyond the containment question.
670 ///
671 /// # Errors
672 /// Fails on a degenerate triangle, which the constructor already forbids but which a
673 /// caller can still produce by mutating the fields.
674 pub fn barycentric(&self, pt: &Pt) -> Outcome<(f64, f64, f64)> {
675 let d = self.signed_area2();
676 if d.abs() <= EPSILON {
677 return Err(err!(
678 "Cannot take barycentric coordinates in a degenerate triangle {:?}.", self;
679 Invalid, Input, Range));
680 }
681 // Each weight is the signed area of the sub-triangle opposite its corner, over the
682 // whole, which is the definition rather than a rearrangement of it.
683 let u = (self.b - *pt).cross(&(self.c - *pt)) / d;
684 let v = (self.c - *pt).cross(&(self.a - *pt)) / d;
685 Ok((u, v, 1.0 - u - v))
686 }
687
688 /// Returns the point with the given barycentric weights, the inverse of
689 /// [`Triangle::barycentric`]. The weights are not required to sum to one, so a caller
690 /// extrapolating deliberately outside the triangle gets the point it asked for.
691 pub fn point_at(&self, u: f64, v: f64, w: f64) -> Pt {
692 Pt::new(
693 self.a.x * u + self.b.x * v + self.c.x * w,
694 self.a.y * u + self.b.y * v + self.c.y * w,
695 )
696 }
697
698 /// Whether `pt` lies inside the triangle or on its boundary, to within `eps` of the
699 /// boundary. The tolerance is in barycentric weight, so it scales with the triangle
700 /// rather than with the coordinate system.
701 pub fn contains(&self, pt: &Pt, eps: f64) -> bool {
702 match self.barycentric(pt) {
703 Ok((u, v, w)) => u >= -eps && v >= -eps && w >= -eps,
704 Err(_) => false,
705 }
706 }
707
708 /// Carries a point from this triangle into `other`, by reading its weights here and
709 /// rebuilding it there. This is the mapping a warp is made of, and it is affine, so a
710 /// straight line inside the triangle stays straight.
711 ///
712 /// # Errors
713 /// Propagates the degeneracy failure of [`Triangle::barycentric`].
714 pub fn map_point(&self, other: &Triangle, pt: &Pt) -> Outcome<Pt> {
715 let (u, v, w) = res!(self.barycentric(pt));
716 Ok(other.point_at(u, v, w))
717 }
718}
719
720// ---------------------------------------------------------------------------------------------
721// Circle
722// ---------------------------------------------------------------------------------------------
723
724/// A circle, stored as a centre and a radius.
725#[derive(Clone, Copy, Debug, PartialEq)]
726pub struct Circle {
727 /// The centre of the circle.
728 pub centre: Pt,
729 /// The radius of the circle.
730 pub radius: f64,
731}
732
733impl Circle {
734 /// Creates a circle from a centre and a radius.
735 ///
736 /// # Errors
737 /// Fails when the radius is not strictly positive (greater than [`EPSILON`]).
738 pub fn new(centre: Pt, radius: f64) -> Outcome<Self> {
739 if radius <= EPSILON {
740 return Err(err!(
741 "Circle radius must be positive, got {}.", radius;
742 Invalid, Input, Range));
743 }
744 Ok(Self { centre, radius })
745 }
746
747 /// Returns `true` when `pt` lies on the circle, within `eps` of the boundary.
748 pub fn contains(&self, pt: &Pt, eps: f64) -> bool {
749 approx_eq(self.centre.distance(pt), self.radius, eps)
750 }
751
752 /// Intersects this circle with `line` (delegates to [`Line::intersect_circle`]).
753 pub fn intersect_line(&self, line: &Line, eps: f64) -> LineCircleIntersect {
754 line.intersect_circle(self, eps)
755 }
756
757 /// Intersects this circle with `other`.
758 ///
759 /// Returns zero points (separate, one wholly inside the other, or concentric), a single
760 /// tangent point, or two points. Two identical circles are reported as coincident, having
761 /// infinitely many common points.
762 pub fn intersect_circle(&self, other: &Circle, eps: f64) -> CircleIntersect {
763 let between = other.centre - self.centre;
764 let d = between.length(); // Distance between centres.
765 let r0 = self.radius;
766 let r1 = other.radius;
767 if d <= eps {
768 // Concentric centres.
769 if approx_eq(r0, r1, eps) {
770 return CircleIntersect::Coincident;
771 }
772 return CircleIntersect::None;
773 }
774 if d > r0 + r1 + eps {
775 return CircleIntersect::None; // Too far apart to meet.
776 }
777 if d < (r0 - r1).abs() - eps {
778 return CircleIntersect::None; // One circle lies wholly inside the other.
779 }
780 // Tangent: externally when d == r0 + r1, internally when d == |r0 - r1|.
781 if approx_eq(d, r0 + r1, eps) || approx_eq(d, (r0 - r1).abs(), eps) {
782 let t = r0 / d;
783 return CircleIntersect::Tangent(self.centre + between * t);
784 }
785 // Two intersection points. `a` is the distance from this centre to the chord midpoint
786 // along the line of centres; `h` is the half-chord perpendicular to it.
787 let a = (d * d + r0 * r0 - r1 * r1) / (2.0 * d);
788 let h = (r0 * r0 - a * a).max(0.0).sqrt();
789 let unit = between / d; // Unit vector along the line of centres.
790 let mid = self.centre + unit * a; // Midpoint of the common chord.
791 let perp = unit.perp() * h; // Half-chord offset.
792 CircleIntersect::Two(mid + perp, mid - perp)
793 }
794}
795
796/// The outcome of intersecting two [`Circle`]s.
797#[derive(Clone, Copy, Debug, PartialEq)]
798pub enum CircleIntersect {
799 /// The circles do not meet (separate, nested, or concentric with unequal radii).
800 None,
801 /// The circles touch at a single point.
802 Tangent(Pt),
803 /// The circles cross at two points.
804 Two(Pt, Pt),
805 /// The circles are identical, meeting at every point.
806 Coincident,
807}
808
809// ---------------------------------------------------------------------------------------------
810// Arc
811// ---------------------------------------------------------------------------------------------
812
813/// An arc: a portion of a circle swept anticlockwise from a start angle to an end angle.
814///
815/// The angles are measured from the positive `x` axis about the circle's centre. The swept
816/// interval runs anticlockwise from `start` to `end`; when `end` precedes `start` the arc
817/// wraps through `2π`.
818#[derive(Clone, Copy, Debug)]
819pub struct Arc {
820 /// The circle the arc lies on.
821 pub circle: Circle,
822 /// The start angle, measured anticlockwise from the positive `x` axis.
823 pub start: Angle,
824 /// The end angle, measured anticlockwise from the positive `x` axis.
825 pub end: Angle,
826}
827
828impl Arc {
829 /// Creates an arc on `circle` swept anticlockwise from `start` to `end`.
830 pub fn new(circle: Circle, start: Angle, end: Angle) -> Self {
831 Self { circle, start, end }
832 }
833
834 /// Returns the point on the circle at the given angle.
835 fn point_at(&self, ang: &Angle) -> Pt {
836 self.circle.centre + Vec2::new(ang.cos(), ang.sin()) * self.circle.radius
837 }
838
839 /// Returns the point at the start of the arc.
840 pub fn start_point(&self) -> Pt {
841 self.point_at(&self.start)
842 }
843
844 /// Returns the point at the end of the arc.
845 pub fn end_point(&self) -> Pt {
846 self.point_at(&self.end)
847 }
848
849 /// Returns the anticlockwise angular sweep of the arc, as an [`Angle`] in `[0, 2π)`.
850 pub fn sweep(&self) -> Angle {
851 Angle::from_radians(self.end.radians() - self.start.radians()).normalised()
852 }
853
854 /// Returns `true` when the given angle lies within the arc's anticlockwise sweep.
855 ///
856 /// The `eps` tolerance, in radians, widens each end of the sweep so an angle sitting
857 /// exactly on an endpoint counts as inside.
858 pub fn contains_angle(&self, ang: &Angle, eps: f64) -> bool {
859 let sweep = self.sweep().radians();
860 // Offset of the query angle from the start, brought into [0, 2π).
861 let offset = Angle::from_radians(ang.radians() - self.start.radians())
862 .normalised()
863 .radians();
864 offset <= sweep + eps || offset >= 2.0 * PI - eps
865 }
866}
867
868#[cfg(test)]
869mod test {
870 use super::*;
871
872 /// The tolerance used to check oracle values in these tests.
873 const T: f64 = 1.0e-9;
874
875 // -- point-on-line incidence --------------------------------------------------------------
876
877 #[test]
878 fn test_point_on_line_00() -> Outcome<()> {
879 // The line y = x, built through the origin and (1, 1).
880 let line = res!(Line::through(Pt::new(0.0, 0.0), Pt::new(1.0, 1.0)));
881 // (2, 2) lies on it.
882 assert!(line.contains(&Pt::new(2.0, 2.0), T));
883 // (2, 3) does not.
884 assert!(!line.contains(&Pt::new(2.0, 3.0), T));
885 Ok(())
886 }
887
888 #[test]
889 fn test_point_on_line_coincident_endpoints() {
890 // A line cannot be built through two coincident points.
891 let res = Line::through(Pt::new(1.0, 1.0), Pt::new(1.0, 1.0));
892 assert!(res.is_err());
893 }
894
895 // -- line-line intersection ---------------------------------------------------------------
896
897 #[test]
898 fn test_line_line_point_00() -> Outcome<()> {
899 // y = x and y = -x + 2 meet at (1, 1).
900 let l1 = res!(Line::through(Pt::new(0.0, 0.0), Pt::new(1.0, 1.0)));
901 let l2 = res!(Line::through(Pt::new(0.0, 2.0), Pt::new(1.0, 1.0)));
902 match l1.intersect_line(&l2, T) {
903 LineIntersect::Point(p) => {
904 assert!(p.approx_eq(&Pt::new(1.0, 1.0), T));
905 },
906 other => return Err(err!("Expected a single point, got {:?}.", other; Test)),
907 }
908 Ok(())
909 }
910
911 #[test]
912 fn test_line_line_parallel() -> Outcome<()> {
913 // y = x and y = x + 1 are parallel and distinct.
914 let l1 = res!(Line::through(Pt::new(0.0, 0.0), Pt::new(1.0, 1.0)));
915 let l2 = res!(Line::through(Pt::new(0.0, 1.0), Pt::new(1.0, 2.0)));
916 assert_eq!(l1.intersect_line(&l2, T), LineIntersect::Parallel);
917 Ok(())
918 }
919
920 #[test]
921 fn test_line_line_coincident() -> Outcome<()> {
922 // y = x described two different ways is the same line.
923 let l1 = res!(Line::through(Pt::new(0.0, 0.0), Pt::new(1.0, 1.0)));
924 let l2 = res!(Line::through(Pt::new(2.0, 2.0), Pt::new(3.0, 3.0)));
925 assert_eq!(l1.intersect_line(&l2, T), LineIntersect::Coincident);
926 Ok(())
927 }
928
929 // -- line-circle intersection -------------------------------------------------------------
930
931 #[test]
932 fn test_line_circle_secant() -> Outcome<()> {
933 // The unit circle and the line x = 0 meet at (0, ±1).
934 let circle = res!(Circle::new(Pt::new(0.0, 0.0), 1.0));
935 let line = res!(Line::through(Pt::new(0.0, 0.0), Pt::new(0.0, 1.0)));
936 match line.intersect_circle(&circle, T) {
937 LineCircleIntersect::Secant(p0, p1) => {
938 // Ordered along +y: (0, -1) then (0, 1).
939 assert!(p0.approx_eq(&Pt::new(0.0, -1.0), T));
940 assert!(p1.approx_eq(&Pt::new(0.0, 1.0), T));
941 },
942 other => return Err(err!("Expected two points, got {:?}.", other; Test)),
943 }
944 Ok(())
945 }
946
947 #[test]
948 fn test_line_circle_tangent() -> Outcome<()> {
949 // The unit circle and the line y = 1 touch at (0, 1).
950 let circle = res!(Circle::new(Pt::new(0.0, 0.0), 1.0));
951 let line = res!(Line::through(Pt::new(0.0, 1.0), Pt::new(1.0, 1.0)));
952 match line.intersect_circle(&circle, T) {
953 LineCircleIntersect::Tangent(p) => {
954 assert!(p.approx_eq(&Pt::new(0.0, 1.0), T));
955 },
956 other => return Err(err!("Expected a tangent point, got {:?}.", other; Test)),
957 }
958 Ok(())
959 }
960
961 #[test]
962 fn test_line_circle_none() -> Outcome<()> {
963 // The unit circle and the line y = 2 do not meet.
964 let circle = res!(Circle::new(Pt::new(0.0, 0.0), 1.0));
965 let line = res!(Line::through(Pt::new(0.0, 2.0), Pt::new(1.0, 2.0)));
966 assert_eq!(line.intersect_circle(&circle, T), LineCircleIntersect::None);
967 Ok(())
968 }
969
970 // -- circle-circle intersection -----------------------------------------------------------
971
972 #[test]
973 fn test_circle_circle_two() -> Outcome<()> {
974 // Unit circles centred at (0, 0) and (1, 0) meet at (0.5, ±√3/2).
975 let c0 = res!(Circle::new(Pt::new(0.0, 0.0), 1.0));
976 let c1 = res!(Circle::new(Pt::new(1.0, 0.0), 1.0));
977 let root3_2 = 3.0_f64.sqrt() / 2.0;
978 match c0.intersect_circle(&c1, T) {
979 CircleIntersect::Two(p0, p1) => {
980 // The two returned points, in either order, are (0.5, ±√3/2).
981 let up = Pt::new(0.5, root3_2);
982 let dn = Pt::new(0.5, -root3_2);
983 let ok = (p0.approx_eq(&up, T) && p1.approx_eq(&dn, T))
984 || (p0.approx_eq(&dn, T) && p1.approx_eq(&up, T));
985 assert!(ok, "Got {:?} and {:?}.", p0, p1);
986 },
987 other => return Err(err!("Expected two points, got {:?}.", other; Test)),
988 }
989 Ok(())
990 }
991
992 #[test]
993 fn test_circle_circle_tangent() -> Outcome<()> {
994 // Unit circles at (0, 0) and (2, 0) touch externally at (1, 0).
995 let c0 = res!(Circle::new(Pt::new(0.0, 0.0), 1.0));
996 let c1 = res!(Circle::new(Pt::new(2.0, 0.0), 1.0));
997 match c0.intersect_circle(&c1, T) {
998 CircleIntersect::Tangent(p) => {
999 assert!(p.approx_eq(&Pt::new(1.0, 0.0), T));
1000 },
1001 other => return Err(err!("Expected a tangent point, got {:?}.", other; Test)),
1002 }
1003 Ok(())
1004 }
1005
1006 #[test]
1007 fn test_circle_circle_none() -> Outcome<()> {
1008 // Unit circles at (0, 0) and (5, 0) are too far apart to meet.
1009 let c0 = res!(Circle::new(Pt::new(0.0, 0.0), 1.0));
1010 let c1 = res!(Circle::new(Pt::new(5.0, 0.0), 1.0));
1011 assert_eq!(c0.intersect_circle(&c1, T), CircleIntersect::None);
1012 Ok(())
1013 }
1014
1015 #[test]
1016 fn test_circle_circle_coincident() -> Outcome<()> {
1017 // The same circle described twice meets everywhere.
1018 let c0 = res!(Circle::new(Pt::new(0.0, 0.0), 1.0));
1019 let c1 = res!(Circle::new(Pt::new(0.0, 0.0), 1.0));
1020 assert_eq!(c0.intersect_circle(&c1, T), CircleIntersect::Coincident);
1021 Ok(())
1022 }
1023
1024 // -- foot of perpendicular ----------------------------------------------------------------
1025
1026 #[test]
1027 fn test_foot_of_perpendicular_00() -> Outcome<()> {
1028 // The foot of the perpendicular from (0, 2) to y = 0 is (0, 0).
1029 let line = res!(Line::through(Pt::new(0.0, 0.0), Pt::new(1.0, 0.0)));
1030 let foot = line.foot_of_perpendicular(&Pt::new(0.0, 2.0));
1031 assert!(foot.approx_eq(&Pt::new(0.0, 0.0), T));
1032 Ok(())
1033 }
1034
1035 #[test]
1036 fn test_foot_of_perpendicular_diagonal() -> Outcome<()> {
1037 // The foot of the perpendicular from (0, 2) onto y = x is (1, 1).
1038 let line = res!(Line::through(Pt::new(0.0, 0.0), Pt::new(1.0, 1.0)));
1039 let foot = line.foot_of_perpendicular(&Pt::new(0.0, 2.0));
1040 assert!(foot.approx_eq(&Pt::new(1.0, 1.0), T));
1041 Ok(())
1042 }
1043
1044 // -- projection onto a line ---------------------------------------------------------------
1045
1046 #[test]
1047 fn test_project_point_onto_line() -> Outcome<()> {
1048 // Projecting (2, 5) onto y = 0 gives (2, 0).
1049 let line = res!(Line::through(Pt::new(0.0, 0.0), Pt::new(1.0, 0.0)));
1050 let p = line.project_point(&Pt::new(2.0, 5.0));
1051 assert!(p.approx_eq(&Pt::new(2.0, 0.0), T));
1052 Ok(())
1053 }
1054
1055 // -- projection onto a segment (clamped) --------------------------------------------------
1056
1057 #[test]
1058 fn test_project_point_onto_segment_interior() -> Outcome<()> {
1059 // Projecting (0.5, 2) onto the segment (0,0)-(1,0) lands at (0.5, 0).
1060 let seg = res!(Segment::new(Pt::new(0.0, 0.0), Pt::new(1.0, 0.0)));
1061 let p = seg.project_point(&Pt::new(0.5, 2.0));
1062 assert!(p.approx_eq(&Pt::new(0.5, 0.0), T));
1063 Ok(())
1064 }
1065
1066 #[test]
1067 fn test_project_point_onto_segment_clamped_ends() -> Outcome<()> {
1068 // Beyond an end, the projection clamps to that endpoint.
1069 let seg = res!(Segment::new(Pt::new(0.0, 0.0), Pt::new(1.0, 0.0)));
1070 let past_b = seg.project_point(&Pt::new(2.0, 2.0));
1071 assert!(past_b.approx_eq(&Pt::new(1.0, 0.0), T));
1072 let past_a = seg.project_point(&Pt::new(-1.0, 2.0));
1073 assert!(past_a.approx_eq(&Pt::new(0.0, 0.0), T));
1074 Ok(())
1075 }
1076
1077 // -- angle between two rays ---------------------------------------------------------------
1078
1079 #[test]
1080 fn test_angle_between_rays_right_angle() -> Outcome<()> {
1081 // The +x and +y rays are a quarter turn apart: π/2.
1082 let rx = res!(Ray::new(Pt::new(0.0, 0.0), Vec2::new(1.0, 0.0)));
1083 let ry = res!(Ray::new(Pt::new(0.0, 0.0), Vec2::new(0.0, 1.0)));
1084 let ang = rx.angle_between(&ry);
1085 assert!(approx_eq(ang.radians(), PI / 2.0, T), "Got {} rad.", ang.radians());
1086 Ok(())
1087 }
1088
1089 #[test]
1090 fn test_angle_between_rays_opposite() -> Outcome<()> {
1091 // The +x and -x rays are a half turn apart: π (the degenerate straight-angle case).
1092 let rx = res!(Ray::new(Pt::new(0.0, 0.0), Vec2::new(1.0, 0.0)));
1093 let rev = res!(Ray::new(Pt::new(0.0, 0.0), Vec2::new(-1.0, 0.0)));
1094 let ang = rx.angle_between(&rev);
1095 assert!(approx_eq(ang.radians(), PI, T), "Got {} rad.", ang.radians());
1096 Ok(())
1097 }
1098
1099 // -- supporting types ---------------------------------------------------------------------
1100
1101 #[test]
1102 fn test_vec2_normalise_zero_fails() {
1103 // A zero vector has no direction and cannot be normalised.
1104 let res = Vec2::new(0.0, 0.0).normalise();
1105 assert!(res.is_err());
1106 }
1107
1108 #[test]
1109 fn test_circle_zero_radius_fails() {
1110 // A circle needs a positive radius.
1111 let res = Circle::new(Pt::new(0.0, 0.0), 0.0);
1112 assert!(res.is_err());
1113 }
1114
1115 #[test]
1116 fn test_angle_degrees_radians() {
1117 // 180 degrees is π radians, and the conversion round-trips.
1118 let a = Angle::from_degrees(180.0);
1119 assert!(approx_eq(a.radians(), PI, T));
1120 assert!(approx_eq(a.degrees(), 180.0, T));
1121 }
1122
1123 // -- triangle: area, barycentric coordinates and the affine map ----------------------------
1124
1125 #[test]
1126 fn test_triangle_area_00() -> Outcome<()> {
1127 // The 3-4-5 right triangle: area 6 by the half-base-times-height oracle.
1128 let t = res!(Triangle::new(Pt::new(0.0, 0.0), Pt::new(4.0, 0.0), Pt::new(0.0, 3.0)));
1129 assert!(approx_eq(t.area(), 6.0, T));
1130 // Anticlockwise in a y-up frame, so the signed area is positive.
1131 assert!(t.signed_area() > 0.0);
1132 // Reversing two corners reverses the winding and nothing else.
1133 let r = res!(Triangle::new(Pt::new(0.0, 0.0), Pt::new(0.0, 3.0), Pt::new(4.0, 0.0)));
1134 assert!(approx_eq(r.signed_area(), -6.0, T));
1135 assert!(approx_eq(r.area(), 6.0, T));
1136 Ok(())
1137 }
1138
1139 #[test]
1140 fn test_triangle_collinear_fails() {
1141 // Three points on the line y = 2x have no interior.
1142 let res = Triangle::new(Pt::new(0.0, 0.0), Pt::new(1.0, 2.0), Pt::new(3.0, 6.0));
1143 assert!(res.is_err());
1144 }
1145
1146 #[test]
1147 fn test_triangle_barycentric_corners_and_centroid() -> Outcome<()> {
1148 let (a, b, c) = (Pt::new(0.0, 0.0), Pt::new(4.0, 0.0), Pt::new(0.0, 3.0));
1149 let t = res!(Triangle::new(a, b, c));
1150 // Each corner carries all of its own weight.
1151 for (p, want) in [(a, (1.0, 0.0, 0.0)), (b, (0.0, 1.0, 0.0)), (c, (0.0, 0.0, 1.0))] {
1152 let (u, v, w) = res!(t.barycentric(&p));
1153 assert!(approx_eq(u, want.0, T) && approx_eq(v, want.1, T) && approx_eq(w, want.2, T),
1154 "corner {:?} gave ({}, {}, {})", p, u, v, w);
1155 }
1156 // The centroid is a third of each.
1157 let (u, v, w) = res!(t.barycentric(&t.centroid()));
1158 assert!(approx_eq(u, 1.0 / 3.0, T) && approx_eq(v, 1.0 / 3.0, T) && approx_eq(w, 1.0 / 3.0, T));
1159 // A midpoint of an edge splits its two corners and excludes the third.
1160 let (u, v, w) = res!(t.barycentric(&a.midpoint(&b)));
1161 assert!(approx_eq(u, 0.5, T) && approx_eq(v, 0.5, T) && approx_eq(w, 0.0, T));
1162 Ok(())
1163 }
1164
1165 #[test]
1166 fn test_triangle_contains_00() -> Outcome<()> {
1167 let t = res!(Triangle::new(Pt::new(0.0, 0.0), Pt::new(4.0, 0.0), Pt::new(0.0, 3.0)));
1168 // Inside, on an edge, on a corner, and outside past the hypotenuse.
1169 assert!(t.contains(&Pt::new(1.0, 1.0), T));
1170 assert!(t.contains(&Pt::new(2.0, 0.0), T));
1171 assert!(t.contains(&Pt::new(0.0, 3.0), T));
1172 assert!(!t.contains(&Pt::new(3.0, 3.0), T));
1173 assert!(!t.contains(&Pt::new(-0.5, 1.0), T));
1174 Ok(())
1175 }
1176
1177 #[test]
1178 fn test_triangle_barycentric_round_trips() -> Outcome<()> {
1179 // Reading a point's weights and rebuilding it returns the same point.
1180 let t = res!(Triangle::new(Pt::new(-2.0, 1.0), Pt::new(5.0, -3.0), Pt::new(1.0, 6.0)));
1181 for p in [Pt::new(1.0, 1.0), Pt::new(0.0, 0.0), Pt::new(4.0, 2.0)] {
1182 let (u, v, w) = res!(t.barycentric(&p));
1183 assert!(approx_eq(u + v + w, 1.0, T), "weights must sum to one");
1184 let q = t.point_at(u, v, w);
1185 assert!(p.approx_eq(&q, 1.0e-9), "{:?} rebuilt as {:?}", p, q);
1186 }
1187 Ok(())
1188 }
1189
1190 #[test]
1191 fn test_triangle_map_point_is_affine() -> Outcome<()> {
1192 // Map the unit right triangle onto one scaled by two in x, three in y and shifted.
1193 // The image of (0.25, 0.5) is therefore (10 + 0.5, 20 + 1.5) by hand.
1194 let src = res!(Triangle::new(Pt::new(0.0, 0.0), Pt::new(1.0, 0.0), Pt::new(0.0, 1.0)));
1195 let dst = res!(Triangle::new(Pt::new(10.0, 20.0), Pt::new(12.0, 20.0), Pt::new(10.0, 23.0)));
1196 let q = res!(src.map_point(&dst, &Pt::new(0.25, 0.5)));
1197 assert!(q.approx_eq(&Pt::new(10.5, 21.5), 1.0e-9), "mapped to {:?}", q);
1198 // Being affine, the image of a midpoint is the midpoint of the images.
1199 let (p0, p1) = (Pt::new(0.1, 0.1), Pt::new(0.6, 0.2));
1200 let (i0, i1) = (res!(src.map_point(&dst, &p0)), res!(src.map_point(&dst, &p1)));
1201 let mid = res!(src.map_point(&dst, &p0.midpoint(&p1)));
1202 assert!(mid.approx_eq(&i0.midpoint(&i1), 1.0e-9));
1203 Ok(())
1204 }
1205
1206 #[test]
1207 fn test_segment_cross_line_reports_the_fraction_along_the_segment() -> Outcome<()> {
1208 // A vertical line at x = 3 cuts the segment from (1,0) to (5,0) at a quarter of
1209 // the way along, by hand: (3 - 1) / (5 - 1).
1210 let seg = res!(Segment::new(Pt::new(1.0, 0.0), Pt::new(5.0, 0.0)));
1211 let line = res!(Line::new(Pt::new(3.0, -7.0), Vec2::new(0.0, 1.0)));
1212 let t = res!(seg.cross_line(&line, 1.0e-9)
1213 .ok_or_else(|| err!("the line crosses the segment"; Test)));
1214 assert!(approx_eq(t, 0.5, 1.0e-9), "crossed at {}", t);
1215 assert!(seg.point_at(t).approx_eq(&Pt::new(3.0, 0.0), 1.0e-9));
1216
1217 // A line beyond the far end misses it, and one along it reports nothing to name.
1218 let past = res!(Line::new(Pt::new(9.0, 0.0), Vec2::new(0.0, 1.0)));
1219 assert_eq!(seg.cross_line(&past, 1.0e-9), None);
1220 let along = res!(Line::new(Pt::new(0.0, 0.0), Vec2::new(1.0, 0.0)));
1221 assert_eq!(seg.cross_line(&along, 1.0e-9), None);
1222
1223 // A crossing exactly at an endpoint counts: dropping it would open a gap in a sweep.
1224 let at_end = res!(Line::new(Pt::new(5.0, 0.0), Vec2::new(0.0, 1.0)));
1225 assert_eq!(seg.cross_line(&at_end, 1.0e-9), Some(1.0));
1226 Ok(())
1227 }
1228
1229 #[test]
1230 fn test_segment_intersect_segment_only_within_both() -> Outcome<()> {
1231 // The diagonals of the unit square meet at its centre.
1232 let d1 = res!(Segment::new(Pt::new(0.0, 0.0), Pt::new(1.0, 1.0)));
1233 let d2 = res!(Segment::new(Pt::new(0.0, 1.0), Pt::new(1.0, 0.0)));
1234 let p = res!(d1.intersect_segment(&d2, 1.0e-9)
1235 .ok_or_else(|| err!("the diagonals meet"; Test)));
1236 assert!(p.approx_eq(&Pt::new(0.5, 0.5), 1.0e-9), "met at {:?}", p);
1237
1238 // Crossing lines whose segments stop short of each other do not meet.
1239 let short = res!(Segment::new(Pt::new(0.0, 1.0), Pt::new(0.2, 0.8)));
1240 assert_eq!(d1.intersect_segment(&short, 1.0e-9), None);
1241
1242 // Parallel and collinear pairs both report nothing.
1243 let par = res!(Segment::new(Pt::new(0.0, 1.0), Pt::new(1.0, 2.0)));
1244 assert_eq!(d1.intersect_segment(&par, 1.0e-9), None);
1245 let over = res!(Segment::new(Pt::new(0.5, 0.5), Pt::new(2.0, 2.0)));
1246 assert_eq!(d1.intersect_segment(&over, 1.0e-9), None);
1247 Ok(())
1248 }
1249}