splashsurf_lib 0.14.1

Library for surface reconstruction of SPH particle data
Documentation 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95

use crate::marching_cubes::marching_cubes_lut::marching_cubes_triangulation_iter;
use crate::marching_cubes::{CellData, MarchingCubesInput};
use crate::mesh::TriMesh3d;
use crate::{Index, Real, profile};
use anyhow::Context;
use log::trace;
use thiserror::Error as ThisError;

/// Error enum for the marching cubes triangulation stage
#[derive(Debug, ThisError)]
pub enum TriangulationError {
    /// Error that occurred while generating the triangle connectivity for a case
    #[error("triangle connectivity error")]
    TriangleConnectivityError(
        #[source]
        #[from]
        anyhow::Error,
    ),
}

/// Converts the marching cubes input cell data into a triangle surface mesh, appends triangles to existing mesh
#[inline(never)]
pub(crate) fn triangulate<I: Index, R: Real>(
    input: MarchingCubesInput<I>,
    mesh: &mut TriMesh3d<R>,
) -> Result<(), TriangulationError> {
    profile!("triangulate");

    let MarchingCubesInput { cell_data } = input;

    trace!(
        "Starting marching cubes triangulation (Input: cell data map with {} cells, surface mesh with {} triangles and {} vertices)",
        cell_data.len(),
        mesh.triangles.len(),
        mesh.vertices.len(),
    );

    // Triangulate affected cells
    for (&_flat_cell_index, cell_data) in &cell_data {
        // TODO: Replace `are_vertices_above_unchecked` with something that can return an error
        for triangle in marching_cubes_triangulation_iter(&cell_data.are_vertices_above_unchecked())
        {
            // TODO: Allow user to set option to skip invalid triangles?
            let global_triangle = get_triangle(cell_data, triangle)
                .map_err(TriangulationError::TriangleConnectivityError)?;
            mesh.triangles.push(global_triangle);
        }
    }

    trace!(
        "Triangulation done. (Output: surface mesh with {} triangles and {} vertices)",
        mesh.triangles.len(),
        mesh.vertices.len()
    );

    Ok(())
}

/// Helper function that extracts vertex indices from the [`CellData`] for the triangle with the given edge indices
fn get_triangle(cell_data: &CellData, edge_indices: [i32; 3]) -> Result<[usize; 3], anyhow::Error> {
    let [edge_idx_0, edge_idx_1, edge_idx_2] = edge_indices;

    Ok([
        cell_data
            .iso_surface_vertices
            .get(edge_idx_0 as usize)
            .with_context(|| "Invalid edge index. This is a bug.")?
            .with_context(|| {
                format!(
                    "Missing iso surface vertex at edge {}. This is a bug.",
                    edge_idx_0
                )
            })?,
        cell_data
            .iso_surface_vertices
            .get(edge_idx_1 as usize)
            .with_context(|| "Invalid edge index. This is a bug.")?
            .with_context(|| {
                format!(
                    "Missing iso surface vertex at edge {}. This is a bug.",
                    edge_idx_1
                )
            })?,
        cell_data
            .iso_surface_vertices
            .get(edge_idx_2 as usize)
            .with_context(|| "Invalid edge index. This is a bug.")?
            .with_context(|| {
                format!(
                    "Missing iso surface vertex at edge {}. This is a bug.",
                    edge_idx_2
                )
            })?,
    ])
}