voronoi_builder.hpp 19 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517
  1. // Boost.Polygon library voronoi_builder.hpp header file
  2. // Copyright Andrii Sydorchuk 2010-2012.
  3. // Distributed under the Boost Software License, Version 1.0.
  4. // (See accompanying file LICENSE_1_0.txt or copy at
  5. // http://www.boost.org/LICENSE_1_0.txt)
  6. // See http://www.boost.org for updates, documentation, and revision history.
  7. #ifndef BOOST_POLYGON_VORONOI_BUILDER
  8. #define BOOST_POLYGON_VORONOI_BUILDER
  9. #include <algorithm>
  10. #include <map>
  11. #include <queue>
  12. #include <utility>
  13. #include <vector>
  14. #include "detail/voronoi_ctypes.hpp"
  15. #include "detail/voronoi_predicates.hpp"
  16. #include "detail/voronoi_structures.hpp"
  17. #include "voronoi_geometry_type.hpp"
  18. namespace boost {
  19. namespace polygon {
  20. // GENERAL INFO:
  21. // The sweepline algorithm implementation to compute Voronoi diagram of
  22. // points and non-intersecting segments (except endpoints).
  23. // Complexity - O(N*logN), memory usage - O(N), where N is the total number
  24. // of input geometries. Input geometries should have integer coordinate type.
  25. //
  26. // IMPLEMENTATION DETAILS:
  27. // Each input point creates one site event. Each input segment creates three
  28. // site events: two for its endpoints and one for the segment itself (this is
  29. // made to simplify output construction). All the site events are constructed
  30. // and sorted at the algorithm initialization step. Priority queue is used to
  31. // dynamically hold circle events. At each step of the algorithm execution the
  32. // leftmost event is retrieved by comparing the current site event and the
  33. // topmost element from the circle event queue. STL map (red-black tree)
  34. // container was chosen to hold state of the beach line. The keys of the map
  35. // correspond to the neighboring sites that form a bisector and values map to
  36. // the corresponding Voronoi edges in the output data structure.
  37. template <typename T,
  38. typename CTT = detail::voronoi_ctype_traits<T>,
  39. typename VP = detail::voronoi_predicates<CTT> >
  40. class voronoi_builder {
  41. public:
  42. typedef typename CTT::int_type int_type;
  43. typedef typename CTT::fpt_type fpt_type;
  44. voronoi_builder() : index_(0) {}
  45. // Each point creates a single site event.
  46. std::size_t insert_point(const int_type& x, const int_type& y) {
  47. site_events_.push_back(site_event_type(x, y));
  48. site_events_.back().initial_index(index_);
  49. site_events_.back().source_category(SOURCE_CATEGORY_SINGLE_POINT);
  50. return index_++;
  51. }
  52. // Each segment creates three site events that correspond to:
  53. // 1) the start point of the segment;
  54. // 2) the end point of the segment;
  55. // 3) the segment itself defined by its start point.
  56. std::size_t insert_segment(
  57. const int_type& x1, const int_type& y1,
  58. const int_type& x2, const int_type& y2) {
  59. // Set up start point site.
  60. point_type p1(x1, y1);
  61. site_events_.push_back(site_event_type(p1));
  62. site_events_.back().initial_index(index_);
  63. site_events_.back().source_category(SOURCE_CATEGORY_SEGMENT_START_POINT);
  64. // Set up end point site.
  65. point_type p2(x2, y2);
  66. site_events_.push_back(site_event_type(p2));
  67. site_events_.back().initial_index(index_);
  68. site_events_.back().source_category(SOURCE_CATEGORY_SEGMENT_END_POINT);
  69. // Set up segment site.
  70. if (point_comparison_(p1, p2)) {
  71. site_events_.push_back(site_event_type(p1, p2));
  72. site_events_.back().source_category(SOURCE_CATEGORY_INITIAL_SEGMENT);
  73. } else {
  74. site_events_.push_back(site_event_type(p2, p1));
  75. site_events_.back().source_category(SOURCE_CATEGORY_REVERSE_SEGMENT);
  76. }
  77. site_events_.back().initial_index(index_);
  78. return index_++;
  79. }
  80. // Run sweepline algorithm and fill output data structure.
  81. template <typename OUTPUT>
  82. void construct(OUTPUT* output) {
  83. // Init structures.
  84. output->_reserve(site_events_.size());
  85. init_sites_queue();
  86. init_beach_line(output);
  87. // The algorithm stops when there are no events to process.
  88. event_comparison_predicate event_comparison;
  89. while (!circle_events_.empty() ||
  90. !(site_event_iterator_ == site_events_.end())) {
  91. if (circle_events_.empty()) {
  92. process_site_event(output);
  93. } else if (site_event_iterator_ == site_events_.end()) {
  94. process_circle_event(output);
  95. } else {
  96. if (event_comparison(*site_event_iterator_,
  97. circle_events_.top().first)) {
  98. process_site_event(output);
  99. } else {
  100. process_circle_event(output);
  101. }
  102. }
  103. while (!circle_events_.empty() &&
  104. !circle_events_.top().first.is_active()) {
  105. circle_events_.pop();
  106. }
  107. }
  108. beach_line_.clear();
  109. // Finish construction.
  110. output->_build();
  111. }
  112. void clear() {
  113. index_ = 0;
  114. site_events_.clear();
  115. }
  116. private:
  117. typedef detail::point_2d<int_type> point_type;
  118. typedef detail::site_event<int_type> site_event_type;
  119. typedef typename std::vector<site_event_type>::const_iterator
  120. site_event_iterator_type;
  121. typedef detail::circle_event<fpt_type> circle_event_type;
  122. typedef typename VP::template point_comparison_predicate<point_type>
  123. point_comparison_predicate;
  124. typedef typename VP::
  125. template event_comparison_predicate<site_event_type, circle_event_type>
  126. event_comparison_predicate;
  127. typedef typename VP::
  128. template circle_formation_predicate<site_event_type, circle_event_type>
  129. circle_formation_predicate_type;
  130. typedef void edge_type;
  131. typedef detail::beach_line_node_key<site_event_type> key_type;
  132. typedef detail::beach_line_node_data<edge_type, circle_event_type>
  133. value_type;
  134. typedef typename VP::template node_comparison_predicate<key_type>
  135. node_comparer_type;
  136. typedef std::map< key_type, value_type, node_comparer_type > beach_line_type;
  137. typedef typename beach_line_type::iterator beach_line_iterator;
  138. typedef std::pair<circle_event_type, beach_line_iterator> event_type;
  139. typedef struct {
  140. bool operator()(const event_type& lhs, const event_type& rhs) const {
  141. return predicate(rhs.first, lhs.first);
  142. }
  143. event_comparison_predicate predicate;
  144. } event_comparison_type;
  145. typedef detail::ordered_queue<event_type, event_comparison_type>
  146. circle_event_queue_type;
  147. typedef std::pair<point_type, beach_line_iterator> end_point_type;
  148. void init_sites_queue() {
  149. // Sort site events.
  150. std::sort(site_events_.begin(), site_events_.end(),
  151. event_comparison_predicate());
  152. // Remove duplicates.
  153. site_events_.erase(std::unique(
  154. site_events_.begin(), site_events_.end()), site_events_.end());
  155. // Index sites.
  156. for (std::size_t cur = 0; cur < site_events_.size(); ++cur) {
  157. site_events_[cur].sorted_index(cur);
  158. }
  159. // Init site iterator.
  160. site_event_iterator_ = site_events_.begin();
  161. }
  162. template <typename OUTPUT>
  163. void init_beach_line(OUTPUT* output) {
  164. if (site_events_.empty())
  165. return;
  166. if (site_events_.size() == 1) {
  167. // Handle single site event case.
  168. output->_process_single_site(site_events_[0]);
  169. ++site_event_iterator_;
  170. } else {
  171. int skip = 0;
  172. while (site_event_iterator_ != site_events_.end() &&
  173. VP::is_vertical(site_event_iterator_->point0(),
  174. site_events_.begin()->point0()) &&
  175. VP::is_vertical(*site_event_iterator_)) {
  176. ++site_event_iterator_;
  177. ++skip;
  178. }
  179. if (skip == 1) {
  180. // Init beach line with the first two sites.
  181. init_beach_line_default(output);
  182. } else {
  183. // Init beach line with collinear vertical sites.
  184. init_beach_line_collinear_sites(output);
  185. }
  186. }
  187. }
  188. // Init beach line with the two first sites.
  189. // The first site is always a point.
  190. template <typename OUTPUT>
  191. void init_beach_line_default(OUTPUT* output) {
  192. // Get the first and the second site event.
  193. site_event_iterator_type it_first = site_events_.begin();
  194. site_event_iterator_type it_second = site_events_.begin();
  195. ++it_second;
  196. insert_new_arc(
  197. *it_first, *it_first, *it_second, beach_line_.end(), output);
  198. // The second site was already processed. Move the iterator.
  199. ++site_event_iterator_;
  200. }
  201. // Init beach line with collinear sites.
  202. template <typename OUTPUT>
  203. void init_beach_line_collinear_sites(OUTPUT* output) {
  204. site_event_iterator_type it_first = site_events_.begin();
  205. site_event_iterator_type it_second = site_events_.begin();
  206. ++it_second;
  207. while (it_second != site_event_iterator_) {
  208. // Create a new beach line node.
  209. key_type new_node(*it_first, *it_second);
  210. // Update the output.
  211. edge_type* edge = output->_insert_new_edge(*it_first, *it_second).first;
  212. // Insert a new bisector into the beach line.
  213. beach_line_.insert(beach_line_.end(),
  214. std::pair<key_type, value_type>(new_node, value_type(edge)));
  215. // Update iterators.
  216. ++it_first;
  217. ++it_second;
  218. }
  219. }
  220. void deactivate_circle_event(value_type* value) {
  221. if (value->circle_event()) {
  222. value->circle_event()->deactivate();
  223. value->circle_event(NULL);
  224. }
  225. }
  226. template <typename OUTPUT>
  227. void process_site_event(OUTPUT* output) {
  228. // Get next site event to process.
  229. site_event_type site_event = *site_event_iterator_;
  230. // Move site iterator.
  231. site_event_iterator_type last = site_event_iterator_ + 1;
  232. // If a new site is an end point of some segment,
  233. // remove temporary nodes from the beach line data structure.
  234. if (!site_event.is_segment()) {
  235. while (!end_points_.empty() &&
  236. end_points_.top().first == site_event.point0()) {
  237. beach_line_iterator b_it = end_points_.top().second;
  238. end_points_.pop();
  239. beach_line_.erase(b_it);
  240. }
  241. } else {
  242. while (last != site_events_.end() &&
  243. last->is_segment() && last->point0() == site_event.point0())
  244. ++last;
  245. }
  246. // Find the node in the binary search tree with left arc
  247. // lying above the new site point.
  248. key_type new_key(*site_event_iterator_);
  249. beach_line_iterator right_it = beach_line_.lower_bound(new_key);
  250. for (; site_event_iterator_ != last; ++site_event_iterator_) {
  251. site_event = *site_event_iterator_;
  252. beach_line_iterator left_it = right_it;
  253. // Do further processing depending on the above node position.
  254. // For any two neighboring nodes the second site of the first node
  255. // is the same as the first site of the second node.
  256. if (right_it == beach_line_.end()) {
  257. // The above arc corresponds to the second arc of the last node.
  258. // Move the iterator to the last node.
  259. --left_it;
  260. // Get the second site of the last node
  261. const site_event_type& site_arc = left_it->first.right_site();
  262. // Insert new nodes into the beach line. Update the output.
  263. right_it = insert_new_arc(
  264. site_arc, site_arc, site_event, right_it, output);
  265. // Add a candidate circle to the circle event queue.
  266. // There could be only one new circle event formed by
  267. // a new bisector and the one on the left.
  268. activate_circle_event(left_it->first.left_site(),
  269. left_it->first.right_site(),
  270. site_event, right_it);
  271. } else if (right_it == beach_line_.begin()) {
  272. // The above arc corresponds to the first site of the first node.
  273. const site_event_type& site_arc = right_it->first.left_site();
  274. // Insert new nodes into the beach line. Update the output.
  275. left_it = insert_new_arc(
  276. site_arc, site_arc, site_event, right_it, output);
  277. // If the site event is a segment, update its direction.
  278. if (site_event.is_segment()) {
  279. site_event.inverse();
  280. }
  281. // Add a candidate circle to the circle event queue.
  282. // There could be only one new circle event formed by
  283. // a new bisector and the one on the right.
  284. activate_circle_event(site_event, right_it->first.left_site(),
  285. right_it->first.right_site(), right_it);
  286. right_it = left_it;
  287. } else {
  288. // The above arc corresponds neither to the first,
  289. // nor to the last site in the beach line.
  290. const site_event_type& site_arc2 = right_it->first.left_site();
  291. const site_event_type& site3 = right_it->first.right_site();
  292. // Remove the candidate circle from the event queue.
  293. deactivate_circle_event(&right_it->second);
  294. --left_it;
  295. const site_event_type& site_arc1 = left_it->first.right_site();
  296. const site_event_type& site1 = left_it->first.left_site();
  297. // Insert new nodes into the beach line. Update the output.
  298. beach_line_iterator new_node_it =
  299. insert_new_arc(site_arc1, site_arc2, site_event, right_it, output);
  300. // Add candidate circles to the circle event queue.
  301. // There could be up to two circle events formed by
  302. // a new bisector and the one on the left or right.
  303. activate_circle_event(site1, site_arc1, site_event, new_node_it);
  304. // If the site event is a segment, update its direction.
  305. if (site_event.is_segment()) {
  306. site_event.inverse();
  307. }
  308. activate_circle_event(site_event, site_arc2, site3, right_it);
  309. right_it = new_node_it;
  310. }
  311. }
  312. }
  313. // In general case circle event is made of the three consecutive sites
  314. // that form two bisectors in the beach line data structure.
  315. // Let circle event sites be A, B, C, two bisectors that define
  316. // circle event are (A, B), (B, C). During circle event processing
  317. // we remove (A, B), (B, C) and insert (A, C). As beach line comparison
  318. // works correctly only if one of the nodes is a new one we remove
  319. // (B, C) bisector and change (A, B) bisector to the (A, C). That's
  320. // why we use const_cast there and take all the responsibility that
  321. // map data structure keeps correct ordering.
  322. template <typename OUTPUT>
  323. void process_circle_event(OUTPUT* output) {
  324. // Get the topmost circle event.
  325. const event_type& e = circle_events_.top();
  326. const circle_event_type& circle_event = e.first;
  327. beach_line_iterator it_first = e.second;
  328. beach_line_iterator it_last = it_first;
  329. // Get the C site.
  330. site_event_type site3 = it_first->first.right_site();
  331. // Get the half-edge corresponding to the second bisector - (B, C).
  332. edge_type* bisector2 = it_first->second.edge();
  333. // Get the half-edge corresponding to the first bisector - (A, B).
  334. --it_first;
  335. edge_type* bisector1 = it_first->second.edge();
  336. // Get the A site.
  337. site_event_type site1 = it_first->first.left_site();
  338. if (!site1.is_segment() && site3.is_segment() &&
  339. site3.point1() == site1.point0()) {
  340. site3.inverse();
  341. }
  342. // Change the (A, B) bisector node to the (A, C) bisector node.
  343. const_cast<key_type&>(it_first->first).right_site(site3);
  344. // Insert the new bisector into the beach line.
  345. it_first->second.edge(output->_insert_new_edge(
  346. site1, site3, circle_event, bisector1, bisector2).first);
  347. // Remove the (B, C) bisector node from the beach line.
  348. beach_line_.erase(it_last);
  349. it_last = it_first;
  350. // Pop the topmost circle event from the event queue.
  351. circle_events_.pop();
  352. // Check new triplets formed by the neighboring arcs
  353. // to the left for potential circle events.
  354. if (it_first != beach_line_.begin()) {
  355. deactivate_circle_event(&it_first->second);
  356. --it_first;
  357. const site_event_type& site_l1 = it_first->first.left_site();
  358. activate_circle_event(site_l1, site1, site3, it_last);
  359. }
  360. // Check the new triplet formed by the neighboring arcs
  361. // to the right for potential circle events.
  362. ++it_last;
  363. if (it_last != beach_line_.end()) {
  364. deactivate_circle_event(&it_last->second);
  365. const site_event_type& site_r1 = it_last->first.right_site();
  366. activate_circle_event(site1, site3, site_r1, it_last);
  367. }
  368. }
  369. // Insert new nodes into the beach line. Update the output.
  370. template <typename OUTPUT>
  371. beach_line_iterator insert_new_arc(
  372. const site_event_type& site_arc1, const site_event_type &site_arc2,
  373. const site_event_type& site_event, beach_line_iterator position,
  374. OUTPUT* output) {
  375. // Create two new bisectors with opposite directions.
  376. key_type new_left_node(site_arc1, site_event);
  377. key_type new_right_node(site_event, site_arc2);
  378. // Set correct orientation for the first site of the second node.
  379. if (site_event.is_segment()) {
  380. new_right_node.left_site().inverse();
  381. }
  382. // Update the output.
  383. std::pair<edge_type*, edge_type*> edges =
  384. output->_insert_new_edge(site_arc2, site_event);
  385. position = beach_line_.insert(position,
  386. typename beach_line_type::value_type(
  387. new_right_node, value_type(edges.second)));
  388. if (site_event.is_segment()) {
  389. // Update the beach line with temporary bisector, that will
  390. // disappear after processing site event corresponding to the
  391. // second endpoint of the segment site.
  392. key_type new_node(site_event, site_event);
  393. new_node.right_site().inverse();
  394. position = beach_line_.insert(position,
  395. typename beach_line_type::value_type(new_node, value_type(NULL)));
  396. // Update the data structure that holds temporary bisectors.
  397. end_points_.push(std::make_pair(site_event.point1(), position));
  398. }
  399. position = beach_line_.insert(position,
  400. typename beach_line_type::value_type(
  401. new_left_node, value_type(edges.first)));
  402. return position;
  403. }
  404. // Add a new circle event to the event queue.
  405. // bisector_node corresponds to the (site2, site3) bisector.
  406. void activate_circle_event(const site_event_type& site1,
  407. const site_event_type& site2,
  408. const site_event_type& site3,
  409. beach_line_iterator bisector_node) {
  410. circle_event_type c_event;
  411. // Check if the three input sites create a circle event.
  412. if (circle_formation_predicate_(site1, site2, site3, c_event)) {
  413. // Add the new circle event to the circle events queue.
  414. // Update bisector's circle event iterator to point to the
  415. // new circle event in the circle event queue.
  416. event_type& e = circle_events_.push(
  417. std::pair<circle_event_type, beach_line_iterator>(
  418. c_event, bisector_node));
  419. bisector_node->second.circle_event(&e.first);
  420. }
  421. }
  422. private:
  423. point_comparison_predicate point_comparison_;
  424. struct end_point_comparison {
  425. bool operator() (const end_point_type& end1,
  426. const end_point_type& end2) const {
  427. return point_comparison(end2.first, end1.first);
  428. }
  429. point_comparison_predicate point_comparison;
  430. };
  431. std::vector<site_event_type> site_events_;
  432. site_event_iterator_type site_event_iterator_;
  433. std::priority_queue< end_point_type, std::vector<end_point_type>,
  434. end_point_comparison > end_points_;
  435. circle_event_queue_type circle_events_;
  436. beach_line_type beach_line_;
  437. circle_formation_predicate_type circle_formation_predicate_;
  438. std::size_t index_;
  439. // Disallow copy constructor and operator=
  440. voronoi_builder(const voronoi_builder&);
  441. void operator=(const voronoi_builder&);
  442. };
  443. typedef voronoi_builder<detail::int32> default_voronoi_builder;
  444. } // polygon
  445. } // boost
  446. #endif // BOOST_POLYGON_VORONOI_BUILDER