Loading...
Searching...
No Matches
Alpha_complex_3d.h
1/* This file is part of the Gudhi Library - https://gudhi.inria.fr/ - which is released under MIT.
2 * See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details.
3 * Author(s): Vincent Rouvreau
4 *
5 * Copyright (C) 2018 Inria
6 *
7 * Modification(s):
8 * - 2019/08 Vincent Rouvreau: Fix issue #10 for CGAL and Eigen3
9 * - YYYY/MM Author: Description of the modification
10 */
11
12#ifndef ALPHA_COMPLEX_3D_H_
13#define ALPHA_COMPLEX_3D_H_
14
15#include <boost/variant.hpp>
16#include <boost/range/size.hpp>
17#include <boost/range/combine.hpp>
18#include <boost/range/adaptor/transformed.hpp>
19
20#include <gudhi/Debug_utils.h>
21#include <gudhi/Alpha_complex_options.h>
22
23#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
24#include <CGAL/Exact_predicates_exact_constructions_kernel.h>
25#include <CGAL/Delaunay_triangulation_3.h>
26#include <CGAL/Periodic_3_Delaunay_triangulation_traits_3.h>
27#include <CGAL/Periodic_3_Delaunay_triangulation_3.h>
28#include <CGAL/Periodic_3_regular_triangulation_traits_3.h>
29#include <CGAL/Periodic_3_regular_triangulation_3.h>
30#include <CGAL/Regular_triangulation_3.h>
31#include <CGAL/Alpha_shape_3.h>
32#include <CGAL/Alpha_shape_cell_base_3.h>
33#include <CGAL/Alpha_shape_vertex_base_3.h>
34
35#include <CGAL/Object.h>
36#include <CGAL/tuple.h>
37#include <CGAL/iterator.h>
38#include <CGAL/version.h> // for CGAL_VERSION_NR
39
40#include <boost/container/static_vector.hpp>
41
42#include <iostream>
43#include <vector>
44#include <unordered_map>
45#include <stdexcept>
46#include <cstddef> // for std::size_t
47#include <memory> // for std::unique_ptr
48#include <type_traits> // for std::conditional and std::enable_if
49#include <limits> // for numeric_limits<>
50
51// Make compilation fail - required for external projects - https://github.com/GUDHI/gudhi-devel/issues/10
52#if CGAL_VERSION_NR < 1041101000
53# error Alpha_complex_3d is only available for CGAL >= 4.11
54#endif
55
56namespace Gudhi {
57
58namespace alpha_complex {
59
60// Value_from_iterator returns the filtration value from an iterator on alpha shapes values
61//
62// FAST SAFE EXACT
63// CGAL::to_double(*iterator) CGAL::to_double(*iterator) CGAL::to_double(iterator->exact())
64
65template <complexity Complexity>
66struct Value_from_iterator {
67 template <typename Iterator>
68 static double perform(Iterator it) {
69 // Default behaviour
70 return CGAL::to_double(*it);
71 }
72};
73
74template <>
75struct Value_from_iterator<complexity::EXACT> {
76 template <typename Iterator>
77 static double perform(Iterator it) {
78 return CGAL::to_double(it->exact());
79 }
80};
81
117template <complexity Complexity = complexity::SAFE, bool Weighted = false, bool Periodic = false>
119 // Epick = Exact_predicates_inexact_constructions_kernel
120 // Epeck = Exact_predicates_exact_constructions_kernel
121 // Exact_alpha_comparison_tag = exact version of CGAL Alpha_shape_3 and of its objects (Alpha_shape_vertex_base_3 and
122 // Alpha_shape_cell_base_3). Not available if weighted or periodic.
123 // Can be CGAL::Tag_false or CGAL::Tag_true. Default is False.
124 // cf. https://doc.cgal.org/latest/Alpha_shapes_3/classCGAL_1_1Alpha__shape__3.html
125 //
126 // We could use Epick + CGAL::Tag_true for not weighted nor periodic, but during benchmark, we found a bug
127 // https://github.com/CGAL/cgal/issues/3460
128 // This is the reason we only use Epick + CGAL::Tag_false, or Epeck
129 //
130 // FAST SAFE EXACT
131 // Epick + CGAL::Tag_false Epeck Epeck
132 using Predicates = typename std::conditional<(Complexity == complexity::FAST),
133 CGAL::Exact_predicates_inexact_constructions_kernel,
134 CGAL::Exact_predicates_exact_constructions_kernel>::type;
135
136 // The other way to do a conditional type. Here there are 3 possibilities
137 template <typename Predicates, bool Weighted_version, bool Periodic_version>
138 struct Kernel_3 {};
139
140 template <typename Predicates, bool Is_periodic>
141 struct Kernel_3<Predicates, Is_periodic, false> {
142 using Kernel = Predicates;
143 };
144
145 template <typename Predicates>
146 struct Kernel_3<Predicates, false, true> {
147 using Kernel = CGAL::Periodic_3_Delaunay_triangulation_traits_3<Predicates>;
148 };
149 template <typename Predicates>
150 struct Kernel_3<Predicates, true, true> {
151 using Kernel = CGAL::Periodic_3_regular_triangulation_traits_3<Predicates>;
152 };
153
154 public:
155 using Kernel = typename Kernel_3<Predicates, Weighted, Periodic>::Kernel;
156
157 private:
158 using TdsVb = typename std::conditional<Periodic, CGAL::Periodic_3_triangulation_ds_vertex_base_3<>,
159 CGAL::Triangulation_ds_vertex_base_3<>>::type;
160
161 using Tvb = typename std::conditional<Weighted, CGAL::Regular_triangulation_vertex_base_3<Kernel, TdsVb>,
162 CGAL::Triangulation_vertex_base_3<Kernel, TdsVb>>::type;
163
164 using Vb = CGAL::Alpha_shape_vertex_base_3<Kernel, Tvb>;
165
166 using TdsCb = typename std::conditional<Periodic, CGAL::Periodic_3_triangulation_ds_cell_base_3<>,
167 CGAL::Triangulation_ds_cell_base_3<>>::type;
168
169 using Tcb = typename std::conditional<Weighted, CGAL::Regular_triangulation_cell_base_3<Kernel, TdsCb>,
170 CGAL::Triangulation_cell_base_3<Kernel, TdsCb>>::type;
171
172 using Cb = CGAL::Alpha_shape_cell_base_3<Kernel, Tcb>;
173 using Tds = CGAL::Triangulation_data_structure_3<Vb, Cb>;
174
175 // The other way to do a conditional type. Here there 4 possibilities, cannot use std::conditional
176 template <typename Kernel, typename Tds, bool Weighted_version, bool Periodic_version>
177 struct Triangulation_3 {};
178
179 template <typename Kernel, typename Tds>
180 struct Triangulation_3<Kernel, Tds, false, false> {
181 using Dt = CGAL::Delaunay_triangulation_3<Kernel, Tds>;
182 using Weighted_point_3 = void;
183 };
184 template <typename Kernel, typename Tds>
185 struct Triangulation_3<Kernel, Tds, true, false> {
186 using Dt = CGAL::Regular_triangulation_3<Kernel, Tds>;
187 using Weighted_point_3 = typename Dt::Weighted_point;
188 };
189 template <typename Kernel, typename Tds>
190 struct Triangulation_3<Kernel, Tds, false, true> {
191 using Dt = CGAL::Periodic_3_Delaunay_triangulation_3<Kernel, Tds>;
192 using Weighted_point_3 = void;
193 };
194 template <typename Kernel, typename Tds>
195 struct Triangulation_3<Kernel, Tds, true, true> {
196 using Dt = CGAL::Periodic_3_regular_triangulation_3<Kernel, Tds>;
197 using Weighted_point_3 = typename Dt::Weighted_point;
198 };
199
207 using Dt = typename Triangulation_3<Kernel, Tds, Weighted, Periodic>::Dt;
208
209 public:
215 using Alpha_shape_3 = CGAL::Alpha_shape_3<Dt>;
216
219 using FT = typename Alpha_shape_3::FT;
220
230 using Bare_point_3 = typename Kernel::Point_3;
231
243 using Weighted_point_3 = typename Triangulation_3<Kernel, Tds, Weighted, Periodic>::Weighted_point_3;
244
248 using Point_3 = typename Alpha_shape_3::Point;
249
250 private:
251 using Dispatch =
252 CGAL::Dispatch_output_iterator<CGAL::cpp11::tuple<CGAL::Object, FT>,
253 CGAL::cpp11::tuple<std::back_insert_iterator<std::vector<CGAL::Object>>,
254 std::back_insert_iterator<std::vector<FT>>>>;
255
256 using Cell_handle = typename Alpha_shape_3::Cell_handle;
257 using Facet = typename Alpha_shape_3::Facet;
258 using Edge = typename Alpha_shape_3::Edge;
259 using Alpha_vertex_handle = typename Alpha_shape_3::Vertex_handle;
260 using Vertex_list = boost::container::static_vector<Alpha_vertex_handle, 4>;
261
262 public:
272 template <typename InputPointRange>
273 Alpha_complex_3d(const InputPointRange& points) {
274 static_assert(!Periodic, "This constructor is not available for periodic versions of Alpha_complex_3d");
275
276 alpha_shape_3_ptr_ = std::make_unique<Alpha_shape_3>(
277 std::begin(points), std::end(points), 0, Alpha_shape_3::GENERAL);
278 }
279
294 template <typename InputPointRange, typename WeightRange>
295 Alpha_complex_3d(const InputPointRange& points, WeightRange weights) {
296 static_assert(Weighted, "This constructor is not available for non-weighted versions of Alpha_complex_3d");
297 static_assert(!Periodic, "This constructor is not available for periodic versions of Alpha_complex_3d");
298 // FIXME: this test is only valid if we have a forward range
299 GUDHI_CHECK(boost::size(weights) == boost::size(points),
300 std::invalid_argument("Points number in range different from weights range number"));
301
302 auto weighted_points_3 = boost::range::combine(points, weights)
303 | boost::adaptors::transformed([](auto const&t){return Weighted_point_3(boost::get<0>(t), boost::get<1>(t));});
304
305 alpha_shape_3_ptr_ = std::make_unique<Alpha_shape_3>(
306 std::begin(weighted_points_3), std::end(weighted_points_3), 0, Alpha_shape_3::GENERAL);
307 }
308
329 template <typename InputPointRange>
330 Alpha_complex_3d(const InputPointRange& points, FT x_min, FT y_min,
331 FT z_min, FT x_max, FT y_max, FT z_max) {
332 static_assert(Periodic, "This constructor is not available for non-periodic versions of Alpha_complex_3d");
333 // Checking if the cuboid is the same in x,y and z direction. If not, CGAL will not process it.
334 GUDHI_CHECK(
335 (x_max - x_min == y_max - y_min) && (x_max - x_min == z_max - z_min) && (z_max - z_min == y_max - y_min),
336 std::invalid_argument("The size of the cuboid in every directions is not the same."));
337
338 // Define the periodic cube
339 Dt pdt(typename Kernel::Iso_cuboid_3(x_min, y_min, z_min, x_max, y_max, z_max));
340 // Heuristic for inserting large point sets (if pts is reasonably large)
341 pdt.insert(std::begin(points), std::end(points), true);
342 // As pdt won't be modified anymore switch to 1-sheeted cover if possible
343 if (!pdt.is_triangulation_in_1_sheet()) {
344 throw std::invalid_argument("Unable to construct a triangulation within a single periodic domain.");
345 }
346 pdt.convert_to_1_sheeted_covering();
347
348 // alpha shape construction from points. CGAL has a strange behavior in REGULARIZED mode. This is the default mode
349 // Maybe need to set it to GENERAL mode
350 alpha_shape_3_ptr_ = std::make_unique<Alpha_shape_3>(pdt, 0, Alpha_shape_3::GENERAL);
351 }
352
377 template <typename InputPointRange, typename WeightRange>
378 Alpha_complex_3d(const InputPointRange& points, WeightRange weights, FT x_min, FT y_min,
379 FT z_min, FT x_max, FT y_max, FT z_max) {
380 static_assert(Weighted, "This constructor is not available for non-weighted versions of Alpha_complex_3d");
381 static_assert(Periodic, "This constructor is not available for non-periodic versions of Alpha_complex_3d");
382 // FIXME: this test is only valid if we have a forward range
383 GUDHI_CHECK(boost::size(weights) == boost::size(points),
384 std::invalid_argument("Points number in range different from weights range number"));
385 // Checking if the cuboid is the same in x,y and z direction. If not, CGAL will not process it.
386 GUDHI_CHECK(
387 (x_max - x_min == y_max - y_min) && (x_max - x_min == z_max - z_min) && (z_max - z_min == y_max - y_min),
388 std::invalid_argument("The size of the cuboid in every directions is not the same."));
389
390#ifdef GUDHI_DEBUG
391 // Defined in GUDHI_DEBUG to avoid unused variable warning for GUDHI_CHECK
392 FT maximal_possible_weight = 0.015625 * (x_max - x_min) * (x_max - x_min);
393#endif
394
395 auto weighted_points_3 = boost::range::combine(points, weights)
396 | boost::adaptors::transformed([=](auto const&t){
397 auto w = boost::get<1>(t);
398 GUDHI_CHECK((w < maximal_possible_weight) && (w >= 0),
399 std::invalid_argument("Invalid weight " + std::to_string(w) +
400 ". Must be non-negative and less than maximal possible weight = 1/64*cuboid length squared."));
401 return Weighted_point_3(boost::get<0>(t), w);
402 });
403
404 // Define the periodic cube
405 Dt pdt(typename Kernel::Iso_cuboid_3(x_min, y_min, z_min, x_max, y_max, z_max));
406 // Heuristic for inserting large point sets (if pts is reasonably large)
407 pdt.insert(std::begin(weighted_points_3), std::end(weighted_points_3), true);
408 // As pdt won't be modified anymore switch to 1-sheeted cover if possible
409 if (!pdt.is_triangulation_in_1_sheet()) {
410 throw std::invalid_argument("Unable to construct a triangulation within a single periodic domain.");
411 }
412 pdt.convert_to_1_sheeted_covering();
413
414 // alpha shape construction from points. CGAL has a strange behavior in REGULARIZED mode. This is the default mode
415 // Maybe need to set it to GENERAL mode
416 alpha_shape_3_ptr_ = std::make_unique<Alpha_shape_3>(pdt, 0, Alpha_shape_3::GENERAL);
417 }
418
433 template <typename SimplicialComplexForAlpha3d,
436 Filtration_value max_alpha_square = std::numeric_limits<Filtration_value>::infinity()) {
437 if (complex.num_vertices() > 0) {
438 std::cerr << "Alpha_complex_3d create_complex - complex is not empty\n";
439 return false; // ----- >>
440 }
441
442 using Complex_vertex_handle = typename SimplicialComplexForAlpha3d::Vertex_handle;
443 using Simplex_tree_vector_vertex = std::vector<Complex_vertex_handle>;
444
445#ifdef DEBUG_TRACES
446 std::size_t count_vertices = 0;
447 std::size_t count_edges = 0;
448 std::size_t count_facets = 0;
449 std::size_t count_cells = 0;
450#endif // DEBUG_TRACES
451 std::vector<CGAL::Object> objects;
452 std::vector<FT> alpha_values;
453
454 Dispatch dispatcher = CGAL::dispatch_output<CGAL::Object, FT>(std::back_inserter(objects),
455 std::back_inserter(alpha_values));
456
457 alpha_shape_3_ptr_->filtration_with_alpha_values(dispatcher);
458#ifdef DEBUG_TRACES
459 std::clog << "filtration_with_alpha_values returns : " << objects.size() << " objects" << std::endl;
460#endif // DEBUG_TRACES
461 if (objects.size() == 0) {
462 std::cerr << "Alpha_complex_3d create_complex - no triangulation as points are on a 2d plane\n";
463 return false; // ----- >>
464 }
465
466 using Alpha_value_iterator = typename std::vector<FT>::const_iterator;
467 Alpha_value_iterator alpha_value_iterator = alpha_values.begin();
468 for (auto object_iterator : objects) {
469 Vertex_list vertex_list;
470
471 // Retrieve Alpha shape vertex list from object
472 if (const Cell_handle* cell = CGAL::object_cast<Cell_handle>(&object_iterator)) {
473 for (auto i = 0; i < 4; i++) {
474#ifdef DEBUG_TRACES
475 std::clog << "from cell[" << i << "] - Point coordinates (" << (*cell)->vertex(i)->point() << ")"
476 << std::endl;
477#endif // DEBUG_TRACES
478 vertex_list.push_back((*cell)->vertex(i));
479 }
480#ifdef DEBUG_TRACES
481 count_cells++;
482#endif // DEBUG_TRACES
483 } else if (const Facet* facet = CGAL::object_cast<Facet>(&object_iterator)) {
484 for (auto i = 0; i < 4; i++) {
485 if ((*facet).second != i) {
486#ifdef DEBUG_TRACES
487 std::clog << "from facet=[" << i << "] - Point coordinates (" << (*facet).first->vertex(i)->point() << ")"
488 << std::endl;
489#endif // DEBUG_TRACES
490 vertex_list.push_back((*facet).first->vertex(i));
491 }
492 }
493#ifdef DEBUG_TRACES
494 count_facets++;
495#endif // DEBUG_TRACES
496 } else if (const Edge* edge = CGAL::object_cast<Edge>(&object_iterator)) {
497 for (auto i : {(*edge).second, (*edge).third}) {
498#ifdef DEBUG_TRACES
499 std::clog << "from edge[" << i << "] - Point coordinates (" << (*edge).first->vertex(i)->point() << ")"
500 << std::endl;
501#endif // DEBUG_TRACES
502 vertex_list.push_back((*edge).first->vertex(i));
503 }
504#ifdef DEBUG_TRACES
505 count_edges++;
506#endif // DEBUG_TRACES
507 } else if (const Alpha_vertex_handle* vertex = CGAL::object_cast<Alpha_vertex_handle>(&object_iterator)) {
508#ifdef DEBUG_TRACES
509 count_vertices++;
510 std::clog << "from vertex - Point coordinates (" << (*vertex)->point() << ")" << std::endl;
511#endif // DEBUG_TRACES
512 vertex_list.push_back((*vertex));
513 }
514 // Construction of the vector of simplex_tree vertex from list of alpha_shapes vertex
515 Simplex_tree_vector_vertex the_simplex;
516 for (auto the_alpha_shape_vertex : vertex_list) {
517 auto the_map_iterator = map_cgal_simplex_tree.find(the_alpha_shape_vertex);
518 if (the_map_iterator == map_cgal_simplex_tree.end()) {
519 // alpha shape not found
520 Complex_vertex_handle vertex = map_cgal_simplex_tree.size();
521#ifdef DEBUG_TRACES
522 std::clog << "Point (" << the_alpha_shape_vertex->point() << ") not found - insert new vertex id " << vertex
523 << std::endl;
524#endif // DEBUG_TRACES
525 the_simplex.push_back(vertex);
526 map_cgal_simplex_tree.emplace(the_alpha_shape_vertex, vertex);
527 } else {
528 // alpha shape found
529 Complex_vertex_handle vertex = the_map_iterator->second;
530#ifdef DEBUG_TRACES
531 std::clog << "Point (" << the_alpha_shape_vertex->point() << ") found as vertex id " << vertex << std::endl;
532#endif // DEBUG_TRACES
533 the_simplex.push_back(vertex);
534 }
535 }
536 // Construction of the simplex_tree
537 Filtration_value filtr = Value_from_iterator<Complexity>::perform(alpha_value_iterator);
538
539#ifdef DEBUG_TRACES
540 std::clog << "filtration = " << filtr << std::endl;
541#endif // DEBUG_TRACES
542 complex.insert_simplex(the_simplex, static_cast<Filtration_value>(filtr));
543 GUDHI_CHECK(alpha_value_iterator != alpha_values.end(), "CGAL provided more simplices than values");
544 ++alpha_value_iterator;
545 }
546
547#ifdef DEBUG_TRACES
548 std::clog << "vertices \t" << count_vertices << std::endl;
549 std::clog << "edges \t\t" << count_edges << std::endl;
550 std::clog << "facets \t\t" << count_facets << std::endl;
551 std::clog << "cells \t\t" << count_cells << std::endl;
552#endif // DEBUG_TRACES
553 // --------------------------------------------------------------------------------------------
554 if (Complexity == complexity::FAST)
555 // As Alpha value is an approximation, we have to make filtration non decreasing while increasing the dimension
556 // Only in FAST version, cf. https://github.com/GUDHI/gudhi-devel/issues/57
558 // Remove all simplices that have a filtration value greater than max_alpha_square
559 complex.prune_above_filtration(max_alpha_square);
560 // --------------------------------------------------------------------------------------------
561 return true;
562 }
563
570 const Point_3& get_point(std::size_t vertex) {
571 if (map_cgal_simplex_tree.size() != cgal_vertex_iterator_vector.size()) {
572 cgal_vertex_iterator_vector.resize(map_cgal_simplex_tree.size());
573 for (auto map_iterator : map_cgal_simplex_tree) {
574 cgal_vertex_iterator_vector[map_iterator.second] = map_iterator.first;
575 }
576
577 }
578 auto cgal_vertex_iterator = cgal_vertex_iterator_vector.at(vertex);
579 return cgal_vertex_iterator->point();
580 }
581
582 private:
583 // use of a unique_ptr on cgal Alpha_shape_3, as copy and default constructor is not available - no need to be freed
584 std::unique_ptr<Alpha_shape_3> alpha_shape_3_ptr_;
585
586 // Map type to switch from CGAL vertex iterator to simplex tree vertex handle.
587 std::unordered_map<Alpha_vertex_handle, std::size_t> map_cgal_simplex_tree;
588 // Vector type to switch from simplex tree vertex handle to CGAL vertex iterator.
589 std::vector<Alpha_vertex_handle> cgal_vertex_iterator_vector;
590};
591
592} // namespace alpha_complex
593
594} // namespace Gudhi
595
596#endif // ALPHA_COMPLEX_3D_H_
Alpha complex data structure for 3d specific case.
Definition: Alpha_complex_3d.h:118
CGAL::Alpha_shape_3< Dt > Alpha_shape_3
The CGAL 3D Alpha Shapes type.
Definition: Alpha_complex_3d.h:215
Alpha_complex_3d(const InputPointRange &points)
Alpha_complex constructor from a list of points.
Definition: Alpha_complex_3d.h:273
bool create_complex(SimplicialComplexForAlpha3d &complex, Filtration_value max_alpha_square=std::numeric_limits< Filtration_value >::infinity())
Inserts all Delaunay triangulation into the simplicial complex. It also computes the filtration value...
Definition: Alpha_complex_3d.h:435
Alpha_complex_3d(const InputPointRange &points, WeightRange weights, FT x_min, FT y_min, FT z_min, FT x_max, FT y_max, FT z_max)
Alpha_complex constructor from a list of points, associated weights and an iso-cuboid coordinates.
Definition: Alpha_complex_3d.h:378
Alpha_complex_3d(const InputPointRange &points, WeightRange weights)
Alpha_complex constructor from a list of points and associated weights.
Definition: Alpha_complex_3d.h:295
typename Kernel::Point_3 Bare_point_3
Gives public access to the Bare_point_3 (bare aka. unweighed) type. Here is a Bare_point_3 constructo...
Definition: Alpha_complex_3d.h:230
typename Triangulation_3< Kernel, Tds, Weighted, Periodic >::Weighted_point_3 Weighted_point_3
Gives public access to the Weighted_point_3 type. A Weighted point can be constructed as follows:
Definition: Alpha_complex_3d.h:243
typename Alpha_shape_3::Point Point_3
Alpha_complex_3d::Point_3 type is either a Alpha_complex_3d::Bare_point_3 (Weighted = false) or a Alp...
Definition: Alpha_complex_3d.h:248
typename Alpha_shape_3::FT FT
The alpha values type. Must be compatible with double.
Definition: Alpha_complex_3d.h:219
Alpha_complex_3d(const InputPointRange &points, FT x_min, FT y_min, FT z_min, FT x_max, FT y_max, FT z_max)
Alpha_complex constructor from a list of points and an iso-cuboid coordinates.
Definition: Alpha_complex_3d.h:330
const Point_3 & get_point(std::size_t vertex)
get_point returns the point corresponding to the vertex given as parameter.
Definition: Alpha_complex_3d.h:570
complexity
Alpha complex complexity template parameter possible values.
Definition: Alpha_complex_options.h:23
Value type for a filtration function on a cell complex.
Definition: FiltrationValue.h:20
The concept SimplicialComplexForAlpha3d describes the requirements for a type to implement a simplici...
Definition: SimplicialComplexForAlpha3d.h:21
void insert_simplex(std::vector< Vertex_handle > const &vertex_range, Filtration_value filtration)
Inserts a simplex from a given simplex (represented by a vector of Vertex_handle) in the simplicial c...
unspecified Filtration_value
Definition: SimplicialComplexForAlpha3d.h:25
unspecified Vertex_handle
Definition: SimplicialComplexForAlpha3d.h:23
void prune_above_filtration(Filtration_value filtration)