diff options
10 files changed, 520 insertions, 508 deletions
diff --git a/src/Alpha_complex/benchmark/Alpha_complex_3d_benchmark.cpp b/src/Alpha_complex/benchmark/Alpha_complex_3d_benchmark.cpp index 3e556cdf..1a33f2b4 100644 --- a/src/Alpha_complex/benchmark/Alpha_complex_3d_benchmark.cpp +++ b/src/Alpha_complex/benchmark/Alpha_complex_3d_benchmark.cpp @@ -18,92 +18,100 @@ std::ofstream results_csv("results.csv"); template <typename Kernel> -void benchmark_points_on_torus_dD(const std::string& msg ) { +void benchmark_points_on_torus_dD(const std::string& msg) { std::cout << "+ " << msg << std::endl; results_csv << "\"" << msg << "\";" << std::endl; - results_csv << "\"nb_points\";" << "\"nb_simplices\";" << "\"alpha_creation_time(sec.)\";" << - "\"simplex_creation_time(sec.)\";" << std::endl; + results_csv << "\"nb_points\";" + << "\"nb_simplices\";" + << "\"alpha_creation_time(sec.)\";" + << "\"simplex_creation_time(sec.)\";" << std::endl; - using K = CGAL::Epick_d< CGAL::Dimension_tag<3> >; - for (int nb_points = 1000; nb_points <= 125000 ; nb_points *= 5) { + using K = CGAL::Epick_d<CGAL::Dimension_tag<3>>; + for (int nb_points = 1000; nb_points <= 125000; nb_points *= 5) { std::cout << " Alpha complex dD on torus with " << nb_points << " points." << std::endl; std::vector<K::Point_d> points_on_torus = Gudhi::generate_points_on_torus_3D<K>(nb_points, 1.0, 0.5); std::vector<typename Kernel::Point_d> points; - for(auto p:points_on_torus) { - points.push_back(typename Kernel::Point_d(p.begin(),p.end())); + for (auto p : points_on_torus) { + points.push_back(typename Kernel::Point_d(p.begin(), p.end())); } Gudhi::Clock ac_create_clock(" benchmark_points_on_torus_dD - Alpha complex 3d creation"); ac_create_clock.begin(); Gudhi::alpha_complex::Alpha_complex<Kernel> alpha_complex_from_points(points); - ac_create_clock.end(); std::cout << ac_create_clock; + ac_create_clock.end(); + std::cout << ac_create_clock; Gudhi::Simplex_tree<> simplex; Gudhi::Clock st_create_clock(" benchmark_points_on_torus_dD - simplex tree creation"); st_create_clock.begin(); alpha_complex_from_points.create_complex(simplex); - st_create_clock.end(); std::cout << st_create_clock; + st_create_clock.end(); + std::cout << st_create_clock; - results_csv << nb_points << ";" << simplex.num_simplices() << ";" << ac_create_clock.num_seconds() << ";" << - st_create_clock.num_seconds() << ";" << std::endl; + results_csv << nb_points << ";" << simplex.num_simplices() << ";" << ac_create_clock.num_seconds() << ";" + << st_create_clock.num_seconds() << ";" << std::endl; std::cout << " benchmark_points_on_torus_dD - nb simplices = " << simplex.num_simplices() << std::endl; } - } template <typename Alpha_complex_3d> -void benchmark_points_on_torus_3D(const std::string& msg ) { - using K = CGAL::Epick_d< CGAL::Dimension_tag<3> >; +void benchmark_points_on_torus_3D(const std::string& msg) { + using K = CGAL::Epick_d<CGAL::Dimension_tag<3>>; std::cout << "+ " << msg << std::endl; results_csv << "\"" << msg << "\";" << std::endl; - results_csv << "\"nb_points\";" << "\"nb_simplices\";" << "\"alpha_creation_time(sec.)\";" << - "\"simplex_creation_time(sec.)\";" << std::endl; + results_csv << "\"nb_points\";" + << "\"nb_simplices\";" + << "\"alpha_creation_time(sec.)\";" + << "\"simplex_creation_time(sec.)\";" << std::endl; - for (int nb_points = 1000; nb_points <= 125000 ; nb_points *= 5) { + for (int nb_points = 1000; nb_points <= 125000; nb_points *= 5) { std::cout << " Alpha complex 3d on torus with " << nb_points << " points." << std::endl; std::vector<K::Point_d> points_on_torus = Gudhi::generate_points_on_torus_3D<K>(nb_points, 1.0, 0.5); std::vector<typename Alpha_complex_3d::Point_3> points; - for(auto p:points_on_torus) { - points.push_back(typename Alpha_complex_3d::Point_3(p[0],p[1],p[2])); + for (auto p : points_on_torus) { + points.push_back(typename Alpha_complex_3d::Point_3(p[0], p[1], p[2])); } Gudhi::Clock ac_create_clock(" benchmark_points_on_torus_3D - Alpha complex 3d creation"); ac_create_clock.begin(); Alpha_complex_3d alpha_complex_from_points(points); - ac_create_clock.end(); std::cout << ac_create_clock; + ac_create_clock.end(); + std::cout << ac_create_clock; Gudhi::Simplex_tree<> simplex; Gudhi::Clock st_create_clock(" benchmark_points_on_torus_3D - simplex tree creation"); st_create_clock.begin(); alpha_complex_from_points.create_complex(simplex); - st_create_clock.end(); std::cout << st_create_clock; + st_create_clock.end(); + std::cout << st_create_clock; - results_csv << nb_points << ";" << simplex.num_simplices() << ";" << ac_create_clock.num_seconds() << ";" << - st_create_clock.num_seconds() << ";" << std::endl; + results_csv << nb_points << ";" << simplex.num_simplices() << ";" << ac_create_clock.num_seconds() << ";" + << st_create_clock.num_seconds() << ";" << std::endl; std::cout << " benchmark_points_on_torus_3D - nb simplices = " << simplex.num_simplices() << std::endl; } - } template <typename Weighted_alpha_complex_3d> -void benchmark_weighted_points_on_torus_3D(const std::string& msg ) { - using K = CGAL::Epick_d< CGAL::Dimension_tag<3> >; +void benchmark_weighted_points_on_torus_3D(const std::string& msg) { + using K = CGAL::Epick_d<CGAL::Dimension_tag<3>>; std::cout << "+ " << msg << std::endl; results_csv << "\"" << msg << "\";" << std::endl; - results_csv << "\"nb_points\";" << "\"nb_simplices\";" << "\"alpha_creation_time(sec.)\";" << - "\"simplex_creation_time(sec.)\";" << std::endl; + results_csv << "\"nb_points\";" + << "\"nb_simplices\";" + << "\"alpha_creation_time(sec.)\";" + << "\"simplex_creation_time(sec.)\";" << std::endl; CGAL::Random random(8); - for (int nb_points = 1000; nb_points <= 125000 ; nb_points *= 5) { + for (int nb_points = 1000; nb_points <= 125000; nb_points *= 5) { std::cout << " Alpha complex 3d on torus with " << nb_points << " points." << std::endl; std::vector<K::Point_d> points_on_torus = Gudhi::generate_points_on_torus_3D<K>(nb_points, 1.0, 0.5); @@ -112,50 +120,52 @@ void benchmark_weighted_points_on_torus_3D(const std::string& msg ) { std::vector<Weighted_point> points; - for(auto p:points_on_torus) { - points.push_back(Weighted_point(Point(p[0],p[1],p[2]), 0.9 + random.get_double(0., 0.01))); + for (auto p : points_on_torus) { + points.push_back(Weighted_point(Point(p[0], p[1], p[2]), 0.9 + random.get_double(0., 0.01))); } Gudhi::Clock ac_create_clock(" benchmark_weighted_points_on_torus_3D - Alpha complex 3d creation"); ac_create_clock.begin(); Weighted_alpha_complex_3d alpha_complex_from_points(points); - ac_create_clock.end(); std::cout << ac_create_clock; + ac_create_clock.end(); + std::cout << ac_create_clock; Gudhi::Simplex_tree<> simplex; Gudhi::Clock st_create_clock(" benchmark_weighted_points_on_torus_3D - simplex tree creation"); st_create_clock.begin(); alpha_complex_from_points.create_complex(simplex); - st_create_clock.end(); std::cout << st_create_clock; + st_create_clock.end(); + std::cout << st_create_clock; - results_csv << nb_points << ";" << simplex.num_simplices() << ";" << ac_create_clock.num_seconds() << ";" << - st_create_clock.num_seconds() << ";" << std::endl; + results_csv << nb_points << ";" << simplex.num_simplices() << ";" << ac_create_clock.num_seconds() << ";" + << st_create_clock.num_seconds() << ";" << std::endl; std::cout << " benchmark_weighted_points_on_torus_3D - nb simplices = " << simplex.num_simplices() << std::endl; } - } template <typename Periodic_alpha_complex_3d> -void benchmark_periodic_points(const std::string& msg ) { +void benchmark_periodic_points(const std::string& msg) { std::cout << "+ " << msg << std::endl; results_csv << "\"" << msg << "\";" << std::endl; - results_csv << "\"nb_points\";" << "\"nb_simplices\";" << "\"alpha_creation_time(sec.)\";" << - "\"simplex_creation_time(sec.)\";" << std::endl; + results_csv << "\"nb_points\";" + << "\"nb_simplices\";" + << "\"alpha_creation_time(sec.)\";" + << "\"simplex_creation_time(sec.)\";" << std::endl; CGAL::Random random(8); - for (double nb_points = 10.; nb_points <= 40. ; nb_points += 10.) { - std::cout << " Periodic alpha complex 3d with " << nb_points*nb_points*nb_points << " points." << std::endl; + for (double nb_points = 10.; nb_points <= 40.; nb_points += 10.) { + std::cout << " Periodic alpha complex 3d with " << nb_points * nb_points * nb_points << " points." << std::endl; using Point = typename Periodic_alpha_complex_3d::Point_3; std::vector<Point> points; for (double i = 0; i < nb_points; i++) { for (double j = 0; j < nb_points; j++) { for (double k = 0; k < nb_points; k++) { - points.push_back(Point(i + random.get_double(0., 0.1), - j + random.get_double(0., 0.1), - k + random.get_double(0., 0.1))); + points.push_back( + Point(i + random.get_double(0., 0.1), j + random.get_double(0., 0.1), k + random.get_double(0., 0.1))); } } } @@ -163,34 +173,38 @@ void benchmark_periodic_points(const std::string& msg ) { Gudhi::Clock ac_create_clock(" benchmark_periodic_points - Alpha complex 3d creation"); ac_create_clock.begin(); Periodic_alpha_complex_3d alpha_complex_from_points(points, 0., 0., 0., nb_points, nb_points, nb_points); - ac_create_clock.end(); std::cout << ac_create_clock; + ac_create_clock.end(); + std::cout << ac_create_clock; Gudhi::Simplex_tree<> simplex; Gudhi::Clock st_create_clock(" benchmark_periodic_points - simplex tree creation"); st_create_clock.begin(); alpha_complex_from_points.create_complex(simplex); - st_create_clock.end(); std::cout << st_create_clock; + st_create_clock.end(); + std::cout << st_create_clock; - results_csv << nb_points*nb_points*nb_points << ";" << simplex.num_simplices() << ";" << - ac_create_clock.num_seconds() << ";" << st_create_clock.num_seconds() << ";" << std::endl; + results_csv << nb_points * nb_points * nb_points << ";" << simplex.num_simplices() << ";" + << ac_create_clock.num_seconds() << ";" << st_create_clock.num_seconds() << ";" << std::endl; std::cout << " benchmark_periodic_points - nb simplices = " << simplex.num_simplices() << std::endl; } - } template <typename Weighted_periodic_alpha_complex_3d> -void benchmark_weighted_periodic_points(const std::string& msg ) { +void benchmark_weighted_periodic_points(const std::string& msg) { std::cout << "+ " << msg << std::endl; results_csv << "\"" << msg << "\";" << std::endl; - results_csv << "\"nb_points\";" << "\"nb_simplices\";" << "\"alpha_creation_time(sec.)\";" << - "\"simplex_creation_time(sec.)\";" << std::endl; + results_csv << "\"nb_points\";" + << "\"nb_simplices\";" + << "\"alpha_creation_time(sec.)\";" + << "\"simplex_creation_time(sec.)\";" << std::endl; CGAL::Random random(8); - for (double nb_points = 10.; nb_points <= 40. ; nb_points += 10.) { - std::cout << " Weighted periodic alpha complex 3d with " << nb_points*nb_points*nb_points << " points." << std::endl; + for (double nb_points = 10.; nb_points <= 40.; nb_points += 10.) { + std::cout << " Weighted periodic alpha complex 3d with " << nb_points * nb_points * nb_points << " points." + << std::endl; using Point = typename Weighted_periodic_alpha_complex_3d::Point_3; using Weighted_point = typename Weighted_periodic_alpha_complex_3d::Triangulation_3::Weighted_point; @@ -199,10 +213,9 @@ void benchmark_weighted_periodic_points(const std::string& msg ) { for (double i = 0; i < nb_points; i++) { for (double j = 0; j < nb_points; j++) { for (double k = 0; k < nb_points; k++) { - points.push_back(Weighted_point(Point(i + random.get_double(0., 0.1), - j + random.get_double(0., 0.1), - k + random.get_double(0., 0.1)), - random.get_double(0., (nb_points*nb_points)/64.))); + points.push_back(Weighted_point( + Point(i + random.get_double(0., 0.1), j + random.get_double(0., 0.1), k + random.get_double(0., 0.1)), + random.get_double(0., (nb_points * nb_points) / 64.))); } } } @@ -210,57 +223,56 @@ void benchmark_weighted_periodic_points(const std::string& msg ) { Gudhi::Clock ac_create_clock(" benchmark_weighted_periodic_points - Alpha complex 3d creation"); ac_create_clock.begin(); Weighted_periodic_alpha_complex_3d alpha_complex_from_points(points, 0., 0., 0., nb_points, nb_points, nb_points); - ac_create_clock.end(); std::cout << ac_create_clock; + ac_create_clock.end(); + std::cout << ac_create_clock; Gudhi::Simplex_tree<> simplex; Gudhi::Clock st_create_clock(" benchmark_weighted_periodic_points - simplex tree creation"); st_create_clock.begin(); alpha_complex_from_points.create_complex(simplex); - st_create_clock.end(); std::cout << st_create_clock; + st_create_clock.end(); + std::cout << st_create_clock; - results_csv << nb_points*nb_points*nb_points << ";" << simplex.num_simplices() << ";" << - ac_create_clock.num_seconds() << ";" << st_create_clock.num_seconds() << ";" << std::endl; + results_csv << nb_points * nb_points * nb_points << ";" << simplex.num_simplices() << ";" + << ac_create_clock.num_seconds() << ";" << st_create_clock.num_seconds() << ";" << std::endl; std::cout << " benchmark_weighted_periodic_points - nb simplices = " << simplex.num_simplices() << std::endl; } - } -int main(int argc, char **argv) { - - benchmark_points_on_torus_dD<CGAL::Epick_d< CGAL::Dimension_tag<3> >>("Fast static dimension version"); - benchmark_points_on_torus_dD<CGAL::Epick_d< CGAL::Dynamic_dimension_tag >>("Fast dynamic dimension version"); - benchmark_points_on_torus_dD<CGAL::Epeck_d< CGAL::Dimension_tag<3> >>("Exact static dimension version"); - benchmark_points_on_torus_dD<CGAL::Epeck_d< CGAL::Dynamic_dimension_tag >>("Exact dynamic dimension version"); - - benchmark_points_on_torus_3D<Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::FAST, - false, false>>("Fast version"); - benchmark_points_on_torus_3D<Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::SAFE, - false, false>>("Safe version"); - benchmark_points_on_torus_3D<Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::EXACT, - false, false>>("Exact version"); - - - benchmark_weighted_points_on_torus_3D<Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::FAST, - true, false>>("Fast version"); - benchmark_weighted_points_on_torus_3D<Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::SAFE, - true, false>>("Safe version"); - benchmark_weighted_points_on_torus_3D<Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::EXACT, - true, false>>("Exact version"); - - benchmark_periodic_points<Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::FAST, - false, true>>("Fast version"); - benchmark_periodic_points<Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::SAFE, - false, true>>("Safe version"); - benchmark_periodic_points<Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::EXACT, - false, true>>("Exact version"); - - benchmark_weighted_periodic_points<Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::FAST, - true, true>>("Fast version"); - benchmark_weighted_periodic_points<Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::SAFE, - true, true>>("Safe version"); - benchmark_weighted_periodic_points<Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::EXACT, - true, true>>("Exact version"); +int main(int argc, char** argv) { + benchmark_points_on_torus_dD<CGAL::Epick_d<CGAL::Dimension_tag<3>>>("Fast static dimension version"); + benchmark_points_on_torus_dD<CGAL::Epick_d<CGAL::Dynamic_dimension_tag>>("Fast dynamic dimension version"); + benchmark_points_on_torus_dD<CGAL::Epeck_d<CGAL::Dimension_tag<3>>>("Exact static dimension version"); + benchmark_points_on_torus_dD<CGAL::Epeck_d<CGAL::Dynamic_dimension_tag>>("Exact dynamic dimension version"); + + benchmark_points_on_torus_3D< + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::FAST, false, false>>("Fast version"); + benchmark_points_on_torus_3D< + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::SAFE, false, false>>("Safe version"); + benchmark_points_on_torus_3D< + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::EXACT, false, false>>("Exact version"); + + benchmark_weighted_points_on_torus_3D< + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::FAST, true, false>>("Fast version"); + benchmark_weighted_points_on_torus_3D< + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::SAFE, true, false>>("Safe version"); + benchmark_weighted_points_on_torus_3D< + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::EXACT, true, false>>("Exact version"); + + benchmark_periodic_points< + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::FAST, false, true>>("Fast version"); + benchmark_periodic_points< + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::SAFE, false, true>>("Safe version"); + benchmark_periodic_points< + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::EXACT, false, true>>("Exact version"); + + benchmark_weighted_periodic_points< + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::FAST, true, true>>("Fast version"); + benchmark_weighted_periodic_points< + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::SAFE, true, true>>("Safe version"); + benchmark_weighted_periodic_points< + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::EXACT, true, true>>("Exact version"); return 0; } diff --git a/src/Alpha_complex/concept/SimplicialComplexForAlpha3d.h b/src/Alpha_complex/concept/SimplicialComplexForAlpha3d.h index f6085a26..7acdf105 100644 --- a/src/Alpha_complex/concept/SimplicialComplexForAlpha3d.h +++ b/src/Alpha_complex/concept/SimplicialComplexForAlpha3d.h @@ -1,5 +1,5 @@ -/* This file is part of the Gudhi Library. The Gudhi library - * (Geometric Understanding in Higher Dimensions) is a generic C++ +/* This file is part of the Gudhi Library. The Gudhi library + * (Geometric Understanding in Higher Dimensions) is a generic C++ * library for computational topology. * * Author(s): Vincent Rouvreau @@ -41,14 +41,13 @@ struct SimplicialComplexForAlpha3d { /** \brief Inserts a simplex from a given simplex (represented by a vector of Vertex_handle) in the * simplicial complex with the given 'filtration' value. */ - void insert_simplex(std::vector<Vertex_handle> const & vertex_range, Filtration_value filtration); + void insert_simplex(std::vector<Vertex_handle> const& vertex_range, Filtration_value filtration); /** Browses the simplicial complex to make the filtration non-decreasing. */ void make_filtration_non_decreasing(); /** Prune the simplicial complex above 'filtration' value given as parameter. */ void prune_above_filtration(Filtration_value filtration); - }; } // namespace alpha_complex diff --git a/src/Alpha_complex/doc/Intro_alpha_complex.h b/src/Alpha_complex/doc/Intro_alpha_complex.h index 82aee275..648fb6d6 100644 --- a/src/Alpha_complex/doc/Intro_alpha_complex.h +++ b/src/Alpha_complex/doc/Intro_alpha_complex.h @@ -29,34 +29,34 @@ namespace Gudhi { namespace alpha_complex { /** \defgroup alpha_complex Alpha complex - * + * * \author Vincent Rouvreau * * @{ - * + * * \section definition Definition - * + * * Alpha_complex is a <a target="_blank" href="https://en.wikipedia.org/wiki/Simplicial_complex">simplicial complex</a> * constructed from the finite cells of a Delaunay Triangulation. - * + * * The filtration value of each simplex is computed as the square of the circumradius of the simplex if the * circumsphere is empty (the simplex is then said to be Gabriel), and as the minimum of the filtration * values of the codimension 1 cofaces that make it not Gabriel otherwise. - * + * * All simplices that have a filtration value strictly greater than a given alpha squared value are not inserted into * the complex. - * + * * \image html "alpha_complex_representation.png" "Alpha-complex representation" - * + * * Alpha_complex is constructing a <a target="_blank" * href="http://doc.cgal.org/latest/Triangulation/index.html#Chapter_Triangulations">Delaunay Triangulation</a> * \cite cgal:hdj-t-15b from <a target="_blank" href="http://www.cgal.org/">CGAL</a> (the Computational Geometry * Algorithms Library \cite cgal:eb-15b) and is able to create a `SimplicialComplexForAlpha`. - * + * * The complex is a template class requiring an Epick_d <a target="_blank" * href="http://doc.cgal.org/latest/Kernel_d/index.html#Chapter_dD_Geometry_Kernel">dD Geometry Kernel</a> * \cite cgal:s-gkd-15b from CGAL as template parameter. - * + * * \remark * - When the simplicial complex is constructed with an infinite value of alpha, the complex is a Delaunay * complex. @@ -65,30 +65,30 @@ namespace alpha_complex { * \ref cech_complex can still make sense in higher dimension precisely because you can bound the radii. * * \section pointsexample Example from points - * + * * This example builds the Delaunay triangulation from the given points in a 2D static kernel, and creates a * `Simplex_tree` with it. - * + * * Then, it is asked to display information about the simplicial complex. - * + * * \include Alpha_complex/Alpha_complex_from_points.cpp - * + * * When launching: - * + * * \code $> ./Alpha_complex_example_from_points * \endcode * * the program output is: - * + * * \include Alpha_complex/alphaoffreader_for_doc_60.txt - * + * * \section createcomplexalgorithm Create complex algorithm - * + * * \subsection datastructure Data structure - * + * * In order to create the simplicial complex, first, it is built from the cells of the Delaunay Triangulation. * The filtration values are set to NaN, which stands for unknown value. - * + * * In example, : * \image html "alpha_complex_doc.png" "Simplicial complex structure construction example" * @@ -118,53 +118,53 @@ namespace alpha_complex { * \f$ * * \subsubsection dimension2 Dimension 2 - * + * * From the example above, it means the algorithm looks into each triangle ([0,1,2], [0,2,4], [1,2,3], ...), * computes the filtration value of the triangle, and then propagates the filtration value as described * here : * \image html "alpha_complex_doc_420.png" "Filtration value propagation example" - * + * * \subsubsection dimension1 Dimension 1 - * + * * Then, the algorithm looks into each edge ([0,1], [0,2], [1,2], ...), * computes the filtration value of the edge (in this case, propagation will have no effect). - * + * * \subsubsection dimension0 Dimension 0 - * + * * Finally, the algorithm looks into each vertex ([0], [1], [2], [3], [4], [5] and [6]) and * sets the filtration value (0 in case of a vertex - propagation will have no effect). - * + * * \subsubsection nondecreasing Non decreasing filtration values - * + * * As the squared radii computed by CGAL are an approximation, it might happen that these alpha squared values do not * quite define a proper filtration (i.e. non-decreasing with respect to inclusion). * We fix that up by calling `SimplicialComplexForAlpha::make_filtration_non_decreasing()`. - * + * * \subsubsection pruneabove Prune above given filtration value - * + * * The simplex tree is pruned from the given maximum alpha squared value (cf. * `SimplicialComplexForAlpha::prune_above_filtration()`). * In the following example, the value is given by the user as argument of the program. - * - * + * + * * \section offexample Example from OFF file - * + * * This example builds the Delaunay triangulation in a dynamic kernel, and initializes the alpha complex with it. - * - * + * + * * Then, it is asked to display information about the alpha complex. - * + * * \include Alpha_complex/Alpha_complex_from_off.cpp - * + * * When launching: - * + * * \code $> ./Alpha_complex_example_from_off ../../data/points/alphacomplexdoc.off 32.0 * \endcode * * the program output is: - * + * * \include Alpha_complex/alphaoffreader_for_doc_32.txt - * + * * * \section weighted3dexample 3d specific example * diff --git a/src/Alpha_complex/example/Alpha_complex_3d_from_points.cpp b/src/Alpha_complex/example/Alpha_complex_3d_from_points.cpp index 5d6e65cb..3acebd2e 100644 --- a/src/Alpha_complex/example/Alpha_complex_3d_from_points.cpp +++ b/src/Alpha_complex/example/Alpha_complex_3d_from_points.cpp @@ -8,7 +8,7 @@ #include <limits> // for numeric limits using Alpha_complex_3d = Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::EXACT, false, false>; -using Point = Alpha_complex_3d::Point_3 ; +using Point = Alpha_complex_3d::Point_3; using Vector_of_points = std::vector<Point>; int main(int argc, char **argv) { @@ -38,9 +38,8 @@ int main(int argc, char **argv) { // ---------------------------------------------------------------------------- // Display information about the alpha complex // ---------------------------------------------------------------------------- - std::cout << "Alpha complex is of dimension " << simplex.dimension() << - " - " << simplex.num_simplices() << " simplices - " << - simplex.num_vertices() << " vertices." << std::endl; + std::cout << "Alpha complex is of dimension " << simplex.dimension() << " - " << simplex.num_simplices() + << " simplices - " << simplex.num_vertices() << " vertices." << std::endl; std::cout << "Iterator on alpha complex simplices in the filtration order, with [filtration value]:" << std::endl; for (auto f_simplex : simplex.filtration_simplex_range()) { @@ -48,7 +47,8 @@ int main(int argc, char **argv) { for (auto vertex : simplex.simplex_vertex_range(f_simplex)) { std::cout << vertex << " "; } - std::cout << ") -> " << "[" << simplex.filtration(f_simplex) << "] "; + std::cout << ") -> " + << "[" << simplex.filtration(f_simplex) << "] "; std::cout << std::endl; } } diff --git a/src/Alpha_complex/example/Weighted_alpha_complex_3d_from_points.cpp b/src/Alpha_complex/example/Weighted_alpha_complex_3d_from_points.cpp index 52c39bf3..68f72f0a 100644 --- a/src/Alpha_complex/example/Weighted_alpha_complex_3d_from_points.cpp +++ b/src/Alpha_complex/example/Weighted_alpha_complex_3d_from_points.cpp @@ -7,7 +7,8 @@ #include <vector> #include <limits> // for numeric limits -using Weighted_alpha_complex_3d = Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::EXACT, true, false>; +using Weighted_alpha_complex_3d = + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::EXACT, true, false>; using Point = Weighted_alpha_complex_3d::Point_3; using Weighted_point = Weighted_alpha_complex_3d::Triangulation_3::Weighted_point; using Vector_of_weighted_points = std::vector<Weighted_point>; @@ -27,8 +28,8 @@ int main(int argc, char **argv) { weighted_points.push_back(Weighted_point(Point(1, -1, -1), 4.)); weighted_points.push_back(Weighted_point(Point(-1, 1, -1), 4.)); weighted_points.push_back(Weighted_point(Point(-1, -1, 1), 4.)); - weighted_points.push_back(Weighted_point(Point(1, 1, 1), 4.)); - weighted_points.push_back(Weighted_point(Point(2, 2, 2), 1.)); + weighted_points.push_back(Weighted_point(Point(1, 1, 1), 4.)); + weighted_points.push_back(Weighted_point(Point(2, 2, 2), 1.)); // ---------------------------------------------------------------------------- // Init of an alpha complex from the list of points @@ -40,9 +41,8 @@ int main(int argc, char **argv) { // ---------------------------------------------------------------------------- // Display information about the alpha complex // ---------------------------------------------------------------------------- - std::cout << "Alpha complex is of dimension " << simplex.dimension() << - " - " << simplex.num_simplices() << " simplices - " << - simplex.num_vertices() << " vertices." << std::endl; + std::cout << "Alpha complex is of dimension " << simplex.dimension() << " - " << simplex.num_simplices() + << " simplices - " << simplex.num_vertices() << " vertices." << std::endl; std::cout << "Iterator on alpha complex simplices in the filtration order, with [filtration value]:" << std::endl; for (auto f_simplex : simplex.filtration_simplex_range()) { @@ -50,7 +50,8 @@ int main(int argc, char **argv) { for (auto vertex : simplex.simplex_vertex_range(f_simplex)) { std::cout << vertex << " "; } - std::cout << ") -> " << "[" << simplex.filtration(f_simplex) << "] "; + std::cout << ") -> " + << "[" << simplex.filtration(f_simplex) << "] "; std::cout << std::endl; } } diff --git a/src/Alpha_complex/include/gudhi/Alpha_complex_3d.h b/src/Alpha_complex/include/gudhi/Alpha_complex_3d.h index 1ba52ad0..0333abbd 100644 --- a/src/Alpha_complex/include/gudhi/Alpha_complex_3d.h +++ b/src/Alpha_complex/include/gudhi/Alpha_complex_3d.h @@ -55,14 +55,12 @@ #include <unordered_map> #include <stdexcept> #include <cstddef> -#include <memory> // for std::unique_ptr +#include <memory> // for std::unique_ptr #include <type_traits> // for std::conditional and std::enable_if - #if CGAL_VERSION_NR < 1041101000 - // Make compilation fail - required for external projects - https://gitlab.inria.fr/GUDHI/gudhi-devel/issues/10 - static_assert(false, - "Alpha_complex_3d is only available for CGAL >= 4.11"); +// Make compilation fail - required for external projects - https://gitlab.inria.fr/GUDHI/gudhi-devel/issues/10 +static_assert(false, "Alpha_complex_3d is only available for CGAL >= 4.11"); #endif namespace Gudhi { @@ -75,39 +73,32 @@ namespace alpha_complex { // *iterator CGAL::to_double(*iterator) CGAL::to_double(iterator->exact()) template <complexity Complexity> -struct Value_from_iterator -{ - template<typename Iterator> - static double perform(Iterator it) - { +struct Value_from_iterator { + template <typename Iterator> + static double perform(Iterator it) { // Default behaviour is to return the value pointed by the given iterator return *it; } }; template <> -struct Value_from_iterator <complexity::SAFE> -{ - template<typename Iterator> - static double perform(Iterator it) - { +struct Value_from_iterator<complexity::SAFE> { + template <typename Iterator> + static double perform(Iterator it) { // In SAFE mode, we are with Epeck or Epick with EXACT value set to CGAL::Tag_true. return CGAL::to_double(*it); } }; template <> -struct Value_from_iterator <complexity::EXACT> -{ - template<typename Iterator> - static double perform(Iterator it) - { +struct Value_from_iterator<complexity::EXACT> { + template <typename Iterator> + static double perform(Iterator it) { // In EXACT mode, we are with Epeck or Epick with EXACT value set to CGAL::Tag_true. return CGAL::to_double(it->exact()); } }; - /** * \class Alpha_complex_3d * \brief Alpha complex data structure for 3d specific case. @@ -143,7 +134,7 @@ struct Value_from_iterator <complexity::EXACT> * 3d Delaunay complex. * */ -template<complexity Complexity = complexity::FAST, bool Weighted = false, bool Periodic = false> +template <complexity Complexity = complexity::FAST, bool Weighted = false, bool Periodic = false> class Alpha_complex_3d { // Epick = Exact_predicates_inexact_constructions_kernel // Epeck = Exact_predicates_exact_constructions_kernel @@ -159,89 +150,85 @@ class Alpha_complex_3d { // // otherwise Epick + CGAL::Tag_false Epeck Epeck using Predicates = typename std::conditional<((!Weighted && !Periodic) || (Complexity == complexity::FAST)), - CGAL::Exact_predicates_inexact_constructions_kernel, - CGAL::Exact_predicates_exact_constructions_kernel>::type; + CGAL::Exact_predicates_inexact_constructions_kernel, + CGAL::Exact_predicates_exact_constructions_kernel>::type; // The other way to do a conditional type. Here there are 3 possibilities - template<typename Predicates, bool Weighted_version, bool Periodic_version> struct Kernel_3 {}; + template <typename Predicates, bool Weighted_version, bool Periodic_version> + struct Kernel_3 {}; - template < typename Predicates > + template <typename Predicates> struct Kernel_3<Predicates, false, false> { using Kernel = Predicates; }; - template < typename Predicates > + template <typename Predicates> struct Kernel_3<Predicates, true, false> { using Kernel = Predicates; }; - template < typename Predicates > + template <typename Predicates> struct Kernel_3<Predicates, false, true> { using Kernel = CGAL::Periodic_3_Delaunay_triangulation_traits_3<Predicates>; }; - template < typename Predicates > + template <typename Predicates> struct Kernel_3<Predicates, true, true> { using Kernel = CGAL::Periodic_3_regular_triangulation_traits_3<Predicates>; }; using Kernel = typename Kernel_3<Predicates, Weighted, Periodic>::Kernel; - using Exact_tag = typename std::conditional<(Complexity == complexity::FAST), - CGAL::Tag_false, - CGAL::Tag_true>::type; + using Exact_tag = typename std::conditional<(Complexity == complexity::FAST), CGAL::Tag_false, CGAL::Tag_true>::type; - using TdsVb = typename std::conditional<Periodic, - CGAL::Periodic_3_triangulation_ds_vertex_base_3<>, - CGAL::Triangulation_ds_vertex_base_3<>>::type; + using TdsVb = typename std::conditional<Periodic, CGAL::Periodic_3_triangulation_ds_vertex_base_3<>, + CGAL::Triangulation_ds_vertex_base_3<>>::type; - using Tvb = typename std::conditional<Weighted, - CGAL::Regular_triangulation_vertex_base_3<Kernel, TdsVb>, - CGAL::Triangulation_vertex_base_3<Kernel, TdsVb>>::type; + using Tvb = typename std::conditional<Weighted, CGAL::Regular_triangulation_vertex_base_3<Kernel, TdsVb>, + CGAL::Triangulation_vertex_base_3<Kernel, TdsVb>>::type; using Vb = CGAL::Alpha_shape_vertex_base_3<Kernel, Tvb, Exact_tag>; - using TdsCb = typename std::conditional<Periodic, - CGAL::Periodic_3_triangulation_ds_cell_base_3<>, - CGAL::Triangulation_ds_cell_base_3<>>::type; + using TdsCb = typename std::conditional<Periodic, CGAL::Periodic_3_triangulation_ds_cell_base_3<>, + CGAL::Triangulation_ds_cell_base_3<>>::type; - using Tcb = typename std::conditional<Weighted, - CGAL::Regular_triangulation_cell_base_3<Kernel, TdsCb>, - CGAL::Triangulation_cell_base_3<Kernel, TdsCb>>::type; + using Tcb = typename std::conditional<Weighted, CGAL::Regular_triangulation_cell_base_3<Kernel, TdsCb>, + CGAL::Triangulation_cell_base_3<Kernel, TdsCb>>::type; using Cb = CGAL::Alpha_shape_cell_base_3<Kernel, Tcb, Exact_tag>; using Tds = CGAL::Triangulation_data_structure_3<Vb, Cb>; // The other way to do a conditional type. Here there 4 possibilities, cannot use std::conditional - template<typename Kernel, typename Tds, bool Weighted_version, bool Periodic_version> struct Triangulation {}; + template <typename Kernel, typename Tds, bool Weighted_version, bool Periodic_version> + struct Triangulation {}; - template < typename Kernel, typename Tds > + template <typename Kernel, typename Tds> struct Triangulation<Kernel, Tds, false, false> { using Triangulation_3 = CGAL::Delaunay_triangulation_3<Kernel, Tds>; }; - template < typename Kernel, typename Tds > + template <typename Kernel, typename Tds> struct Triangulation<Kernel, Tds, true, false> { using Triangulation_3 = CGAL::Regular_triangulation_3<Kernel, Tds>; }; - template < typename Kernel, typename Tds > + template <typename Kernel, typename Tds> struct Triangulation<Kernel, Tds, false, true> { using Triangulation_3 = CGAL::Periodic_3_Delaunay_triangulation_3<Kernel, Tds>; }; - template < typename Kernel, typename Tds > + template <typename Kernel, typename Tds> struct Triangulation<Kernel, Tds, true, true> { using Triangulation_3 = CGAL::Periodic_3_regular_triangulation_3<Kernel, Tds>; }; -public: + public: using Triangulation_3 = typename Triangulation<Kernel, Tds, Weighted, Periodic>::Triangulation_3; using Alpha_shape_3 = CGAL::Alpha_shape_3<Triangulation_3, Exact_tag>; using Point_3 = typename Kernel::Point_3; -private: + private: using Alpha_value_type = typename Alpha_shape_3::FT; using Dispatch = - CGAL::Dispatch_output_iterator<CGAL::cpp11::tuple<CGAL::Object, Alpha_value_type>, - CGAL::cpp11::tuple<std::back_insert_iterator<std::vector<CGAL::Object> >, - std::back_insert_iterator<std::vector<Alpha_value_type> > > >; + CGAL::Dispatch_output_iterator<CGAL::cpp11::tuple<CGAL::Object, Alpha_value_type>, + CGAL::cpp11::tuple<std::back_insert_iterator<std::vector<CGAL::Object>>, + std::back_insert_iterator<std::vector<Alpha_value_type>>>>; using Cell_handle = typename Alpha_shape_3::Cell_handle; using Facet = typename Alpha_shape_3::Facet; @@ -253,46 +240,43 @@ private: using Vertex_list = std::vector<Alpha_vertex_handle>; #endif -public: + public: /** \brief Alpha_complex constructor from a list of points. - * - * @param[in] points Range of points to triangulate. Points must be in `Alpha_complex_3d::Point_3` or - * `Alpha_complex_3d::Triangulation_3::Weighted_point`. - * - * @pre Available if Alpha_complex_3d is not Periodic. - * - * The type InputPointRange must be a range for which std::begin and std::end return input iterators on a - * `Alpha_complex_3d::Point_3` or a `Alpha_complex_3d::Triangulation_3::Weighted_point`. - */ - template<typename InputPointRange > + * + * @param[in] points Range of points to triangulate. Points must be in `Alpha_complex_3d::Point_3` or + * `Alpha_complex_3d::Triangulation_3::Weighted_point`. + * + * @pre Available if Alpha_complex_3d is not Periodic. + * + * The type InputPointRange must be a range for which std::begin and std::end return input iterators on a + * `Alpha_complex_3d::Point_3` or a `Alpha_complex_3d::Triangulation_3::Weighted_point`. + */ + template <typename InputPointRange> Alpha_complex_3d(const InputPointRange& points) { - static_assert(!Periodic, - "This constructor is not available for periodic versions of Alpha_complex_3d"); + static_assert(!Periodic, "This constructor is not available for periodic versions of Alpha_complex_3d"); - alpha_shape_3_ptr_ = std::unique_ptr<Alpha_shape_3>(new Alpha_shape_3(std::begin(points), std::end(points), 0, - Alpha_shape_3::GENERAL)); + alpha_shape_3_ptr_ = std::unique_ptr<Alpha_shape_3>( + new Alpha_shape_3(std::begin(points), std::end(points), 0, Alpha_shape_3::GENERAL)); } /** \brief Alpha_complex constructor from a list of points and associated weights. - * - * @exception std::invalid_argument In debug mode, if points and weights do not have the same size. - * - * @param[in] points Range of points to triangulate. Points must be in `Alpha_complex_3d::Point_3` - * @param[in] weights Range of weights on points. Weights shall be in `Alpha_complex_3d::Alpha_shape_3::FT` - * - * @pre Available if Alpha_complex_3d is Weighted and not Periodic. - * - * The type InputPointRange must be a range for which std::begin and - * std::end return input iterators on a `Alpha_complex_3d::Point_3`. - * The type WeightRange must be a range for which std::begin and - * std::end return an input iterator on a `Alpha_complex_3d::Alpha_shape_3::FT`. - */ - template<typename InputPointRange , typename WeightRange> + * + * @exception std::invalid_argument In debug mode, if points and weights do not have the same size. + * + * @param[in] points Range of points to triangulate. Points must be in `Alpha_complex_3d::Point_3` + * @param[in] weights Range of weights on points. Weights shall be in `Alpha_complex_3d::Alpha_shape_3::FT` + * + * @pre Available if Alpha_complex_3d is Weighted and not Periodic. + * + * The type InputPointRange must be a range for which std::begin and + * std::end return input iterators on a `Alpha_complex_3d::Point_3`. + * The type WeightRange must be a range for which std::begin and + * std::end return an input iterator on a `Alpha_complex_3d::Alpha_shape_3::FT`. + */ + template <typename InputPointRange, typename WeightRange> Alpha_complex_3d(const InputPointRange& points, WeightRange weights) { - static_assert(Weighted, - "This constructor is not available for non-weighted versions of Alpha_complex_3d"); - static_assert(!Periodic, - "This constructor is not available for periodic versions of Alpha_complex_3d"); + static_assert(Weighted, "This constructor is not available for non-weighted versions of Alpha_complex_3d"); + static_assert(!Periodic, "This constructor is not available for periodic versions of Alpha_complex_3d"); GUDHI_CHECK((weights.size() == points.size()), std::invalid_argument("Points number in range different from weights range number")); @@ -306,44 +290,39 @@ public: index++; } - alpha_shape_3_ptr_ = std::unique_ptr<Alpha_shape_3>(new Alpha_shape_3(std::begin(weighted_points_3), - std::end(weighted_points_3), - 0, - Alpha_shape_3::GENERAL)); + alpha_shape_3_ptr_ = std::unique_ptr<Alpha_shape_3>( + new Alpha_shape_3(std::begin(weighted_points_3), std::end(weighted_points_3), 0, Alpha_shape_3::GENERAL)); } /** \brief Alpha_complex constructor from a list of points and an iso-cuboid coordinates. - * - * @exception std::invalid_argument In debug mode, if the size of the cuboid in every directions is not the same. - * - * @param[in] points Range of points to triangulate. Points must be in `Alpha_complex_3d::Point_3` or - * `Alpha_complex_3d::Triangulation_3::Weighted_point`. - * @param[in] x_min Iso-oriented cuboid x_min. - * @param[in] y_min Iso-oriented cuboid y_min. - * @param[in] z_min Iso-oriented cuboid z_min. - * @param[in] x_max Iso-oriented cuboid x_max. - * @param[in] y_max Iso-oriented cuboid y_max. - * @param[in] z_max Iso-oriented cuboid z_max. - * - * @pre Available if Alpha_complex_3d is Periodic. - * - * The type InputPointRange must be a range for which std::begin and std::end return input iterators on a - * `Alpha_complex_3d::Point_3` or a `Alpha_complex_3d::Triangulation_3::Weighted_point`. - * - * @note In weighted version, please check weights are greater than zero, and lower than 1/64*cuboid length - * squared. - */ - template<typename InputPointRange> - Alpha_complex_3d(const InputPointRange& points, - Alpha_value_type x_min, Alpha_value_type y_min, Alpha_value_type z_min, - Alpha_value_type x_max, Alpha_value_type y_max, Alpha_value_type z_max) { - static_assert(Periodic, - "This constructor is not available for non-periodic versions of Alpha_complex_3d"); + * + * @exception std::invalid_argument In debug mode, if the size of the cuboid in every directions is not the same. + * + * @param[in] points Range of points to triangulate. Points must be in `Alpha_complex_3d::Point_3` or + * `Alpha_complex_3d::Triangulation_3::Weighted_point`. + * @param[in] x_min Iso-oriented cuboid x_min. + * @param[in] y_min Iso-oriented cuboid y_min. + * @param[in] z_min Iso-oriented cuboid z_min. + * @param[in] x_max Iso-oriented cuboid x_max. + * @param[in] y_max Iso-oriented cuboid y_max. + * @param[in] z_max Iso-oriented cuboid z_max. + * + * @pre Available if Alpha_complex_3d is Periodic. + * + * The type InputPointRange must be a range for which std::begin and std::end return input iterators on a + * `Alpha_complex_3d::Point_3` or a `Alpha_complex_3d::Triangulation_3::Weighted_point`. + * + * @note In weighted version, please check weights are greater than zero, and lower than 1/64*cuboid length + * squared. + */ + template <typename InputPointRange> + Alpha_complex_3d(const InputPointRange& points, Alpha_value_type x_min, Alpha_value_type y_min, + Alpha_value_type z_min, Alpha_value_type x_max, Alpha_value_type y_max, Alpha_value_type z_max) { + static_assert(Periodic, "This constructor is not available for non-periodic versions of Alpha_complex_3d"); // Checking if the cuboid is the same in x,y and z direction. If not, CGAL will not process it. - GUDHI_CHECK((x_max - x_min == y_max - y_min) && - (x_max - x_min == z_max - z_min) && - (z_max - z_min == y_max - y_min), - std::invalid_argument("The size of the cuboid in every directions is not the same.")); + GUDHI_CHECK( + (x_max - x_min == y_max - y_min) && (x_max - x_min == z_max - z_min) && (z_max - z_min == y_max - y_min), + std::invalid_argument("The size of the cuboid in every directions is not the same.")); // Define the periodic cube Triangulation_3 pdt(typename Kernel::Iso_cuboid_3(x_min, y_min, z_min, x_max, y_max, z_max)); @@ -357,49 +336,44 @@ public: // alpha shape construction from points. CGAL has a strange behavior in REGULARIZED mode. This is the default mode // Maybe need to set it to GENERAL mode - alpha_shape_3_ptr_ = std::unique_ptr<Alpha_shape_3>(new Alpha_shape_3(pdt, 0, - Alpha_shape_3::GENERAL)); + alpha_shape_3_ptr_ = std::unique_ptr<Alpha_shape_3>(new Alpha_shape_3(pdt, 0, Alpha_shape_3::GENERAL)); } /** \brief Alpha_complex constructor from a list of points, associated weights and an iso-cuboid coordinates. - * - * @exception std::invalid_argument In debug mode, if points and weights do not have the same size. - * @exception std::invalid_argument In debug mode, if the size of the cuboid in every directions is not the same. - * @exception std::invalid_argument In debug mode, if a weight is negative, zero, or greater than 1/64*cuboid length - * squared. - * - * @param[in] points Range of points to triangulate. Points must be in `Alpha_complex_3d::Point_3` - * @param[in] weights Range of weights on points. Weights shall be in `Alpha_complex_3d::Alpha_shape_3::FT` - * @param[in] x_min Iso-oriented cuboid x_min. - * @param[in] y_min Iso-oriented cuboid y_min. - * @param[in] z_min Iso-oriented cuboid z_min. - * @param[in] x_max Iso-oriented cuboid x_max. - * @param[in] y_max Iso-oriented cuboid y_max. - * @param[in] z_max Iso-oriented cuboid z_max. - * - * @pre Available if Alpha_complex_3d is Weighted and Periodic. - * - * The type InputPointRange must be a range for which std::begin and - * std::end return input iterators on a `Alpha_complex_3d::Point_3`. - * The type WeightRange must be a range for which std::begin and - * std::end return an input iterator on a `Alpha_complex_3d::Alpha_shape_3::FT`. - * The type of x_min, y_min, z_min, x_max, y_max and z_max is `Alpha_complex_3d::Alpha_shape_3::FT`. - */ - template<typename InputPointRange , typename WeightRange> - Alpha_complex_3d(const InputPointRange& points, WeightRange weights, - Alpha_value_type x_min, Alpha_value_type y_min, Alpha_value_type z_min, - Alpha_value_type x_max, Alpha_value_type y_max, Alpha_value_type z_max) { - static_assert(Weighted, - "This constructor is not available for non-weighted versions of Alpha_complex_3d"); - static_assert(Periodic, - "This constructor is not available for non-periodic versions of Alpha_complex_3d"); + * + * @exception std::invalid_argument In debug mode, if points and weights do not have the same size. + * @exception std::invalid_argument In debug mode, if the size of the cuboid in every directions is not the same. + * @exception std::invalid_argument In debug mode, if a weight is negative, zero, or greater than 1/64*cuboid length + * squared. + * + * @param[in] points Range of points to triangulate. Points must be in `Alpha_complex_3d::Point_3` + * @param[in] weights Range of weights on points. Weights shall be in `Alpha_complex_3d::Alpha_shape_3::FT` + * @param[in] x_min Iso-oriented cuboid x_min. + * @param[in] y_min Iso-oriented cuboid y_min. + * @param[in] z_min Iso-oriented cuboid z_min. + * @param[in] x_max Iso-oriented cuboid x_max. + * @param[in] y_max Iso-oriented cuboid y_max. + * @param[in] z_max Iso-oriented cuboid z_max. + * + * @pre Available if Alpha_complex_3d is Weighted and Periodic. + * + * The type InputPointRange must be a range for which std::begin and + * std::end return input iterators on a `Alpha_complex_3d::Point_3`. + * The type WeightRange must be a range for which std::begin and + * std::end return an input iterator on a `Alpha_complex_3d::Alpha_shape_3::FT`. + * The type of x_min, y_min, z_min, x_max, y_max and z_max is `Alpha_complex_3d::Alpha_shape_3::FT`. + */ + template <typename InputPointRange, typename WeightRange> + Alpha_complex_3d(const InputPointRange& points, WeightRange weights, Alpha_value_type x_min, Alpha_value_type y_min, + Alpha_value_type z_min, Alpha_value_type x_max, Alpha_value_type y_max, Alpha_value_type z_max) { + static_assert(Weighted, "This constructor is not available for non-weighted versions of Alpha_complex_3d"); + static_assert(Periodic, "This constructor is not available for non-periodic versions of Alpha_complex_3d"); GUDHI_CHECK((weights.size() == points.size()), std::invalid_argument("Points number in range different from weights range number")); // Checking if the cuboid is the same in x,y and z direction. If not, CGAL will not process it. - GUDHI_CHECK((x_max - x_min == y_max - y_min) && - (x_max - x_min == z_max - z_min) && - (z_max - z_min == y_max - y_min), - std::invalid_argument("The size of the cuboid in every directions is not the same.")); + GUDHI_CHECK( + (x_max - x_min == y_max - y_min) && (x_max - x_min == z_max - z_min) && (z_max - z_min == y_max - y_min), + std::invalid_argument("The size of the cuboid in every directions is not the same.")); using Weighted_point_3 = typename Triangulation_3::Weighted_point; std::vector<Weighted_point_3> weighted_points_3; @@ -433,38 +407,36 @@ public: // alpha shape construction from points. CGAL has a strange behavior in REGULARIZED mode. This is the default mode // Maybe need to set it to GENERAL mode - alpha_shape_3_ptr_ = std::unique_ptr<Alpha_shape_3>(new Alpha_shape_3(pdt, 0, - Alpha_shape_3::GENERAL)); + alpha_shape_3_ptr_ = std::unique_ptr<Alpha_shape_3>(new Alpha_shape_3(pdt, 0, Alpha_shape_3::GENERAL)); } /** \brief Inserts all Delaunay triangulation into the simplicial complex. - * It also computes the filtration values accordingly to the \ref createcomplexalgorithm - * - * \tparam SimplicialComplexForAlpha3d must meet `SimplicialComplexForAlpha3d` concept. - * - * @param[in] complex SimplicialComplexForAlpha3d to be created. - * @param[in] max_alpha_square maximum for alpha square value. Default value is +\f$\infty\f$. - * - * @return true if creation succeeds, false otherwise. - * - * @pre The simplicial complex must be empty (no vertices) - * - * Initialization can be launched once. - * - */ - template <typename SimplicialComplexForAlpha3d, typename Filtration_value = - typename SimplicialComplexForAlpha3d::Filtration_value> - bool create_complex(SimplicialComplexForAlpha3d& complex, Filtration_value max_alpha_square = - std::numeric_limits<Filtration_value>::infinity()) { + * It also computes the filtration values accordingly to the \ref createcomplexalgorithm + * + * \tparam SimplicialComplexForAlpha3d must meet `SimplicialComplexForAlpha3d` concept. + * + * @param[in] complex SimplicialComplexForAlpha3d to be created. + * @param[in] max_alpha_square maximum for alpha square value. Default value is +\f$\infty\f$. + * + * @return true if creation succeeds, false otherwise. + * + * @pre The simplicial complex must be empty (no vertices) + * + * Initialization can be launched once. + * + */ + template <typename SimplicialComplexForAlpha3d, + typename Filtration_value = typename SimplicialComplexForAlpha3d::Filtration_value> + bool create_complex(SimplicialComplexForAlpha3d& complex, + Filtration_value max_alpha_square = std::numeric_limits<Filtration_value>::infinity()) { if (complex.num_vertices() > 0) { std::cerr << "Alpha_complex_3d create_complex - complex is not empty\n"; return false; // ----- >> } - //using Filtration_value = typename SimplicialComplexForAlpha3d::Filtration_value; + // using Filtration_value = typename SimplicialComplexForAlpha3d::Filtration_value; using Complex_vertex_handle = typename SimplicialComplexForAlpha3d::Vertex_handle; - using Alpha_shape_simplex_tree_map = std::unordered_map<Alpha_vertex_handle, - Complex_vertex_handle>; + using Alpha_shape_simplex_tree_map = std::unordered_map<Alpha_vertex_handle, Complex_vertex_handle>; using Simplex_tree_vector_vertex = std::vector<Complex_vertex_handle>; #ifdef DEBUG_TRACES @@ -483,7 +455,7 @@ public: #ifdef DEBUG_TRACES std::cout << "filtration_with_alpha_values returns : " << objects.size() << " objects" << std::endl; #endif // DEBUG_TRACES - + Alpha_shape_simplex_tree_map map_cgal_simplex_tree; using Alpha_value_iterator = typename std::vector<Alpha_value_type>::const_iterator; Alpha_value_iterator alpha_value_iterator = alpha_values.begin(); @@ -491,7 +463,7 @@ public: Vertex_list vertex_list; // Retrieve Alpha shape vertex list from object - if (const Cell_handle *cell = CGAL::object_cast<Cell_handle>(&object_iterator)) { + if (const Cell_handle* cell = CGAL::object_cast<Cell_handle>(&object_iterator)) { for (auto i = 0; i < 4; i++) { #ifdef DEBUG_TRACES std::cout << "from cell[" << i << "]=" << (*cell)->vertex(i)->point() << std::endl; @@ -501,29 +473,29 @@ public: #ifdef DEBUG_TRACES count_cells++; #endif // DEBUG_TRACES - } else if (const Facet *facet = CGAL::object_cast<Facet>(&object_iterator)) { - for (auto i = 0; i < 4; i++) { - if ((*facet).second != i) { + } else if (const Facet* facet = CGAL::object_cast<Facet>(&object_iterator)) { + for (auto i = 0; i < 4; i++) { + if ((*facet).second != i) { #ifdef DEBUG_TRACES - std::cout << "from facet=[" << i << "]" << (*facet).first->vertex(i)->point() << std::endl; + std::cout << "from facet=[" << i << "]" << (*facet).first->vertex(i)->point() << std::endl; #endif // DEBUG_TRACES - vertex_list.push_back((*facet).first->vertex(i)); - } + vertex_list.push_back((*facet).first->vertex(i)); } + } #ifdef DEBUG_TRACES count_facets++; #endif // DEBUG_TRACES - } else if (const Edge *edge = CGAL::object_cast<Edge>(&object_iterator)) { - for (auto i : {(*edge).second, (*edge).third}) { + } else if (const Edge* edge = CGAL::object_cast<Edge>(&object_iterator)) { + for (auto i : {(*edge).second, (*edge).third}) { #ifdef DEBUG_TRACES - std::cout << "from edge[" << i << "]=" << (*edge).first->vertex(i)->point() << std::endl; + std::cout << "from edge[" << i << "]=" << (*edge).first->vertex(i)->point() << std::endl; #endif // DEBUG_TRACES - vertex_list.push_back((*edge).first->vertex(i)); - } + vertex_list.push_back((*edge).first->vertex(i)); + } #ifdef DEBUG_TRACES count_edges++; #endif // DEBUG_TRACES - } else if (const Alpha_vertex_handle *vertex = CGAL::object_cast<Alpha_vertex_handle>(&object_iterator)) { + } else if (const Alpha_vertex_handle* vertex = CGAL::object_cast<Alpha_vertex_handle>(&object_iterator)) { #ifdef DEBUG_TRACES count_vertices++; std::cout << "from vertex=" << (*vertex)->point() << std::endl; @@ -577,10 +549,9 @@ public: return true; } -private: + private: // use of a unique_ptr on cgal Alpha_shape_3, as copy and default constructor is not available - no need to be freed std::unique_ptr<Alpha_shape_3> alpha_shape_3_ptr_; - }; } // namespace alpha_complex diff --git a/src/Alpha_complex/include/gudhi/Alpha_complex_options.h b/src/Alpha_complex/include/gudhi/Alpha_complex_options.h index 29eb514a..7a555fa1 100644 --- a/src/Alpha_complex/include/gudhi/Alpha_complex_options.h +++ b/src/Alpha_complex/include/gudhi/Alpha_complex_options.h @@ -23,7 +23,6 @@ #ifndef ALPHA_COMPLEX_OPTIONS_H_ #define ALPHA_COMPLEX_OPTIONS_H_ - namespace Gudhi { namespace alpha_complex { @@ -33,11 +32,10 @@ namespace alpha_complex { * * \ingroup alpha_complex */ -enum class complexity: char -{ - FAST='f', ///< Fast version. - SAFE='s', ///< Safe version. - EXACT='e', ///< Exact version. +enum class complexity : char { + FAST = 'f', ///< Fast version. + SAFE = 's', ///< Safe version. + EXACT = 'e', ///< Exact version. }; } // namespace alpha_complex diff --git a/src/Alpha_complex/test/Alpha_complex_3d_unit_test.cpp b/src/Alpha_complex/test/Alpha_complex_3d_unit_test.cpp index a32e88a6..9e071195 100644 --- a/src/Alpha_complex/test/Alpha_complex_3d_unit_test.cpp +++ b/src/Alpha_complex/test/Alpha_complex_3d_unit_test.cpp @@ -42,21 +42,33 @@ #include <CGAL/Random.h> #include <CGAL/point_generators_3.h> -using Fast_alpha_complex_3d = Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::FAST, false, false>; -using Safe_alpha_complex_3d = Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::SAFE, false, false>; -using Exact_alpha_complex_3d = Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::EXACT, false, false>; - -using Fast_weighted_alpha_complex_3d = Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::FAST, true, false>; -using Safe_weighted_alpha_complex_3d = Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::SAFE, true, false>; -using Exact_weighted_alpha_complex_3d = Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::EXACT, true, false>; - -using Fast_periodic_alpha_complex_3d = Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::FAST, false, true>; -using Safe_periodic_alpha_complex_3d = Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::SAFE, false, true>; -using Exact_periodic_alpha_complex_3d = Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::EXACT, false, true>; - -using Fast_weighted_periodic_alpha_complex_3d = Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::FAST, true, true>; -using Safe_weighted_periodic_alpha_complex_3d = Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::SAFE, true, true>; -using Exact_weighted_periodic_alpha_complex_3d = Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::EXACT, true, true>; +using Fast_alpha_complex_3d = + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::FAST, false, false>; +using Safe_alpha_complex_3d = + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::SAFE, false, false>; +using Exact_alpha_complex_3d = + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::EXACT, false, false>; + +using Fast_weighted_alpha_complex_3d = + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::FAST, true, false>; +using Safe_weighted_alpha_complex_3d = + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::SAFE, true, false>; +using Exact_weighted_alpha_complex_3d = + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::EXACT, true, false>; + +using Fast_periodic_alpha_complex_3d = + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::FAST, false, true>; +using Safe_periodic_alpha_complex_3d = + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::SAFE, false, true>; +using Exact_periodic_alpha_complex_3d = + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::EXACT, false, true>; + +using Fast_weighted_periodic_alpha_complex_3d = + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::FAST, true, true>; +using Safe_weighted_periodic_alpha_complex_3d = + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::SAFE, true, true>; +using Exact_weighted_periodic_alpha_complex_3d = + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::EXACT, true, true>; BOOST_AUTO_TEST_CASE(Alpha_complex_3d_from_points) { // ----------------- @@ -89,18 +101,18 @@ BOOST_AUTO_TEST_CASE(Alpha_complex_3d_from_points) { // --------------------- // Compare both versions // --------------------- - std::cout << "Exact Alpha complex 3d is of dimension " << exact_stree.dimension() - << " - Non exact is " << stree.dimension() << std::endl; + std::cout << "Exact Alpha complex 3d is of dimension " << exact_stree.dimension() << " - Non exact is " + << stree.dimension() << std::endl; BOOST_CHECK(exact_stree.dimension() == stree.dimension()); - std::cout << "Exact Alpha complex 3d num_simplices " << exact_stree.num_simplices() - << " - Non exact is " << stree.num_simplices() << std::endl; + std::cout << "Exact Alpha complex 3d num_simplices " << exact_stree.num_simplices() << " - Non exact is " + << stree.num_simplices() << std::endl; BOOST_CHECK(exact_stree.num_simplices() == stree.num_simplices()); - std::cout << "Exact Alpha complex 3d num_vertices " << exact_stree.num_vertices() - << " - Non exact is " << stree.num_vertices() << std::endl; + std::cout << "Exact Alpha complex 3d num_vertices " << exact_stree.num_vertices() << " - Non exact is " + << stree.num_vertices() << std::endl; BOOST_CHECK(exact_stree.num_vertices() == stree.num_vertices()); auto sh = stree.filtration_simplex_range().begin(); - while(sh != stree.filtration_simplex_range().end()) { + while (sh != stree.filtration_simplex_range().end()) { std::vector<int> simplex; std::vector<int> exact_simplex; std::cout << "Non-exact ( "; @@ -108,7 +120,8 @@ BOOST_AUTO_TEST_CASE(Alpha_complex_3d_from_points) { simplex.push_back(vertex); std::cout << vertex << " "; } - std::cout << ") -> " << "[" << stree.filtration(*sh) << "] "; + std::cout << ") -> " + << "[" << stree.filtration(*sh) << "] "; std::cout << std::endl; // Find it in the exact structure @@ -133,18 +146,18 @@ BOOST_AUTO_TEST_CASE(Alpha_complex_3d_from_points) { // --------------------- // Compare both versions // --------------------- - std::cout << "Exact Alpha complex 3d is of dimension " << safe_stree.dimension() - << " - Non exact is " << stree.dimension() << std::endl; + std::cout << "Exact Alpha complex 3d is of dimension " << safe_stree.dimension() << " - Non exact is " + << stree.dimension() << std::endl; BOOST_CHECK(safe_stree.dimension() == stree.dimension()); - std::cout << "Exact Alpha complex 3d num_simplices " << safe_stree.num_simplices() - << " - Non exact is " << stree.num_simplices() << std::endl; + std::cout << "Exact Alpha complex 3d num_simplices " << safe_stree.num_simplices() << " - Non exact is " + << stree.num_simplices() << std::endl; BOOST_CHECK(safe_stree.num_simplices() == stree.num_simplices()); - std::cout << "Exact Alpha complex 3d num_vertices " << safe_stree.num_vertices() - << " - Non exact is " << stree.num_vertices() << std::endl; + std::cout << "Exact Alpha complex 3d num_vertices " << safe_stree.num_vertices() << " - Non exact is " + << stree.num_vertices() << std::endl; BOOST_CHECK(safe_stree.num_vertices() == stree.num_vertices()); auto safe_sh = stree.filtration_simplex_range().begin(); - while(safe_sh != stree.filtration_simplex_range().end()) { + while (safe_sh != stree.filtration_simplex_range().end()) { std::vector<int> simplex; std::vector<int> exact_simplex; #ifdef DEBUG_TRACES @@ -157,7 +170,8 @@ BOOST_AUTO_TEST_CASE(Alpha_complex_3d_from_points) { #endif } #ifdef DEBUG_TRACES - std::cout << ") -> " << "[" << stree.filtration(*safe_sh) << "] "; + std::cout << ") -> " + << "[" << stree.filtration(*safe_sh) << "] "; std::cout << std::endl; #endif @@ -172,9 +186,9 @@ BOOST_AUTO_TEST_CASE(Alpha_complex_3d_from_points) { } } -typedef boost::mpl::list<Fast_weighted_alpha_complex_3d, - Safe_weighted_alpha_complex_3d, - Exact_weighted_alpha_complex_3d> weighted_variants_type_list; +typedef boost::mpl::list<Fast_weighted_alpha_complex_3d, Safe_weighted_alpha_complex_3d, + Exact_weighted_alpha_complex_3d> + weighted_variants_type_list; #ifdef GUDHI_DEBUG BOOST_AUTO_TEST_CASE_TEMPLATE(Alpha_complex_weighted_throw, Weighted_alpha_complex_3d, weighted_variants_type_list) { @@ -191,7 +205,7 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(Alpha_complex_weighted_throw, Weighted_alpha_compl std::vector<double> weights = {0.01, 0.005, 0.006, 0.01, 0.009, 0.001}; std::cout << "Check exception throw in debug mode" << std::endl; - BOOST_CHECK_THROW (Weighted_alpha_complex_3d wac(w_points, weights), std::invalid_argument); + BOOST_CHECK_THROW(Weighted_alpha_complex_3d wac(w_points, weights), std::invalid_argument); } #endif @@ -218,7 +232,7 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(Alpha_complex_weighted, Weighted_alpha_complex_3d, std::vector<Weighted_point_3> weighted_points; - for (std::size_t i=0; i < w_points.size(); i++) { + for (std::size_t i = 0; i < w_points.size(); i++) { weighted_points.push_back(Weighted_point_3(w_points[i], weights[i])); } Weighted_alpha_complex_3d alpha_complex_w_p(weighted_points); @@ -229,18 +243,18 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(Alpha_complex_weighted, Weighted_alpha_complex_3d, // --------------------- // Compare both versions // --------------------- - std::cout << "Weighted alpha complex 3d is of dimension " << stree_bis.dimension() - << " - versus " << stree.dimension() << std::endl; + std::cout << "Weighted alpha complex 3d is of dimension " << stree_bis.dimension() << " - versus " + << stree.dimension() << std::endl; BOOST_CHECK(stree_bis.dimension() == stree.dimension()); - std::cout << "Weighted alpha complex 3d num_simplices " << stree_bis.num_simplices() - << " - versus " << stree.num_simplices() << std::endl; + std::cout << "Weighted alpha complex 3d num_simplices " << stree_bis.num_simplices() << " - versus " + << stree.num_simplices() << std::endl; BOOST_CHECK(stree_bis.num_simplices() == stree.num_simplices()); - std::cout << "Weighted alpha complex 3d num_vertices " << stree_bis.num_vertices() - << " - versus " << stree.num_vertices() << std::endl; + std::cout << "Weighted alpha complex 3d num_vertices " << stree_bis.num_vertices() << " - versus " + << stree.num_vertices() << std::endl; BOOST_CHECK(stree_bis.num_vertices() == stree.num_vertices()); auto sh = stree.filtration_simplex_range().begin(); - while(sh != stree.filtration_simplex_range().end()) { + while (sh != stree.filtration_simplex_range().end()) { std::vector<int> simplex; std::vector<int> exact_simplex; #ifdef DEBUG_TRACES @@ -253,7 +267,8 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(Alpha_complex_weighted, Weighted_alpha_complex_3d, #endif } #ifdef DEBUG_TRACES - std::cout << ") -> " << "[" << stree.filtration(*sh) << "] "; + std::cout << ") -> " + << "[" << stree.filtration(*sh) << "] "; std::cout << std::endl; #endif @@ -266,13 +281,12 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(Alpha_complex_weighted, Weighted_alpha_complex_3d, ++sh; } - } #ifdef GUDHI_DEBUG -typedef boost::mpl::list<Fast_periodic_alpha_complex_3d, - Safe_periodic_alpha_complex_3d, - Exact_periodic_alpha_complex_3d> periodic_variants_type_list; +typedef boost::mpl::list<Fast_periodic_alpha_complex_3d, Safe_periodic_alpha_complex_3d, + Exact_periodic_alpha_complex_3d> + periodic_variants_type_list; BOOST_AUTO_TEST_CASE_TEMPLATE(Alpha_complex_periodic_throw, Periodic_alpha_complex_3d, periodic_variants_type_list) { std::cout << "Periodic alpha complex 3d exception throw" << std::endl; @@ -284,18 +298,18 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(Alpha_complex_periodic_throw, Periodic_alpha_compl std::cout << "Check exception throw in debug mode" << std::endl; // Check it throws an exception when the cuboid is not iso - BOOST_CHECK_THROW (Periodic_alpha_complex_3d periodic_alpha_complex(p_points, 0., 0., 0., 0.9, 1., 1.), - std::invalid_argument); - BOOST_CHECK_THROW (Periodic_alpha_complex_3d periodic_alpha_complex(p_points, 0., 0., 0., 1., 0.9, 1.), - std::invalid_argument); - BOOST_CHECK_THROW (Periodic_alpha_complex_3d periodic_alpha_complex(p_points, 0., 0., 0., 1., 1., 0.9), - std::invalid_argument); - BOOST_CHECK_THROW (Periodic_alpha_complex_3d periodic_alpha_complex(p_points, 0., 0., 0., 1.1, 1., 1.), - std::invalid_argument); - BOOST_CHECK_THROW (Periodic_alpha_complex_3d periodic_alpha_complex(p_points, 0., 0., 0., 1., 1.1, 1.), - std::invalid_argument); - BOOST_CHECK_THROW (Periodic_alpha_complex_3d periodic_alpha_complex(p_points, 0., 0., 0., 1., 1., 1.1), - std::invalid_argument); + BOOST_CHECK_THROW(Periodic_alpha_complex_3d periodic_alpha_complex(p_points, 0., 0., 0., 0.9, 1., 1.), + std::invalid_argument); + BOOST_CHECK_THROW(Periodic_alpha_complex_3d periodic_alpha_complex(p_points, 0., 0., 0., 1., 0.9, 1.), + std::invalid_argument); + BOOST_CHECK_THROW(Periodic_alpha_complex_3d periodic_alpha_complex(p_points, 0., 0., 0., 1., 1., 0.9), + std::invalid_argument); + BOOST_CHECK_THROW(Periodic_alpha_complex_3d periodic_alpha_complex(p_points, 0., 0., 0., 1.1, 1., 1.), + std::invalid_argument); + BOOST_CHECK_THROW(Periodic_alpha_complex_3d periodic_alpha_complex(p_points, 0., 0., 0., 1., 1.1, 1.), + std::invalid_argument); + BOOST_CHECK_THROW(Periodic_alpha_complex_3d periodic_alpha_complex(p_points, 0., 0., 0., 1., 1., 1.1), + std::invalid_argument); } #endif @@ -310,7 +324,7 @@ BOOST_AUTO_TEST_CASE(Alpha_complex_periodic) { CGAL::Random_points_in_cube_3<Fast_periodic_alpha_complex_3d::Point_3, Creator> in_cube(1, random); std::vector<Fast_periodic_alpha_complex_3d::Point_3> p_points; - for (int i=0 ; i < 50 ; i++) { + for (int i = 0; i < 50; i++) { Fast_periodic_alpha_complex_3d::Point_3 p = *in_cube++; p_points.push_back(p); } @@ -327,7 +341,7 @@ BOOST_AUTO_TEST_CASE(Alpha_complex_periodic) { std::vector<Exact_periodic_alpha_complex_3d::Point_3> e_p_points; - for (auto p: p_points) { + for (auto p : p_points) { e_p_points.push_back(Exact_periodic_alpha_complex_3d::Point_3(p[0], p[1], p[2])); } @@ -339,27 +353,25 @@ BOOST_AUTO_TEST_CASE(Alpha_complex_periodic) { // --------------------- // Compare both versions // --------------------- - std::cout << "Exact periodic alpha complex 3d is of dimension " << exact_stree.dimension() - << " - Non exact is " << stree.dimension() << std::endl; + std::cout << "Exact periodic alpha complex 3d is of dimension " << exact_stree.dimension() << " - Non exact is " + << stree.dimension() << std::endl; BOOST_CHECK(exact_stree.dimension() == stree.dimension()); - std::cout << "Exact periodic alpha complex 3d num_simplices " << exact_stree.num_simplices() - << " - Non exact is " << stree.num_simplices() << std::endl; + std::cout << "Exact periodic alpha complex 3d num_simplices " << exact_stree.num_simplices() << " - Non exact is " + << stree.num_simplices() << std::endl; BOOST_CHECK(exact_stree.num_simplices() == stree.num_simplices()); - std::cout << "Exact periodic alpha complex 3d num_vertices " << exact_stree.num_vertices() - << " - Non exact is " << stree.num_vertices() << std::endl; + std::cout << "Exact periodic alpha complex 3d num_vertices " << exact_stree.num_vertices() << " - Non exact is " + << stree.num_vertices() << std::endl; BOOST_CHECK(exact_stree.num_vertices() == stree.num_vertices()); - // We cannot compare as objects from dispatcher on the alpha shape is not deterministic. // cf. https://github.com/CGAL/cgal/issues/3346 auto sh = stree.filtration_simplex_range().begin(); auto sh_exact = exact_stree.filtration_simplex_range().begin(); - while(sh != stree.filtration_simplex_range().end() || sh_exact != exact_stree.filtration_simplex_range().end()) { + while (sh != stree.filtration_simplex_range().end() || sh_exact != exact_stree.filtration_simplex_range().end()) { GUDHI_TEST_FLOAT_EQUALITY_CHECK(stree.filtration(*sh), exact_stree.filtration(*sh_exact), 1e-14); - std::vector<int> vh(stree.simplex_vertex_range(*sh).begin(), - stree.simplex_vertex_range(*sh).end()); + std::vector<int> vh(stree.simplex_vertex_range(*sh).begin(), stree.simplex_vertex_range(*sh).end()); std::vector<int> exact_vh(exact_stree.simplex_vertex_range(*sh_exact).begin(), exact_stree.simplex_vertex_range(*sh_exact).end()); @@ -371,7 +383,6 @@ BOOST_AUTO_TEST_CASE(Alpha_complex_periodic) { BOOST_CHECK(sh == stree.filtration_simplex_range().end()); BOOST_CHECK(sh_exact == exact_stree.filtration_simplex_range().end()); - // ---------------------- // Safe periodic version // ---------------------- @@ -379,7 +390,7 @@ BOOST_AUTO_TEST_CASE(Alpha_complex_periodic) { std::vector<Safe_periodic_alpha_complex_3d::Point_3> s_p_points; - for (auto p: p_points) { + for (auto p : p_points) { s_p_points.push_back(Safe_periodic_alpha_complex_3d::Point_3(p[0], p[1], p[2])); } @@ -396,11 +407,10 @@ BOOST_AUTO_TEST_CASE(Alpha_complex_periodic) { sh = stree.filtration_simplex_range().begin(); auto sh_safe = safe_stree.filtration_simplex_range().begin(); - while(sh != stree.filtration_simplex_range().end() || sh_safe != safe_stree.filtration_simplex_range().end()) { + while (sh != stree.filtration_simplex_range().end() || sh_safe != safe_stree.filtration_simplex_range().end()) { GUDHI_TEST_FLOAT_EQUALITY_CHECK(stree.filtration(*sh), safe_stree.filtration(*sh_safe), 1e-14); - std::vector<int> vh(stree.simplex_vertex_range(*sh).begin(), - stree.simplex_vertex_range(*sh).end()); + std::vector<int> vh(stree.simplex_vertex_range(*sh).begin(), stree.simplex_vertex_range(*sh).end()); std::vector<int> safe_vh(safe_stree.simplex_vertex_range(*sh_safe).begin(), safe_stree.simplex_vertex_range(*sh_safe).end()); @@ -411,11 +421,11 @@ BOOST_AUTO_TEST_CASE(Alpha_complex_periodic) { BOOST_CHECK(sh == stree.filtration_simplex_range().end()); BOOST_CHECK(sh_safe == safe_stree.filtration_simplex_range().end()); - } typedef boost::mpl::list<Fast_weighted_periodic_alpha_complex_3d, Exact_weighted_periodic_alpha_complex_3d, - Safe_weighted_periodic_alpha_complex_3d> wp_variants_type_list; + Safe_weighted_periodic_alpha_complex_3d> + wp_variants_type_list; #ifdef GUDHI_DEBUG BOOST_AUTO_TEST_CASE_TEMPLATE(Alpha_complex_weighted_periodic_throw, Weighted_periodic_alpha_complex_3d, @@ -427,7 +437,7 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(Alpha_complex_weighted_periodic_throw, Weighted_pe CGAL::Random_points_in_cube_3<Weighted_periodic_alpha_complex_3d::Point_3, Creator> in_cube(1, random); std::vector<Weighted_periodic_alpha_complex_3d::Point_3> wp_points; - for (int i=0 ; i < 50 ; i++) { + for (int i = 0; i < 50; i++) { Weighted_periodic_alpha_complex_3d::Point_3 p = *in_cube++; wp_points.push_back(p); } @@ -439,44 +449,50 @@ BOOST_AUTO_TEST_CASE_TEMPLATE(Alpha_complex_weighted_periodic_throw, Weighted_pe std::cout << "Cuboid is not iso exception" << std::endl; // Check it throws an exception when the cuboid is not iso - BOOST_CHECK_THROW (Weighted_periodic_alpha_complex_3d wp_alpha_complex(wp_points, p_weights, -1., -1., -1., 0.9, 1., 1.), - std::invalid_argument); - BOOST_CHECK_THROW (Weighted_periodic_alpha_complex_3d wp_alpha_complex(wp_points, p_weights, -1., -1., -1., 1., 0.9, 1.), - std::invalid_argument); - BOOST_CHECK_THROW (Weighted_periodic_alpha_complex_3d wp_alpha_complex(wp_points, p_weights, -1., -1., -1., 1., 1., 0.9), - std::invalid_argument); - BOOST_CHECK_THROW (Weighted_periodic_alpha_complex_3d wp_alpha_complex(wp_points, p_weights, -1., -1., -1., 1.1, 1., 1.), - std::invalid_argument); - BOOST_CHECK_THROW (Weighted_periodic_alpha_complex_3d wp_alpha_complex(wp_points, p_weights, -1., -1., -1., 1., 1.1, 1.), - std::invalid_argument); - BOOST_CHECK_THROW (Weighted_periodic_alpha_complex_3d wp_alpha_complex(wp_points, p_weights, -1., -1., -1., 1., 1., 1.1), - std::invalid_argument); + BOOST_CHECK_THROW( + Weighted_periodic_alpha_complex_3d wp_alpha_complex(wp_points, p_weights, -1., -1., -1., 0.9, 1., 1.), + std::invalid_argument); + BOOST_CHECK_THROW( + Weighted_periodic_alpha_complex_3d wp_alpha_complex(wp_points, p_weights, -1., -1., -1., 1., 0.9, 1.), + std::invalid_argument); + BOOST_CHECK_THROW( + Weighted_periodic_alpha_complex_3d wp_alpha_complex(wp_points, p_weights, -1., -1., -1., 1., 1., 0.9), + std::invalid_argument); + BOOST_CHECK_THROW( + Weighted_periodic_alpha_complex_3d wp_alpha_complex(wp_points, p_weights, -1., -1., -1., 1.1, 1., 1.), + std::invalid_argument); + BOOST_CHECK_THROW( + Weighted_periodic_alpha_complex_3d wp_alpha_complex(wp_points, p_weights, -1., -1., -1., 1., 1.1, 1.), + std::invalid_argument); + BOOST_CHECK_THROW( + Weighted_periodic_alpha_complex_3d wp_alpha_complex(wp_points, p_weights, -1., -1., -1., 1., 1., 1.1), + std::invalid_argument); std::cout << "0 <= point.weight() < 1/64 * domain_size * domain_size exception" << std::endl; // Weights must be in range ]0, 1/64 = 0.015625[ double temp = p_weights[25]; p_weights[25] = 1.0; - BOOST_CHECK_THROW (Weighted_periodic_alpha_complex_3d wp_alpha_complex(wp_points, p_weights, 0., 0., 0., 1., 1., 1.), - std::invalid_argument); + BOOST_CHECK_THROW(Weighted_periodic_alpha_complex_3d wp_alpha_complex(wp_points, p_weights, 0., 0., 0., 1., 1., 1.), + std::invalid_argument); // Weights must be in range ]0, 1/64 = 0.015625[ p_weights[25] = temp; temp = p_weights[14]; p_weights[14] = -1e-10; - BOOST_CHECK_THROW (Weighted_periodic_alpha_complex_3d wp_alpha_complex(wp_points, p_weights, 0., 0., 0., 1., 1., 1.), - std::invalid_argument); + BOOST_CHECK_THROW(Weighted_periodic_alpha_complex_3d wp_alpha_complex(wp_points, p_weights, 0., 0., 0., 1., 1., 1.), + std::invalid_argument); p_weights[14] = temp; std::cout << "wp_points and p_weights size exception" << std::endl; // Weights and points must have the same size // + 1 p_weights.push_back(1e-10); - BOOST_CHECK_THROW (Weighted_periodic_alpha_complex_3d wp_alpha_complex(wp_points, p_weights, 0., 0., 0., 1., 1., 1.), - std::invalid_argument); + BOOST_CHECK_THROW(Weighted_periodic_alpha_complex_3d wp_alpha_complex(wp_points, p_weights, 0., 0., 0., 1., 1., 1.), + std::invalid_argument); // - 1 p_weights.pop_back(); p_weights.pop_back(); - BOOST_CHECK_THROW (Weighted_periodic_alpha_complex_3d wp_alpha_complex(wp_points, p_weights, 0., 0., 0., 1., 1., 1.), - std::invalid_argument); + BOOST_CHECK_THROW(Weighted_periodic_alpha_complex_3d wp_alpha_complex(wp_points, p_weights, 0., 0., 0., 1., 1., 1.), + std::invalid_argument); } #endif @@ -491,7 +507,7 @@ BOOST_AUTO_TEST_CASE(Alpha_complex_weighted_periodic) { CGAL::Random_points_in_cube_3<Fast_weighted_periodic_alpha_complex_3d::Point_3, Creator> in_cube(1, random); std::vector<Fast_weighted_periodic_alpha_complex_3d::Point_3> p_points; - for (int i=0 ; i < 50 ; i++) { + for (int i = 0; i < 50; i++) { Fast_weighted_periodic_alpha_complex_3d::Point_3 p = *in_cube++; p_points.push_back(p); } @@ -513,7 +529,7 @@ BOOST_AUTO_TEST_CASE(Alpha_complex_weighted_periodic) { std::vector<Exact_weighted_periodic_alpha_complex_3d::Point_3> e_p_points; - for (auto p: p_points) { + for (auto p : p_points) { e_p_points.push_back(Exact_weighted_periodic_alpha_complex_3d::Point_3(p[0], p[1], p[2])); } @@ -535,17 +551,15 @@ BOOST_AUTO_TEST_CASE(Alpha_complex_weighted_periodic) { << " - Non exact is " << stree.num_vertices() << std::endl; BOOST_CHECK(exact_stree.num_vertices() == stree.num_vertices()); - // We cannot compare as objects from dispatcher on the alpha shape is not deterministic. // cf. https://github.com/CGAL/cgal/issues/3346 auto sh = stree.filtration_simplex_range().begin(); auto sh_exact = exact_stree.filtration_simplex_range().begin(); - while(sh != stree.filtration_simplex_range().end() || sh_exact != exact_stree.filtration_simplex_range().end()) { + while (sh != stree.filtration_simplex_range().end() || sh_exact != exact_stree.filtration_simplex_range().end()) { GUDHI_TEST_FLOAT_EQUALITY_CHECK(stree.filtration(*sh), exact_stree.filtration(*sh_exact), 1e-14); - std::vector<int> vh(stree.simplex_vertex_range(*sh).begin(), - stree.simplex_vertex_range(*sh).end()); + std::vector<int> vh(stree.simplex_vertex_range(*sh).begin(), stree.simplex_vertex_range(*sh).end()); std::vector<int> exact_vh(exact_stree.simplex_vertex_range(*sh_exact).begin(), exact_stree.simplex_vertex_range(*sh_exact).end()); @@ -557,7 +571,6 @@ BOOST_AUTO_TEST_CASE(Alpha_complex_weighted_periodic) { BOOST_CHECK(sh == stree.filtration_simplex_range().end()); BOOST_CHECK(sh_exact == exact_stree.filtration_simplex_range().end()); - // ---------------------- // Safe weighted periodic version // ---------------------- @@ -565,7 +578,7 @@ BOOST_AUTO_TEST_CASE(Alpha_complex_weighted_periodic) { std::vector<Safe_weighted_periodic_alpha_complex_3d::Point_3> s_p_points; - for (auto p: p_points) { + for (auto p : p_points) { s_p_points.push_back(Safe_weighted_periodic_alpha_complex_3d::Point_3(p[0], p[1], p[2])); } @@ -582,11 +595,10 @@ BOOST_AUTO_TEST_CASE(Alpha_complex_weighted_periodic) { sh = stree.filtration_simplex_range().begin(); auto sh_safe = safe_stree.filtration_simplex_range().begin(); - while(sh != stree.filtration_simplex_range().end() || sh_safe != safe_stree.filtration_simplex_range().end()) { + while (sh != stree.filtration_simplex_range().end() || sh_safe != safe_stree.filtration_simplex_range().end()) { GUDHI_TEST_FLOAT_EQUALITY_CHECK(stree.filtration(*sh), safe_stree.filtration(*sh_safe), 1e-14); - std::vector<int> vh(stree.simplex_vertex_range(*sh).begin(), - stree.simplex_vertex_range(*sh).end()); + std::vector<int> vh(stree.simplex_vertex_range(*sh).begin(), stree.simplex_vertex_range(*sh).end()); std::vector<int> safe_vh(safe_stree.simplex_vertex_range(*sh_safe).begin(), safe_stree.simplex_vertex_range(*sh_safe).end()); diff --git a/src/Alpha_complex/utilities/CMakeLists.txt b/src/Alpha_complex/utilities/CMakeLists.txt index 80444de8..65ca1624 100644 --- a/src/Alpha_complex/utilities/CMakeLists.txt +++ b/src/Alpha_complex/utilities/CMakeLists.txt @@ -30,9 +30,15 @@ if(CGAL_FOUND) "${CMAKE_SOURCE_DIR}/data/points/tore3D_300.off" "-p" "2" "-m" "0.45" "-o" "exact.pers" "-e") + add_test(NAME Alpha_complex_utilities_safe_alpha_complex_3d COMMAND $<TARGET_FILE:alpha_complex_3d_persistence> + "${CMAKE_SOURCE_DIR}/data/points/tore3D_300.off" + "-p" "2" "-m" "0.45" "-o" "safe.pers" "-s") + if (DIFF_PATH) add_test(Alpha_complex_utilities_diff_alpha_complex_3d ${DIFF_PATH} "exact.pers" "alpha.pers") + add_test(Alpha_complex_utilities_diff_alpha_complex_3d ${DIFF_PATH} + "safe.pers" "alpha.pers") endif() add_test(NAME Alpha_complex_utilities_periodic_alpha_complex_3d_persistence COMMAND $<TARGET_FILE:alpha_complex_3d_persistence> diff --git a/src/Alpha_complex/utilities/alpha_complex_3d_persistence.cpp b/src/Alpha_complex/utilities/alpha_complex_3d_persistence.cpp index c2b49fed..d14ba375 100644 --- a/src/Alpha_complex/utilities/alpha_complex_3d_persistence.cpp +++ b/src/Alpha_complex/utilities/alpha_complex_3d_persistence.cpp @@ -38,11 +38,12 @@ using Filtration_value = Simplex_tree::Filtration_value; using Persistent_cohomology = Gudhi::persistent_cohomology::Persistent_cohomology<Simplex_tree, Gudhi::persistent_cohomology::Field_Zp>; -void program_options(int argc, char *argv[], std::string &off_file_points, bool& EXACT, std::string &weight_file, - std::string &cuboid_file, std::string &output_file_diag, Filtration_value &alpha_square_max_value, - int &coeff_field_characteristic, Filtration_value &min_persistence); +void program_options(int argc, char *argv[], std::string &off_file_points, bool &exact, bool &safe, + std::string &weight_file, std::string &cuboid_file, std::string &output_file_diag, + Filtration_value &alpha_square_max_value, int &coeff_field_characteristic, + Filtration_value &min_persistence); -bool read_weight_file(const std::string& weight_file, std::vector<double>& weights) { +bool read_weight_file(const std::string &weight_file, std::vector<double> &weights) { // Read weights information from file std::ifstream weights_ifstr(weight_file); if (weights_ifstr.good()) { @@ -57,8 +58,8 @@ bool read_weight_file(const std::string& weight_file, std::vector<double>& weigh return true; } -bool read_cuboid_file(const std::string& cuboid_file, double& x_min, double& y_min, double& z_min, - double& x_max, double& y_max, double& z_max) { +bool read_cuboid_file(const std::string &cuboid_file, double &x_min, double &y_min, double &z_min, double &x_max, + double &y_max, double &z_max) { // Read weights information from file std::ifstream iso_cuboid_str(cuboid_file); if (iso_cuboid_str.is_open()) { @@ -71,8 +72,8 @@ bool read_cuboid_file(const std::string& cuboid_file, double& x_min, double& y_m return true; } -template<typename AlphaComplex3d> -std::vector<typename AlphaComplex3d::Point_3> read_off(const std::string& off_file_points) { +template <typename AlphaComplex3d> +std::vector<typename AlphaComplex3d::Point_3> read_off(const std::string &off_file_points) { // Read the OFF file (input file name given as parameter) and triangulate points Gudhi::Points_3D_off_reader<typename AlphaComplex3d::Point_3> off_reader(off_file_points); // Check the read operation was correct @@ -92,10 +93,11 @@ int main(int argc, char **argv) { int coeff_field_characteristic = 0; Filtration_value min_persistence = 0.; bool exact_version = false; + bool safe_version = false; bool weighted_version = false; bool periodic_version = false; - program_options(argc, argv, off_file_points, exact_version, weight_file, cuboid_file, output_file_diag, + program_options(argc, argv, off_file_points, exact_version, safe_version, weight_file, cuboid_file, output_file_diag, alpha_square_max_value, coeff_field_characteristic, min_persistence); std::vector<double> weights; @@ -107,7 +109,7 @@ int main(int argc, char **argv) { weighted_version = true; } - double x_min=0., y_min=0., z_min=0., x_max=0., y_max=0., z_max=0.; + double x_min = 0., y_min = 0., z_min = 0., x_max = 0., y_max = 0., z_max = 0.; std::ifstream iso_cuboid_str(argv[3]); if (cuboid_file != std::string()) { if (!read_cuboid_file(cuboid_file, x_min, y_min, z_min, x_max, y_max, z_max)) { @@ -119,37 +121,44 @@ int main(int argc, char **argv) { Gudhi::alpha_complex::complexity complexity = Gudhi::alpha_complex::complexity::FAST; if (exact_version) { + if (safe_version) { + std::cerr << "You cannot set the exact and the safe version." << std::endl; + exit(-1); + } complexity = Gudhi::alpha_complex::complexity::EXACT; } + if (safe_version) { + complexity = Gudhi::alpha_complex::complexity::SAFE; + } Simplex_tree simplex_tree; - switch(complexity) { + switch (complexity) { case Gudhi::alpha_complex::complexity::FAST: if (weighted_version) { if (periodic_version) { - using Alpha_complex_3d = Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::FAST, - true, true>; + using Alpha_complex_3d = + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::FAST, true, true>; auto points = read_off<Alpha_complex_3d>(off_file_points); Alpha_complex_3d alpha_complex(points, weights, x_min, y_min, z_min, x_max, y_max, z_max); alpha_complex.create_complex(simplex_tree, alpha_square_max_value); } else { - using Alpha_complex_3d = Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::FAST, - true, false>; + using Alpha_complex_3d = + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::FAST, true, false>; auto points = read_off<Alpha_complex_3d>(off_file_points); Alpha_complex_3d alpha_complex(points, weights); alpha_complex.create_complex(simplex_tree, alpha_square_max_value); } } else { if (periodic_version) { - using Alpha_complex_3d = Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::FAST, - false, true>; + using Alpha_complex_3d = + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::FAST, false, true>; auto points = read_off<Alpha_complex_3d>(off_file_points); Alpha_complex_3d alpha_complex(points, x_min, y_min, z_min, x_max, y_max, z_max); alpha_complex.create_complex(simplex_tree, alpha_square_max_value); } else { - using Alpha_complex_3d = Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::FAST, - false, false>; + using Alpha_complex_3d = + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::FAST, false, false>; auto points = read_off<Alpha_complex_3d>(off_file_points); Alpha_complex_3d alpha_complex(points); alpha_complex.create_complex(simplex_tree, alpha_square_max_value); @@ -159,28 +168,28 @@ int main(int argc, char **argv) { case Gudhi::alpha_complex::complexity::EXACT: if (weighted_version) { if (periodic_version) { - using Alpha_complex_3d = Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::EXACT, - true, true>; + using Alpha_complex_3d = + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::EXACT, true, true>; auto points = read_off<Alpha_complex_3d>(off_file_points); Alpha_complex_3d alpha_complex(points, weights, x_min, y_min, z_min, x_max, y_max, z_max); alpha_complex.create_complex(simplex_tree, alpha_square_max_value); } else { - using Alpha_complex_3d = Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::EXACT, - true, false>; + using Alpha_complex_3d = + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::EXACT, true, false>; auto points = read_off<Alpha_complex_3d>(off_file_points); Alpha_complex_3d alpha_complex(points, weights); alpha_complex.create_complex(simplex_tree, alpha_square_max_value); } } else { if (periodic_version) { - using Alpha_complex_3d = Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::EXACT, - false, true>; + using Alpha_complex_3d = + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::EXACT, false, true>; auto points = read_off<Alpha_complex_3d>(off_file_points); Alpha_complex_3d alpha_complex(points, x_min, y_min, z_min, x_max, y_max, z_max); alpha_complex.create_complex(simplex_tree, alpha_square_max_value); } else { - using Alpha_complex_3d = Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::EXACT, - false, false>; + using Alpha_complex_3d = + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::EXACT, false, false>; auto points = read_off<Alpha_complex_3d>(off_file_points); Alpha_complex_3d alpha_complex(points); alpha_complex.create_complex(simplex_tree, alpha_square_max_value); @@ -190,28 +199,28 @@ int main(int argc, char **argv) { case Gudhi::alpha_complex::complexity::SAFE: if (weighted_version) { if (periodic_version) { - using Alpha_complex_3d = Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::SAFE, - true, true>; + using Alpha_complex_3d = + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::SAFE, true, true>; auto points = read_off<Alpha_complex_3d>(off_file_points); Alpha_complex_3d alpha_complex(points, weights, x_min, y_min, z_min, x_max, y_max, z_max); alpha_complex.create_complex(simplex_tree, alpha_square_max_value); } else { - using Alpha_complex_3d = Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::SAFE, - true, false>; + using Alpha_complex_3d = + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::SAFE, true, false>; auto points = read_off<Alpha_complex_3d>(off_file_points); Alpha_complex_3d alpha_complex(points, weights); alpha_complex.create_complex(simplex_tree, alpha_square_max_value); } } else { if (periodic_version) { - using Alpha_complex_3d = Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::SAFE, - false, true>; + using Alpha_complex_3d = + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::SAFE, false, true>; auto points = read_off<Alpha_complex_3d>(off_file_points); Alpha_complex_3d alpha_complex(points, x_min, y_min, z_min, x_max, y_max, z_max); alpha_complex.create_complex(simplex_tree, alpha_square_max_value); } else { - using Alpha_complex_3d = Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::SAFE, - false, false>; + using Alpha_complex_3d = + Gudhi::alpha_complex::Alpha_complex_3d<Gudhi::alpha_complex::complexity::SAFE, false, false>; auto points = read_off<Alpha_complex_3d>(off_file_points); Alpha_complex_3d alpha_complex(points); alpha_complex.create_complex(simplex_tree, alpha_square_max_value); @@ -248,9 +257,10 @@ int main(int argc, char **argv) { return 0; } -void program_options(int argc, char *argv[], std::string &off_file_points, bool& EXACT, std::string &weight_file, - std::string &cuboid_file, std::string &output_file_diag, Filtration_value &alpha_square_max_value, - int &coeff_field_characteristic, Filtration_value &min_persistence) { +void program_options(int argc, char *argv[], std::string &off_file_points, bool &exact, bool &safe, + std::string &weight_file, std::string &cuboid_file, std::string &output_file_diag, + Filtration_value &alpha_square_max_value, int &coeff_field_characteristic, + Filtration_value &min_persistence) { namespace po = boost::program_options; po::options_description hidden("Hidden options"); hidden.add_options()("input-file", po::value<std::string>(&off_file_points), @@ -258,15 +268,18 @@ void program_options(int argc, char *argv[], std::string &off_file_points, bool& po::options_description visible("Allowed options", 100); visible.add_options()("help,h", "produce help message")( - "EXACT,e", po::bool_switch(&EXACT), - "To activate EXACT version of Alpha complex 3d (default is false, not available for weighted and/or periodic)")( + "exact,e", po::bool_switch(&exact), + "To activate exact version of Alpha complex 3d (default is false, not available if safe is set)")( + "safe,s", po::bool_switch(&safe), + "To activate safe version of Alpha complex 3d (default is false, not available if exact is set)")( "weight-file,w", po::value<std::string>(&weight_file)->default_value(std::string()), "Name of file containing a point weights. Format is one weight per line:\n W1\n ...\n Wn ")( "cuboid-file,c", po::value<std::string>(&cuboid_file), "Name of file describing the periodic domain. Format is:\n min_hx min_hy min_hz\n max_hx max_hy max_hz")( "output-file,o", po::value<std::string>(&output_file_diag)->default_value(std::string()), "Name of file in which the persistence diagram is written. Default print in std::cout")( - "max-alpha-square-value,r", po::value<Filtration_value>(&alpha_square_max_value) + "max-alpha-square-value,r", + po::value<Filtration_value>(&alpha_square_max_value) ->default_value(std::numeric_limits<Filtration_value>::infinity()), "Maximal alpha square value for the Alpha complex construction.")( "field-charac,p", po::value<int>(&coeff_field_characteristic)->default_value(11), @@ -289,7 +302,7 @@ void program_options(int argc, char *argv[], std::string &off_file_points, bool& std::cout << std::endl; std::cout << "Compute the persistent homology with coefficient field Z/pZ \n"; std::cout << "of a 3D Alpha complex defined on a set of input points.\n"; - std::cout << "3D Alpha complex can be EXACT or weighted and/or periodic\n\n"; + std::cout << "3D Alpha complex can be exact or safe, weighted and/or periodic\n\n"; std::cout << "The output diagram contains one bar per line, written with the convention: \n"; std::cout << " p dim b d \n"; std::cout << "where dim is the dimension of the homological feature,\n"; |