diff --git a/.gitignore b/.gitignore index a25c4ac1e74..3d8f09650c8 100644 --- a/.gitignore +++ b/.gitignore @@ -446,6 +446,8 @@ Surface_reconstruction_3/doc_doxygen Surface_reconstruction_3/doc_html Surface_reconstruction_3/doc_pdf Surface_reconstruction_3/doc_ps +Surface_reconstruction_3/examples/Surface_reconstruction_3/polygirl_apss.off +Surface_reconstruction_3/examples/Surface_reconstruction_3/polygirl_poisson.off Surface_reconstruction_3/test/Surface_reconstruction_3/*.kdev* Surface_reconstruction_3/test/Surface_reconstruction_3/*.ncb Surface_reconstruction_3/test/Surface_reconstruction_3/*.suo diff --git a/Surface_reconstruction_3/demo/Surface_reconstruction_3/poisson/PoissonDoc.cpp b/Surface_reconstruction_3/demo/Surface_reconstruction_3/poisson/PoissonDoc.cpp index 5b54277bb1b..d6e21723e02 100644 --- a/Surface_reconstruction_3/demo/Surface_reconstruction_3/poisson/PoissonDoc.cpp +++ b/Surface_reconstruction_3/demo/Surface_reconstruction_3/poisson/PoissonDoc.cpp @@ -20,6 +20,7 @@ #include #include #include +#include // This package #include @@ -326,7 +327,7 @@ BOOL CPoissonDoc::OnOpenDocument(LPCTSTR lpszPathName) m_points.invalidate_bounding_box(); m_edit_mode = POINT_SET; - status_message("Load point set...done (%lf s)", task_timer.time()); + status_message("Load point set...done (%.2lf s)", task_timer.time()); update_status(); UpdateAllViews(NULL); return TRUE; @@ -477,8 +478,12 @@ void CPoissonDoc::update_status() CString selected_points; selected_points.Format("%d selected",m_points.nb_selected_points()); + long memory = CGAL::Memory_sizer().virtual_size(); + // write message to cerr - std::cerr << "=> " << points << " (" << selected_points << ")" << std::endl; + std::cerr << "=> " << points << " (" << selected_points << "), " + << (memory>>20) << " Mb allocated" + << std::endl; // Update status bar pStatus->SetPaneText(1,points); @@ -495,8 +500,12 @@ void CPoissonDoc::update_status() CString tets; tets.Format("%d tets",m_poisson_dt->number_of_cells()); + long memory = CGAL::Memory_sizer().virtual_size(); + // write message to cerr - std::cerr << "=> " << vertices << " and " << tets << std::endl; + std::cerr << "=> " << vertices << ", " << tets << ", " + << (memory>>20) << " Mb allocated" + << std::endl; // Update status bar pStatus->SetPaneText(1,vertices); @@ -597,7 +606,7 @@ void CPoissonDoc::OnAlgorithmsEstimateNormalsByPCA() m_points.normals_begin(), m_number_of_neighbours); - status_message("Estimate Normals Direction by PCA...done (%lf s)", task_timer.time()); + status_message("Estimate Normals Direction by PCA...done (%.2lf s)", task_timer.time()); update_status(); UpdateAllViews(NULL); EndWaitCursor(); @@ -619,7 +628,7 @@ void CPoissonDoc::OnAlgorithmsEstimateNormalsByJetFitting() m_points.normals_begin(), m_number_of_neighbours); - status_message("Estimate Normals Direction by Jet Fitting...done (%lf s)", task_timer.time()); + status_message("Estimate Normals Direction by Jet Fitting...done (%.2lf s)", task_timer.time()); update_status(); UpdateAllViews(NULL); EndWaitCursor(); @@ -649,7 +658,7 @@ void CPoissonDoc::OnAlgorithmsOrientNormalsWithMST() if ( ! m_points[i].normal().is_oriented() ) m_points.select(&m_points[i]); - status_message("Orient Normals with MST...done (%lf s)", task_timer.time()); + status_message("Orient Normals with MST...done (%.2lf s)", task_timer.time()); update_status(); UpdateAllViews(NULL); EndWaitCursor(); @@ -684,7 +693,7 @@ void CPoissonDoc::OnAlgorithmsOrientNormalsWrtCameras() if (m_points[i].normal().get_vector() * normals_copy[i].get_vector() < 0) m_points.select(&m_points[i]); - status_message("Orient Normals wrt Cameras...done (%lf s)", task_timer.time()); + status_message("Orient Normals wrt Cameras...done (%.2lf s)", task_timer.time()); update_status(); UpdateAllViews(NULL); EndWaitCursor(); @@ -720,7 +729,7 @@ void CPoissonDoc::OnCreatePoissonTriangulation() m_triangulation_refined = false; // Need to apply Delaunay refinement m_poisson_solved = false; // Need to solve Poisson equation - status_message("Create Poisson Triangulation...done (%lf s)", task_timer.time()); + status_message("Create Poisson Triangulation...done (%.2lf s)", task_timer.time()); update_status(); UpdateAllViews(NULL); EndWaitCursor(); @@ -752,7 +761,8 @@ void CPoissonDoc::OnReconstructionDelaunayRefinement() unsigned int nb_vertices_added = m_poisson_function->delaunay_refinement(quality,max_vertices,enlarge_ratio,50000); m_triangulation_refined = true; - status_message("Delaunay refinement...done (%lf s, %d vertices inserted)", task_timer.time(),nb_vertices_added); + status_message("Delaunay refinement...done (%.2lf s, %d vertices inserted)", + task_timer.time(), nb_vertices_added); update_status(); UpdateAllViews(NULL); EndWaitCursor(); @@ -779,7 +789,8 @@ void CPoissonDoc::OnAlgorithmsRefineInShell() unsigned int nb_vertices_added = m_poisson_function->delaunay_refinement_shell(m_dr_shell_size,m_dr_sizing,m_dr_max_vertices); m_triangulation_refined = true; - status_message("Delaunay refinement in surface shell...done (%lf s, %d vertices inserted)", task_timer.time(),nb_vertices_added); + status_message("Delaunay refinement in surface shell...done (%.2lf s, %d vertices inserted)", + task_timer.time(), nb_vertices_added); update_status(); UpdateAllViews(NULL); EndWaitCursor(); @@ -803,7 +814,7 @@ void CPoissonDoc::OnAlgorithmsExtrapolatenormals() m_poisson_function->extrapolate_normals(); - status_message("Extrapolate the normals field...done (%lf s)", task_timer.time()); + status_message("Extrapolate the normals field...done (%.2lf s)", task_timer.time()); update_status(); UpdateAllViews(NULL); EndWaitCursor(); @@ -835,11 +846,10 @@ void CPoissonDoc::OnReconstructionPoisson() m_poisson_function->set_contouring_value(m_poisson_function->median_value_at_input_vertices()); m_contouring_value = 0.0; - double total_duration = task_timer.time(); if (!m_poisson_solved) status_message("Solve Poisson equation...solver failed"); else - status_message("Solve Poisson equation...done (%lf s)", total_duration); + status_message("Solve Poisson equation...done (%.2lf s)", task_timer.time()); update_status(); UpdateAllViews(NULL); EndWaitCursor(); @@ -865,11 +875,10 @@ void CPoissonDoc::OnReconstructionPoissonNormalized() m_poisson_function->set_contouring_value(m_poisson_function->median_value_at_input_vertices()); m_contouring_value = 0.0; - double total_duration = task_timer.time(); if (!m_poisson_solved) status_message("Solve Poisson equation...solver failed"); else - status_message("Solve Poisson equation...done (%lf s)", total_duration); + status_message("Solve Poisson equation...done (%.2lf s)", task_timer.time()); update_status(); UpdateAllViews(NULL); EndWaitCursor(); @@ -926,10 +935,13 @@ void CPoissonDoc::OnReconstructionPoissonSurfaceMeshing() m_sm_radius*size, // upper bound of Delaunay balls radii m_sm_distance_poisson*size); // upper bound of distance to surface -std::cerr << "Implicit_surface_3(dichotomy error="<set_contouring_value(0.0); // Print status - status_message("Surface meshing...done (%d vertices, %lf s)", + status_message("Surface meshing...done (%d output vertices, %.2lf s)", m_surface_mesher_dt.number_of_vertices(), task_timer.time()); update_status(); UpdateAllViews(NULL); @@ -972,7 +984,8 @@ void CPoissonDoc::OnAlgorithmsMarchingTetContouring() int nb = m_poisson_dt->marching_tet(std::back_inserter(triangles), m_contouring_value); m_contour.insert(m_contour.end(), triangles.begin(), triangles.end()); - status_message("Marching tet contouring...done (%d triangles, %lf s)",nb, task_timer.time()); + status_message("Marching tet contouring...done (%d triangles, %.2lf s)", + nb, task_timer.time()); update_status(); UpdateAllViews(NULL); EndWaitCursor(); @@ -1027,7 +1040,7 @@ void CPoissonDoc::OnAlgorithmsSmoothUsingJetFitting() m_points[i].set_position(output[i]); m_points.invalidate_bounding_box(); - status_message("Smooth Point Set...done (%lf s)", task_timer.time()); + status_message("Smooth Point Set...done (%.2lf s)", task_timer.time()); update_status(); UpdateAllViews(NULL); EndWaitCursor(); @@ -1103,7 +1116,7 @@ void CPoissonDoc::OnOneStepPoissonReconstruction() OnReconstructionPoisson(); OnReconstructionPoissonSurfaceMeshing(); - status_message("1-step Poisson reconstruction...done (%lf s)", task_timer.time()); + status_message("1-step Poisson reconstruction...done (%.2lf s)", task_timer.time()); update_status(); UpdateAllViews(NULL); EndWaitCursor(); @@ -1127,7 +1140,7 @@ void CPoissonDoc::OnOneStepPoissonReconstructionWithNormalizedDivergence() OnReconstructionPoissonNormalized(); OnReconstructionPoissonSurfaceMeshing(); - status_message("1-step Poisson reconstruction with normalized divergence...done (%lf s)", task_timer.time()); + status_message("1-step Poisson reconstruction with normalized divergence...done (%.2lf s)", task_timer.time()); update_status(); UpdateAllViews(NULL); EndWaitCursor(); @@ -1147,7 +1160,7 @@ void CPoissonDoc::OnAlgorithmsOutliersRemovalWrtCamerasConeAngle() m_min_cameras_cone_angle*M_PI/180.0); m_points.select(first_iterator_to_remove, m_points.end(), true); - status_message("Remove outliers / cameras cone's angle is low...done (%lf s)", task_timer.time()); + status_message("Remove outliers / cameras cone's angle is low...done (%.2lf s)", task_timer.time()); update_status(); UpdateAllViews(NULL); EndWaitCursor(); @@ -1177,7 +1190,7 @@ void CPoissonDoc::OnOutliersRemovalWrtAvgKnnSqDist() m_threshold_percent_avg_knn_sq_dst); m_points.select(first_iterator_to_remove, m_points.end(), true); - status_message("Remove outliers wrt average squared distance to knn...done (%lf s)", task_timer.time()); + status_message("Remove outliers wrt average squared distance to knn...done (%.2lf s)", task_timer.time()); update_status(); UpdateAllViews(NULL); EndWaitCursor(); @@ -1223,6 +1236,7 @@ void CPoissonDoc::OnReconstructionApssReconstruction() m_surface.clear(); // Create implicit function + std::cerr << " APSS_implicit_function(knn="< - - diff --git a/Surface_reconstruction_3/examples/Surface_reconstruction_3/APSS_reconstruction.cpp b/Surface_reconstruction_3/examples/Surface_reconstruction_3/APSS_reconstruction.cpp index 0ad15fd67b3..d1a042c179a 100644 --- a/Surface_reconstruction_3/examples/Surface_reconstruction_3/APSS_reconstruction.cpp +++ b/Surface_reconstruction_3/examples/Surface_reconstruction_3/APSS_reconstruction.cpp @@ -18,6 +18,7 @@ #include // include basic.h before testing #defines #include #include +#include #include // Surface mesher @@ -65,43 +66,18 @@ typedef CGAL::Surface_mesh_complex_2_in_triangulation_3 C2t3; typedef CGAL::Implicit_surface_3 Surface_3; -// ---------------------------------------------------------------------------- -// Private functions -// ---------------------------------------------------------------------------- - -//// Scale point set to [-1,1]^3 -//void reshape(PointList& pwns) -//{ -// Point cmin = pwns[0]; -// Point cmax = pwns[0]; -// for (unsigned int i=1 ; idiag.y() ? (diag.x()>diag.z() ? diag.x() : diag.z()) : (diag.y()>diag.z() ? diag.y() : diag.z())); -// for (unsigned int i=0 ; i Radius upper bound (default = 0.1 * point set radius)\n"; - std::cerr << " -sm_distance Distance upper bound (default = 0.005 * point set radius)\n"; - std::cerr << " -k Number of neighbors (default = 7)\n"; - std::cerr << " - should be greater than 7,\n"; - std::cerr << " - high numbers lead to smoother surfaces.\n"; - return(EXIT_FAILURE); + std::cerr << " -sm_radius Radius upper bound (default=0.1 * point set radius)\n"; + std::cerr << " -sm_distance Distance upper bound (default=0.005 * point set radius)\n"; + std::cerr << " -k Number of neighbors (default=7)\n"; + std::cerr << " - should be greater than 7,\n"; + std::cerr << " - high numbers lead to smoother surfaces.\n"; + return EXIT_FAILURE; } // Default APSS options @@ -221,21 +197,30 @@ int main(int argc, char * argv[]) } //*************************************** - // Create implicit function + // Compute implicit function //*************************************** + std::cerr << "Compute implicit function...\n"; + CGAL::Timer task_timer; task_timer.start(); - //reshape(pwns); // Scale point set to [-1,1]^3 - // Create implicit function + std::cerr << " APSS_implicit_function(knn="<>20) << " Mb allocated" + << std::endl; + task_timer.reset(); + //*************************************** // Surface mesh generation //*************************************** + std::cerr << "Surface meshing...\n"; STr tr; // 3D-Delaunay triangulation C2t3 c2t3 (tr); // 2D-complex in 3D-Delaunay triangulation @@ -264,19 +249,23 @@ int main(int argc, char * argv[]) CGAL::Surface_mesh_default_criteria_3 criteria(sm_angle, // lower bound of facets angles (degrees) sm_radius*size, // upper bound of Delaunay balls radii sm_distance*size); // upper bound of distance to surface - -std::cerr << "APSS_implicit_function(knn="<>20) << " Mb allocated" << std::endl; task_timer.reset(); @@ -287,6 +276,6 @@ std::cerr << "make_surface_mesh(sphere={center=("< #include +#include #include // Surface mesher @@ -73,12 +74,12 @@ typedef CGAL::Implicit_surface_3 Surface_3; int main(int argc, char * argv[]) { - std::cerr << "RECONSTRUCTION" << std::endl; - std::cerr << "Poisson Delaunay Reconstruction method" << std::endl; + std::cerr << "RECONSTRUCTION" << std::endl; + std::cerr << "Poisson Delaunay Reconstruction method" << std::endl; - //*************************************** - // decode parameters - //*************************************** + //*************************************** + // decode parameters + //*************************************** if (argc<3) { @@ -86,9 +87,9 @@ int main(int argc, char * argv[]) std::cerr << "Input file formats are .off and .xyz.\n"; std::cerr << "Output file format is .off.\n"; std::cerr << "Options:\n"; - std::cerr << " -sm_radius Radius upper bound (default = 0.1 * point set radius)\n"; - std::cerr << " -sm_distance Distance upper bound (default = 0.002 * point set radius)\n"; - return(EXIT_FAILURE); + std::cerr << " -sm_radius Radius upper bound (default=0.1 * point set radius)\n"; + std::cerr << " -sm_distance Distance upper bound (default=0.002 * point set radius)\n"; + return EXIT_FAILURE; } // Default Surface Mesher options @@ -197,6 +198,8 @@ int main(int argc, char * argv[]) // Compute implicit function //*************************************** + std::cerr << "Compute implicit function...\n"; + CGAL::Timer task_timer; task_timer.start(); // Create implicit function @@ -211,15 +214,20 @@ int main(int argc, char * argv[]) } // Print status + long memory = CGAL::Memory_sizer().virtual_size(); int nb_vertices2 = dt.number_of_vertices(); - std::cerr << "Solve Poisson equation: " - << "(added " << nb_vertices2-nb_vertices << " vertices)" - << std::endl; + std::cerr << "Compute implicit function: " << "added " << nb_vertices2-nb_vertices << " Steiner points, " + << task_timer.time() << " seconds, " + << (memory>>20) << " Mb allocated" + << std::endl; + task_timer.reset(); //*************************************** // Surface mesh generation //*************************************** + std::cerr << "Surface meshing...\n"; + STr tr; // 3D-Delaunay triangulation C2t3 c2t3 (tr); // 2D-complex in 3D-Delaunay triangulation @@ -248,18 +256,23 @@ int main(int argc, char * argv[]) CGAL::Surface_mesh_default_criteria_3 criteria(sm_angle, // lower bound of facets angles (degrees) sm_radius*size, // upper bound of Delaunay balls radii sm_distance*size); // upper bound of distance to surface - -std::cerr << "Implicit_surface_3(dichotomy error="<>20) << " Mb allocated" << std::endl; task_timer.reset(); @@ -270,7 +283,7 @@ std::cerr << "make_surface_mesh(sphere={center=("< #include #include +#include +#include #include #include CGAL_BEGIN_NAMESPACE + +// Traces? +//#define CGAL_TRACE printf + +#ifndef CGAL_TRACE + #define CGAL_TRACE if (false) printf +#endif + + /// Estimate normal direction using jet fitting /// on the K nearest neighbors. /// @@ -66,12 +77,11 @@ estimate_normal_jet_fitting_3(const typename Kernel::Point_3& query, ///< 3D poi typedef typename CGAL::Monge_via_jet_fitting Monge_jet_fitting; typedef typename Monge_jet_fitting::Monge_form Monge_form; - // gather set of (KNN+1) neighboring points + // Gather set of (KNN+1) neighboring points. + // Perform KNN+1 queries (as in point set, the query point is + // output first). Search may be aborted when KNN is greater + // than number of input points. std::vector points; - - // performs KNN + 1 queries (if unique the query point is - // output first). search may be aborted when KNN is greater - // than number of input points Neighbor_search search(tree,query,KNN+1); Search_iterator search_iterator = search.begin(); unsigned int i; @@ -120,6 +130,8 @@ estimate_normals_jet_fitting_3(InputIterator first, ///< input points const Kernel& /*kernel*/, const unsigned int degre_fitting = 2) { + CGAL_TRACE("Call estimate_normals_jet_fitting_3()\n"); + // value_type_traits is a workaround as back_insert_iterator's value_type is void typedef typename value_type_traits::type Normal; @@ -127,7 +139,6 @@ estimate_normals_jet_fitting_3(InputIterator first, ///< input points typedef typename CGAL::Search_traits_3 Tree_traits; typedef typename CGAL::Orthogonal_k_neighbor_search Neighbor_search; typedef typename Neighbor_search::Tree Tree; - typedef typename Neighbor_search::iterator Search_iterator; // precondition: at least one element in the container. // to fix: should have at least three distinct points @@ -137,9 +148,15 @@ estimate_normals_jet_fitting_3(InputIterator first, ///< input points // precondition: at least 2 nearest neighbors CGAL_surface_reconstruction_precondition(KNN >= 2); + long memory = CGAL::Memory_sizer().virtual_size(); CGAL_TRACE(" %ld Mb allocated\n", memory>>20); + CGAL_TRACE(" Create KD-tree\n"); + // instanciate a KD-tree search Tree tree(first,beyond); + /*long*/ memory = CGAL::Memory_sizer().virtual_size(); CGAL_TRACE(" %ld Mb allocated\n", memory>>20); + CGAL_TRACE(" Compute normals\n"); + // iterate over input points, compute and output normal // vectors (already normalized) InputIterator it; @@ -148,6 +165,10 @@ estimate_normals_jet_fitting_3(InputIterator first, ///< input points *normals = estimate_normal_jet_fitting_3(*it,tree,KNN,degre_fitting); normals++; } + + /*long*/ memory = CGAL::Memory_sizer().virtual_size(); CGAL_TRACE(" %ld Mb allocated\n", memory>>20); + CGAL_TRACE("End of estimate_normals_jet_fitting_3()\n"); + return normals; } @@ -178,6 +199,9 @@ estimate_normals_jet_fitting_3(InputIterator first, ///< input points } +// Avoid clash with other header +#undef CGAL_TRACE + CGAL_END_NAMESPACE #endif // CGAL_ESTIMATE_NORMALS_JET_FITTING_3_H diff --git a/Surface_reconstruction_3/include/CGAL/estimate_normals_pca_3.h b/Surface_reconstruction_3/include/CGAL/estimate_normals_pca_3.h index 9e76353cc33..037a86afa8b 100644 --- a/Surface_reconstruction_3/include/CGAL/estimate_normals_pca_3.h +++ b/Surface_reconstruction_3/include/CGAL/estimate_normals_pca_3.h @@ -26,12 +26,22 @@ #include #include #include +#include #include #include CGAL_BEGIN_NAMESPACE + +// Traces? +//#define CGAL_TRACE printf + +#ifndef CGAL_TRACE + #define CGAL_TRACE if (false) printf +#endif + + /// Estimate normal direction using linear least /// squares fitting of a plane on the K nearest neighbors. /// @@ -63,12 +73,11 @@ estimate_normal_pca_3(const typename Kernel::Point_3& query, ///< 3D point whose typedef typename CGAL::Orthogonal_k_neighbor_search Neighbor_search; typedef typename Neighbor_search::iterator Search_iterator; - // gather set of (KNN+1) neighboring points - std::list points; - - // performs KNN + 1 queries (if unique the query point is - // output first). search may be aborted when KNN is greater - // than number of input points + // Gather set of (KNN+1) neighboring points. + // Perform KNN+1 queries (as in point set, the query point is + // output first). Search may be aborted when KNN is greater + // than number of input points. + std::vector points; Neighbor_search search(tree,query,KNN+1); Search_iterator search_iterator = search.begin(); unsigned int i; @@ -83,11 +92,7 @@ estimate_normal_pca_3(const typename Kernel::Point_3& query, ///< 3D point whose // performs plane fitting by point-based PCA Plane plane; -#ifndef CGAL_DIMENSION_H // if CGAL 3.3.1 - linear_least_squares_fitting_3(points.begin(),points.end(),plane); -#else // if CGAL >= 3.4 linear_least_squares_fitting_3(points.begin(),points.end(),plane,Dimension_tag<0>()); -#endif // output normal vector (already normalized by PCA) return OrientedNormal_3(plane.orthogonal_vector(), @@ -118,6 +123,8 @@ estimate_normals_pca_3(InputIterator first, ///< input points const unsigned int KNN, ///< number of neighbors const Kernel& /*kernel*/) { + CGAL_TRACE("Call estimate_normals_pca_3()\n"); + // value_type_traits is a workaround as back_insert_iterator's value_type is void typedef typename value_type_traits::type Normal; @@ -134,9 +141,15 @@ estimate_normals_pca_3(InputIterator first, ///< input points // precondition: at least 2 nearest neighbors CGAL_surface_reconstruction_precondition(KNN >= 2); + long memory = CGAL::Memory_sizer().virtual_size(); CGAL_TRACE(" %ld Mb allocated\n", memory>>20); + CGAL_TRACE(" Create KD-tree\n"); + // instanciate a KD-tree search Tree tree(first,beyond); + /*long*/ memory = CGAL::Memory_sizer().virtual_size(); CGAL_TRACE(" %ld Mb allocated\n", memory>>20); + CGAL_TRACE(" Compute normals\n"); + // iterate over input points, compute and output normal // vectors (already normalized) InputIterator it; @@ -145,6 +158,10 @@ estimate_normals_pca_3(InputIterator first, ///< input points *normals = estimate_normal_pca_3(*it,tree,KNN); normals++; } + + /*long*/ memory = CGAL::Memory_sizer().virtual_size(); CGAL_TRACE(" %ld Mb allocated\n", memory>>20); + CGAL_TRACE("End of estimate_normals_pca_3()\n"); + return normals; } @@ -174,6 +191,9 @@ estimate_normals_pca_3(InputIterator first, ///< input points } +// Avoid clash with other header +#undef CGAL_TRACE + CGAL_END_NAMESPACE #endif // CGAL_ESTIMATE_NORMALS_PCA_3_H diff --git a/Surface_reconstruction_3/include/CGAL/remove_outliers_wrt_avg_knn_sq_distance_3.h b/Surface_reconstruction_3/include/CGAL/remove_outliers_wrt_avg_knn_sq_distance_3.h index 3c74afc8641..b86d0eb28c5 100644 --- a/Surface_reconstruction_3/include/CGAL/remove_outliers_wrt_avg_knn_sq_distance_3.h +++ b/Surface_reconstruction_3/include/CGAL/remove_outliers_wrt_avg_knn_sq_distance_3.h @@ -57,10 +57,11 @@ compute_avg_knn_sq_distance_3( typedef typename Neighbor_search::iterator Search_iterator; // Gather set of (KNN+1) neighboring points. - // Perform KNN + 1 queries (if in point set, the query point is + // Perform KNN+1 queries (if in point set, the query point is // output first). Search may be aborted when KNN is greater // than number of input points. std::vector points; + points.reserve(KNN+1); Neighbor_search search(tree,query,KNN+1); Search_iterator search_iterator = search.begin(); unsigned int i; diff --git a/Surface_reconstruction_3/test/Surface_reconstruction_3/APSS_reconstruction_test.cpp b/Surface_reconstruction_3/test/Surface_reconstruction_3/APSS_reconstruction_test.cpp index a5dab275dd2..96349d443dc 100644 --- a/Surface_reconstruction_3/test/Surface_reconstruction_3/APSS_reconstruction_test.cpp +++ b/Surface_reconstruction_3/test/Surface_reconstruction_3/APSS_reconstruction_test.cpp @@ -7,8 +7,8 @@ //---------------------------------------------------------- // Test the APSS reconstruction method -// No output -// Input files are .off and .xyz +// Input file formats are .off and .xyz. +// No output. //---------------------------------------------------------- // APSS_reconstruction_test mesh1.off point_set2.xyz... @@ -18,6 +18,7 @@ #include // include basic.h before testing #defines #include #include +#include #include // Surface mesher @@ -64,31 +65,6 @@ typedef CGAL::Surface_mesh_complex_2_in_triangulation_3 C2t3; typedef CGAL::Implicit_surface_3 Surface_3; -// ---------------------------------------------------------------------------- -// Private functions -// ---------------------------------------------------------------------------- - -//// Scale point set to [-1,1]^3 -//void reshape(PointList& pwns) -//{ -// Point cmin = pwns[0]; -// Point cmax = pwns[0]; -// for (unsigned int i=1 ; idiag.y() ? (diag.x()>diag.z() ? diag.x() : diag.z()) : (diag.y()>diag.z() ? diag.y() : diag.z())); -// for (unsigned int i=0 ; i>20) << " Mb allocated" + << std::endl; + task_timer.reset(); + //*************************************** // Surface mesh generation //*************************************** + std::cerr << "Surface meshing...\n"; + STr tr; // 3D-Delaunay triangulation C2t3 c2t3 (tr); // 2D-complex in 3D-Delaunay triangulation @@ -256,19 +247,23 @@ int main(int argc, char * argv[]) CGAL::Surface_mesh_default_criteria_3 criteria(sm_angle, // lower bound of facets angles (degrees) sm_radius*size, // upper bound of Delaunay balls radii sm_distance*size); // upper bound of distance to surface - -std::cerr << "APSS_implicit_function(knn="<>20) << " Mb allocated" << std::endl; task_timer.reset(); diff --git a/Surface_reconstruction_3/test/Surface_reconstruction_3/VC/APSS_reconstruction_test_vc80.vcproj b/Surface_reconstruction_3/test/Surface_reconstruction_3/VC/APSS_reconstruction_test_vc80.vcproj index b3de620d2b9..26a0c857647 100644 --- a/Surface_reconstruction_3/test/Surface_reconstruction_3/VC/APSS_reconstruction_test_vc80.vcproj +++ b/Surface_reconstruction_3/test/Surface_reconstruction_3/VC/APSS_reconstruction_test_vc80.vcproj @@ -41,12 +41,13 @@ Name="VCCLCompilerTool" AdditionalOptions="/Zm900" Optimization="0" - AdditionalIncludeDirectories="..\include;..\..\..\include;..\..\..\..\Spatial_searching\include;"$(CGALROOT)\include\CGAL\config\msvc";"$(CGALROOT)\include";"$(CGALROOT)\auxiliary\gmp\include";"$(BOOSTROOT)"" + AdditionalIncludeDirectories="..\include;..\..\..\include;..\..\..\..\Spatial_searching\include;..\..\..\..\Jet_fitting_3\include;"$(CGALROOT)\include\CGAL\config\msvc";"$(CGALROOT)\include";"$(CGALROOT)\auxiliary\gmp\include";"$(BOOSTROOT)"" PreprocessorDefinitions="WIN32;_SCL_SECURE_NO_DEPRECATE;_CRT_SECURE_NO_DEPRECATE;CGAL_USE_GMP;CGAL_SURFACE_RECONSTRUCTION_CHECK_EXPENSIVE" MinimalRebuild="true" BasicRuntimeChecks="3" RuntimeLibrary="3" RuntimeTypeInfo="true" + ProgramDataBaseFileName="$(IntDir)\$(ProjectName).pdb" BrowseInformation="1" WarningLevel="2" Detect64BitPortabilityProblems="true" @@ -64,6 +65,7 @@ /> -#include // This package #include diff --git a/Surface_reconstruction_3/test/Surface_reconstruction_3/normal_estimation_test.cpp b/Surface_reconstruction_3/test/Surface_reconstruction_3/normal_estimation_test.cpp index 888e210f616..f1e4103b75c 100644 --- a/Surface_reconstruction_3/test/Surface_reconstruction_3/normal_estimation_test.cpp +++ b/Surface_reconstruction_3/test/Surface_reconstruction_3/normal_estimation_test.cpp @@ -1,5 +1,6 @@ // normal_estimation_test.cpp + // ---------------------------------------------------------------------------- // USAGE EXAMPLES // ---------------------------------------------------------------------------- @@ -20,6 +21,7 @@ // CGAL #include #include +#include #include // This package @@ -57,18 +59,30 @@ void test_pca(const std::vector& points, // input point set std::vector& normals, // computed normals unsigned int k) // number of neighbors { - std::cerr << " Estimate normals using KNN and point-based PCA..."; + std::cerr << "Estimate normals using KNN and point-based PCA..." << std::endl; + CGAL::Timer task_timer; task_timer.start(); + CGAL::estimate_normals_pca_3(points.begin(),points.end(),std::back_inserter(normals),k); - std::cerr << "ok" << std::endl; + + long memory = CGAL::Memory_sizer().virtual_size(); + std::cerr << "ok: " << task_timer.time() << " seconds, " + << (memory>>20) << " Mb allocated" + << std::endl; } void test_jet_fitting(const std::vector& points, // input point set - std::vector& normals, // computed normals - unsigned int k) // number of neighbors) + std::vector& normals, // computed normals + unsigned int k) // number of neighbors) { - std::cerr << " Estimate normals using KNN and jet fitting..."; + std::cerr << "Estimate normals using KNN and jet fitting..." << std::endl; + CGAL::Timer task_timer; task_timer.start(); + CGAL::estimate_normals_jet_fitting_3(points.begin(),points.end(),std::back_inserter(normals),k); - std::cerr << "ok" << std::endl; + + long memory = CGAL::Memory_sizer().virtual_size(); + std::cerr << "ok: " << task_timer.time() << " seconds, " + << (memory>>20) << " Mb allocated" + << std::endl; } void test_orient_normals_MST( @@ -76,7 +90,8 @@ void test_orient_normals_MST( std::vector& normals, // normals to orient unsigned int k) // number of neighbors { - std::cerr << " Orient normals using a minimum spanning tree..."; + std::cerr << "Orient normals using a minimum spanning tree..." << std::endl; + CGAL::Timer task_timer; task_timer.start(); // orient_normals_minimum_spanning_tree_3() requires an iterator over points // + property maps to access each point's index, position and normal. @@ -89,7 +104,11 @@ void test_orient_normals_MST( boost::make_iterator_property_map(normals.begin(), index_id), // index -> normal prop. map k); - std::cerr << "ok" << std::endl; + + long memory = CGAL::Memory_sizer().virtual_size(); + std::cerr << "ok: " << task_timer.time() << " seconds, " + << (memory>>20) << " Mb allocated" + << std::endl; } @@ -101,9 +120,11 @@ int main(int argc, char * argv[]) { std::cerr << "Normal estimation test" << std::endl; std::cerr << "No output" << std::endl; + +char* leak1 = new char[100]; // decode parameters - if(argc < 2) + if (argc-1 == 0) { std::cerr << "Usage: " << argv[0] << " file1.xyz file2.xyz ..." << std::endl; return EXIT_FAILURE; @@ -119,9 +140,11 @@ int main(int argc, char * argv[]) { std::cerr << std::endl; +char* leak2 = new char[20]; + // Load point set std::vector points; - std::cerr << " Open " << argv[i] << " for reading..."; + std::cerr << "Open " << argv[i] << " for reading..." << std::endl; if(CGAL::surface_reconstruction_read_xyz(argv[i], std::back_inserter(points), false /*skip normals*/)) @@ -135,7 +158,7 @@ int main(int argc, char * argv[]) } else { - std::cerr << " FATAL ERROR: cannot read file " << argv[i] << std::endl; + std::cerr << "FATAL ERROR: cannot read file " << argv[i] << std::endl; accumulated_fatal_err = EXIT_FAILURE; } } // for each input file diff --git a/Surface_reconstruction_3/test/Surface_reconstruction_3/outliers_removal_test.cpp b/Surface_reconstruction_3/test/Surface_reconstruction_3/outliers_removal_test.cpp index 810dcdd16b3..ec7bdef3072 100644 --- a/Surface_reconstruction_3/test/Surface_reconstruction_3/outliers_removal_test.cpp +++ b/Surface_reconstruction_3/test/Surface_reconstruction_3/outliers_removal_test.cpp @@ -16,7 +16,6 @@ // CGAL #include -#include #include // This package diff --git a/Surface_reconstruction_3/test/Surface_reconstruction_3/poisson_reconstruction_test.cpp b/Surface_reconstruction_3/test/Surface_reconstruction_3/poisson_reconstruction_test.cpp index 537a7e3801f..940887c2120 100644 --- a/Surface_reconstruction_3/test/Surface_reconstruction_3/poisson_reconstruction_test.cpp +++ b/Surface_reconstruction_3/test/Surface_reconstruction_3/poisson_reconstruction_test.cpp @@ -7,8 +7,8 @@ //---------------------------------------------------------- // Test the Poisson Delaunay Reconstruction method -// No output -// Input files are .off and .xyz +// Input file formats are .off and .xyz. +// No output. //---------------------------------------------------------- // poisson_reconstruction_test mesh1.off point_set2.xyz... @@ -21,6 +21,7 @@ // CGAL #include #include +#include #include // Surface mesher @@ -74,7 +75,6 @@ int main(int argc, char * argv[]) { std::cerr << "RECONSTRUCTION" << std::endl; std::cerr << "Test the Poisson Delaunay Reconstruction method" << std::endl; - std::cerr << "No output" << std::endl; //*************************************** // decode parameters @@ -83,33 +83,38 @@ int main(int argc, char * argv[]) if (argc-1 == 0) { std::cerr << "Usage: " << argv[0] << " mesh1.off point_set2.xyz ..." << std::endl; - return(EXIT_FAILURE); + std::cerr << "Input file formats are .off and .xyz.\n"; + std::cerr << "No output." << std::endl; + return EXIT_FAILURE; } - // Surface Mesher options - FT sm_angle = 20.0; // theorical guaranty if angle >= 30, but slower - FT sm_radius = 0.1; // as suggested by LR - FT sm_distance = 0.002; // AG: was 0.005 (upper bound of distance to surface/Poisson) - FT sm_error_bound = 1e-3; + // Surface Mesher options + FT sm_angle = 20.0; // theorical guaranty if angle >= 30, but slower + FT sm_radius = 0.1; // as suggested by LR + FT sm_distance = 0.002; // AG: was 0.005 (upper bound of distance to surface/Poisson) + FT sm_error_bound = 1e-3; // Accumulated errors int accumulated_fatal_err = EXIT_SUCCESS; + //*************************************** // Process each input file + //*************************************** + for (int arg_index = 1; arg_index <= argc-1; arg_index++) { std::cerr << std::endl; + //*************************************** + // Load mesh/point set + //*************************************** + // File name is: std::string input_filename = argv[arg_index]; // get extension std::string extension = input_filename.substr(input_filename.find_last_of('.')); - //*************************************** - // Load mesh/point set - //*************************************** - Dt3 dt; if (extension == ".off" || extension == ".OFF") @@ -194,6 +199,8 @@ int main(int argc, char * argv[]) // Compute implicit function //*************************************** + std::cerr << "Compute implicit function...\n"; + CGAL::Timer task_timer; task_timer.start(); // Create implicit function @@ -203,21 +210,26 @@ int main(int argc, char * argv[]) /// at each vertex of the triangulation if ( ! poisson_function.compute_implicit_function() ) { - std::cerr << "FATAL ERROR: cannot solve Poisson equation" << std::endl; + std::cerr << "FATAL ERROR: cannot compute implicit function" << std::endl; accumulated_fatal_err = EXIT_FAILURE; continue; } // Print status + long memory = CGAL::Memory_sizer().virtual_size(); int nb_vertices2 = dt.number_of_vertices(); - std::cerr << "Solve Poisson equation: " - << "(added " << nb_vertices2-nb_vertices << " vertices)" - << std::endl; + std::cerr << "Compute implicit function: " << "added " << nb_vertices2-nb_vertices << " Steiner points, " + << task_timer.time() << " seconds, " + << (memory>>20) << " Mb allocated" + << std::endl; + task_timer.reset(); //*************************************** // Surface mesh generation //*************************************** + std::cerr << "Surface meshing...\n"; + STr tr; // 3D-Delaunay triangulation C2t3 c2t3 (tr); // 2D-complex in 3D-Delaunay triangulation @@ -247,18 +259,23 @@ int main(int argc, char * argv[]) CGAL::Surface_mesh_default_criteria_3 criteria(sm_angle, // lower bound of facets angles (degrees) sm_radius*size, // upper bound of Delaunay balls radii sm_distance*size); // upper bound of distance to surface - -std::cerr << "Implicit_surface_3(dichotomy error="<>20) << " Mb allocated" << std::endl; task_timer.reset(); diff --git a/Surface_reconstruction_3/test/Surface_reconstruction_3/smoothing_test.cpp b/Surface_reconstruction_3/test/Surface_reconstruction_3/smoothing_test.cpp index 882003e5e75..0d3165176b8 100644 --- a/Surface_reconstruction_3/test/Surface_reconstruction_3/smoothing_test.cpp +++ b/Surface_reconstruction_3/test/Surface_reconstruction_3/smoothing_test.cpp @@ -19,7 +19,6 @@ // CGAL #include -#include // This package #include