diff --git a/AABB_tree/demo/AABB_tree/MainWindow.cpp b/AABB_tree/demo/AABB_tree/MainWindow.cpp index 98563109327..513606d295e 100644 --- a/AABB_tree/demo/AABB_tree/MainWindow.cpp +++ b/AABB_tree/demo/AABB_tree/MainWindow.cpp @@ -283,3 +283,24 @@ void MainWindow::on_actionClear_distance_function_triggered() { m_pScene->clear_distance_function(); } + +void MainWindow::on_actionRefine_bisection_triggered() +{ + bool ok; + const double max_len = + QInputDialog::getDouble(NULL, "Max edge len", + "Max edge len:",0.1,0.001,100.0,9,&ok); + if(!ok) + return; + + QApplication::setOverrideCursor(Qt::WaitCursor); + m_pScene->refine_bisection(max_len * max_len); + QApplication::restoreOverrideCursor(); +} + +void MainWindow::on_actionBench_memory_triggered() +{ + QApplication::setOverrideCursor(Qt::WaitCursor); + m_pScene->bench_memory(); + QApplication::restoreOverrideCursor(); +} \ No newline at end of file diff --git a/AABB_tree/demo/AABB_tree/MainWindow.h b/AABB_tree/demo/AABB_tree/MainWindow.h index 20e186d426e..875fffce9a9 100644 --- a/AABB_tree/demo/AABB_tree/MainWindow.h +++ b/AABB_tree/demo/AABB_tree/MainWindow.h @@ -50,12 +50,14 @@ public: void on_actionEdge_points_triggered(); void on_actionInside_points_triggered(); void on_actionBoundary_points_triggered(); + void on_actionRefine_bisection_triggered(); void on_actionBoundary_segments_triggered(); void on_actionSigned_distance_function_to_facets_triggered(); void on_actionUnsigned_distance_function_to_edges_triggered(); void on_actionUnsigned_distance_function_to_facets_triggered(); // benchmark menu + void on_actionBench_memory_triggered(); void on_actionBench_distances_triggered(); void on_actionBench_intersections_triggered(); diff --git a/AABB_tree/demo/AABB_tree/MainWindow.ui b/AABB_tree/demo/AABB_tree/MainWindow.ui index 721616a2890..7ef86f82521 100644 --- a/AABB_tree/demo/AABB_tree/MainWindow.ui +++ b/AABB_tree/demo/AABB_tree/MainWindow.ui @@ -60,6 +60,12 @@ Algorithms + + + Refine + + + @@ -68,11 +74,14 @@ + + Benchmark + @@ -186,6 +195,21 @@ Clear distance function + + + Longest edge bisection + + + + + Loop subdivision + + + + + Memory + + diff --git a/AABB_tree/demo/AABB_tree/Refiner.h b/AABB_tree/demo/AABB_tree/Refiner.h new file mode 100644 index 00000000000..73e0619b6ed --- /dev/null +++ b/AABB_tree/demo/AABB_tree/Refiner.h @@ -0,0 +1,197 @@ +#ifndef _REFINER_H_ +#define _REFINER_H_ + +#include +#include +#include + +template +class CEdge +{ +public: + typedef typename Kernel::FT FT; + typedef typename Polyhedron::Halfedge_handle Halfedge_handle; + +private: + FT m_sqlen; + Halfedge_handle m_he; + +public: + // life cycle + CEdge(const Halfedge_handle& he) + { + m_sqlen = squared_len(he); + m_he = he; + } + CEdge(const CEdge& edge) + : m_sqlen(edge.sqlen()), + m_he(edge.he()) + { + } + ~CEdge() {} + +public: + // squared length of an edge + static FT squared_len(Halfedge_handle he) + { + return CGAL::squared_distance(he->vertex()->point(), + he->opposite()->vertex()->point()); + } + +public: + // accessors + FT& sqlen() { return m_sqlen; } + const FT sqlen() const { return m_sqlen; } + + Halfedge_handle he() { return m_he; } + const Halfedge_handle he() const { return m_he; } +}; + +// functor for priority queue +template +struct less // read more priority +{ + bool operator()(const Edge& e1, + const Edge& e2) const + { + return e1.sqlen() < e2.sqlen(); + } +}; + +template +class Refiner +{ + // types + typedef typename Kernel::FT FT; + typedef typename CEdge Edge; + typedef typename Polyhedron::Halfedge_handle Halfedge_handle; + typedef typename Polyhedron::Edge_iterator Edge_iterator; + typedef typename std::priority_queue, + less > PQueue; + // data + PQueue m_queue; + Polyhedron* m_pMesh; + +public : + // life cycle + Refiner(Polyhedron* pMesh) + { + m_pMesh = pMesh; + } + ~Refiner() {} + +public : + + void fill_queue(const FT& max_sqlen) + { + for(Edge_iterator he = m_pMesh->edges_begin(); + he != m_pMesh->edges_end(); + he++) + if(Edge::squared_len(he) > max_sqlen) + m_queue.push(Edge(he)); + } + + void fill_queue() + { + for(Edge_iterator he = m_pMesh->edges_begin(); + he != m_pMesh->edges_end(); + he++) + m_queue.push(Edge(he)); + } + + // run + void run_nb_splits(const unsigned int nb_splits) + { + // fill queue + fill_queue(); + + unsigned int nb = 0; + while(nb < nb_splits) + { + // get a copy of the candidate edge + Edge edge = m_queue.top(); + m_queue.pop(); + + Halfedge_handle he = edge.he(); + // update point + Halfedge_handle hnew = m_pMesh->split_edge(he); + hnew->vertex()->point() = CGAL::midpoint(he->vertex()->point(), + he->opposite()->vertex()->point()); + + // hit has been split into two edges + m_queue.push(Edge(hnew)); + m_queue.push(Edge(he)); + + // split facet if possible + if(!hnew->is_border()) + { + m_pMesh->split_facet(hnew,hnew->next()->next()); + m_queue.push(Edge(hnew->next())); + } + + // split facet if possible + if(!hnew->opposite()->is_border()) + { + m_pMesh->split_facet(hnew->opposite()->next(), + hnew->opposite()->next()->next()->next()); + m_queue.push(Edge(hnew->opposite()->prev())); + } + + nb++; + } // end while + } // end run + + + + // run + unsigned int operator()(const FT& max_sqlen) + { + // fill queue + fill_queue(max_sqlen); + + unsigned int nb_split = 0; + while(!m_queue.empty()) + { + // get a copy of the candidate edge + Edge edge = m_queue.top(); + m_queue.pop(); + + Halfedge_handle he = edge.he(); + FT sqlen = Edge::squared_len(he); + if(sqlen > max_sqlen) + { + // update point + Halfedge_handle hnew = m_pMesh->split_edge(he); + hnew->vertex()->point() = CGAL::midpoint(he->vertex()->point(), + he->opposite()->vertex()->point()); + + // hit has been split into two edges + m_queue.push(Edge(hnew)); + m_queue.push(Edge(he)); + + // split facet if possible + if(!hnew->is_border()) + { + m_pMesh->split_facet(hnew,hnew->next()->next()); + m_queue.push(Edge(hnew->next())); + } + + // split facet if possible + if(!hnew->opposite()->is_border()) + { + m_pMesh->split_facet(hnew->opposite()->next(), + hnew->opposite()->next()->next()->next()); + m_queue.push(Edge(hnew->opposite()->prev())); + } + + nb_split++; + } // end if(sqlen > max_sqlen) + } // end while(!m_queue.empty()) + return nb_split; + } // end run +}; + +#endif // _REFINER_H_ + + diff --git a/AABB_tree/demo/AABB_tree/Scene.cpp b/AABB_tree/demo/AABB_tree/Scene.cpp index 25d0461bc50..8a20980be2e 100644 --- a/AABB_tree/demo/AABB_tree/Scene.cpp +++ b/AABB_tree/demo/AABB_tree/Scene.cpp @@ -9,6 +9,7 @@ #include #include +#include "Refiner.h" #include "render_edges.h" #include @@ -556,9 +557,6 @@ unsigned_distance : m_max_distance_function; m_signed_distance_function = true; } - - - void Scene::toggle_view_poyhedron() { m_view_polyhedron = !m_view_polyhedron; @@ -579,4 +577,10 @@ void Scene::toggle_view_distance_function() m_view_distance_function = !m_view_distance_function; } - +void Scene::refine_bisection(const FT max_sqlen) +{ + std::cout << "Refine through recursive longest edge bisection..."; + Refiner refiner(m_pPolyhedron); + refiner(max_sqlen); + std::cout << "done (" << m_pPolyhedron->size_of_facets() << " facets)" << std::endl; +} diff --git a/AABB_tree/demo/AABB_tree/Scene.h b/AABB_tree/demo/AABB_tree/Scene.h index c6c5f29e197..09e4fb04709 100644 --- a/AABB_tree/demo/AABB_tree/Scene.h +++ b/AABB_tree/demo/AABB_tree/Scene.h @@ -38,7 +38,6 @@ public: public: void draw(); - Polyhedron* polyhedron() const; typedef CGAL::Bbox_3 Bbox; Bbox bbox() { return Bbox(); } @@ -66,6 +65,7 @@ public: void generate_inside_points(const unsigned int nb_points); void generate_boundary_points(const unsigned int nb_points); void generate_boundary_segments(const unsigned int nb_slices); + void refine_bisection(const FT max_sqlen); // distance functions void signed_distance_function(); @@ -96,6 +96,8 @@ public: NB_INTERSECTIONS, ALL_INTERSECTIONS, ALL_INTERSECTED_PRIMITIVES}; + void bench_memory(); + unsigned int nb_digits(const unsigned int value); void benchmark_intersections(const int duration); void bench_intersection_rays(Facet_tree& tree,const int function,const int duration); void bench_intersection_lines(Facet_tree& tree,const int function,const int duration); diff --git a/AABB_tree/demo/AABB_tree/benchmarks.cpp b/AABB_tree/demo/AABB_tree/benchmarks.cpp index baa2c8c3612..817e1d91dfb 100644 --- a/AABB_tree/demo/AABB_tree/benchmarks.cpp +++ b/AABB_tree/demo/AABB_tree/benchmarks.cpp @@ -1,5 +1,6 @@ #include "Scene.h" #include +#include void Scene::benchmark_intersections(const int duration) { @@ -35,6 +36,41 @@ void Scene::benchmark_distances(const int duration) bench_closest_point_and_primitive(tree,duration); } +unsigned int Scene::nb_digits(unsigned int value) +{ + unsigned int nb_digits = 0; + while(value > 0) + { + nb_digits++; + value /= 10; + } + return nb_digits; +} + +// bench memory against number of facets in the tree +// the tree is reconstructed each time in the mesh +// refinement loop +void Scene::bench_memory() +{ + std::cout << std::endl << "Benchmark memory" << std::endl; + std::cout << "#Facets MBytes" << std::endl; + while(m_pPolyhedron->size_of_facets() < 2000000) + { + // refine mesh at increasing speed + Refiner refiner(m_pPolyhedron); + unsigned int digits = nb_digits(m_pPolyhedron->size_of_facets()); + unsigned int nb_splits = 0.2 * std::pow(10.0,(double)digits - 1.0); + refiner.run_nb_splits(nb_splits); + + // constructs tree and measure memory before then after + long before = CGAL::Memory_sizer().virtual_size(); + Facet_tree tree(m_pPolyhedron->facets_begin(),m_pPolyhedron->facets_end()); + long after = CGAL::Memory_sizer().virtual_size(); + double memory = (double)(after - before) / 1048576.0; // in MBytes + std::cout << m_pPolyhedron->size_of_facets() << " " << memory << std::endl; + } +} + ////////////////////////////////////////////////////////////// // BENCH INTERSECTIONS //////////////////////////////////////////////////////////////