PGROUTING  2.5
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
alpha_driver.cpp File Reference
#include "./alpha_driver.h"
#include <CGAL/Simple_cartesian.h>
#include <CGAL/Filtered_kernel.h>
#include <CGAL/algorithm.h>
#include <CGAL/Polygon_2.h>
#include <CGAL/Delaunay_triangulation_2.h>
#include <CGAL/Triangulation_2.h>
#include <CGAL/Triangulation_hierarchy_vertex_base_2.h>
#include <CGAL/Triangulation_hierarchy_2.h>
#include <CGAL/Triangulation_face_base_2.h>
#include <CGAL/Triangulation_euclidean_traits_2.h>
#include <CGAL/Alpha_shape_2.h>
#include <CGAL/Alpha_shape_face_base_2.h>
#include <CGAL/Alpha_shape_vertex_base_2.h>
#include <vector>
#include <list>
#include <cmath>
#include <utility>
#include <algorithm>
#include <set>
#include "cpp_common/pgr_alloc.hpp"
Include dependency graph for alpha_driver.cpp:

Go to the source code of this file.

Typedefs

typedef
CGAL::Alpha_shape_face_base_2
< K, Tf
Af
 
typedef
Alpha_shape_2::Alpha_iterator 
Alpha_iterator
 
typedef CGAL::Alpha_shape_2< HtAlpha_shape_2
 
typedef
Alpha_shape_2::Alpha_shape_edges_iterator 
Alpha_shape_edges_iterator
 
typedef
CGAL::Triangulation_hierarchy_vertex_base_2
< Avb
Av
 
typedef
CGAL::Alpha_shape_vertex_base_2
< K
Avb
 
typedef double coord_type
 
typedef
CGAL::Delaunay_triangulation_2
< K, Tds
Dt
 
typedef
Alpha_shape_2::Face_circulator 
Face_circulator
 
typedef
CGAL::Triangulation_hierarchy_2
< Dt
Ht
 
typedef CGAL::Filtered_kernel< SCK
 
typedef K::Point_2 Point
 
typedef CGAL::Polygon_2< KPolygon_2
 
typedef CGAL::Simple_cartesian
< coord_type
SC
 
typedef K::Segment_2 Segment
 
typedef
CGAL::Triangulation_default_data_structure_2
< K, Av, Af
Tds
 
typedef
CGAL::Triangulation_face_base_2
< K
Tf
 
typedef K::Vector_2 Vector
 
typedef
Alpha_shape_2::Vertex_circulator 
Vertex_circulator
 

Functions

template<class OutputIterator >
void alpha_edges (const Alpha_shape_2 &A, OutputIterator out)
 
int alpha_shape (vertex_t *vertices, size_t count, double alpha, vertex_t **res, size_t *res_count, char **err_msg)
 
void find_next_edge (Segment s, std::vector< Segment > &segments, std::set< int > &unusedIndexes, std::vector< Polygon_2 > &rings)
 
double get_angle (Point p, Point q, Point r)
 

Variables

size_t prev_size = 0
 

Typedef Documentation

typedef CGAL::Alpha_shape_face_base_2<K, Tf> Af

Definition at line 77 of file alpha_driver.cpp.

typedef Alpha_shape_2::Alpha_iterator Alpha_iterator

Definition at line 87 of file alpha_driver.cpp.

typedef CGAL::Alpha_shape_2<Ht> Alpha_shape_2

Definition at line 82 of file alpha_driver.cpp.

typedef Alpha_shape_2::Alpha_shape_edges_iterator Alpha_shape_edges_iterator

Definition at line 88 of file alpha_driver.cpp.

typedef CGAL::Triangulation_hierarchy_vertex_base_2<Avb> Av

Definition at line 74 of file alpha_driver.cpp.

typedef CGAL::Alpha_shape_vertex_base_2<K> Avb

Definition at line 73 of file alpha_driver.cpp.

typedef double coord_type

Definition at line 64 of file alpha_driver.cpp.

typedef CGAL::Delaunay_triangulation_2<K, Tds> Dt

Definition at line 80 of file alpha_driver.cpp.

typedef Alpha_shape_2::Face_circulator Face_circulator

Definition at line 84 of file alpha_driver.cpp.

typedef CGAL::Triangulation_hierarchy_2<Dt> Ht

Definition at line 81 of file alpha_driver.cpp.

typedef CGAL::Filtered_kernel<SC> K

Definition at line 67 of file alpha_driver.cpp.

typedef K::Point_2 Point

Definition at line 68 of file alpha_driver.cpp.

typedef CGAL::Polygon_2<K> Polygon_2

Definition at line 71 of file alpha_driver.cpp.

typedef CGAL::Simple_cartesian<coord_type> SC

Definition at line 66 of file alpha_driver.cpp.

typedef K::Segment_2 Segment

Definition at line 69 of file alpha_driver.cpp.

typedef CGAL::Triangulation_default_data_structure_2<K, Av, Af> Tds

Definition at line 79 of file alpha_driver.cpp.

typedef CGAL::Triangulation_face_base_2<K> Tf

Definition at line 76 of file alpha_driver.cpp.

typedef K::Vector_2 Vector

Definition at line 70 of file alpha_driver.cpp.

typedef Alpha_shape_2::Vertex_circulator Vertex_circulator

Definition at line 85 of file alpha_driver.cpp.

Function Documentation

template<class OutputIterator >
void alpha_edges ( const Alpha_shape_2 A,
OutputIterator  out 
)

Definition at line 161 of file alpha_driver.cpp.

Referenced by alpha_shape().

162  {
163  for (Alpha_shape_edges_iterator it = A.alpha_shape_edges_begin();
164  it != A.alpha_shape_edges_end();
165  ++it) {
166  *out++ = A.segment(*it);
167  }
168 }
Alpha_shape_2::Alpha_shape_edges_iterator Alpha_shape_edges_iterator

Here is the caller graph for this function:

int alpha_shape ( vertex_t vertices,
size_t  count,
double  alpha,
vertex_t **  res,
size_t *  res_count,
char **  err_msg 
)

Definition at line 171 of file alpha_driver.cpp.

References alpha_edges(), find_next_edge(), and pgr_alloc().

Referenced by compute_alpha_shape().

172  {
173  try {
174  std::list<Point> points;
175  {
176  std::vector<Point> pv;
177 
178  for (std::size_t j = 0; j < count; ++j) {
179  Point p(vertices[j].x, vertices[j].y);
180  pv.push_back(p);
181  }
182 
183  std::sort(pv.begin(), pv.end(),
184  [](const Point &e1, const Point &e2)->bool {
185  return e2.y() < e1.y();
186  });
187  std::stable_sort(pv.begin(), pv.end(),
188  [](const Point &e1, const Point &e2)->bool {
189  return e2.x() < e1.x();
190  });
191  pv.erase(std::unique(pv.begin(), pv.end()), pv.end());
192  if (pv.size() != count && pv.size() < 3) {
193  *err_msg = strdup("After eliminating duplicated points, less than 3 points remain!!. Alpha shape calculation needs at least 3 vertices.");
194  return -1;
195  }
196  points.insert(points.begin(), pv.begin(), pv.end());
197  }
198 
199  Alpha_shape_2 A(points.begin(), points.end(),
200  coord_type(10000),
201  Alpha_shape_2::REGULARIZED);
202 
203  std::vector<Segment> segments;
204  // std::vector<Segment> result;
205 
206  // Alpha_shape_2::Alpha_shape_vertices_iterator vit;
207  // Alpha_shape_2::Vertex_handle vertex;
208  // Alpha_shape_2::Alpha_shape_edges_iterator eit;
209  // Alpha_shape_2::Edge edge;
210  // Alpha_shape_2::Face_iterator fit;
211  // Alpha_shape_2::Face_handle face;
212 
213  if (alpha <= 0.0) {
214  alpha = *A.find_optimal_alpha(1);
215  }
216  A.set_alpha(alpha);
217 
218  alpha_edges(A, std::back_inserter(segments));
219 
220  // Segment s = segments.at(0);
221  // find_next_edge(s, segments, result);
222  if (segments.empty()) {
223  *res = NULL;
224  *res_count = 0;
225  } else {
226  std::set<int> unusedIndexes;
227  for (unsigned int i = 0; i < segments.size(); i++) {
228  unusedIndexes.insert(i);
229  }
230 
231  std::vector<Polygon_2> rings;
232  Polygon_2 ring;
233  ring.push_back(segments.at(0).source());
234  rings.push_back(ring);
235  unusedIndexes.erase(0);
236  find_next_edge(segments.at(0), segments, unusedIndexes, rings);
237 
238  size_t result_count = 0;
239  for (unsigned int i = 0; i < rings.size(); i++) {
240  Polygon_2 ring = rings.at(i);
241  result_count += ring.size();
242  }
243  result_count += rings.size() - 1;
244  *res = pgr_alloc(result_count, (*res));
245  *res_count = result_count;
246 
247  int idx = 0;
248  for (unsigned int i = 0; i < rings.size(); i++) {
249  if (i > 0) {
250  (*res)[idx].x = DBL_MAX;
251  (*res)[idx].y = DBL_MAX;
252  idx++;
253  }
254  Polygon_2 ring = rings.at(i);
255  for (unsigned int j = 0; j < ring.size(); j++) {
256  Point point = ring.vertex(j);
257  (*res)[idx].x = point.x();
258  (*res)[idx].y = point.y();
259  idx++;
260  }
261  }
262  }
263  *err_msg = NULL;
264 
265  return EXIT_SUCCESS;
266  } catch ( ... ) {
267  *err_msg = strdup("Caught unknown exception!");
268  }
269  return -1;
270 }
K::Point_2 Point
void alpha_edges(const Alpha_shape_2 &A, OutputIterator out)
CGAL::Alpha_shape_2< Ht > Alpha_shape_2
double coord_type
CGAL::Polygon_2< K > Polygon_2
T * pgr_alloc(std::size_t size, T *ptr)
allocates memory
Definition: pgr_alloc.hpp:64
void find_next_edge(Segment s, std::vector< Segment > &segments, std::set< int > &unusedIndexes, std::vector< Polygon_2 > &rings)

Here is the call graph for this function:

Here is the caller graph for this function:

void find_next_edge ( Segment  s,
std::vector< Segment > &  segments,
std::set< int > &  unusedIndexes,
std::vector< Polygon_2 > &  rings 
)

Definition at line 106 of file alpha_driver.cpp.

References get_angle(), and prev_size.

Referenced by alpha_shape().

107  {
108  if (unusedIndexes.empty()
109  || prev_size == unusedIndexes.size()) {
110  return;
111  }
112 
113  prev_size = unusedIndexes.size();
114 
115  Point start = s.source();
116  Point end = s.target();
117  rings.back().push_back(end);
118 
119  std::vector<int> nextIndexes;
120  for (unsigned int i = 0; i < segments.size(); i++) {
121  if (unusedIndexes.find(i) != unusedIndexes.end()) {
122  Point source = segments.at(i).source();
123  if (source == end) {
124  nextIndexes.push_back(i);
125  }
126  }
127  }
128  if (nextIndexes.size() == 1) {
129  int i = nextIndexes.at(0);
130  unusedIndexes.erase(i);
131  find_next_edge(segments.at(i), segments, unusedIndexes, rings);
132  } else if (nextIndexes.size() > 1) {
133  std::vector< std::pair<double, int> > nextAngles;
134  for (unsigned int i = 0; i < nextIndexes.size(); i++) {
135  int j = nextIndexes.at(i);
136  Point target = segments.at(j).target();
137  double angle = get_angle(start, end, target);
138  nextAngles.push_back(std::pair<double, int>(angle, j));
139  }
140  std::sort(nextAngles.begin(), nextAngles.end());
141  int i = nextAngles.begin()->second;
142  unusedIndexes.erase(i);
143  find_next_edge(segments.at(i), segments, unusedIndexes, rings);
144  }
145 
146  if (!unusedIndexes.empty()) {
147  for (unsigned int i = 0; i < segments.size(); i++) {
148  if (unusedIndexes.find(i) != unusedIndexes.end()) {
149  Polygon_2 ring;
150  ring.push_back(segments.at(i).source());
151  rings.push_back(ring);
152  unusedIndexes.erase(i);
153  find_next_edge(segments.at(i), segments, unusedIndexes, rings);
154  }
155  }
156  }
157 }
size_t prev_size
K::Point_2 Point
double get_angle(Point p, Point q, Point r)
CGAL::Polygon_2< K > Polygon_2
void find_next_edge(Segment s, std::vector< Segment > &segments, std::set< int > &unusedIndexes, std::vector< Polygon_2 > &rings)

Here is the call graph for this function:

Here is the caller graph for this function:

double get_angle ( Point  p,
Point  q,
Point  r 
)

Definition at line 92 of file alpha_driver.cpp.

Referenced by find_next_edge().

92  {
93  double m_pi(3.14159265358979323846);
94  Vector v1(q, p);
95  Vector v2(q, r);
96  double cross = v1.x() * v2.y() - v1.y() * v2.x();
97  double dot = v1.x() * v2.x() + v1.y() * v2.y();
98  double angle = atan2(cross, dot);
99  if (angle < 0.0) {
100  angle += 2 * m_pi;
101  }
102  return angle;
103 }
K::Vector_2 Vector

Here is the caller graph for this function:

Variable Documentation

size_t prev_size = 0

Definition at line 105 of file alpha_driver.cpp.

Referenced by find_next_edge().