Panda3D
Loading...
Searching...
No Matches
triangulator.cxx
Go to the documentation of this file.
1/**
2 * PANDA 3D SOFTWARE
3 * Copyright (c) Carnegie Mellon University. All rights reserved.
4 *
5 * All use of this software is subject to the terms of the revised BSD
6 * license. You should have received a copy of this license along
7 * with this source code in a file named "LICENSE."
8 *
9 * @file triangulator.cxx
10 * @author drose
11 * @date 2007-01-18
12 */
13
14#include "triangulator.h"
15#include "randomizer.h"
16
17/**
18 *
19 */
20Triangulator::
21Triangulator() {
22}
23
24/**
25 * Removes all vertices and polygon specifications from the Triangulator, and
26 * prepares it to start over.
27 */
29clear() {
30 _vertices.clear();
32}
33
34/**
35 * Adds a new vertex to the vertex pool. Returns the vertex index number.
36 */
38add_vertex(const LPoint2d &point) {
39 int index = (int)_vertices.size();
40 _vertices.push_back(point);
41 return index;
42}
43
44/**
45 * Removes the current polygon definition (and its set of holes), but does not
46 * clear the vertex pool.
47 */
50 _polygon.clear();
51 _holes.clear();
52}
53
54/**
55 * Adds the next consecutive vertex of the polygon. This vertex should index
56 * into the vertex pool established by repeated calls to add_vertex().
57 *
58 * The vertices may be listed in either clockwise or counterclockwise order.
59 * Vertices should not be repeated. In particular, do not repeat the first
60 * vertex at the end.
61 */
63add_polygon_vertex(int index) {
64 _polygon.push_back(index);
65}
66
67/**
68 * Finishes the previous hole, if any, and prepares to add a new hole.
69 */
71begin_hole() {
72 _holes.push_back(vector_int());
73}
74
75/**
76 * Adds the next consecutive vertex of the current hole. This vertex should
77 * index into the vertex pool established by repeated calls to add_vertex().
78 *
79 * The vertices may be listed in either clockwise or counterclockwise order.
80 * Vertices should not be repeated.
81 */
83add_hole_vertex(int index) {
84 nassertv(!_holes.empty());
85 _holes.back().push_back(index);
86}
87
88/**
89 * Does the work of triangulating the specified polygon. After this call, you
90 * may retrieve the new triangles one at a time by iterating through
91 * get_triangle_v0/1/2().
92 */
95 _result.clear();
96
97 // Make sure our index numbers are reasonable.
98 cleanup_polygon_indices(_polygon);
99 Holes::iterator hi;
100 for (hi = _holes.begin(); hi != _holes.end(); ++hi) {
101 cleanup_polygon_indices(*hi);
102 }
103
104 if (_polygon.size() < 3) {
105 // Degenerate case.
106 return;
107 }
108
109 // Set up the list of segments.
110 seg.clear();
111 seg.push_back(segment_t()); // we don't use the first entry.
112 make_segment(_polygon, true);
113
114 for (hi = _holes.begin(); hi != _holes.end(); ++hi) {
115 if ((*hi).size() >= 3) {
116 make_segment(*hi, false);
117 }
118 }
119
120 // Shuffle the segment index.
121 int num_segments = (int)seg.size() - 1;
122 permute.reserve(num_segments);
123 int i;
124 for (i = 0; i < num_segments; ++i) {
125 permute.push_back(i + 1);
126 }
127
128 // Actually, I'm not sure why we should shuffle the index. That makes the
129 // result non-deterministic, and isn't one order--for instance, the initial
130 // order--as good as any other?
131 /*
132 Randomizer randomizer;
133 for (i = 0; i < num_segments; ++i) {
134 int j = randomizer.random_int(num_segments);
135 nassertv(j >= 0 && j < num_segments);
136 int t = permute[i];
137 permute[i] = permute[j];
138 permute[j] = t;
139 }
140 */
141 choose_idx = 0;
142
143 /*
144 // cerr << "got " << num_segments << " segments\n";
145 for (i = 1; i < (int)seg.size(); ++i) {
146 segment_t &s = seg[i];
147 printf(" %d. (%g %g), (%g %g)\n", i, s.v0.x, s.v0.y, s.v1.x, s.v1.y);
148 printf(" root0 = %d, root1 = %d\n", s.root0, s.root1);
149 printf(" next = %d, prev = %d\n", s.next, s.prev);
150 }
151 */
152
153 while (construct_trapezoids(num_segments) != 0) {
154 // If there's an error, re-shuffle the index and try again.
155 Randomizer randomizer;
156 for (i = 0; i < num_segments; ++i) {
157 int j = randomizer.random_int(num_segments);
158 nassertv(j >= 0 && j < num_segments);
159 int t = permute[i];
160 permute[i] = permute[j];
161 permute[j] = t;
162 }
163 choose_idx = 0;
164
165 /*
166 // cerr << "got " << num_segments << " segments\n";
167 for (i = 1; i < (int)seg.size(); ++i) {
168 segment_t &s = seg[i];
169 printf(" %d. (%g %g), (%g %g)\n", i, s.v0.x, s.v0.y, s.v1.x, s.v1.y);
170 printf(" root0 = %d, root1 = %d\n", s.root0, s.root1);
171 printf(" next = %d, prev = %d\n", s.next, s.prev);
172 }
173 */
174 }
175
176 /*
177 // cerr << "got " << tr.size() - 1 << " trapezoids\n";
178 for (i = 1; i < (int)tr.size(); ++i) {
179 trap_t &t = tr[i];
180 // cerr << " " << i << ". state = " << t.state << "\n"; cerr << " lseg
181 // = " << t.lseg << " rseg = " << t.rseg << "\n"; cerr << " hi = " <<
182 // t.hi.x << " " << t.hi.y << " lo = " << t.lo.x << " " << t.lo.y << "\n";
183 }
184 */
185
186 int nmonpoly = monotonate_trapezoids(num_segments);
187
188 // cerr << "got " << nmonpoly << " monotone polygons\n";
189
190 triangulate_monotone_polygons(num_segments, nmonpoly);
191
192 /*
193 Result::iterator ri;
194 for (ri = _result.begin(); ri != _result.end(); ++ri) {
195 // cerr << "tri: " << (*ri)._v0 << " " << (*ri)._v1 << " " << (*ri)._v2 <<
196 // "\n";
197 }
198 */
199}
200
201/**
202 * Returns the number of triangles generated by the previous call to
203 * triangulate().
204 */
206get_num_triangles() const {
207 return _result.size();
208}
209
210/**
211 * Returns vertex 0 of the nth triangle generated by the previous call to
212 * triangulate().
213 *
214 * This is a zero-based index into the vertices added by repeated calls to
215 * add_vertex().
216 */
218get_triangle_v0(int n) const {
219 nassertr(n >= 0 && n < (int)_result.size(), -1);
220 return _result[n]._v0;
221}
222
223/**
224 * Returns vertex 1 of the nth triangle generated by the previous call to
225 * triangulate().
226 *
227 * This is a zero-based index into the vertices added by repeated calls to
228 * add_vertex().
229 */
231get_triangle_v1(int n) const {
232 nassertr(n >= 0 && n < (int)_result.size(), -1);
233 return _result[n]._v1;
234}
235
236/**
237 * Returns vertex 2 of the nth triangle generated by the previous call to
238 * triangulate().
239 *
240 * This is a zero-based index into the vertices added by repeated calls to
241 * add_vertex().
242 */
244get_triangle_v2(int n) const {
245 nassertr(n >= 0 && n < (int)_result.size(), -1);
246 return _result[n]._v2;
247}
248
249/**
250 * Removes any invalid index numbers from the list.
251 */
252void Triangulator::
253cleanup_polygon_indices(vector_int &polygon) {
254 // First, check for index bounds.
255 size_t pi = 0;
256 while (pi < polygon.size()) {
257 if (polygon[pi] >= 0 && (size_t)polygon[pi] < _vertices.size()) {
258 // This vertex is OK.
259 ++pi;
260 } else {
261 // This index is out-of-bounds; remove it.
262 polygon.erase(polygon.begin() + pi);
263 }
264 }
265
266 // Now, remove any consecutive repeated vertices.
267 pi = 1;
268 while (pi < polygon.size()) {
269 if (_vertices[polygon[pi]] != _vertices[polygon[pi - 1]]) {
270 // This vertex is OK.
271 ++pi;
272 } else {
273 // This vertex repeats the previous one; remove it.
274 polygon.erase(polygon.begin() + pi);
275 }
276 }
277
278 if (polygon.size() > 1 && _vertices[polygon.back()] == _vertices[polygon.front()]) {
279 // The last vertex repeats the first one; remove it.
280 polygon.pop_back();
281 }
282
283 // Another pass over the polygons, this time removing any "tails".
284 while (polygon.size() >= 3) {
285 bool removed_any = false;
286
287 int prevprev = polygon[polygon.size() - 2];
288 int prev = polygon[polygon.size() - 1];
289
290 for (size_t i = 0; i < polygon.size(); ++i) {
291 int cur = polygon[i];
292 if (_vertices[prevprev] == _vertices[cur]) {
293 // Cut off the tail.
294 removed_any = true;
295 polygon.erase(polygon.begin() + i);
296 if (i == 0) {
297 polygon.pop_back();
298 } else {
299 polygon.erase(polygon.begin() + i - 1);
300 }
301 break;
302 }
303
304 prevprev = prev;
305 prev = cur;
306 }
307
308 // This might have been the tip of a longer tail, so if we removed
309 // something, go again.
310 if (!removed_any) {
311 break;
312 }
313 }
314}
315
316
317// The remainder of the code in this file is adapted more or less from the C
318// code published with the referenced paper.
319
320#define T_X 1
321#define T_Y 2
322#define T_SINK 3
323
324
325#define FIRSTPT 1 /* checking whether pt. is inserted */
326#define LASTPT 2
327
328
329#define REALLY_BIG 1<<30
330#define C_EPS 1.0e-7 /* tolerance value: Used for making */
331 /* all decisions about collinearity or */
332 /* left/right of segment. Decrease */
333 /* this value if the input points are */
334 /* spaced very close together */
335
336
337#define S_LEFT 1 /* for merge-direction */
338#define S_RIGHT 2
339
340
341#define ST_VALID 1 /* for trapezium state */
342#define ST_INVALID 2
343
344
345#define SP_SIMPLE_LRUP 1 /* for splitting trapezoids */
346#define SP_SIMPLE_LRDN 2
347#define SP_2UP_2DN 3
348#define SP_2UP_LEFT 4
349#define SP_2UP_RIGHT 5
350#define SP_2DN_LEFT 6
351#define SP_2DN_RIGHT 7
352#define SP_NOSPLIT -1
353
354#define TR_FROM_UP 1 /* for traverse-direction */
355#define TR_FROM_DN 2
356
357#define TRI_LHS 1
358#define TRI_RHS 2
359
360
361#define CROSS(v0, v1, v2) (((v1).x - (v0).x)*((v2).y - (v0).y) - \
362 ((v1).y - (v0).y)*((v2).x - (v0).x))
363
364#define DOT(v0, v1) ((v0).x * (v1).x + (v0).y * (v1).y)
365
366#define FP_EQUAL(s, t) (fabs(s - t) <= C_EPS)
367
368
369#define CROSS_SINE(v0, v1) ((v0).x * (v1).y - (v1).x * (v0).y)
370#define LENGTH(v0) (sqrt((v0).x * (v0).x + (v0).y * (v0).y))
371
372
373/**
374 * Returns true if the list of vertices is counter-clockwise, false if it is
375 * clockwise.
376 */
377bool Triangulator::
378check_left_winding(const vector_int &range) const {
379 // We do this by computing the polygon's signed area. If it comes out
380 // negative, the polygon is right-winding.
381
382 double area = 0.0;
383 size_t j = range.size() - 1;
384 for (size_t i = 0; i < range.size(); ++i) {
385 const LPoint2d &p0 = _vertices[range[j]];
386 const LPoint2d &p1 = _vertices[range[i]];
387 area += p0[0] * p1[1] - p0[1] * p1[0];
388 j = i;
389 }
390
391 return area >= 0.0;
392}
393
394/**
395 * Converts a linear list of integer vertices to a list of segment_t. If
396 * want_left_winding is true, the list is reversed if necessary to make it
397 * left-winding; otherwise, it is reversed to make it right-winding.
398 */
399void Triangulator::
400make_segment(const vector_int &range, bool want_left_winding) {
401 int num_points = (int)range.size();
402 nassertv(num_points >= 2);
403
404 int first = (int)seg.size();
405 int last = first + num_points - 1;
406
407 if (want_left_winding == check_left_winding(range)) {
408 // Keep it in its natural order.
409 int first = (int)seg.size();
410 int last = first + num_points - 1;
411
412 seg.push_back(segment_t(this, range[0], range[1],
413 last, first + 1));
414
415 for (int i = 1; i < num_points - 1; ++i) {
416 seg.push_back(segment_t(this, range[i], range[i + 1],
417 first + i - 1, first + i + 1));
418 }
419
420 seg.push_back(segment_t(this, range[num_points - 1], range[0],
421 last - 1, first));
422
423 } else {
424 // Reverse it.
425 seg.push_back(segment_t(this, range[0], range[num_points - 1],
426 last, first + 1));
427
428 for (int i = 1; i < num_points - 1; ++i) {
429 seg.push_back(segment_t(this, range[num_points - i], range[num_points - i - 1],
430 first + i - 1, first + i + 1));
431 }
432
433 seg.push_back(segment_t(this, range[1], range[0],
434 last - 1, first));
435
436 }
437}
438
439/* Return the next segment in the generated random ordering of all the */
440/* segments in S */
441int Triangulator::
442choose_segment() {
443 nassertr(choose_idx < (int)permute.size(), 0);
444 // segment_t &s = seg[permute[choose_idx]]; cerr << "choose_segment " <<
445 // permute[choose_idx] << ": " << s.v0.x << ", " << s.v0.y << " to " <<
446 // s.v1.x << ", " << s.v1.y << "\n";
447 return permute[choose_idx++];
448}
449
450double Triangulator::
451math_log2(double v) {
452 static const double log2 = log(2.0);
453 return log(v) / log2;
454}
455
456/* Get log*n for given n */
457int Triangulator::
458math_logstar_n(int n) {
459 int i;
460 double v;
461
462 for (i = 0, v = (double) n; v >= 1; i++)
463 v = math_log2(v);
464
465 return (i - 1);
466}
467
468
469int Triangulator::
470math_N(int n, int h) {
471 int i;
472 double v;
473
474 for (i = 0, v = (int) n; i < h; i++)
475 v = math_log2(v);
476
477 return (int) ceil((double) 1.0*n/v);
478}
479
480
481/* Return a new node to be added into the query tree */
482int Triangulator::newnode() {
483 int index = (int)qs.size();
484 qs.push_back(node_t());
485 // cerr << "creating new node " << index << "\n";
486 return index;
487}
488
489/* Return a free trapezoid */
490int Triangulator::newtrap() {
491 int tr_idx = (int)tr.size();
492 tr.push_back(trap_t());
493 tr[tr_idx].lseg = -1;
494 tr[tr_idx].rseg = -1;
495 tr[tr_idx].state = ST_VALID;
496 // cerr << "creating new trapezoid " << tr_idx << "\n";
497 return tr_idx;
498}
499
500
501/* Return the maximum of the two points into the yval structure */
502int Triangulator::_max(point_t *yval, point_t *v0, point_t *v1) {
503 if (v0->y > v1->y + C_EPS)
504 *yval = *v0;
505 else if (FP_EQUAL(v0->y, v1->y))
506 {
507 if (v0->x > v1->x + C_EPS)
508 *yval = *v0;
509 else
510 *yval = *v1;
511 }
512 else
513 *yval = *v1;
514
515 return 0;
516}
517
518
519/* Return the minimum of the two points into the yval structure */
520int Triangulator::_min(point_t *yval, point_t *v0, point_t *v1) {
521 if (v0->y < v1->y - C_EPS)
522 *yval = *v0;
523 else if (FP_EQUAL(v0->y, v1->y))
524 {
525 if (v0->x < v1->x)
526 *yval = *v0;
527 else
528 *yval = *v1;
529 }
530 else
531 *yval = *v1;
532
533 return 0;
534}
535
536
537int Triangulator::
538_greater_than(point_t *v0, point_t *v1) {
539 if (v0->y > v1->y + C_EPS)
540 return true;
541 else if (v0->y < v1->y - C_EPS)
542 return false;
543 else
544 return (v0->x > v1->x);
545}
546
547
548int Triangulator::
549_equal_to(point_t *v0, point_t *v1) {
550 return (FP_EQUAL(v0->y, v1->y) && FP_EQUAL(v0->x, v1->x));
551}
552
553int Triangulator::
554_greater_than_equal_to(point_t *v0, point_t *v1) {
555 if (v0->y > v1->y + C_EPS)
556 return true;
557 else if (v0->y < v1->y - C_EPS)
558 return false;
559 else
560 return (v0->x >= v1->x);
561}
562
563int Triangulator::
564_less_than(point_t *v0, point_t *v1) {
565 if (v0->y < v1->y - C_EPS)
566 return true;
567 else if (v0->y > v1->y + C_EPS)
568 return false;
569 else
570 return (v0->x < v1->x);
571}
572
573
574/* Initilialise the query structure (Q) and the trapezoid table (T)
575 * when the first segment is added to start the trapezoidation. The
576 * query-tree starts out with 4 trapezoids, one S-node and 2 Y-nodes
577 *
578 * 4
579 * -----------------------------------
580 * \
581 * 1 \ 2
582 * \
583 * -----------------------------------
584 * 3
585 */
586
587int Triangulator::
588init_query_structure(int segnum) {
589 int i1, i2, i3, i4, i5, i6, i7, root;
590 int t1, t2, t3, t4;
591 segment_t *s = &seg[segnum];
592
593 tr.clear();
594 qs.clear();
595
596 // We don't use the first elements.
597 tr.push_back(trap_t());
598 qs.push_back(node_t());
599
600 i1 = newnode();
601 qs[i1].nodetype = T_Y;
602 _max(&qs[i1].yval, &s->v0, &s->v1); /* root */
603 root = i1;
604
605 i2 = newnode();
606 qs[i1].right = i2;
607 qs[i2].nodetype = T_SINK;
608 qs[i2].parent = i1;
609
610 i3 = newnode();
611 qs[i1].left = i3;
612
613 qs[i3].nodetype = T_Y;
614 _min(&qs[i3].yval, &s->v0, &s->v1); /* root */
615 qs[i3].parent = i1;
616
617 i4 = newnode();
618 qs[i3].left = i4;
619 qs[i4].nodetype = T_SINK;
620 qs[i4].parent = i3;
621
622 i5 = newnode();
623 qs[i3].right = i5;
624 qs[i5].nodetype = T_X;
625 qs[i5].segnum = segnum;
626 qs[i5].parent = i3;
627
628 i6 = newnode();
629 qs[i5].left = i6;
630 qs[i6].nodetype = T_SINK;
631 qs[i6].parent = i5;
632
633 i7 = newnode();
634 qs[i5].right = i7;
635 qs[i7].nodetype = T_SINK;
636 qs[i7].parent = i5;
637
638 t1 = newtrap(); /* middle left */
639 t2 = newtrap(); /* middle right */
640 t3 = newtrap(); /* bottom-most */
641 t4 = newtrap(); /* topmost */
642
643 tr[t4].lo = qs[i1].yval;
644 tr[t2].hi = qs[i1].yval;
645 tr[t1].hi = qs[i1].yval;
646 tr[t3].hi = qs[i3].yval;
647 tr[t2].lo = qs[i3].yval;
648 tr[t1].lo = qs[i3].yval;
649 tr[t4].hi.y = (double) (REALLY_BIG);
650 tr[t4].hi.x = (double) (REALLY_BIG);
651 tr[t3].lo.y = (double) -1* (REALLY_BIG);
652 tr[t3].lo.x = (double) -1* (REALLY_BIG);
653 tr[t2].lseg = segnum;
654 tr[t1].rseg = segnum;
655 tr[t2].u0 = t4;
656 tr[t1].u0 = t4;
657 tr[t2].d0 = t3;
658 tr[t1].d0 = t3;
659 tr[t3].u0 = t1;
660 tr[t4].d0 = t1;
661 tr[t3].u1 = t2;
662 tr[t4].d1 = t2;
663
664 tr[t1].sink = i6;
665 tr[t2].sink = i7;
666 tr[t3].sink = i4;
667 tr[t4].sink = i2;
668
669 tr[t2].state = ST_VALID;
670 tr[t1].state = ST_VALID;
671 tr[t4].state = ST_VALID;
672 tr[t3].state = ST_VALID;
673
674 qs[i2].trnum = t4;
675 qs[i4].trnum = t3;
676 qs[i6].trnum = t1;
677 qs[i7].trnum = t2;
678
679 s->is_inserted = true;
680 return root;
681}
682
683
684/* Retun true if the vertex v is to the left of line segment no.
685 * segnum. Takes care of the degenerate cases when both the vertices
686 * have the same y--cood, etc.
687 */
688
689int Triangulator::
690is_left_of(int segnum, point_t *v) {
691 segment_t *s = &seg[segnum];
692 double area;
693
694 if (_greater_than(&s->v1, &s->v0)) /* seg. going upwards */
695 {
696 if (FP_EQUAL(s->v1.y, v->y))
697 {
698 if (v->x < s->v1.x)
699 area = 1.0;
700 else
701 area = -1.0;
702 }
703 else if (FP_EQUAL(s->v0.y, v->y))
704 {
705 if (v->x < s->v0.x)
706 area = 1.0;
707 else
708 area = -1.0;
709 }
710 else
711 area = CROSS(s->v0, s->v1, (*v));
712 }
713 else /* v0 > v1 */
714 {
715 if (FP_EQUAL(s->v1.y, v->y))
716 {
717 if (v->x < s->v1.x)
718 area = 1.0;
719 else
720 area = -1.0;
721 }
722 else if (FP_EQUAL(s->v0.y, v->y))
723 {
724 if (v->x < s->v0.x)
725 area = 1.0;
726 else
727 area = -1.0;
728 }
729 else
730 area = CROSS(s->v1, s->v0, (*v));
731 }
732
733 if (area > 0.0)
734 return true;
735 else
736 return false;
737}
738
739
740
741/* Returns true if the corresponding endpoint of the given segment is */
742/* already inserted into the segment tree. Use the simple test of */
743/* whether the segment which shares this endpoint is already inserted */
744
745int Triangulator::
746inserted(int segnum, int whichpt) {
747 if (whichpt == FIRSTPT)
748 return seg[seg[segnum].prev].is_inserted;
749 else
750 return seg[seg[segnum].next].is_inserted;
751}
752
753/* This is query routine which determines which trapezoid does the
754 * point v lie in. The return value is the trapezoid number.
755 */
756
757int Triangulator::
758locate_endpoint(point_t *v, point_t *vo, int r) {
759 // cerr << "locate_endpoint(" << v->x << " " << v->y << ", " << vo->x << " "
760 // << vo->y << ", " << r << ")\n";
761 node_t *rptr = &qs[r];
762
763 switch (rptr->nodetype)
764 {
765 case T_SINK:
766 return rptr->trnum;
767
768 case T_Y:
769 if (_greater_than(v, &rptr->yval)) /* above */
770 return locate_endpoint(v, vo, rptr->right);
771 else if (_equal_to(v, &rptr->yval)) /* the point is already */
772 { /* inserted. */
773 if (_greater_than(vo, &rptr->yval)) /* above */
774 return locate_endpoint(v, vo, rptr->right);
775 else
776 return locate_endpoint(v, vo, rptr->left); /* below */
777 }
778 else
779 return locate_endpoint(v, vo, rptr->left); /* below */
780
781 case T_X:
782 if (_equal_to(v, &seg[rptr->segnum].v0) ||
783 _equal_to(v, &seg[rptr->segnum].v1))
784 {
785 if (FP_EQUAL(v->y, vo->y)) /* horizontal segment */
786 {
787 if (vo->x < v->x)
788 return locate_endpoint(v, vo, rptr->left); /* left */
789 else
790 return locate_endpoint(v, vo, rptr->right); /* right */
791 }
792
793 else if (is_left_of(rptr->segnum, vo))
794 return locate_endpoint(v, vo, rptr->left); /* left */
795 else
796 return locate_endpoint(v, vo, rptr->right); /* right */
797 }
798 else if (is_left_of(rptr->segnum, v))
799 return locate_endpoint(v, vo, rptr->left); /* left */
800 else
801 return locate_endpoint(v, vo, rptr->right); /* right */
802
803 default:
804 fprintf(stderr, "Haggu !!!!!\n");
805 nassertr(false, -1);
806 return -1;
807 }
808}
809
810
811/* Thread in the segment into the existing trapezoidation. The
812 * limiting trapezoids are given by tfirst and tlast (which are the
813 * trapezoids containing the two endpoints of the segment. Merges all
814 * possible trapezoids which flank this segment and have been recently
815 * divided because of its insertion
816 */
817
818int Triangulator::
819merge_trapezoids(int segnum, int tfirst, int tlast, int side) {
820 int t, tnext, cond;
821 int ptnext;
822
823 // cerr << "merge_trapezoids(" << segnum << ", " << tfirst << ", " << tlast
824 // << ", " << side << ")\n";
825
826 /* First merge polys on the LHS */
827 t = tfirst;
828 while ((t > 0) && _greater_than_equal_to(&tr[t].lo, &tr[tlast].lo))
829 {
830 if (side == S_LEFT)
831 cond = ((((tnext = tr[t].d0) > 0) && (tr[tnext].rseg == segnum)) ||
832 (((tnext = tr[t].d1) > 0) && (tr[tnext].rseg == segnum)));
833 else
834 cond = ((((tnext = tr[t].d0) > 0) && (tr[tnext].lseg == segnum)) ||
835 (((tnext = tr[t].d1) > 0) && (tr[tnext].lseg == segnum)));
836
837 if (cond)
838 {
839 if ((tr[t].lseg == tr[tnext].lseg) &&
840 (tr[t].rseg == tr[tnext].rseg)) /* good neighbours */
841 { /* merge them */
842 /* Use the upper node as the new node i.e. t */
843
844 ptnext = qs[tr[tnext].sink].parent;
845
846 if (qs[ptnext].left == tr[tnext].sink)
847 qs[ptnext].left = tr[t].sink;
848 else
849 qs[ptnext].right = tr[t].sink; /* redirect parent */
850
851
852 /* Change the upper neighbours of the lower trapezoids */
853
854 if ((tr[t].d0 = tr[tnext].d0) > 0) {
855 if (tr[tr[t].d0].u0 == tnext) {
856 tr[tr[t].d0].u0 = t;
857 } else if (tr[tr[t].d0].u1 == tnext) {
858 tr[tr[t].d0].u1 = t;
859 }
860 }
861
862 if ((tr[t].d1 = tr[tnext].d1) > 0) {
863 if (tr[tr[t].d1].u0 == tnext) {
864 tr[tr[t].d1].u0 = t;
865 } else if (tr[tr[t].d1].u1 == tnext) {
866 tr[tr[t].d1].u1 = t;
867 }
868 }
869
870 tr[t].lo = tr[tnext].lo;
871 tr[tnext].state = ST_INVALID; /* invalidate the lower */
872 /* trapezium */
873 }
874 else /* not good neighbours */
875 t = tnext;
876 }
877 else /* do not satisfy the outer if */
878 t = tnext;
879
880 } /* end-while */
881
882 return 0;
883}
884
885
886/* Add in the new segment into the trapezoidation and update Q and T
887 * structures. First locate the two endpoints of the segment in the
888 * Q-structure. Then start from the topmost trapezoid and go down to
889 * the lower trapezoid dividing all the trapezoids in between .
890 */
891
892int Triangulator::
893add_segment(int segnum) {
894 // cerr << "add_segment(" << segnum << ")\n";
895
896 segment_t s;
897 // segment_t *so = &seg[segnum];
898 int tu, tl, sk, tfirst, tlast; //, tnext;
899 int tfirstr = 0, tlastr = 0, tfirstl = 0, tlastl = 0;
900 int i1, i2, t, tn; // t1, t2,
901 point_t tpt;
902 int tribot = 0, is_swapped = 0;
903 int tmptriseg;
904
905 s = seg[segnum];
906 if (_greater_than(&s.v1, &s.v0)) /* Get higher vertex in v0 */
907 {
908 int tmp;
909 tpt = s.v0;
910 s.v0 = s.v1;
911 s.v1 = tpt;
912 tmp = s.root0;
913 s.root0 = s.root1;
914 s.root1 = tmp;
915 is_swapped = true;
916 }
917
918 if ((is_swapped) ? !inserted(segnum, LASTPT) :
919 !inserted(segnum, FIRSTPT)) /* insert v0 in the tree */
920 {
921 int tmp_d;
922
923 tu = locate_endpoint(&s.v0, &s.v1, s.root0);
924 tl = newtrap(); /* tl is the new lower trapezoid */
925 tr[tl].state = ST_VALID;
926 tr[tl] = tr[tu];
927 tr[tl].hi.y = s.v0.y;
928 tr[tu].lo.y = s.v0.y;
929 tr[tl].hi.x = s.v0.x;
930 tr[tu].lo.x = s.v0.x;
931 tr[tu].d0 = tl;
932 tr[tu].d1 = 0;
933 tr[tl].u0 = tu;
934 tr[tl].u1 = 0;
935
936 if (((tmp_d = tr[tl].d0) > 0) && (tr[tmp_d].u0 == tu))
937 tr[tmp_d].u0 = tl;
938 if (((tmp_d = tr[tl].d0) > 0) && (tr[tmp_d].u1 == tu))
939 tr[tmp_d].u1 = tl;
940
941 if (((tmp_d = tr[tl].d1) > 0) && (tr[tmp_d].u0 == tu))
942 tr[tmp_d].u0 = tl;
943 if (((tmp_d = tr[tl].d1) > 0) && (tr[tmp_d].u1 == tu))
944 tr[tmp_d].u1 = tl;
945
946 /* Now update the query structure and obtain the sinks for the */
947 /* two trapezoids */
948
949 i1 = newnode(); /* Upper trapezoid sink */
950 i2 = newnode(); /* Lower trapezoid sink */
951 sk = tr[tu].sink;
952
953 qs[sk].nodetype = T_Y;
954 qs[sk].yval = s.v0;
955 qs[sk].segnum = segnum; /* not really reqd ... maybe later */
956 qs[sk].left = i2;
957 qs[sk].right = i1;
958
959 qs[i1].nodetype = T_SINK;
960 qs[i1].trnum = tu;
961 qs[i1].parent = sk;
962
963 qs[i2].nodetype = T_SINK;
964 qs[i2].trnum = tl;
965 qs[i2].parent = sk;
966
967 tr[tu].sink = i1;
968 tr[tl].sink = i2;
969 tfirst = tl;
970 }
971 else /* v0 already present */
972 { /* Get the topmost intersecting trapezoid */
973 tfirst = locate_endpoint(&s.v0, &s.v1, s.root0);
974 }
975
976
977 if ((is_swapped) ? !inserted(segnum, FIRSTPT) :
978 !inserted(segnum, LASTPT)) /* insert v1 in the tree */
979 {
980 int tmp_d;
981
982 tu = locate_endpoint(&s.v1, &s.v0, s.root1);
983
984 tl = newtrap(); /* tl is the new lower trapezoid */
985 tr[tl].state = ST_VALID;
986 tr[tl] = tr[tu];
987 tr[tl].hi.y = s.v1.y;
988 tr[tu].lo.y = s.v1.y;
989 tr[tl].hi.x = s.v1.x;
990 tr[tu].lo.x = s.v1.x;
991 tr[tu].d0 = tl;
992 tr[tu].d1 = 0;
993 tr[tl].u0 = tu;
994 tr[tl].u1 = 0;
995
996 if (((tmp_d = tr[tl].d0) > 0) && (tr[tmp_d].u0 == tu))
997 tr[tmp_d].u0 = tl;
998 if (((tmp_d = tr[tl].d0) > 0) && (tr[tmp_d].u1 == tu))
999 tr[tmp_d].u1 = tl;
1000
1001 if (((tmp_d = tr[tl].d1) > 0) && (tr[tmp_d].u0 == tu))
1002 tr[tmp_d].u0 = tl;
1003 if (((tmp_d = tr[tl].d1) > 0) && (tr[tmp_d].u1 == tu))
1004 tr[tmp_d].u1 = tl;
1005
1006 /* Now update the query structure and obtain the sinks for the */
1007 /* two trapezoids */
1008
1009 i1 = newnode(); /* Upper trapezoid sink */
1010 i2 = newnode(); /* Lower trapezoid sink */
1011 sk = tr[tu].sink;
1012
1013 qs[sk].nodetype = T_Y;
1014 qs[sk].yval = s.v1;
1015 qs[sk].segnum = segnum; /* not really reqd ... maybe later */
1016 qs[sk].left = i2;
1017 qs[sk].right = i1;
1018
1019 qs[i1].nodetype = T_SINK;
1020 qs[i1].trnum = tu;
1021 qs[i1].parent = sk;
1022
1023 qs[i2].nodetype = T_SINK;
1024 qs[i2].trnum = tl;
1025 qs[i2].parent = sk;
1026
1027 tr[tu].sink = i1;
1028 tr[tl].sink = i2;
1029 tlast = tu;
1030 }
1031 else /* v1 already present */
1032 { /* Get the lowermost intersecting trapezoid */
1033 tlast = locate_endpoint(&s.v1, &s.v0, s.root1);
1034 tribot = 1;
1035 }
1036
1037 /* Thread the segment into the query tree creating a new X-node */
1038 /* First, split all the trapezoids which are intersected by s into */
1039 /* two */
1040
1041 t = tfirst; /* topmost trapezoid */
1042
1043 while ((t > 0) &&
1044 _greater_than_equal_to(&tr[t].lo, &tr[tlast].lo))
1045 /* traverse from top to bot */
1046 {
1047 int t_sav, tn_sav;
1048 sk = tr[t].sink;
1049 i1 = newnode(); /* left trapezoid sink */
1050 i2 = newnode(); /* right trapezoid sink */
1051
1052 qs[sk].nodetype = T_X;
1053 qs[sk].segnum = segnum;
1054 qs[sk].left = i1;
1055 qs[sk].right = i2;
1056
1057 qs[i1].nodetype = T_SINK; /* left trapezoid (use existing one) */
1058 qs[i1].trnum = t;
1059 qs[i1].parent = sk;
1060
1061 qs[i2].nodetype = T_SINK; /* right trapezoid (allocate new) */
1062 tn = newtrap();
1063 qs[i2].trnum = tn;
1064 tr[tn].state = ST_VALID;
1065 qs[i2].parent = sk;
1066
1067 if (t == tfirst)
1068 tfirstr = tn;
1069 if (_equal_to(&tr[t].lo, &tr[tlast].lo))
1070 tlastr = tn;
1071
1072 tr[tn] = tr[t];
1073 tr[t].sink = i1;
1074 tr[tn].sink = i2;
1075 t_sav = t;
1076 tn_sav = tn;
1077
1078 /* error */
1079
1080 if ((tr[t].d0 <= 0) && (tr[t].d1 <= 0)) /* case cannot arise */
1081 {
1082 /* Actually, this case does sometimes arise. Huh. */
1083 fprintf(stderr, "add_segment: error\n");
1084 return 1;
1085 }
1086
1087 /* only one trapezoid below. partition t into two and make the */
1088 /* two resulting trapezoids t and tn as the upper neighbours of */
1089 /* the sole lower trapezoid */
1090
1091 else if ((tr[t].d0 > 0) && (tr[t].d1 <= 0))
1092 { /* Only one trapezoid below */
1093 if ((tr[t].u0 > 0) && (tr[t].u1 > 0))
1094 { /* continuation of a chain from abv. */
1095 if (tr[t].usave > 0) /* three upper neighbours */
1096 {
1097 if (tr[t].uside == S_LEFT)
1098 {
1099 tr[tn].u0 = tr[t].u1;
1100 tr[t].u1 = -1;
1101 tr[tn].u1 = tr[t].usave;
1102
1103 tr[tr[t].u0].d0 = t;
1104 tr[tr[tn].u0].d0 = tn;
1105 tr[tr[tn].u1].d0 = tn;
1106 }
1107 else /* intersects in the right */
1108 {
1109 tr[tn].u1 = -1;
1110 tr[tn].u0 = tr[t].u1;
1111 tr[t].u1 = tr[t].u0;
1112 tr[t].u0 = tr[t].usave;
1113
1114 tr[tr[t].u0].d0 = t;
1115 tr[tr[t].u1].d0 = t;
1116 tr[tr[tn].u0].d0 = tn;
1117 }
1118 tr[tn].usave = 0;
1119 tr[t].usave = 0;
1120 }
1121 else /* No usave.... simple case */
1122 {
1123 tr[tn].u0 = tr[t].u1;
1124 tr[tn].u1 = -1;
1125 tr[t].u1 = -1;
1126 tr[tr[tn].u0].d0 = tn;
1127 }
1128 }
1129 else
1130 { /* fresh seg. or upward cusp */
1131 int tmp_u = tr[t].u0;
1132 int td0, td1;
1133 if (((td0 = tr[tmp_u].d0) > 0) &&
1134 ((td1 = tr[tmp_u].d1) > 0))
1135 { /* upward cusp */
1136 if ((tr[td0].rseg > 0) &&
1137 !is_left_of(tr[td0].rseg, &s.v1))
1138 {
1139 tr[tn].u1 = -1;
1140 tr[t].u1 = -1;
1141 tr[t].u0 = -1;
1142 tr[tr[tn].u0].d1 = tn;
1143 }
1144 else /* cusp going leftwards */
1145 {
1146 tr[t].u1 = -1;
1147 tr[tn].u1 = -1;
1148 tr[tn].u0 = -1;
1149 tr[tr[t].u0].d0 = t;
1150 }
1151 }
1152 else /* fresh segment */
1153 {
1154 tr[tr[t].u0].d0 = t;
1155 tr[tr[t].u0].d1 = tn;
1156 }
1157 }
1158
1159 if (FP_EQUAL(tr[t].lo.y, tr[tlast].lo.y) &&
1160 FP_EQUAL(tr[t].lo.x, tr[tlast].lo.x) && tribot)
1161 { /* bottom forms a triangle */
1162
1163 if (is_swapped)
1164 tmptriseg = seg[segnum].prev;
1165 else
1166 tmptriseg = seg[segnum].next;
1167
1168 if ((tmptriseg > 0) && is_left_of(tmptriseg, &s.v0))
1169 {
1170 /* L-R downward cusp */
1171 tr[tr[t].d0].u0 = t;
1172 tr[tn].d1 = -1;
1173 tr[tn].d0 = -1;
1174 }
1175 else
1176 {
1177 /* R-L downward cusp */
1178 tr[tr[tn].d0].u1 = tn;
1179 tr[t].d1 = -1;
1180 tr[t].d0 = -1;
1181 }
1182 }
1183 else
1184 {
1185 if ((tr[tr[t].d0].u0 > 0) && (tr[tr[t].d0].u1 > 0))
1186 {
1187 if (tr[tr[t].d0].u0 == t) /* passes thru LHS */
1188 {
1189 tr[tr[t].d0].usave = tr[tr[t].d0].u1;
1190 tr[tr[t].d0].uside = S_LEFT;
1191 }
1192 else
1193 {
1194 tr[tr[t].d0].usave = tr[tr[t].d0].u0;
1195 tr[tr[t].d0].uside = S_RIGHT;
1196 }
1197 }
1198 tr[tr[t].d0].u0 = t;
1199 tr[tr[t].d0].u1 = tn;
1200 }
1201
1202 t = tr[t].d0;
1203 }
1204
1205
1206 else if ((tr[t].d0 <= 0) && (tr[t].d1 > 0))
1207 { /* Only one trapezoid below */
1208 if ((tr[t].u0 > 0) && (tr[t].u1 > 0))
1209 { /* continuation of a chain from abv. */
1210 if (tr[t].usave > 0) /* three upper neighbours */
1211 {
1212 if (tr[t].uside == S_LEFT)
1213 {
1214 tr[tn].u0 = tr[t].u1;
1215 tr[t].u1 = -1;
1216 tr[tn].u1 = tr[t].usave;
1217
1218 tr[tr[t].u0].d0 = t;
1219 tr[tr[tn].u0].d0 = tn;
1220 tr[tr[tn].u1].d0 = tn;
1221 }
1222 else /* intersects in the right */
1223 {
1224 tr[tn].u1 = -1;
1225 tr[tn].u0 = tr[t].u1;
1226 tr[t].u1 = tr[t].u0;
1227 tr[t].u0 = tr[t].usave;
1228
1229 tr[tr[t].u0].d0 = t;
1230 tr[tr[t].u1].d0 = t;
1231 tr[tr[tn].u0].d0 = tn;
1232 }
1233
1234 tr[tn].usave = 0;
1235 tr[t].usave = 0;
1236 }
1237 else /* No usave.... simple case */
1238 {
1239 tr[tn].u0 = tr[t].u1;
1240 tr[tn].u1 = -1;
1241 tr[t].u1 = -1;
1242 tr[tr[tn].u0].d0 = tn;
1243 }
1244 }
1245 else
1246 { /* fresh seg. or upward cusp */
1247 int tmp_u = tr[t].u0;
1248 int td0, td1;
1249 if (((td0 = tr[tmp_u].d0) > 0) &&
1250 ((td1 = tr[tmp_u].d1) > 0))
1251 { /* upward cusp */
1252 if ((tr[td0].rseg > 0) &&
1253 !is_left_of(tr[td0].rseg, &s.v1))
1254 {
1255 tr[tn].u1 = -1;
1256 tr[t].u1 = -1;
1257 tr[t].u0 = -1;
1258 tr[tr[tn].u0].d1 = tn;
1259 }
1260 else
1261 {
1262 tr[t].u1 = -1;
1263 tr[tn].u1 = -1;
1264 tr[tn].u0 = -1;
1265 tr[tr[t].u0].d0 = t;
1266 }
1267 }
1268 else /* fresh segment */
1269 {
1270 tr[tr[t].u0].d0 = t;
1271 tr[tr[t].u0].d1 = tn;
1272 }
1273 }
1274
1275 if (FP_EQUAL(tr[t].lo.y, tr[tlast].lo.y) &&
1276 FP_EQUAL(tr[t].lo.x, tr[tlast].lo.x) && tribot)
1277 { /* bottom forms a triangle */
1278 if (is_swapped)
1279 tmptriseg = seg[segnum].prev;
1280 else
1281 tmptriseg = seg[segnum].next;
1282
1283 if ((tmptriseg > 0) && is_left_of(tmptriseg, &s.v0))
1284 {
1285 /* L-R downward cusp */
1286 tr[tr[t].d1].u0 = t;
1287 tr[tn].d1 = -1;
1288 tr[tn].d0 = -1;
1289 }
1290 else
1291 {
1292 /* R-L downward cusp */
1293 tr[tr[tn].d1].u1 = tn;
1294 tr[t].d1 = -1;
1295 tr[t].d0 = -1;
1296 }
1297 }
1298 else
1299 {
1300 if ((tr[tr[t].d1].u0 > 0) && (tr[tr[t].d1].u1 > 0))
1301 {
1302 if (tr[tr[t].d1].u0 == t) /* passes thru LHS */
1303 {
1304 tr[tr[t].d1].usave = tr[tr[t].d1].u1;
1305 tr[tr[t].d1].uside = S_LEFT;
1306 }
1307 else
1308 {
1309 tr[tr[t].d1].usave = tr[tr[t].d1].u0;
1310 tr[tr[t].d1].uside = S_RIGHT;
1311 }
1312 }
1313 tr[tr[t].d1].u0 = t;
1314 tr[tr[t].d1].u1 = tn;
1315 }
1316
1317 t = tr[t].d1;
1318 }
1319
1320 /* two trapezoids below. Find out which one is intersected by */
1321 /* this segment and proceed down that one */
1322
1323 else
1324 {
1325 // int tmpseg = tr[tr[t].d0].rseg;
1326 double y0, yt;
1327 point_t tmppt;
1328 int tnext, i_d0;
1329
1330 i_d0 = false;
1331 if (FP_EQUAL(tr[t].lo.y, s.v0.y))
1332 {
1333 if (tr[t].lo.x > s.v0.x)
1334 i_d0 = true;
1335 }
1336 else
1337 {
1338 y0 = tr[t].lo.y;
1339 tmppt.y = y0;
1340 yt = (y0 - s.v0.y)/(s.v1.y - s.v0.y);
1341 tmppt.x = s.v0.x + yt * (s.v1.x - s.v0.x);
1342
1343 if (_less_than(&tmppt, &tr[t].lo))
1344 i_d0 = true;
1345 }
1346
1347 /* check continuity from the top so that the lower-neighbour */
1348 /* values are properly filled for the upper trapezoid */
1349
1350 if ((tr[t].u0 > 0) && (tr[t].u1 > 0))
1351 { /* continuation of a chain from abv. */
1352 if (tr[t].usave > 0) /* three upper neighbours */
1353 {
1354 if (tr[t].uside == S_LEFT)
1355 {
1356 tr[tn].u0 = tr[t].u1;
1357 tr[t].u1 = -1;
1358 tr[tn].u1 = tr[t].usave;
1359
1360 tr[tr[t].u0].d0 = t;
1361 tr[tr[tn].u0].d0 = tn;
1362 tr[tr[tn].u1].d0 = tn;
1363 }
1364 else /* intersects in the right */
1365 {
1366 tr[tn].u1 = -1;
1367 tr[tn].u0 = tr[t].u1;
1368 tr[t].u1 = tr[t].u0;
1369 tr[t].u0 = tr[t].usave;
1370
1371 tr[tr[t].u0].d0 = t;
1372 tr[tr[t].u1].d0 = t;
1373 tr[tr[tn].u0].d0 = tn;
1374 }
1375
1376 tr[tn].usave = 0;
1377 tr[t].usave = 0;
1378 }
1379 else /* No usave.... simple case */
1380 {
1381 tr[tn].u0 = tr[t].u1;
1382 tr[tn].u1 = -1;
1383 tr[t].u1 = -1;
1384 tr[tr[tn].u0].d0 = tn;
1385 }
1386 }
1387 else
1388 { /* fresh seg. or upward cusp */
1389 int tmp_u = tr[t].u0;
1390 int td0, td1;
1391 if (((td0 = tr[tmp_u].d0) > 0) &&
1392 ((td1 = tr[tmp_u].d1) > 0))
1393 { /* upward cusp */
1394 if ((tr[td0].rseg > 0) &&
1395 !is_left_of(tr[td0].rseg, &s.v1))
1396 {
1397 tr[tn].u1 = -1;
1398 tr[t].u1 = -1;
1399 tr[t].u0 = -1;
1400 tr[tr[tn].u0].d1 = tn;
1401 }
1402 else
1403 {
1404 tr[t].u1 = -1;
1405 tr[tn].u1 = -1;
1406 tr[tn].u0 = -1;
1407 tr[tr[t].u0].d0 = t;
1408 }
1409 }
1410 else /* fresh segment */
1411 {
1412 tr[tr[t].u0].d0 = t;
1413 tr[tr[t].u0].d1 = tn;
1414 }
1415 }
1416
1417 if (FP_EQUAL(tr[t].lo.y, tr[tlast].lo.y) &&
1418 FP_EQUAL(tr[t].lo.x, tr[tlast].lo.x) && tribot)
1419 {
1420 /* this case arises only at the lowest trapezoid.. i.e.
1421 tlast, if the lower endpoint of the segment is
1422 already inserted in the structure */
1423
1424 tr[tr[t].d0].u0 = t;
1425 tr[tr[t].d0].u1 = -1;
1426 tr[tr[t].d1].u0 = tn;
1427 tr[tr[t].d1].u1 = -1;
1428
1429 tr[tn].d0 = tr[t].d1;
1430 tr[tn].d1 = -1;
1431 tr[t].d1 = -1;
1432
1433 tnext = tr[t].d1;
1434 }
1435 else if (i_d0)
1436 /* intersecting d0 */
1437 {
1438 tr[tr[t].d0].u0 = t;
1439 tr[tr[t].d0].u1 = tn;
1440 tr[tr[t].d1].u0 = tn;
1441 tr[tr[t].d1].u1 = -1;
1442
1443 /* new code to determine the bottom neighbours of the */
1444 /* newly partitioned trapezoid */
1445
1446 tr[t].d1 = -1;
1447
1448 tnext = tr[t].d0;
1449 }
1450 else /* intersecting d1 */
1451 {
1452 tr[tr[t].d0].u0 = t;
1453 tr[tr[t].d0].u1 = -1;
1454 tr[tr[t].d1].u0 = t;
1455 tr[tr[t].d1].u1 = tn;
1456
1457 /* new code to determine the bottom neighbours of the */
1458 /* newly partitioned trapezoid */
1459
1460 tr[tn].d0 = tr[t].d1;
1461 tr[tn].d1 = -1;
1462
1463 tnext = tr[t].d1;
1464 }
1465
1466 t = tnext;
1467 }
1468
1469 tr[tn_sav].lseg = segnum;
1470 tr[t_sav].rseg = segnum;
1471 } /* end-while */
1472
1473 /* Now combine those trapezoids which share common segments. We can */
1474 /* use the pointers to the parent to connect these together. This */
1475 /* works only because all these new trapezoids have been formed */
1476 /* due to splitting by the segment, and hence have only one parent */
1477
1478 tfirstl = tfirst;
1479 tlastl = tlast;
1480 merge_trapezoids(segnum, tfirstl, tlastl, S_LEFT);
1481 merge_trapezoids(segnum, tfirstr, tlastr, S_RIGHT);
1482
1483 seg[segnum].is_inserted = true;
1484 return 0;
1485}
1486
1487
1488/* Update the roots stored for each of the endpoints of the segment.
1489 * This is done to speed up the location-query for the endpoint when
1490 * the segment is inserted into the trapezoidation subsequently
1491 */
1492int Triangulator::
1493find_new_roots(int segnum) {
1494 // cerr << "find_new_roots(" << segnum << ")\n";
1495 segment_t *s = &seg[segnum];
1496
1497 if (s->is_inserted)
1498 return 0;
1499
1500 s->root0 = locate_endpoint(&s->v0, &s->v1, s->root0);
1501 s->root0 = tr[s->root0].sink;
1502
1503 s->root1 = locate_endpoint(&s->v1, &s->v0, s->root1);
1504 s->root1 = tr[s->root1].sink;
1505 return 0;
1506}
1507
1508
1509/* Main routine to perform trapezoidation */
1510int Triangulator::
1511construct_trapezoids(int nseg) {
1512 // cerr << "construct_trapezoids(" << nseg << ")\n";
1513 int i;
1514 int root, h;
1515
1516 /* Add the first segment and get the query structure and trapezoid */
1517 /* list initialised */
1518
1519 root = init_query_structure(choose_segment());
1520
1521 for (i = 1; i <= nseg; i++) {
1522 seg[i].root1 = root;
1523 seg[i].root0 = root;
1524 }
1525
1526 for (h = 1; h <= math_logstar_n(nseg); h++)
1527 {
1528 for (i = math_N(nseg, h -1) + 1; i <= math_N(nseg, h); i++) {
1529 if (add_segment(choose_segment()) != 0) {
1530 // error in add_segment.
1531 return 1;
1532 }
1533 }
1534
1535 /* Find a new root for each of the segment endpoints */
1536 for (i = 1; i <= nseg; i++)
1537 find_new_roots(i);
1538 }
1539
1540 for (i = math_N(nseg, math_logstar_n(nseg)) + 1; i <= nseg; i++)
1541 add_segment(choose_segment());
1542
1543 return 0;
1544}
1545
1546
1547
1548/* Function returns true if the trapezoid lies inside the polygon */
1549int Triangulator::
1550inside_polygon(trap_t *t) {
1551 int rseg = t->rseg;
1552
1553 if (t->state == ST_INVALID)
1554 return 0;
1555
1556 if ((t->lseg <= 0) || (t->rseg <= 0))
1557 return 0;
1558
1559 if (((t->u0 <= 0) && (t->u1 <= 0)) ||
1560 ((t->d0 <= 0) && (t->d1 <= 0))) /* triangle */
1561 return (_greater_than(&seg[rseg].v1, &seg[rseg].v0));
1562
1563 return 0;
1564}
1565
1566
1567/* return a new mon structure from the table */
1568int Triangulator::
1569newmon() {
1570 int index = (int)mon.size();
1571 mon.push_back(0);
1572 // cerr << "newmon " << index << "\n";
1573 return index;
1574}
1575
1576
1577/* return a new chain element from the table */
1578int Triangulator::
1579new_chain_element() {
1580 int index = (int)mchain.size();
1581 mchain.push_back(monchain_t());
1582 // cerr << "new_chain_element " << index << "\n";
1583 return index;
1584}
1585
1586
1587double Triangulator::
1588get_angle(point_t *vp0, point_t *vpnext, point_t *vp1) {
1589 point_t v0, v1;
1590
1591 v0.x = vpnext->x - vp0->x;
1592 v0.y = vpnext->y - vp0->y;
1593
1594 v1.x = vp1->x - vp0->x;
1595 v1.y = vp1->y - vp0->y;
1596
1597 if (CROSS_SINE(v0, v1) >= 0) /* sine is positive */
1598 return DOT(v0, v1)/LENGTH(v0)/LENGTH(v1);
1599 else
1600 return (-1.0 * DOT(v0, v1)/LENGTH(v0)/LENGTH(v1) - 2);
1601}
1602
1603
1604/* (v0, v1) is the new diagonal to be added to the polygon. Find which */
1605/* chain to use and return the positions of v0 and v1 in p and q */
1606int Triangulator::
1607get_vertex_positions(int v0, int v1, int *ip, int *iq) {
1608 vertexchain_t *vp0, *vp1;
1609 int i;
1610 double angle, temp;
1611 int tp = 0, tq = 0;
1612
1613 vp0 = &vert[v0];
1614 vp1 = &vert[v1];
1615
1616 /* p is identified as follows. Scan from (v0, v1) rightwards till */
1617 /* you hit the first segment starting from v0. That chain is the */
1618 /* chain of our interest */
1619
1620 angle = -4.0;
1621 for (i = 0; i < 4; i++)
1622 {
1623 if (vp0->vnext[i] <= 0)
1624 continue;
1625 if ((temp = get_angle(&vp0->pt, &(vert[vp0->vnext[i]].pt),
1626 &vp1->pt)) > angle)
1627 {
1628 angle = temp;
1629 tp = i;
1630 }
1631 }
1632
1633 *ip = tp;
1634
1635 /* Do similar actions for q */
1636
1637 angle = -4.0;
1638 for (i = 0; i < 4; i++)
1639 {
1640 if (vp1->vnext[i] <= 0)
1641 continue;
1642 if ((temp = get_angle(&vp1->pt, &(vert[vp1->vnext[i]].pt),
1643 &vp0->pt)) > angle)
1644 {
1645 angle = temp;
1646 tq = i;
1647 }
1648 }
1649
1650 *iq = tq;
1651
1652 return 0;
1653}
1654
1655
1656/* v0 and v1 are specified in anti-clockwise order with respect to
1657 * the current monotone polygon mcur. Split the current polygon into
1658 * two polygons using the diagonal (v0, v1)
1659 */
1660int Triangulator::
1661make_new_monotone_poly(int mcur, int v0, int v1) {
1662 int p, q, ip, iq;
1663 int mnew = newmon();
1664 int i, j, nf0, nf1;
1665 vertexchain_t *vp0, *vp1;
1666
1667 if (v0 <= 0 || v1 <= 0) {
1668 return -1;
1669 }
1670
1671 vp0 = &vert[v0];
1672 vp1 = &vert[v1];
1673
1674 get_vertex_positions(v0, v1, &ip, &iq);
1675
1676 p = vp0->vpos[ip];
1677 q = vp1->vpos[iq];
1678
1679 /* At this stage, we have got the positions of v0 and v1 in the */
1680 /* desired chain. Now modify the linked lists */
1681
1682 i = new_chain_element(); /* for the new list */
1683 j = new_chain_element();
1684
1685 mchain[i].vnum = v0;
1686 mchain[j].vnum = v1;
1687
1688 mchain[i].next = mchain[p].next;
1689 mchain[mchain[p].next].prev = i;
1690 mchain[i].prev = j;
1691 mchain[j].next = i;
1692 mchain[j].prev = mchain[q].prev;
1693 mchain[mchain[q].prev].next = j;
1694
1695 mchain[p].next = q;
1696 mchain[q].prev = p;
1697
1698 nf0 = vp0->nextfree;
1699 nf1 = vp1->nextfree;
1700
1701 vp0->vnext[ip] = v1;
1702
1703 vp0->vpos[nf0] = i;
1704 vp0->vnext[nf0] = mchain[mchain[i].next].vnum;
1705 vp1->vpos[nf1] = j;
1706 vp1->vnext[nf1] = v0;
1707
1708 vp0->nextfree++;
1709 vp1->nextfree++;
1710
1711#ifdef DEBUG
1712 fprintf(stderr, "make_poly: mcur = %d, (v0, v1) = (%d, %d)\n",
1713 mcur, v0, v1);
1714 fprintf(stderr, "next posns = (p, q) = (%d, %d)\n", p, q);
1715#endif
1716
1717 mon[mcur] = p;
1718 mon[mnew] = i;
1719 return mnew;
1720}
1721
1722/* Main routine to get monotone polygons from the trapezoidation of
1723 * the polygon.
1724 */
1725
1726int Triangulator::
1727monotonate_trapezoids(int n) {
1728 int i;
1729 int tr_start;
1730
1731 vert.clear();
1732 visited.clear();
1733 mchain.clear();
1734 mon.clear();
1735
1736 vert.insert(vert.begin(), n + 1, vertexchain_t());
1737 mchain.insert(mchain.begin(), n + 1, monchain_t());
1738
1739 visited.insert(visited.begin(), tr.size(), 0);
1740
1741 /* First locate a trapezoid which lies inside the polygon */
1742 /* and which is triangular */
1743 for (i = 1; i < (int)tr.size(); i++)
1744 if (inside_polygon(&tr[i]))
1745 break;
1746 if (i >= (int)tr.size()) {
1747 // No valid trapezoids.
1748 return 0;
1749 }
1750 // printf("start = %d\n", i);
1751 tr_start = i;
1752
1753 /* Initialise the mon data-structure and start spanning all the */
1754 /* trapezoids within the polygon */
1755
1756#if 0
1757 for (i = 1; i <= n; i++)
1758 {
1759 mchain[i].prev = i - 1;
1760 mchain[i].next = i + 1;
1761 mchain[i].vnum = i;
1762 vert[i].pt = seg[i].v0;
1763 vert[i].vnext[0] = i + 1; /* next vertex */
1764 vert[i].vpos[0] = i; /* locn. of next vertex */
1765 vert[i].nextfree = 1;
1766 vert[i].user_i = seg[i].v0_i;
1767 }
1768 mchain[1].prev = n;
1769 mchain[n].next = 1;
1770 vert[n].vnext[0] = 1;
1771 vert[n].vpos[0] = n;
1772 mon.push_back(1); /* position of any vertex in the first */
1773 /* chain */
1774
1775#else
1776
1777 for (i = 1; i <= n; i++)
1778 {
1779 mchain[i].prev = seg[i].prev;
1780 mchain[i].next = seg[i].next;
1781 mchain[i].vnum = i;
1782 vert[i].pt = seg[i].v0;
1783 vert[i].vnext[0] = seg[i].next; /* next vertex */
1784 vert[i].vpos[0] = i; /* locn. of next vertex */
1785 vert[i].nextfree = 1;
1786 vert[i].user_i = seg[i].v0_i;
1787 }
1788
1789 mon.push_back(1); /* position of any vertex in the first */
1790 /* chain */
1791
1792#endif
1793
1794 /* traverse the polygon */
1795 if (tr[tr_start].u0 > 0)
1796 traverse_polygon(0, tr_start, tr[tr_start].u0, TR_FROM_UP);
1797 else if (tr[tr_start].d0 > 0)
1798 traverse_polygon(0, tr_start, tr[tr_start].d0, TR_FROM_DN);
1799
1800 /* return the number of polygons created */
1801 return newmon();
1802}
1803
1804
1805/* recursively visit all the trapezoids */
1806int Triangulator::
1807traverse_polygon(int mcur, int trnum, int from, int dir) {
1808 // printf("traverse_polygon(%d, %d, %d, %d)\n", mcur, trnum, from, dir);
1809
1810 if (mcur < 0 || trnum <= 0)
1811 return 0;
1812
1813 if (visited[trnum])
1814 return 0;
1815
1816 trap_t *t = &tr[trnum];
1817 // int howsplit;
1818 int mnew;
1819 int v0, v1; //, v0next, v1next;
1820 int retval = 0; //, tmp;
1821
1822 // printf("visited size = %d, visited[trnum] = %d\n", visited.size(),
1823 // visited[trnum]);
1824
1825 visited[trnum] = true;
1826
1827 /* We have much more information available here. */
1828 /* rseg: goes upwards */
1829 /* lseg: goes downwards */
1830
1831 /* Initially assume that dir = TR_FROM_DN (from the left) */
1832 /* Switch v0 and v1 if necessary afterwards */
1833
1834
1835 /* special cases for triangles with cusps at the opposite ends. */
1836 /* take care of this first */
1837 if ((t->u0 <= 0) && (t->u1 <= 0))
1838 {
1839 if ((t->d0 > 0) && (t->d1 > 0)) /* downward opening triangle */
1840 {
1841 v0 = tr[t->d1].lseg;
1842 v1 = t->lseg;
1843 if (from == t->d1)
1844 {
1845 mnew = make_new_monotone_poly(mcur, v1, v0);
1846 traverse_polygon(mcur, t->d1, trnum, TR_FROM_UP);
1847 traverse_polygon(mnew, t->d0, trnum, TR_FROM_UP);
1848 }
1849 else
1850 {
1851 mnew = make_new_monotone_poly(mcur, v0, v1);
1852 traverse_polygon(mcur, t->d0, trnum, TR_FROM_UP);
1853 traverse_polygon(mnew, t->d1, trnum, TR_FROM_UP);
1854 }
1855 }
1856 else
1857 {
1858 retval = SP_NOSPLIT; /* Just traverse all neighbours */
1859 traverse_polygon(mcur, t->u0, trnum, TR_FROM_DN);
1860 traverse_polygon(mcur, t->u1, trnum, TR_FROM_DN);
1861 traverse_polygon(mcur, t->d0, trnum, TR_FROM_UP);
1862 traverse_polygon(mcur, t->d1, trnum, TR_FROM_UP);
1863 }
1864 }
1865
1866 else if ((t->d0 <= 0) && (t->d1 <= 0))
1867 {
1868 if ((t->u0 > 0) && (t->u1 > 0)) /* upward opening triangle */
1869 {
1870 v0 = t->rseg;
1871 v1 = tr[t->u0].rseg;
1872 if (from == t->u1)
1873 {
1874 mnew = make_new_monotone_poly(mcur, v1, v0);
1875 traverse_polygon(mcur, t->u1, trnum, TR_FROM_DN);
1876 traverse_polygon(mnew, t->u0, trnum, TR_FROM_DN);
1877 }
1878 else
1879 {
1880 mnew = make_new_monotone_poly(mcur, v0, v1);
1881 traverse_polygon(mcur, t->u0, trnum, TR_FROM_DN);
1882 traverse_polygon(mnew, t->u1, trnum, TR_FROM_DN);
1883 }
1884 }
1885 else
1886 {
1887 retval = SP_NOSPLIT; /* Just traverse all neighbours */
1888 traverse_polygon(mcur, t->u0, trnum, TR_FROM_DN);
1889 traverse_polygon(mcur, t->u1, trnum, TR_FROM_DN);
1890 traverse_polygon(mcur, t->d0, trnum, TR_FROM_UP);
1891 traverse_polygon(mcur, t->d1, trnum, TR_FROM_UP);
1892 }
1893 }
1894
1895 else if ((t->u0 > 0) && (t->u1 > 0))
1896 {
1897 if ((t->d0 > 0) && (t->d1 > 0)) /* downward + upward cusps */
1898 {
1899 v0 = tr[t->d1].lseg;
1900 v1 = tr[t->u0].rseg;
1901 retval = SP_2UP_2DN;
1902 if (((dir == TR_FROM_DN) && (t->d1 == from)) ||
1903 ((dir == TR_FROM_UP) && (t->u1 == from)))
1904 {
1905 mnew = make_new_monotone_poly(mcur, v1, v0);
1906 traverse_polygon(mcur, t->u1, trnum, TR_FROM_DN);
1907 traverse_polygon(mcur, t->d1, trnum, TR_FROM_UP);
1908 traverse_polygon(mnew, t->u0, trnum, TR_FROM_DN);
1909 traverse_polygon(mnew, t->d0, trnum, TR_FROM_UP);
1910 }
1911 else
1912 {
1913 mnew = make_new_monotone_poly(mcur, v0, v1);
1914 traverse_polygon(mcur, t->u0, trnum, TR_FROM_DN);
1915 traverse_polygon(mcur, t->d0, trnum, TR_FROM_UP);
1916 traverse_polygon(mnew, t->u1, trnum, TR_FROM_DN);
1917 traverse_polygon(mnew, t->d1, trnum, TR_FROM_UP);
1918 }
1919 }
1920 else /* only downward cusp */
1921 {
1922 if (_equal_to(&t->lo, &seg[t->lseg].v1))
1923 {
1924 v0 = tr[t->u0].rseg;
1925 v1 = seg[t->lseg].next;
1926
1927 retval = SP_2UP_LEFT;
1928 if ((dir == TR_FROM_UP) && (t->u0 == from))
1929 {
1930 mnew = make_new_monotone_poly(mcur, v1, v0);
1931 traverse_polygon(mcur, t->u0, trnum, TR_FROM_DN);
1932 traverse_polygon(mnew, t->d0, trnum, TR_FROM_UP);
1933 traverse_polygon(mnew, t->u1, trnum, TR_FROM_DN);
1934 traverse_polygon(mnew, t->d1, trnum, TR_FROM_UP);
1935 }
1936 else
1937 {
1938 mnew = make_new_monotone_poly(mcur, v0, v1);
1939 traverse_polygon(mcur, t->u1, trnum, TR_FROM_DN);
1940 traverse_polygon(mcur, t->d0, trnum, TR_FROM_UP);
1941 traverse_polygon(mcur, t->d1, trnum, TR_FROM_UP);
1942 traverse_polygon(mnew, t->u0, trnum, TR_FROM_DN);
1943 }
1944 }
1945 else
1946 {
1947 v0 = t->rseg;
1948 v1 = tr[t->u0].rseg;
1949 retval = SP_2UP_RIGHT;
1950 if ((dir == TR_FROM_UP) && (t->u1 == from))
1951 {
1952 mnew = make_new_monotone_poly(mcur, v1, v0);
1953 traverse_polygon(mcur, t->u1, trnum, TR_FROM_DN);
1954 traverse_polygon(mnew, t->d1, trnum, TR_FROM_UP);
1955 traverse_polygon(mnew, t->d0, trnum, TR_FROM_UP);
1956 traverse_polygon(mnew, t->u0, trnum, TR_FROM_DN);
1957 }
1958 else
1959 {
1960 mnew = make_new_monotone_poly(mcur, v0, v1);
1961 traverse_polygon(mcur, t->u0, trnum, TR_FROM_DN);
1962 traverse_polygon(mcur, t->d0, trnum, TR_FROM_UP);
1963 traverse_polygon(mcur, t->d1, trnum, TR_FROM_UP);
1964 traverse_polygon(mnew, t->u1, trnum, TR_FROM_DN);
1965 }
1966 }
1967 }
1968 }
1969 else if ((t->u0 > 0) || (t->u1 > 0)) /* no downward cusp */
1970 {
1971 if ((t->d0 > 0) && (t->d1 > 0)) /* only upward cusp */
1972 {
1973 if (_equal_to(&t->hi, &seg[t->lseg].v0))
1974 {
1975 v0 = tr[t->d1].lseg;
1976 v1 = t->lseg;
1977 retval = SP_2DN_LEFT;
1978 if (!((dir == TR_FROM_DN) && (t->d0 == from)))
1979 {
1980 mnew = make_new_monotone_poly(mcur, v1, v0);
1981 traverse_polygon(mcur, t->u1, trnum, TR_FROM_DN);
1982 traverse_polygon(mcur, t->d1, trnum, TR_FROM_UP);
1983 traverse_polygon(mcur, t->u0, trnum, TR_FROM_DN);
1984 traverse_polygon(mnew, t->d0, trnum, TR_FROM_UP);
1985 }
1986 else
1987 {
1988 mnew = make_new_monotone_poly(mcur, v0, v1);
1989 traverse_polygon(mcur, t->d0, trnum, TR_FROM_UP);
1990 traverse_polygon(mnew, t->u0, trnum, TR_FROM_DN);
1991 traverse_polygon(mnew, t->u1, trnum, TR_FROM_DN);
1992 traverse_polygon(mnew, t->d1, trnum, TR_FROM_UP);
1993 }
1994 }
1995 else
1996 {
1997 v0 = tr[t->d1].lseg;
1998 v1 = seg[t->rseg].next;
1999
2000 retval = SP_2DN_RIGHT;
2001 if ((dir == TR_FROM_DN) && (t->d1 == from))
2002 {
2003 mnew = make_new_monotone_poly(mcur, v1, v0);
2004 traverse_polygon(mcur, t->d1, trnum, TR_FROM_UP);
2005 traverse_polygon(mnew, t->u1, trnum, TR_FROM_DN);
2006 traverse_polygon(mnew, t->u0, trnum, TR_FROM_DN);
2007 traverse_polygon(mnew, t->d0, trnum, TR_FROM_UP);
2008 }
2009 else
2010 {
2011 mnew = make_new_monotone_poly(mcur, v0, v1);
2012 traverse_polygon(mcur, t->u0, trnum, TR_FROM_DN);
2013 traverse_polygon(mcur, t->d0, trnum, TR_FROM_UP);
2014 traverse_polygon(mcur, t->u1, trnum, TR_FROM_DN);
2015 traverse_polygon(mnew, t->d1, trnum, TR_FROM_UP);
2016 }
2017 }
2018 }
2019 else /* no cusp */
2020 {
2021 if (_equal_to(&t->hi, &seg[t->lseg].v0) &&
2022 _equal_to(&t->lo, &seg[t->rseg].v0))
2023 {
2024 v0 = t->rseg;
2025 v1 = t->lseg;
2026 retval = SP_SIMPLE_LRDN;
2027 if (dir == TR_FROM_UP)
2028 {
2029 mnew = make_new_monotone_poly(mcur, v1, v0);
2030 traverse_polygon(mcur, t->u0, trnum, TR_FROM_DN);
2031 traverse_polygon(mcur, t->u1, trnum, TR_FROM_DN);
2032 traverse_polygon(mnew, t->d1, trnum, TR_FROM_UP);
2033 traverse_polygon(mnew, t->d0, trnum, TR_FROM_UP);
2034 }
2035 else
2036 {
2037 mnew = make_new_monotone_poly(mcur, v0, v1);
2038 traverse_polygon(mcur, t->d1, trnum, TR_FROM_UP);
2039 traverse_polygon(mcur, t->d0, trnum, TR_FROM_UP);
2040 traverse_polygon(mnew, t->u0, trnum, TR_FROM_DN);
2041 traverse_polygon(mnew, t->u1, trnum, TR_FROM_DN);
2042 }
2043 }
2044 else if (_equal_to(&t->hi, &seg[t->rseg].v1) &&
2045 _equal_to(&t->lo, &seg[t->lseg].v1))
2046 {
2047 v0 = seg[t->rseg].next;
2048 v1 = seg[t->lseg].next;
2049
2050 retval = SP_SIMPLE_LRUP;
2051 if (dir == TR_FROM_UP)
2052 {
2053 mnew = make_new_monotone_poly(mcur, v1, v0);
2054 traverse_polygon(mcur, t->u0, trnum, TR_FROM_DN);
2055 traverse_polygon(mcur, t->u1, trnum, TR_FROM_DN);
2056 traverse_polygon(mnew, t->d1, trnum, TR_FROM_UP);
2057 traverse_polygon(mnew, t->d0, trnum, TR_FROM_UP);
2058 }
2059 else
2060 {
2061 mnew = make_new_monotone_poly(mcur, v0, v1);
2062 traverse_polygon(mcur, t->d1, trnum, TR_FROM_UP);
2063 traverse_polygon(mcur, t->d0, trnum, TR_FROM_UP);
2064 traverse_polygon(mnew, t->u0, trnum, TR_FROM_DN);
2065 traverse_polygon(mnew, t->u1, trnum, TR_FROM_DN);
2066 }
2067 }
2068 else /* no split possible */
2069 {
2070 retval = SP_NOSPLIT;
2071 traverse_polygon(mcur, t->u0, trnum, TR_FROM_DN);
2072 traverse_polygon(mcur, t->d0, trnum, TR_FROM_UP);
2073 traverse_polygon(mcur, t->u1, trnum, TR_FROM_DN);
2074 traverse_polygon(mcur, t->d1, trnum, TR_FROM_UP);
2075 }
2076 }
2077 }
2078
2079 return retval;
2080}
2081
2082
2083/* For each monotone polygon, find the ymax and ymin (to determine the */
2084/* two y-monotone chains) and pass on this monotone polygon for greedy */
2085/* triangulation. */
2086/* Take care not to triangulate duplicate monotone polygons */
2087
2088void Triangulator::
2089triangulate_monotone_polygons(int nvert, int nmonpoly) {
2090 int i;
2091 point_t ymax, ymin;
2092 int p, vfirst, posmax, posmin, v;
2093 int vcount, processed;
2094
2095 #ifdef DEBUG
2096 for (i = 0; i < nmonpoly; i++)
2097 {
2098 fprintf(stderr, "\n\nPolygon %d: ", i);
2099 vfirst = mchain[mon[i]].vnum;
2100 p = mchain[mon[i]].next;
2101 fprintf (stderr, "%d ", mchain[mon[i]].vnum);
2102 while (mchain[p].vnum != vfirst)
2103 {
2104 fprintf(stderr, "%d ", mchain[p].vnum);
2105 if (mchain[p].vnum == -1) {
2106 fprintf(stderr, " xx");
2107 break;
2108 }
2109 p = mchain[p].next;
2110 }
2111 }
2112 fprintf(stderr, "\n");
2113 #endif
2114
2115 for (i = 0; i < nmonpoly; i++)
2116 {
2117 vcount = 1;
2118 processed = false;
2119 vfirst = mchain[mon[i]].vnum;
2120 if (vfirst <= 0) {
2121 return;
2122 }
2123 ymin = vert[vfirst].pt;
2124 ymax = ymin;
2125 posmin = mon[i];
2126 posmax = posmin;
2127 mchain[mon[i]].marked = true;
2128 p = mchain[mon[i]].next;
2129 while ((v = mchain[p].vnum) != vfirst)
2130 {
2131 if (v <= 0) {
2132 return;
2133 }
2134 if (mchain[p].marked)
2135 {
2136 processed = true;
2137 break; /* break from while */
2138 }
2139 else
2140 mchain[p].marked = true;
2141
2142 if (_greater_than(&vert[v].pt, &ymax))
2143 {
2144 ymax = vert[v].pt;
2145 posmax = p;
2146 }
2147 if (_less_than(&vert[v].pt, &ymin))
2148 {
2149 ymin = vert[v].pt;
2150 posmin = p;
2151 }
2152 p = mchain[p].next;
2153 vcount++;
2154 }
2155
2156 if (processed) /* Go to next polygon */
2157 continue;
2158
2159 if (vcount == 3) /* already a triangle */
2160 {
2161 _result.push_back(Triangle(this, mchain[p].vnum,
2162 mchain[mchain[p].next].vnum,
2163 mchain[mchain[p].prev].vnum));
2164 }
2165 else /* triangulate the polygon */
2166 {
2167 v = mchain[mchain[posmax].next].vnum;
2168 if (_equal_to(&vert[v].pt, &ymin))
2169 { /* LHS is a single line */
2170 triangulate_single_polygon(nvert, posmax, TRI_LHS);
2171 }
2172 else
2173 triangulate_single_polygon(nvert, posmax, TRI_RHS);
2174 }
2175 }
2176}
2177
2178
2179/* A greedy corner-cutting algorithm to triangulate a y-monotone
2180 * polygon in O(n) time.
2181 * Joseph O-Rourke, Computational Geometry in C.
2182 */
2183void Triangulator::
2184triangulate_single_polygon(int nvert, int posmax, int side) {
2185 int v;
2186 vector_int rc; /* reflex chain */
2187 int ri;
2188 int endv, tmp, vpos;
2189
2190 // cerr << "triangulate_single_polygon(" << nvert << ", " << posmax << ", "
2191 // << side << ")\n";
2192
2193 if (side == TRI_RHS) /* RHS segment is a single segment */
2194 {
2195 rc.push_back(mchain[posmax].vnum);
2196 tmp = mchain[posmax].next;
2197 rc.push_back(mchain[tmp].vnum);
2198 ri = 1;
2199
2200 vpos = mchain[tmp].next;
2201 v = mchain[vpos].vnum;
2202
2203 if ((endv = mchain[mchain[posmax].prev].vnum) == 0)
2204 endv = nvert;
2205 }
2206 else /* LHS is a single segment */
2207 {
2208 tmp = mchain[posmax].next;
2209 rc.push_back(mchain[tmp].vnum);
2210 tmp = mchain[tmp].next;
2211 rc.push_back(mchain[tmp].vnum);
2212 ri = 1;
2213
2214 vpos = mchain[tmp].next;
2215 v = mchain[vpos].vnum;
2216
2217 endv = mchain[posmax].vnum;
2218 }
2219
2220 int num_triangles = 0;
2221
2222 while ((v != endv) || (ri > 1))
2223 {
2224 // cerr << " v = " << v << " ri = " << ri << " rc = " << rc.size() << "
2225 // _result = " << _result.size() << "\n";
2226 if (v <= 0) {
2227 // Something went wrong.
2228 return;
2229 }
2230 if (ri > 0) /* reflex chain is non-empty */
2231 {
2232 PN_stdfloat crossResult = CROSS(vert[v].pt, vert[rc[ri - 1]].pt,
2233 vert[rc[ri]].pt);
2234 if ( crossResult >= 0 ) /* could be convex corner or straight */
2235 {
2236 if (crossResult > 0) { /* convex corner: cut it off */
2237 _result.push_back(Triangle(this, rc[ri - 1], rc[ri], v));
2238 if (++num_triangles >= nvert - 2) {
2239 // We can't generate more than this number of triangles.
2240 return;
2241 }
2242 }
2243 /* else : perfectly straight, will be abandoned anyway */
2244 ri--;
2245 rc.pop_back();
2246 nassertv(ri + 1 == (int)rc.size());
2247 }
2248 else /* concave, add v to the chain */
2249 {
2250 ri++;
2251 rc.push_back(v);
2252 nassertv(ri + 1 == (int)rc.size());
2253 vpos = mchain[vpos].next;
2254 v = mchain[vpos].vnum;
2255 }
2256 }
2257 else /* reflex-chain empty: add v to the */
2258 { /* reflex chain and advance it */
2259 ri++;
2260 rc.push_back(v);
2261 nassertv(ri + 1 == (int)rc.size());
2262 vpos = mchain[vpos].next;
2263 v = mchain[vpos].vnum;
2264 }
2265 } /* end-while */
2266
2267 /* reached the bottom vertex. Add in the triangle formed */
2268 _result.push_back(Triangle(this, rc[ri - 1], rc[ri], v));
2269 ri--;
2270}
A handy class to return random numbers.
Definition randomizer.h:26
int random_int(int range)
Returns a random integer in the range [0, range).
Definition randomizer.I:45
int get_triangle_v0(int n) const
Returns vertex 0 of the nth triangle generated by the previous call to triangulate().
void clear()
Removes all vertices and polygon specifications from the Triangulator, and prepares it to start over.
void add_polygon_vertex(int index)
Adds the next consecutive vertex of the polygon.
void clear_polygon()
Removes the current polygon definition (and its set of holes), but does not clear the vertex pool.
int get_triangle_v1(int n) const
Returns vertex 1 of the nth triangle generated by the previous call to triangulate().
int get_num_triangles() const
Returns the number of triangles generated by the previous call to triangulate().
int add_vertex(const LPoint2d &point)
Adds a new vertex to the vertex pool.
void add_hole_vertex(int index)
Adds the next consecutive vertex of the current hole.
void triangulate()
Does the work of triangulating the specified polygon.
void begin_hole()
Finishes the previous hole, if any, and prepares to add a new hole.
int get_triangle_v2(int n) const
Returns vertex 2 of the nth triangle generated by the previous call to triangulate().
PANDA 3D SOFTWARE Copyright (c) Carnegie Mellon University.
PANDA 3D SOFTWARE Copyright (c) Carnegie Mellon University.