remesh 0.0.5

// SPDX-License-Identifier: MIT OR Apache-2.0
// Copyright (c) 2025 lacklustr@protonmail.com https://github.com/eadf

mod dihedral_angle;

#[cfg(test)]
mod tests;

use crate::bvh::StaticBVH;
use crate::common::VertexIndex;
use crate::common::macros::{integrity_assert, integrity_assert_eq, integrity_println};
use crate::common::remesh_error::RemeshError;
use crate::corner_table::{CornerIndex, CornerTable};
use crate::isotropic_remesh::{FlipStrategy, RemeshParams};
use crate::prelude::{CollapseStrategy, SplitStrategy};
use crate::util::index_pool::{Allocation, SifoPool};
use std::fmt::Debug;
use std::marker::PhantomData;
use vector_traits::glam::DVec3;
use vector_traits::num_traits::AsPrimitive;
use vector_traits::prelude::{GenericVector3, SimdUpgradable};
use vob::Vob;

/// The isotropic remesh builder
pub struct IsotropicRemesh<S, V, const ENABLE_UNSAFE: bool = false>
where
    S: crate::common::sealed::ScalarType,
    V: Copy + From<[S; 3]> + Into<[S; 3]> + Sync + 'static,
{
    pub(super) vertices: Vec<S::Vec3>, // S::Vec3 is one of glam::Vec3 or glam::DVec3
    pub(super) indices: Vec<VertexIndex>,
    pub(super) max_index: u32,
    pub(super) params: RemeshParams<S>,
    pub(super) _pd: PhantomData<V>,
}

impl<S, V, const ENABLE_UNSAFE: bool> IsotropicRemesh<S, V, ENABLE_UNSAFE>
where
    S: crate::common::sealed::ScalarType,
    f64: AsPrimitive<S>,
    V: Debug + Copy + From<[S; 3]> + Into<[S; 3]> + Sync + 'static,
{
    /// Run the isotropic remeshing algorithm without final mesh compression
    pub(crate) fn new_algo(self) -> Result<IsotropicRemeshAlgo<S, V, ENABLE_UNSAFE>, RemeshError> {
        let mut vertices = self.vertices;
        let corner_table = if self.params.fix_non_manifold {
            if let Some(vertex_limit) = self.params.print_stats {
                IsotropicRemeshAlgo::<S, V, ENABLE_UNSAFE>::print_input_status(
                    vertex_limit,
                    &vertices,
                    &self.indices,
                );
            }

            CornerTable::from_non_manifold_triangles(self.indices, self.max_index, &mut vertices)?
        } else {
            CornerTable::from_manifold_triangles(self.indices, self.max_index)?
        };
        let bvh = if self.params.smooth_weight.is_some() {
            StaticBVH::new(corner_table.vertex_of_corners(), &vertices)
        } else {
            StaticBVH::default()
        };
        let algo = IsotropicRemeshAlgo::<S, V, ENABLE_UNSAFE>::new(
            vertices,
            corner_table,
            bvh,
            self.params,
        );

        Ok(algo)
    }

    /// Run the isotropic remeshing algorithm
    pub fn run<I>(self, iterations: I) -> Result<(Vec<V>, Vec<u32>), RemeshError>
    where
        I: TryInto<u32>,
        <I as TryInto<u32>>::Error: Debug,
    {
        let mut algo = self.new_algo()?;
        algo.run_internal(iterations)?;
        algo.compress_mesh()
    }

    #[inline(always)]
    pub(crate) fn to_glam(v: V) -> S::Vec3 {
        v.into().into()
    }
}

/// The isotropic remesh algorithm instance
pub struct IsotropicRemeshAlgo<S, V, const ENABLE_UNSAFE: bool = false>
where
    S: crate::common::sealed::ScalarType,
    V: Copy + From<[S; 3]> + Into<[S; 3]> + Sync + 'static,
{
    pub(super) vertices: Vec<S::Vec3>, // S::Vec3 is one of glam::Vec3 or glam::DVec3
    pub(super) vertices_buffer: Vec<S::Vec3>, // smoothing will write to this and then swap
    pub(super) bvh: StaticBVH<S>,
    pub(super) vertex_pool: SifoPool<VertexIndex, ENABLE_UNSAFE>,
    pub(super) corner_table: CornerTable<ENABLE_UNSAFE>,
    pub(super) params: RemeshParams<S>,
    /// Indicator to tell if a vertex was already touched in this iteration and phase
    pub(super) dirty_vertices: DirtyVertices<ENABLE_UNSAFE>,
    _pd: PhantomData<V>,
}

impl<S, V, const ENABLE_UNSAFE: bool> IsotropicRemeshAlgo<S, V, ENABLE_UNSAFE>
where
    S: crate::common::sealed::ScalarType,
    f64: AsPrimitive<S>,
    V: Debug + Copy + From<[S; 3]> + Into<[S; 3]> + Sync + 'static,
{
    pub(super) fn new(
        vertices: Vec<S::Vec3>,
        corner_table: CornerTable<ENABLE_UNSAFE>,
        bvh: StaticBVH<S>,
        params: RemeshParams<S>,
    ) -> Self {
        let vertices_len = vertices.len();
        let vertices_buffer = vec![S::Vec3::ZERO; vertices_len];
        Self {
            vertices,
            vertices_buffer,
            bvh,
            corner_table,
            vertex_pool: SifoPool::<VertexIndex, ENABLE_UNSAFE>::new(vertices_len as u32),
            params,
            dirty_vertices: DirtyVertices::with_capacity(vertices_len),
            _pd: PhantomData,
        }
    }

    /// Run the isotropic remeshing algorithm without final compression
    pub(crate) fn run_internal<I>(&mut self, iterations: I) -> Result<(), RemeshError>
    where
        I: TryInto<u32>,
        <I as TryInto<u32>>::Error: Debug,
    {
        self.params.iterations = self.check_iterations(iterations)?;
        self.params.validate()?;

        self.params.print_essentials();
        if let Some(vertex_limit) = self.params.print_stats {
            Self::print_input_status(
                vertex_limit,
                &self.vertices,
                self.corner_table.vertex_of_corners(),
            );
        }
        for _i in 0..self.params.iterations {
            integrity_println!(
                "#######   running remesh iteration {:?}. Vertices:{} allocated vertices:{} triangles:{}",
                _i,
                self.vertex_pool.active_count(),
                self.vertices.len(),
                self.corner_table.active_triangles()
            );
            integrity_assert_eq!(self.vertex_pool.total_count() as usize, self.vertices.len());

            if !self.remesh_iteration() {
                integrity_println!("terminating with no changes after {} iterations ", _i + 1);
                break;
            }
        }
        Ok(())
    }

    /// Perform one iteration of remeshing
    fn remesh_iteration(&mut self) -> bool {
        let mut did_something = false;

        #[cfg(feature = "integrity_check")]
        self.check_mesh_integrity("start of iteration").unwrap();

        match self.params.split_strategy {
            SplitStrategy::DihedralAngle => {
                // 1. Split long edges
                did_something |= self.split_edges_dihedral(self.params.split_threshold_sq);
                #[cfg(feature = "integrity_check")]
                self.check_mesh_integrity("after split_edges_dihedral()")
                    .unwrap();
            }
            SplitStrategy::Disabled => {}
            _ => unimplemented!(),
        }

        match self.params.collapse_strategy {
            CollapseStrategy::DihedralAngle => {
                // 2. Collapse short edges using dihedral angle criteria
                did_something |= self.collapse_edges_dihedral(self.params.collapse_threshold_sq);
                #[cfg(feature = "integrity_check")]
                self.check_mesh_integrity("after collapse_edges_dihedral()")
                    .unwrap();
            }
            CollapseStrategy::Qem => {
                // 2. Collapse short edges using QEM
                did_something |= self.collapse_edges_qem(self.params.collapse_threshold_sq);
                #[cfg(feature = "integrity_check")]
                self.check_mesh_integrity("after collapse_edges_qem()")
                    .unwrap();
            }
            CollapseStrategy::Disabled => {}
            _ => unimplemented!(),
        }

        if let Some(defrag_ratio) = self.params.corner_table_defragmentation_ratio {
            let deleted_triangles: f64 = self.corner_table.deleted_triangle_count().as_();
            let total_triangle_count: f64 = self.corner_table.total_triangle_count().as_();

            if (deleted_triangles / total_triangle_count).as_() > defrag_ratio {
                // 3. defragment the corner table
                self.corner_table.defragment_triangles();
                #[cfg(feature = "integrity_check")]
                self.check_mesh_integrity("after defragment_triangles()")
                    .unwrap();
            }
        }

        match self.params.flip_strategy {
            FlipStrategy::Disabled => {}
            FlipStrategy::Valence => {
                // 4. Edge flips to improve triangle quality
                did_something |= self.flip_edges_dihedral();

                #[cfg(feature = "integrity_check")]
                self.check_mesh_integrity("after flip_edges_dihedral()")
                    .unwrap();
            }
            FlipStrategy::WeightedQuality {
                quality_threshold: quality_weight,
            } => {
                // 4. Edge flips to improve triangle quality
                did_something |= self.flip_edges_aspect_ratio(quality_weight);

                #[cfg(feature = "integrity_check")]
                self.check_mesh_integrity("after flip_edges_aspect_ratio()")
                    .unwrap();
            }
        }

        if let Some(weight) = self.params.smooth_weight {
            // 5. Vertex smoothing (tangential relaxation)
            did_something |= self.smooth_project_vertices(
                weight,
                self.params.max_projection_distance_sq,
                self.params.smooth_normal_threshold,
            );
            #[cfg(feature = "integrity_check")]
            self.check_mesh_integrity("after smooth_vertices()")
                .unwrap();
        }
        did_something
    }

    #[inline(always)]
    pub(crate) fn from_glam(v: S::Vec3) -> V {
        v.into().into()
    }

    #[inline(always)]
    /// Calculate edge length squared between two vertices
    pub(crate) fn edge_length_sq(&self, v1_idx: VertexIndex, v2_idx: VertexIndex) -> S {
        let v1 = self.vertex(v1_idx);
        let v2 = self.vertex(v2_idx);

        (v2 - v1).magnitude_sq()
    }

    #[inline(always)]
    pub(crate) fn vertex_of_corner(&self, ci: CornerIndex) -> S::Vec3 {
        debug_assert!(ci.is_valid());
        let vi = self.corner_table.vertex(ci);
        debug_assert!(vi.is_valid());
        self.vertex(vi)
    }

    #[inline(always)]
    pub(crate) fn vertex_of_corner_f64(&self, ci: CornerIndex) -> DVec3 {
        self.vertex_f64(self.corner_table.vertex(ci))
    }

    #[inline(always)]
    pub(crate) fn is_vertex_deleted(&self, vertex: VertexIndex) -> bool {
        !vertex.is_valid()
            || if cfg!(any(feature = "integrity_check", debug_assertions)) {
                let slow = !self.vertex_pool.is_used(vertex);
                let fast = !self.corner_table.corner(vertex).is_valid();
                assert_eq!(slow, fast);
                fast
            } else {
                !self.corner_table.corner(vertex).is_valid()
            }
    }

    #[inline(always)]
    /// Allocates a vertex slot and returns its index.
    ///
    /// Reuses an available slot from the pool if possible,
    /// otherwise creates a new slot. The provided vertex data
    /// is stored in the allocated slot.
    ///
    /// # Returns
    /// [`VertexIndex`] identifying the allocated vertex slot
    pub(crate) fn add_vertex(&mut self, v: S::Vec3Simd) -> VertexIndex {
        match self.vertex_pool.pop() {
            Allocation::New(index) => {
                self.vertices.push(S::Vec3::from_simd(v));
                self.corner_table.handle_new_vertex(index);
                integrity_println!("Added new vertex {index:?}:{:?}", v.into());
                index
            }
            Allocation::Reused(index) => {
                self.vertices[index.0 as usize] = S::Vec3::from_simd(v);
                self.corner_table.handle_reused_vertex(index);
                integrity_println!("Reused vertex {index:?}:{:?}", v.into());
                index
            }
        }
    }

    #[inline(always)]
    pub(crate) fn delete_vertex(&mut self, vertex: VertexIndex) {
        integrity_assert!(vertex.is_valid());
        integrity_assert!(
            !self.is_vertex_deleted(vertex),
            "vertex {vertex:?} already deleted"
        );
        self.vertex_pool.push(vertex);
        self.corner_table.delete_vertex(vertex);
    }

    #[inline(always)]
    pub(crate) fn vertex(&self, index: VertexIndex) -> S::Vec3 {
        #[cfg(feature = "integrity_check")]
        assert!(index.is_valid() && self.vertex_pool.is_used(index));

        debug_assert!(index.is_valid(), "{index:?}");
        if ENABLE_UNSAFE {
            unsafe { *self.vertices.get_unchecked(index.0 as usize) }
        } else {
            self.vertices[index.0 as usize]
        }
    }

    #[inline(always)]
    pub(crate) fn vertex_f64(&self, index: VertexIndex) -> DVec3 {
        let v = self.vertex(index);
        DVec3::new(v[0].as_(), v[1].as_(), v[2].as_())
    }

    /*#[allow(dead_code)]
    #[inline(always)]
    pub(crate) fn get_vertex_nf64(&self, index: VertexIndex) -> nalgebra::Vector3<f64> {
        let v = self.get_vertex(index);
        nalgebra::Vector3::new(v[0].as_(), v[1].as_(), v[2].as_())
    }*/

    /// prints input data info
    pub(crate) fn print_input_status(
        vertex_limit: usize,
        vertices: &[S::Vec3],
        indices: &[VertexIndex],
    ) {
        if vertices.len() <= vertex_limit {
            println!(
                "//******** vertices:{} indices:{}",
                vertices.len(),
                indices.len()
            );
            println!(
                "let vertices = {:?};",
                vertices
                    .iter()
                    .map(|v| (*v).into())
                    .collect::<Vec<[S; 3]>>(),
            );
            println!(
                "let indices = {:?};",
                indices.iter().map(|i| i.0).collect::<Vec<u32>>(),
            );
        } else {
            println!("********");
            println!("vertices:{} indices:{}", vertices.len(), indices.len());
        }
    }
}

pub(crate) struct DirtyVertices<const ENABLE_UNSAFE: bool>(Vob);

impl<const ENABLE_UNSAFE: bool> DirtyVertices<ENABLE_UNSAFE> {
    pub(crate) fn with_capacity(capacity: usize) -> Self {
        Self(Vob::from_elem(false, capacity))
    }

    #[inline(always)]
    pub(crate) fn mark_dirty(&mut self, vertex: VertexIndex) {
        integrity_assert!(vertex.is_valid());
        if ENABLE_UNSAFE {
            unsafe {
                let _ = self.0.set_unchecked(vertex.0 as usize, true);
            }
        } else {
            let _ = self.0.set(vertex.0 as usize, true);
        }
    }

    #[inline(always)]
    pub(crate) fn is_dirty(&self, vertex: VertexIndex) -> bool {
        integrity_assert!(vertex.is_valid());
        if ENABLE_UNSAFE {
            unsafe { self.0.get_unchecked(vertex.0 as usize) }
        } else {
            self.0.get(vertex.0 as usize).unwrap()
        }
    }

    #[inline(always)]
    pub(crate) fn clear(&mut self) {
        self.0.set_all(false);
    }

    pub(crate) fn prepare(&mut self, vertex_count: usize) {
        self.clear();
        // Resize if needed
        if self.0.len() < vertex_count {
            self.0.resize(vertex_count, false);
        }
    }
}