Skip to main content

ifc_lite_geometry/
void_analysis.rs

1// This Source Code Form is subject to the terms of the Mozilla Public
2// License, v. 2.0. If a copy of the MPL was not distributed with this
3// file, You can obtain one at https://mozilla.org/MPL/2.0/.
4
5//! Void Analysis Module
6//!
7//! Analyzes void geometry to determine optimal processing strategy:
8//! - Coplanar voids can be subtracted in 2D at the profile level
9//! - Non-planar voids require full 3D CSG operations
10//!
11//! This module implements the detection and projection logic needed to
12//! efficiently route voids to the appropriate processing path.
13
14use crate::bool2d::{ensure_ccw, is_valid_contour};
15use crate::mesh::Mesh;
16use crate::profile::VoidInfo;
17use nalgebra::{Matrix4, Point2, Point3, Vector3};
18use rustc_hash::FxHashMap;
19
20/// Default epsilon for planarity detection
21const DEFAULT_PLANARITY_EPSILON: f64 = 0.02;
22
23/// Maximum depth tolerance for considering a void as "through"
24const THROUGH_VOID_TOLERANCE: f64 = 0.01;
25
26/// Classification of a void relative to a host extrusion
27#[derive(Debug, Clone)]
28pub enum VoidClassification {
29    /// Void is coplanar with the profile plane and can be subtracted in 2D
30    Coplanar {
31        /// The 2D footprint of the void in profile space
32        profile_hole: Vec<Point2<f64>>,
33        /// Start depth in extrusion space (0.0 = bottom cap)
34        depth_start: f64,
35        /// End depth in extrusion space
36        depth_end: f64,
37        /// Whether the void extends through the full depth
38        is_through: bool,
39    },
40    /// Void is at an angle to the extrusion direction - requires 3D CSG
41    NonPlanar {
42        /// The void mesh for 3D CSG subtraction
43        mesh: Mesh,
44    },
45    /// Void doesn't intersect the host geometry - can be skipped
46    NonIntersecting,
47}
48
49/// Analyzer for void geometry classification
50pub struct VoidAnalyzer {
51    /// Epsilon for planarity detection (dot product threshold)
52    planarity_epsilon: f64,
53    /// Whether to use adaptive epsilon refinement
54    adaptive_epsilon: bool,
55}
56
57impl VoidAnalyzer {
58    /// Create a new void analyzer with default settings
59    pub fn new() -> Self {
60        Self {
61            planarity_epsilon: DEFAULT_PLANARITY_EPSILON,
62            adaptive_epsilon: true,
63        }
64    }
65
66    /// Create a void analyzer with custom planarity epsilon
67    pub fn with_epsilon(epsilon: f64) -> Self {
68        Self {
69            planarity_epsilon: epsilon,
70            adaptive_epsilon: false,
71        }
72    }
73
74    /// Classify a void relative to host extrusion parameters
75    ///
76    /// # Arguments
77    /// * `void_mesh` - The mesh geometry of the void/opening
78    /// * `profile_transform` - Transform from profile space to world space
79    /// * `extrusion_direction` - Direction of extrusion (normalized)
80    /// * `extrusion_depth` - Total depth of extrusion
81    ///
82    /// # Returns
83    /// Classification indicating how the void should be processed
84    pub fn classify_void(
85        &self,
86        void_mesh: &Mesh,
87        profile_transform: &Matrix4<f64>,
88        extrusion_direction: &Vector3<f64>,
89        extrusion_depth: f64,
90    ) -> VoidClassification {
91        if void_mesh.is_empty() {
92            return VoidClassification::NonIntersecting;
93        }
94
95        // Get profile plane normal (Z-axis in profile space transformed to world)
96        let profile_normal = transform_direction(profile_transform, &Vector3::new(0.0, 0.0, 1.0));
97
98        // Check if void is coplanar with profile plane
99        let is_coplanar = self.check_coplanarity(void_mesh, &profile_normal, extrusion_direction);
100
101        if !is_coplanar {
102            return VoidClassification::NonPlanar {
103                mesh: void_mesh.clone(),
104            };
105        }
106
107        // Extract 2D footprint
108        let inverse_transform = match profile_transform.try_inverse() {
109            Some(inv) => inv,
110            None => {
111                return VoidClassification::NonPlanar {
112                    mesh: void_mesh.clone(),
113                }
114            }
115        };
116
117        match self.extract_footprint(void_mesh, &inverse_transform, extrusion_direction) {
118            Some((footprint, depth_start, depth_end)) => {
119                // Check if footprint is valid
120                if !is_valid_contour(&footprint) {
121                    return VoidClassification::NonPlanar {
122                        mesh: void_mesh.clone(),
123                    };
124                }
125
126                // Check if it's a through void
127                let is_through = depth_start <= THROUGH_VOID_TOLERANCE
128                    && depth_end >= extrusion_depth - THROUGH_VOID_TOLERANCE;
129
130                VoidClassification::Coplanar {
131                    profile_hole: footprint,
132                    depth_start,
133                    depth_end,
134                    is_through,
135                }
136            }
137            None => VoidClassification::NonPlanar {
138                mesh: void_mesh.clone(),
139            },
140        }
141    }
142
143    /// Check if void geometry is coplanar with the profile plane
144    fn check_coplanarity(
145        &self,
146        void_mesh: &Mesh,
147        profile_normal: &Vector3<f64>,
148        extrusion_direction: &Vector3<f64>,
149    ) -> bool {
150        // Get dominant faces of the void mesh
151        let face_normals = self.extract_face_normals(void_mesh);
152
153        if face_normals.is_empty() {
154            return false;
155        }
156
157        // A void is coplanar if its dominant face normals are either:
158        // 1. Parallel to the profile normal (top/bottom faces of void)
159        // 2. Perpendicular to the extrusion direction (side faces)
160
161        let mut has_parallel_face = false;
162        let mut epsilon = self.planarity_epsilon;
163
164        // Adaptive epsilon: try progressively smaller values
165        let epsilons = if self.adaptive_epsilon {
166            vec![0.02, 0.01, 0.005, 0.001]
167        } else {
168            vec![epsilon]
169        };
170
171        for eps in epsilons {
172            epsilon = eps;
173            has_parallel_face = false;
174
175            for normal in &face_normals {
176                let dot_profile = normal.dot(profile_normal).abs();
177                let dot_extrusion = normal.dot(extrusion_direction).abs();
178
179                // Face is parallel to profile plane (top/bottom of void)
180                if dot_profile > 1.0 - epsilon {
181                    has_parallel_face = true;
182                    break;
183                }
184
185                // Face is perpendicular to extrusion (side of void, parallel to extrusion)
186                if dot_extrusion < epsilon {
187                    has_parallel_face = true;
188                    break;
189                }
190            }
191
192            if has_parallel_face {
193                break;
194            }
195        }
196
197        has_parallel_face
198    }
199
200    /// Extract face normals from mesh, grouped by direction
201    fn extract_face_normals(&self, mesh: &Mesh) -> Vec<Vector3<f64>> {
202        let mut normal_groups: FxHashMap<(i32, i32, i32), (Vector3<f64>, usize)> =
203            FxHashMap::default();
204
205        let quantize_normal = |n: &Vector3<f64>| -> (i32, i32, i32) {
206            (
207                (n.x * 100.0).round() as i32,
208                (n.y * 100.0).round() as i32,
209                (n.z * 100.0).round() as i32,
210            )
211        };
212
213        for i in (0..mesh.indices.len()).step_by(3) {
214            let i0 = mesh.indices[i] as usize;
215            let i1 = mesh.indices[i + 1] as usize;
216            let i2 = mesh.indices[i + 2] as usize;
217
218            if i0 * 3 + 2 >= mesh.positions.len()
219                || i1 * 3 + 2 >= mesh.positions.len()
220                || i2 * 3 + 2 >= mesh.positions.len()
221            {
222                continue;
223            }
224
225            let v0 = Point3::new(
226                mesh.positions[i0 * 3] as f64,
227                mesh.positions[i0 * 3 + 1] as f64,
228                mesh.positions[i0 * 3 + 2] as f64,
229            );
230            let v1 = Point3::new(
231                mesh.positions[i1 * 3] as f64,
232                mesh.positions[i1 * 3 + 1] as f64,
233                mesh.positions[i1 * 3 + 2] as f64,
234            );
235            let v2 = Point3::new(
236                mesh.positions[i2 * 3] as f64,
237                mesh.positions[i2 * 3 + 1] as f64,
238                mesh.positions[i2 * 3 + 2] as f64,
239            );
240
241            let edge1 = v1 - v0;
242            let edge2 = v2 - v0;
243            let normal = match edge1.cross(&edge2).try_normalize(1e-10) {
244                Some(n) => n,
245                None => continue,
246            };
247
248            let key = quantize_normal(&normal);
249            normal_groups
250                .entry(key)
251                .and_modify(|(sum, count)| {
252                    *sum += normal;
253                    *count += 1;
254                })
255                .or_insert((normal, 1));
256        }
257
258        // Return normalized group averages
259        normal_groups
260            .values()
261            .filter_map(|(sum, count)| {
262                if *count > 0 {
263                    (sum / *count as f64).try_normalize(1e-10)
264                } else {
265                    None
266                }
267            })
268            .collect()
269    }
270
271    /// Extract 2D footprint from void mesh
272    ///
273    /// Projects the void mesh onto the profile plane and extracts the
274    /// boundary contour.
275    fn extract_footprint(
276        &self,
277        void_mesh: &Mesh,
278        inverse_profile_transform: &Matrix4<f64>,
279        _extrusion_direction: &Vector3<f64>,
280    ) -> Option<(Vec<Point2<f64>>, f64, f64)> {
281        if void_mesh.is_empty() {
282            return None;
283        }
284
285        // Transform all vertices to profile space
286        let mut min_z = f64::MAX;
287        let mut max_z = f64::MIN;
288        let mut projected_points: Vec<Point2<f64>> = Vec::new();
289
290        for i in (0..void_mesh.positions.len()).step_by(3) {
291            let world_point = Point3::new(
292                void_mesh.positions[i] as f64,
293                void_mesh.positions[i + 1] as f64,
294                void_mesh.positions[i + 2] as f64,
295            );
296
297            // Transform to profile space
298            let profile_point = inverse_profile_transform.transform_point(&world_point);
299
300            // Track Z range for depth
301            min_z = min_z.min(profile_point.z);
302            max_z = max_z.max(profile_point.z);
303
304            // Project to 2D (XY plane in profile space)
305            projected_points.push(Point2::new(profile_point.x, profile_point.y));
306        }
307
308        if projected_points.len() < 3 {
309            return None;
310        }
311
312        // Extract convex hull or boundary of projected points
313        let hull = self.compute_convex_hull(&projected_points);
314
315        if hull.len() < 3 {
316            return None;
317        }
318
319        // Ensure counter-clockwise winding
320        let footprint = ensure_ccw(&hull);
321
322        // Clamp depth values
323        let depth_start = min_z.max(0.0);
324        let depth_end = max_z;
325
326        Some((footprint, depth_start, depth_end))
327    }
328
329    /// Compute convex hull of 2D points using Graham scan
330    fn compute_convex_hull(&self, points: &[Point2<f64>]) -> Vec<Point2<f64>> {
331        if points.len() < 3 {
332            return points.to_vec();
333        }
334
335        // Find bottom-most point (lowest Y, then leftmost X)
336        let mut start_idx = 0;
337        for (i, p) in points.iter().enumerate() {
338            if p.y < points[start_idx].y
339                || (p.y == points[start_idx].y && p.x < points[start_idx].x)
340            {
341                start_idx = i;
342            }
343        }
344
345        let start = points[start_idx];
346
347        // Sort points by polar angle with respect to start
348        let mut sorted: Vec<Point2<f64>> = points
349            .iter()
350            .filter(|p| **p != start)
351            .cloned()
352            .collect();
353
354        sorted.sort_by(|a, b| {
355            let angle_a = (a.y - start.y).atan2(a.x - start.x);
356            let angle_b = (b.y - start.y).atan2(b.x - start.x);
357            angle_a.total_cmp(&angle_b)
358        });
359
360        // Graham scan
361        let mut hull = vec![start];
362
363        for p in sorted {
364            while hull.len() > 1 {
365                let top = hull[hull.len() - 1];
366                let second = hull[hull.len() - 2];
367
368                // Cross product to check turn direction
369                let cross = (top.x - second.x) * (p.y - second.y)
370                    - (top.y - second.y) * (p.x - second.x);
371
372                if cross <= 0.0 {
373                    hull.pop();
374                } else {
375                    break;
376                }
377            }
378            hull.push(p);
379        }
380
381        hull
382    }
383
384    /// Compute depth range of void in extrusion space
385    pub fn compute_depth_range(
386        &self,
387        void_mesh: &Mesh,
388        profile_origin: &Point3<f64>,
389        extrusion_direction: &Vector3<f64>,
390    ) -> (f64, f64) {
391        if void_mesh.is_empty() {
392            return (0.0, 0.0);
393        }
394
395        let mut min_depth = f64::MAX;
396        let mut max_depth = f64::MIN;
397
398        for i in (0..void_mesh.positions.len()).step_by(3) {
399            let point = Point3::new(
400                void_mesh.positions[i] as f64,
401                void_mesh.positions[i + 1] as f64,
402                void_mesh.positions[i + 2] as f64,
403            );
404
405            // Project onto extrusion direction
406            let relative = point - profile_origin;
407            let depth = relative.dot(extrusion_direction);
408
409            min_depth = min_depth.min(depth);
410            max_depth = max_depth.max(depth);
411        }
412
413        (min_depth.max(0.0), max_depth)
414    }
415}
416
417impl Default for VoidAnalyzer {
418    fn default() -> Self {
419        Self::new()
420    }
421}
422
423/// Transform a direction vector (ignoring translation)
424fn transform_direction(transform: &Matrix4<f64>, direction: &Vector3<f64>) -> Vector3<f64> {
425    let transformed = transform.transform_vector(direction);
426    transformed.normalize()
427}
428
429/// Batch classify multiple voids for a single host
430pub fn classify_voids_batch(
431    void_meshes: &[Mesh],
432    profile_transform: &Matrix4<f64>,
433    extrusion_direction: &Vector3<f64>,
434    extrusion_depth: f64,
435) -> Vec<VoidClassification> {
436    let analyzer = VoidAnalyzer::new();
437
438    void_meshes
439        .iter()
440        .map(|mesh| {
441            analyzer.classify_void(mesh, profile_transform, extrusion_direction, extrusion_depth)
442        })
443        .collect()
444}
445
446/// Extract coplanar voids from a batch classification result
447pub fn extract_coplanar_voids(classifications: &[VoidClassification]) -> Vec<VoidInfo> {
448    classifications
449        .iter()
450        .filter_map(|c| match c {
451            VoidClassification::Coplanar {
452                profile_hole,
453                depth_start,
454                depth_end,
455                is_through,
456            } => Some(VoidInfo {
457                contour: profile_hole.clone(),
458                depth_start: *depth_start,
459                depth_end: *depth_end,
460                is_through: *is_through,
461            }),
462            _ => None,
463        })
464        .collect()
465}
466
467/// Extract non-planar void meshes from a batch classification result
468pub fn extract_nonplanar_voids(classifications: Vec<VoidClassification>) -> Vec<Mesh> {
469    classifications
470        .into_iter()
471        .filter_map(|c| match c {
472            VoidClassification::NonPlanar { mesh } => Some(mesh),
473            _ => None,
474        })
475        .collect()
476}
477
478#[cfg(test)]
479mod tests {
480    use super::*;
481
482    fn create_box_mesh(min: Point3<f64>, max: Point3<f64>) -> Mesh {
483        let mut mesh = Mesh::with_capacity(8, 36);
484
485        // 8 vertices of a box
486        let vertices = [
487            Point3::new(min.x, min.y, min.z),
488            Point3::new(max.x, min.y, min.z),
489            Point3::new(max.x, max.y, min.z),
490            Point3::new(min.x, max.y, min.z),
491            Point3::new(min.x, min.y, max.z),
492            Point3::new(max.x, min.y, max.z),
493            Point3::new(max.x, max.y, max.z),
494            Point3::new(min.x, max.y, max.z),
495        ];
496
497        // Add all vertices with dummy normals
498        for v in &vertices {
499            mesh.add_vertex(*v, Vector3::new(0.0, 0.0, 1.0));
500        }
501
502        // 12 triangles (2 per face)
503        // Front face
504        mesh.add_triangle(0, 1, 2);
505        mesh.add_triangle(0, 2, 3);
506        // Back face
507        mesh.add_triangle(4, 6, 5);
508        mesh.add_triangle(4, 7, 6);
509        // Left face
510        mesh.add_triangle(0, 3, 7);
511        mesh.add_triangle(0, 7, 4);
512        // Right face
513        mesh.add_triangle(1, 5, 6);
514        mesh.add_triangle(1, 6, 2);
515        // Bottom face
516        mesh.add_triangle(0, 4, 5);
517        mesh.add_triangle(0, 5, 1);
518        // Top face
519        mesh.add_triangle(3, 2, 6);
520        mesh.add_triangle(3, 6, 7);
521
522        mesh
523    }
524
525    #[test]
526    fn test_void_analyzer_coplanar() {
527        let analyzer = VoidAnalyzer::new();
528
529        // Create a vertical box void (aligned with Z-axis extrusion)
530        let void_mesh = create_box_mesh(
531            Point3::new(2.0, 2.0, 0.0),
532            Point3::new(4.0, 4.0, 10.0),
533        );
534
535        let profile_transform = Matrix4::identity();
536        let extrusion_direction = Vector3::new(0.0, 0.0, 1.0);
537        let extrusion_depth = 10.0;
538
539        let classification = analyzer.classify_void(
540            &void_mesh,
541            &profile_transform,
542            &extrusion_direction,
543            extrusion_depth,
544        );
545
546        match classification {
547            VoidClassification::Coplanar { is_through, .. } => {
548                assert!(is_through);
549            }
550            _ => panic!("Expected Coplanar classification"),
551        }
552    }
553
554    #[test]
555    fn test_void_analyzer_partial_depth() {
556        let analyzer = VoidAnalyzer::new();
557
558        // Create a box void that only goes halfway through
559        let void_mesh = create_box_mesh(
560            Point3::new(2.0, 2.0, 2.0),
561            Point3::new(4.0, 4.0, 8.0),
562        );
563
564        let profile_transform = Matrix4::identity();
565        let extrusion_direction = Vector3::new(0.0, 0.0, 1.0);
566        let extrusion_depth = 10.0;
567
568        let classification = analyzer.classify_void(
569            &void_mesh,
570            &profile_transform,
571            &extrusion_direction,
572            extrusion_depth,
573        );
574
575        match classification {
576            VoidClassification::Coplanar {
577                depth_start,
578                depth_end,
579                is_through,
580                ..
581            } => {
582                assert!(!is_through);
583                assert!(depth_start >= 1.9 && depth_start <= 2.1);
584                assert!(depth_end >= 7.9 && depth_end <= 8.1);
585            }
586            _ => panic!("Expected Coplanar classification"),
587        }
588    }
589
590    #[test]
591    fn test_extract_coplanar_voids() {
592        let classifications = vec![
593            VoidClassification::Coplanar {
594                profile_hole: vec![
595                    Point2::new(0.0, 0.0),
596                    Point2::new(1.0, 0.0),
597                    Point2::new(1.0, 1.0),
598                    Point2::new(0.0, 1.0),
599                ],
600                depth_start: 0.0,
601                depth_end: 10.0,
602                is_through: true,
603            },
604            VoidClassification::NonPlanar {
605                mesh: Mesh::new(),
606            },
607            VoidClassification::NonIntersecting,
608        ];
609
610        let coplanar = extract_coplanar_voids(&classifications);
611        assert_eq!(coplanar.len(), 1);
612        assert!(coplanar[0].is_through);
613    }
614
615    #[test]
616    fn test_compute_convex_hull() {
617        let analyzer = VoidAnalyzer::new();
618
619        let points = vec![
620            Point2::new(0.0, 0.0),
621            Point2::new(1.0, 0.0),
622            Point2::new(0.5, 0.5), // Interior point
623            Point2::new(1.0, 1.0),
624            Point2::new(0.0, 1.0),
625        ];
626
627        let hull = analyzer.compute_convex_hull(&points);
628
629        // Hull should have 4 points (excluding interior)
630        assert_eq!(hull.len(), 4);
631    }
632}