Merge pull request #1150 from Geode-solutions/feat/geometry_adding_degeneration_functions

Feat/geometry_adding_degeneration_functions

Author: MelchiorSchuh

include/geode/geometry/basic_objects/polygon.hpp

@@ -77,6 +77,7 @@ namespace geode
         [[nodiscard]] BoundingBox< dimension > bounding_box() const;
         [[nodiscard]] std::vector< std::array< index_t, 3 > >
             triangulate() const;
+        [[nodiscard]] double minimum_height() const;
         [[nodiscard]] bool is_degenerated() const;
         [[nodiscard]] std::string string() const;
 

include/geode/mesh/core/solid_mesh.hpp

@@ -592,6 +592,9 @@ namespace geode
         [[nodiscard]] std::optional< PolyhedronVertex >
             polyhedron_around_vertex( index_t vertex_id ) const;
 
+        [[nodiscard]] double polyhedron_minimum_height(
+            index_t polyhedron_id ) const;
+
     public:
         void associate_polyhedron_vertex_to_vertex(
             const PolyhedronVertex& polyhedron_vertex,

include/geode/mesh/core/surface_mesh.hpp

@@ -442,6 +442,8 @@ namespace geode
         [[nodiscard]] std::optional< PolygonVertex > polygon_around_vertex(
             index_t vertex_id ) const;
 
+        [[nodiscard]] double polygon_minimum_height( index_t polygon_id ) const;
+
     public:
         void associate_polygon_vertex_to_vertex(
             const PolygonVertex& polygon_vertex,

src/geode/geometry/basic_objects/polygon.cpp

@@ -235,43 +235,47 @@ namespace geode
     }
 
     template < typename PointType, index_t dimension >
-    bool GenericPolygon< PointType, dimension >::is_degenerated() const
+    double GenericPolygon< PointType, dimension >::minimum_height() const
     {
+        const auto nb_vertices = vertices_.size();
         double max_length{ 0. };
-        index_t max_length_edge{ 0 };
-        for( const auto e : Range{ nb_vertices() } )
+        geode::index_t max_length_edge{ 0 };
+        for( const auto e : geode::Range{ nb_vertices } )
         {
-            const Point< dimension >& point0 = vertices_[e];
-            const Point< dimension >& point1 =
-                vertices_[e == nb_vertices() - 1 ? 0 : e + 1];
-            const auto cur_length = point_point_distance( point0, point1 );
+            const geode::Point< dimension >& point0 = vertices_[e];
+            const geode::Point< dimension >& point1 =
+                vertices_[e == nb_vertices - 1 ? 0 : e + 1];
+            const auto cur_length =
+                geode::point_point_distance( point0, point1 );
             if( cur_length > max_length )
             {
                 max_length = cur_length;
                 max_length_edge = e;
             }
         }
-        if( max_length < GLOBAL_EPSILON )
-        {
-            return true;
-        }
         const auto next =
-            max_length_edge + 1 == nb_vertices() ? 0 : max_length_edge + 1;
-        const InfiniteLine< dimension > line{ Segment< dimension >{
-            vertices_[max_length_edge], vertices_[next] } };
-        for( const auto v : Range{ nb_vertices() } )
+            max_length_edge + 1 == nb_vertices ? 0 : max_length_edge + 1;
+        const geode::InfiniteLine< dimension > line{
+            geode::Segment< dimension >{
+                vertices_[max_length_edge], vertices_[next] }
+        };
+        auto opposite_vertex = 0;
+        for( const auto vertex : geode::Range{ nb_vertices } )
         {
-            if( v == max_length_edge || v == next )
+            if( vertex == max_length_edge || vertex == next )
             {
                 continue;
             }
-            const Point< dimension >& point = vertices_[v];
-            if( point_line_distance( point, line ) > GLOBAL_EPSILON )
-            {
-                return false;
-            }
+            opposite_vertex = vertex;
         }
-        return true;
+        const geode::Point< dimension >& point = vertices_[opposite_vertex];
+        return geode::point_line_distance( point, line );
+    }
+
+    template < typename PointType, index_t dimension >
+    bool GenericPolygon< PointType, dimension >::is_degenerated() const
+    {
+        return minimum_height() < GLOBAL_EPSILON;
     }
 
     template < typename PointType, index_t dimension >

src/geode/mesh/core/solid_mesh.cpp

@@ -458,7 +458,7 @@ namespace geode
         }
 
         std::optional< PolyhedronVertex > polyhedron_around_vertex(
-            const index_t vertex_id ) const
+            index_t vertex_id ) const
         {
             const auto& value = polyhedron_around_vertex_->value( vertex_id );
             if( value.polyhedron_id != NO_ID )
@@ -468,6 +468,46 @@ namespace geode
             return std::nullopt;
         }
 
+        double polyhedron_minimum_height(
+            const geode::SolidMesh< dimension >& mesh,
+            const index_t polyhedron_id ) const
+        {
+            double max_area{ 0. };
+            local_index_t max_area_facet{ 0 };
+            for( const auto f :
+                LRange{ mesh.nb_polyhedron_facets( polyhedron_id ) } )
+            {
+                const auto cur_area =
+                    mesh.polyhedron_facet_area( { polyhedron_id, f } );
+                if( cur_area > max_area )
+                {
+                    max_area = cur_area;
+                    max_area_facet = f;
+                }
+            }
+            const auto vertices = mesh.polyhedron_vertices( polyhedron_id );
+            const auto normal = mesh.polyhedron_facet_normal(
+                { polyhedron_id, max_area_facet } );
+            if( !normal )
+            {
+                return true;
+            }
+            const auto facet_vertices = mesh.polyhedron_facet_vertices(
+                { polyhedron_id, max_area_facet } );
+            Plane plane{ normal.value(), mesh.point( facet_vertices[0] ) };
+            auto opposite_vertex{ 0 };
+            for( const auto vertex_id : vertices )
+            {
+                if( absl::c_contains( facet_vertices, vertex_id ) )
+                {
+                    continue;
+                }
+                opposite_vertex = vertex_id;
+            }
+            return std::get< 0 >(
+                point_plane_distance( mesh.point( opposite_vertex ), plane ) );
+        }
+
         void reset_polyhedra_around_vertex( index_t vertex_id )
         {
             polyhedra_around_vertex_->modify_value(
@@ -960,41 +1000,7 @@ namespace geode
     bool SolidMesh< dimension >::is_polyhedron_degenerated(
         index_t polyhedron_id ) const
     {
-        double max_area{ 0. };
-        local_index_t max_area_facet{ 0 };
-        for( const auto f : LRange{ nb_polyhedron_facets( polyhedron_id ) } )
-        {
-            const auto cur_area = polyhedron_facet_area( { polyhedron_id, f } );
-            if( cur_area > max_area )
-            {
-                max_area = cur_area;
-                max_area_facet = f;
-            }
-        }
-        const auto vertices = polyhedron_vertices( polyhedron_id );
-        const auto normal =
-            polyhedron_facet_normal( { polyhedron_id, max_area_facet } );
-        if( !normal )
-        {
-            return true;
-        }
-        const auto facet_vertices =
-            polyhedron_facet_vertices( { polyhedron_id, max_area_facet } );
-        Plane plane{ normal.value(), this->point( facet_vertices[0] ) };
-        for( const auto vertex_id : vertices )
-        {
-            if( absl::c_contains( facet_vertices, vertex_id ) )
-            {
-                continue;
-            }
-            if( std::get< 0 >(
-                    point_plane_distance( this->point( vertex_id ), plane ) )
-                > GLOBAL_EPSILON )
-            {
-                return false;
-            }
-        }
-        return true;
+        return polyhedron_minimum_height( polyhedron_id ) < GLOBAL_EPSILON;
     }
 
     template < index_t dimension >
@@ -1358,6 +1364,13 @@ namespace geode
         return get_polyhedron_around_vertex( vertex_id );
     }
 
+    template < index_t dimension >
+    double SolidMesh< dimension >::polyhedron_minimum_height(
+        const index_t polyhedron_id ) const
+    {
+        return impl_->polyhedron_minimum_height( *this, polyhedron_id );
+    }
+
     template < index_t dimension >
     std::optional< PolyhedronVertex >
         SolidMesh< dimension >::get_polyhedron_around_vertex(

src/geode/mesh/core/surface_mesh.cpp

@@ -285,6 +285,14 @@ namespace geode
             return result;
         }
 
+        double polygon_minimum_height(
+            const SurfaceMesh< dimension >& mesh, index_t polygon_id ) const
+        {
+            OPENGEODE_EXCEPTION( polygon_id < mesh.nb_polygons(),
+                "[Impl::polygon_minimum_height] Wrong polygon id" );
+            return mesh.polygon( polygon_id ).minimum_height();
+        }
+
         std::optional< PolygonVertex > polygon_around_vertex(
             const index_t vertex_id ) const
         {
@@ -594,6 +602,13 @@ namespace geode
         return get_polygon_around_vertex( vertex_id );
     }
 
+    template < index_t dimension >
+    double SurfaceMesh< dimension >::polygon_minimum_height(
+        index_t polygon_id ) const
+    {
+        return impl_->polygon_minimum_height( *this, polygon_id );
+    }
+
     template < index_t dimension >
     std::optional< PolygonVertex >
         SurfaceMesh< dimension >::get_polygon_around_vertex(

src/geode/mesh/helpers/triangulated_surface_scalar_function.cpp

@@ -25,6 +25,7 @@
 
 #include <geode/basic/attribute_manager.hpp>
 #include <geode/basic/pimpl_impl.hpp>
+#include <geode/basic/variable_attribute.hpp>
 
 #include <geode/geometry/barycentric_coordinates.hpp>
 #include <geode/geometry/basic_objects/triangle.hpp>