CGAL users discussion list

[Sebastien Loriot <>]

Several remarks:

You are using float as NT. If the mesh was created with double, you are 
doing a potentially harmful rounding.

You are using self-intersections with a Kernel with inexact predicates 
so the output is not guaranteed to be valid. Using 
CGAL::Exact_predicates_inexact_constructions_kernel is recommended.

You are getting an assertion because faces are indeed removed but
the mesh connectivity is not updated. It basically simply marks the 
faces as removed. You need to use Euler operations for that:
https://doc.cgal.org/latest/BGL/group__PkgBGLEulerOperations.html#gacfae7ff8e782da55b941e4487e86c738

            if (!mesh.is_removed(inter_pair.first))
              CGAL::Euler::remove_face( 
halfedge(inter_pair.first,mesh), mesh );
            if (!mesh.is_removed(inter_pair.second))
              CGAL::Euler::remove_face( halfedge(inter_pair.second, 
mesh), mesh );

Then you will see that you'll get the following assertion:
terminate called after throwing an instance of 'CGAL::Assertion_exception'
  what():  CGAL ERROR: assertion violation!
Expr: abad1 != NULL_VECTOR
File: 
/home/sloriot/CGAL/git/master/Kernel_23/include/CGAL/Kernel/function_objects.h
Line: 251
Explanation: ab1 and ad1 are collinear


That is because of the inexact kernel used.

About face removal, I advice you to use consider the set of faces to 
remove as a selection and use that function:
https://doc.cgal.org/latest/BGL/group__PkgBGLSelectionFct.html#ga6f8d338df28c24bbf93677f24f7bcbbe
to avoid creating new non-manifold vertices.

A function like this one
https://doc.cgal.org/latest/Polygon_mesh_processing/group__keep__connected__components__grp.html#ga17cb532198888791f2ceb00666f4c1d4

could then be used to remove all selected faces (considering them as one 
connected component).

HTH

Sebastien.


On 12/2/21 6:07 AM, Yaoyu Hu ( via cgal-discuss 
Mailing List) wrote:
> Hi,
>
> I am trying to fix and clear some issues from a surface mesh. I do the 
> following three steps in sequence:
> (1) Duplicate non-manifold vertices.
> (2) Remove self-intersections. And,
> (3) Remove isolated vertices.
> (4) Remove small connected components.
>
> I always get an assertion error in (4). In my code, (1), (3), and (4) 
> are directly calling CGAL functions to do the job. For (2) I have a 
> simple procedure where I first identify all the self-intersection pairs 
> and then remove all the faces that have self-intersection. I have set up 
> a sample code and a sample mesh. The code is also available as a GitHub 
> Gist
> https://gist.github.com/huyaoyu/c10ff67a2fe67c141af77142e475a83d 
> <https://gist.github.com/huyaoyu/c10ff67a2fe67c141af77142e475a83d>
>
> The sample mesh is saved in my personal Google Drive.
> https://drive.google.com/drive/folders/1AmQONiss-q26cgs4wNBp3fMSlVmo09Ej?usp=sharing 
> <https://drive.google.com/drive/folders/1AmQONiss-q26cgs4wNBp3fMSlVmo09Ej?usp=sharing>
>
> I am using CGAL 5.3 on Ubuntu 20.04.
>
> === Code begins. ===
>
> #include <fstream>
> #include <iostream>
> #include <sstream>
> #include <vector>
>
> #include <CGAL/Simple_cartesian.h>
> #include <CGAL/Surface_mesh.h>
>
> #include <CGAL/Surface_mesh/IO/PLY.h>
>
> #include <CGAL/Polygon_mesh_processing/connected_components.h>
> #include <CGAL/Polygon_mesh_processing/repair.h>
> #include <CGAL/Polygon_mesh_processing/self_intersections.h>
>
> namespace PMP = CGAL::Polygon_mesh_processing;
>
> typedef CGAL::Simple_cartesian<float> Kernel_t;
> typedef Kernel_t::Point_3  Point_t;
>
> typedef CGAL::Surface_mesh<Point_t> Mesh_t;
>
> typedef Kernel_t::Compare_dihedral_angle_3 Compare_dihedral_angle_3;
>
> template < typename MeshT >
> void read_mesh(const std::string& fn, MeshT& mesh) {
>      std::ifstream ifs(fn, std::ios::binary);
>
>      if ( !ifs ) {
>          std::stringstream ss;
>          ss << "Cannot open " << fn;
>          throw std::runtime_error(ss.str());
>      }
>
>      std::string comments;
>      CGAL::IO::read_PLY(ifs, mesh, comments);
> }
>
> template < typename MeshT >
> void duplicate_non_manifold_vertices(MeshT& mesh) {
>      PMP::duplicate_non_manifold_vertices(mesh);
> }
>
> template < typename MeshT >
> int remove_self_intersections( MeshT& mesh ) {
>      typedef typename boost::graph_traits<MeshT>::face_descriptor Face_t;
>
>      std::vector< std::pair< Face_t, Face_t > > intersections;
>      PMP::self_intersections<CGAL::Parallel_if_available_tag>(
>          faces(mesh),
>          mesh,
>          std::back_inserter(intersections)
>      );
>
>      const int n_intersections = intersections.size();
>      if ( n_intersections > 0 ) {
>          for ( const auto& inter_pair : intersections ) {
>              mesh.remove_face( inter_pair.first );
>              mesh.remove_face( inter_pair.second );
>          }
>
>          mesh.collect_garbage();
>      }
>
>      return n_intersections;
> }
>
> template < typename MeshT >
> int remove_isolated_vertices( MeshT& mesh ) {
>      const int N = PMP::remove_isolated_vertices(mesh);
>      return N;
> }
>
> // For computing the connected components.
> template < typename MeshT >
> struct ConnectionConstraint : public boost::put_get_helper< bool, 
> ConnectionConstraint<MeshT> > {
>      typedef typename boost::graph_traits<MeshT>::edge_descriptor EdgeT;
>
>      typedef boost::readable_property_map_tag category;
>      typedef bool                             value_type;
>      typedef bool                             reference;
>      typedef EdgeT                            key_type;
>
>      ConnectionConstraint()
>      : m_mesh(nullptr) {}
>
>      ConnectionConstraint( MeshT& mesh, double bound )
>      : m_mesh(&mesh), m_bound(bound) {}
>
>      bool operator[] ( EdgeT edge ) const {
>          const MeshT& mesh = *m_mesh;
>          return m_compare(
>              mesh.point( source( edge, mesh ) ),
>              mesh.point( target( edge, mesh ) ),
>              mesh.point( target( next( halfedge( edge, mesh ), mesh ), 
> mesh ) ),
>              mesh.point( target( next( opposite( halfedge( edge, mesh ), 
> mesh), mesh ) ,mesh ) ),
>              m_bound ) == CGAL::SMALLER;
>      }
>
>      const MeshT* m_mesh;
>      Compare_dihedral_angle_3 m_compare;
>      double m_bound;
> };
>
> template < typename MeshT >
> void remove_small_connected_components(
>      MeshT& mesh,
>      int limit_num=10,
>      double limit_di_angle_rad=2.36 /* 3/4 pi */ ) {
>
>      const double bound = std::cos( limit_di_angle_rad );
>
>      PMP::keep_large_connected_components(
>          mesh,
>          limit_num,
>          PMP::parameters::edge_is_constrained_map(
>              ConnectionConstraint<MeshT>(mesh, bound) )
>      );
> }
>
> int main(int argc, char** argv) {
>      const std::string in_mesh_fn =
>          "please download from 
> https://drive.google.com/drive/folders/1AmQONiss-q26cgs4wNBp3fMSlVmo09Ej?usp=sharing 
> <https://drive.google.com/drive/folders/1AmQONiss-q26cgs4wNBp3fMSlVmo09Ej?usp=sharing>";
>
>      // Read the mesh.
>      Mesh_t mesh;
>      read_mesh( in_mesh_fn, mesh );
>
>      // Duplicate non-manifold_vertices.
>      duplicate_non_manifold_vertices(mesh);
>
>      // Remove self-intersection faces.
>      remove_self_intersections(mesh);
>
>      // Remove isolated vertices.
>      remove_isolated_vertices(mesh);
>
>      // Remove small connected components.
>      remove_small_connected_components(mesh);
>
>      return 0;
> }
>
> === Code ends. ===
>
> Step (2) is the remove_self_intersections() function. The assertion 
> error shows in my terminal is
>
> ===  Terminal output begins. ===
>
> terminate called after throwing an instance of 'CGAL::Assertion_exception'
>    what():  CGAL ERROR: assertion violation!
> Expr: _idx < data_.size()
> File: /usr/local/include/CGAL/Surface_mesh/Properties.h
> Line: 193
>
> === Terminal output ends. ===
>
> I have also saved the backtrace from GDB in the Google Drive folder 
> (a txt and a screenshot). please download from 
> https://drive.google.com/drive/folders/1AmQONiss-q26cgs4wNBp3fMSlVmo09Ej?usp=sharing 
> <https://drive.google.com/drive/folders/1AmQONiss-q26cgs4wNBp3fMSlVmo09Ej?usp=sharing>
>
>  From the backtrace, it looks like when iterating through the faces of 
> the mesh for computing the connected components, a face descriptor is 
> out of the valid range.
>
> Another thing that I found during the debugging is that, if you check 
> the output of mesh.has_garbage() after calling 
> PMP::remove_isolated_vertices() in my remove_isolated_vertices() 
> function, mesh.has_garbage() is true. I suppose  
> PMP::remove_isolated_vertices() should not leave the mesh in a state 
> where it has uncollected garbage.
>
> Any help is appreciated!
>
> Thank you!
>
> Yaoyu
>
>
> --
> You are currently subscribed to cgal-discuss.
> To unsubscribe or access the archives, go to
> https://sympa.inria.fr/sympa/info/cgal-discuss