diff options
author | mcarrier <mcarrier@636b058d-ea47-450e-bf9e-a15bfbe3eedb> | 2018-02-15 16:41:04 +0000 |
---|---|---|
committer | mcarrier <mcarrier@636b058d-ea47-450e-bf9e-a15bfbe3eedb> | 2018-02-15 16:41:04 +0000 |
commit | 69c683e663329d8410ca77c371f877bcc3bef906 (patch) | |
tree | fb18c914cf4056881b2b31875eb6b44e5ce23895 /src/Alpha_complex/utilities | |
parent | be131d6f74a9264e15a0b1c1e72fa8967c4518bd (diff) | |
parent | 265484997185f3bf900744406206a2d64ca0a20d (diff) |
integrated kernel code in pers_representation branch
git-svn-id: svn+ssh://scm.gforge.inria.fr/svnroot/gudhi/branches/kernels@3249 636b058d-ea47-450e-bf9e-a15bfbe3eedb
Former-commit-id: 58e77263a0da3674e4699cef832b6d357dcf12e2
Diffstat (limited to 'src/Alpha_complex/utilities')
-rw-r--r-- | src/Alpha_complex/utilities/CMakeLists.txt | 2 | ||||
-rw-r--r-- | src/Alpha_complex/utilities/alpha_complex_3d_helper.h | 8 | ||||
-rw-r--r-- | src/Alpha_complex/utilities/alpha_complex_3d_persistence.cpp | 59 | ||||
-rw-r--r-- | src/Alpha_complex/utilities/alphacomplex.md (renamed from src/Alpha_complex/utilities/README) | 173 | ||||
-rw-r--r-- | src/Alpha_complex/utilities/exact_alpha_complex_3d_persistence.cpp | 43 | ||||
-rw-r--r-- | src/Alpha_complex/utilities/periodic_alpha_complex_3d_persistence.cpp | 54 | ||||
-rw-r--r-- | src/Alpha_complex/utilities/weighted_alpha_complex_3d_persistence.cpp | 51 | ||||
-rw-r--r-- | src/Alpha_complex/utilities/weighted_periodic_alpha_complex_3d_persistence.cpp | 95 |
8 files changed, 206 insertions, 279 deletions
diff --git a/src/Alpha_complex/utilities/CMakeLists.txt b/src/Alpha_complex/utilities/CMakeLists.txt index 79d9e7dd..a2dfac20 100644 --- a/src/Alpha_complex/utilities/CMakeLists.txt +++ b/src/Alpha_complex/utilities/CMakeLists.txt @@ -54,7 +54,7 @@ if(CGAL_FOUND) target_link_libraries(weighted_periodic_alpha_complex_3d_persistence ${TBB_LIBRARIES}) endif(TBB_FOUND) - add_test(NAME Persistent_cohomology_example_weigted_periodic_alpha_complex_3d COMMAND $<TARGET_FILE:weighted_periodic_alpha_complex_3d_persistence> + add_test(NAME Alpha_complex_utilities_weigted_periodic_alpha_complex_3d COMMAND $<TARGET_FILE:weighted_periodic_alpha_complex_3d_persistence> "${CMAKE_SOURCE_DIR}/data/points/grid_10_10_10_in_0_1.off" "${CMAKE_SOURCE_DIR}/data/points/grid_10_10_10_in_0_1.weights" "${CMAKE_SOURCE_DIR}/data/points/iso_cuboid_3_in_0_1.txt" "3" "1.0") diff --git a/src/Alpha_complex/utilities/alpha_complex_3d_helper.h b/src/Alpha_complex/utilities/alpha_complex_3d_helper.h index 6b3b7d5d..a59f0654 100644 --- a/src/Alpha_complex/utilities/alpha_complex_3d_helper.h +++ b/src/Alpha_complex/utilities/alpha_complex_3d_helper.h @@ -52,13 +52,11 @@ Vertex_list from_facet(const Facet& fct) { template <class Vertex_list, class Edge_3> Vertex_list from_edge(const Edge_3& edg) { Vertex_list the_list; - for (auto i = 0; i < 4; i++) { - if ((edg.second == i) || (edg.third == i)) { + for (auto i : {edg.second, edg.third}) { #ifdef DEBUG_TRACES - std::cout << "from edge[" << i << "]=" << edg.first->vertex(i)->point() << std::endl; + std::cout << "from edge[" << i << "]=" << edg.first->vertex(i)->point() << std::endl; #endif // DEBUG_TRACES - the_list.push_back(edg.first->vertex(i)); - } + the_list.push_back(edg.first->vertex(i)); } return the_list; } diff --git a/src/Alpha_complex/utilities/alpha_complex_3d_persistence.cpp b/src/Alpha_complex/utilities/alpha_complex_3d_persistence.cpp index 0a021a0f..8ef5ffb2 100644 --- a/src/Alpha_complex/utilities/alpha_complex_3d_persistence.cpp +++ b/src/Alpha_complex/utilities/alpha_complex_3d_persistence.cpp @@ -20,9 +20,14 @@ * along with this program. If not, see <http://www.gnu.org/licenses/>. */ +#include <boost/version.hpp> #include <boost/program_options.hpp> #include <boost/variant.hpp> +#if BOOST_VERSION >= 105400 +#include <boost/container/static_vector.hpp> +#endif + #include <gudhi/Simplex_tree.h> #include <gudhi/Persistent_cohomology.h> #include <gudhi/Points_3D_off_io.h> @@ -38,7 +43,6 @@ #include <tuple> #include <map> #include <utility> -#include <list> #include <vector> #include <cstdlib> @@ -66,14 +70,18 @@ using Cell_handle = Alpha_shape_3::Cell_handle; using Facet = Alpha_shape_3::Facet; using Edge_3 = Alpha_shape_3::Edge; using Vertex_handle = Alpha_shape_3::Vertex_handle; -using Vertex_list = std::list<Alpha_shape_3::Vertex_handle>; + +#if BOOST_VERSION >= 105400 +using Vertex_list = boost::container::static_vector<Alpha_shape_3::Vertex_handle, 4>; +#else +using Vertex_list = std::vector<Alpha_shape_3::Vertex_handle>; +#endif // gudhi type definition using ST = Gudhi::Simplex_tree<Gudhi::Simplex_tree_options_fast_persistence>; using Filtration_value = ST::Filtration_value; using Simplex_tree_vertex = ST::Vertex_handle; using Alpha_shape_simplex_tree_map = std::map<Alpha_shape_3::Vertex_handle, Simplex_tree_vertex>; -using Alpha_shape_simplex_tree_pair = std::pair<Alpha_shape_3::Vertex_handle, Simplex_tree_vertex>; using Simplex_tree_vector_vertex = std::vector<Simplex_tree_vertex>; using Persistent_cohomology = Gudhi::persistent_cohomology::Persistent_cohomology<ST, Gudhi::persistent_cohomology::Field_Zp>; @@ -97,7 +105,7 @@ int main(int argc, char **argv) { exit(-1); } - // Retrieve the triangulation + // Retrieve the points std::vector<Point_3> lp = off_reader.get_point_cloud(); // alpha shape construction from points. CGAL has a strange behavior in REGULARIZED mode. @@ -128,37 +136,23 @@ int main(int argc, char **argv) { ST simplex_tree; Alpha_shape_simplex_tree_map map_cgal_simplex_tree; std::vector<Alpha_value_type>::iterator the_alpha_value_iterator = the_alpha_values.begin(); - int dim_max = 0; - Filtration_value filtration_max = 0.0; for (auto object_iterator : the_objects) { // 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)) { vertex_list = from_cell<Vertex_list, Cell_handle>(*cell); count_cells++; - if (dim_max < 3) { - // Cell is of dim 3 - dim_max = 3; - } - } else if (const Facet* facet = CGAL::object_cast<Facet>(&object_iterator)) { + } else if (const Facet *facet = CGAL::object_cast<Facet>(&object_iterator)) { vertex_list = from_facet<Vertex_list, Facet>(*facet); count_facets++; - if (dim_max < 2) { - // Facet is of dim 2 - dim_max = 2; - } - } else if (const Edge_3* edge = CGAL::object_cast<Edge_3>(&object_iterator)) { + } else if (const Edge_3 *edge = CGAL::object_cast<Edge_3>(&object_iterator)) { vertex_list = from_edge<Vertex_list, Edge_3>(*edge); count_edges++; - if (dim_max < 1) { - // Edge_3 is of dim 1 - dim_max = 1; - } - } else if (const Vertex_handle* vertex = CGAL::object_cast<Vertex_handle>(&object_iterator)) { + } else if (const Vertex_handle *vertex = CGAL::object_cast<Vertex_handle>(&object_iterator)) { count_vertices++; vertex_list = from_vertex<Vertex_list, Vertex_handle>(*vertex); } // Construction of the vector of simplex_tree vertex from list of alpha_shapes vertex - Simplex_tree_vector_vertex the_simplex_tree; + Simplex_tree_vector_vertex the_simplex; for (auto the_alpha_shape_vertex : vertex_list) { Alpha_shape_simplex_tree_map::iterator the_map_iterator = map_cgal_simplex_tree.find(the_alpha_shape_vertex); if (the_map_iterator == map_cgal_simplex_tree.end()) { @@ -167,30 +161,25 @@ int main(int argc, char **argv) { #ifdef DEBUG_TRACES std::cout << "vertex [" << the_alpha_shape_vertex->point() << "] not found - insert " << vertex << std::endl; #endif // DEBUG_TRACES - the_simplex_tree.push_back(vertex); - map_cgal_simplex_tree.insert(Alpha_shape_simplex_tree_pair(the_alpha_shape_vertex, vertex)); + the_simplex.push_back(vertex); + map_cgal_simplex_tree.emplace(the_alpha_shape_vertex, vertex); } else { // alpha shape found Simplex_tree_vertex vertex = the_map_iterator->second; #ifdef DEBUG_TRACES std::cout << "vertex [" << the_alpha_shape_vertex->point() << "] found in " << vertex << std::endl; #endif // DEBUG_TRACES - the_simplex_tree.push_back(vertex); + the_simplex.push_back(vertex); } } // Construction of the simplex_tree - Filtration_value filtr = /*std::sqrt*/(*the_alpha_value_iterator); + Filtration_value filtr = /*std::sqrt*/ (*the_alpha_value_iterator); #ifdef DEBUG_TRACES std::cout << "filtration = " << filtr << std::endl; #endif // DEBUG_TRACES - if (filtr > filtration_max) { - filtration_max = filtr; - } - simplex_tree.insert_simplex(the_simplex_tree, filtr); - if (the_alpha_value_iterator != the_alpha_values.end()) - ++the_alpha_value_iterator; - else - std::cout << "This shall not happen" << std::endl; + simplex_tree.insert_simplex(the_simplex, filtr); + GUDHI_CHECK(the_alpha_value_iterator != the_alpha_values.end(), "CGAL provided more simplices than values"); + ++the_alpha_value_iterator; } #ifdef DEBUG_TRACES diff --git a/src/Alpha_complex/utilities/README b/src/Alpha_complex/utilities/alphacomplex.md index 1cd2ca95..aace85d3 100644 --- a/src/Alpha_complex/utilities/README +++ b/src/Alpha_complex/utilities/alphacomplex.md @@ -1,14 +1,26 @@ -# Alpha_complex #
-## `alpha_complex_3d_persistence` ##
-This program computes the persistent homology with coefficient field Z/pZ of the 3D alpha complex built from a 3D point cloud. The output diagram contains one bar per line, written with the convention:
-`p dim birth death`
+# Alpha complex #
-where `dim` is the dimension of the homological feature, `birth` and `death` are respectively the birth and death of the feature, and `p` is the characteristic of the field *Z/pZ* used for homology coefficients (`p` must be a prime number).
+
+## alpha_complex_persistence ##
+
+This program computes the persistent homology with coefficient field Z/pZ of the dD alpha complex built from a dD point cloud.
+The output diagram contains one bar per line, written with the convention:
+
+```
+ p dim birth death
+```
+
+where `dim` is the dimension of the homological feature, `birth` and `death` are respectively the birth and death of the feature,
+and `p` is the characteristic of the field *Z/pZ* used for homology coefficients (`p` must be a prime number).
**Usage**
-`alpha_complex_3d_persistence [options] <input OFF file>`
+
+```
+ alpha_complex_persistence [options] <input OFF file>
+```
+
where
`<input OFF file>` is the path to the input point cloud in [nOFF ASCII format](http://www.geomview.org/docs/html/OFF.html).
@@ -16,47 +28,37 @@ where * `-h [ --help ]` Produce help message
* `-o [ --output-file ]` Name of file in which the persistence diagram is written. Default print in standard output.
-* `-p [ --field-charac ]` (default=11) Characteristic p of the coefficient field Z/pZ for computing homology.
+* `-r [ --max-alpha-square-value ]` (default = inf) Maximal alpha square value for the Alpha complex construction.
+* `-p [ --field-charac ]` (default = 11) Characteristic p of the coefficient field Z/pZ for computing homology.
* `-m [ --min-persistence ]` (default = 0) Minimal lifetime of homology feature to be recorded. Enter a negative value to see zero length intervals.
**Example**
-`alpha_complex_3d_persistence ../../data/points/tore3D_300.off -p 2 -m 0.45`
-outputs:
```
-Simplex_tree dim: 3
-2 0 0 inf
-2 1 0.0682162 1.0001
-2 1 0.0934117 1.00003
-2 2 0.56444 1.03938
+ alpha_complex_persistence -r 32 -p 2 -m 0.45 ../../data/points/tore3D_300.off
```
-Here we retrieve expected Betti numbers on a tore 3D:
-```
-Betti numbers[0] = 1
-Betti numbers[1] = 2
-Betti numbers[2] = 1
-```
+N.B.:
-N.B.:
-* `alpha_complex_3d_persistence` only accepts OFF files in dimension 3.
* Filtration values are alpha square values.
+## alpha_complex_3d_persistence ##
+This program computes the persistent homology with coefficient field Z/pZ of the 3D alpha complex built from a 3D point cloud. The output diagram contains one bar per line, written with the convention:
-## `exact_alpha_complex_3d_persistence` ##
-Same as `alpha_complex_3d_persistence`, but using exact computation. It is slower, but it is necessary when points are on a grid for instance.
+```
+p dim birth death
+```
+where `dim` is the dimension of the homological feature, `birth` and `death` are respectively the birth and death of the feature, and `p` is the characteristic of the field *Z/pZ* used for homology coefficients (`p` must be a prime number).
+**Usage**
-## `weighted_alpha_complex_3d_persistence` ##
-Same as `alpha_complex_3d_persistence`, but using weighted points.
+```
+ alpha_complex_3d_persistence [options] <input OFF file>
+```
-**Usage**
-`weighted_alpha_complex_3d_persistence [options] <input OFF file> <weights input file>`
-where
-`<input OFF file>` is the path to the input point cloud in [nOFF ASCII format](http://www.geomview.org/docs/html/OFF.html).
-`<input weights file>` is the path to the file containing the weights of the points (one value per line).
+where `<input OFF file>` is the path to the input point cloud in [nOFF ASCII format](http://www.geomview.org/docs/html/OFF.html).
**Allowed options**
@@ -66,112 +68,91 @@ where * `-m [ --min-persistence ]` (default = 0) Minimal lifetime of homology feature to be recorded. Enter a negative value to see zero length intervals.
**Example**
-`weighted_alpha_complex_3d_persistence ../../data/points/tore3D_300.off ../../data/points/tore3D_300.weights -p 2 -m 0.45`
-outputs:
```
-Simplex_tree dim: 3
-2 0 -1 inf
-2 1 -0.931784 0.000103311
-2 1 -0.906588 2.60165e-05
-2 2 -0.43556 0.0393798
+alpha_complex_3d_persistence ../../data/points/tore3D_300.off -p 2 -m 0.45
```
N.B.:
-* Weights values are explained on CGAL [Alpha shape](https://doc.cgal.org/latest/Alpha_shapes_3/index.html#title0)
-and [Regular triangulation](https://doc.cgal.org/latest/Triangulation_3/index.html#Triangulation3secclassRegulartriangulation) documentation.
+
+* `alpha_complex_3d_persistence` only accepts OFF files in dimension 3.
* Filtration values are alpha square values.
-## `periodic_alpha_complex_3d_persistence` ##
-Same as `alpha_complex_3d_persistence`, but using periodic alpha shape 3d.
-Refer to the [CGAL's 3D Periodic Triangulations User Manual](https://doc.cgal.org/latest/Periodic_3_triangulation_3/index.html) for more details.
+## exact_alpha_complex_3d_persistence ##
+
+Same as `alpha_complex_3d_persistence`, but using exact computation.
+It is slower, but it is necessary when points are on a grid for instance.
+
+
+
+## weighted_alpha_complex_3d_persistence ##
+
+Same as `alpha_complex_3d_persistence`, but using weighted points.
**Usage**
-`periodic_alpha_complex_3d_persistence [options] <input OFF file> <cuboid file>`
+
+```
+ weighted_alpha_complex_3d_persistence [options] <input OFF file> <weights input file>
+```
+
where
-`<input OFF file>` is the path to the input point cloud in [nOFF ASCII format](http://www.geomview.org/docs/html/OFF.html).
-`<cuboid file>` is the path to the file describing the periodic domain. It must be in the format described [here](http://gudhi.gforge.inria.fr/doc/latest/fileformats.html#FileFormatsIsoCuboid).
+
+* `<input OFF file>` is the path to the input point cloud in [nOFF ASCII format](http://www.geomview.org/docs/html/OFF.html).
+* `<input weights file>` is the path to the file containing the weights of the points (one value per line).
**Allowed options**
* `-h [ --help ]` Produce help message
* `-o [ --output-file ]` Name of file in which the persistence diagram is written. Default print in standard output.
* `-p [ --field-charac ]` (default=11) Characteristic p of the coefficient field Z/pZ for computing homology.
-* `-m [ --min-persistence ]` (default = 0) Minimal lifetime of homology feature to be recorded. Enter a negative value to see zero length intervals
-
+* `-m [ --min-persistence ]` (default = 0) Minimal lifetime of homology feature to be recorded. Enter a negative value to see zero length intervals.
**Example**
-`periodic_alpha_complex_3d_persistence ../../data/points/grid_10_10_10_in_0_1.off ../../data/points/iso_cuboid_3_in_0_1.txt -p 3 -m 1.0`
-
-outputs:
-```
-Periodic Delaunay computed.
-Simplex_tree dim: 3
-3 0 0 inf
-3 1 0.0025 inf
-3 1 0.0025 inf
-3 1 0.0025 inf
-3 2 0.005 inf
-3 2 0.005 inf
-3 2 0.005 inf
-3 3 0.0075 inf
-```
-Here we retrieve expected Betti numbers on an 3D iso-oriented cuboids:
```
-Betti numbers[0] = 1
-Betti numbers[1] = 3
-Betti numbers[2] = 3
-Betti numbers[3] = 1
+ weighted_alpha_complex_3d_persistence ../../data/points/tore3D_300.off ../../data/points/tore3D_300.weights -p 2 -m 0.45
```
-N.B.:
-* Cuboid file must be in the format described [here](http://gudhi.gforge.inria.fr/doc/latest/fileformats.html#FileFormatsIsoCuboid).
-* Filtration values are alpha square values.
+N.B.:
+* Weights values are explained on CGAL [Alpha shape](https://doc.cgal.org/latest/Alpha_shapes_3/index.html#title0)
+and [Regular triangulation](https://doc.cgal.org/latest/Triangulation_3/index.html#Triangulation3secclassRegulartriangulation) documentation.
+* Filtration values are alpha square values.
-## `alpha_complex_persistence` ##
-This program computes the persistent homology with coefficient field Z/pZ of the dD alpha complex built from a dD point cloud. The output diagram contains one bar per line, written with the convention:
+## periodic_alpha_complex_3d_persistence ##
+Same as `alpha_complex_3d_persistence`, but using periodic alpha shape 3d.
+Refer to the [CGAL's 3D Periodic Triangulations User Manual](https://doc.cgal.org/latest/Periodic_3_triangulation_3/index.html) for more details.
-`p dim birth death`
+**Usage**
-where `dim` is the dimension of the homological feature, `birth` and `death` are respectively the birth and death of the feature, and `p` is the characteristic of the field *Z/pZ* used for homology coefficients (`p` must be a prime number).
+```
+ periodic_alpha_complex_3d_persistence [options] <input OFF file> <cuboid file>
+```
-**Usage**
-`alpha_complex_persistence [options] <input OFF file>`
where
-`<input OFF file>` is the path to the input point cloud in [nOFF ASCII format](http://www.geomview.org/docs/html/OFF.html).
+
+* `<input OFF file>` is the path to the input point cloud in [nOFF ASCII format](http://www.geomview.org/docs/html/OFF.html).
+* `<cuboid file>` is the path to the file describing the periodic domain. It must be in the format described
+[here](/doc/latest/fileformats.html#FileFormatsIsoCuboid).
**Allowed options**
* `-h [ --help ]` Produce help message
* `-o [ --output-file ]` Name of file in which the persistence diagram is written. Default print in standard output.
-* `-r [ --max-alpha-square-value ]` (default = inf) Maximal alpha square value for the Alpha complex construction.
-* `-p [ --field-charac ]` (default = 11) Characteristic p of the coefficient field Z/pZ for computing homology.
-* `-m [ --min-persistence ]` (default = 0) Minimal lifetime of homology feature to be recorded. Enter a negative value to see zero length intervals.
+* `-p [ --field-charac ]` (default=11) Characteristic p of the coefficient field Z/pZ for computing homology.
+* `-m [ --min-persistence ]` (default = 0) Minimal lifetime of homology feature to be recorded. Enter a negative value to see zero length intervals
-**Example**
-`alpha_complex_persistence -r 32 -p 2 -m 0.45 ../../data/points/tore3D_300.off`
-outputs:
-```
-Alpha complex is of dimension 3 - 9273 simplices - 300 vertices.
-Simplex_tree dim: 3
-2 0 0 inf
-2 1 0.0682162 1.0001
-2 1 0.0934117 1.00003
-2 2 0.56444 1.03938
-```
+**Example**
-Here we retrieve expected Betti numbers on a tore 3D:
```
-Betti numbers[0] = 1
-Betti numbers[1] = 2
-Betti numbers[2] = 1
+periodic_alpha_complex_3d_persistence ../../data/points/grid_10_10_10_in_0_1.off ../../data/points/iso_cuboid_3_in_0_1.txt -p 3 -m 1.0
```
N.B.:
+
+* Cuboid file must be in the format described [here](/doc/latest/fileformats.html#FileFormatsIsoCuboid).
* Filtration values are alpha square values.
diff --git a/src/Alpha_complex/utilities/exact_alpha_complex_3d_persistence.cpp b/src/Alpha_complex/utilities/exact_alpha_complex_3d_persistence.cpp index 9a266418..cceac46e 100644 --- a/src/Alpha_complex/utilities/exact_alpha_complex_3d_persistence.cpp +++ b/src/Alpha_complex/utilities/exact_alpha_complex_3d_persistence.cpp @@ -38,7 +38,6 @@ #include <tuple> #include <map> #include <utility> -#include <list> #include <vector> #include <cstdlib> @@ -67,14 +66,13 @@ using Cell_handle = Alpha_shape_3::Cell_handle; using Facet = Alpha_shape_3::Facet; using Edge_3 = Alpha_shape_3::Edge; using Vertex_handle = Alpha_shape_3::Vertex_handle; -using Vertex_list = std::list<Vertex_handle>; +using Vertex_list = std::vector<Vertex_handle>; // gudhi type definition using ST = Gudhi::Simplex_tree<Gudhi::Simplex_tree_options_fast_persistence>; using Filtration_value = ST::Filtration_value; using Simplex_tree_vertex = ST::Vertex_handle; using Alpha_shape_simplex_tree_map = std::map<Alpha_shape_3::Vertex_handle, Simplex_tree_vertex>; -using Alpha_shape_simplex_tree_pair = std::pair<Alpha_shape_3::Vertex_handle, Simplex_tree_vertex>; using Simplex_tree_vector_vertex = std::vector<Simplex_tree_vertex>; using Persistent_cohomology = Gudhi::persistent_cohomology::Persistent_cohomology<ST, Gudhi::persistent_cohomology::Field_Zp>; @@ -98,7 +96,7 @@ int main(int argc, char **argv) { exit(-1); } - // Retrieve the triangulation + // Retrieve the points std::vector<Point_3> lp = off_reader.get_point_cloud(); // alpha shape construction from points. CGAL has a strange behavior in REGULARIZED mode. @@ -129,37 +127,23 @@ int main(int argc, char **argv) { ST simplex_tree; Alpha_shape_simplex_tree_map map_cgal_simplex_tree; std::vector<Alpha_value_type>::iterator the_alpha_value_iterator = the_alpha_values.begin(); - int dim_max = 0; - Filtration_value filtration_max = 0.0; for (auto object_iterator : the_objects) { // 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)) { vertex_list = from_cell<Vertex_list, Cell_handle>(*cell); count_cells++; - if (dim_max < 3) { - // Cell is of dim 3 - dim_max = 3; - } - } else if (const Facet* facet = CGAL::object_cast<Facet>(&object_iterator)) { + } else if (const Facet *facet = CGAL::object_cast<Facet>(&object_iterator)) { vertex_list = from_facet<Vertex_list, Facet>(*facet); count_facets++; - if (dim_max < 2) { - // Facet is of dim 2 - dim_max = 2; - } - } else if (const Edge_3* edge = CGAL::object_cast<Edge_3>(&object_iterator)) { + } else if (const Edge_3 *edge = CGAL::object_cast<Edge_3>(&object_iterator)) { vertex_list = from_edge<Vertex_list, Edge_3>(*edge); count_edges++; - if (dim_max < 1) { - // Edge_3 is of dim 1 - dim_max = 1; - } - } else if (const Vertex_handle* vertex = CGAL::object_cast<Vertex_handle>(&object_iterator)) { + } else if (const Vertex_handle *vertex = CGAL::object_cast<Vertex_handle>(&object_iterator)) { count_vertices++; vertex_list = from_vertex<Vertex_list, Vertex_handle>(*vertex); } // Construction of the vector of simplex_tree vertex from list of alpha_shapes vertex - Simplex_tree_vector_vertex the_simplex_tree; + Simplex_tree_vector_vertex the_simplex; for (auto the_alpha_shape_vertex : vertex_list) { Alpha_shape_simplex_tree_map::iterator the_map_iterator = map_cgal_simplex_tree.find(the_alpha_shape_vertex); if (the_map_iterator == map_cgal_simplex_tree.end()) { @@ -168,27 +152,24 @@ int main(int argc, char **argv) { #ifdef DEBUG_TRACES std::cout << "vertex [" << the_alpha_shape_vertex->point() << "] not found - insert " << vertex << std::endl; #endif // DEBUG_TRACES - the_simplex_tree.push_back(vertex); - map_cgal_simplex_tree.insert(Alpha_shape_simplex_tree_pair(the_alpha_shape_vertex, vertex)); + the_simplex.push_back(vertex); + map_cgal_simplex_tree.emplace(the_alpha_shape_vertex, vertex); } else { // alpha shape found Simplex_tree_vertex vertex = the_map_iterator->second; #ifdef DEBUG_TRACES std::cout << "vertex [" << the_alpha_shape_vertex->point() << "] found in " << vertex << std::endl; #endif // DEBUG_TRACES - the_simplex_tree.push_back(vertex); + the_simplex.push_back(vertex); } } // Construction of the simplex_tree // you can also use the_alpha_value_iterator->exact() - Filtration_value filtr = /*std::sqrt*/CGAL::to_double(the_alpha_value_iterator->exact()); + Filtration_value filtr = /*std::sqrt*/ CGAL::to_double(the_alpha_value_iterator->exact()); #ifdef DEBUG_TRACES std::cout << "filtration = " << filtr << std::endl; #endif // DEBUG_TRACES - if (filtr > filtration_max) { - filtration_max = filtr; - } - simplex_tree.insert_simplex(the_simplex_tree, filtr); + simplex_tree.insert_simplex(the_simplex, filtr); if (the_alpha_value_iterator != the_alpha_values.end()) ++the_alpha_value_iterator; else diff --git a/src/Alpha_complex/utilities/periodic_alpha_complex_3d_persistence.cpp b/src/Alpha_complex/utilities/periodic_alpha_complex_3d_persistence.cpp index 186a58f8..188cf604 100644 --- a/src/Alpha_complex/utilities/periodic_alpha_complex_3d_persistence.cpp +++ b/src/Alpha_complex/utilities/periodic_alpha_complex_3d_persistence.cpp @@ -3,6 +3,7 @@ * library for computational topology. * * Author(s): Vincent Rouvreau + * Pawel Dlotko - 2017 - Swansea University, UK * * Copyright (C) 2014 INRIA * @@ -39,7 +40,6 @@ #include <tuple> #include <map> #include <utility> -#include <list> #include <vector> #include <cstdlib> @@ -72,14 +72,13 @@ using Cell_handle = Alpha_shape_3::Cell_handle; using Facet = Alpha_shape_3::Facet; using Edge_3 = Alpha_shape_3::Edge; using Vertex_handle = Alpha_shape_3::Vertex_handle; -using Vertex_list = std::list<Alpha_shape_3::Vertex_handle>; +using Vertex_list = std::vector<Alpha_shape_3::Vertex_handle>; // gudhi type definition using ST = Gudhi::Simplex_tree<Gudhi::Simplex_tree_options_fast_persistence>; using Filtration_value = ST::Filtration_value; using Simplex_tree_vertex = ST::Vertex_handle; using Alpha_shape_simplex_tree_map = std::map<Alpha_shape_3::Vertex_handle, Simplex_tree_vertex>; -using Alpha_shape_simplex_tree_pair = std::pair<Alpha_shape_3::Vertex_handle, Simplex_tree_vertex>; using Simplex_tree_vector_vertex = std::vector<Simplex_tree_vertex>; using Persistent_cohomology = Gudhi::persistent_cohomology::Persistent_cohomology<ST, Gudhi::persistent_cohomology::Field_Zp>; @@ -114,8 +113,13 @@ int main(int argc, char **argv) { std::cerr << "Unable to read file " << cuboid_file << std::endl; exit(-1); } + // Checking if the cuboid is the same in x,y and z direction. If not, CGAL will not process it. + if ((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::cerr << "The size of the cuboid in every directions is not the same." << std::endl; + exit(-1); + } - // Retrieve the triangulation + // Retrieve the points std::vector<Point_3> lp = off_reader.get_point_cloud(); // Define the periodic cube @@ -123,7 +127,12 @@ int main(int argc, char **argv) { // Heuristic for inserting large point sets (if pts is reasonably large) pdt.insert(lp.begin(), lp.end(), true); // As pdt won't be modified anymore switch to 1-sheeted cover if possible - if (pdt.is_triangulation_in_1_sheet()) pdt.convert_to_1_sheeted_covering(); + if (pdt.is_triangulation_in_1_sheet()) { + pdt.convert_to_1_sheeted_covering(); + } else { + std::cerr << "ERROR: we were not able to construct a triangulation within a single periodic domain." << std::endl; + exit(-1); + } std::cout << "Periodic Delaunay computed." << std::endl; // alpha shape construction from points. CGAL has a strange behavior in REGULARIZED mode. This is the default mode @@ -152,37 +161,23 @@ int main(int argc, char **argv) { ST simplex_tree; Alpha_shape_simplex_tree_map map_cgal_simplex_tree; std::vector<Alpha_value_type>::iterator the_alpha_value_iterator = the_alpha_values.begin(); - int dim_max = 0; - Filtration_value filtration_max = 0.0; for (auto object_iterator : the_objects) { // 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)) { vertex_list = from_cell<Vertex_list, Cell_handle>(*cell); count_cells++; - if (dim_max < 3) { - // Cell is of dim 3 - dim_max = 3; - } - } else if (const Facet* facet = CGAL::object_cast<Facet>(&object_iterator)) { + } else if (const Facet *facet = CGAL::object_cast<Facet>(&object_iterator)) { vertex_list = from_facet<Vertex_list, Facet>(*facet); count_facets++; - if (dim_max < 2) { - // Facet is of dim 2 - dim_max = 2; - } - } else if (const Edge_3* edge = CGAL::object_cast<Edge_3>(&object_iterator)) { + } else if (const Edge_3 *edge = CGAL::object_cast<Edge_3>(&object_iterator)) { vertex_list = from_edge<Vertex_list, Edge_3>(*edge); count_edges++; - if (dim_max < 1) { - // Edge_3 is of dim 1 - dim_max = 1; - } - } else if (const Vertex_handle* vertex = CGAL::object_cast<Vertex_handle>(&object_iterator)) { + } else if (const Vertex_handle *vertex = CGAL::object_cast<Vertex_handle>(&object_iterator)) { count_vertices++; vertex_list = from_vertex<Vertex_list, Vertex_handle>(*vertex); } // Construction of the vector of simplex_tree vertex from list of alpha_shapes vertex - Simplex_tree_vector_vertex the_simplex_tree; + Simplex_tree_vector_vertex the_simplex; for (auto the_alpha_shape_vertex : vertex_list) { Alpha_shape_simplex_tree_map::iterator the_map_iterator = map_cgal_simplex_tree.find(the_alpha_shape_vertex); if (the_map_iterator == map_cgal_simplex_tree.end()) { @@ -191,15 +186,15 @@ int main(int argc, char **argv) { #ifdef DEBUG_TRACES std::cout << "vertex [" << the_alpha_shape_vertex->point() << "] not found - insert " << vertex << std::endl; #endif // DEBUG_TRACES - the_simplex_tree.push_back(vertex); - map_cgal_simplex_tree.insert(Alpha_shape_simplex_tree_pair(the_alpha_shape_vertex, vertex)); + the_simplex.push_back(vertex); + map_cgal_simplex_tree.emplace(the_alpha_shape_vertex, vertex); } else { // alpha shape found Simplex_tree_vertex vertex = the_map_iterator->second; #ifdef DEBUG_TRACES std::cout << "vertex [" << the_alpha_shape_vertex->point() << "] found in " << vertex << std::endl; #endif // DEBUG_TRACES - the_simplex_tree.push_back(vertex); + the_simplex.push_back(vertex); } } // Construction of the simplex_tree @@ -207,10 +202,7 @@ int main(int argc, char **argv) { #ifdef DEBUG_TRACES std::cout << "filtration = " << filtr << std::endl; #endif // DEBUG_TRACES - if (filtr > filtration_max) { - filtration_max = filtr; - } - simplex_tree.insert_simplex(the_simplex_tree, filtr); + simplex_tree.insert_simplex(the_simplex, filtr); if (the_alpha_value_iterator != the_alpha_values.end()) ++the_alpha_value_iterator; else diff --git a/src/Alpha_complex/utilities/weighted_alpha_complex_3d_persistence.cpp b/src/Alpha_complex/utilities/weighted_alpha_complex_3d_persistence.cpp index 0e73a99b..93be8a05 100644 --- a/src/Alpha_complex/utilities/weighted_alpha_complex_3d_persistence.cpp +++ b/src/Alpha_complex/utilities/weighted_alpha_complex_3d_persistence.cpp @@ -44,7 +44,6 @@ #include <tuple> #include <map> #include <utility> -#include <list> #include <vector> #include <cstdlib> @@ -66,13 +65,13 @@ using Point_3 = Gt::Bare_point; using Weighted_point_3 = Gt::Weighted_point; // For CGAL >= 4.11 -#else // CGAL_VERSION_NR < 1041100000 +#else // CGAL_VERSION_NR < 1041100000 using Rvb = CGAL::Regular_triangulation_vertex_base_3<Kernel>; -using Vb = CGAL::Alpha_shape_vertex_base_3<Kernel,Rvb>; +using Vb = CGAL::Alpha_shape_vertex_base_3<Kernel, Rvb>; using Rcb = CGAL::Regular_triangulation_cell_base_3<Kernel>; -using Cb = CGAL::Alpha_shape_cell_base_3<Kernel,Rcb>; -using Tds = CGAL::Triangulation_data_structure_3<Vb,Cb>; -using Triangulation_3 = CGAL::Regular_triangulation_3<Kernel,Tds>; +using Cb = CGAL::Alpha_shape_cell_base_3<Kernel, Rcb>; +using Tds = CGAL::Triangulation_data_structure_3<Vb, Cb>; +using Triangulation_3 = CGAL::Regular_triangulation_3<Kernel, Tds>; // From file type definition using Point_3 = Triangulation_3::Bare_point; @@ -92,14 +91,13 @@ using Cell_handle = Alpha_shape_3::Cell_handle; using Facet = Alpha_shape_3::Facet; using Edge_3 = Alpha_shape_3::Edge; using Vertex_handle = Alpha_shape_3::Vertex_handle; -using Vertex_list = std::list<Alpha_shape_3::Vertex_handle>; +using Vertex_list = std::vector<Alpha_shape_3::Vertex_handle>; // gudhi type definition using ST = Gudhi::Simplex_tree<Gudhi::Simplex_tree_options_fast_persistence>; using Filtration_value = ST::Filtration_value; using Simplex_tree_vertex = ST::Vertex_handle; using Alpha_shape_simplex_tree_map = std::map<Alpha_shape_3::Vertex_handle, Simplex_tree_vertex>; -using Alpha_shape_simplex_tree_pair = std::pair<Alpha_shape_3::Vertex_handle, Simplex_tree_vertex>; using Simplex_tree_vector_vertex = std::vector<Simplex_tree_vertex>; using Persistent_cohomology = Gudhi::persistent_cohomology::Persistent_cohomology<ST, Gudhi::persistent_cohomology::Field_Zp>; @@ -125,7 +123,7 @@ int main(int argc, char **argv) { exit(-1); } - // Retrieve the triangulation + // Retrieve the points std::vector<Point_3> lp = off_reader.get_point_cloud(); // Read weights information from file @@ -177,37 +175,23 @@ int main(int argc, char **argv) { ST simplex_tree; Alpha_shape_simplex_tree_map map_cgal_simplex_tree; std::vector<Alpha_value_type>::iterator the_alpha_value_iterator = the_alpha_values.begin(); - int dim_max = 0; - Filtration_value filtration_max = 0.0; for (auto object_iterator : the_objects) { // 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)) { vertex_list = from_cell<Vertex_list, Cell_handle>(*cell); count_cells++; - if (dim_max < 3) { - // Cell is of dim 3 - dim_max = 3; - } - } else if (const Facet* facet = CGAL::object_cast<Facet>(&object_iterator)) { + } else if (const Facet *facet = CGAL::object_cast<Facet>(&object_iterator)) { vertex_list = from_facet<Vertex_list, Facet>(*facet); count_facets++; - if (dim_max < 2) { - // Facet is of dim 2 - dim_max = 2; - } - } else if (const Edge_3* edge = CGAL::object_cast<Edge_3>(&object_iterator)) { + } else if (const Edge_3 *edge = CGAL::object_cast<Edge_3>(&object_iterator)) { vertex_list = from_edge<Vertex_list, Edge_3>(*edge); count_edges++; - if (dim_max < 1) { - // Edge_3 is of dim 1 - dim_max = 1; - } - } else if (const Vertex_handle* vertex = CGAL::object_cast<Vertex_handle>(&object_iterator)) { + } else if (const Vertex_handle *vertex = CGAL::object_cast<Vertex_handle>(&object_iterator)) { count_vertices++; vertex_list = from_vertex<Vertex_list, Vertex_handle>(*vertex); } // Construction of the vector of simplex_tree vertex from list of alpha_shapes vertex - Simplex_tree_vector_vertex the_simplex_tree; + Simplex_tree_vector_vertex the_simplex; for (auto the_alpha_shape_vertex : vertex_list) { Alpha_shape_simplex_tree_map::iterator the_map_iterator = map_cgal_simplex_tree.find(the_alpha_shape_vertex); if (the_map_iterator == map_cgal_simplex_tree.end()) { @@ -216,15 +200,15 @@ int main(int argc, char **argv) { #ifdef DEBUG_TRACES std::cout << "vertex [" << the_alpha_shape_vertex->point() << "] not found - insert " << vertex << std::endl; #endif // DEBUG_TRACES - the_simplex_tree.push_back(vertex); - map_cgal_simplex_tree.insert(Alpha_shape_simplex_tree_pair(the_alpha_shape_vertex, vertex)); + the_simplex.push_back(vertex); + map_cgal_simplex_tree.emplace(the_alpha_shape_vertex, vertex); } else { // alpha shape found Simplex_tree_vertex vertex = the_map_iterator->second; #ifdef DEBUG_TRACES std::cout << "vertex [" << the_alpha_shape_vertex->point() << "] found in " << vertex << std::endl; #endif // DEBUG_TRACES - the_simplex_tree.push_back(vertex); + the_simplex.push_back(vertex); } } // Construction of the simplex_tree @@ -232,10 +216,7 @@ int main(int argc, char **argv) { #ifdef DEBUG_TRACES std::cout << "filtration = " << filtr << std::endl; #endif // DEBUG_TRACES - if (filtr > filtration_max) { - filtration_max = filtr; - } - simplex_tree.insert_simplex(the_simplex_tree, filtr); + simplex_tree.insert_simplex(the_simplex, filtr); if (the_alpha_value_iterator != the_alpha_values.end()) ++the_alpha_value_iterator; else diff --git a/src/Alpha_complex/utilities/weighted_periodic_alpha_complex_3d_persistence.cpp b/src/Alpha_complex/utilities/weighted_periodic_alpha_complex_3d_persistence.cpp index 13634ff7..5321bb0a 100644 --- a/src/Alpha_complex/utilities/weighted_periodic_alpha_complex_3d_persistence.cpp +++ b/src/Alpha_complex/utilities/weighted_periodic_alpha_complex_3d_persistence.cpp @@ -3,6 +3,7 @@ * library for computational topology. * * Author(s): Vincent Rouvreau + * Pawel Dlotko - 2017 - Swansea University, UK * * Copyright (C) 2014 INRIA * @@ -38,7 +39,6 @@ #include <tuple> #include <map> #include <utility> -#include <list> #include <vector> #include <cstdlib> @@ -50,14 +50,14 @@ using PK = CGAL::Periodic_3_regular_triangulation_traits_3<Kernel>; // Vertex type using DsVb = CGAL::Periodic_3_triangulation_ds_vertex_base_3<>; -using Vb = CGAL::Regular_triangulation_vertex_base_3<PK,DsVb>; -using AsVb = CGAL::Alpha_shape_vertex_base_3<PK,Vb>; +using Vb = CGAL::Regular_triangulation_vertex_base_3<PK, DsVb>; +using AsVb = CGAL::Alpha_shape_vertex_base_3<PK, Vb>; // Cell type using DsCb = CGAL::Periodic_3_triangulation_ds_cell_base_3<>; -using Cb = CGAL::Regular_triangulation_cell_base_3<PK,DsCb>; -using AsCb = CGAL::Alpha_shape_cell_base_3<PK,Cb>; -using Tds = CGAL::Triangulation_data_structure_3<AsVb,AsCb>; -using P3RT3 = CGAL::Periodic_3_regular_triangulation_3<PK,Tds>; +using Cb = CGAL::Regular_triangulation_cell_base_3<PK, DsCb>; +using AsCb = CGAL::Alpha_shape_cell_base_3<PK, Cb>; +using Tds = CGAL::Triangulation_data_structure_3<AsVb, AsCb>; +using P3RT3 = CGAL::Periodic_3_regular_triangulation_3<PK, Tds>; using Alpha_shape_3 = CGAL::Alpha_shape_3<P3RT3>; using Point_3 = P3RT3::Bare_point; @@ -74,14 +74,13 @@ using Cell_handle = Alpha_shape_3::Cell_handle; using Facet = Alpha_shape_3::Facet; using Edge_3 = Alpha_shape_3::Edge; using Vertex_handle = Alpha_shape_3::Vertex_handle; -using Vertex_list = std::list<Alpha_shape_3::Vertex_handle>; +using Vertex_list = std::vector<Alpha_shape_3::Vertex_handle>; // gudhi type definition using ST = Gudhi::Simplex_tree<Gudhi::Simplex_tree_options_fast_persistence>; using Filtration_value = ST::Filtration_value; using Simplex_tree_vertex = ST::Vertex_handle; using Alpha_shape_simplex_tree_map = std::map<Alpha_shape_3::Vertex_handle, Simplex_tree_vertex>; -using Alpha_shape_simplex_tree_pair = std::pair<Alpha_shape_3::Vertex_handle, Simplex_tree_vertex>; using Simplex_tree_vector_vertex = std::vector<Simplex_tree_vertex>; using Persistent_cohomology = Gudhi::persistent_cohomology::Persistent_cohomology<ST, Gudhi::persistent_cohomology::Field_Zp>; @@ -112,18 +111,46 @@ int main(int argc, char* const argv[]) { usage(argv[0]); } - // Retrieve the triangulation + // Retrieve the points std::vector<Point_3> lp = off_reader.get_point_cloud(); + // Read iso_cuboid_3 information from file + std::ifstream iso_cuboid_str(argv[3]); + double x_min, y_min, z_min, x_max, y_max, z_max; + if (iso_cuboid_str.is_open()) { + if (!(iso_cuboid_str >> x_min >> y_min >> z_min >> x_max >> y_max >> z_max)) { + std::cerr << argv[3] << " - Bad file format." << std::endl; + usage(argv[0]); + } + + } else { + std::cerr << "Unable to read file " << argv[3] << std::endl; + usage(argv[0]); + } + // Checking if the cuboid is the same in x,y and z direction. If not, CGAL will not process it. + if ((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::cerr << "The size of the cuboid in every directions is not the same." << std::endl; + exit(-1); + } + + double maximal_possible_weight = 0.015625 * (x_max - x_min) * (x_max - x_min); + // Read weights information from file std::ifstream weights_ifstr(argv[2]); std::vector<Weighted_point_3> wp; - if (weights_ifstr.good()) { + if (weights_ifstr.is_open()) { double weight = 0.0; std::size_t index = 0; wp.reserve(lp.size()); // Attempt read the weight in a double format, return false if it fails while ((weights_ifstr >> weight) && (index < lp.size())) { + if ((weight >= maximal_possible_weight) || (weight < 0)) { + std::cerr << "At line " << (index + 1) << ", the weight (" << weight + << ") is negative or more than or equal to maximal possible weight (" << maximal_possible_weight + << ") = 1/64*cuboid length squared, which is not an acceptable input." << std::endl; + exit(-1); + } + wp.push_back(Weighted_point_3(lp[index], weight)); index++; } @@ -136,23 +163,18 @@ int main(int argc, char* const argv[]) { usage(argv[0]); } - // Read iso_cuboid_3 information from file - std::ifstream iso_cuboid_str(argv[3]); - double x_min, y_min, z_min, x_max, y_max, z_max; - if (iso_cuboid_str.good()) { - iso_cuboid_str >> x_min >> y_min >> z_min >> x_max >> y_max >> z_max; - } else { - std::cerr << "Unable to read file " << argv[3] << std::endl; - usage(argv[0]); - } - // Define the periodic cube P3RT3 prt(PK::Iso_cuboid_3(x_min, y_min, z_min, x_max, y_max, z_max)); // Heuristic for inserting large point sets (if pts is reasonably large) prt.insert(wp.begin(), wp.end(), true); // As prt won't be modified anymore switch to 1-sheeted cover if possible - if (prt.is_triangulation_in_1_sheet()) prt.convert_to_1_sheeted_covering(); - std::cout << "Periodic Delaunay computed." << std::endl; + if (prt.is_triangulation_in_1_sheet()) { + prt.convert_to_1_sheeted_covering(); + } else { + std::cerr << "ERROR: we were not able to construct a triangulation within a single periodic domain." << std::endl; + exit(-1); + } + std::cout << "Weighted Periodic Delaunay computed." << std::endl; // 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 @@ -180,37 +202,23 @@ int main(int argc, char* const argv[]) { ST simplex_tree; Alpha_shape_simplex_tree_map map_cgal_simplex_tree; std::vector<Alpha_value_type>::iterator the_alpha_value_iterator = the_alpha_values.begin(); - int dim_max = 0; - Filtration_value filtration_max = 0.0; for (auto object_iterator : the_objects) { // Retrieve Alpha shape vertex list from object if (const Cell_handle* cell = CGAL::object_cast<Cell_handle>(&object_iterator)) { vertex_list = from_cell<Vertex_list, Cell_handle>(*cell); count_cells++; - if (dim_max < 3) { - // Cell is of dim 3 - dim_max = 3; - } } else if (const Facet* facet = CGAL::object_cast<Facet>(&object_iterator)) { vertex_list = from_facet<Vertex_list, Facet>(*facet); count_facets++; - if (dim_max < 2) { - // Facet is of dim 2 - dim_max = 2; - } } else if (const Edge_3* edge = CGAL::object_cast<Edge_3>(&object_iterator)) { vertex_list = from_edge<Vertex_list, Edge_3>(*edge); count_edges++; - if (dim_max < 1) { - // Edge_3 is of dim 1 - dim_max = 1; - } } else if (const Vertex_handle* vertex = CGAL::object_cast<Vertex_handle>(&object_iterator)) { count_vertices++; vertex_list = from_vertex<Vertex_list, Vertex_handle>(*vertex); } // Construction of the vector of simplex_tree vertex from list of alpha_shapes vertex - Simplex_tree_vector_vertex the_simplex_tree; + Simplex_tree_vector_vertex the_simplex; for (auto the_alpha_shape_vertex : vertex_list) { Alpha_shape_simplex_tree_map::iterator the_map_iterator = map_cgal_simplex_tree.find(the_alpha_shape_vertex); if (the_map_iterator == map_cgal_simplex_tree.end()) { @@ -219,15 +227,15 @@ int main(int argc, char* const argv[]) { #ifdef DEBUG_TRACES std::cout << "vertex [" << the_alpha_shape_vertex->point() << "] not found - insert " << vertex << std::endl; #endif // DEBUG_TRACES - the_simplex_tree.push_back(vertex); - map_cgal_simplex_tree.insert(Alpha_shape_simplex_tree_pair(the_alpha_shape_vertex, vertex)); + the_simplex.push_back(vertex); + map_cgal_simplex_tree.emplace(the_alpha_shape_vertex, vertex); } else { // alpha shape found Simplex_tree_vertex vertex = the_map_iterator->second; #ifdef DEBUG_TRACES std::cout << "vertex [" << the_alpha_shape_vertex->point() << "] found in " << vertex << std::endl; #endif // DEBUG_TRACES - the_simplex_tree.push_back(vertex); + the_simplex.push_back(vertex); } } // Construction of the simplex_tree @@ -235,10 +243,7 @@ int main(int argc, char* const argv[]) { #ifdef DEBUG_TRACES std::cout << "filtration = " << filtr << std::endl; #endif // DEBUG_TRACES - if (filtr > filtration_max) { - filtration_max = filtr; - } - simplex_tree.insert_simplex(the_simplex_tree, filtr); + simplex_tree.insert_simplex(the_simplex, filtr); if (the_alpha_value_iterator != the_alpha_values.end()) ++the_alpha_value_iterator; else |