Branch data Line data Source code
1 : : /*
2 : : * Poly2Tri Copyright (c) 2009-2018, Poly2Tri Contributors
3 : : * https://github.com/jhasse/poly2tri
4 : : *
5 : : * All rights reserved.
6 : : *
7 : : * Redistribution and use in source and binary forms, with or without modification,
8 : : * are permitted provided that the following conditions are met:
9 : : *
10 : : * * Redistributions of source code must retain the above copyright notice,
11 : : * this list of conditions and the following disclaimer.
12 : : * * Redistributions in binary form must reproduce the above copyright notice,
13 : : * this list of conditions and the following disclaimer in the documentation
14 : : * and/or other materials provided with the distribution.
15 : : * * Neither the name of Poly2Tri nor the names of its contributors may be
16 : : * used to endorse or promote products derived from this software without specific
17 : : * prior written permission.
18 : : *
19 : : * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
20 : : * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
21 : : * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
22 : : * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR
23 : : * CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
24 : : * EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
25 : : * PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
26 : : * PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
27 : : * LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
28 : : * NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
29 : : * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
30 : : */
31 : : #include "sweep.h"
32 : : #include "sweep_context.h"
33 : : #include "advancing_front.h"
34 : : #include "../common/utils.h"
35 : :
36 : : #include <cassert>
37 : : #include <stdexcept>
38 : :
39 : : namespace p2t {
40 : :
41 : : // Triangulate simple polygon with holes
42 : 13 : void Sweep::Triangulate(SweepContext& tcx)
43 : : {
44 : 13 : tcx.InitTriangulation();
45 : 13 : tcx.CreateAdvancingFront();
46 : : // Sweep points; build mesh
47 : 13 : SweepPoints(tcx);
48 : : // Clean up
49 : 13 : FinalizationPolygon(tcx);
50 : 13 : }
51 : :
52 : 13 : void Sweep::SweepPoints(SweepContext& tcx)
53 : : {
54 : 91 : for (size_t i = 1; i < tcx.point_count(); i++) {
55 : 78 : Point& point = *tcx.GetPoint(i);
56 : 78 : Node* node = &PointEvent(tcx, point);
57 : 169 : for (unsigned int i = 0; i < point.edge_list.size(); i++) {
58 : 91 : EdgeEvent(tcx, point.edge_list[i], node);
59 : 91 : }
60 : 78 : }
61 : 13 : }
62 : :
63 : 13 : void Sweep::FinalizationPolygon(SweepContext& tcx)
64 : : {
65 : : // Get an Internal triangle to start with
66 : 13 : Triangle* t = tcx.front()->head()->next->triangle;
67 : 13 : Point* p = tcx.front()->head()->next->point;
68 : 13 : while (!t->GetConstrainedEdgeCW(*p)) {
69 : 0 : t = t->NeighborCCW(*p);
70 : : }
71 : :
72 : : // Collect interior triangles constrained by edges
73 : 13 : tcx.MeshClean(*t);
74 : 13 : }
75 : :
76 : 78 : Node& Sweep::PointEvent(SweepContext& tcx, Point& point)
77 : : {
78 : 78 : Node& node = tcx.LocateNode(point);
79 : 78 : Node& new_node = NewFrontTriangle(tcx, point, node);
80 : :
81 : : // Only need to check +epsilon since point never have smaller
82 : : // x value than node due to how we fetch nodes from the front
83 : 78 : if (point.x <= node.point->x + EPSILON) {
84 : 39 : Fill(tcx, node);
85 : 39 : }
86 : :
87 : : //tcx.AddNode(new_node);
88 : :
89 : 78 : FillAdvancingFront(tcx, new_node);
90 : 78 : return new_node;
91 : : }
92 : :
93 : 91 : void Sweep::EdgeEvent(SweepContext& tcx, Edge* edge, Node* node)
94 : : {
95 : 91 : tcx.edge_event.constrained_edge = edge;
96 : 91 : tcx.edge_event.right = (edge->p->x > edge->q->x);
97 : :
98 : 91 : if (IsEdgeSideOfTriangle(*node->triangle, *edge->p, *edge->q)) {
99 : 39 : return;
100 : : }
101 : :
102 : : // For now we will do all needed filling
103 : : // TODO: integrate with flip process might give some better performance
104 : : // but for now this avoid the issue with cases that needs both flips and fills
105 : 52 : FillEdgeEvent(tcx, edge, node);
106 : 52 : EdgeEvent(tcx, *edge->p, *edge->q, node->triangle, *edge->q);
107 : 91 : }
108 : :
109 : 117 : void Sweep::EdgeEvent(SweepContext& tcx, Point& ep, Point& eq, Triangle* triangle, Point& point)
110 : : {
111 : 117 : if (IsEdgeSideOfTriangle(*triangle, ep, eq)) {
112 : 52 : return;
113 : : }
114 : :
115 : 65 : Point* p1 = triangle->PointCCW(point);
116 : 65 : Orientation o1 = Orient2d(eq, *p1, ep);
117 : 65 : if (o1 == COLLINEAR) {
118 : 0 : if( triangle->Contains(&eq, p1)) {
119 : 0 : triangle->MarkConstrainedEdge(&eq, p1 );
120 : : // We are modifying the constraint maybe it would be better to
121 : : // not change the given constraint and just keep a variable for the new constraint
122 : 0 : tcx.edge_event.constrained_edge->q = p1;
123 : 0 : triangle = &triangle->NeighborAcross(point);
124 : 0 : EdgeEvent( tcx, ep, *p1, triangle, *p1 );
125 : 0 : } else {
126 : 0 : std::runtime_error("EdgeEvent - collinear points not supported");
127 : 0 : assert(0);
128 : : }
129 : 0 : return;
130 : : }
131 : :
132 : 65 : Point* p2 = triangle->PointCW(point);
133 : 65 : Orientation o2 = Orient2d(eq, *p2, ep);
134 : 65 : if (o2 == COLLINEAR) {
135 : 0 : if( triangle->Contains(&eq, p2)) {
136 : 0 : triangle->MarkConstrainedEdge(&eq, p2 );
137 : : // We are modifying the constraint maybe it would be better to
138 : : // not change the given constraint and just keep a variable for the new constraint
139 : 0 : tcx.edge_event.constrained_edge->q = p2;
140 : 0 : triangle = &triangle->NeighborAcross(point);
141 : 0 : EdgeEvent( tcx, ep, *p2, triangle, *p2 );
142 : 0 : } else {
143 : 0 : std::runtime_error("EdgeEvent - collinear points not supported");
144 : 0 : assert(0);
145 : : }
146 : 0 : return;
147 : : }
148 : :
149 : 65 : if (o1 == o2) {
150 : : // Need to decide if we are rotating CW or CCW to get to a triangle
151 : : // that will cross edge
152 : 65 : if (o1 == CW) {
153 : 65 : triangle = triangle->NeighborCCW(point);
154 : 65 : } else{
155 : 0 : triangle = triangle->NeighborCW(point);
156 : : }
157 : 65 : EdgeEvent(tcx, ep, eq, triangle, point);
158 : 65 : } else {
159 : : // This triangle crosses constraint so lets flippin start!
160 : 0 : FlipEdgeEvent(tcx, ep, eq, triangle, point);
161 : : }
162 : 117 : }
163 : :
164 : 208 : bool Sweep::IsEdgeSideOfTriangle(Triangle& triangle, Point& ep, Point& eq)
165 : : {
166 : 208 : const int index = triangle.EdgeIndex(&ep, &eq);
167 : :
168 : 208 : if (index != -1) {
169 : 91 : triangle.MarkConstrainedEdge(index);
170 : 91 : Triangle* t = triangle.GetNeighbor(index);
171 : 91 : if (t) {
172 : 39 : t->MarkConstrainedEdge(&ep, &eq);
173 : 39 : }
174 : 91 : return true;
175 : : }
176 : 117 : return false;
177 : 208 : }
178 : :
179 : 78 : Node& Sweep::NewFrontTriangle(SweepContext& tcx, Point& point, Node& node)
180 : : {
181 : 78 : Triangle* triangle = new Triangle(point, *node.point, *node.next->point);
182 : :
183 : 78 : triangle->MarkNeighbor(*node.triangle);
184 : 78 : tcx.AddToMap(triangle);
185 : :
186 : 78 : Node* new_node = new Node(point);
187 : 78 : nodes_.push_back(new_node);
188 : :
189 : 78 : new_node->next = node.next;
190 : 78 : new_node->prev = &node;
191 : 78 : node.next->prev = new_node;
192 : 78 : node.next = new_node;
193 : :
194 : 78 : if (!Legalize(tcx, *triangle)) {
195 : 65 : tcx.MapTriangleToNodes(*triangle);
196 : 65 : }
197 : :
198 : 78 : return *new_node;
199 : 0 : }
200 : :
201 : 65 : void Sweep::Fill(SweepContext& tcx, Node& node)
202 : : {
203 : 65 : Triangle* triangle = new Triangle(*node.prev->point, *node.point, *node.next->point);
204 : :
205 : : // TODO: should copy the constrained_edge value from neighbor triangles
206 : : // for now constrained_edge values are copied during the legalize
207 : 65 : triangle->MarkNeighbor(*node.prev->triangle);
208 : 65 : triangle->MarkNeighbor(*node.triangle);
209 : :
210 : 65 : tcx.AddToMap(triangle);
211 : :
212 : : // Update the advancing front
213 : 65 : node.prev->next = node.next;
214 : 65 : node.next->prev = node.prev;
215 : :
216 : : // If it was legalized the triangle has already been mapped
217 : 65 : if (!Legalize(tcx, *triangle)) {
218 : 65 : tcx.MapTriangleToNodes(*triangle);
219 : 65 : }
220 : :
221 : 65 : }
222 : :
223 : 78 : void Sweep::FillAdvancingFront(SweepContext& tcx, Node& n)
224 : : {
225 : :
226 : : // Fill right holes
227 : 78 : Node* node = n.next;
228 : :
229 : 91 : while (node->next) {
230 : : // if HoleAngle exceeds 90 degrees then break.
231 : 65 : if (LargeHole_DontFill(node)) break;
232 : 13 : Fill(tcx, *node);
233 : 13 : node = node->next;
234 : : }
235 : :
236 : : // Fill left holes
237 : 78 : node = n.prev;
238 : :
239 : 78 : while (node->prev) {
240 : : // if HoleAngle exceeds 90 degrees then break.
241 : 65 : if (LargeHole_DontFill(node)) break;
242 : 0 : Fill(tcx, *node);
243 : 0 : node = node->prev;
244 : : }
245 : :
246 : : // Fill right basins
247 : 78 : if (n.next && n.next->next) {
248 : 52 : const double angle = BasinAngle(n);
249 : 52 : if (angle < PI_3div4) {
250 : 13 : FillBasin(tcx, n);
251 : 13 : }
252 : 52 : }
253 : 78 : }
254 : :
255 : : // True if HoleAngle exceeds 90 degrees.
256 : 130 : bool Sweep::LargeHole_DontFill(const Node* node) const {
257 : :
258 : 130 : const Node* nextNode = node->next;
259 : 130 : const Node* prevNode = node->prev;
260 : 130 : if (!AngleExceeds90Degrees(node->point, nextNode->point, prevNode->point))
261 : 13 : return false;
262 : :
263 : : // Check additional points on front.
264 : 117 : const Node* next2Node = nextNode->next;
265 : : // "..Plus.." because only want angles on same side as point being added.
266 : 117 : if ((next2Node != NULL) && !AngleExceedsPlus90DegreesOrIsNegative(node->point, next2Node->point, prevNode->point))
267 : 0 : return false;
268 : :
269 : 117 : const Node* prev2Node = prevNode->prev;
270 : : // "..Plus.." because only want angles on same side as point being added.
271 : 117 : if ((prev2Node != NULL) && !AngleExceedsPlus90DegreesOrIsNegative(node->point, nextNode->point, prev2Node->point))
272 : 0 : return false;
273 : :
274 : 117 : return true;
275 : 130 : }
276 : :
277 : 130 : bool Sweep::AngleExceeds90Degrees(const Point* origin, const Point* pa, const Point* pb) const {
278 : 130 : const double angle = Angle(origin, pa, pb);
279 : 130 : return ((angle > PI_div2) || (angle < -PI_div2));
280 : : }
281 : :
282 : 169 : bool Sweep::AngleExceedsPlus90DegreesOrIsNegative(const Point* origin, const Point* pa, const Point* pb) const {
283 : 169 : const double angle = Angle(origin, pa, pb);
284 : 169 : return (angle > PI_div2) || (angle < 0);
285 : : }
286 : :
287 : 299 : double Sweep::Angle(const Point* origin, const Point* pa, const Point* pb) const {
288 : : /* Complex plane
289 : : * ab = cosA +i*sinA
290 : : * ab = (ax + ay*i)(bx + by*i) = (ax*bx + ay*by) + i(ax*by-ay*bx)
291 : : * atan2(y,x) computes the principal value of the argument function
292 : : * applied to the complex number x+iy
293 : : * Where x = ax*bx + ay*by
294 : : * y = ax*by - ay*bx
295 : : */
296 : 299 : const double px = origin->x;
297 : 299 : const double py = origin->y;
298 : 299 : const double ax = pa->x- px;
299 : 299 : const double ay = pa->y - py;
300 : 299 : const double bx = pb->x - px;
301 : 299 : const double by = pb->y - py;
302 : 299 : const double x = ax * by - ay * bx;
303 : 299 : const double y = ax * bx + ay * by;
304 : 299 : return atan2(x, y);
305 : : }
306 : :
307 : 52 : double Sweep::BasinAngle(const Node& node) const
308 : : {
309 : 52 : const double ax = node.point->x - node.next->next->point->x;
310 : 52 : const double ay = node.point->y - node.next->next->point->y;
311 : 52 : return atan2(ay, ax);
312 : : }
313 : :
314 : 0 : double Sweep::HoleAngle(const Node& node) const
315 : : {
316 : : /* Complex plane
317 : : * ab = cosA +i*sinA
318 : : * ab = (ax + ay*i)(bx + by*i) = (ax*bx + ay*by) + i(ax*by-ay*bx)
319 : : * atan2(y,x) computes the principal value of the argument function
320 : : * applied to the complex number x+iy
321 : : * Where x = ax*bx + ay*by
322 : : * y = ax*by - ay*bx
323 : : */
324 : 0 : const double ax = node.next->point->x - node.point->x;
325 : 0 : const double ay = node.next->point->y - node.point->y;
326 : 0 : const double bx = node.prev->point->x - node.point->x;
327 : 0 : const double by = node.prev->point->y - node.point->y;
328 : 0 : return atan2(ax * by - ay * bx, ax * bx + ay * by);
329 : : }
330 : :
331 : 169 : bool Sweep::Legalize(SweepContext& tcx, Triangle& t)
332 : : {
333 : : // To legalize a triangle we start by finding if any of the three edges
334 : : // violate the Delaunay condition
335 : 637 : for (int i = 0; i < 3; i++) {
336 : 481 : if (t.delaunay_edge[i])
337 : 26 : continue;
338 : :
339 : 455 : Triangle* ot = t.GetNeighbor(i);
340 : :
341 : 455 : if (ot) {
342 : 221 : Point* p = t.GetPoint(i);
343 : 221 : Point* op = ot->OppositePoint(t, *p);
344 : 221 : int oi = ot->Index(op);
345 : :
346 : : // If this is a Constrained Edge or a Delaunay Edge(only during recursive legalization)
347 : : // then we should not try to legalize
348 : 221 : if (ot->constrained_edge[oi] || ot->delaunay_edge[oi]) {
349 : 52 : t.constrained_edge[i] = ot->constrained_edge[oi];
350 : 52 : continue;
351 : : }
352 : :
353 : 169 : bool inside = Incircle(*p, *t.PointCCW(*p), *t.PointCW(*p), *op);
354 : :
355 : 169 : if (inside) {
356 : : // Lets mark this shared edge as Delaunay
357 : 13 : t.delaunay_edge[i] = true;
358 : 13 : ot->delaunay_edge[oi] = true;
359 : :
360 : : // Lets rotate shared edge one vertex CW to legalize it
361 : 13 : RotateTrianglePair(t, *p, *ot, *op);
362 : :
363 : : // We now got one valid Delaunay Edge shared by two triangles
364 : : // This gives us 4 new edges to check for Delaunay
365 : :
366 : : // Make sure that triangle to node mapping is done only one time for a specific triangle
367 : 13 : bool not_legalized = !Legalize(tcx, t);
368 : 13 : if (not_legalized) {
369 : 13 : tcx.MapTriangleToNodes(t);
370 : 13 : }
371 : :
372 : 13 : not_legalized = !Legalize(tcx, *ot);
373 : 13 : if (not_legalized)
374 : 13 : tcx.MapTriangleToNodes(*ot);
375 : :
376 : : // Reset the Delaunay edges, since they only are valid Delaunay edges
377 : : // until we add a new triangle or point.
378 : : // XXX: need to think about this. Can these edges be tried after we
379 : : // return to previous recursive level?
380 : 13 : t.delaunay_edge[i] = false;
381 : 13 : ot->delaunay_edge[oi] = false;
382 : :
383 : : // If triangle have been legalized no need to check the other edges since
384 : : // the recursive legalization will handles those so we can end here.
385 : 13 : return true;
386 : : }
387 : 156 : }
388 : 390 : }
389 : 156 : return false;
390 : 169 : }
391 : :
392 : 169 : bool Sweep::Incircle(const Point& pa, const Point& pb, const Point& pc, const Point& pd) const
393 : : {
394 : 169 : const double adx = pa.x - pd.x;
395 : 169 : const double ady = pa.y - pd.y;
396 : 169 : const double bdx = pb.x - pd.x;
397 : 169 : const double bdy = pb.y - pd.y;
398 : :
399 : 169 : const double adxbdy = adx * bdy;
400 : 169 : const double bdxady = bdx * ady;
401 : 169 : const double oabd = adxbdy - bdxady;
402 : :
403 : 169 : if (oabd <= 0)
404 : 65 : return false;
405 : :
406 : 104 : const double cdx = pc.x - pd.x;
407 : 104 : const double cdy = pc.y - pd.y;
408 : :
409 : 104 : const double cdxady = cdx * ady;
410 : 104 : const double adxcdy = adx * cdy;
411 : 104 : const double ocad = cdxady - adxcdy;
412 : :
413 : 104 : if (ocad <= 0)
414 : 52 : return false;
415 : :
416 : 52 : const double bdxcdy = bdx * cdy;
417 : 52 : const double cdxbdy = cdx * bdy;
418 : :
419 : 52 : const double alift = adx * adx + ady * ady;
420 : 52 : const double blift = bdx * bdx + bdy * bdy;
421 : 52 : const double clift = cdx * cdx + cdy * cdy;
422 : :
423 : 52 : const double det = alift * (bdxcdy - cdxbdy) + blift * ocad + clift * oabd;
424 : :
425 : 52 : return det > 0;
426 : 169 : }
427 : :
428 : 13 : void Sweep::RotateTrianglePair(Triangle& t, Point& p, Triangle& ot, Point& op) const
429 : : {
430 : : Triangle* n1, *n2, *n3, *n4;
431 : 13 : n1 = t.NeighborCCW(p);
432 : 13 : n2 = t.NeighborCW(p);
433 : 13 : n3 = ot.NeighborCCW(op);
434 : 13 : n4 = ot.NeighborCW(op);
435 : :
436 : : bool ce1, ce2, ce3, ce4;
437 : 13 : ce1 = t.GetConstrainedEdgeCCW(p);
438 : 13 : ce2 = t.GetConstrainedEdgeCW(p);
439 : 13 : ce3 = ot.GetConstrainedEdgeCCW(op);
440 : 13 : ce4 = ot.GetConstrainedEdgeCW(op);
441 : :
442 : : bool de1, de2, de3, de4;
443 : 13 : de1 = t.GetDelunayEdgeCCW(p);
444 : 13 : de2 = t.GetDelunayEdgeCW(p);
445 : 13 : de3 = ot.GetDelunayEdgeCCW(op);
446 : 13 : de4 = ot.GetDelunayEdgeCW(op);
447 : :
448 : 13 : t.Legalize(p, op);
449 : 13 : ot.Legalize(op, p);
450 : :
451 : : // Remap delaunay_edge
452 : 13 : ot.SetDelunayEdgeCCW(p, de1);
453 : 13 : t.SetDelunayEdgeCW(p, de2);
454 : 13 : t.SetDelunayEdgeCCW(op, de3);
455 : 13 : ot.SetDelunayEdgeCW(op, de4);
456 : :
457 : : // Remap constrained_edge
458 : 13 : ot.SetConstrainedEdgeCCW(p, ce1);
459 : 13 : t.SetConstrainedEdgeCW(p, ce2);
460 : 13 : t.SetConstrainedEdgeCCW(op, ce3);
461 : 13 : ot.SetConstrainedEdgeCW(op, ce4);
462 : :
463 : : // Remap neighbors
464 : : // XXX: might optimize the markNeighbor by keeping track of
465 : : // what side should be assigned to what neighbor after the
466 : : // rotation. Now mark neighbor does lots of testing to find
467 : : // the right side.
468 : 13 : t.ClearNeighbors();
469 : 13 : ot.ClearNeighbors();
470 : 13 : if (n1) ot.MarkNeighbor(*n1);
471 : 13 : if (n2) t.MarkNeighbor(*n2);
472 : 13 : if (n3) t.MarkNeighbor(*n3);
473 : 13 : if (n4) ot.MarkNeighbor(*n4);
474 : 13 : t.MarkNeighbor(ot);
475 : 13 : }
476 : :
477 : 13 : void Sweep::FillBasin(SweepContext& tcx, Node& node)
478 : : {
479 : 13 : if (Orient2d(*node.point, *node.next->point, *node.next->next->point) == CCW) {
480 : 13 : tcx.basin.left_node = node.next->next;
481 : 13 : } else {
482 : 0 : tcx.basin.left_node = node.next;
483 : : }
484 : :
485 : : // Find the bottom and right node
486 : 13 : tcx.basin.bottom_node = tcx.basin.left_node;
487 : 13 : while (tcx.basin.bottom_node->next
488 : 13 : && tcx.basin.bottom_node->point->y >= tcx.basin.bottom_node->next->point->y) {
489 : 0 : tcx.basin.bottom_node = tcx.basin.bottom_node->next;
490 : : }
491 : 13 : if (tcx.basin.bottom_node == tcx.basin.left_node) {
492 : : // No valid basin
493 : 13 : return;
494 : : }
495 : :
496 : 0 : tcx.basin.right_node = tcx.basin.bottom_node;
497 : 0 : while (tcx.basin.right_node->next
498 : 0 : && tcx.basin.right_node->point->y < tcx.basin.right_node->next->point->y) {
499 : 0 : tcx.basin.right_node = tcx.basin.right_node->next;
500 : : }
501 : 0 : if (tcx.basin.right_node == tcx.basin.bottom_node) {
502 : : // No valid basins
503 : 0 : return;
504 : : }
505 : :
506 : 0 : tcx.basin.width = tcx.basin.right_node->point->x - tcx.basin.left_node->point->x;
507 : 0 : tcx.basin.left_highest = tcx.basin.left_node->point->y > tcx.basin.right_node->point->y;
508 : :
509 : 0 : FillBasinReq(tcx, tcx.basin.bottom_node);
510 : 13 : }
511 : :
512 : 0 : void Sweep::FillBasinReq(SweepContext& tcx, Node* node)
513 : : {
514 : : // if shallow stop filling
515 : 0 : if (IsShallow(tcx, *node)) {
516 : 0 : return;
517 : : }
518 : :
519 : 0 : Fill(tcx, *node);
520 : :
521 : 0 : if (node->prev == tcx.basin.left_node && node->next == tcx.basin.right_node) {
522 : 0 : return;
523 : 0 : } else if (node->prev == tcx.basin.left_node) {
524 : 0 : Orientation o = Orient2d(*node->point, *node->next->point, *node->next->next->point);
525 : 0 : if (o == CW) {
526 : 0 : return;
527 : : }
528 : 0 : node = node->next;
529 : 0 : } else if (node->next == tcx.basin.right_node) {
530 : 0 : Orientation o = Orient2d(*node->point, *node->prev->point, *node->prev->prev->point);
531 : 0 : if (o == CCW) {
532 : 0 : return;
533 : : }
534 : 0 : node = node->prev;
535 : 0 : } else {
536 : : // Continue with the neighbor node with lowest Y value
537 : 0 : if (node->prev->point->y < node->next->point->y) {
538 : 0 : node = node->prev;
539 : 0 : } else {
540 : 0 : node = node->next;
541 : : }
542 : : }
543 : :
544 : 0 : FillBasinReq(tcx, node);
545 : 0 : }
546 : :
547 : 0 : bool Sweep::IsShallow(SweepContext& tcx, Node& node)
548 : : {
549 : : double height;
550 : :
551 : 0 : if (tcx.basin.left_highest) {
552 : 0 : height = tcx.basin.left_node->point->y - node.point->y;
553 : 0 : } else {
554 : 0 : height = tcx.basin.right_node->point->y - node.point->y;
555 : : }
556 : :
557 : : // if shallow stop filling
558 : 0 : if (tcx.basin.width > height) {
559 : 0 : return true;
560 : : }
561 : 0 : return false;
562 : 0 : }
563 : :
564 : 52 : void Sweep::FillEdgeEvent(SweepContext& tcx, Edge* edge, Node* node)
565 : : {
566 : 52 : if (tcx.edge_event.right) {
567 : 0 : FillRightAboveEdgeEvent(tcx, edge, node);
568 : 0 : } else {
569 : 52 : FillLeftAboveEdgeEvent(tcx, edge, node);
570 : : }
571 : 52 : }
572 : :
573 : 0 : void Sweep::FillRightAboveEdgeEvent(SweepContext& tcx, Edge* edge, Node* node)
574 : : {
575 : 0 : while (node->next->point->x < edge->p->x) {
576 : : // Check if next node is below the edge
577 : 0 : if (Orient2d(*edge->q, *node->next->point, *edge->p) == CCW) {
578 : 0 : FillRightBelowEdgeEvent(tcx, edge, *node);
579 : 0 : } else {
580 : 0 : node = node->next;
581 : : }
582 : : }
583 : 0 : }
584 : :
585 : 0 : void Sweep::FillRightBelowEdgeEvent(SweepContext& tcx, Edge* edge, Node& node)
586 : : {
587 : 0 : if (node.point->x < edge->p->x) {
588 : 0 : if (Orient2d(*node.point, *node.next->point, *node.next->next->point) == CCW) {
589 : : // Concave
590 : 0 : FillRightConcaveEdgeEvent(tcx, edge, node);
591 : 0 : } else{
592 : : // Convex
593 : 0 : FillRightConvexEdgeEvent(tcx, edge, node);
594 : : // Retry this one
595 : 0 : FillRightBelowEdgeEvent(tcx, edge, node);
596 : : }
597 : 0 : }
598 : 0 : }
599 : :
600 : 0 : void Sweep::FillRightConcaveEdgeEvent(SweepContext& tcx, Edge* edge, Node& node)
601 : : {
602 : 0 : Fill(tcx, *node.next);
603 : 0 : if (node.next->point != edge->p) {
604 : : // Next above or below edge?
605 : 0 : if (Orient2d(*edge->q, *node.next->point, *edge->p) == CCW) {
606 : : // Below
607 : 0 : if (Orient2d(*node.point, *node.next->point, *node.next->next->point) == CCW) {
608 : : // Next is concave
609 : 0 : FillRightConcaveEdgeEvent(tcx, edge, node);
610 : 0 : } else {
611 : : // Next is convex
612 : : }
613 : 0 : }
614 : 0 : }
615 : :
616 : 0 : }
617 : :
618 : 0 : void Sweep::FillRightConvexEdgeEvent(SweepContext& tcx, Edge* edge, Node& node)
619 : : {
620 : : // Next concave or convex?
621 : 0 : if (Orient2d(*node.next->point, *node.next->next->point, *node.next->next->next->point) == CCW) {
622 : : // Concave
623 : 0 : FillRightConcaveEdgeEvent(tcx, edge, *node.next);
624 : 0 : } else{
625 : : // Convex
626 : : // Next above or below edge?
627 : 0 : if (Orient2d(*edge->q, *node.next->next->point, *edge->p) == CCW) {
628 : : // Below
629 : 0 : FillRightConvexEdgeEvent(tcx, edge, *node.next);
630 : 0 : } else{
631 : : // Above
632 : : }
633 : : }
634 : 0 : }
635 : :
636 : 52 : void Sweep::FillLeftAboveEdgeEvent(SweepContext& tcx, Edge* edge, Node* node)
637 : : {
638 : 65 : while (node->prev->point->x > edge->p->x) {
639 : : // Check if next node is below the edge
640 : 13 : if (Orient2d(*edge->q, *node->prev->point, *edge->p) == CW) {
641 : 13 : FillLeftBelowEdgeEvent(tcx, edge, *node);
642 : 13 : } else {
643 : 0 : node = node->prev;
644 : : }
645 : : }
646 : 52 : }
647 : :
648 : 13 : void Sweep::FillLeftBelowEdgeEvent(SweepContext& tcx, Edge* edge, Node& node)
649 : : {
650 : 13 : if (node.point->x > edge->p->x) {
651 : 13 : if (Orient2d(*node.point, *node.prev->point, *node.prev->prev->point) == CW) {
652 : : // Concave
653 : 13 : FillLeftConcaveEdgeEvent(tcx, edge, node);
654 : 13 : } else {
655 : : // Convex
656 : 0 : FillLeftConvexEdgeEvent(tcx, edge, node);
657 : : // Retry this one
658 : 0 : FillLeftBelowEdgeEvent(tcx, edge, node);
659 : : }
660 : 13 : }
661 : 13 : }
662 : :
663 : 0 : void Sweep::FillLeftConvexEdgeEvent(SweepContext& tcx, Edge* edge, Node& node)
664 : : {
665 : : // Next concave or convex?
666 : 0 : if (Orient2d(*node.prev->point, *node.prev->prev->point, *node.prev->prev->prev->point) == CW) {
667 : : // Concave
668 : 0 : FillLeftConcaveEdgeEvent(tcx, edge, *node.prev);
669 : 0 : } else{
670 : : // Convex
671 : : // Next above or below edge?
672 : 0 : if (Orient2d(*edge->q, *node.prev->prev->point, *edge->p) == CW) {
673 : : // Below
674 : 0 : FillLeftConvexEdgeEvent(tcx, edge, *node.prev);
675 : 0 : } else{
676 : : // Above
677 : : }
678 : : }
679 : 0 : }
680 : :
681 : 13 : void Sweep::FillLeftConcaveEdgeEvent(SweepContext& tcx, Edge* edge, Node& node)
682 : : {
683 : 13 : Fill(tcx, *node.prev);
684 : 13 : if (node.prev->point != edge->p) {
685 : : // Next above or below edge?
686 : 0 : if (Orient2d(*edge->q, *node.prev->point, *edge->p) == CW) {
687 : : // Below
688 : 0 : if (Orient2d(*node.point, *node.prev->point, *node.prev->prev->point) == CW) {
689 : : // Next is concave
690 : 0 : FillLeftConcaveEdgeEvent(tcx, edge, node);
691 : 0 : } else{
692 : : // Next is convex
693 : : }
694 : 0 : }
695 : 0 : }
696 : :
697 : 13 : }
698 : :
699 : 0 : void Sweep::FlipEdgeEvent(SweepContext& tcx, Point& ep, Point& eq, Triangle* t, Point& p)
700 : : {
701 : 0 : Triangle& ot = t->NeighborAcross(p);
702 : 0 : Point& op = *ot.OppositePoint(*t, p);
703 : :
704 : 0 : if (InScanArea(p, *t->PointCCW(p), *t->PointCW(p), op)) {
705 : : // Lets rotate shared edge one vertex CW
706 : 0 : RotateTrianglePair(*t, p, ot, op);
707 : 0 : tcx.MapTriangleToNodes(*t);
708 : 0 : tcx.MapTriangleToNodes(ot);
709 : :
710 : 0 : if (p == eq && op == ep) {
711 : 0 : if (eq == *tcx.edge_event.constrained_edge->q && ep == *tcx.edge_event.constrained_edge->p) {
712 : 0 : t->MarkConstrainedEdge(&ep, &eq);
713 : 0 : ot.MarkConstrainedEdge(&ep, &eq);
714 : 0 : Legalize(tcx, *t);
715 : 0 : Legalize(tcx, ot);
716 : 0 : } else {
717 : : // XXX: I think one of the triangles should be legalized here?
718 : : }
719 : 0 : } else {
720 : 0 : Orientation o = Orient2d(eq, op, ep);
721 : 0 : t = &NextFlipTriangle(tcx, (int)o, *t, ot, p, op);
722 : 0 : FlipEdgeEvent(tcx, ep, eq, t, p);
723 : : }
724 : 0 : } else {
725 : 0 : Point& newP = NextFlipPoint(ep, eq, ot, op);
726 : 0 : FlipScanEdgeEvent(tcx, ep, eq, *t, ot, newP);
727 : 0 : EdgeEvent(tcx, ep, eq, t, p);
728 : : }
729 : 0 : }
730 : :
731 : 0 : Triangle& Sweep::NextFlipTriangle(SweepContext& tcx, int o, Triangle& t, Triangle& ot, Point& p, Point& op)
732 : : {
733 : 0 : if (o == CCW) {
734 : : // ot is not crossing edge after flip
735 : 0 : int edge_index = ot.EdgeIndex(&p, &op);
736 : 0 : ot.delaunay_edge[edge_index] = true;
737 : 0 : Legalize(tcx, ot);
738 : 0 : ot.ClearDelunayEdges();
739 : 0 : return t;
740 : : }
741 : :
742 : : // t is not crossing edge after flip
743 : 0 : int edge_index = t.EdgeIndex(&p, &op);
744 : :
745 : 0 : t.delaunay_edge[edge_index] = true;
746 : 0 : Legalize(tcx, t);
747 : 0 : t.ClearDelunayEdges();
748 : 0 : return ot;
749 : 0 : }
750 : :
751 : 0 : Point& Sweep::NextFlipPoint(Point& ep, Point& eq, Triangle& ot, Point& op)
752 : : {
753 : 0 : Orientation o2d = Orient2d(eq, op, ep);
754 : 0 : if (o2d == CW) {
755 : : // Right
756 : 0 : return *ot.PointCCW(op);
757 : 0 : } else if (o2d == CCW) {
758 : : // Left
759 : 0 : return *ot.PointCW(op);
760 : : }
761 : 0 : throw std::runtime_error("[Unsupported] Opposing point on constrained edge");
762 : 0 : }
763 : :
764 : 0 : void Sweep::FlipScanEdgeEvent(SweepContext& tcx, Point& ep, Point& eq, Triangle& flip_triangle,
765 : : Triangle& t, Point& p)
766 : : {
767 : 0 : Triangle& ot = t.NeighborAcross(p);
768 : 0 : Point& op = *ot.OppositePoint(t, p);
769 : :
770 : 0 : if (InScanArea(eq, *flip_triangle.PointCCW(eq), *flip_triangle.PointCW(eq), op)) {
771 : : // flip with new edge op->eq
772 : 0 : FlipEdgeEvent(tcx, eq, op, &ot, op);
773 : : // TODO: Actually I just figured out that it should be possible to
774 : : // improve this by getting the next ot and op before the the above
775 : : // flip and continue the flipScanEdgeEvent here
776 : : // set new ot and op here and loop back to inScanArea test
777 : : // also need to set a new flip_triangle first
778 : : // Turns out at first glance that this is somewhat complicated
779 : : // so it will have to wait.
780 : 0 : } else{
781 : 0 : Point& newP = NextFlipPoint(ep, eq, ot, op);
782 : 0 : FlipScanEdgeEvent(tcx, ep, eq, flip_triangle, ot, newP);
783 : : }
784 : 0 : }
785 : :
786 : 13 : Sweep::~Sweep() {
787 : :
788 : : // Clean up memory
789 : 91 : for(size_t i = 0; i < nodes_.size(); i++) {
790 : 78 : delete nodes_[i];
791 : 78 : }
792 : :
793 : 13 : }
794 : :
795 : : }
796 : :
|