diff options
Diffstat (limited to 'src/Cech_complex/include')
-rw-r--r-- | src/Cech_complex/include/gudhi/Cech_complex.h | 10 | ||||
-rw-r--r-- | src/Cech_complex/include/gudhi/Cech_complex_blocker.h | 24 | ||||
-rw-r--r-- | src/Cech_complex/include/gudhi/Sphere_circumradius.h | 15 |
3 files changed, 31 insertions, 18 deletions
diff --git a/src/Cech_complex/include/gudhi/Cech_complex.h b/src/Cech_complex/include/gudhi/Cech_complex.h index fc39f75b..08b7a72f 100644 --- a/src/Cech_complex/include/gudhi/Cech_complex.h +++ b/src/Cech_complex/include/gudhi/Cech_complex.h @@ -30,7 +30,7 @@ namespace cech_complex { * \ingroup cech_complex * * \details - * Cech complex is a simplicial complex constructed from a proximity graph, where the set of all simplices is filtered + * Cech complex is a simplicial complex where the set of all simplices is filtered * by the radius of their minimal enclosing ball and bounded by the given max_radius. * * \tparam Kernel CGAL kernel: either Epick_d or Epeck_d. @@ -70,7 +70,7 @@ class Cech_complex { point_cloud_.assign(std::begin(points), std::end(points)); cech_skeleton_graph_ = Gudhi::compute_proximity_graph<SimplicialComplexForCechComplex>( - point_cloud_, max_radius_, Sphere_circumradius<Kernel>()); + point_cloud_, max_radius_, Sphere_circumradius<Kernel, Filtration_value>()); } /** \brief Initializes the simplicial complex from the proximity graph and expands it until a given maximal @@ -78,17 +78,19 @@ class Cech_complex { * * @param[in] complex SimplicialComplexForCech to be created. * @param[in] dim_max graph expansion until this given maximal dimension. + * @param[in] exact Exact filtration values computation. Not exact if `Kernel` is not <a target="_blank" + * href="https://doc.cgal.org/latest/Kernel_d/structCGAL_1_1Epeck__d.html">CGAL::Epeck_d</a>. * @exception std::invalid_argument In debug mode, if `complex.num_vertices()` does not return 0. * */ - void create_complex(SimplicialComplexForCechComplex& complex, int dim_max) { + void create_complex(SimplicialComplexForCechComplex& complex, int dim_max, const bool exact = false) { GUDHI_CHECK(complex.num_vertices() == 0, std::invalid_argument("Cech_complex::create_complex - simplicial complex is not empty")); // insert the proximity graph in the simplicial complex complex.insert_graph(cech_skeleton_graph_); // expand the graph until dimension dim_max - complex.expansion_with_blockers(dim_max, cech_blocker(&complex, this)); + complex.expansion_with_blockers(dim_max, cech_blocker(&complex, this, exact)); } /** @return max_radius value given at construction. */ diff --git a/src/Cech_complex/include/gudhi/Cech_complex_blocker.h b/src/Cech_complex/include/gudhi/Cech_complex_blocker.h index 9917999f..1bb205b3 100644 --- a/src/Cech_complex/include/gudhi/Cech_complex_blocker.h +++ b/src/Cech_complex/include/gudhi/Cech_complex_blocker.h @@ -12,6 +12,7 @@ #define CECH_COMPLEX_BLOCKER_H_ #include <CGAL/NT_converter.h> // for casting from FT to Filtration_value +#include <CGAL/Lazy_exact_nt.h> // for CGAL::exact #include <iostream> #include <vector> @@ -84,9 +85,9 @@ class Cech_blocker { Point_cloud face_points; for (auto vertex : sc_ptr_->simplex_vertex_range(face_opposite_vertex.first)) { face_points.push_back(cc_ptr_->get_point(vertex)); - #ifdef DEBUG_TRACES - std::clog << "#(" << vertex << ")#"; - #endif // DEBUG_TRACES +#ifdef DEBUG_TRACES + std::clog << "#(" << vertex << ")#"; +#endif // DEBUG_TRACES } sph = get_sphere(face_points.cbegin(), face_points.cend()); // Put edge sphere in cache @@ -97,10 +98,13 @@ class Cech_blocker { } // Check if the minimal enclosing ball of current face contains the extra point/opposite vertex if (kernel_.squared_distance_d_object()(sph.first, cc_ptr_->get_point(face_opposite_vertex.second)) <= sph.second) { - #ifdef DEBUG_TRACES - std::clog << "center: " << sph.first << ", radius: " << radius << std::endl; - #endif // DEBUG_TRACES +#ifdef DEBUG_TRACES + std::clog << "center: " << sph.first << ", radius: " << radius << std::endl; +#endif // DEBUG_TRACES is_min_enclos_ball = true; +#if CGAL_VERSION_NR >= 1050000000 + if(exact_) CGAL::exact(sph.second); +#endif radius = std::sqrt(cast_to_fv(sph.second)); sc_ptr_->assign_key(sh, cc_ptr_->get_cache().size()); cc_ptr_->get_cache().push_back(sph); @@ -114,10 +118,13 @@ class Cech_blocker { points.push_back(cc_ptr_->get_point(vertex)); } Sphere sph = get_sphere(points.cbegin(), points.cend()); +#if CGAL_VERSION_NR >= 1050000000 + if(exact_) CGAL::exact(sph.second); +#endif radius = std::sqrt(cast_to_fv(sph.second)); sc_ptr_->assign_key(sh, cc_ptr_->get_cache().size()); - cc_ptr_->get_cache().push_back(sph); + cc_ptr_->get_cache().push_back(std::move(sph)); } #ifdef DEBUG_TRACES @@ -128,12 +135,13 @@ class Cech_blocker { } /** \internal \brief Čech complex blocker constructor. */ - Cech_blocker(SimplicialComplexForCech* sc_ptr, Cech_complex* cc_ptr) : sc_ptr_(sc_ptr), cc_ptr_(cc_ptr) {} + Cech_blocker(SimplicialComplexForCech* sc_ptr, Cech_complex* cc_ptr, const bool exact) : sc_ptr_(sc_ptr), cc_ptr_(cc_ptr), exact_(exact) {} private: SimplicialComplexForCech* sc_ptr_; Cech_complex* cc_ptr_; Kernel kernel_; + const bool exact_; }; } // namespace cech_complex diff --git a/src/Cech_complex/include/gudhi/Sphere_circumradius.h b/src/Cech_complex/include/gudhi/Sphere_circumradius.h index b0d9f7cc..790f6950 100644 --- a/src/Cech_complex/include/gudhi/Sphere_circumradius.h +++ b/src/Cech_complex/include/gudhi/Sphere_circumradius.h @@ -11,7 +11,7 @@ #ifndef SPHERE_CIRCUMRADIUS_H_ #define SPHERE_CIRCUMRADIUS_H_ -#include <CGAL/Epeck_d.h> // for #include <CGAL/NewKernel_d/KernelD_converter.h> +#include <CGAL/Epick_d.h> // for #include <CGAL/NT_converter.h> which is not working/compiling alone #include <cmath> // for std::sqrt #include <vector> @@ -22,14 +22,17 @@ namespace cech_complex { /** \private @brief Compute the circumradius of the sphere passing through points given by a range of coordinates. * The points are assumed to have the same dimension. */ -template<typename Kernel> +template<typename Kernel, typename Filtration_value> class Sphere_circumradius { private: Kernel kernel_; public: + using FT = typename Kernel::FT; using Point = typename Kernel::Point_d; using Point_cloud = typename std::vector<Point>; + CGAL::NT_converter<FT, Filtration_value> cast_to_fv; + /** \brief Circumradius of sphere passing through two points using CGAL. * * @param[in] point_1 @@ -38,8 +41,8 @@ class Sphere_circumradius { * \tparam Point must be a Kernel::Point_d from CGAL. * */ - double operator()(const Point& point_1, const Point& point_2) const { - return std::sqrt(CGAL::to_double(kernel_.squared_distance_d_object()(point_1, point_2))) / 2.; + Filtration_value operator()(const Point& point_1, const Point& point_2) const { + return std::sqrt(cast_to_fv(kernel_.squared_distance_d_object()(point_1, point_2))) / 2.; } /** \brief Circumradius of sphere passing through point cloud using CGAL. @@ -49,8 +52,8 @@ class Sphere_circumradius { * \tparam Point_cloud must be a range of Kernel::Point_d points from CGAL. * */ - double operator()(const Point_cloud& point_cloud) const { - return std::sqrt(CGAL::to_double(kernel_.compute_squared_radius_d_object()(point_cloud.begin(), point_cloud.end()))); + Filtration_value operator()(const Point_cloud& point_cloud) const { + return std::sqrt(cast_to_fv(kernel_.compute_squared_radius_d_object()(point_cloud.begin(), point_cloud.end()))); } }; |