math_convex_hull/
lib.rs

1//! 3D Convex Hull and Computational Geometry Library
2//!
3//! This library implements the Quickhull algorithm for computing convex hulls
4//! in 3D space.
5//!
6//! Based on the C implementation by Leo McCormack and the MATLAB Computational
7//! Geometry Toolbox.
8//!
9//! # 3D Convex Hull Example
10//! ```
11//! use math_convex_hull::{ConvexHull3D, Vertex};
12//!
13//! let vertices = vec![
14//!     Vertex::new(0.0, 0.0, 0.0),
15//!     Vertex::new(1.0, 0.0, 0.0),
16//!     Vertex::new(0.0, 1.0, 0.0),
17//!     Vertex::new(0.0, 0.0, 1.0),
18//! ];
19//!
20//! let hull = ConvexHull3D::build(&vertices).unwrap();
21//! println!("Number of faces: {}", hull.num_faces());
22//! ```
23
24mod export;
25mod geometry;
26mod quickhull;
27mod types;
28
29// Make testdata publicly available for tests
30pub mod testdata;
31
32// 3D types and functions
33pub use export::{export_html, export_obj};
34pub use types::{ConvexHull3D, Face, Vertex};
35
36/// Error types for convex hull operations
37#[derive(Debug, thiserror::Error)]
38pub enum ConvexHullError {
39    #[error("Not enough vertices to form a hull (minimum 4 required)")]
40    InsufficientVertices,
41
42    #[error("Vertices are coplanar or collinear")]
43    DegenerateConfiguration,
44
45    #[error("Maximum iterations exceeded")]
46    MaxIterationsExceeded,
47
48    #[error("IO error: {0}")]
49    IoError(#[from] std::io::Error),
50
51    #[error("Invalid face: {0}")]
52    InvalidFace(String),
53}
54
55pub type Result<T> = std::result::Result<T, ConvexHullError>;
56
57/// Numerical tolerance for floating-point comparisons
58/// Used throughout the library for:
59/// - Distance calculations
60/// - Determinant checks
61/// - Degeneracy detection
62pub(crate) const EPSILON: f64 = 1e-10;
63
64/// Compute a scale-aware relative epsilon based on vertex coordinates
65pub(crate) fn compute_relative_epsilon(vertices: &[Vertex]) -> f64 {
66    if vertices.is_empty() {
67        return EPSILON;
68    }
69    let max_coord = vertices.iter().fold(0.0f64, |max, v| {
70        max.max(v.x.abs()).max(v.y.abs()).max(v.z.abs())
71    });
72    // Use relative epsilon based on coordinate scale
73    // For coordinates around 1.0, this gives ~1e-10
74    // For coordinates around 1e8, this gives ~1e3
75    (max_coord * EPSILON).max(EPSILON)
76}
77
78/// Deduplicate vertices, keeping the first occurrence of each unique point
79pub(crate) fn deduplicate_vertices(vertices: &[Vertex], epsilon: f64) -> Vec<Vertex> {
80    let mut unique_points: Vec<Vertex> = Vec::with_capacity(vertices.len());
81    let eps_sq = epsilon * epsilon;
82
83    for point in vertices {
84        let mut is_duplicate = false;
85        for existing in &unique_points {
86            let dx = point.x - existing.x;
87            let dy = point.y - existing.y;
88            let dz = point.z - existing.z;
89            if dx * dx + dy * dy + dz * dz < eps_sq {
90                is_duplicate = true;
91                break;
92            }
93        }
94        if !is_duplicate {
95            unique_points.push(*point);
96        }
97    }
98
99    unique_points
100}