summaryrefslogtreecommitdiff
path: root/src/common/include
diff options
context:
space:
mode:
Diffstat (limited to 'src/common/include')
-rw-r--r--src/common/include/gudhi/Clock.h77
-rw-r--r--src/common/include/gudhi/Debug_utils.h45
-rw-r--r--src/common/include/gudhi/Null_output_iterator.h36
-rw-r--r--src/common/include/gudhi/Off_reader.h173
-rw-r--r--src/common/include/gudhi/Point.h157
-rw-r--r--src/common/include/gudhi/Points_3D_off_io.h190
-rw-r--r--src/common/include/gudhi/Points_off_io.h171
-rw-r--r--src/common/include/gudhi/Simple_object_pool.h69
-rw-r--r--src/common/include/gudhi/Unitary_tests_utils.h28
-rw-r--r--src/common/include/gudhi/allocator.h43
-rw-r--r--src/common/include/gudhi/console_color.h85
-rw-r--r--src/common/include/gudhi/distance_functions.h111
-rw-r--r--src/common/include/gudhi/graph_simplicial_complex.h99
-rw-r--r--src/common/include/gudhi/random_point_generators.h502
-rw-r--r--src/common/include/gudhi/reader_utils.h357
-rw-r--r--src/common/include/gudhi/writing_persistence_to_file.h105
16 files changed, 2248 insertions, 0 deletions
diff --git a/src/common/include/gudhi/Clock.h b/src/common/include/gudhi/Clock.h
new file mode 100644
index 00000000..00ab2f27
--- /dev/null
+++ b/src/common/include/gudhi/Clock.h
@@ -0,0 +1,77 @@
+/* This file is part of the Gudhi Library - https://gudhi.inria.fr/ - which is released under MIT.
+ * See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details.
+ * Author(s): David Salinas
+ *
+ * Copyright (C) 2014 Inria
+ *
+ * Modification(s):
+ * - YYYY/MM Author: Description of the modification
+ */
+
+#ifndef CLOCK_H_
+#define CLOCK_H_
+
+#include <iostream>
+#include <string>
+#include <chrono>
+
+namespace Gudhi {
+
+class Clock {
+ public:
+ // Construct and start the timer
+ Clock(const std::string& msg_ = std::string())
+ : startTime(std::chrono::system_clock::now()),
+ end_called(false),
+ msg(msg_) { }
+
+ // Restart the timer
+ void begin() const {
+ end_called = false;
+ startTime = std::chrono::system_clock::now();
+ }
+
+ // Stop the timer
+ void end() const {
+ end_called = true;
+ endTime = std::chrono::system_clock::now();
+ }
+
+ std::string message() const {
+ return msg;
+ }
+
+ // Print current value to std::cout
+ void print() const {
+ std::cout << *this << std::endl;
+ }
+
+ friend std::ostream& operator<<(std::ostream& stream, const Clock& clock) {
+ if (!clock.msg.empty())
+ stream << clock.msg << ": ";
+
+ stream << clock.num_seconds() << "s\n";
+ return stream;
+ }
+
+ // Get the number of seconds between the timer start and:
+ // - the last call of end() if it was called
+ // - or now otherwise. In this case, the timer is not stopped.
+ double num_seconds() const {
+ if (!end_called) {
+ auto end = std::chrono::system_clock::now();
+ return std::chrono::duration_cast<std::chrono::milliseconds>(end-startTime).count() / 1000.;
+ } else {
+ return std::chrono::duration_cast<std::chrono::milliseconds>(endTime-startTime).count() / 1000.;
+ }
+ }
+
+ private:
+ mutable std::chrono::time_point<std::chrono::system_clock> startTime, endTime;
+ mutable bool end_called;
+ std::string msg;
+};
+
+} // namespace Gudhi
+
+#endif // CLOCK_H_
diff --git a/src/common/include/gudhi/Debug_utils.h b/src/common/include/gudhi/Debug_utils.h
new file mode 100644
index 00000000..38abc06d
--- /dev/null
+++ b/src/common/include/gudhi/Debug_utils.h
@@ -0,0 +1,45 @@
+/* This file is part of the Gudhi Library - https://gudhi.inria.fr/ - which is released under MIT.
+ * See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details.
+ * Author(s): David Salinas
+ *
+ * Copyright (C) 2014 Inria
+ *
+ * Modification(s):
+ * - YYYY/MM Author: Description of the modification
+ */
+#ifndef DEBUG_UTILS_H_
+#define DEBUG_UTILS_H_
+
+#include <iostream>
+
+#ifndef NDEBUG
+ // GUDHI_DEBUG is the Gudhi official flag for debug mode.
+ #define GUDHI_DEBUG
+#endif
+
+// GUDHI_CHECK throw an exception if expression is false in debug mode, but does nothing in release mode
+// Could assert in release mode, but cmake sets NDEBUG (for "NO DEBUG") in this mode, means assert does nothing.
+#ifdef GUDHI_DEBUG
+ #define GUDHI_CHECK(expression, excpt) ((expression) ? (void) 0 : (throw excpt))
+ #define GUDHI_CHECK_code(CODE) CODE
+#else
+ #define GUDHI_CHECK(expression, excpt) (void) 0
+ #define GUDHI_CHECK_code(CODE)
+#endif
+
+#define PRINT(a) std::cerr << #a << ": " << (a) << " (DISP)" << std::endl
+
+// #define DBG_VERBOSE
+#ifdef DBG_VERBOSE
+ #define DBG(a) std::cout << "DBG: " << (a) << std::endl
+ #define DBGMSG(a, b) std::cout << "DBG: " << a << b << std::endl
+ #define DBGVALUE(a) std::cout << "DBG: " << #a << ": " << a << std::endl
+ #define DBGCONT(a) std::cout << "DBG: container " << #a << " -> "; for (auto x : a) std::cout << x << ","; std::cout << std::endl
+#else
+ #define DBG(a) (void) 0
+ #define DBGMSG(a, b) (void) 0
+ #define DBGVALUE(a) (void) 0
+ #define DBGCONT(a) (void) 0
+#endif
+
+#endif // DEBUG_UTILS_H_
diff --git a/src/common/include/gudhi/Null_output_iterator.h b/src/common/include/gudhi/Null_output_iterator.h
new file mode 100644
index 00000000..3d03bca6
--- /dev/null
+++ b/src/common/include/gudhi/Null_output_iterator.h
@@ -0,0 +1,36 @@
+/* This file is part of the Gudhi Library - https://gudhi.inria.fr/ - which is released under MIT.
+ * See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details.
+ * Author(s): Marc Glisse
+ *
+ * Copyright (C) 2017 Inria
+ *
+ * Modification(s):
+ * - YYYY/MM Author: Description of the modification
+ */
+
+#ifndef NULL_OUTPUT_ITERATOR_H_
+#define NULL_OUTPUT_ITERATOR_H_
+
+#include <iterator>
+
+namespace Gudhi {
+
+/** An output iterator that ignores whatever it is given. */
+struct Null_output_iterator {
+ typedef std::output_iterator_tag iterator_category;
+ typedef void value_type;
+ typedef void difference_type;
+ typedef void pointer;
+ typedef void reference;
+
+ Null_output_iterator& operator++() {return *this;}
+ Null_output_iterator operator++(int) {return *this;}
+ struct proxy {
+ template<class T>
+ proxy& operator=(T&&){return *this;}
+ };
+ proxy operator*()const{return {};}
+};
+} // namespace Gudhi
+
+#endif // NULL_OUTPUT_ITERATOR_H_
diff --git a/src/common/include/gudhi/Off_reader.h b/src/common/include/gudhi/Off_reader.h
new file mode 100644
index 00000000..aaff95b8
--- /dev/null
+++ b/src/common/include/gudhi/Off_reader.h
@@ -0,0 +1,173 @@
+/* This file is part of the Gudhi Library - https://gudhi.inria.fr/ - which is released under MIT.
+ * See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details.
+ * Author(s): David Salinas
+ *
+ * Copyright (C) 2014 Inria
+ *
+ * Modification(s):
+ * - YYYY/MM Author: Description of the modification
+ */
+
+
+#ifndef OFF_READER_H_
+#define OFF_READER_H_
+
+
+#include <sstream>
+#include <iostream>
+#include <iterator>
+#include <string>
+#include <vector>
+#include <fstream>
+
+namespace Gudhi {
+
+/** \brief OFF file reader top class visitor.
+ *
+ * OFF file must be conform to \ref FileFormatsOFF
+ */
+class Off_reader {
+ public:
+ Off_reader(std::ifstream& stream) : stream_(stream) { }
+
+ ~Off_reader() {
+ stream_.close();
+ }
+
+ /** \brief
+ * Read an OFF file and calls the following methods :
+ *
+ * <CODE>void init(int dim,int num_vertices,int num_faces,int num_edges); // from file header - num_edges may not be set
+ *
+ * void point(const std::vector<double>& point); // for each point read
+ *
+ * void maximal_face(const std::list<int>& face); // for each face read
+ *
+ * void done(); // upon file read is finished</CODE>
+ *
+ * of the visitor when reading a point or a maximal face. Edges are not taken into account.
+ */
+ template<typename OffVisitor>
+ bool read(OffVisitor& off_visitor) {
+ bool success_read_off_preambule = read_off_preambule(off_visitor);
+ if (!success_read_off_preambule) {
+ std::cerr << "could not read off preambule\n";
+ return false;
+ }
+
+ bool success_read_off_points = read_off_points(off_visitor);
+ if (!success_read_off_points) {
+ std::cerr << "could not read off points\n";
+ return false;
+ }
+
+ bool success_read_off_faces = read_off_faces(off_visitor);
+ if (!success_read_off_faces) {
+ std::cerr << "could not read off faces\n";
+ return false;
+ }
+
+ off_visitor.done();
+ return success_read_off_preambule && success_read_off_points && success_read_off_faces;
+ }
+
+ private:
+ std::ifstream& stream_;
+
+ struct Off_info {
+ int dim;
+ int num_vertices;
+ int num_edges;
+ int num_faces;
+ };
+
+ Off_info off_info_;
+
+ template<typename OffVisitor>
+ bool read_off_preambule(OffVisitor& off_visitor) {
+ std::string line;
+ if (!goto_next_uncomment_line(line)) return false;
+
+ bool is_off_file = (line.find("OFF") != std::string::npos);
+ bool is_noff_file = (line.find("nOFF") != std::string::npos);
+
+
+
+ if (!is_off_file && !is_noff_file) {
+ std::cerr << line << std::endl;
+ std::cerr << "missing off header\n";
+ return false;
+ }
+
+ if (is_noff_file) {
+ // Should be on a separate line, but we accept it on the same line as the number of vertices
+ stream_ >> off_info_.dim;
+ } else {
+ off_info_.dim = 3;
+ }
+
+ if (!goto_next_uncomment_line(line)) return false;
+ std::istringstream iss(line);
+ if (!(iss >> off_info_.num_vertices >> off_info_.num_faces >> off_info_.num_edges)) {
+ std::cerr << "incorrect number of vertices/faces/edges\n";
+ return false;
+ }
+ off_visitor.init(off_info_.dim, off_info_.num_vertices, off_info_.num_faces, off_info_.num_edges);
+
+ return true;
+ }
+
+ bool goto_next_uncomment_line(std::string& uncomment_line) {
+ do {
+ // skip whitespace, including empty lines
+ if (!std::ifstream::sentry(stream_)) return false;
+ std::getline(stream_, uncomment_line);
+ } while (uncomment_line[0] == '#');
+ return static_cast<bool>(stream_);
+ }
+
+ template<typename OffVisitor>
+ bool read_off_points(OffVisitor& visitor) {
+ int num_vertices_to_read = off_info_.num_vertices;
+ while (num_vertices_to_read--) {
+ std::string line;
+ if (!goto_next_uncomment_line(line)) return false;
+ std::vector<double> point;
+ std::istringstream iss(line);
+ point.assign(std::istream_iterator<double>(iss), std::istream_iterator<double>());
+ // if(point.size() != off_info_.dim) return false;
+ visitor.point(point);
+ }
+ return true;
+ }
+
+ template<typename OffVisitor>
+ bool read_off_faces(OffVisitor& visitor) {
+ std::string line;
+ while (goto_next_uncomment_line(line)) {
+ std::istringstream iss(line);
+ int num_face_vertices;
+ iss >> num_face_vertices;
+ std::vector<int> face;
+ face.assign(std::istream_iterator<int>(iss), std::istream_iterator<int>());
+ // if (face.size() != (off_info_.dim + 1)) return false;
+ visitor.maximal_face(face);
+ }
+ return true;
+ }
+};
+
+template<typename OFFVisitor>
+void read_off(const std::string& name_file_off, OFFVisitor& vis) {
+ std::ifstream stream(name_file_off);
+ if (!stream.is_open()) {
+ std::cerr << "could not open file \n";
+ } else {
+ Off_reader off_reader(stream);
+ off_reader.read(vis);
+ }
+}
+
+} // namespace Gudhi
+
+#endif // OFF_READER_H_
diff --git a/src/common/include/gudhi/Point.h b/src/common/include/gudhi/Point.h
new file mode 100644
index 00000000..e85277e9
--- /dev/null
+++ b/src/common/include/gudhi/Point.h
@@ -0,0 +1,157 @@
+/* This file is part of the Gudhi Library - https://gudhi.inria.fr/ - which is released under MIT.
+ * See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details.
+ * Author(s): David Salinas
+ *
+ * Copyright (C) 2014 Inria
+ *
+ * Modification(s):
+ * - YYYY/MM Author: Description of the modification
+ */
+
+#ifndef POINT_H_
+#define POINT_H_
+
+#include <cmath>
+#include <vector>
+#include <cassert>
+#include <cstddef>
+#include <initializer_list>
+
+class Point_d {
+ public:
+ Point_d(size_t dim = 3) : coords_(dim, 0) { }
+
+ Point_d(const Point_d& other) : coords_(other.coords_) { }
+
+ Point_d(const std::initializer_list<double>& list) : coords_(list) { }
+
+ template<typename CoordsIt>
+ Point_d(CoordsIt begin, CoordsIt end) : coords_(begin, end) { }
+
+ size_t dimension() const {
+ return coords_.size();
+ }
+
+ double x() const {
+ return coords_[0];
+ }
+
+ double y() const {
+ return coords_[1];
+ }
+
+ double z() const {
+ return coords_[2];
+ }
+
+ double& x() {
+ return coords_[0];
+ }
+
+ double& y() {
+ return coords_[1];
+ }
+
+ double& z() {
+ return coords_[2];
+ }
+
+ std::vector<double>::const_iterator begin() const {
+ return coords_.begin();
+ }
+
+ std::vector<double>::const_iterator end() const {
+ return coords_.end();
+ }
+
+ double& operator[](unsigned i) {
+ return coords_[i];
+ }
+
+ const double& operator[](unsigned i) const {
+ return coords_[i];
+ }
+
+ double squared_norm() const {
+ double res = 0;
+ for (auto x : coords_)
+ res += x * x;
+ return res;
+ }
+
+ friend double squared_dist(const Point_d& p1, const Point_d& p2) {
+ assert(p1.dimension() == p2.dimension());
+ double res = 0;
+ for (unsigned i = 0; i < p1.coords_.size(); ++i)
+ res += (p1[i] - p2[i])*(p1[i] - p2[i]);
+ return res;
+ }
+
+ /**
+ * dot product
+ */
+ double operator*(const Point_d& other) const {
+ assert(dimension() == other.dimension());
+ double res = 0;
+ for (unsigned i = 0; i < coords_.size(); ++i)
+ res += coords_[i] * other[i];
+ return res;
+ }
+
+ /**
+ * only if points have dimension 3
+ */
+ Point_d cross_product(const Point_d& other) {
+ assert(dimension() == 3 && other.dimension() == 3);
+ Point_d res(3);
+ res[0] = (*this)[1] * other[2] - (*this)[2] * other[1];
+ res[1] = (*this)[2] * other[0] - (*this)[0] * other[2];
+ res[2] = (*this)[0] * other[1] - (*this)[1] * other[0];
+ return res;
+ }
+
+ Point_d operator+(const Point_d& other) const {
+ assert(dimension() == other.dimension());
+ Point_d res(dimension());
+ for (unsigned i = 0; i < coords_.size(); ++i)
+ res[i] = (*this)[i] + other[i];
+ return res;
+ }
+
+ Point_d operator*(double lambda) const {
+ Point_d res(dimension());
+ for (unsigned i = 0; i < coords_.size(); ++i)
+ res[i] = (*this)[i] * lambda;
+ return res;
+ }
+
+ Point_d operator/(double lambda) const {
+ Point_d res(dimension());
+ for (unsigned i = 0; i < coords_.size(); ++i)
+ res[i] = (*this)[i] / lambda;
+ return res;
+ }
+
+ Point_d operator-(const Point_d& other) const {
+ assert(dimension() == other.dimension());
+ Point_d res(dimension());
+ for (unsigned i = 0; i < coords_.size(); ++i)
+ res[i] = (*this)[i] - other[i];
+ return res;
+ }
+
+ friend Point_d unit_normal(const Point_d& p1, const Point_d& p2, const Point_d& p3) {
+ assert(p1.dimension() == 3);
+ assert(p2.dimension() == 3);
+ assert(p3.dimension() == 3);
+ Point_d p1p2 = p2 - p1;
+ Point_d p1p3 = p3 - p1;
+ Point_d res(p1p2.cross_product(p1p3));
+ return res / std::sqrt(res.squared_norm());
+ }
+
+ private:
+ std::vector<double> coords_;
+};
+
+#endif // POINT_H_
diff --git a/src/common/include/gudhi/Points_3D_off_io.h b/src/common/include/gudhi/Points_3D_off_io.h
new file mode 100644
index 00000000..2d110af3
--- /dev/null
+++ b/src/common/include/gudhi/Points_3D_off_io.h
@@ -0,0 +1,190 @@
+/* This file is part of the Gudhi Library - https://gudhi.inria.fr/ - which is released under MIT.
+ * See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details.
+ * Author(s): Vincent Rouvreau
+ *
+ * Copyright (C) 2015 Inria
+ *
+ * Modification(s):
+ * - YYYY/MM Author: Description of the modification
+ */
+#ifndef POINTS_3D_OFF_IO_H_
+#define POINTS_3D_OFF_IO_H_
+
+#include <gudhi/Off_reader.h>
+
+#include <string>
+#include <vector>
+#include <fstream>
+#include <map>
+
+namespace Gudhi {
+
+/**
+ * @brief OFF file visitor implementation according to Off_reader in order to read points from an OFF file.
+ */
+template<typename Point_3>
+class Points_3D_off_visitor_reader {
+ private:
+ std::vector<Point_3> point_cloud_;
+ bool valid_;
+
+ public:
+ /** @brief Off_reader visitor init implementation.
+ *
+ * The init parameters are set from OFF file header.
+ * Dimension value is required and the value must be 3.
+ *
+ * @param[in] dim space dimension of vertices.
+ * @param[in] num_vertices number of vertices in the OFF file (not used).
+ * @param[in] num_faces number of faces in the OFF file (not used).
+ * @param[in] num_edges number of edges in the OFF file (not used).
+ */
+ void init(int dim, int num_vertices, int num_faces, int num_edges) {
+#ifdef DEBUG_TRACES
+ std::cout << "Points_3D_off_visitor_reader::init - dim=" << dim << " - num_vertices=" <<
+ num_vertices << " - num_faces=" << num_faces << " - num_edges=" << num_edges << std::endl;
+#endif // DEBUG_TRACES
+ if (dim == 3) {
+ valid_ = true;
+ } else {
+ valid_ = false;
+ std::cerr << "Points_3D_off_reader::Points_3D_off_reader cannot read OFF files in dimension " << dim << "\n";
+ }
+
+ if (num_faces > 0) {
+ std::cerr << "Points_3D_off_visitor_reader::init faces are not taken into account from OFF file for Points.\n";
+ }
+ if (num_edges > 0) {
+ std::cerr << "Points_3D_off_visitor_reader::init edges are not taken into account from OFF file for Points.\n";
+ }
+ }
+
+ /** @brief Off_reader visitor point implementation.
+ *
+ * The point function is called on each vertex line from OFF file.
+ * This function inserts the vertex in the vector of points.
+ *
+ * @param[in] point vector of vertex coordinates.
+ *
+ * @details
+ * Point_3 must have a constructor with the following form:
+ *
+ * @code template<class InputIterator > Point_3::Point_3(double x, double y, double z) @endcode
+ */
+ void point(const std::vector<double>& point) {
+ if (valid_) {
+#ifdef DEBUG_TRACES
+ std::cout << "Points_3D_off_visitor_reader::point ";
+ for (auto coordinate : point) {
+ std::cout << coordinate << " | ";
+ }
+ std::cout << std::endl;
+#endif // DEBUG_TRACES
+ // Fill the point cloud
+ point_cloud_.push_back(Point_3(point[0], point[1], point[2]));
+ }
+ }
+
+ // Off_reader visitor maximal_face implementation - Only points are read
+
+ void maximal_face(const std::vector<int>& face) { }
+
+ // Off_reader visitor done implementation - Only points are read
+
+ void done() { }
+
+ /** @brief Point cloud getter.
+ *
+ * @return The point cloud.
+ */
+ const std::vector<Point_3>& get_point_cloud() const {
+ return point_cloud_;
+ }
+
+ /** @brief Returns if the OFF file read operation was successful or not.
+ *
+ * @return OFF file read status.
+ */
+ bool is_valid() const {
+ return valid_;
+ }
+};
+
+/**
+ * \@brief OFF file reader implementation in order to read dimension 3 points from an OFF file.
+ *
+ * @details
+ * This class is using the Points_3D_off_visitor_reader to visit the OFF file according to Off_reader.
+ *
+ * Point_3 must have a constructor with the following form:
+ *
+ * @code template<class InputIterator > Point_3::Point_3(double x, double y, double z) @endcode
+ *
+ * @section point3doffioexample Example
+ *
+ * This example loads points from an OFF file and builds a vector of CGAL points in dimension 3.
+ * Then, it is asked to display the points.
+ *
+ * @include common/example_CGAL_3D_points_off_reader.cpp
+ *
+ * When launching:
+ *
+ * @code $> ./cgal3Doffreader ../../data/points/tore3D_300.off
+ * @endcode
+ *
+ * the program output is:
+ *
+ * @include common/cgal3Doffreader_result.txt
+ */
+template<typename Point_3>
+class Points_3D_off_reader {
+ public:
+ /** @brief Reads the OFF file and constructs a vector of points from the points
+ * that are in the OFF file.
+ *
+ * @param[in] name_file OFF file to read.
+ *
+ * @post Check with is_valid() function to see if read operation was successful.
+ */
+ Points_3D_off_reader(const std::string& name_file)
+ : valid_(false) {
+ std::ifstream stream(name_file);
+ if (stream.is_open()) {
+ Off_reader off_reader(stream);
+ Points_3D_off_visitor_reader<Point_3> off_visitor;
+ valid_ = off_reader.read(off_visitor);
+ valid_ = valid_ && off_visitor.is_valid();
+ if (valid_) {
+ point_cloud = off_visitor.get_point_cloud();
+ }
+ } else {
+ std::cerr << "Points_3D_off_reader::Points_3D_off_reader could not open file " << name_file << "\n";
+ }
+ }
+
+ /** @brief Returns if the OFF file read operation was successful or not.
+ *
+ * @return OFF file read status.
+ */
+ bool is_valid() const {
+ return valid_;
+ }
+
+ /** @brief Point cloud getter.
+ *
+ * @return point_cloud.
+ */
+ const std::vector<Point_3>& get_point_cloud() const {
+ return point_cloud;
+ }
+
+ private:
+ /** @brief point_cloud.*/
+ std::vector<Point_3> point_cloud;
+ /** @brief OFF file read status.*/
+ bool valid_;
+};
+
+} // namespace Gudhi
+
+#endif // POINTS_3D_OFF_IO_H_
diff --git a/src/common/include/gudhi/Points_off_io.h b/src/common/include/gudhi/Points_off_io.h
new file mode 100644
index 00000000..99371d56
--- /dev/null
+++ b/src/common/include/gudhi/Points_off_io.h
@@ -0,0 +1,171 @@
+/* This file is part of the Gudhi Library - https://gudhi.inria.fr/ - which is released under MIT.
+ * See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details.
+ * Author(s): Vincent Rouvreau
+ *
+ * Copyright (C) 2015 Inria
+ *
+ * Modification(s):
+ * - YYYY/MM Author: Description of the modification
+ */
+#ifndef POINTS_OFF_IO_H_
+#define POINTS_OFF_IO_H_
+
+#include <gudhi/Off_reader.h>
+
+#include <string>
+#include <vector>
+#include <fstream>
+#include <map>
+
+namespace Gudhi {
+
+/**
+ * \brief OFF file visitor implementation according to Off_reader in order to read points from an OFF file.
+ */
+template<typename Point_d>
+class Points_off_visitor_reader {
+ private:
+ std::vector<Point_d> point_cloud;
+
+ public:
+ /** \brief Off_reader visitor init implementation.
+ *
+ * The init parameters are set from OFF file header.
+ * Dimension value is required in order to construct a vector of points.
+ *
+ * @param[in] dim space dimension of vertices.
+ * @param[in] num_vertices number of vertices in the OFF file (not used).
+ * @param[in] num_faces number of faces in the OFF file (not used).
+ * @param[in] num_edges number of edges in the OFF file (not used).
+ */
+ void init(int dim, int num_vertices, int num_faces, int num_edges) {
+#ifdef DEBUG_TRACES
+ std::cout << "Points_off_visitor_reader::init - dim=" << dim << " - num_vertices=" <<
+ num_vertices << " - num_faces=" << num_faces << " - num_edges=" << num_edges << std::endl;
+#endif // DEBUG_TRACES
+ if (num_faces > 0) {
+ std::cerr << "Points_off_visitor_reader::init faces are not taken into account from OFF file for Points.\n";
+ }
+ if (num_edges > 0) {
+ std::cerr << "Points_off_visitor_reader::init edges are not taken into account from OFF file for Points.\n";
+ }
+ }
+
+ /** @brief Off_reader visitor point implementation.
+ *
+ * The point function is called on each vertex line from OFF file.
+ * This function inserts the vertex in the vector of points.
+ *
+ * @param[in] point vector of vertex coordinates.
+ *
+ * @details
+ * Point_d must have a constructor with the following form:
+ *
+ * @code template<class InputIterator > Point_d::Point_d(InputIterator first, InputIterator last) @endcode
+ *
+ */
+ void point(const std::vector<double>& point) {
+#ifdef DEBUG_TRACES
+ std::cout << "Points_off_visitor_reader::point ";
+ for (auto coordinate : point) {
+ std::cout << coordinate << " | ";
+ }
+ std::cout << std::endl;
+#endif // DEBUG_TRACES
+ // Fill the point cloud
+ point_cloud.push_back(Point_d(point.begin(), point.end()));
+ }
+
+ // Off_reader visitor maximal_face implementation - Only points are read
+ void maximal_face(const std::vector<int>& face) { }
+
+ // Off_reader visitor done implementation - Only points are read
+ void done() { }
+
+ /** \brief Point cloud getter.
+ *
+ * @return point_cloud.
+ */
+ const std::vector<Point_d>& get_point_cloud() const {
+ return point_cloud;
+ }
+};
+
+/**
+ * \brief OFF file reader implementation in order to read points from an OFF file.
+ *
+ * This class is using the Points_off_visitor_reader to visit the OFF file according to Off_reader.
+ *
+ * Point_d must have a constructor with the following form:
+ *
+ * \code template<class InputIterator > Point_d::Point_d(int d, InputIterator first, InputIterator last) \endcode
+ *
+ * where d is the point dimension.
+ *
+ * \section pointoffioexample Example
+ *
+ * This example loads points from an OFF file and builds a vector of points (vector of double).
+ * Then, it is asked to display the points.
+ *
+ * \include common/example_vector_double_points_off_reader.cpp
+ *
+ * When launching:
+ *
+ * \code $> ./vector_double_off_reader ../../data/points/alphacomplexdoc.off
+ * \endcode
+ *
+ * the program outputs a file ../../data/points/alphacomplexdoc.off.txt:
+ *
+ * \include common/vectordoubleoffreader_result.txt
+ */
+template<typename Point_d>
+class Points_off_reader {
+ public:
+ /** \brief Reads the OFF file and constructs a vector of points from the points
+ * that are in the OFF file.
+ *
+ * @param[in] name_file OFF file to read.
+ *
+ * \post Check with is_valid() function to see if read operation was successful.
+ */
+ Points_off_reader(const std::string& name_file)
+ : valid_(false) {
+ std::ifstream stream(name_file);
+ if (stream.is_open()) {
+ Off_reader off_reader(stream);
+ Points_off_visitor_reader<Point_d> off_visitor;
+ valid_ = off_reader.read(off_visitor);
+ if (valid_) {
+ point_cloud = off_visitor.get_point_cloud();
+ }
+ } else {
+ std::cerr << "Points_off_reader::Points_off_reader could not open file " << name_file << "\n";
+ }
+ }
+
+ /** \brief Returns if the OFF file read operation was successful or not.
+ *
+ * @return OFF file read status.
+ */
+ bool is_valid() const {
+ return valid_;
+ }
+
+ /** \brief Point cloud getter.
+ *
+ * @return point_cloud.
+ */
+ const std::vector<Point_d>& get_point_cloud() const {
+ return point_cloud;
+ }
+
+ private:
+ /** \brief point_cloud.*/
+ std::vector<Point_d> point_cloud;
+ /** \brief OFF file read status.*/
+ bool valid_;
+};
+
+} // namespace Gudhi
+
+#endif // POINTS_OFF_IO_H_
diff --git a/src/common/include/gudhi/Simple_object_pool.h b/src/common/include/gudhi/Simple_object_pool.h
new file mode 100644
index 00000000..d1482b44
--- /dev/null
+++ b/src/common/include/gudhi/Simple_object_pool.h
@@ -0,0 +1,69 @@
+/* This file is part of the Gudhi Library - https://gudhi.inria.fr/ - which is released under MIT.
+ * See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details.
+ * Author(s): Marc Glisse
+ *
+ * Copyright (C) 2015 Inria
+ *
+ * Modification(s):
+ * - YYYY/MM Author: Description of the modification
+ */
+
+#ifndef SIMPLE_OBJECT_POOL_H_
+#define SIMPLE_OBJECT_POOL_H_
+
+#include <boost/pool/pool.hpp>
+#include <utility>
+
+namespace Gudhi {
+
+/** \private
+ * This is a simpler version of boost::object_pool, that requires
+ * that users explicitly destroy all objects. This lets the
+ * performance scale much better, see
+ * https://svn.boost.org/trac/boost/ticket/3789 .
+ */
+template <class T>
+class Simple_object_pool : protected boost::pool<boost::default_user_allocator_malloc_free> {
+ protected:
+ typedef boost::pool<boost::default_user_allocator_malloc_free> Base;
+ typedef T* pointer;
+
+ Base& base() {
+ return *this;
+ }
+
+ Base const& base()const {
+ return *this;
+ }
+
+ public:
+ typedef T element_type;
+ typedef boost::default_user_allocator_malloc_free user_allocator;
+ typedef typename Base::size_type size_type;
+ typedef typename Base::difference_type difference_type;
+
+ template<class...U>
+ Simple_object_pool(U&&...u) : Base(sizeof (T), std::forward<U>(u)...) { }
+
+ template<class...U>
+ pointer construct(U&&...u) {
+ void* p = base().malloc BOOST_PREVENT_MACRO_SUBSTITUTION();
+ assert(p);
+ try {
+ new(p) T(std::forward<U>(u)...);
+ } catch (...) {
+ base().free BOOST_PREVENT_MACRO_SUBSTITUTION(p);
+ throw;
+ }
+ return static_cast<pointer> (p);
+ }
+
+ void destroy(pointer p) {
+ p->~T();
+ base().free BOOST_PREVENT_MACRO_SUBSTITUTION(p);
+ }
+};
+
+} // namespace Gudhi
+
+#endif // SIMPLE_OBJECT_POOL_H_
diff --git a/src/common/include/gudhi/Unitary_tests_utils.h b/src/common/include/gudhi/Unitary_tests_utils.h
new file mode 100644
index 00000000..4ad4dae8
--- /dev/null
+++ b/src/common/include/gudhi/Unitary_tests_utils.h
@@ -0,0 +1,28 @@
+/* This file is part of the Gudhi Library - https://gudhi.inria.fr/ - which is released under MIT.
+ * See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details.
+ * Author(s): Vincent Rouvreau
+ *
+ * Copyright (C) 2017 Inria
+ *
+ * Modification(s):
+ * - YYYY/MM Author: Description of the modification
+ */
+#ifndef UNITARY_TESTS_UTILS_H_
+#define UNITARY_TESTS_UTILS_H_
+
+#include <boost/test/unit_test.hpp>
+
+#include <iostream>
+#include <limits> // for std::numeric_limits<>
+
+template<typename FloatingType >
+void GUDHI_TEST_FLOAT_EQUALITY_CHECK(FloatingType a, FloatingType b,
+ FloatingType epsilon = std::numeric_limits<FloatingType>::epsilon()) {
+#ifdef DEBUG_TRACES
+ std::cout << "GUDHI_TEST_FLOAT_EQUALITY_CHECK - " << a << " versus " << b
+ << " | diff = " << std::fabs(a - b) << " - epsilon = " << epsilon << std::endl;
+#endif
+ BOOST_CHECK(std::fabs(a - b) <= epsilon);
+}
+
+#endif // UNITARY_TESTS_UTILS_H_
diff --git a/src/common/include/gudhi/allocator.h b/src/common/include/gudhi/allocator.h
new file mode 100644
index 00000000..b7ccd180
--- /dev/null
+++ b/src/common/include/gudhi/allocator.h
@@ -0,0 +1,43 @@
+/* This file is part of the Gudhi Library - https://gudhi.inria.fr/ - which is released under MIT.
+ * See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details.
+ * Author(s): Marc Glisse
+ *
+ * Copyright (C) 2015 Inria
+ *
+ * Modification(s):
+ * - YYYY/MM Author: Description of the modification
+ */
+
+#ifndef ALLOCATOR_H_
+#define ALLOCATOR_H_
+
+#include <memory>
+#include <utility>
+
+namespace Gudhi {
+
+/** \private
+ * An allocator that can be used to build an uninitialized vector.
+ */
+template <class T, class Base = std::allocator<T>>
+struct no_init_allocator : Base {
+ typedef std::allocator_traits<Base> Base_traits;
+ template <class U> struct rebind {
+ typedef no_init_allocator<U, typename Base_traits::template rebind_alloc<U>> other;
+ };
+
+ // Inherit constructors.
+ using Base::Base;
+
+ // Do nothing: that's the whole point!
+ template<class P>
+ void construct(P*) noexcept {}
+
+ template<class P, class...U> void construct(P*p, U&&...u) {
+ Base_traits::construct(*(Base*)this, p, std::forward<U>(u)...);
+ }
+};
+
+} // namespace Gudhi
+
+#endif // ALLOCATOR_H_
diff --git a/src/common/include/gudhi/console_color.h b/src/common/include/gudhi/console_color.h
new file mode 100644
index 00000000..f9167119
--- /dev/null
+++ b/src/common/include/gudhi/console_color.h
@@ -0,0 +1,85 @@
+/* This file is part of the Gudhi Library - https://gudhi.inria.fr/ - which is released under MIT.
+ * See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details.
+ * Author(s): Clement Jamin
+ *
+ * Copyright (C) 2016 Inria
+ *
+ * Modification(s):
+ * - YYYY/MM Author: Description of the modification
+ */
+
+#ifndef CONSOLE_COLOR_H_
+#define CONSOLE_COLOR_H_
+
+#include <iostream>
+
+#if defined(WIN32)
+#include <windows.h>
+#endif
+
+inline std::ostream& blue(std::ostream &s) {
+#if defined(WIN32)
+ HANDLE hStdout = GetStdHandle(STD_OUTPUT_HANDLE);
+ SetConsoleTextAttribute(hStdout,
+ FOREGROUND_BLUE | FOREGROUND_GREEN | FOREGROUND_INTENSITY);
+#else
+ s << "\x1b[0;34m";
+#endif
+ return s;
+}
+
+inline std::ostream& red(std::ostream &s) {
+#if defined(WIN32)
+ HANDLE hStdout = GetStdHandle(STD_OUTPUT_HANDLE);
+ SetConsoleTextAttribute(hStdout, FOREGROUND_RED | FOREGROUND_INTENSITY);
+#else
+ s << "\x1b[0;31m";
+#endif
+ return s;
+}
+
+inline std::ostream& green(std::ostream &s) {
+#if defined(WIN32)
+ HANDLE hStdout = GetStdHandle(STD_OUTPUT_HANDLE);
+ SetConsoleTextAttribute(hStdout, FOREGROUND_GREEN | FOREGROUND_INTENSITY);
+#else
+ s << "\x1b[0;32m";
+#endif
+ return s;
+}
+
+inline std::ostream& yellow(std::ostream &s) {
+#if defined(WIN32)
+ HANDLE hStdout = GetStdHandle(STD_OUTPUT_HANDLE);
+ SetConsoleTextAttribute(hStdout,
+ FOREGROUND_GREEN | FOREGROUND_RED | FOREGROUND_INTENSITY);
+#else
+ s << "\x1b[0;33m";
+#endif
+ return s;
+}
+
+inline std::ostream& white(std::ostream &s) {
+#if defined(WIN32)
+ HANDLE hStdout = GetStdHandle(STD_OUTPUT_HANDLE);
+ SetConsoleTextAttribute(hStdout,
+ FOREGROUND_RED | FOREGROUND_GREEN | FOREGROUND_BLUE);
+#else
+ s << "\x1b[0;37m";
+#endif
+ return s;
+}
+
+inline std::ostream& black_on_white(std::ostream &s) {
+#if defined(WIN32)
+ HANDLE hStdout = GetStdHandle(STD_OUTPUT_HANDLE);
+ SetConsoleTextAttribute(hStdout,
+ BACKGROUND_RED | BACKGROUND_GREEN | BACKGROUND_BLUE);
+#else
+ s << "\x1b[0;33m";
+#endif
+ return s;
+}
+
+
+#endif // CONSOLE_COLOR_H_
diff --git a/src/common/include/gudhi/distance_functions.h b/src/common/include/gudhi/distance_functions.h
new file mode 100644
index 00000000..94cf9ccc
--- /dev/null
+++ b/src/common/include/gudhi/distance_functions.h
@@ -0,0 +1,111 @@
+/* This file is part of the Gudhi Library - https://gudhi.inria.fr/ - which is released under MIT.
+ * See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details.
+ * Author(s): Clément Maria
+ *
+ * Copyright (C) 2014 Inria
+ *
+ * Modification(s):
+ * - YYYY/MM Author: Description of the modification
+ */
+
+#ifndef DISTANCE_FUNCTIONS_H_
+#define DISTANCE_FUNCTIONS_H_
+
+#include <gudhi/Debug_utils.h>
+
+#include <gudhi/Miniball.hpp>
+
+#include <boost/range/metafunctions.hpp>
+#include <boost/range/size.hpp>
+
+#include <cmath> // for std::sqrt
+#include <type_traits> // for std::decay
+#include <iterator> // for std::begin, std::end
+#include <utility>
+
+namespace Gudhi {
+
+/** @file
+ * @brief Global distance functions
+ */
+
+/** @brief Compute the Euclidean distance between two Points given by a range of coordinates. The points are assumed to
+ * have the same dimension. */
+class Euclidean_distance {
+ public:
+ // boost::range_value is not SFINAE-friendly so we cannot use it in the return type
+ template< typename Point >
+ typename std::iterator_traits<typename boost::range_iterator<Point>::type>::value_type
+ operator()(const Point& p1, const Point& p2) const {
+ auto it1 = std::begin(p1);
+ auto it2 = std::begin(p2);
+ typedef typename boost::range_value<Point>::type NT;
+ NT dist = 0;
+ for (; it1 != std::end(p1); ++it1, ++it2) {
+ GUDHI_CHECK(it2 != std::end(p2), "inconsistent point dimensions");
+ NT tmp = *it1 - *it2;
+ dist += tmp*tmp;
+ }
+ GUDHI_CHECK(it2 == std::end(p2), "inconsistent point dimensions");
+ using std::sqrt;
+ return sqrt(dist);
+ }
+ template< typename T >
+ T operator() (const std::pair< T, T >& f, const std::pair< T, T >& s) const {
+ T dx = f.first - s.first;
+ T dy = f.second - s.second;
+ using std::sqrt;
+ return sqrt(dx*dx + dy*dy);
+ }
+};
+
+/** @brief Compute the radius of the minimal enclosing ball between Points given by a range of coordinates.
+ * The points are assumed to have the same dimension. */
+class Minimal_enclosing_ball_radius {
+ public:
+ /** \brief Minimal_enclosing_ball_radius from two points.
+ *
+ * @param[in] point_1 First point.
+ * @param[in] point_2 second point.
+ * @return The minimal enclosing ball radius for the two points (aka. Euclidean distance / 2.).
+ *
+ * \tparam Point must be a range of Cartesian coordinates.
+ *
+ */
+ template< typename Point >
+ typename std::iterator_traits<typename boost::range_iterator<Point>::type>::value_type
+ operator()(const Point& point_1, const Point& point_2) const {
+ return Euclidean_distance()(point_1, point_2) / 2.;
+ }
+ /** \brief Minimal_enclosing_ball_radius from a point cloud.
+ *
+ * @param[in] point_cloud The points.
+ * @return The minimal enclosing ball radius for the points.
+ *
+ * \tparam Point_cloud must be a range of points with Cartesian coordinates.
+ * Point_cloud is a range over a range of Coordinate.
+ *
+ */
+ template< typename Point_cloud,
+ typename Point_iterator = typename boost::range_const_iterator<Point_cloud>::type,
+ typename Point = typename std::iterator_traits<Point_iterator>::value_type,
+ typename Coordinate_iterator = typename boost::range_const_iterator<Point>::type,
+ typename Coordinate = typename std::iterator_traits<Coordinate_iterator>::value_type>
+ Coordinate
+ operator()(const Point_cloud& point_cloud) const {
+ using Min_sphere = Miniball::Miniball<Miniball::CoordAccessor<Point_iterator, Coordinate_iterator>>;
+
+ Min_sphere ms(boost::size(*point_cloud.begin()), point_cloud.begin(), point_cloud.end());
+#ifdef DEBUG_TRACES
+ std::cout << "Minimal_enclosing_ball_radius = " << std::sqrt(ms.squared_radius()) << " | nb points = "
+ << boost::size(point_cloud) << " | dimension = "
+ << boost::size(*point_cloud.begin()) << std::endl;
+#endif // DEBUG_TRACES
+
+ return std::sqrt(ms.squared_radius());
+ }
+};
+
+} // namespace Gudhi
+
+#endif // DISTANCE_FUNCTIONS_H_
diff --git a/src/common/include/gudhi/graph_simplicial_complex.h b/src/common/include/gudhi/graph_simplicial_complex.h
new file mode 100644
index 00000000..b8508697
--- /dev/null
+++ b/src/common/include/gudhi/graph_simplicial_complex.h
@@ -0,0 +1,99 @@
+/* This file is part of the Gudhi Library - https://gudhi.inria.fr/ - which is released under MIT.
+ * See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details.
+ * Author(s): Clément Maria
+ *
+ * Copyright (C) 2014 Inria
+ *
+ * Modification(s):
+ * - YYYY/MM Author: Description of the modification
+ */
+
+#ifndef GRAPH_SIMPLICIAL_COMPLEX_H_
+#define GRAPH_SIMPLICIAL_COMPLEX_H_
+
+#include <boost/graph/adjacency_list.hpp>
+
+#include <utility> // for pair<>
+#include <vector>
+#include <map>
+#include <tuple> // for std::tie
+
+namespace Gudhi {
+
+/* Edge tag for Boost PropertyGraph. */
+struct edge_filtration_t {
+ typedef boost::edge_property_tag kind;
+};
+
+/* Vertex tag for Boost PropertyGraph. */
+struct vertex_filtration_t {
+ typedef boost::vertex_property_tag kind;
+};
+
+/** \brief Proximity_graph contains the vertices and edges with their filtration values in order to store the result
+ * of `Gudhi::compute_proximity_graph` function.
+ *
+ * \tparam SimplicialComplexForProximityGraph furnishes `Filtration_value` type definition.
+ *
+ */
+template <typename SimplicialComplexForProximityGraph>
+using Proximity_graph = typename boost::adjacency_list < boost::vecS, boost::vecS, boost::directedS
+, boost::property < vertex_filtration_t, typename SimplicialComplexForProximityGraph::Filtration_value >
+, boost::property < edge_filtration_t, typename SimplicialComplexForProximityGraph::Filtration_value >>;
+
+/** \brief Computes the proximity graph of the points.
+ *
+ * If points contains n elements, the proximity graph is the graph with n vertices, and an edge [u,v] iff the
+ * distance function between points u and v is smaller than threshold.
+ *
+ * \tparam ForwardPointRange furnishes `.begin()` and `.end()` methods.
+ *
+ * \tparam Distance furnishes `operator()(const Point& p1, const Point& p2)`, where
+ * `Point` is a point from the `ForwardPointRange`, and that returns a `Filtration_value`.
+ */
+template< typename SimplicialComplexForProximityGraph
+ , typename ForwardPointRange
+ , typename Distance >
+Proximity_graph<SimplicialComplexForProximityGraph> compute_proximity_graph(
+ const ForwardPointRange& points,
+ typename SimplicialComplexForProximityGraph::Filtration_value threshold,
+ Distance distance) {
+ using Vertex_handle = typename SimplicialComplexForProximityGraph::Vertex_handle;
+ using Filtration_value = typename SimplicialComplexForProximityGraph::Filtration_value;
+
+ std::vector<std::pair< Vertex_handle, Vertex_handle >> edges;
+ std::vector< Filtration_value > edges_fil;
+ std::map< Vertex_handle, Filtration_value > vertices;
+
+ Vertex_handle idx_u, idx_v;
+ Filtration_value fil;
+ idx_u = 0;
+ for (auto it_u = points.begin(); it_u != points.end(); ++it_u) {
+ idx_v = idx_u + 1;
+ for (auto it_v = it_u + 1; it_v != points.end(); ++it_v, ++idx_v) {
+ fil = distance(*it_u, *it_v);
+ if (fil <= threshold) {
+ edges.emplace_back(idx_u, idx_v);
+ edges_fil.push_back(fil);
+ }
+ }
+ ++idx_u;
+ }
+
+ // Points are labeled from 0 to idx_u-1
+ Proximity_graph<SimplicialComplexForProximityGraph> skel_graph(edges.begin(), edges.end(), edges_fil.begin(), idx_u);
+
+ auto vertex_prop = boost::get(vertex_filtration_t(), skel_graph);
+
+ typename boost::graph_traits<Proximity_graph<SimplicialComplexForProximityGraph>>::vertex_iterator vi, vi_end;
+ for (std::tie(vi, vi_end) = boost::vertices(skel_graph);
+ vi != vi_end; ++vi) {
+ boost::put(vertex_prop, *vi, 0.);
+ }
+
+ return skel_graph;
+}
+
+} // namespace Gudhi
+
+#endif // GRAPH_SIMPLICIAL_COMPLEX_H_
diff --git a/src/common/include/gudhi/random_point_generators.h b/src/common/include/gudhi/random_point_generators.h
new file mode 100644
index 00000000..fb69f832
--- /dev/null
+++ b/src/common/include/gudhi/random_point_generators.h
@@ -0,0 +1,502 @@
+/* This file is part of the Gudhi Library - https://gudhi.inria.fr/ - which is released under MIT.
+ * See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details.
+ * Author(s): Clement Jamin
+ *
+ * Copyright (C) 2016 Inria
+ *
+ * Modification(s):
+ * - 2019/08 Vincent Rouvreau: Fix issue #10 for CGAL
+ * - YYYY/MM Author: Description of the modification
+ */
+
+#ifndef RANDOM_POINT_GENERATORS_H_
+#define RANDOM_POINT_GENERATORS_H_
+
+#include <CGAL/number_utils.h>
+#include <CGAL/Random.h>
+#include <CGAL/point_generators_d.h>
+#include <CGAL/version.h> // for CGAL_VERSION_NR
+
+#include <vector> // for vector<>
+
+// Make compilation fail - required for external projects - https://github.com/GUDHI/gudhi-devel/issues/10
+#if CGAL_VERSION_NR < 1041101000
+# error Alpha_complex_3d is only available for CGAL >= 4.11
+#endif
+
+namespace Gudhi {
+
+///////////////////////////////////////////////////////////////////////////////
+// Note: All these functions have been tested with the CGAL::Epick_d kernel
+///////////////////////////////////////////////////////////////////////////////
+
+// construct_point: dim 2
+
+template <typename Kernel>
+typename Kernel::Point_d construct_point(const Kernel &k,
+ typename Kernel::FT x1, typename Kernel::FT x2) {
+ typename Kernel::FT tab[2];
+ tab[0] = x1;
+ tab[1] = x2;
+ return k.construct_point_d_object()(2, &tab[0], &tab[2]);
+}
+
+// construct_point: dim 3
+
+template <typename Kernel>
+typename Kernel::Point_d construct_point(const Kernel &k,
+ typename Kernel::FT x1, typename Kernel::FT x2, typename Kernel::FT x3) {
+ typename Kernel::FT tab[3];
+ tab[0] = x1;
+ tab[1] = x2;
+ tab[2] = x3;
+ return k.construct_point_d_object()(3, &tab[0], &tab[3]);
+}
+
+// construct_point: dim 4
+
+template <typename Kernel>
+typename Kernel::Point_d construct_point(const Kernel &k,
+ typename Kernel::FT x1, typename Kernel::FT x2, typename Kernel::FT x3,
+ typename Kernel::FT x4) {
+ typename Kernel::FT tab[4];
+ tab[0] = x1;
+ tab[1] = x2;
+ tab[2] = x3;
+ tab[3] = x4;
+ return k.construct_point_d_object()(4, &tab[0], &tab[4]);
+}
+
+// construct_point: dim 5
+
+template <typename Kernel>
+typename Kernel::Point_d construct_point(const Kernel &k,
+ typename Kernel::FT x1, typename Kernel::FT x2, typename Kernel::FT x3,
+ typename Kernel::FT x4, typename Kernel::FT x5) {
+ typename Kernel::FT tab[5];
+ tab[0] = x1;
+ tab[1] = x2;
+ tab[2] = x3;
+ tab[3] = x4;
+ tab[4] = x5;
+ return k.construct_point_d_object()(5, &tab[0], &tab[5]);
+}
+
+// construct_point: dim 6
+
+template <typename Kernel>
+typename Kernel::Point_d construct_point(const Kernel &k,
+ typename Kernel::FT x1, typename Kernel::FT x2, typename Kernel::FT x3,
+ typename Kernel::FT x4, typename Kernel::FT x5, typename Kernel::FT x6) {
+ typename Kernel::FT tab[6];
+ tab[0] = x1;
+ tab[1] = x2;
+ tab[2] = x3;
+ tab[3] = x4;
+ tab[4] = x5;
+ tab[5] = x6;
+ return k.construct_point_d_object()(6, &tab[0], &tab[6]);
+}
+
+template <typename Kernel>
+std::vector<typename Kernel::Point_d> generate_points_on_plane(std::size_t num_points, int intrinsic_dim,
+ int ambient_dim,
+ double coord_min = -5., double coord_max = 5.) {
+ typedef typename Kernel::Point_d Point;
+ typedef typename Kernel::FT FT;
+ Kernel k;
+ CGAL::Random rng;
+ std::vector<Point> points;
+ points.reserve(num_points);
+ for (std::size_t i = 0; i < num_points;) {
+ std::vector<FT> pt(ambient_dim, FT(0));
+ for (int j = 0; j < intrinsic_dim; ++j)
+ pt[j] = rng.get_double(coord_min, coord_max);
+
+ Point p = k.construct_point_d_object()(ambient_dim, pt.begin(), pt.end());
+ points.push_back(p);
+ ++i;
+ }
+ return points;
+}
+
+template <typename Kernel>
+std::vector<typename Kernel::Point_d> generate_points_on_moment_curve(std::size_t num_points, int dim,
+ typename Kernel::FT min_x,
+ typename Kernel::FT max_x) {
+ typedef typename Kernel::Point_d Point;
+ typedef typename Kernel::FT FT;
+ Kernel k;
+ CGAL::Random rng;
+ std::vector<Point> points;
+ points.reserve(num_points);
+ for (std::size_t i = 0; i < num_points;) {
+ FT x = rng.get_double(min_x, max_x);
+ std::vector<FT> coords;
+ coords.reserve(dim);
+ for (int p = 1; p <= dim; ++p)
+ coords.push_back(std::pow(CGAL::to_double(x), p));
+ Point p = k.construct_point_d_object()(
+ dim, coords.begin(), coords.end());
+ points.push_back(p);
+ ++i;
+ }
+ return points;
+}
+
+
+// R = big radius, r = small radius
+template <typename Kernel/*, typename TC_basis*/>
+std::vector<typename Kernel::Point_d> generate_points_on_torus_3D(std::size_t num_points, double R, double r,
+ bool uniform = false) {
+ typedef typename Kernel::Point_d Point;
+ typedef typename Kernel::FT FT;
+ Kernel k;
+ CGAL::Random rng;
+
+ // if uniform
+ std::size_t num_lines = (std::size_t)sqrt(num_points);
+
+ std::vector<Point> points;
+ points.reserve(num_points);
+ for (std::size_t i = 0; i < num_points;) {
+ FT u, v;
+ if (uniform) {
+ std::size_t k1 = i / num_lines;
+ std::size_t k2 = i % num_lines;
+ u = 6.2832 * k1 / num_lines;
+ v = 6.2832 * k2 / num_lines;
+ } else {
+ u = rng.get_double(0, 6.2832);
+ v = rng.get_double(0, 6.2832);
+ }
+ Point p = construct_point(k,
+ (R + r * std::cos(u)) * std::cos(v),
+ (R + r * std::cos(u)) * std::sin(v),
+ r * std::sin(u));
+ points.push_back(p);
+ ++i;
+ }
+ return points;
+}
+
+// "Private" function used by generate_points_on_torus_d
+template <typename Kernel, typename OutputIterator>
+static void generate_uniform_points_on_torus_d(const Kernel &k, int dim, std::size_t num_slices,
+ OutputIterator out,
+ double radius_noise_percentage = 0.,
+ std::vector<typename Kernel::FT> current_point =
+ std::vector<typename Kernel::FT>()) {
+ CGAL::Random rng;
+ int point_size = static_cast<int>(current_point.size());
+ if (point_size == 2 * dim) {
+ *out++ = k.construct_point_d_object()(point_size, current_point.begin(), current_point.end());
+ } else {
+ for (std::size_t slice_idx = 0; slice_idx < num_slices; ++slice_idx) {
+ double radius_noise_ratio = 1.;
+ if (radius_noise_percentage > 0.) {
+ radius_noise_ratio = rng.get_double(
+ (100. - radius_noise_percentage) / 100.,
+ (100. + radius_noise_percentage) / 100.);
+ }
+ std::vector<typename Kernel::FT> cp2 = current_point;
+ double alpha = 6.2832 * slice_idx / num_slices;
+ cp2.push_back(radius_noise_ratio * std::cos(alpha));
+ cp2.push_back(radius_noise_ratio * std::sin(alpha));
+ generate_uniform_points_on_torus_d(
+ k, dim, num_slices, out, radius_noise_percentage, cp2);
+ }
+ }
+}
+
+template <typename Kernel>
+std::vector<typename Kernel::Point_d> generate_points_on_torus_d(std::size_t num_points, int dim, bool uniform = false,
+ double radius_noise_percentage = 0.) {
+ typedef typename Kernel::Point_d Point;
+ typedef typename Kernel::FT FT;
+ Kernel k;
+ CGAL::Random rng;
+
+ std::vector<Point> points;
+ points.reserve(num_points);
+ if (uniform) {
+ std::size_t num_slices = (std::size_t)std::pow(num_points, 1. / dim);
+ generate_uniform_points_on_torus_d(
+ k, dim, num_slices, std::back_inserter(points), radius_noise_percentage);
+ } else {
+ for (std::size_t i = 0; i < num_points;) {
+ double radius_noise_ratio = 1.;
+ if (radius_noise_percentage > 0.) {
+ radius_noise_ratio = rng.get_double(
+ (100. - radius_noise_percentage) / 100.,
+ (100. + radius_noise_percentage) / 100.);
+ }
+ std::vector<typename Kernel::FT> pt;
+ pt.reserve(dim * 2);
+ for (int curdim = 0; curdim < dim; ++curdim) {
+ FT alpha = rng.get_double(0, 6.2832);
+ pt.push_back(radius_noise_ratio * std::cos(alpha));
+ pt.push_back(radius_noise_ratio * std::sin(alpha));
+ }
+
+ Point p = k.construct_point_d_object()(pt.begin(), pt.end());
+ points.push_back(p);
+ ++i;
+ }
+ }
+ return points;
+}
+
+template <typename Kernel>
+std::vector<typename Kernel::Point_d> generate_points_on_sphere_d(std::size_t num_points, int dim, double radius,
+ double radius_noise_percentage = 0.) {
+ typedef typename Kernel::Point_d Point;
+ Kernel k;
+ CGAL::Random rng;
+ CGAL::Random_points_on_sphere_d<Point> generator(dim, radius);
+ std::vector<Point> points;
+ points.reserve(num_points);
+ for (std::size_t i = 0; i < num_points;) {
+ Point p = *generator++;
+ if (radius_noise_percentage > 0.) {
+ double radius_noise_ratio = rng.get_double(
+ (100. - radius_noise_percentage) / 100.,
+ (100. + radius_noise_percentage) / 100.);
+
+ typename Kernel::Point_to_vector_d k_pt_to_vec =
+ k.point_to_vector_d_object();
+ typename Kernel::Vector_to_point_d k_vec_to_pt =
+ k.vector_to_point_d_object();
+ typename Kernel::Scaled_vector_d k_scaled_vec =
+ k.scaled_vector_d_object();
+ p = k_vec_to_pt(k_scaled_vec(k_pt_to_vec(p), radius_noise_ratio));
+ }
+ points.push_back(p);
+ ++i;
+ }
+ return points;
+}
+
+template <typename Kernel>
+std::vector<typename Kernel::Point_d> generate_points_in_ball_d(std::size_t num_points, int dim, double radius) {
+ typedef typename Kernel::Point_d Point;
+ Kernel k;
+ CGAL::Random rng;
+ CGAL::Random_points_in_ball_d<Point> generator(dim, radius);
+ std::vector<Point> points;
+ points.reserve(num_points);
+ for (std::size_t i = 0; i < num_points;) {
+ Point p = *generator++;
+ points.push_back(p);
+ ++i;
+ }
+ return points;
+}
+
+template <typename Kernel>
+std::vector<typename Kernel::Point_d> generate_points_in_cube_d(std::size_t num_points, int dim, double radius) {
+ typedef typename Kernel::Point_d Point;
+ Kernel k;
+ CGAL::Random rng;
+ CGAL::Random_points_in_cube_d<Point> generator(dim, radius);
+ std::vector<Point> points;
+ points.reserve(num_points);
+ for (std::size_t i = 0; i < num_points;) {
+ Point p = *generator++;
+ points.push_back(p);
+ ++i;
+ }
+ return points;
+}
+
+template <typename Kernel>
+std::vector<typename Kernel::Point_d> generate_points_on_two_spheres_d(std::size_t num_points, int dim, double radius,
+ double distance_between_centers,
+ double radius_noise_percentage = 0.) {
+ typedef typename Kernel::FT FT;
+ typedef typename Kernel::Point_d Point;
+ typedef typename Kernel::Vector_d Vector;
+ Kernel k;
+ CGAL::Random rng;
+ CGAL::Random_points_on_sphere_d<Point> generator(dim, radius);
+ std::vector<Point> points;
+ points.reserve(num_points);
+
+ std::vector<FT> t(dim, FT(0));
+ t[0] = distance_between_centers;
+ Vector c1_to_c2(t.begin(), t.end());
+
+ for (std::size_t i = 0; i < num_points;) {
+ Point p = *generator++;
+ if (radius_noise_percentage > 0.) {
+ double radius_noise_ratio = rng.get_double(
+ (100. - radius_noise_percentage) / 100.,
+ (100. + radius_noise_percentage) / 100.);
+
+ typename Kernel::Point_to_vector_d k_pt_to_vec =
+ k.point_to_vector_d_object();
+ typename Kernel::Vector_to_point_d k_vec_to_pt =
+ k.vector_to_point_d_object();
+ typename Kernel::Scaled_vector_d k_scaled_vec =
+ k.scaled_vector_d_object();
+ p = k_vec_to_pt(k_scaled_vec(k_pt_to_vec(p), radius_noise_ratio));
+ }
+
+ typename Kernel::Translated_point_d k_transl =
+ k.translated_point_d_object();
+ Point p2 = k_transl(p, c1_to_c2);
+ points.push_back(p);
+ points.push_back(p2);
+ i += 2;
+ }
+ return points;
+}
+
+// Product of a 3-sphere and a circle => d = 3 / D = 5
+
+template <typename Kernel>
+std::vector<typename Kernel::Point_d> generate_points_on_3sphere_and_circle(std::size_t num_points,
+ double sphere_radius) {
+ typedef typename Kernel::FT FT;
+ typedef typename Kernel::Point_d Point;
+ Kernel k;
+ CGAL::Random rng;
+ CGAL::Random_points_on_sphere_d<Point> generator(3, sphere_radius);
+ std::vector<Point> points;
+ points.reserve(num_points);
+
+ typename Kernel::Compute_coordinate_d k_coord =
+ k.compute_coordinate_d_object();
+ for (std::size_t i = 0; i < num_points;) {
+ Point p_sphere = *generator++; // First 3 coords
+
+ FT alpha = rng.get_double(0, 6.2832);
+ std::vector<FT> pt(5);
+ pt[0] = k_coord(p_sphere, 0);
+ pt[1] = k_coord(p_sphere, 1);
+ pt[2] = k_coord(p_sphere, 2);
+ pt[3] = std::cos(alpha);
+ pt[4] = std::sin(alpha);
+ Point p(pt.begin(), pt.end());
+ points.push_back(p);
+ ++i;
+ }
+ return points;
+}
+
+// a = big radius, b = small radius
+template <typename Kernel>
+std::vector<typename Kernel::Point_d> generate_points_on_klein_bottle_3D(std::size_t num_points, double a, double b,
+ bool uniform = false) {
+ typedef typename Kernel::Point_d Point;
+ typedef typename Kernel::FT FT;
+ Kernel k;
+ CGAL::Random rng;
+
+ // if uniform
+ std::size_t num_lines = (std::size_t)sqrt(num_points);
+
+ std::vector<Point> points;
+ points.reserve(num_points);
+ for (std::size_t i = 0; i < num_points;) {
+ FT u, v;
+ if (uniform) {
+ std::size_t k1 = i / num_lines;
+ std::size_t k2 = i % num_lines;
+ u = 6.2832 * k1 / num_lines;
+ v = 6.2832 * k2 / num_lines;
+ } else {
+ u = rng.get_double(0, 6.2832);
+ v = rng.get_double(0, 6.2832);
+ }
+ double tmp = cos(u / 2) * sin(v) - sin(u / 2) * sin(2. * v);
+ Point p = construct_point(k,
+ (a + b * tmp) * cos(u),
+ (a + b * tmp) * sin(u),
+ b * (sin(u / 2) * sin(v) + cos(u / 2) * sin(2. * v)));
+ points.push_back(p);
+ ++i;
+ }
+ return points;
+}
+
+// a = big radius, b = small radius
+template <typename Kernel>
+std::vector<typename Kernel::Point_d> generate_points_on_klein_bottle_4D(std::size_t num_points, double a, double b,
+ double noise = 0., bool uniform = false) {
+ typedef typename Kernel::Point_d Point;
+ typedef typename Kernel::FT FT;
+ Kernel k;
+ CGAL::Random rng;
+
+ // if uniform
+ std::size_t num_lines = (std::size_t)sqrt(num_points);
+
+ std::vector<Point> points;
+ points.reserve(num_points);
+ for (std::size_t i = 0; i < num_points;) {
+ FT u, v;
+ if (uniform) {
+ std::size_t k1 = i / num_lines;
+ std::size_t k2 = i % num_lines;
+ u = 6.2832 * k1 / num_lines;
+ v = 6.2832 * k2 / num_lines;
+ } else {
+ u = rng.get_double(0, 6.2832);
+ v = rng.get_double(0, 6.2832);
+ }
+ Point p = construct_point(k,
+ (a + b * cos(v)) * cos(u) + (noise == 0. ? 0. : rng.get_double(0, noise)),
+ (a + b * cos(v)) * sin(u) + (noise == 0. ? 0. : rng.get_double(0, noise)),
+ b * sin(v) * cos(u / 2) + (noise == 0. ? 0. : rng.get_double(0, noise)),
+ b * sin(v) * sin(u / 2) + (noise == 0. ? 0. : rng.get_double(0, noise)));
+ points.push_back(p);
+ ++i;
+ }
+ return points;
+}
+
+
+// a = big radius, b = small radius
+
+template <typename Kernel>
+std::vector<typename Kernel::Point_d>
+generate_points_on_klein_bottle_variant_5D(
+ std::size_t num_points, double a, double b, bool uniform = false) {
+ typedef typename Kernel::Point_d Point;
+ typedef typename Kernel::FT FT;
+ Kernel k;
+ CGAL::Random rng;
+
+ // if uniform
+ std::size_t num_lines = (std::size_t)sqrt(num_points);
+
+ std::vector<Point> points;
+ points.reserve(num_points);
+ for (std::size_t i = 0; i < num_points;) {
+ FT u, v;
+ if (uniform) {
+ std::size_t k1 = i / num_lines;
+ std::size_t k2 = i % num_lines;
+ u = 6.2832 * k1 / num_lines;
+ v = 6.2832 * k2 / num_lines;
+ } else {
+ u = rng.get_double(0, 6.2832);
+ v = rng.get_double(0, 6.2832);
+ }
+ FT x1 = (a + b * cos(v)) * cos(u);
+ FT x2 = (a + b * cos(v)) * sin(u);
+ FT x3 = b * sin(v) * cos(u / 2);
+ FT x4 = b * sin(v) * sin(u / 2);
+ FT x5 = x1 + x2 + x3 + x4;
+
+ Point p = construct_point(k, x1, x2, x3, x4, x5);
+ points.push_back(p);
+ ++i;
+ }
+ return points;
+}
+
+} // namespace Gudhi
+
+#endif // RANDOM_POINT_GENERATORS_H_
diff --git a/src/common/include/gudhi/reader_utils.h b/src/common/include/gudhi/reader_utils.h
new file mode 100644
index 00000000..98335552
--- /dev/null
+++ b/src/common/include/gudhi/reader_utils.h
@@ -0,0 +1,357 @@
+/* This file is part of the Gudhi Library - https://gudhi.inria.fr/ - which is released under MIT.
+ * See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details.
+ * Author(s): Clement Maria, Pawel Dlotko, Clement Jamin
+ *
+ * Copyright (C) 2014 Inria
+ *
+ * Modification(s):
+ * - YYYY/MM Author: Description of the modification
+ */
+
+#ifndef READER_UTILS_H_
+#define READER_UTILS_H_
+
+#include <gudhi/graph_simplicial_complex.h>
+#include <gudhi/Debug_utils.h>
+
+#include <boost/function_output_iterator.hpp>
+#include <boost/graph/adjacency_list.hpp>
+
+#include <iostream>
+#include <fstream>
+#include <map>
+#include <limits> // for numeric_limits
+#include <string>
+#include <vector>
+#include <utility> // for pair
+#include <tuple> // for std::make_tuple
+
+namespace Gudhi {
+
+// Keep this file tag for Doxygen to parse the code, otherwise, functions are not documented.
+// It is required for global functions and variables.
+
+/** @file
+ * @brief This file includes common file reader for GUDHI
+ */
+
+/**
+ * @brief Read a set of points to turn it into a vector< vector<double> > by filling points.
+ *
+ * File format: 1 point per line<br>
+ * X11 X12 ... X1d<br>
+ * X21 X22 ... X2d<br>
+ * etc<br>
+ */
+inline void read_points(std::string file_name, std::vector<std::vector<double>>& points) {
+ std::ifstream in_file(file_name.c_str(), std::ios::in);
+ if (!in_file.is_open()) {
+ std::cerr << "Unable to open file " << file_name << std::endl;
+ return;
+ }
+
+ std::string line;
+ double x;
+ while (getline(in_file, line)) {
+ std::vector<double> point;
+ std::istringstream iss(line);
+ while (iss >> x) {
+ point.push_back(x);
+ }
+ // Check for empty lines
+ if (!point.empty()) points.push_back(point);
+ }
+ in_file.close();
+}
+
+/**
+ * @brief Read a graph from a file.
+ *
+ * \tparam Graph_t Type for the return graph. Must be constructible from iterators on pairs of Vertex_handle
+ * \tparam Filtration_value Type for the value of the read filtration
+ * \tparam Vertex_handle Type for the value of the read vertices
+ *
+ * File format: 1 simplex per line<br>
+ * Dim1 X11 X12 ... X1d Fil1<br>
+ * Dim2 X21 X22 ... X2d Fil2<br>
+ * etc<br>
+ *
+ * The vertices must be labeled from 0 to n-1.
+ * Every simplex must appear exactly once.
+ * Simplices of dimension more than 1 are ignored.
+ */
+template <typename Graph_t, typename Filtration_value, typename Vertex_handle>
+Graph_t read_graph(std::string file_name) {
+ std::ifstream in_(file_name.c_str(), std::ios::in);
+ if (!in_.is_open()) {
+ std::string error_str("read_graph - Unable to open file ");
+ error_str.append(file_name);
+ std::cerr << error_str << std::endl;
+ throw std::invalid_argument(error_str);
+ }
+
+ typedef std::pair<Vertex_handle, Vertex_handle> Edge_t;
+ std::vector<Edge_t> edges;
+ std::vector<Filtration_value> edges_fil;
+ std::map<Vertex_handle, Filtration_value> vertices;
+
+ std::string line;
+ int dim;
+ Vertex_handle u, v, max_h = -1;
+ Filtration_value fil;
+ while (getline(in_, line)) {
+ std::istringstream iss(line);
+ while (iss >> dim) {
+ switch (dim) {
+ case 0: {
+ iss >> u;
+ iss >> fil;
+ vertices[u] = fil;
+ if (max_h < u) {
+ max_h = u;
+ }
+ break;
+ }
+ case 1: {
+ iss >> u;
+ iss >> v;
+ iss >> fil;
+ edges.push_back(Edge_t(u, v));
+ edges_fil.push_back(fil);
+ break;
+ }
+ default: { break; }
+ }
+ }
+ }
+ in_.close();
+
+ if ((size_t)(max_h + 1) != vertices.size()) {
+ std::cerr << "Error: vertices must be labeled from 0 to n-1 \n";
+ }
+
+ Graph_t skel_graph(edges.begin(), edges.end(), edges_fil.begin(), vertices.size());
+ auto vertex_prop = boost::get(vertex_filtration_t(), skel_graph);
+
+ typename boost::graph_traits<Graph_t>::vertex_iterator vi, vi_end;
+ auto v_it = vertices.begin();
+ for (std::tie(vi, vi_end) = boost::vertices(skel_graph); vi != vi_end; ++vi, ++v_it) {
+ boost::put(vertex_prop, *vi, v_it->second);
+ }
+
+ return skel_graph;
+}
+
+/**
+ * @brief Read a face from a file.
+ *
+ * File format: 1 simplex per line<br>
+ * Dim1 X11 X12 ... X1d Fil1<br>
+ * Dim2 X21 X22 ... X2d Fil2<br>
+ * etc<br>
+ *
+ * The vertices must be labeled from 0 to n-1.
+ * Every simplex must appear exactly once.
+ * Simplices of dimension more than 1 are ignored.
+ */
+template <typename Vertex_handle, typename Filtration_value>
+bool read_simplex(std::istream& in_, std::vector<Vertex_handle>& simplex, Filtration_value& fil) {
+ int dim = 0;
+ if (!(in_ >> dim)) return false;
+ Vertex_handle v;
+ for (int i = 0; i < dim + 1; ++i) {
+ if (!(in_ >> v)) return false;
+ simplex.push_back(v);
+ }
+ if (!(in_ >> fil)) return false;
+ in_.ignore((std::numeric_limits<std::streamsize>::max)(), '\n'); // ignore until the carriage return
+ return true;
+}
+
+/**
+ * @brief Read a hasse simplex from a file.
+ *
+ * File format: 1 simplex per line<br>
+ * Dim1 k11 k12 ... k1Dim1 Fil1<br>
+ * Dim2 k21 k22 ... k2Dim2 Fil2<br>
+ * etc<br>
+ *
+ * The key of a simplex is its position in the filtration order and also the number of its row in the file.
+ * Dimi ki1 ki2 ... kiDimi Fili means that the ith simplex in the filtration has dimension Dimi, filtration value
+ * fil1 and simplices with key ki1 ... kiDimi in its boundary.*/
+template <typename Simplex_key, typename Filtration_value>
+bool read_hasse_simplex(std::istream& in_, std::vector<Simplex_key>& boundary, Filtration_value& fil) {
+ int dim;
+ if (!(in_ >> dim)) return false;
+ if (dim == 0) {
+ in_ >> fil;
+ return true;
+ }
+ Simplex_key key;
+ for (int i = 0; i < dim + 1; ++i) {
+ in_ >> key;
+ boundary.push_back(key);
+ }
+ in_ >> fil;
+ return true;
+}
+
+/**
+ * @brief Read a lower triangular distance matrix from a csv file. We assume that the .csv store the whole
+ * (square) matrix.
+ *
+ * @author Pawel Dlotko
+ *
+ * Square matrix file format:<br>
+ * 0;D12;...;D1j<br>
+ * D21;0;...;D2j<br>
+ * ...<br>
+ * Dj1;Dj2;...;0<br>
+ *
+ * lower matrix file format:<br>
+ * 0<br>
+ * D21;<br>
+ * D31;D32;<br>
+ * ...<br>
+ * Dj1;Dj2;...;Dj(j-1);<br>
+ *
+ **/
+template <typename Filtration_value>
+std::vector<std::vector<Filtration_value>> read_lower_triangular_matrix_from_csv_file(const std::string& filename,
+ const char separator = ';') {
+#ifdef DEBUG_TRACES
+ std::cout << "Using procedure read_lower_triangular_matrix_from_csv_file \n";
+#endif // DEBUG_TRACES
+ std::vector<std::vector<Filtration_value>> result;
+ std::ifstream in;
+ in.open(filename.c_str());
+ if (!in.is_open()) {
+ return result;
+ }
+
+ std::string line;
+
+ // the first line is emtpy, so we ignore it:
+ std::getline(in, line);
+ std::vector<Filtration_value> values_in_this_line;
+ result.push_back(values_in_this_line);
+
+ int number_of_line = 0;
+
+ // first, read the file line by line to a string:
+ while (std::getline(in, line)) {
+ // if line is empty, break
+ if (line.size() == 0) break;
+
+ // if the last element of a string is comma:
+ if (line[line.size() - 1] == separator) {
+ // then shrink the string by one
+ line.pop_back();
+ }
+
+ // replace all commas with spaces
+ std::replace(line.begin(), line.end(), separator, ' ');
+
+ // put the new line to a stream
+ std::istringstream iss(line);
+ // and now read the doubles.
+
+ int number_of_entry = 0;
+ std::vector<Filtration_value> values_in_this_line;
+ while (iss.good()) {
+ double entry;
+ iss >> entry;
+ if (number_of_entry <= number_of_line) {
+ values_in_this_line.push_back(entry);
+ }
+ ++number_of_entry;
+ }
+ if (!values_in_this_line.empty()) result.push_back(values_in_this_line);
+ ++number_of_line;
+ }
+ in.close();
+
+#ifdef DEBUG_TRACES
+ std::cerr << "Here is the matrix we read : \n";
+ for (size_t i = 0; i != result.size(); ++i) {
+ for (size_t j = 0; j != result[i].size(); ++j) {
+ std::cerr << result[i][j] << " ";
+ }
+ std::cerr << std::endl;
+ }
+#endif // DEBUG_TRACES
+
+ return result;
+} // read_lower_triangular_matrix_from_csv_file
+
+/**
+Reads a file containing persistence intervals.
+Each line might contain 2, 3 or 4 values: [[field] dimension] birth death
+The output iterator `out` is used this way: `*out++ = std::make_tuple(dim, birth, death);`
+where `dim` is an `int`, `birth` a `double`, and `death` a `double`.
+Note: the function does not check that birth <= death.
+**/
+template <typename OutputIterator>
+void read_persistence_intervals_and_dimension(std::string const& filename, OutputIterator out) {
+ std::ifstream in(filename);
+ if (!in.is_open()) {
+ std::string error_str("read_persistence_intervals_and_dimension - Unable to open file ");
+ error_str.append(filename);
+ std::cerr << error_str << std::endl;
+ throw std::invalid_argument(error_str);
+ }
+
+ while (!in.eof()) {
+ std::string line;
+ getline(in, line);
+ if (line.length() != 0 && line[0] != '#') {
+ double numbers[4];
+ int n = sscanf(line.c_str(), "%lf %lf %lf %lf", &numbers[0], &numbers[1], &numbers[2], &numbers[3]);
+ if (n >= 2) {
+ int dim = (n >= 3 ? static_cast<int>(numbers[n - 3]) : -1);
+ *out++ = std::make_tuple(dim, numbers[n - 2], numbers[n - 1]);
+ }
+ }
+ }
+}
+
+/**
+Reads a file containing persistence intervals.
+Each line might contain 2, 3 or 4 values: [[field] dimension] birth death
+The return value is an `std::map<dim, std::vector<std::pair<birth, death>>>`
+where `dim` is an `int`, `birth` a `double`, and `death` a `double`.
+Note: the function does not check that birth <= death.
+**/
+inline std::map<int, std::vector<std::pair<double, double>>> read_persistence_intervals_grouped_by_dimension(
+ std::string const& filename) {
+ std::map<int, std::vector<std::pair<double, double>>> ret;
+ read_persistence_intervals_and_dimension(
+ filename, boost::make_function_output_iterator([&ret](std::tuple<int, double, double> t) {
+ ret[get<0>(t)].push_back(std::make_pair(get<1>(t), get<2>(t)));
+ }));
+ return ret;
+}
+
+/**
+Reads a file containing persistence intervals.
+Each line might contain 2, 3 or 4 values: [[field] dimension] birth death
+If `only_this_dim` = -1, dimension is ignored and all lines are returned.
+If `only_this_dim` is >= 0, only the lines where dimension = `only_this_dim`
+(or where dimension is not specified) are returned.
+The return value is an `std::vector<std::pair<birth, death>>`
+where `dim` is an `int`, `birth` a `double`, and `death` a `double`.
+Note: the function does not check that birth <= death.
+**/
+inline std::vector<std::pair<double, double>> read_persistence_intervals_in_dimension(std::string const& filename,
+ int only_this_dim = -1) {
+ std::vector<std::pair<double, double>> ret;
+ read_persistence_intervals_and_dimension(
+ filename, boost::make_function_output_iterator([only_this_dim, &ret](std::tuple<int, double, double> t) {
+ if (only_this_dim == get<0>(t) || only_this_dim == -1) ret.emplace_back(get<1>(t), get<2>(t));
+ }));
+ return ret;
+}
+
+} // namespace Gudhi
+
+#endif // READER_UTILS_H_
diff --git a/src/common/include/gudhi/writing_persistence_to_file.h b/src/common/include/gudhi/writing_persistence_to_file.h
new file mode 100644
index 00000000..2e36b831
--- /dev/null
+++ b/src/common/include/gudhi/writing_persistence_to_file.h
@@ -0,0 +1,105 @@
+/* This file is part of the Gudhi Library - https://gudhi.inria.fr/ - which is released under MIT.
+ * See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details.
+ * Author(s): Pawel Dlotko
+ *
+ * Copyright (C) 2017 Swansea University, UK
+ *
+ * Modification(s):
+ * - YYYY/MM Author: Description of the modification
+ */
+
+#ifndef WRITING_PERSISTENCE_TO_FILE_H_
+#define WRITING_PERSISTENCE_TO_FILE_H_
+
+#include <iostream>
+#include <string>
+#include <limits>
+
+namespace Gudhi {
+
+/**
+* This is a class to store persistence intervals. Its main purpose is to
+* exchange data in between different packages and provide unified way
+* of writing a collection of persistence intervals to file.
+**/
+template <typename Filtration_type, typename Coefficient_field>
+class Persistence_interval_common {
+ public:
+ /**
+ * Constructor taking as an input birth and death of the pair.
+ **/
+ Persistence_interval_common(Filtration_type birth, Filtration_type death)
+ : birth_(birth),
+ death_(death),
+ dimension_(std::numeric_limits<unsigned>::max()),
+ arith_element_(std::numeric_limits<Coefficient_field>::max()) {}
+
+ /**
+ * Constructor taking as an input birth, death and dimension of the pair.
+ **/
+ Persistence_interval_common(Filtration_type birth, Filtration_type death, unsigned dim)
+ : birth_(birth), death_(death), dimension_(dim), arith_element_(std::numeric_limits<Coefficient_field>::max()) {}
+
+ /**
+* Constructor taking as an input birth, death, dimension of the pair as well
+* as the number p such that this interval is present over Z_p field.
+**/
+ Persistence_interval_common(Filtration_type birth, Filtration_type death, unsigned dim, Coefficient_field field)
+ : birth_(birth), death_(death), dimension_(dim), arith_element_(field) {}
+
+ /**
+ * Operator to compare two persistence pairs. During the comparision all the
+ * fields: birth, death, dimensiona and arith_element_ are taken into account
+ * and they all have to be equal for two pairs to be equal.
+ **/
+ bool operator==(const Persistence_interval_common& i2) const {
+ return ((this->birth_ == i2.birth_) && (this->death_ == i2.death_) && (this->dimension_ == i2.dimension_) &&
+ (this->arith_element_ == i2.arith_element_));
+ }
+
+ /**
+ * Check if two persistence paris are not equal.
+ **/
+ bool operator!=(const Persistence_interval_common& i2) const { return (!((*this) == i2)); }
+
+ /**
+ * Operator to compare objects of a type Persistence_interval_common.
+ * One intervals is smaller than the other if it has lower persistence.
+ * Note that this operator do not take Arith_element into account when doing comparisions.
+ **/
+ bool operator<(const Persistence_interval_common& i2) const {
+ return fabs(this->death_ - this->birth_) < fabs(i2.death_ - i2.birth_);
+ }
+
+ friend std::ostream& operator<<(std::ostream& out, const Persistence_interval_common& it) {
+ if (it.arith_element_ != std::numeric_limits<Coefficient_field>::max()) {
+ out << it.arith_element_ << " ";
+ }
+ if (it.dimension_ != std::numeric_limits<unsigned>::max()) {
+ out << it.dimension_ << " ";
+ }
+ out << it.birth_ << " " << it.death_ << " ";
+ return out;
+ }
+
+ private:
+ Filtration_type birth_;
+ Filtration_type death_;
+ unsigned dimension_;
+ Coefficient_field arith_element_;
+};
+
+/**
+ * This function write a vector<Persistence_interval_common> to a stream
+**/
+template <typename Persistence_interval_range>
+void write_persistence_intervals_to_stream(const Persistence_interval_range& intervals,
+ std::ostream& out = std::cout) {
+ for (auto interval : intervals) {
+ out << interval << "\n";
+ }
+}
+
+} // namespace Gudhi
+
+#endif // WRITING_PERSISTENCE_TO_FILE_H_
generated by cgit v1.2.3 (git 2.39.1) at 2024-05-17 09:04:44 +0000