FreeWRL / FreeX3D 4.3.0
monoChain.cc
1/*
2** License Applicability. Except to the extent portions of this file are
3** made subject to an alternative license as permitted in the SGI Free
4** Software License B, Version 1.1 (the "License"), the contents of this
5** file are subject only to the provisions of the License. You may not use
6** this file except in compliance with the License. You may obtain a copy
7** of the License at Silicon Graphics, Inc., attn: Legal Services, 1600
8** Amphitheatre Parkway, Mountain View, CA 94043-1351, or at:
9**
10** http://oss.sgi.com/projects/FreeB
11**
12** Note that, as provided in the License, the Software is distributed on an
13** "AS IS" basis, with ALL EXPRESS AND IMPLIED WARRANTIES AND CONDITIONS
14** DISCLAIMED, INCLUDING, WITHOUT LIMITATION, ANY IMPLIED WARRANTIES AND
15** CONDITIONS OF MERCHANTABILITY, SATISFACTORY QUALITY, FITNESS FOR A
16** PARTICULAR PURPOSE, AND NON-INFRINGEMENT.
17**
18** Original Code. The Original Code is: OpenGL Sample Implementation,
19** Version 1.2.1, released January 26, 2000, developed by Silicon Graphics,
20** Inc. The Original Code is Copyright (c) 1991-2000 Silicon Graphics, Inc.
21** Copyright in any portions created by third parties is as indicated
22** elsewhere herein. All Rights Reserved.
23**
24** Additional Notice Provisions: The application programming interfaces
25** established by SGI in conjunction with the Original Code are The
26** OpenGL(R) Graphics System: A Specification (Version 1.2.1), released
27** April 1, 1999; The OpenGL(R) Graphics System Utility Library (Version
28** 1.3), released November 4, 1998; and OpenGL(R) Graphics with the X
29** Window System(R) (Version 1.3), released October 19, 1998. This software
30** was created using the OpenGL(R) version 1.2.1 Sample Implementation
31** published by SGI, but has not been independently verified as being
32** compliant with the OpenGL(R) version 1.2.1 Specification.
33**
34*/
35/*
36*/
37
38#include "gluos.h"
39#include <stdlib.h>
40#include <stdio.h>
41#include <libnurbs2.h>
42//#include <GL/gl.h>
43
44#include "glimports.h"
45#include "zlassert.h"
46
47#include "monoChain.h"
48#include "quicksort.h"
49#include "searchTree.h"
50#include "polyUtil.h"
51
52#ifndef max
53#define max(a,b) ((a>b)? a:b)
54#endif
55#ifndef min
56#define min(a,b) ((a>b)? b:a)
57#endif
58
59extern Int isCusp(directedLine *v);
60extern Int deleteRepeatDiagonals(Int num_diagonals, directedLine** diagonal_vertices, directedLine** new_vertices);
61
62//for debug purpose only
63#if 0 // UNUSED
64static void drawDiagonals(Int num_diagonals, directedLine** diagonal_vertices)
65{
66 Int i;
67 for(i=0; i<num_diagonals; i++)
68 {
69 glBegin(GL_LINE);
70 glVertex2fv(diagonal_vertices[2*i]->head());
71 glVertex2fv(diagonal_vertices[2*i+1]->head());
72 glEnd();
73 }
74}
75#endif
76
77/*given (x_1, y_1) and (x_2, y_2), and y
78 *return x such that (x,y) is on the line
79 */
80inline Real intersectHoriz(Real x1, Real y1, Real x2, Real y2, Real y)
81{
82 return ((y2==y1)? (x1+x2)*0.5 : x1 + ((y-y1)/(y2-y1)) * (x2-x1));
83}
84
85//compare the heads of the two chains
86static int compChainHeadInY(monoChain* mc1, monoChain* mc2)
87{
88 return compV2InY(mc1->getHead()->head(), mc2->getHead()->head());
89}
90
91monoChain::monoChain(directedLine* cHead, directedLine* cTail)
92{
93 chainHead = cHead;
94 chainTail = cTail;
95 next = this;
96 prev = this;
97
98 nextPolygon = NULL;
99
100 //compute bounding box
101 directedLine* temp;
102 minX = maxX = chainTail->head()[0];
103 minY = maxY = chainTail->head()[1];
104
105 for(temp=chainHead; temp!=cTail; temp = temp->getNext())
106 {
107 if(temp->head()[0] < minX)
108 minX = temp->head()[0];
109 if(temp->head()[0] > maxX)
110 maxX = temp->head()[0];
111
112 if(temp->head()[1] < minY)
113 minY = temp->head()[1];
114 if(temp->head()[1] > maxY)
115 maxY = temp->head()[1];
116 }
117
118 //check whether the chain is increasing or decreasing
119 if(chainHead->compInY(chainTail) <0)
120 isIncrease = 1;
121 else
122 isIncrease = 0;
123
124 //initilize currrent, this is used for accelerating search
125 if(isIncrease)
126 current = chainHead;
127 else
128 current = chainTail;
129
130 isKey = 0;
131}
132
133//insert a new line between prev and this
134void monoChain::insert(monoChain* nc)
135{
136 nc->next = this;
137 nc->prev = prev;
138 prev->next = nc;
139 prev = nc;
140}
141
142void monoChain::deleteLoop()
143{
144 monoChain *temp, *tempNext;
145 prev->next = NULL;
146 for(temp=this; temp != NULL; temp = tempNext)
147 {
148 tempNext = temp->next;
149 delete temp;
150 }
151}
152
153void monoChain::deleteLoopList()
154{
155 monoChain *temp, *tempNext;
156 for(temp=this; temp != NULL; temp = tempNext)
157 {
158 tempNext = temp->nextPolygon;
159 temp->deleteLoop();
160 }
161}
162
163Int monoChain::toArraySingleLoop(monoChain** array, Int index)
164{
165 monoChain *temp;
166 array[index++] = this;
167 for(temp = next; temp != this; temp = temp->next)
168 {
169 array[index++] = temp;
170 }
171 return index;
172}
173
174monoChain** monoChain::toArrayAllLoops(Int& num_chains)
175{
176 num_chains = numChainsAllLoops();
177 monoChain **ret = (monoChain**) malloc(sizeof(monoChain*) * num_chains);
178 assert(ret);
179 monoChain *temp;
180 Int index = 0;
181 for(temp = this; temp != NULL; temp=temp->nextPolygon){
182 index = temp->toArraySingleLoop(ret, index);
183 }
184 return ret;
185}
186
187Int monoChain::numChainsSingleLoop()
188{
189 Int ret=0;
190 monoChain* temp;
191 if(next == this) return 1;
192 ret = 1;
193 for(temp=next; temp != this; temp = temp->next)
194 ret++;
195 return ret;
196}
197
198Int monoChain::numChainsAllLoops()
199{
200 Int ret=0;
201 monoChain *temp;
202 for(temp =this; temp != NULL; temp = temp->nextPolygon)
203 ret += temp->numChainsSingleLoop();
204 return ret;
205}
206
207//update 'current'
208Real monoChain::chainIntersectHoriz(Real y)
209{
210 directedLine* temp;
211 if(isIncrease)
212 {
213 for(temp= current; temp != chainTail; temp = temp->getNext())
214 {
215 if(temp->head()[1] > y)
216 break;
217 }
218 current = temp->getPrev();
219 }
220 else
221 {
222 for(temp = current; temp != chainHead; temp = temp->getPrev())
223 {
224 if(temp->head()[1] > y)
225 break;
226 }
227 current = temp->getNext();
228 }
229 return intersectHoriz(current->head()[0], current->head()[1], current->tail()[0], current->tail()[1], y);
230}
231
232monoChain* directedLineLoopToMonoChainLoop(directedLine* loop)
233{
234 directedLine *temp;
235 monoChain *ret=NULL;
236
237 //find the first cusp
238 directedLine *prevCusp=NULL;
239 directedLine *firstCusp;
240
241 if(isCusp(loop))
242 prevCusp = loop;
243 else
244 {
245 for(temp = loop->getNext(); temp != loop; temp = temp->getNext())
246 if(isCusp(temp))
247 break;
248 prevCusp = temp;
249 }
250 firstCusp = prevCusp;
251//printf("first cusp is (%f,%f), (%f,%f), (%f,%f)\n", prevCusp->getPrev()->head()[0], prevCusp->getPrev()->head()[1], prevCusp->head()[0], prevCusp->head()[1], prevCusp->tail()[0], prevCusp->tail()[1]);
252
253 for(temp = prevCusp->getNext(); temp != loop; temp = temp->getNext())
254 {
255 if(isCusp(temp))
256 {
257//printf("the cusp is (%f,%f), (%f,%f), (%f,%f)\n", temp->getPrev()->head()[0], temp->getPrev()->head()[1], temp->head()[0], temp->head()[1], temp->tail()[0], temp->tail()[1]);
258 if(ret == NULL)
259 {
260 ret = new monoChain(prevCusp, temp);
261 }
262 else
263 ret->insert(new monoChain(prevCusp, temp));
264 prevCusp = temp;
265 }
266 }
267 ret->insert(new monoChain(prevCusp, firstCusp));
268
269 return ret;
270}
271
272monoChain* directedLineLoopListToMonoChainLoopList(directedLine* list)
273{
274 directedLine* temp;
275 monoChain* mc;
276 monoChain* mcEnd;
277 mc = directedLineLoopToMonoChainLoop(list);
278 mcEnd = mc;
279 for(temp = list->getNextPolygon(); temp != NULL; temp = temp->getNextPolygon())
280 {
281 monoChain *newLoop = directedLineLoopToMonoChainLoop(temp);
282 mcEnd->setNextPolygon(newLoop);
283 mcEnd = newLoop;
284 }
285 return mc;
286}
287
288/*compare two edges of a polygon.
289 *edge A < edge B if there is a horizontal line so that the intersection
290 *with A is to the left of the intersection with B.
291 *This function is used in sweepY for the dynamic search tree insertion to
292 *order the edges.
293 * Implementation: (x_1,y_1) and (x_2, y_2)
294 */
295static Int compEdges(directedLine *e1, directedLine *e2)
296{
297 Real* head1 = e1->head();
298 Real* tail1 = e1->tail();
299 Real* head2 = e2->head();
300 Real* tail2 = e2->tail();
301/*
302 Real h10 = head1[0];
303 Real h11 = head1[1];
304 Real t10 = tail1[0];
305 Real t11 = tail1[1];
306 Real h20 = head2[0];
307 Real h21 = head2[1];
308 Real t20 = tail2[0];
309 Real t21 = tail2[1];
310*/
311 Real e1_Ymax, e1_Ymin, e2_Ymax, e2_Ymin;
312/*
313 if(h11>t11) {
314 e1_Ymax= h11;
315 e1_Ymin= t11;
316 }
317 else{
318 e1_Ymax = t11;
319 e1_Ymin = h11;
320 }
321
322 if(h21>t21) {
323 e2_Ymax= h21;
324 e2_Ymin= t21;
325 }
326 else{
327 e2_Ymax = t21;
328 e2_Ymin = h21;
329 }
330*/
331
332 if(head1[1]>tail1[1]) {
333 e1_Ymax= head1[1];
334 e1_Ymin= tail1[1];
335 }
336 else{
337 e1_Ymax = tail1[1];
338 e1_Ymin = head1[1];
339 }
340
341 if(head2[1]>tail2[1]) {
342 e2_Ymax= head2[1];
343 e2_Ymin= tail2[1];
344 }
345 else{
346 e2_Ymax = tail2[1];
347 e2_Ymin = head2[1];
348 }
349
350
351 /*Real e1_Ymax = max(head1[1], tail1[1]);*/ /*max(e1->head()[1], e1->tail()[1]);*/
352 /*Real e1_Ymin = min(head1[1], tail1[1]);*/ /*min(e1->head()[1], e1->tail()[1]);*/
353 /*Real e2_Ymax = max(head2[1], tail2[1]);*/ /*max(e2->head()[1], e2->tail()[1]);*/
354 /*Real e2_Ymin = min(head2[1], tail2[1]);*/ /*min(e2->head()[1], e2->tail()[1]);*/
355
356 Real Ymax = min(e1_Ymax, e2_Ymax);
357 Real Ymin = max(e1_Ymin, e2_Ymin);
358
359 Real y = 0.5*(Ymax + Ymin);
360
361/* Real x1 = intersectHoriz(e1->head()[0], e1->head()[1], e1->tail()[0], e1->tail()[1], y);
362 Real x2 = intersectHoriz(e2->head()[0], e2->head()[1], e2->tail()[0], e2->tail()[1], y);
363*/
364/*
365 Real x1 = intersectHoriz(h10, h11, t10, t11, y);
366 Real x2 = intersectHoriz(h20, h21, t20, t21, y);
367*/
368 Real x1 = intersectHoriz(head1[0], head1[1], tail1[0], tail1[1], y);
369 Real x2 = intersectHoriz(head2[0], head2[1], tail2[0], tail2[1], y);
370
371 if(x1<= x2) return -1;
372 else return 1;
373}
374
375Int compChains(monoChain* mc1, monoChain* mc2)
376{
377 Real y;
378 assert(mc1->isKey || mc2->isKey);
379 if(mc1->isKey)
380 y = mc1->keyY;
381 else
382 y = mc2->keyY;
383 directedLine *d1 = mc1->find(y);
384 directedLine *d2 = mc2->find(y);
385 mc2->find(y);
386// Real x1 = mc1->chainIntersectHoriz(y);
387// Real x2 = mc2->chainIntersectHoriz(y);
388 return compEdges(d1, d2);
389}
390
391//this function modifies current for efficiency
392directedLine* monoChain::find(Real y)
393{
394 directedLine *ret;
395 directedLine *temp;
396 assert(current->head()[1] <= y);
397 if(isIncrease)
398 {
399 assert(chainTail->head()[1] >=y);
400 for(temp=current; temp!=chainTail; temp = temp->getNext())
401 {
402 if(temp->head()[1] > y)
403 break;
404 }
405 current = temp->getPrev();
406 ret = current;
407 }
408 else
409 {
410 for(temp=current; temp != chainHead; temp = temp->getPrev())
411 {
412 if(temp->head()[1] > y)
413 break;
414 }
415 current = temp->getNext();
416 ret = temp;
417 }
418 return ret;
419}
420
421void monoChain::printOneChain()
422{
423 directedLine* temp;
424 for(temp = chainHead; temp != chainTail; temp = temp->getNext())
425 {
426 printf("(%f,%f) ", temp->head()[0], temp->head()[1]);
427 }
428 printf("(%f,%f) \n", chainTail->head()[0], chainTail->head()[1]);
429}
430
431void monoChain::printChainLoop()
432{
433 monoChain* temp;
434 this->printOneChain();
435 for(temp = next; temp != this; temp = temp->next)
436 {
437 temp->printOneChain();
438 }
439 printf("\n");
440}
441
442void monoChain::printAllLoops()
443{
444 monoChain* temp;
445 for(temp=this; temp != NULL; temp = temp->nextPolygon)
446 temp->printChainLoop();
447}
448
449//return 1 if error occures
450Int MC_sweepY(Int nVertices, monoChain** sortedVertices, sweepRange** ret_ranges)
451{
452 Int i;
453 Real keyY;
454 Int errOccur=0;
455//printf("enter MC_sweepY\n");
456//printf("nVertices=%i\n", nVertices);
457 /*for each vertex in the sorted list, update the binary search tree.
458 *and store the range information for each vertex.
459 */
460 treeNode* searchTree = NULL;
461//printf("nVertices=%i\n", nVertices);
462 for(i=0; i<nVertices; i++)
463 {
464 monoChain* vert = sortedVertices[i];
465 keyY = vert->getHead()->head()[1]; //the sweep line
466 directedLine *dline = vert->getHead();
467 directedLine *dlinePrev = dline->getPrev();
468 if(isBelow(dline, dline) && isBelow(dline, dlinePrev))
469 {
470//printf("case 1\n");
471 //this<v and prev < v
472 //delete both edges
473 vert->isKey = 1;
474 vert->keyY = keyY;
475 treeNode* thisNode = TreeNodeFind(searchTree, vert, (Int (*) (void *, void *))compChains);
476 vert->isKey = 0;
477
478 vert->getPrev()->isKey = 1;
479 vert->getPrev()->keyY = keyY;
480 treeNode* prevNode = TreeNodeFind(searchTree, vert->getPrev(), (Int (*) (void *, void *))compChains);
481 vert->getPrev()->isKey = 0;
482
483 if(cuspType(dline) == 1)//interior cusp
484 {
485
486 treeNode* leftEdge = TreeNodePredecessor(prevNode);
487 treeNode* rightEdge = TreeNodeSuccessor(thisNode);
488 if(leftEdge == NULL || rightEdge == NULL)
489 {
490 errOccur = 1;
491 goto JUMP_HERE;
492 }
493
494 directedLine* leftEdgeDline = ((monoChain* ) leftEdge->key)->find(keyY);
495
496
497
498 directedLine* rightEdgeDline = ((monoChain* ) rightEdge->key)->find(keyY);
499
500 ret_ranges[i] = sweepRangeMake(leftEdgeDline, 1, rightEdgeDline, 1);
501 }
502 else /*exterior cusp*/
503 {
504 ret_ranges[i] = sweepRangeMake( dline, 1, dlinePrev, 1);
505 }
506
507 searchTree = TreeNodeDeleteSingleNode(searchTree, thisNode);
508 searchTree = TreeNodeDeleteSingleNode(searchTree, prevNode);
509
510 }
511 else if(isAbove(dline, dline) && isAbove(dline, dlinePrev))
512 {
513//printf("case 2\n");
514 //insert both edges
515 treeNode* thisNode = TreeNodeMake(vert);
516 treeNode* prevNode = TreeNodeMake(vert->getPrev());
517
518 vert->isKey = 1;
519 vert->keyY = keyY;
520 searchTree = TreeNodeInsert(searchTree, thisNode, (Int (*) (void *, void *))compChains);
521 vert->isKey = 0;
522
523 vert->getPrev()->isKey = 1;
524 vert->getPrev()->keyY = keyY;
525 searchTree = TreeNodeInsert(searchTree, prevNode, (Int (*) (void *, void *))compChains);
526 vert->getPrev()->isKey = 0;
527
528 if(cuspType(dline) == 1) //interior cusp
529 {
530//printf("cuspType is 1\n");
531 treeNode* leftEdge = TreeNodePredecessor(thisNode);
532 treeNode* rightEdge = TreeNodeSuccessor(prevNode);
533 if(leftEdge == NULL || rightEdge == NULL)
534 {
535 errOccur = 1;
536 goto JUMP_HERE;
537 }
538//printf("leftEdge is %i, rightEdge is %i\n", leftEdge, rightEdge);
539 directedLine* leftEdgeDline = ((monoChain*) leftEdge->key)->find(keyY);
540 directedLine* rightEdgeDline = ((monoChain*) rightEdge->key)->find(keyY);
541 ret_ranges[i] = sweepRangeMake( leftEdgeDline, 1, rightEdgeDline, 1);
542 }
543 else //exterior cusp
544 {
545//printf("cuspType is not 1\n");
546 ret_ranges[i] = sweepRangeMake(dlinePrev, 1, dline, 1);
547 }
548 }
549 else
550 {
551//printf("%i,%i\n", isAbove(dline, dline), isAbove(dline, dlinePrev));
552 errOccur = 1;
553 goto JUMP_HERE;
554
555 fprintf(stderr, "error in MC_sweepY\n");
556 exit(1);
557 }
558 }
559
560 JUMP_HERE:
561 //finally clean up space: delete the search tree
562 TreeNodeDeleteWholeTree(searchTree);
563 return errOccur;
564}
565
566void MC_findDiagonals(Int total_num_edges, monoChain** sortedVertices,
567 sweepRange** ranges, Int& num_diagonals,
568 directedLine** diagonal_vertices)
569{
570 Int i,j,k;
571 k=0;
572 //reset 'current' of all the monoChains
573 for(i=0; i<total_num_edges; i++)
574 sortedVertices[i]->resetCurrent();
575
576 for(i=0; i<total_num_edges; i++)
577 {
578 directedLine* vert = sortedVertices[i]->getHead();
579 directedLine* thisEdge = vert;
580 directedLine* prevEdge = vert->getPrev();
581 if(isBelow(vert, thisEdge) && isBelow(vert, prevEdge) && compEdges(prevEdge, thisEdge)<0)
582 {
583 //this is an upward interior cusp
584 diagonal_vertices[k++] = vert;
585
586 directedLine* leftEdge = ranges[i]->left;
587 directedLine* rightEdge = ranges[i]->right;
588
589 directedLine* leftVert = leftEdge;
590 directedLine* rightVert = rightEdge->getNext();
591 assert(leftVert->head()[1] >= vert->head()[1]);
592 assert(rightVert->head()[1] >= vert->head()[1]);
593 directedLine* minVert = (leftVert->head()[1] <= rightVert->head()[1])?leftVert:rightVert;
594 Int found = 0;
595 for(j=i+1; j<total_num_edges; j++)
596 {
597 if(sortedVertices[j]->getHead()->head()[1] > minVert->head()[1])
598 break;
599
600 if(sweepRangeEqual(ranges[i], ranges[j]))
601 {
602 found = 1;
603 break;
604 }
605 }
606
607 if(found)
608 diagonal_vertices[k++] = sortedVertices[j]->getHead();
609 else
610 diagonal_vertices[k++] = minVert;
611 }
612 else if(isAbove(vert, thisEdge) && isAbove(vert, prevEdge) && compEdges(prevEdge, thisEdge)>0)
613 {
614 //downward interior cusp
615 diagonal_vertices[k++] = vert;
616 directedLine* leftEdge = ranges[i]->left;
617 directedLine* rightEdge = ranges[i]->right;
618 directedLine* leftVert = leftEdge->getNext();
619 directedLine* rightVert = rightEdge;
620 assert(leftVert->head()[1] <= vert->head()[1]);
621 assert(rightVert->head()[1] <= vert->head()[1]);
622 directedLine* maxVert = (leftVert->head()[1] > rightVert->head()[1])? leftVert:rightVert;
623 Int found=0;
624 for(j=i-1; j>=0; j--)
625 {
626 if(sortedVertices[j]->getHead()->head()[1] < maxVert->head()[1])
627 break;
628 if(sweepRangeEqual(ranges[i], ranges[j]))
629 {
630 found = 1;
631 break;
632 }
633 }
634 if(found)
635 diagonal_vertices[k++] = sortedVertices[j]->getHead();
636 else
637 diagonal_vertices[k++] = maxVert;
638 }
639 }
640 num_diagonals = k/2;
641}
642
643
644
645
646directedLine* MC_partitionY(directedLine *polygons, sampledLine **retSampledLines)
647{
648//printf("enter mc_partitionY\n");
649 Int total_num_chains = 0;
650 monoChain* loopList = directedLineLoopListToMonoChainLoopList(polygons);
651 monoChain** array = loopList->toArrayAllLoops(total_num_chains);
652
653 if(total_num_chains<=2) //there is just one single monotone polygon
654 {
655 loopList->deleteLoopList();
656 free(array);
657 *retSampledLines = NULL;
658 return polygons;
659 }
660
661//loopList->printAllLoops();
662//printf("total_num_chains=%i\n", total_num_chains);
663 quicksort( (void**)array, 0, total_num_chains-1, (Int (*)(void*, void*))compChainHeadInY);
664//printf("after quicksort\n");
665
666 sweepRange** ranges = (sweepRange**)malloc(sizeof(sweepRange*) * (total_num_chains));
667 assert(ranges);
668
669 if(MC_sweepY(total_num_chains, array, ranges))
670 {
671 loopList->deleteLoopList();
672 free(array);
673 *retSampledLines = NULL;
674 return NULL;
675 }
676//printf("after MC_sweepY\n");
677
678
679 Int num_diagonals;
680 /*number diagonals is < total_num_edges*total_num_edges*/
681 directedLine** diagonal_vertices = (directedLine**) malloc(sizeof(directedLine*) * total_num_chains*2/*total_num_edges*/);
682 assert(diagonal_vertices);
683
684//printf("before call MC_findDiagonales\n");
685
686 MC_findDiagonals(total_num_chains, array, ranges, num_diagonals, diagonal_vertices);
687//printf("after call MC_findDia, num_diagnla=%i\n", num_diagonals);
688
689 directedLine* ret_polygons = polygons;
690 sampledLine* newSampledLines = NULL;
691 Int i,k;
692
693 num_diagonals=deleteRepeatDiagonals(num_diagonals, diagonal_vertices, diagonal_vertices);
694
695
696
697//drawDiagonals(num_diagonals, diagonal_vertices);
698//printf("diagoanls are \n");
699//for(i=0; i<num_diagonals; i++)
700// {
701// printf("(%f,%f)\n", diagonal_vertices[2*i]->head()[0], diagonal_vertices[2*i]->head()[1]);
702// printf("**(%f,%f)\n", diagonal_vertices[2*i+1]->head()[0], diagonal_vertices[2*i+1]->head()[1]);
703// }
704
705 Int *removedDiagonals=(Int*)malloc(sizeof(Int) * num_diagonals);
706 for(i=0; i<num_diagonals; i++)
707 removedDiagonals[i] = 0;
708// printf("first pass\n");
709
710
711 for(i=0,k=0; i<num_diagonals; i++,k+=2)
712 {
713
714
715 directedLine* v1=diagonal_vertices[k];
716 directedLine* v2=diagonal_vertices[k+1];
717 directedLine* ret_p1;
718 directedLine* ret_p2;
719
720 /*we ahve to determine whether v1 and v2 belong to the same polygon before
721 *their structure are modified by connectDiagonal().
722 */
723/*
724 directedLine *root1 = v1->findRoot();
725 directedLine *root2 = v2->findRoot();
726 assert(root1);
727 assert(root2);
728*/
729
730directedLine* root1 = v1->rootLinkFindRoot();
731directedLine* root2 = v2->rootLinkFindRoot();
732
733 if(root1 != root2)
734 {
735
736 removedDiagonals[i] = 1;
737 sampledLine* generatedLine;
738
739
740
741 v1->connectDiagonal(v1,v2, &ret_p1, &ret_p2, &generatedLine, ret_polygons);
742
743
744
745 newSampledLines = generatedLine->insert(newSampledLines);
746/*
747 ret_polygons = ret_polygons->cutoffPolygon(root1);
748
749 ret_polygons = ret_polygons->cutoffPolygon(root2);
750 ret_polygons = ret_p1->insertPolygon(ret_polygons);
751root1->rootLinkSet(ret_p1);
752root2->rootLinkSet(ret_p1);
753ret_p1->rootLinkSet(NULL);
754ret_p2->rootLinkSet(ret_p1);
755*/
756 ret_polygons = ret_polygons->cutoffPolygon(root2);
757
758
759
760root2->rootLinkSet(root1);
761ret_p1->rootLinkSet(root1);
762ret_p2->rootLinkSet(root1);
763
764 /*now that we have connected the diagonal v1 and v2,
765 *we have to check those unprocessed diagonals which
766 *have v1 or v2 as an end point. Notice that the head of v1
767 *has the same coodinates as the head of v2->prev, and the head of
768 *v2 has the same coordinate as the head of v1->prev.
769 *Suppose these is a diagonal (v1, x). If (v1,x) is still a valid
770 *diagonal, then x should be on the left hand side of the directed line: *v1->prev->head -- v1->head -- v1->tail. Otherwise, (v1,x) should be
771 *replaced by (v2->prev, x), that is, x is on the left of
772 * v2->prev->prev->head, v2->prev->head, v2->prev->tail.
773 */
774 Int ii, kk;
775 for(ii=0, kk=0; ii<num_diagonals; ii++, kk+=2)
776 if( removedDiagonals[ii]==0)
777 {
778 directedLine* d1=diagonal_vertices[kk];
779 directedLine* d2=diagonal_vertices[kk+1];
780 /*check d1, and replace diagonal_vertices[kk] if necessary*/
781 if(d1 == v1) {
782 /*check if d2 is to left of v1->prev->head:v1->head:v1->tail*/
783 if(! pointLeft2Lines(v1->getPrev()->head(),
784 v1->head(), v1->tail(), d2->head()))
785 {
786/*
787 assert(pointLeft2Lines(v2->getPrev()->getPrev()->head(),
788 v2->getPrev()->head(),
789 v2->getPrev()->tail(), d2->head()));
790*/
791 diagonal_vertices[kk] = v2->getPrev();
792 }
793 }
794 if(d1 == v2) {
795 /*check if d2 is to left of v2->prev->head:v2->head:v2->tail*/
796 if(! pointLeft2Lines(v2->getPrev()->head(),
797 v2->head(), v2->tail(), d2->head()))
798 {
799/*
800 assert(pointLeft2Lines(v1->getPrev()->getPrev()->head(),
801 v1->getPrev()->head(),
802 v1->getPrev()->tail(), d2->head()));
803*/
804 diagonal_vertices[kk] = v1->getPrev();
805 }
806 }
807 /*check d2 and replace diagonal_vertices[k+1] if necessary*/
808 if(d2 == v1) {
809 /*check if d1 is to left of v1->prev->head:v1->head:v1->tail*/
810 if(! pointLeft2Lines(v1->getPrev()->head(),
811 v1->head(), v1->tail(), d1->head()))
812 {
813/* assert(pointLeft2Lines(v2->getPrev()->getPrev()->head(),
814 v2->getPrev()->head(),
815 v2->getPrev()->tail(), d1->head()));
816*/
817 diagonal_vertices[kk+1] = v2->getPrev();
818 }
819 }
820 if(d2 == v2) {
821 /*check if d1 is to left of v2->prev->head:v2->head:v2->tail*/
822 if(! pointLeft2Lines(v2->getPrev()->head(),
823 v2->head(), v2->tail(), d1->head()))
824 {
825/* assert(pointLeft2Lines(v1->getPrev()->getPrev()->head(),
826 v1->getPrev()->head(),
827 v1->getPrev()->tail(), d1->head()));
828*/
829 diagonal_vertices[kk+1] = v1->getPrev();
830 }
831 }
832 }
833}/*end if (root1 not equal to root 2)*/
834}
835
836 /*second pass, now all diagoals should belong to the same polygon*/
837//printf("second pass: \n");
838
839// for(i=0; i<num_diagonals; i++)
840// printf("%i ", removedDiagonals[i]);
841
842
843 for(i=0,k=0; i<num_diagonals; i++, k += 2)
844 if(removedDiagonals[i] == 0)
845 {
846
847
848 directedLine* v1=diagonal_vertices[k];
849 directedLine* v2=diagonal_vertices[k+1];
850
851
852
853 directedLine* ret_p1;
854 directedLine* ret_p2;
855
856 /*we ahve to determine whether v1 and v2 belong to the same polygon before
857 *their structure are modified by connectDiagonal().
858 */
859 directedLine *root1 = v1->findRoot();
860/*
861 directedLine *root2 = v2->findRoot();
862
863
864
865 assert(root1);
866 assert(root2);
867 assert(root1 == root2);
868 */
869 sampledLine* generatedLine;
870
871
872
873 v1->connectDiagonal(v1,v2, &ret_p1, &ret_p2, &generatedLine, ret_polygons);
874 newSampledLines = generatedLine->insert(newSampledLines);
875
876 ret_polygons = ret_polygons->cutoffPolygon(root1);
877
878 ret_polygons = ret_p1->insertPolygon(ret_polygons);
879
880 ret_polygons = ret_p2->insertPolygon(ret_polygons);
881
882
883
884 for(Int j=i+1; j<num_diagonals; j++)
885 {
886 if(removedDiagonals[j] ==0)
887 {
888
889 directedLine* temp1=diagonal_vertices[2*j];
890 directedLine* temp2=diagonal_vertices[2*j+1];
891 if(temp1==v1 || temp1==v2 || temp2==v1 || temp2==v2)
892 if(! temp1->samePolygon(temp1, temp2))
893 {
894 /*if temp1 and temp2 are in different polygons,
895 *then one of them must be v1 or v2.
896 */
897
898
899
900 assert(temp1==v1 || temp1 == v2 || temp2==v1 || temp2 ==v2);
901 if(temp1==v1)
902 {
903 diagonal_vertices[2*j] = v2->getPrev();
904 }
905 if(temp2==v1)
906 {
907 diagonal_vertices[2*j+1] = v2->getPrev();
908 }
909 if(temp1==v2)
910 {
911 diagonal_vertices[2*j] = v1->getPrev();
912 }
913 if(temp2==v2)
914 {
915 diagonal_vertices[2*j+1] = v1->getPrev();
916 }
917 }
918 }
919 }
920
921 }
922
923
924 //clean up
925 loopList->deleteLoopList();
926 free(array);
927 free(ranges);
928 free(diagonal_vertices);
929 free(removedDiagonals);
930
931 *retSampledLines = newSampledLines;
932 return ret_polygons;
933}
934
935