1 #pragma once
2 
3 #define _USE_MATH_DEFINES
4 #include <cmath>
5 
6 #include <algorithm>
7 #include <list>
8 #include <map>
9 #include <set>
10 #include <unordered_map>
11 #include <utility>
12 
13 #include <mapbox/geometry/wagyu/config.hpp>
14 #include <mapbox/geometry/wagyu/ring.hpp>
15 #include <mapbox/geometry/wagyu/ring_util.hpp>
16 
17 #ifdef DEBUG
18 #include <iostream>
19 #include <sstream>
20 #endif
21 
22 namespace mapbox {
23 namespace geometry {
24 namespace wagyu {
25 
26 template <typename T>
27 struct point_ptr_pair {
28     point_ptr<T> op1;
29     point_ptr<T> op2;
30 
point_ptr_pairmapbox::geometry::wagyu::point_ptr_pair31     constexpr point_ptr_pair(point_ptr<T> o1, point_ptr<T> o2) : op1(o1), op2(o2) {
32     }
33 
34     point_ptr_pair(point_ptr_pair<T> const& p) = default;
35 
point_ptr_pairmapbox::geometry::wagyu::point_ptr_pair36     point_ptr_pair(point_ptr_pair<T>&& p) : op1(std::move(p.op1)), op2(std::move(p.op2)) {
37     }
38 
operator =mapbox::geometry::wagyu::point_ptr_pair39     point_ptr_pair& operator=(point_ptr_pair<T>&& p) {
40         op1 = std::move(p.op1);
41         op2 = std::move(p.op2);
42         return *this;
43     }
44 };
45 
46 #ifdef DEBUG
47 
48 template <class charT, class traits, typename T>
49 inline std::basic_ostream<charT, traits>&
operator <<(std::basic_ostream<charT,traits> & out,const std::unordered_multimap<ring_ptr<T>,point_ptr_pair<T>> & dupe_ring)50 operator<<(std::basic_ostream<charT, traits>& out,
51            const std::unordered_multimap<ring_ptr<T>, point_ptr_pair<T>>& dupe_ring) {
52 
53     out << " BEGIN CONNECTIONS: " << std::endl;
54     for (auto& r : dupe_ring) {
55         out << "  Ring: ";
56         if (r.second.op1->ring) {
57             out << r.second.op1->ring->ring_index;
58         } else {
59             out << "---";
60         }
61         out << " to ";
62         if (r.second.op2->ring) {
63             out << r.second.op2->ring->ring_index;
64         } else {
65             out << "---";
66         }
67         out << "  ( at " << r.second.op1->x << ", " << r.second.op1->y << " )";
68         out << "  Ring1 ( ";
69         if (r.second.op1->ring) {
70             out << "area: " << r.second.op1->ring->area << " parent: ";
71             if (r.second.op1->ring->parent) {
72                 out << r.second.op1->ring->parent->ring_index;
73             } else {
74                 out << "---";
75             }
76         } else {
77             out << "---";
78         }
79         out << " )";
80         out << "  Ring2 ( ";
81         if (r.second.op2->ring) {
82             out << "area: " << r.second.op2->ring->area << " parent: ";
83             if (r.second.op2->ring->parent) {
84                 out << r.second.op2->ring->parent->ring_index;
85             } else {
86                 out << "---";
87             }
88         } else {
89             out << "---";
90         }
91         out << " )";
92         out << std::endl;
93     }
94     out << " END CONNECTIONS: " << std::endl;
95     return out;
96 }
97 
98 #endif
99 
100 template <typename T>
find_intersect_loop(std::unordered_multimap<ring_ptr<T>,point_ptr_pair<T>> & dupe_ring,std::list<std::pair<ring_ptr<T>,point_ptr_pair<T>>> & iList,ring_ptr<T> ring_parent,ring_ptr<T> ring_origin,ring_ptr<T> ring_search,std::set<ring_ptr<T>> & visited,point_ptr<T> orig_pt,point_ptr<T> prev_pt,ring_manager<T> & rings)101 bool find_intersect_loop(std::unordered_multimap<ring_ptr<T>, point_ptr_pair<T>>& dupe_ring,
102                          std::list<std::pair<ring_ptr<T>, point_ptr_pair<T>>>& iList,
103                          ring_ptr<T> ring_parent,
104                          ring_ptr<T> ring_origin,
105                          ring_ptr<T> ring_search,
106                          std::set<ring_ptr<T>>& visited,
107                          point_ptr<T> orig_pt,
108                          point_ptr<T> prev_pt,
109                          ring_manager<T>& rings) {
110     {
111         auto range = dupe_ring.equal_range(ring_search);
112         // Check for direct connection
113         for (auto& it = range.first; it != range.second;) {
114             ring_ptr<T> it_ring1 = it->second.op1->ring;
115             ring_ptr<T> it_ring2 = it->second.op2->ring;
116             if (!it_ring1 || !it_ring2 || it_ring1 != ring_search ||
117                 (!it_ring1->is_hole() && !it_ring2->is_hole())) {
118                 it = dupe_ring.erase(it);
119                 continue;
120             }
121             if (it_ring2 == ring_origin &&
122                 (ring_parent == it_ring2 || ring_parent == it_ring2->parent) &&
123                 *prev_pt != *it->second.op2 && *orig_pt != *it->second.op2) {
124                 iList.emplace_front(ring_search, it->second);
125                 return true;
126             }
127             ++it;
128         }
129     }
130     auto range = dupe_ring.equal_range(ring_search);
131     visited.insert(ring_search);
132     // Check for connection through chain of other intersections
133     for (auto& it = range.first;
134          it != range.second && it != dupe_ring.end() && it->first == ring_search; ++it) {
135         ring_ptr<T> it_ring = it->second.op2->ring;
136         if (visited.count(it_ring) > 0 || it_ring == nullptr ||
137             (ring_parent != it_ring && ring_parent != it_ring->parent) ||
138             value_is_zero(it_ring->area()) || *prev_pt == *it->second.op2) {
139             continue;
140         }
141         if (find_intersect_loop(dupe_ring, iList, ring_parent, ring_origin, it_ring, visited,
142                                 orig_pt, it->second.op2, rings)) {
143             iList.emplace_front(ring_search, it->second);
144             return true;
145         }
146     }
147     return false;
148 }
149 
150 template <typename T>
151 struct point_ptr_cmp {
operator ()mapbox::geometry::wagyu::point_ptr_cmp152     inline bool operator()(point_ptr<T> op1, point_ptr<T> op2) {
153         if (op1->y != op2->y) {
154             return (op1->y > op2->y);
155         } else if (op1->x != op2->x) {
156             return (op1->x < op2->x);
157         } else {
158             std::size_t depth_1 = ring_depth(op1->ring);
159             std::size_t depth_2 = ring_depth(op2->ring);
160             return depth_1 > depth_2;
161         }
162     }
163 };
164 
165 template <typename T>
correct_orientations(ring_manager<T> & manager)166 void correct_orientations(ring_manager<T>& manager) {
167     for (auto& r : manager.rings) {
168         if (!r.points) {
169             continue;
170         }
171         r.recalculate_stats();
172         if (r.size() < 3) {
173             remove_ring_and_points(&r, manager, false);
174             continue;
175         }
176         if (ring_is_hole(&r) != r.is_hole()) {
177             reverse_ring(r.points);
178             r.recalculate_stats();
179         }
180     }
181 }
182 
183 template <typename T>
sort_ring_points(ring_ptr<T> r)184 point_vector<T> sort_ring_points(ring_ptr<T> r) {
185     point_vector<T> sorted_points;
186     point_ptr<T> point_itr = r->points;
187     point_ptr<T> last_point = point_itr->prev;
188     while (point_itr != last_point) {
189         sorted_points.push_back(point_itr);
190         point_itr = point_itr->next;
191     }
192     sorted_points.push_back(last_point);
193     std::stable_sort(sorted_points.begin(), sorted_points.end(),
194                      [](point_ptr<T> const& pt1, point_ptr<T> const& pt2) {
195                          if (pt1->y != pt2->y) {
196                              return (pt1->y > pt2->y);
197                          }
198                          return (pt1->x < pt2->x);
199                      });
200     return sorted_points;
201 }
202 
203 template <typename T>
correct_self_intersection(point_ptr<T> pt1,point_ptr<T> pt2,ring_manager<T> & manager)204 ring_ptr<T> correct_self_intersection(point_ptr<T> pt1,
205                                       point_ptr<T> pt2,
206                                       ring_manager<T>& manager) {
207     if (pt1->ring != pt2->ring) {
208         return static_cast<ring_ptr<T>>(nullptr);
209     }
210 
211     ring_ptr<T> ring = pt1->ring;
212 
213     // split the polygon into two ...
214     point_ptr<T> pt3 = pt1->prev;
215     point_ptr<T> pt4 = pt2->prev;
216     pt1->prev = pt4;
217     pt4->next = pt1;
218     pt2->prev = pt3;
219     pt3->next = pt2;
220 
221     ring_ptr<T> new_ring = create_new_ring(manager);
222     std::size_t size_1 = 0;
223     std::size_t size_2 = 0;
224     mapbox::geometry::box<T> box1({ 0, 0 }, { 0, 0 });
225     mapbox::geometry::box<T> box2({ 0, 0 }, { 0, 0 });
226     double area_1 = area_from_point(pt1, size_1, box1);
227     double area_2 = area_from_point(pt2, size_2, box2);
228 
229     if (std::fabs(area_1) > std::fabs(area_2)) {
230         ring->points = pt1;
231         ring->set_stats(area_1, size_1, box1);
232         new_ring->points = pt2;
233         new_ring->set_stats(area_2, size_2, box2);
234     } else {
235         ring->points = pt2;
236         ring->set_stats(area_2, size_2, box2);
237         new_ring->points = pt1;
238         new_ring->set_stats(area_1, size_1, box1);
239     }
240     update_points_ring(new_ring);
241     return new_ring;
242 }
243 
244 template <typename T>
correct_repeated_points(ring_manager<T> & manager,ring_vector<T> & new_rings,point_vector_itr<T> const & begin,point_vector_itr<T> const & end)245 void correct_repeated_points(ring_manager<T>& manager,
246                              ring_vector<T>& new_rings,
247                              point_vector_itr<T> const& begin,
248                              point_vector_itr<T> const& end) {
249     for (auto itr1 = begin; itr1 != end; ++itr1) {
250         if ((*itr1)->ring == nullptr) {
251             continue;
252         }
253         for (auto itr2 = std::next(itr1); itr2 != end; ++itr2) {
254             if ((*itr2)->ring == nullptr) {
255                 continue;
256             }
257             ring_ptr<T> new_ring = correct_self_intersection(*itr1, *itr2, manager);
258             if (new_ring != nullptr) {
259                 new_rings.push_back(new_ring);
260             }
261         }
262     }
263 }
264 
265 template <typename T>
find_and_correct_repeated_points(ring_ptr<T> r,ring_manager<T> & manager,ring_vector<T> & new_rings)266 void find_and_correct_repeated_points(ring_ptr<T> r,
267                                       ring_manager<T>& manager,
268                                       ring_vector<T>& new_rings) {
269     auto sorted_points = sort_ring_points(r);
270     // Find sets of repeated points
271     std::size_t count = 0;
272     auto prev_itr = sorted_points.begin();
273     auto itr = std::next(prev_itr);
274     while (itr != sorted_points.end()) {
275         if (*(*prev_itr) == *(*(itr))) {
276             ++count;
277             ++prev_itr;
278             ++itr;
279             if (itr != sorted_points.end()) {
280                 continue;
281             } else {
282                 ++prev_itr;
283             }
284         } else {
285             ++prev_itr;
286             ++itr;
287         }
288         if (count == 0) {
289             continue;
290         }
291         auto first = prev_itr;
292         std::advance(first, -(static_cast<int>(count) + 1));
293         correct_repeated_points(manager, new_rings, first, prev_itr);
294         count = 0;
295     }
296 }
297 
298 template <typename T>
reassign_children_if_necessary(ring_ptr<T> new_ring,ring_ptr<T> sibling_ring,ring_manager<T> & manager,ring_vector<T> & new_rings)299 void reassign_children_if_necessary(ring_ptr<T> new_ring,
300                                     ring_ptr<T> sibling_ring,
301                                     ring_manager<T>& manager,
302                                     ring_vector<T>& new_rings) {
303     auto& children = sibling_ring == nullptr ? manager.children : sibling_ring->children;
304     for (auto c : children) {
305         if (c == nullptr) {
306             continue;
307         }
308         if (std::find(new_rings.begin(), new_rings.end(), c) != new_rings.end()) {
309             continue;
310         }
311         if (poly2_contains_poly1(c, new_ring)) {
312             reassign_as_child(c, new_ring, manager);
313         }
314     }
315 }
316 
317 template <typename T>
find_parent_in_tree(ring_ptr<T> r,ring_ptr<T> possible_parent,ring_manager<T> & manager)318 bool find_parent_in_tree(ring_ptr<T> r, ring_ptr<T> possible_parent, ring_manager<T>& manager) {
319     // Before starting this we are assuming that possible_parent
320     // and r have opposite signs of their areas
321 
322     // First we must search all grandchildren
323     for (auto c : possible_parent->children) {
324         if (c == nullptr) {
325             continue;
326         }
327         for (auto gc : c->children) {
328             if (gc == nullptr) {
329                 continue;
330             }
331             if (find_parent_in_tree(r, gc, manager)) {
332                 return true;
333             }
334         }
335     }
336 
337     if (poly2_contains_poly1(r, possible_parent)) {
338         reassign_as_child(r, possible_parent, manager);
339         return true;
340     }
341     return false;
342 }
343 
344 template <typename T>
assign_new_ring_parents(ring_manager<T> & manager,ring_ptr<T> original_ring,ring_vector<T> & new_rings)345 void assign_new_ring_parents(ring_manager<T>& manager,
346                              ring_ptr<T> original_ring,
347                              ring_vector<T>& new_rings) {
348 
349     // First lets remove any rings that have zero area
350     // or have no points
351     new_rings.erase(std::remove_if(new_rings.begin(), new_rings.end(),
352                                    [](ring_ptr<T> const& r) {
353                                        if (r->points == nullptr) {
354                                            return true;
355                                        }
356                                        return value_is_zero(r->area());
357                                    }),
358                     new_rings.end());
359 
360     if (new_rings.empty()) {
361         // No new rings created simply return;
362         return;
363     }
364 
365     // We should not have to re-assign the parent of the original ring
366     // because we always maintained the largest ring during splitting
367     // on repeated points.
368 
369     double original_ring_area = original_ring->area();
370     bool original_positive = original_ring_area > 0.0;
371 
372     // If there is only one new ring the logic is very simple and we
373     // do not have to check which ring contains, we only need to compare
374     // the areas of the original ring and that of the new ring.
375     if (new_rings.size() == 1) {
376         double new_ring_area = new_rings.front()->area();
377         bool new_positive = new_ring_area > 0.0;
378         if (original_positive == new_positive) {
379             // The rings should be siblings
380             assign_as_child(new_rings.front(), original_ring->parent, manager);
381             reassign_children_if_necessary(new_rings.front(), original_ring, manager, new_rings);
382         } else {
383             // The new ring is a child of original ring
384             // Check the
385             assign_as_child(new_rings.front(), original_ring, manager);
386             reassign_children_if_necessary(new_rings.front(), original_ring->parent, manager,
387                                            new_rings);
388         }
389         return;
390     }
391 
392     // Now we want to sort rings from the largest in absolute area to the smallest
393     // as we will assign the rings with the largest areas first
394     std::stable_sort(new_rings.begin(), new_rings.end(),
395                      [](ring_ptr<T> const& r1, ring_ptr<T> const& r2) {
396                          return std::fabs(r1->area()) > std::fabs(r2->area());
397                      });
398 
399     for (auto r_itr = new_rings.begin(); r_itr != new_rings.end(); ++r_itr) {
400         double new_ring_area = (*r_itr)->area();
401         bool new_positive = new_ring_area > 0.0;
402         bool same_orientation = new_positive == original_positive;
403         bool found = false;
404         // First lets check the trees of any new_rings that might have
405         // been assigned as siblings to the original ring.
406         for (auto s_itr = new_rings.begin(); s_itr != r_itr; ++s_itr) {
407             if ((*s_itr)->parent != original_ring->parent) {
408                 continue;
409             }
410             if (same_orientation) {
411                 for (auto s_child : (*s_itr)->children) {
412                     if (s_child == nullptr) {
413                         continue;
414                     }
415                     if (find_parent_in_tree(*r_itr, s_child, manager)) {
416                         reassign_children_if_necessary(*r_itr, original_ring, manager, new_rings);
417                         found = true;
418                         break;
419                     }
420                 }
421             } else {
422                 if (find_parent_in_tree(*r_itr, *s_itr, manager)) {
423                     reassign_children_if_necessary(*r_itr, original_ring->parent, manager,
424                                                    new_rings);
425                     found = true;
426                 }
427             }
428             if (found) {
429                 break;
430             }
431         }
432 
433         if (found) {
434             continue;
435         }
436 
437         // Next lets check the tree of the original_ring
438         if (same_orientation) {
439             for (auto o_child : original_ring->children) {
440                 if (o_child == nullptr) {
441                     continue;
442                 }
443                 if (find_parent_in_tree(*r_itr, o_child, manager)) {
444                     reassign_children_if_necessary(*r_itr, original_ring, manager, new_rings);
445                     found = true;
446                     break;
447                 }
448             }
449             if (!found) {
450                 // If we didn't find any parent and the same orientation
451                 // then it must be a sibling of the original ring
452                 assign_as_child(*r_itr, original_ring->parent, manager);
453                 reassign_children_if_necessary(*r_itr, original_ring, manager, new_rings);
454             }
455         } else {
456             if (find_parent_in_tree(*r_itr, original_ring, manager)) {
457                 reassign_children_if_necessary(*r_itr, original_ring->parent, manager, new_rings);
458             } else {
459                 throw std::runtime_error("Unable to find a proper parent ring");
460             }
461         }
462     }
463 }
464 
465 template <typename T>
correct_ring_self_intersections(ring_manager<T> & manager,ring_ptr<T> r,bool correct_tree)466 bool correct_ring_self_intersections(ring_manager<T>& manager, ring_ptr<T> r, bool correct_tree) {
467 
468     if (r->corrected || !r->points) {
469         return false;
470     }
471 
472     ring_vector<T> new_rings;
473 
474     find_and_correct_repeated_points(r, manager, new_rings);
475 
476     if (correct_tree) {
477         assign_new_ring_parents(manager, r, new_rings);
478     }
479 
480     r->corrected = true;
481     return true;
482 }
483 
484 template <typename T>
process_single_intersection(std::unordered_multimap<ring_ptr<T>,point_ptr_pair<T>> & connection_map,point_ptr<T> op_j,point_ptr<T> op_k,ring_manager<T> & manager)485 void process_single_intersection(
486     std::unordered_multimap<ring_ptr<T>, point_ptr_pair<T>>& connection_map,
487     point_ptr<T> op_j,
488     point_ptr<T> op_k,
489     ring_manager<T>& manager) {
490     ring_ptr<T> ring_j = op_j->ring;
491     ring_ptr<T> ring_k = op_k->ring;
492     if (ring_j == ring_k) {
493         return;
494     }
495 
496     if (!ring_j->is_hole() && !ring_k->is_hole()) {
497         // Both are not holes, return nothing to do.
498         return;
499     }
500 
501     ring_ptr<T> ring_origin;
502     ring_ptr<T> ring_search;
503     ring_ptr<T> ring_parent;
504     point_ptr<T> op_origin_1;
505     point_ptr<T> op_origin_2;
506     if (!ring_j->is_hole()) {
507         ring_origin = ring_j;
508         ring_parent = ring_origin;
509         ring_search = ring_k;
510         op_origin_1 = op_j;
511         op_origin_2 = op_k;
512     } else if (!ring_k->is_hole()) {
513         ring_origin = ring_k;
514         ring_parent = ring_origin;
515         ring_search = ring_j;
516         op_origin_1 = op_k;
517         op_origin_2 = op_j;
518 
519     } else {
520         // both are holes
521         // Order doesn't matter
522         ring_origin = ring_j;
523         ring_parent = ring_origin->parent;
524         ring_search = ring_k;
525         op_origin_1 = op_j;
526         op_origin_2 = op_k;
527     }
528     if (ring_parent != ring_search->parent) {
529         // The two holes do not have the same parent, do not add them
530         // simply return!
531         return;
532     }
533     bool found = false;
534     std::list<std::pair<ring_ptr<T>, point_ptr_pair<T>>> iList;
535     {
536         auto range = connection_map.equal_range(ring_search);
537         // Check for direct connection
538         for (auto& it = range.first; it != range.second;) {
539             if (!it->second.op1->ring) {
540                 it = connection_map.erase(it);
541                 continue;
542             }
543             if (!it->second.op2->ring) {
544                 it = connection_map.erase(it);
545                 continue;
546             }
547             ring_ptr<T> it_ring2 = it->second.op2->ring;
548             if (it_ring2 == ring_origin) {
549                 found = true;
550                 if (*op_origin_1 != *(it->second.op2)) {
551                     iList.emplace_back(ring_search, it->second);
552                     break;
553                 }
554             }
555             ++it;
556         }
557     }
558     if (iList.empty()) {
559         auto range = connection_map.equal_range(ring_search);
560         std::set<ring_ptr<T>> visited;
561         visited.insert(ring_search);
562         // Check for connection through chain of other intersections
563         for (auto& it = range.first;
564              it != range.second && it != connection_map.end() && it->first == ring_search; ++it) {
565             ring_ptr<T> it_ring = it->second.op2->ring;
566             if (it_ring != ring_search && *op_origin_2 != *it->second.op2 && it_ring != nullptr &&
567                 (ring_parent == it_ring || ring_parent == it_ring->parent) &&
568                 !value_is_zero(it_ring->area()) &&
569                 find_intersect_loop(connection_map, iList, ring_parent, ring_origin, it_ring,
570                                     visited, op_origin_2, it->second.op2, manager)) {
571                 found = true;
572                 iList.emplace_front(ring_search, it->second);
573                 break;
574             }
575         }
576     }
577     if (!found) {
578         point_ptr_pair<T> intPt_origin = { op_origin_1, op_origin_2 };
579         point_ptr_pair<T> intPt_search = { op_origin_2, op_origin_1 };
580         connection_map.emplace(ring_origin, std::move(intPt_origin));
581         connection_map.emplace(ring_search, std::move(intPt_search));
582         return;
583     }
584 
585     if (iList.empty()) {
586         // The situation where both origin and search are holes might have a missing
587         // search condition, we must check if a new pair must be added.
588         bool missing = true;
589         auto rng = connection_map.equal_range(ring_origin);
590         // Check for direct connection
591         for (auto& it = rng.first; it != rng.second; ++it) {
592             ring_ptr<T> it_ring2 = it->second.op2->ring;
593             if (it_ring2 == ring_search) {
594                 missing = false;
595             }
596         }
597         if (missing) {
598             point_ptr_pair<T> intPt_origin = { op_origin_1, op_origin_2 };
599             connection_map.emplace(ring_origin, std::move(intPt_origin));
600         }
601         return;
602     }
603     if (ring_origin->is_hole()) {
604         for (auto& iRing : iList) {
605             ring_ptr<T> ring_itr = iRing.first;
606             if (!ring_itr->is_hole()) {
607                 // Make the hole the origin!
608                 point_ptr<T> op1 = op_origin_1;
609                 op_origin_1 = iRing.second.op1;
610                 iRing.second.op1 = op1;
611                 point_ptr<T> op2 = op_origin_2;
612                 op_origin_2 = iRing.second.op2;
613                 iRing.second.op2 = op2;
614                 iRing.first = ring_origin;
615                 ring_origin = ring_itr;
616                 ring_parent = ring_origin;
617                 break;
618             }
619         }
620     }
621     bool origin_is_hole = ring_origin->is_hole();
622 
623     // Switch
624     point_ptr<T> op_origin_1_next = op_origin_1->next;
625     point_ptr<T> op_origin_2_next = op_origin_2->next;
626     op_origin_1->next = op_origin_2_next;
627     op_origin_2->next = op_origin_1_next;
628     op_origin_1_next->prev = op_origin_2;
629     op_origin_2_next->prev = op_origin_1;
630 
631     for (auto& iRing : iList) {
632         point_ptr<T> op_search_1 = iRing.second.op1;
633         point_ptr<T> op_search_2 = iRing.second.op2;
634         point_ptr<T> op_search_1_next = op_search_1->next;
635         point_ptr<T> op_search_2_next = op_search_2->next;
636         op_search_1->next = op_search_2_next;
637         op_search_2->next = op_search_1_next;
638         op_search_1_next->prev = op_search_2;
639         op_search_2_next->prev = op_search_1;
640     }
641 
642     ring_ptr<T> ring_new = create_new_ring(manager);
643     ring_origin->corrected = false;
644     std::size_t size_1 = 0;
645     std::size_t size_2 = 0;
646     mapbox::geometry::box<T> box1({ 0, 0 }, { 0, 0 });
647     mapbox::geometry::box<T> box2({ 0, 0 }, { 0, 0 });
648     double area_1 = area_from_point(op_origin_1, size_1, box1);
649     double area_2 = area_from_point(op_origin_2, size_2, box2);
650     if (origin_is_hole && ((area_1 < 0.0))) {
651         ring_origin->points = op_origin_1;
652         ring_origin->set_stats(area_1, size_1, box1);
653         ring_new->points = op_origin_2;
654         ring_new->set_stats(area_2, size_2, box2);
655     } else {
656         ring_origin->points = op_origin_2;
657         ring_origin->set_stats(area_2, size_2, box2);
658         ring_new->points = op_origin_1;
659         ring_new->set_stats(area_1, size_1, box1);
660     }
661 
662     update_points_ring(ring_origin);
663     update_points_ring(ring_new);
664 
665     ring_origin->bottom_point = nullptr;
666 
667     for (auto& iRing : iList) {
668         ring_ptr<T> ring_itr = iRing.first;
669         ring_itr->bottom_point = nullptr;
670         if (origin_is_hole) {
671             ring1_replaces_ring2(ring_origin, ring_itr, manager);
672         } else {
673             ring1_replaces_ring2(ring_origin->parent, ring_itr, manager);
674         }
675     }
676     if (origin_is_hole) {
677         assign_as_child(ring_new, ring_origin, manager);
678         // The parent ring in this situation might need to give up children
679         // to the new ring.
680         for (auto c : ring_parent->children) {
681             if (c == nullptr) {
682                 continue;
683             }
684             if (poly2_contains_poly1(c, ring_new)) {
685                 reassign_as_child(c, ring_new, manager);
686             }
687         }
688     } else {
689         // The new ring and the origin ring need to be siblings
690         // however some children ring from the ring origin might
691         // need to be re-assigned to the new ring
692         assign_as_sibling(ring_new, ring_origin, manager);
693         for (auto c : ring_origin->children) {
694             if (c == nullptr) {
695                 continue;
696             }
697             if (poly2_contains_poly1(c, ring_new)) {
698                 reassign_as_child(c, ring_new, manager);
699             }
700         }
701     }
702 
703     std::list<std::pair<ring_ptr<T>, point_ptr_pair<T>>> move_list;
704 
705     for (auto& iRing : iList) {
706         auto range_itr = connection_map.equal_range(iRing.first);
707         if (range_itr.first != range_itr.second) {
708             for (auto& it = range_itr.first; it != range_itr.second; ++it) {
709                 ring_ptr<T> it_ring = it->second.op1->ring;
710                 ring_ptr<T> it_ring2 = it->second.op2->ring;
711                 if (it_ring == nullptr || it_ring2 == nullptr || it_ring == it_ring2) {
712                     continue;
713                 }
714                 if (it_ring->is_hole() || it_ring2->is_hole()) {
715                     move_list.emplace_back(it_ring, it->second);
716                 }
717             }
718             connection_map.erase(iRing.first);
719         }
720     }
721 
722     auto range_itr = connection_map.equal_range(ring_origin);
723     for (auto& it = range_itr.first; it != range_itr.second;) {
724         ring_ptr<T> it_ring = it->second.op1->ring;
725         ring_ptr<T> it_ring2 = it->second.op2->ring;
726         if (it_ring == nullptr || it_ring2 == nullptr || it_ring == it_ring2) {
727             it = connection_map.erase(it);
728             continue;
729         }
730         if (it_ring != ring_origin) {
731             if (it_ring->is_hole() || it_ring2->is_hole()) {
732                 move_list.emplace_back(it_ring, it->second);
733             }
734             it = connection_map.erase(it);
735         } else {
736             if (it_ring->is_hole() || it_ring2->is_hole()) {
737                 ++it;
738             } else {
739                 it = connection_map.erase(it);
740             }
741         }
742     }
743 
744     if (!move_list.empty()) {
745         connection_map.insert(move_list.begin(), move_list.end());
746     }
747 
748     return;
749 }
750 
751 template <typename T>
correct_chained_repeats(ring_manager<T> & manager,std::unordered_multimap<ring_ptr<T>,point_ptr_pair<T>> & connection_map,point_vector_itr<T> const & begin,point_vector_itr<T> const & end)752 void correct_chained_repeats(
753     ring_manager<T>& manager,
754     std::unordered_multimap<ring_ptr<T>, point_ptr_pair<T>>& connection_map,
755     point_vector_itr<T> const& begin,
756     point_vector_itr<T> const& end) {
757     for (auto itr1 = begin; itr1 != end; ++itr1) {
758         if ((*itr1)->ring == nullptr) {
759             continue;
760         }
761         for (auto itr2 = std::next(itr1); itr2 != end; ++itr2) {
762             if ((*itr2)->ring == nullptr) {
763                 continue;
764             }
765             process_single_intersection(connection_map, *itr1, *itr2, manager);
766         }
767     }
768 }
769 
770 template <typename T>
correct_chained_rings(ring_manager<T> & manager)771 void correct_chained_rings(ring_manager<T>& manager) {
772 
773     if (manager.all_points.size() < 2) {
774         return;
775     }
776     // Setup connection map which is a map of rings and their
777     // connection point pairs with other rings.
778     std::unordered_multimap<ring_ptr<T>, point_ptr_pair<T>> connection_map;
779     connection_map.reserve(manager.rings.size());
780 
781     // Now lets find and process any points
782     // that overlap -- we should have solved
783     // all situations where these points
784     // would be self intersections of a ring with
785     // earlier processing so this should just be
786     // points where different rings are touching.
787     std::size_t count = 0;
788     auto prev_itr = manager.all_points.begin();
789     auto itr = std::next(prev_itr);
790     while (itr != manager.all_points.end()) {
791         if (*(*prev_itr) == *(*(itr))) {
792             ++count;
793             ++prev_itr;
794             ++itr;
795             if (itr != manager.all_points.end()) {
796                 continue;
797             } else {
798                 ++prev_itr;
799             }
800         } else {
801             ++prev_itr;
802             ++itr;
803         }
804         if (count == 0) {
805             continue;
806         }
807         auto first = prev_itr;
808         std::advance(first, -(static_cast<int>(count) + 1));
809         correct_chained_repeats(manager, connection_map, first, prev_itr);
810         count = 0;
811     }
812 }
813 
814 template <typename T>
sort_rings_largest_to_smallest(ring_manager<T> & manager)815 ring_vector<T> sort_rings_largest_to_smallest(ring_manager<T>& manager) {
816     ring_vector<T> sorted_rings;
817     sorted_rings.reserve(manager.rings.size());
818     for (auto& r : manager.rings) {
819         sorted_rings.push_back(&r);
820     }
821     std::stable_sort(sorted_rings.begin(), sorted_rings.end(),
822                      [](ring_ptr<T> const& r1, ring_ptr<T> const& r2) {
823                          if (!r1->points || !r2->points) {
824                              return r1->points != nullptr;
825                          }
826                          return std::fabs(r1->area()) > std::fabs(r2->area());
827                      });
828     return sorted_rings;
829 }
830 
831 template <typename T>
sort_rings_smallest_to_largest(ring_manager<T> & manager)832 ring_vector<T> sort_rings_smallest_to_largest(ring_manager<T>& manager) {
833     ring_vector<T> sorted_rings;
834     sorted_rings.reserve(manager.rings.size());
835     for (auto& r : manager.rings) {
836         sorted_rings.push_back(&r);
837     }
838     std::stable_sort(sorted_rings.begin(), sorted_rings.end(),
839                      [](ring_ptr<T> const& r1, ring_ptr<T> const& r2) {
840                          if (!r1->points || !r2->points) {
841                              return r1->points != nullptr;
842                          }
843                          return std::fabs(r1->area()) < std::fabs(r2->area());
844                      });
845     return sorted_rings;
846 }
847 
848 template <typename T>
849 struct collinear_path {
850     // Collinear edges are in opposite directions
851     // such that start_1 is at the same position
852     // of end_2 and vise versa. Start to end is
853     // always forward (next) on path.
854     point_ptr<T> start_1;
855     point_ptr<T> end_1;
856     point_ptr<T> start_2;
857     point_ptr<T> end_2;
858 };
859 
860 template <typename T>
861 struct collinear_result {
862     point_ptr<T> pt1;
863     point_ptr<T> pt2;
864 };
865 
866 template <typename T>
fix_collinear_path(collinear_path<T> & path)867 collinear_result<T> fix_collinear_path(collinear_path<T>& path) {
868     // Left and right are just the opposite ends of the
869     // collinear path, they may not be actually left
870     // and right of each other.
871     // The left end is start_1 and end_2
872     // The right end is start_2 and end_1
873 
874     // NOTE spike detection is checking that the
875     // pointers are the same values, not position!
876     // additionally duplicate points we will treat
877     // if they are a spike left.
878     bool spike_left = (path.start_1 == path.end_2);
879     bool spike_right = (path.start_2 == path.end_1);
880 
881     if (spike_left && spike_right) {
882         // If both ends are spikes we should simply
883         // delete all the points. (they should be in a loop)
884         point_ptr<T> itr = path.start_1;
885         while (itr != nullptr) {
886             itr->prev->next = nullptr;
887             itr->prev = nullptr;
888             itr->ring = nullptr;
889             itr = itr->next;
890         }
891         return { nullptr, nullptr };
892     } else if (spike_left) {
893         point_ptr<T> prev = path.start_2->prev;
894         point_ptr<T> itr = path.start_2;
895         while (itr != path.end_1) {
896             itr->prev->next = nullptr;
897             itr->prev = nullptr;
898             itr->ring = nullptr;
899             itr = itr->next;
900         }
901         prev->next = path.end_1;
902         path.end_1->prev = prev;
903         return { path.end_1, nullptr };
904     } else if (spike_right) {
905         point_ptr<T> prev = path.start_1->prev;
906         point_ptr<T> itr = path.start_1;
907         while (itr != path.end_2) {
908             itr->prev->next = nullptr;
909             itr->prev = nullptr;
910             itr->ring = nullptr;
911             itr = itr->next;
912         }
913         prev->next = path.end_2;
914         path.end_2->prev = prev;
915         return { path.end_2, nullptr };
916     } else {
917         point_ptr<T> prev_1 = path.start_1->prev;
918         point_ptr<T> prev_2 = path.start_2->prev;
919         point_ptr<T> itr = path.start_1;
920         do {
921             itr->prev->next = nullptr;
922             itr->prev = nullptr;
923             itr->ring = nullptr;
924             itr = itr->next;
925         } while (itr != path.end_1 && itr != nullptr);
926         itr = path.start_2;
927         do {
928             itr->prev->next = nullptr;
929             itr->prev = nullptr;
930             itr->ring = nullptr;
931             itr = itr->next;
932         } while (itr != path.end_2 && itr != nullptr);
933         if (path.start_1 == path.end_1 && path.start_2 == path.end_2) {
934             return { nullptr, nullptr };
935         } else if (path.start_1 == path.end_1) {
936             prev_2->next = path.end_2;
937             path.end_2->prev = prev_2;
938             return { path.end_2, nullptr };
939         } else if (path.start_2 == path.end_2) {
940             prev_1->next = path.end_1;
941             path.end_1->prev = prev_1;
942             return { path.end_1, nullptr };
943         } else {
944             prev_1->next = path.end_2;
945             path.end_2->prev = prev_1;
946             prev_2->next = path.end_1;
947             path.end_1->prev = prev_2;
948             return { path.end_1, path.end_2 };
949         }
950     }
951 }
952 
953 template <typename T>
find_start_and_end_of_collinear_edges(point_ptr<T> pt_a,point_ptr<T> pt_b)954 collinear_path<T> find_start_and_end_of_collinear_edges(point_ptr<T> pt_a, point_ptr<T> pt_b) {
955     // Search backward on A, forwards on B first
956     bool same_ring = pt_a->ring == pt_b->ring;
957     point_ptr<T> back = pt_a;
958     point_ptr<T> forward = pt_b;
959     bool first = true;
960     do {
961         while (*(back->prev) == *back && back != forward) {
962             back = back->prev;
963             if (back == pt_a) {
964                 break;
965             }
966         }
967         if (back == forward) {
968             back = back->prev;
969             forward = forward->next;
970             break;
971         }
972         while (*(forward->next) == *forward && back != forward) {
973             forward = forward->next;
974             if (forward == pt_b) {
975                 break;
976             }
977         }
978         if (!first && (back == pt_a || forward == pt_b)) {
979             break;
980         }
981         if (back == forward) {
982             back = back->prev;
983             forward = forward->next;
984             break;
985         }
986         back = back->prev;
987         forward = forward->next;
988         first = false;
989     } while (*back == *forward);
990     point_ptr<T> start_a = back->next;
991     // If there are repeated points at the diverge we want to select
992     // only the first of those repeated points.
993     while (!same_ring && *start_a == *start_a->next && start_a != pt_a) {
994         start_a = start_a->next;
995     }
996     point_ptr<T> end_b = forward->prev;
997     while (!same_ring && *end_b == *end_b->prev && end_b != pt_b) {
998         end_b = end_b->prev;
999     }
1000     // Search backward on B, forwards on A next
1001     back = pt_b;
1002     forward = pt_a;
1003     first = true;
1004     do {
1005         while (*(back->prev) == *back && back != forward) {
1006             back = back->prev;
1007             if (back == pt_b) {
1008                 break;
1009             }
1010         }
1011         if (back == forward) {
1012             back = back->prev;
1013             forward = forward->next;
1014             break;
1015         }
1016         while (*(forward->next) == *forward && back != forward) {
1017             forward = forward->next;
1018             if (forward == pt_a) {
1019                 break;
1020             }
1021         }
1022         if (!first && (back == pt_b || forward == pt_a)) {
1023             break;
1024         }
1025         if (back == forward || (!first && (back == end_b || forward == start_a))) {
1026             back = back->prev;
1027             forward = forward->next;
1028             break;
1029         }
1030         back = back->prev;
1031         forward = forward->next;
1032         first = false;
1033     } while (*back == *forward);
1034     point_ptr<T> start_b = back->next;
1035     while (!same_ring && *start_b == *start_b->next && start_b != pt_b) {
1036         start_b = start_b->next;
1037     }
1038     point_ptr<T> end_a = forward->prev;
1039     while (!same_ring && *end_a == *end_a->prev && end_a != pt_a) {
1040         end_a = end_a->prev;
1041     }
1042     return { start_a, end_a, start_b, end_b };
1043 }
1044 
1045 template <typename T>
has_collinear_edge(point_ptr<T> pt_a,point_ptr<T> pt_b)1046 bool has_collinear_edge(point_ptr<T> pt_a, point_ptr<T> pt_b) {
1047     // It is assumed pt_a and pt_b are at the same location.
1048     return (*pt_a->next == *pt_b->prev || *pt_b->next == *pt_a->prev);
1049 }
1050 
1051 template <typename T>
process_collinear_edges_same_ring(point_ptr<T> pt_a,point_ptr<T> pt_b,ring_manager<T> & manager)1052 void process_collinear_edges_same_ring(point_ptr<T> pt_a,
1053                                        point_ptr<T> pt_b,
1054                                        ring_manager<T>& manager) {
1055     ring_ptr<T> original_ring = pt_a->ring;
1056     // As they are the same ring that are forming a collinear edge
1057     // we should expect the creation of two different rings.
1058     auto path = find_start_and_end_of_collinear_edges(pt_a, pt_b);
1059     auto results = fix_collinear_path(path);
1060     if (results.pt1 == nullptr) {
1061         // If pt1 is a nullptr, both values
1062         // are nullptrs. This mean the ring was
1063         // removed during this processing.
1064         remove_ring(original_ring, manager, false);
1065     } else if (results.pt2 == nullptr) {
1066         // If a single point is only returned, we simply removed a spike.
1067         // In this case, we don't need to worry about parent or children
1068         // and we simply need to set the points and clear the area
1069         original_ring->points = results.pt1;
1070         original_ring->recalculate_stats();
1071     } else {
1072         // If we have two seperate points, the ring has split into
1073         // two different rings.
1074         ring_ptr<T> ring_new = create_new_ring(manager);
1075         ring_new->points = results.pt2;
1076         ring_new->recalculate_stats();
1077         update_points_ring(ring_new);
1078         original_ring->points = results.pt1;
1079         original_ring->recalculate_stats();
1080     }
1081 }
1082 
1083 template <typename T>
process_collinear_edges_different_rings(point_ptr<T> pt_a,point_ptr<T> pt_b,ring_manager<T> & manager)1084 void process_collinear_edges_different_rings(point_ptr<T> pt_a,
1085                                              point_ptr<T> pt_b,
1086                                              ring_manager<T>& manager) {
1087     ring_ptr<T> ring_a = pt_a->ring;
1088     ring_ptr<T> ring_b = pt_b->ring;
1089     bool ring_a_larger = std::fabs(ring_a->area()) > std::fabs(ring_b->area());
1090     auto path = find_start_and_end_of_collinear_edges(pt_a, pt_b);
1091     // This should result in two rings becoming one.
1092     auto results = fix_collinear_path(path);
1093     if (results.pt1 == nullptr) {
1094         remove_ring(ring_a, manager, false);
1095         remove_ring(ring_b, manager, false);
1096         return;
1097     }
1098     // Rings should merge into a single ring of the same orientation.
1099     // Therefore, we we will need to replace one ring with the other
1100     ring_ptr<T> merged_ring = ring_a_larger ? ring_a : ring_b;
1101     ring_ptr<T> deleted_ring = ring_a_larger ? ring_b : ring_a;
1102 
1103     merged_ring->points = results.pt1;
1104     update_points_ring(merged_ring);
1105     merged_ring->recalculate_stats();
1106     if (merged_ring->size() < 3) {
1107         remove_ring_and_points(merged_ring, manager, false);
1108     }
1109     remove_ring(deleted_ring, manager, false);
1110 }
1111 
1112 template <typename T>
remove_duplicate_points(point_ptr<T> pt_a,point_ptr<T> pt_b,ring_manager<T> & manager)1113 bool remove_duplicate_points(point_ptr<T> pt_a, point_ptr<T> pt_b, ring_manager<T>& manager) {
1114     if (pt_a->ring == pt_b->ring) {
1115         if (pt_a->next == pt_b) {
1116             pt_a->next = pt_b->next;
1117             pt_a->next->prev = pt_a;
1118             pt_b->next = nullptr;
1119             pt_b->prev = nullptr;
1120             pt_b->ring = nullptr;
1121             if (pt_a->ring->points == pt_b) {
1122                 pt_a->ring->points = pt_a;
1123             }
1124             return true;
1125         } else if (pt_b->next == pt_a) {
1126             pt_a->prev = pt_b->prev;
1127             pt_a->prev->next = pt_a;
1128             pt_b->next = nullptr;
1129             pt_b->prev = nullptr;
1130             pt_b->ring = nullptr;
1131             if (pt_a->ring->points == pt_b) {
1132                 pt_a->ring->points = pt_a;
1133             }
1134             return true;
1135         }
1136     }
1137     while (*pt_a->next == *pt_a && pt_a->next != pt_a) {
1138         point_ptr<T> remove = pt_a->next;
1139         pt_a->next = remove->next;
1140         pt_a->next->prev = pt_a;
1141         remove->next = nullptr;
1142         remove->prev = nullptr;
1143         remove->ring = nullptr;
1144         if (pt_a->ring->points == remove) {
1145             pt_a->ring->points = pt_a;
1146         }
1147     }
1148     while (*pt_a->prev == *pt_a && pt_a->prev != pt_a) {
1149         point_ptr<T> remove = pt_a->prev;
1150         pt_a->prev = remove->prev;
1151         pt_a->prev->next = pt_a;
1152         remove->next = nullptr;
1153         remove->prev = nullptr;
1154         remove->ring = nullptr;
1155         if (pt_a->ring->points == remove) {
1156             pt_a->ring->points = pt_a;
1157         }
1158     }
1159     if (pt_a->next == pt_a) {
1160         remove_ring_and_points(pt_a->ring, manager, false);
1161         return true;
1162     }
1163     if (pt_b->ring == nullptr) {
1164         return true;
1165     }
1166     while (*pt_b->next == *pt_b && pt_b->next != pt_b) {
1167         point_ptr<T> remove = pt_b->next;
1168         pt_b->next = remove->next;
1169         pt_b->next->prev = pt_b;
1170         remove->next = nullptr;
1171         remove->prev = nullptr;
1172         remove->ring = nullptr;
1173         if (pt_b->ring->points == remove) {
1174             pt_b->ring->points = pt_b;
1175         }
1176     }
1177     while (*pt_b->prev == *pt_b && pt_b->prev != pt_b) {
1178         point_ptr<T> remove = pt_b->prev;
1179         pt_b->prev = remove->prev;
1180         pt_b->prev->next = pt_b;
1181         remove->next = nullptr;
1182         remove->prev = nullptr;
1183         remove->ring = nullptr;
1184         if (pt_b->ring->points == remove) {
1185             pt_b->ring->points = pt_b;
1186         }
1187     }
1188     if (pt_b->next == pt_b) {
1189         remove_ring_and_points(pt_b->ring, manager, false);
1190         return true;
1191     }
1192     if (pt_a->ring == nullptr) {
1193         return true;
1194     }
1195     return false;
1196 }
1197 
1198 template <typename T>
process_collinear_edges(point_ptr<T> pt_a,point_ptr<T> pt_b,ring_manager<T> & manager)1199 bool process_collinear_edges(point_ptr<T> pt_a, point_ptr<T> pt_b, ring_manager<T>& manager) {
1200     // Neither point assigned to a ring (deleted points)
1201     if (!pt_a->ring || !pt_b->ring) {
1202         return false;
1203     }
1204 
1205     if (remove_duplicate_points(pt_a, pt_b, manager)) {
1206         return true;
1207     }
1208 
1209     if (!has_collinear_edge(pt_a, pt_b)) {
1210         if (pt_a->ring == pt_b->ring) {
1211             correct_self_intersection(pt_a, pt_b, manager);
1212             return true;
1213         }
1214         return false;
1215     }
1216 
1217     if (pt_a->ring == pt_b->ring) {
1218         process_collinear_edges_same_ring(pt_a, pt_b, manager);
1219     } else {
1220         process_collinear_edges_different_rings(pt_a, pt_b, manager);
1221     }
1222     return true;
1223 }
1224 
1225 template <typename T>
correct_collinear_repeats(ring_manager<T> & manager,point_vector_itr<T> const & begin,point_vector_itr<T> const & end)1226 void correct_collinear_repeats(ring_manager<T>& manager,
1227                                point_vector_itr<T> const& begin,
1228                                point_vector_itr<T> const& end) {
1229     for (auto itr1 = begin; itr1 != end; ++itr1) {
1230         if ((*itr1)->ring == nullptr) {
1231             continue;
1232         }
1233         for (auto itr2 = begin; itr2 != end;) {
1234             if ((*itr1)->ring == nullptr) {
1235                 break;
1236             }
1237             if ((*itr2)->ring == nullptr || *itr2 == *itr1) {
1238                 ++itr2;
1239                 continue;
1240             }
1241             if (process_collinear_edges(*itr1, *itr2, manager)) {
1242                 itr2 = begin;
1243             } else {
1244                 ++itr2;
1245             }
1246         }
1247     }
1248 }
1249 
1250 template <typename T>
correct_collinear_edges(ring_manager<T> & manager)1251 void correct_collinear_edges(ring_manager<T>& manager) {
1252 
1253     if (manager.all_points.size() < 2) {
1254         return;
1255     }
1256     std::size_t count = 0;
1257     auto prev_itr = manager.all_points.begin();
1258     auto itr = std::next(prev_itr);
1259     while (itr != manager.all_points.end()) {
1260         if (*(*prev_itr) == *(*(itr))) {
1261             ++count;
1262             ++prev_itr;
1263             ++itr;
1264             if (itr != manager.all_points.end()) {
1265                 continue;
1266             } else {
1267                 ++prev_itr;
1268             }
1269         } else {
1270             ++prev_itr;
1271             ++itr;
1272         }
1273         if (count == 0) {
1274             continue;
1275         }
1276         auto first = prev_itr;
1277         std::advance(first, -(static_cast<int>(count) + 1));
1278         correct_collinear_repeats(manager, first, prev_itr);
1279         count = 0;
1280     }
1281 }
1282 
1283 template <typename T>
correct_tree(ring_manager<T> & manager)1284 void correct_tree(ring_manager<T>& manager) {
1285 
1286     // It is possible that vatti resulted in some parent child
1287     // relationships that are not quite right, therefore, we
1288     // need to just rebuild the entire tree of rings
1289 
1290     // First lets process the rings in order of size from largest
1291     // area to smallest, we know right away that no smaller ring could ever
1292     // contain a larger ring so we can use this to our advantage
1293     // as we iterate over the rings.
1294     using rev_itr = typename ring_vector<T>::reverse_iterator;
1295     ring_vector<T> sorted_rings = sort_rings_largest_to_smallest(manager);
1296     for (auto itr = sorted_rings.begin(); itr != sorted_rings.end(); ++itr) {
1297         if ((*itr)->points == nullptr) {
1298             continue;
1299         }
1300         if ((*itr)->size() < 3 || value_is_zero((*itr)->area())) {
1301             remove_ring_and_points(*itr, manager, false);
1302             continue;
1303         }
1304         (*itr)->corrected = true;
1305         bool found = false;
1306         // Search in reverse from the current iterator back to the begining
1307         // to see if any of those rings might be its parent.
1308         for (auto r_itr = rev_itr(itr); r_itr != sorted_rings.rend(); ++r_itr) {
1309             // If orientations are not different, this can't be its parent.
1310             if ((*r_itr)->is_hole() == (*itr)->is_hole()) {
1311                 continue;
1312             }
1313             if (poly2_contains_poly1(*itr, *r_itr)) {
1314                 reassign_as_child(*itr, *r_itr, manager);
1315                 found = true;
1316                 break;
1317             }
1318         }
1319         if (!found) {
1320             if ((*itr)->is_hole()) {
1321                 throw std::runtime_error("Could not properly place hole to a parent.");
1322             } else {
1323                 // Assign to base of tree by passing nullptr
1324                 reassign_as_child(*itr, static_cast<ring_ptr<T>>(nullptr), manager);
1325             }
1326         }
1327     }
1328 }
1329 
1330 template <typename T>
correct_self_intersections(ring_manager<T> & manager,bool correct_tree)1331 bool correct_self_intersections(ring_manager<T>& manager, bool correct_tree) {
1332     bool fixed_intersections = false;
1333     auto sorted_rings = sort_rings_smallest_to_largest(manager);
1334     for (auto const& r : sorted_rings) {
1335         if (correct_ring_self_intersections(manager, r, correct_tree)) {
1336             fixed_intersections = true;
1337         }
1338     }
1339     return fixed_intersections;
1340 }
1341 
1342 template <typename T>
correct_topology(ring_manager<T> & manager)1343 void correct_topology(ring_manager<T>& manager) {
1344 
1345     // Sort all the points, this will be used for the locating of chained rings
1346     // and the collinear edges and only needs to be done once.
1347     std::stable_sort(manager.all_points.begin(), manager.all_points.end(), point_ptr_cmp<T>());
1348 
1349     // Initially the orientations of the rings
1350     // could be incorrect, we need to adjust them
1351     correct_orientations(manager);
1352 
1353     // We should only have to fix collinear edges once.
1354     // During this we also correct self intersections
1355     correct_collinear_edges(manager);
1356 
1357     correct_self_intersections(manager, false);
1358 
1359     correct_tree(manager);
1360 
1361     bool fixed_intersections = true;
1362     while (fixed_intersections) {
1363         correct_chained_rings(manager);
1364         fixed_intersections = correct_self_intersections(manager, true);
1365     }
1366 }
1367 }
1368 }
1369 }
1370