summaryrefslogtreecommitdiff
path: root/src/GudhUI/utils
diff options
context:
space:
mode:
Diffstat (limited to 'src/GudhUI/utils')
-rw-r--r--src/GudhUI/utils/Bar_code_persistence.h100
-rw-r--r--src/GudhUI/utils/Critical_points.h120
-rw-r--r--src/GudhUI/utils/Edge_collapsor.h84
-rw-r--r--src/GudhUI/utils/Edge_contractor.h84
-rw-r--r--src/GudhUI/utils/Furthest_point_epsilon_net.h119
-rw-r--r--src/GudhUI/utils/Is_manifold.h88
-rw-r--r--src/GudhUI/utils/K_nearest_builder.h77
-rw-r--r--src/GudhUI/utils/Lloyd_builder.h78
-rw-r--r--src/GudhUI/utils/MClock.h57
-rw-r--r--src/GudhUI/utils/Persistence_compute.h81
-rw-r--r--src/GudhUI/utils/Rips_builder.h56
-rw-r--r--src/GudhUI/utils/UI_utils.h32
-rw-r--r--src/GudhUI/utils/Vertex_collapsor.h76
13 files changed, 1052 insertions, 0 deletions
diff --git a/src/GudhUI/utils/Bar_code_persistence.h b/src/GudhUI/utils/Bar_code_persistence.h
new file mode 100644
index 00000000..cd9b009f
--- /dev/null
+++ b/src/GudhUI/utils/Bar_code_persistence.h
@@ -0,0 +1,100 @@
+/* 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
+ */
+
+#include <math.h> // isfinite
+
+#include <QtGui/QApplication>
+
+#include <QGraphicsView>
+#include <QGraphicsScene>
+#include <QPointF>
+#include <QVector>
+#include <QGraphicsTextItem>
+
+#include <iostream>
+#include <vector>
+#include <limits> // NaN, infinity
+#include <utility> // for pair
+#include <string>
+
+#ifndef UTILS_BAR_CODE_PERSISTENCE_H_
+#define UTILS_BAR_CODE_PERSISTENCE_H_
+
+class Bar_code_persistence {
+ private:
+ typedef std::vector<std::pair<double, double>> Persistence;
+ Persistence persistence_vector;
+ double min_birth;
+ double max_death;
+
+ public:
+ Bar_code_persistence()
+ : min_birth(std::numeric_limits<double>::quiet_NaN()),
+ max_death(std::numeric_limits<double>::quiet_NaN()) { }
+
+ void insert(double birth, double death) {
+ persistence_vector.push_back(std::make_pair(birth, death));
+ if (std::isfinite(birth)) {
+ if ((birth < min_birth) || (std::isnan(min_birth)))
+ min_birth = birth;
+ if ((birth > max_death) || (std::isnan(max_death)))
+ max_death = birth;
+ }
+ if (std::isfinite(death))
+ if ((death > max_death) || (std::isnan(max_death)))
+ max_death = death;
+ }
+
+ void show(const std::string& window_title) {
+ // Create a view, put a scene in it
+ QGraphicsView * view = new QGraphicsView();
+ QGraphicsScene * scene = new QGraphicsScene();
+ view->setScene(scene);
+ double ratio = 600.0 / (max_death - min_birth);
+ // std::cout << "min_birth=" << min_birth << " - max_death=" << max_death << " - ratio=" << ratio << std::endl;
+
+ double height = 0.0, birth = 0.0, death = 0.0;
+ int pers_num = 1;
+ for (auto& persistence : persistence_vector) {
+ height = 5.0 * pers_num;
+ // std::cout << "[" << pers_num << "] birth=" << persistence.first << " - death=" << persistence.second << std::endl;
+ if (std::isfinite(persistence.first))
+ birth = ((persistence.first - min_birth) * ratio) + 50.0;
+ else
+ birth = 0.0;
+
+ if (std::isfinite(persistence.second))
+ death = ((persistence.second - min_birth) * ratio) + 50.0;
+ else
+ death = 700.0;
+
+ scene->addLine(birth, height, death, height, QPen(Qt::blue, 2));
+ pers_num++;
+ }
+ height += 10.0;
+ // scale line
+ scene->addLine(0, height, 700.0, height, QPen(Qt::black, 1));
+ int modulo = 0;
+ for (double scale = 50.0; scale < 700.0; scale += 50.0) {
+ modulo++;
+ // scale small dash
+ scene->addLine(scale, height - 3.0, scale, height + 3.0, QPen(Qt::black, 1));
+ // scale text
+ QString scale_value = QString::number(((scale - 50.0) / ratio) + min_birth);
+ QGraphicsTextItem* dimText = scene->addText(scale_value, QFont("Helvetica", 8));
+ dimText->setPos(scale - (3.0 * scale_value.size()), height + 9.0 * (modulo % 2));
+ }
+ view->setWindowTitle(window_title.c_str());
+ // Show the view
+ view->show();
+ }
+};
+
+#endif // UTILS_BAR_CODE_PERSISTENCE_H_
diff --git a/src/GudhUI/utils/Critical_points.h b/src/GudhUI/utils/Critical_points.h
new file mode 100644
index 00000000..32fcf32e
--- /dev/null
+++ b/src/GudhUI/utils/Critical_points.h
@@ -0,0 +1,120 @@
+/* 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 UTILS_CRITICAL_POINTS_H_
+#define UTILS_CRITICAL_POINTS_H_
+
+#include <deque>
+#include <utility> // for pair<>
+#include <algorithm> // for sort
+
+#include "utils/Edge_contractor.h"
+
+/**
+ * Iteratively tries to anticollapse smallest edge non added so far.
+ * If its link is contractible then no topological change and else possible topological change.
+ *
+ * todo do a sparsification with some parameter eps while growing
+ */
+template<typename SkBlComplex> class Critical_points {
+ private:
+ SkBlComplex filled_complex_;
+ const SkBlComplex& input_complex_;
+ double max_length_;
+ std::ostream& stream_;
+
+ public:
+ typedef typename SkBlComplex::Vertex_handle Vertex_handle;
+ typedef typename SkBlComplex::Edge_handle Edge_handle;
+ typedef typename std::pair<Vertex_handle, Vertex_handle> Edge;
+
+ /**
+ * @brief check all pair of points with length smaller than max_length
+ */
+ Critical_points(const SkBlComplex& input_complex, std::ostream& stream, double max_length) :
+ input_complex_(input_complex), max_length_(max_length), stream_(stream) {
+ std::deque<Edge> edges;
+ auto vertices = input_complex.vertex_range();
+ for (auto p = vertices.begin(); p != vertices.end(); ++p) {
+ filled_complex_.add_vertex(input_complex.point(*p));
+ for (auto q = p; ++q != vertices.end(); /**/)
+ if (squared_eucl_distance(input_complex.point(*p), input_complex.point(*q)) < max_length_ * max_length_)
+ edges.emplace_back(*p, *q);
+ }
+
+ std::sort(edges.begin(), edges.end(),
+ [&](Edge e1, Edge e2) {
+ return squared_edge_length(e1) < squared_edge_length(e2);
+ });
+
+ anti_collapse_edges(edges);
+ }
+
+ private:
+ double squared_eucl_distance(const Point& p1, const Point& p2) const {
+ return Geometry_trait::Squared_distance_d()(p1, p2);
+ }
+
+ void anti_collapse_edges(const std::deque<Edge>& edges) {
+ unsigned pos = 0;
+ for (Edge e : edges) {
+ std::cout << "edge " << pos++ << "/" << edges.size() << "\n";
+ auto eh = filled_complex_.add_edge_without_blockers(e.first, e.second);
+ int is_contractible(is_link_reducible(eh));
+
+ switch (is_contractible) {
+ case 0:
+ stream_ << "alpha=" << std::sqrt(squared_edge_length(e)) << " topological change" << std::endl;
+ break;
+ case 2:
+ stream_ << "alpha=" << std::sqrt(squared_edge_length(e)) << " maybe a topological change" << std::endl;
+ break;
+ default:
+ break;
+ }
+ }
+ }
+
+ // 0 -> not
+ // 1 -> yes
+ // 2 -> maybe
+
+ int is_link_reducible(Edge_handle e) {
+ auto link = filled_complex_.link(e);
+
+ if (link.empty())
+ return 0;
+
+ Edge_contractor<Complex> contractor(link, link.num_vertices() - 1);
+ (void)contractor;
+
+ if (link.num_connected_components() > 1)
+ // one than more CC -> not contractible
+ return 0;
+
+ if (link.num_vertices() == 1)
+ // reduced to one point -> contractible
+ return 1;
+ else
+ // we dont know
+ return 2;
+ }
+
+ double squared_edge_length(Edge_handle e) const {
+ return squared_eucl_distance(input_complex_.point(input_complex_.first_vertex(e)),
+ input_complex_.point(input_complex_.second_vertex(e)));
+ }
+
+ double squared_edge_length(Edge e) const {
+ return squared_eucl_distance(input_complex_.point(e.first), input_complex_.point(e.second));
+ }
+};
+
+#endif // UTILS_CRITICAL_POINTS_H_
diff --git a/src/GudhUI/utils/Edge_collapsor.h b/src/GudhUI/utils/Edge_collapsor.h
new file mode 100644
index 00000000..89e032f0
--- /dev/null
+++ b/src/GudhUI/utils/Edge_collapsor.h
@@ -0,0 +1,84 @@
+/* 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 UTILS_EDGE_COLLAPSOR_H_
+#define UTILS_EDGE_COLLAPSOR_H_
+
+#include <list>
+#include "utils/Edge_contractor.h"
+#include "utils/UI_utils.h"
+
+/**
+ * Iteratively puts every vertex at the center of its neighbors
+ */
+template<typename SkBlComplex> class Edge_collapsor {
+ private:
+ SkBlComplex& complex_;
+ unsigned num_collapses_;
+
+ public:
+ typedef typename SkBlComplex::Vertex_handle Vertex_handle;
+ typedef typename SkBlComplex::Edge_handle Edge_handle;
+
+ /**
+ * @brief Collapse num_collapses edges. If num_collapses<0 then it collapses all possible edges.
+ * Largest edges are collapsed first.
+ *
+ */
+ Edge_collapsor(SkBlComplex& complex, unsigned num_collapses) :
+ complex_(complex), num_collapses_(num_collapses) {
+ std::list<Edge_handle> edges;
+ edges.insert(edges.begin(), complex_.edge_range().begin(), complex_.edge_range().end());
+
+ edges.sort(
+ [&](Edge_handle e1, Edge_handle e2) {
+ return squared_edge_length(e1) < squared_edge_length(e2);
+ });
+
+ collapse_edges(edges);
+ }
+
+ private:
+ void collapse_edges(std::list<Edge_handle>& edges) {
+ while (!edges.empty() && num_collapses_--) {
+ Edge_handle current_edge = edges.front();
+ edges.pop_front();
+ if (is_link_reducible(current_edge))
+ complex_.remove_edge(current_edge);
+ }
+ }
+
+ bool is_link_reducible(Edge_handle e) {
+ auto link = complex_.link(e);
+
+ if (link.empty())
+ return false;
+
+ if (link.is_cone())
+ return true;
+
+ if (link.num_connected_components() > 1)
+ return false;
+
+ Edge_contractor<SkBlComplex> contractor(link, link.num_vertices() - 1);
+
+ return (link.num_vertices() == 1);
+ }
+
+ double squared_edge_length(Edge_handle e) const {
+ return squared_eucl_distance(complex_.point(complex_.first_vertex(e)), complex_.point(complex_.second_vertex(e)));
+ }
+
+ double squared_eucl_distance(const Point& p1, const Point& p2) const {
+ return Geometry_trait::Squared_distance_d()(p1, p2);
+ }
+};
+
+#endif // UTILS_EDGE_COLLAPSOR_H_
diff --git a/src/GudhUI/utils/Edge_contractor.h b/src/GudhUI/utils/Edge_contractor.h
new file mode 100644
index 00000000..0707b186
--- /dev/null
+++ b/src/GudhUI/utils/Edge_contractor.h
@@ -0,0 +1,84 @@
+/* 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 UTILS_EDGE_CONTRACTOR_H_
+#define UTILS_EDGE_CONTRACTOR_H_
+
+
+#include <gudhi/Skeleton_blocker_contractor.h>
+#include <gudhi/Contraction/Edge_profile.h>
+#include <gudhi/Contraction/policies/Cost_policy.h>
+
+#include <vector>
+
+/**
+ * Iteratively puts every vertex at the center of its neighbors
+ */
+template<typename SkBlComplex> class Edge_contractor {
+ private:
+ SkBlComplex& complex_;
+ unsigned num_contractions_;
+
+ /**
+ * @brief return a cost corresponding to the squared length of the edge
+ */
+ template< typename EdgeProfile> class Length_cost : public contraction::Cost_policy<EdgeProfile> {
+ public:
+ typedef typename contraction::Cost_policy<EdgeProfile>::Cost_type Cost_type;
+ typedef typename EdgeProfile::Point Point;
+
+ Cost_type operator()(const EdgeProfile& profile, const boost::optional<Point>& placement) const override {
+ Cost_type res;
+ if (!placement)
+ return res;
+ return Geometry_trait::Squared_distance_d()(profile.p0(), profile.p1()); // not working??
+ }
+ };
+
+ /**
+ * @brief return a cost corresponding to the squared length of the edge
+ */
+ template< typename EdgeProfile> class Middle_placement : public contraction::Placement_policy<EdgeProfile> {
+ public:
+ typedef typename contraction::Placement_policy<EdgeProfile>::Placement_type Placement_type;
+ typedef typename EdgeProfile::Point Point;
+
+ Placement_type operator()(const EdgeProfile& profile) const override {
+ std::vector<double> mid_coords(profile.p0().dimension(), 0);
+ for (int i = 0; i < profile.p0().dimension(); ++i) {
+ mid_coords[i] = (profile.p0()[i] + profile.p1()[i]) / 2.;
+ }
+ return Point(profile.p0().dimension(), mid_coords.begin(), mid_coords.end());
+ }
+ };
+
+ public:
+ typedef typename SkBlComplex::Vertex_handle Vertex_handle;
+ typedef typename SkBlComplex::Edge_handle Edge_handle;
+
+ /**
+ * @brief Modify complex to be the expansion of the k-nearest neighbor
+ * symetric graph.
+ */
+ Edge_contractor(SkBlComplex& complex, unsigned num_contractions) :
+ complex_(complex), num_contractions_(num_contractions) {
+ typedef typename contraction::Edge_profile<Complex> Profile;
+ num_contractions = (std::min)(static_cast<int>(num_contractions), static_cast<int>(complex_.num_vertices() - 1));
+ typedef typename contraction::Skeleton_blocker_contractor<Complex> Contractor;
+ Contractor contractor(complex_,
+ new Length_cost<contraction::Edge_profile < Complex >> (),
+ new Middle_placement<contraction::Edge_profile < Complex >> (),
+ contraction::make_link_valid_contraction<Profile>(),
+ contraction::make_remove_popable_blockers_visitor<Profile>());
+ contractor.contract_edges(num_contractions);
+ }
+};
+
+#endif // UTILS_EDGE_CONTRACTOR_H_
diff --git a/src/GudhUI/utils/Furthest_point_epsilon_net.h b/src/GudhUI/utils/Furthest_point_epsilon_net.h
new file mode 100644
index 00000000..6eb71071
--- /dev/null
+++ b/src/GudhUI/utils/Furthest_point_epsilon_net.h
@@ -0,0 +1,119 @@
+/* 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 UTILS_FURTHEST_POINT_EPSILON_NET_H_
+#define UTILS_FURTHEST_POINT_EPSILON_NET_H_
+
+#include <vector>
+#include <algorithm> // for sort
+
+#include "utils/UI_utils.h"
+
+/**
+ * Computes an epsilon net with furthest point strategy.
+ */
+template<typename SkBlComplex> class Furthest_point_epsilon_net {
+ private:
+ SkBlComplex& complex_;
+ typedef typename SkBlComplex::Vertex_handle Vertex_handle;
+ typedef typename SkBlComplex::Edge_handle Edge_handle;
+
+ /**
+ * Let V be the set of vertices.
+ * Initially v0 is one arbitrarly vertex and the set V0 is {v0}.
+ * Then Vk is computed as follows.
+ * First we compute the vertex pk that is the furthest from Vk
+ * then Vk = Vk \cup pk.
+ * The radius of pk is its distance to Vk and its meeting vertex
+ * is the vertex of Vk for which this distance is achieved.
+ */
+ struct Net_filtration_vertex {
+ Vertex_handle vertex_handle;
+ Vertex_handle meeting_vertex;
+ double radius;
+
+ Net_filtration_vertex(Vertex_handle vertex_handle_,
+ Vertex_handle meeting_vertex_,
+ double radius_) :
+ vertex_handle(vertex_handle_), meeting_vertex(meeting_vertex_), radius(radius_) { }
+
+ bool operator<(const Net_filtration_vertex& other) const {
+ return radius < other.radius;
+ }
+ };
+
+ public:
+ std::vector<Net_filtration_vertex> net_filtration_;
+
+ /**
+ * @brief Modify complex to be the expansion of the k-nearest neighbor
+ * symetric graph.
+ */
+ Furthest_point_epsilon_net(SkBlComplex& complex) :
+ complex_(complex) {
+ if (!complex.empty()) {
+ init_filtration();
+ for (std::size_t k = 2; k < net_filtration_.size(); ++k) {
+ update_radius_value(k);
+ }
+ }
+ }
+
+ // xxx does not work if complex not full
+
+ double radius(Vertex_handle v) {
+ return net_filtration_[v.vertex].radius;
+ }
+
+ private:
+ void init_filtration() {
+ Vertex_handle v0 = *(complex_.vertex_range().begin());
+ net_filtration_.reserve(complex_.num_vertices());
+ for (auto v : complex_.vertex_range()) {
+ if (v != v0)
+ net_filtration_.push_back(Net_filtration_vertex(v,
+ Vertex_handle(-1),
+ squared_eucl_distance(v, v0)));
+ }
+ net_filtration_.push_back(Net_filtration_vertex(v0, Vertex_handle(-1), 1e10));
+ auto n = net_filtration_.size();
+ sort_filtration(n - 1);
+ }
+
+ void update_radius_value(int k) {
+ int n = net_filtration_.size();
+ int index_to_update = n - k;
+ for (int i = 0; i < index_to_update; ++i) {
+ net_filtration_[i].radius = (std::min)(net_filtration_[i].radius,
+ squared_eucl_distance(net_filtration_[i].vertex_handle,
+ net_filtration_[index_to_update].vertex_handle));
+ }
+ sort_filtration(n - k);
+ }
+
+ /**
+ * sort all i first elements.
+ */
+ void sort_filtration(int i) {
+ std::sort(net_filtration_.begin(), net_filtration_.begin() + i);
+ }
+
+ double squared_eucl_distance(Vertex_handle v1, Vertex_handle v2) const {
+ return std::sqrt(Geometry_trait::Squared_distance_d()(complex_.point(v1), complex_.point(v2)));
+ }
+
+ void print_filtration() const {
+ for (auto v : net_filtration_) {
+ std::cerr << "v=" << v.vertex_handle << "-> d=" << v.radius << std::endl;
+ }
+ }
+};
+
+#endif // UTILS_FURTHEST_POINT_EPSILON_NET_H_
diff --git a/src/GudhUI/utils/Is_manifold.h b/src/GudhUI/utils/Is_manifold.h
new file mode 100644
index 00000000..2eb79d21
--- /dev/null
+++ b/src/GudhUI/utils/Is_manifold.h
@@ -0,0 +1,88 @@
+/* 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 UTILS_IS_MANIFOLD_H_
+#define UTILS_IS_MANIFOLD_H_
+
+#include "utils/UI_utils.h"
+#include "utils/Edge_contractor.h"
+
+/**
+ * Iteratively tries to anticollapse smallest edge non added so far.
+ * If its link is contractible then no topological change and else possible topological change.
+ *
+ * todo do a sparsification with some parameter eps while growing
+ */
+template<typename SkBlComplex> class Is_manifold {
+ private:
+ const SkBlComplex& input_complex_;
+ typedef typename SkBlComplex::Vertex_handle Vertex_handle;
+
+ public:
+ /*
+ * return dim the maximum dimension around one simplex and res which is true if the complex is a manifold.
+ * If the complex has dimension <= 3 then if res is false, the complex is not a manifold.
+ * For d-manifold with d>=4, res may be false while the complex is a manifold.
+ */
+ Is_manifold(const SkBlComplex& input_complex, unsigned& dim, bool & res) : input_complex_(input_complex) {
+ res = true;
+ dim = -1;
+ if (!input_complex_.empty()) {
+ for (auto v : input_complex_.vertex_range()) {
+ dim = local_dimension(v);
+ break;
+ }
+ // check that the link of every vertex is a dim-1 sphere
+ for (auto v : input_complex_.vertex_range()) {
+ if (!is_k_sphere(v, dim - 1)) {
+ res = false;
+ break;
+ }
+ }
+ }
+ }
+
+ private:
+ unsigned local_dimension(Vertex_handle v) {
+ unsigned dim = 0;
+ for (const auto& s : input_complex_.star_simplex_range(v))
+ dim = (std::max)(dim, (unsigned) s.dimension());
+ return dim;
+ }
+
+ bool is_k_sphere(Vertex_handle v, int k) {
+ auto link = input_complex_.link(v);
+ Edge_contractor<Complex> contractor(link, link.num_vertices() - 1);
+ (void)contractor;
+ return (is_sphere_simplex(link) == k);
+ }
+
+ // A minimal sphere is a complex that contains vertices v1...vn and all faces
+ // made upon this set except the face {v1,...,vn}
+ // return -2 if not a minimal sphere
+ // and d otherwise if complex is a d minimal sphere
+
+ template<typename SubComplex>
+ int is_sphere_simplex(const SubComplex& complex) {
+ if (complex.empty()) return -1;
+ if (complex.num_blockers() != 1) return -2;
+
+ // necessary and sufficient condition : there exists a unique blocker that passes through all vertices
+ auto first_blocker = *(complex.const_blocker_range().begin());
+
+ if (first_blocker->dimension() + 1 != complex.num_vertices())
+ return -2;
+ else
+ return (first_blocker->dimension() - 1);
+ }
+};
+
+#endif // UTILS_IS_MANIFOLD_H_
diff --git a/src/GudhUI/utils/K_nearest_builder.h b/src/GudhUI/utils/K_nearest_builder.h
new file mode 100644
index 00000000..34483e58
--- /dev/null
+++ b/src/GudhUI/utils/K_nearest_builder.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 UTILS_K_NEAREST_BUILDER_H_
+#define UTILS_K_NEAREST_BUILDER_H_
+
+#include <CGAL/Euclidean_distance.h>
+#include <CGAL/Orthogonal_k_neighbor_search.h>
+#include <CGAL/Search_traits_d.h>
+#include <CGAL/Search_traits_adapter.h>
+#include <CGAL/property_map.h>
+
+#include <unordered_map>
+#include <list>
+#include <utility>
+
+#include "utils/UI_utils.h"
+#include "model/Complex_typedefs.h"
+
+template<typename SkBlComplex> class K_nearest_builder {
+ private:
+ typedef Geometry_trait Kernel;
+ typedef Point Point_d;
+ typedef std::pair<Point_d, unsigned> Point_d_with_id;
+ typedef CGAL::Search_traits_d<Kernel> Traits_base;
+ typedef CGAL::Search_traits_adapter<Point_d_with_id, CGAL::First_of_pair_property_map<Point_d_with_id>,
+ Traits_base> Traits;
+ typedef CGAL::Orthogonal_k_neighbor_search<Traits> Neighbor_search;
+ typedef Neighbor_search::Tree Tree;
+ typedef Neighbor_search::Distance Distance;
+
+ SkBlComplex& complex_;
+
+ public:
+ /**
+ * @brief Modify complex to be the expansion of the k-nearest neighbor
+ * symetric graph.
+ */
+ K_nearest_builder(SkBlComplex& complex, unsigned k) : complex_(complex) {
+ complex.keep_only_vertices();
+ compute_edges(k);
+ }
+
+ private:
+ double squared_eucl_distance(const Point& p1, const Point& p2) const {
+ return Geometry_trait::Squared_distance_d()(p1, p2);
+ }
+
+ void compute_edges(unsigned k) {
+ std::list<Point_d_with_id> points_with_id;
+ for (auto v : complex_.vertex_range()) {
+ Point_d_with_id point_v(complex_.point(v), v.vertex);
+ points_with_id.push_back(point_v);
+ }
+
+ Tree tree(points_with_id.begin(), points_with_id.end());
+
+ typedef typename SkBlComplex::Vertex_handle Vertex_handle;
+ for (auto p : complex_.vertex_range()) {
+ Neighbor_search search(tree, complex_.point(p), k + 1);
+ for (auto it = ++search.begin(); it != search.end(); ++it) {
+ Vertex_handle q(std::get<1>(it->first));
+ if (p != q && complex_.contains_vertex(p) && complex_.contains_vertex(q))
+ complex_.add_edge_without_blockers(p, q);
+ }
+ }
+ }
+};
+
+#endif // UTILS_K_NEAREST_BUILDER_H_
diff --git a/src/GudhUI/utils/Lloyd_builder.h b/src/GudhUI/utils/Lloyd_builder.h
new file mode 100644
index 00000000..c042564f
--- /dev/null
+++ b/src/GudhUI/utils/Lloyd_builder.h
@@ -0,0 +1,78 @@
+/* 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 UTILS_LLOYD_BUILDER_H_
+#define UTILS_LLOYD_BUILDER_H_
+
+#include <vector>
+#include <list>
+
+/**
+ * Iteratively puts every vertex at the center of its neighbors
+ */
+template<typename SkBlComplex> class Lloyd_builder {
+ private:
+ SkBlComplex& complex_;
+ int dim;
+
+ public:
+ typedef typename SkBlComplex::Vertex_handle Vertex_handle;
+
+ /**
+ * @brief Modify complex to be the expansion of the k-nearest neighbor
+ * symetric graph.
+ */
+ Lloyd_builder(SkBlComplex& complex, unsigned num_iterations) : complex_(complex), dim(-1) {
+ if (!complex_.empty()) {
+ dim = get_dimension();
+ while (num_iterations--) {
+ std::list<Point> new_points;
+ for (auto v : complex.vertex_range())
+ new_points.push_back(barycenter_neighbors(v));
+
+ auto new_points_it = new_points.begin();
+ for (auto v : complex.vertex_range())
+ complex_.point(v) = *(new_points_it++);
+ }
+ }
+ }
+
+ private:
+ int get_dimension() {
+ assert(!complex_.empty());
+ for (auto v : complex_.vertex_range())
+ return complex_.point(v).dimension();
+ return -1;
+ }
+
+ Point barycenter_neighbors(Vertex_handle v) const {
+ if (complex_.degree(v) == 0)
+ return complex_.point(v);
+
+ std::vector<double> res(dim, 0);
+ unsigned num_points = 0;
+ for (auto nv : complex_.vertex_range(v)) {
+ ++num_points;
+ const Point& point = complex_.point(nv);
+ assert(point.dimension() == dim);
+ for (int i = 0; i < point.dimension(); ++i)
+ res[i] += point[i];
+ }
+ for (auto& x : res)
+ x /= num_points;
+ return Point(dim, res.begin(), res.end());
+ }
+
+ double squared_eucl_distance(const Point& p1, const Point& p2) const {
+ return Geometry_trait::Squared_distance_d()(p1, p2);
+ }
+};
+
+#endif // UTILS_LLOYD_BUILDER_H_
diff --git a/src/GudhUI/utils/MClock.h b/src/GudhUI/utils/MClock.h
new file mode 100644
index 00000000..54b28211
--- /dev/null
+++ b/src/GudhUI/utils/MClock.h
@@ -0,0 +1,57 @@
+/* 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 UTILS_MCLOCK_H_
+#define UTILS_MCLOCK_H_
+
+#include <sys/time.h>
+
+class MClock {
+ public:
+ MClock() {
+ end_called = false;
+ begin();
+ }
+
+ void begin() const {
+ end_called = false;
+ gettimeofday(&startTime, NULL);
+ }
+
+ void end() const {
+ end_called = true;
+ gettimeofday(&endTime, NULL);
+ }
+
+ friend std::ostream& operator<<(std::ostream& stream, const MClock& clock) {
+ // if(!clock.end_called) clock.end();
+ if (!clock.end_called) {
+ stream << "end not called";
+ } else {
+ long totalTime = (clock.endTime.tv_sec - clock.startTime.tv_sec) * 1000000L;
+ totalTime += (clock.endTime.tv_usec - clock.startTime.tv_usec);
+ stream << clock.num_seconds() << "s";
+ }
+ return stream;
+ }
+
+ double num_seconds() const {
+ if (!end_called) end();
+ long totalTime = (endTime.tv_sec - startTime.tv_sec) * 1000000L;
+ totalTime += (endTime.tv_usec - startTime.tv_usec);
+ return (totalTime / 1000L) / static_cast<double>(1000);
+ }
+
+ private:
+ mutable struct timeval startTime, endTime;
+ mutable bool end_called;
+};
+
+#endif // UTILS_MCLOCK_H_
diff --git a/src/GudhUI/utils/Persistence_compute.h b/src/GudhUI/utils/Persistence_compute.h
new file mode 100644
index 00000000..874f31cf
--- /dev/null
+++ b/src/GudhUI/utils/Persistence_compute.h
@@ -0,0 +1,81 @@
+/* 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 UTILS_PERSISTENCE_COMPUTE_H_
+#define UTILS_PERSISTENCE_COMPUTE_H_
+
+#include <gudhi/graph_simplicial_complex.h>
+#include <gudhi/Simplex_tree.h>
+#include <gudhi/distance_functions.h>
+#include <gudhi/Persistent_cohomology.h>
+#include <gudhi/Rips_complex.h>
+
+#include <vector>
+
+struct Persistence_params {
+ int p;
+ double threshold;
+ int max_dim;
+ double min_pers;
+
+ Persistence_params(int p_, double th_, int max_dim_ = 10, double min_pers_ = 0)
+ : p(p_), threshold(th_), max_dim(max_dim_), min_pers(min_pers_) { }
+};
+
+/**
+ * Show persistence into output stream
+ */
+template<typename SkBlComplex> class Persistence_compute {
+ public:
+ typedef typename SkBlComplex::Vertex_handle Vertex_handle;
+ typedef typename SkBlComplex::Edge_handle Edge_handle;
+
+ /**
+ * @brief Compute persistence
+ * parameters :
+ * unsigned dim_max
+ * double threshold
+ * int p for coefficient Z_p
+ */
+ Persistence_compute(SkBlComplex& complex, std::ostream& stream, const Persistence_params& params) {
+ // for now everything is copied, todo boost adapt iterators to points of SkBlComplex instead of copying to an
+ // initial vector
+ typedef std::vector<double> Point_t;
+ std::vector< Point_t > points;
+ points.reserve(complex.num_vertices());
+ for (auto v : complex.vertex_range()) {
+ const auto & pt = complex.point(v);
+ Point_t pt_to_add(pt.cartesian_begin(), pt.cartesian_end());
+ points.emplace_back(std::move(pt_to_add));
+ }
+
+ using Simplex_tree = Gudhi::Simplex_tree<>;
+ using Filtration_value = Simplex_tree::Filtration_value;
+ using Rips_complex = Gudhi::rips_complex::Rips_complex<Filtration_value>;
+ using Field_Zp = Gudhi::persistent_cohomology::Field_Zp;
+ using Persistent_cohomology = Gudhi::persistent_cohomology::Persistent_cohomology<Simplex_tree, Field_Zp>;
+
+ Rips_complex rips_complex(points, params.threshold, Euclidean_distance());
+
+ Simplex_tree st;
+ rips_complex.create_complex(st, params.max_dim);
+ Persistent_cohomology pcoh(st);
+ // initializes the coefficient field for homology
+ pcoh.init_coefficients(params.p);
+ // put params.min_pers
+ pcoh.compute_persistent_cohomology(params.min_pers);
+ stream << "persistence: \n";
+ stream << "p dimension birth death: \n";
+ pcoh.output_diagram(stream);
+ }
+};
+
+#endif // UTILS_PERSISTENCE_COMPUTE_H_
diff --git a/src/GudhUI/utils/Rips_builder.h b/src/GudhUI/utils/Rips_builder.h
new file mode 100644
index 00000000..aba1a8e4
--- /dev/null
+++ b/src/GudhUI/utils/Rips_builder.h
@@ -0,0 +1,56 @@
+/* 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 UTILS_RIPS_BUILDER_H_
+#define UTILS_RIPS_BUILDER_H_
+
+#include <boost/iterator/iterator_facade.hpp>
+
+#include <CGAL/Euclidean_distance.h>
+#include <CGAL/Orthogonal_k_neighbor_search.h>
+#include <CGAL/Search_traits_d.h>
+
+#include "utils/UI_utils.h"
+#include "model/Complex_typedefs.h"
+
+template<typename SkBlComplex> class Rips_builder {
+ private:
+ SkBlComplex& complex_;
+
+ public:
+ /**
+ * @brief Modify complex to be the Rips complex
+ * of its points with offset alpha.
+ */
+ Rips_builder(SkBlComplex& complex, double alpha) : complex_(complex) {
+ complex.keep_only_vertices();
+ if (alpha <= 0) return;
+ compute_edges(alpha);
+ }
+
+ private:
+ double squared_eucl_distance(const Point& p1, const Point& p2) const {
+ return Geometry_trait::Squared_distance_d()(p1, p2);
+ }
+
+ void compute_edges(double alpha) {
+ auto vertices = complex_.vertex_range();
+ for (auto p = vertices.begin(); p != vertices.end(); ++p) {
+ std::cout << *p << " ";
+ std::cout.flush();
+ for (auto q = p; ++q != vertices.end(); /**/)
+ if (squared_eucl_distance(complex_.point(*p), complex_.point(*q)) < 4 * alpha * alpha)
+ complex_.add_edge_without_blockers(*p, *q);
+ }
+ std::cout << std::endl;
+ }
+};
+
+#endif // UTILS_RIPS_BUILDER_H_
diff --git a/src/GudhUI/utils/UI_utils.h b/src/GudhUI/utils/UI_utils.h
new file mode 100644
index 00000000..571a3cf2
--- /dev/null
+++ b/src/GudhUI/utils/UI_utils.h
@@ -0,0 +1,32 @@
+/* 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 UTILS_UI_UTILS_H_
+#define UTILS_UI_UTILS_H_
+
+#define UIDBG_VERBOSE
+
+#ifdef UIDBG_VERBOSE
+#define UIDBG(a) std::cerr << "UIDBG: " << (a) << std::endl
+#define UIDBGMSG(a, b) std::cerr << "UIDBG: " << a << b << std::endl
+#define UIDBGVALUE(a) std::cerr << "UIDBG: " << #a << ": " << a << std::endl
+#define UIDBGCONT(a) std::cerr << "UIDBG: container " << #a << " -> "; for (auto x : a) std::cerr << x << ","; std::cerr << std::endl }
+#else
+// #define DBG(a) a
+// #define DBGMSG(a, b) b
+// #define DBGVALUE(a) a
+// #define DBGCONT(a) a
+#define UIDBG(a)
+#define UIDBGMSG(a, b)
+#define UIDBGVALUE(a)
+#define UIDBGCONT(a)
+#endif
+
+#endif // UTILS_UI_UTILS_H_
diff --git a/src/GudhUI/utils/Vertex_collapsor.h b/src/GudhUI/utils/Vertex_collapsor.h
new file mode 100644
index 00000000..030e4bb0
--- /dev/null
+++ b/src/GudhUI/utils/Vertex_collapsor.h
@@ -0,0 +1,76 @@
+/* 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 UTILS_VERTEX_COLLAPSOR_H_
+#define UTILS_VERTEX_COLLAPSOR_H_
+
+#include <list>
+
+#include "utils/Edge_contractor.h"
+#include "utils/Furthest_point_epsilon_net.h"
+#include "utils/UI_utils.h"
+
+/**
+ * Iteratively puts every vertex at the center of its neighbors
+ */
+template<typename SkBlComplex> class Vertex_collapsor {
+ private:
+ SkBlComplex& complex_;
+ size_t num_collapses_;
+
+ public:
+ typedef typename SkBlComplex::Vertex_handle Vertex_handle;
+ typedef typename SkBlComplex::Edge_handle Edge_handle;
+
+ /**
+ * @brief Modify complex to be the expansion of the k-nearest neighbor
+ * symetric graph.
+ */
+ Vertex_collapsor(SkBlComplex& complex, size_t num_collapses) :
+ complex_(complex), num_collapses_(num_collapses) {
+ // std::list<Vertex_handle> vertices;
+ // vertices.insert(vertices.begin(),complex_.vertex_range().begin(),complex_.vertex_range().end());
+ // UIDBG("Collapse vertices");
+ // collapse_vertices(vertices);
+
+ std::list<Vertex_handle> vertices;
+
+ UIDBG("Compute eps net");
+ Furthest_point_epsilon_net<Complex> eps_net(complex_);
+
+ for (auto vh : eps_net.net_filtration_)
+ vertices.push_back(vh.vertex_handle);
+
+ UIDBG("Collapse vertices");
+ collapse_vertices(vertices);
+ }
+
+ private:
+ void collapse_vertices(std::list<Vertex_handle>& vertices) {
+ while (!vertices.empty() && num_collapses_--) {
+ Vertex_handle current_vertex = vertices.front();
+ vertices.pop_front();
+ if (is_link_reducible(current_vertex))
+ complex_.remove_vertex(current_vertex);
+ }
+ }
+
+ bool is_link_reducible(Vertex_handle v) {
+ auto link = complex_.link(v);
+ if (link.empty()) return false;
+ if (link.is_cone()) return true;
+ if (link.num_connected_components() > 1) return false;
+ Edge_contractor<Complex> contractor(link, link.num_vertices() - 1);
+ (void)contractor;
+ return (link.num_vertices() == 1);
+ }
+};
+
+#endif // UTILS_VERTEX_COLLAPSOR_H_
generated by cgit v1.2.3 (git 2.39.1) at 2024-05-17 08:17:10 +0000