CGAL users discussion list

["Sebastien Loriot (GeometryFactory)" <>]

TroyMcLure wrote:
> i did many attempts (both on linux and windows) and my conclusion is: or cgal
> mesh simplification is done really bad , or i am doing mistakes...
> Please consider the following example. I want to do a surface simplification
> with Edge_length_cost ( the cost is the squared distance from 2 points) and
> Midpoint_placement ( the point is choosed half distance, in the exact middle
> of the edge). I want to do my own stop predicate, that stops when the
> collapsing cost is greater then a certain value.
>
> so this is the code of my example (part of this is from the cgal examples):
>
> #include <iostream>
> #include <fstream>
>
> #include <CGAL/Simple_cartesian.h>
> #include <CGAL/Polyhedron_3.h>
> #include <CGAL/IO/Polyhedron_iostream.h>
>
> // Adaptor for Polyhedron_3
> #include <CGAL/Surface_mesh_simplification/HalfedgeGraph_Polyhedron_3.h>
>
> // Simplification function
> #include <CGAL/Surface_mesh_simplification/edge_collapse.h>
>
> // Extended polyhedron items which include an id() field
> #include <CGAL/Polyhedron_items_with_id_3.h>
>
> // Non-default cost and placement policies
> #include
> <CGAL/Surface_mesh_simplification/Policies/Edge_collapse/Midpoint_and_length.h> 
>
> #include
> <CGAL/Surface_mesh_simplification/Policies/Edge_collapse/Edge_profile.h>
>
> typedef CGAL::Simple_cartesian<double> Kernel;
>
> typedef Kernel::Point_3 Point ;
>
> //
> // Setup an enriched polyhedron type which stores an id() field in the items
> //
> typedef CGAL::Polyhedron_3<Kernel,CGAL::Polyhedron_items_with_id_3> Surface; 
>
> typedef Surface::Halfedge_handle Halfedge_handle ;
> typedef Surface::Vertex_handle   Vertex_handle ;
>
> namespace SMS = CGAL::Surface_mesh_simplification ;
>
> typedef SMS::Edge_profile<Surface> Profile ;
>
> typedef  boost::graph_traits<Surface>::edges_size_type size_type ;
>
> typedef Kernel::FT FT;
>
> typedef bool (*Max_cost_predicate)( Kernel::FT const& , Profile const& ,
> size_type , size_type);
>       
> //Custom stop predicate, returning when the current cost is greater then
> 0.01
>
>   bool my_stop_predicate( FT const&  aCurrentCost
>                  , Profile const& //aEdgeProfile
>                  , size_type      //aInitialCount
>                  , size_type      //aCurrentCount
>                  )  
>   {
>     
> 	  bool quitting= aCurrentCost>0.01;
> 	  
> 	  if (quitting) std::cerr<<"Quitting, cause current cost
> ("<<aCurrentCost<<") is greater then target cost(0.01)"<<std::endl;
>         return quitting;
>
>   }
>
>
> int main( int argc, char** argv ) 
> {
>   Surface surface; 
>   
>   std::ifstream is(argv[1]) ; is >> surface ;
>
>   // The items in this polyhedron have an "id()" field 
>   // which the default index maps used in the algorithm
>   // need to get the index of a vertex/edge.
>   // However, the Polyhedron_3 class doesn't assign any value to
>   // this id(), so we must do it here:
>   int index = 0 ;
>   
>   for( Surface::Halfedge_iterator eb = surface.halfedges_begin()
>      , ee = surface.halfedges_end()
>      ; eb != ee
>      ; ++ eb
>      ) 
>     eb->id() = index++;
>
>   index = 0 ;
>   for( Surface::Vertex_iterator vb = surface.vertices_begin()
>      , ve = surface.vertices_end()
>      ; vb != ve
>      ; ++ vb
>      ) 
>     vb->id() = index++;
>     
>   int r = SMS::edge_collapse<Surface,Max_cost_predicate>
>            (surface
>            ,my_stop_predicate
>            ,CGAL::get_cost     (SMS::Edge_length_cost  <Surface>())
>                  .get_placement(SMS::Midpoint_placement<Surface>())
>             );
>     
>   std::cout << "Finished...\n" << r << " edges removed.\n" 
>             << (surface.size_of_halfedges()/2) << " final edges.\n" ;
>
>   Surface::Halfedge_iterator itB= surface.halfedges_begin();
>   Surface::Halfedge_iterator itE= surface.halfedges_end();
>
>   int halfedges_with_a_greter_cost=0;
>   int halfedges_with_a_lower_cost=0;
>     for (;itB!=itE;itB++){
>
>                 Surface::Halfedge_handle hhCurr= itB;
>                 FT error=
> CGAL::squared_distance(hhCurr->vertex()->point(),hhCurr->prev()->vertex()->point());
>                 if (error>0.01) halfedges_with_a_greter_cost++;
>                 else halfedges_with_a_lower_cost++;
>         }
>   
> 	std::cout<<"halfedges_with_a_greter_cost =
> "<<halfedges_with_a_greter_cost<<std::endl;
>         std::cout<<"halfedges_with_a_lower_cost =
> "<<halfedges_with_a_lower_cost<<std::endl;
>
>         std::ofstream os( argc > 2 ? argv[2] : "out.off" ) ; os << surface ;
>   
>   return 0 ;      
> }
>
>
> /*
> output is:
>
> C:\Users\TroyMcLure\Desktop\cgalRel22\build\Debug>esempio zampa2.off
> Quitting, cause current cost (0.0111412) is greater then target cost(0.01)
> Finished...
> 28339 edges removed.
> 1000 final edges.
> halfedges_with_a_greter_cost = 1688
> halfedges_with_a_lower_cost = 312
>
> */
>
> // EOF //
>
> so the output of the program tells this: the algorithm quits cause it
> selects an edge with a greater cost (0.0111412) then the target cost (0.01).
> But there are 312/2 edges with a cost lower then 0.01 in the final mesh! 
> all this 312/2 edges have no topologial issues ( i checked ). So what is it?
> surface mesh simplification does not work, or i am doing mistakes? ... and
> where?
From the user manual:
"Not all edges selected for processing are collapsed. A processed edge
can be discarded right away, without being collapsed, if it doesn't
satisfy certain topological and geometric conditions."

Have you checked that geometric condition was validated?

Looking at void EdgeCollapse<M,SP,VIM,EIM,EBM,CF,PF,V>::Loop() (l 151)
in include/CGAL/Surface_mesh_simplification/Detail/Edge_collapse_impl.h
you will see that if an edge is topologivally valid but not collapsed
it is because collapsing it is geometrically invalid.

Is_collapse_geometrically_valid is defined in the same file l. 567
the comments in the code says:
"A collase is geometrically valid if, in the resulting local mesh no two
 adjacent triangles form an internal dihedral angle greater than a fixed
 threshold (i.e. triangles do not "fold" into each other)"


HTH

S.