Critical_points.h
1 /* This file is part of the Gudhi Library - https://gudhi.inria.fr/ - which is released under MIT.
2  * See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details.
3  * Author(s): David Salinas
4  *
5  * Copyright (C) 2014 Inria
6  *
7  * Modification(s):
8  * - YYYY/MM Author: Description of the modification
9  */
10 
11 #ifndef UTILS_CRITICAL_POINTS_H_
12 #define UTILS_CRITICAL_POINTS_H_
13 
14 #include <deque>
15 #include <utility> // for pair<>
16 #include <algorithm> // for sort
17 
18 #include "utils/Edge_contractor.h"
19 
26 template<typename SkBlComplex> class Critical_points {
27  private:
28  SkBlComplex filled_complex_;
29  const SkBlComplex& input_complex_;
30  double max_length_;
31  std::ostream& stream_;
32 
33  public:
34  typedef typename SkBlComplex::Vertex_handle Vertex_handle;
35  typedef typename SkBlComplex::Edge_handle Edge_handle;
36  typedef typename std::pair<Vertex_handle, Vertex_handle> Edge;
37 
41  Critical_points(const SkBlComplex& input_complex, std::ostream& stream, double max_length) :
42  input_complex_(input_complex), max_length_(max_length), stream_(stream) {
43  std::deque<Edge> edges;
44  auto vertices = input_complex.vertex_range();
45  for (auto p = vertices.begin(); p != vertices.end(); ++p) {
46  filled_complex_.add_vertex(input_complex.point(*p));
47  for (auto q = p; ++q != vertices.end(); )
48  if (squared_eucl_distance(input_complex.point(*p), input_complex.point(*q)) < max_length_ * max_length_)
49  edges.emplace_back(*p, *q);
50  }
51 
52  std::sort(edges.begin(), edges.end(),
53  [&](Edge e1, Edge e2) {
54  return squared_edge_length(e1) < squared_edge_length(e2);
55  });
56 
57  anti_collapse_edges(edges);
58  }
59 
60  private:
61  double squared_eucl_distance(const Point& p1, const Point& p2) const {
62  return Geometry_trait::Squared_distance_d()(p1, p2);
63  }
64 
65  void anti_collapse_edges(const std::deque<Edge>& edges) {
66  unsigned pos = 0;
67  for (Edge e : edges) {
68  std::clog << "edge " << pos++ << "/" << edges.size() << "\n";
69  auto eh = filled_complex_.add_edge_without_blockers(e.first, e.second);
70  int is_contractible(is_link_reducible(eh));
71 
72  switch (is_contractible) {
73  case 0:
74  stream_ << "alpha=" << std::sqrt(squared_edge_length(e)) << " topological change" << std::endl;
75  break;
76  case 2:
77  stream_ << "alpha=" << std::sqrt(squared_edge_length(e)) << " maybe a topological change" << std::endl;
78  break;
79  default:
80  break;
81  }
82  }
83  }
84 
85  // 0 -> not
86  // 1 -> yes
87  // 2 -> maybe
88 
89  int is_link_reducible(Edge_handle e) {
90  auto link = filled_complex_.link(e);
91 
92  if (link.empty())
93  return 0;
94 
95  Edge_contractor<Complex> contractor(link, link.num_vertices() - 1);
96  (void)contractor;
97 
98  if (link.num_connected_components() > 1)
99  // one than more CC -> not contractible
100  return 0;
101 
102  if (link.num_vertices() == 1)
103  // reduced to one point -> contractible
104  return 1;
105  else
106  // we dont know
107  return 2;
108  }
109 
110  double squared_edge_length(Edge_handle e) const {
111  return squared_eucl_distance(input_complex_.point(input_complex_.first_vertex(e)),
112  input_complex_.point(input_complex_.second_vertex(e)));
113  }
114 
115  double squared_edge_length(Edge e) const {
116  return squared_eucl_distance(input_complex_.point(e.first), input_complex_.point(e.second));
117  }
118 };
119 
120 #endif // UTILS_CRITICAL_POINTS_H_
Definition: Critical_points.h:26
Critical_points(const SkBlComplex &input_complex, std::ostream &stream, double max_length)
check all pair of points with length smaller than max_length
Definition: Critical_points.h:41
Definition: Edge_contractor.h:24
GUDHIdev  Version 3.5.0  - C++ library for Topological Data Analysis (TDA) and Higher Dimensional Geometry Understanding.  - Copyright : MIT Generated on Sun May 1 2022 09:19:32 for GUDHIdev by Doxygen 1.9.1