Skip to main content

axiolid_inspect/
sheet.rs

1//! Two faces of one mesh that face each other within a linear distance.
2//!
3//! A boolean leaves such a thin sheet where an operand's face lands a
4//! rounding residue off the face it should have met: the result is a Solid
5//! whose lid is thinner than its coordinates can carry (#400). The query is
6//! the mesh's own, so it needs no second operand.
7
8use axiolid_core::{Point3, Vec3};
9use axiolid_mesh::TriMesh;
10use axiolid_spatial::SpatialItem;
11
12use crate::clearance::{build_index, triangle_at, triangle_distance};
13use crate::coplanar::corner_distances;
14
15/// A pair of triangles of one mesh that face each other within a linear
16/// distance.
17#[non_exhaustive]
18#[derive(Debug, Clone, Copy, PartialEq)]
19pub struct ThinSheet {
20    /// The lower triangle index.
21    pub first: usize,
22    /// The higher triangle index.
23    pub second: usize,
24    /// The largest distance from a corner of either triangle to the other's
25    /// plane.
26    pub offset: f64,
27}
28
29/// The first pair of triangles of `mesh` that form a sheet thinner than
30/// `linear`, in ascending `(first, second)` order with `first < second`.
31///
32/// A pair qualifies when both triangles are non-degenerate, their unit
33/// normals point opposite ways, the distance between them is at most
34/// `linear`, every corner of each lies within `linear` of the other's plane,
35/// and, projected onto the first's plane, they overlap in a polygon whose
36/// doubled area exceeds `linear` times its perimeter: contact along a single
37/// edge is not a sheet. Exactly coplanar pairs count. `None` when no pair
38/// qualifies, and when `linear` is negative or not finite.
39///
40/// # Panics
41///
42/// Never. A malformed index buffer yields `None` rather than a panic.
43pub fn thin_sheet(mesh: &TriMesh, linear: f64) -> Option<ThinSheet> {
44    if !linear.is_finite() || linear < 0.0 {
45        return None;
46    }
47    let index = build_index(mesh, linear);
48    let mut candidates = Vec::new();
49    for first in 0..mesh.indices.len() / 3 {
50        let Some(first_triangle) = triangle_at(mesh, first) else {
51            continue;
52        };
53        let probe = crate::clearance::grow(crate::clearance::bounds_of(&first_triangle), linear);
54        candidates.clear();
55        index.query_aabb(&probe, &mut candidates);
56        let mut seconds: Vec<usize> = candidates
57            .iter()
58            .filter_map(|&candidate| {
59                index
60                    .item(candidate)
61                    .map(|item: &SpatialItem<u32>| item.key as usize)
62            })
63            .filter(|&second| second > first)
64            .collect();
65        seconds.sort_unstable();
66        for second in seconds {
67            let Some(second_triangle) = triangle_at(mesh, second) else {
68                continue;
69            };
70            if let Some(offset) = facing_within(&first_triangle, &second_triangle, linear) {
71                return Some(ThinSheet {
72                    first,
73                    second,
74                    offset,
75                });
76            }
77        }
78    }
79    None
80}
81
82/// The pair's offset when it qualifies.
83fn facing_within(a: &[Point3; 3], b: &[Point3; 3], linear: f64) -> Option<f64> {
84    let to_a = corner_distances(a, b)?;
85    let to_b = corner_distances(b, a)?;
86    let offset = to_a.max(to_b);
87    if offset > linear {
88        return None;
89    }
90    let na = (a[1] - a[0]).cross(a[2] - a[0]).normalize();
91    let nb = (b[1] - b[0]).cross(b[2] - b[0]).normalize();
92    if na.dot(nb) >= 0.0 || triangle_distance(a, b) > linear {
93        return None;
94    }
95    let (area2, perimeter) = overlap(a, b, na);
96    (area2 > linear * perimeter).then_some(offset)
97}
98
99/// Twice the area and the perimeter of the polygon in which `b`, projected
100/// onto the plane of `a` (unit normal `n`), overlaps `a`.
101fn overlap(a: &[Point3; 3], b: &[Point3; 3], n: Vec3) -> (f64, f64) {
102    let u = (a[1] - a[0]).normalize();
103    let v = n.cross(u);
104    let flat = |p: Point3| [(p - a[0]).dot(u), (p - a[0]).dot(v)];
105    let clip = a.map(flat);
106    let mut polygon: Vec<[f64; 2]> = b.iter().map(|&p| flat(p)).collect();
107    // `b` faces away from `a`: its projection winds the other way.
108    if signed_area2(&polygon) < 0.0 {
109        polygon.reverse();
110    }
111    for i in 0..3 {
112        let (p, q) = (clip[i], clip[(i + 1) % 3]);
113        let side = |x: [f64; 2]| (q[0] - p[0]) * (x[1] - p[1]) - (q[1] - p[1]) * (x[0] - p[0]);
114        let mut kept = Vec::with_capacity(polygon.len() + 1);
115        for j in 0..polygon.len() {
116            let (s, e) = (polygon[j], polygon[(j + 1) % polygon.len()]);
117            let (ds, de) = (side(s), side(e));
118            if ds >= 0.0 {
119                kept.push(s);
120            }
121            if (ds > 0.0 && de < 0.0) || (ds < 0.0 && de > 0.0) {
122                let t = ds / (ds - de);
123                kept.push([s[0] + t * (e[0] - s[0]), s[1] + t * (e[1] - s[1])]);
124            }
125        }
126        polygon = kept;
127        if polygon.len() < 3 {
128            return (0.0, 0.0);
129        }
130    }
131    let perimeter = (0..polygon.len())
132        .map(|i| {
133            let (s, e) = (polygon[i], polygon[(i + 1) % polygon.len()]);
134            (e[0] - s[0]).hypot(e[1] - s[1])
135        })
136        .sum();
137    (signed_area2(&polygon).abs(), perimeter)
138}
139
140fn signed_area2(polygon: &[[f64; 2]]) -> f64 {
141    (0..polygon.len())
142        .map(|i| {
143            let (s, e) = (polygon[i], polygon[(i + 1) % polygon.len()]);
144            s[0] * e[1] - e[0] * s[1]
145        })
146        .sum()
147}