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 <stdlib.h>
39 #include <stdio.h>
40 //#include <time.h>
41 
42 #include "zlassert.h"
43 #include "partitionY.h"
44 #include "searchTree.h"
45 #include "quicksort.h"
46 #include "polyUtil.h"
47 
48 
49 #define max(a,b) ((a>b)? a:b)
50 #define min(a,b) ((a>b)? b:a)
51 
52 
53 /*retrurn
54  *-1: if A < B (Ya<Yb) || (Ya==Yb)
55  * 0: if A == B
56  * 1: if A>B
57  */
compVertInY(Real A[2],Real B[2])58 static Int compVertInY(Real A[2], Real B[2])
59 {
60   if( (A[1] < B[1]) || (A[1]==B[1] && A[0]<B[0]))
61     return -1;
62   else if
63     ( A[1] == B[1] && A[0] == B[0]) return 0;
64   else
65     return 1;
66 }
67 
68 /*v is a vertex: the head of en edge,
69  *e is an edge,
70  *return 1 if e is below v: assume v1 and v2 are the two endpoints of e:
71  * v1<= v, v2<=v.
72  */
isBelow(directedLine * v,directedLine * e)73 Int isBelow(directedLine *v, directedLine *e)
74 {
75   Real* vert = v->head();
76   if(   compVertInY(e->head(), vert) != 1
77      && compVertInY(e->tail(), vert) != 1
78      )
79     return 1;
80   else
81     return 0;
82 }
83 
84 /*v is a vertex: the head of en edge,
85  *e is an edge,
86  *return 1 if e is below v: assume v1 and v2 are the two endpoints of e:
87  * v1>= v, v2>=v.
88  */
isAbove(directedLine * v,directedLine * e)89 Int isAbove(directedLine *v, directedLine *e)
90 {
91   Real* vert = v->head();
92   if(   compVertInY(e->head(), vert) != -1
93      && compVertInY(e->tail(), vert) != -1
94      )
95     return 1;
96   else
97     return 0;
98 }
99 
isCusp(directedLine * v)100 Int isCusp(directedLine *v)
101 {
102   Real *A=v->getPrev()->head();
103   Real *B=v->head();
104   Real *C=v->tail();
105   if(A[1] < B[1] && B[1] < C[1])
106     return 0;
107   else if(A[1] > B[1] && B[1] > C[1])
108     return 0;
109   else if(A[1] < B[1] && C[1] < B[1])
110     return 1;
111   else if(A[1] > B[1] && C[1] > B[1])
112     return 1;
113 
114   if((isAbove(v, v) && isAbove(v, v->getPrev())) ||
115      (isBelow(v, v) && isBelow(v, v->getPrev())))
116     return 1;
117   else
118     return 0;
119 }
120 
121 /*crossproduct is strictly less than 0*/
isReflex(directedLine * v)122 Int isReflex(directedLine *v)
123 {
124   Real* A = v->getPrev()->head();
125   Real* B = v->head();
126   Real* C = v->tail();
127   Real Bx,By, Cx, Cy;
128   Bx = B[0] - A[0];
129   By = B[1] - A[1];
130   Cx = C[0] - A[0];
131   Cy = C[1] - A[1];
132 
133   if(Bx*Cy - Cx*By < 0) return 1;
134   else return 0;
135 }
136 
137  /*return
138  *0: not-cusp
139  *1: interior cusp
140  *2: exterior cusp
141  */
cuspType(directedLine * v)142 Int cuspType(directedLine *v)
143 {
144   if(! isCusp(v)) return 0;
145   else if(isReflex(v)) return 1;
146   else
147     return 2;
148 }
149 
sweepRangeMake(directedLine * left,Int leftType,directedLine * right,Int rightType)150 sweepRange* sweepRangeMake(directedLine* left, Int leftType,
151 			   directedLine* right, Int rightType)
152 {
153   sweepRange* ret = (sweepRange*)malloc(sizeof(sweepRange));
154   assert(ret);
155   ret->left = left;
156   ret->leftType = leftType;
157   ret->right = right;
158   ret->rightType = rightType;
159   return ret;
160 }
161 
sweepRangeDelete(sweepRange * range)162 void sweepRangeDelete(sweepRange* range)
163 {
164   free(range);
165 }
166 
sweepRangeEqual(sweepRange * src1,sweepRange * src2)167 Int sweepRangeEqual(sweepRange* src1, sweepRange* src2)
168 {
169   Int leftEqual;
170   Int rightEqual;
171 
172 
173   /*The case when both are vertices should not happen*/
174   assert(! (src1->leftType == 0 && src2->leftType == 0));
175   if(src1->leftType == 0 && src2->leftType == 1){
176     if(src1->left == src2->left ||
177        src1->left->getPrev() == src2->left
178        )
179       leftEqual = 1;
180     else
181       leftEqual = 0;
182   }
183   else if(src1->leftType == 1 && src2->leftType == 1){
184     if(src1->left == src2->left)
185       leftEqual = 1;
186     else
187       leftEqual = 0;
188   }
189   else /*src1->leftType == 1 && src2->leftType == 0*/{
190     if(src1->left == src2->left ||
191        src1->left == src2->left->getPrev()
192        )
193       leftEqual = 1;
194     else
195       leftEqual = 0;
196   }
197 
198   /*the same thing for right*/
199   /*The case when both are vertices should not happen*/
200   assert(! (src1->rightType == 0 && src2->rightType == 0));
201   if(src1->rightType == 0 && src2->rightType == 1){
202     if(src1->right == src2->right ||
203        src1->right->getPrev() == src2->right
204        )
205       rightEqual = 1;
206     else
207       rightEqual = 0;
208   }
209   else if(src1->rightType == 1 && src2->rightType == 1){
210     if(src1->right == src2->right)
211       rightEqual = 1;
212     else
213       rightEqual = 0;
214   }
215   else /*src1->rightType == 1 && src2->rightType == 0*/{
216     if(src1->right == src2->right ||
217        src1->right == src2->right->getPrev()
218        )
219       rightEqual = 1;
220     else
221       rightEqual = 0;
222   }
223 
224   return (leftEqual == 1 || rightEqual == 1);
225 }
226 
227 /*given (x_1, y_1) and (x_2, y_2), and y
228  *return x such that (x,y) is on the line
229  */
intersectHoriz(Real x1,Real y1,Real x2,Real y2,Real y)230 inline/*static*/ Real intersectHoriz(Real x1, Real y1, Real x2, Real y2, Real y)
231 {
232   return ((y2==y1)? (x1+x2)*Real(0.5) : x1 + ((y-y1)/(y2-y1)) * (x2-x1));
233 /*
234   if(y2 == y1) return (x1+x2)*0.5;
235   else return x1 + ((y-y1)/(y2-y1)) * (x2-x1);
236 */
237 }
238 
239 /*compare two edges of a polygon.
240  *edge A < edge B if there is a horizontal line so that the intersection
241  *with A is to the left of the intersection with B.
242  *This function is used in sweepY for the dynamic search tree insertion to
243  *order the edges.
244  * Implementation: (x_1,y_1) and (x_2, y_2)
245  */
compEdges(directedLine * e1,directedLine * e2)246 static Int compEdges(directedLine *e1, directedLine *e2)
247 {
248   Real* head1 = e1->head();
249   Real* tail1 = e1->tail();
250   Real* head2 = e2->head();
251   Real* tail2 = e2->tail();
252 /*
253   Real h10 = head1[0];
254   Real h11 = head1[1];
255   Real t10 = tail1[0];
256   Real t11 = tail1[1];
257   Real h20 = head2[0];
258   Real h21 = head2[1];
259   Real t20 = tail2[0];
260   Real t21 = tail2[1];
261 */
262   Real e1_Ymax, e1_Ymin, e2_Ymax, e2_Ymin;
263 /*
264   if(h11>t11) {
265     e1_Ymax= h11;
266     e1_Ymin= t11;
267   }
268   else{
269     e1_Ymax = t11;
270     e1_Ymin = h11;
271   }
272 
273   if(h21>t21) {
274     e2_Ymax= h21;
275     e2_Ymin= t21;
276   }
277   else{
278     e2_Ymax = t21;
279     e2_Ymin = h21;
280   }
281 */
282 
283   if(head1[1]>tail1[1]) {
284     e1_Ymax= head1[1];
285     e1_Ymin= tail1[1];
286   }
287   else{
288     e1_Ymax = tail1[1];
289     e1_Ymin = head1[1];
290   }
291 
292   if(head2[1]>tail2[1]) {
293     e2_Ymax= head2[1];
294     e2_Ymin= tail2[1];
295   }
296   else{
297     e2_Ymax = tail2[1];
298     e2_Ymin = head2[1];
299   }
300 
301 
302   /*Real e1_Ymax = max(head1[1], tail1[1]);*/ /*max(e1->head()[1], e1->tail()[1]);*/
303   /*Real e1_Ymin = min(head1[1], tail1[1]);*/ /*min(e1->head()[1], e1->tail()[1]);*/
304   /*Real e2_Ymax = max(head2[1], tail2[1]);*/ /*max(e2->head()[1], e2->tail()[1]);*/
305   /*Real e2_Ymin = min(head2[1], tail2[1]);*/ /*min(e2->head()[1], e2->tail()[1]);*/
306 
307   Real Ymax = min(e1_Ymax, e2_Ymax);
308   Real Ymin = max(e1_Ymin, e2_Ymin);
309 
310   Real y = Real(0.5)*(Ymax + Ymin);
311 
312 /*  Real x1 = intersectHoriz(e1->head()[0], e1->head()[1], e1->tail()[0], e1->tail()[1], y);
313   Real x2 = intersectHoriz(e2->head()[0], e2->head()[1], e2->tail()[0], e2->tail()[1], y);
314 */
315 /*
316   Real x1 = intersectHoriz(h10, h11, t10, t11, y);
317   Real x2 = intersectHoriz(h20, h21, t20, t21, y);
318 */
319   Real x1 = intersectHoriz(head1[0], head1[1], tail1[0], tail1[1], y);
320   Real x2 = intersectHoriz(head2[0], head2[1], tail2[0], tail2[1], y);
321 
322   if(x1<= x2) return -1;
323   else return 1;
324 }
325 
326 /*used by sort precedures
327  */
compInY(directedLine * v1,directedLine * v2)328 static Int compInY(directedLine* v1, directedLine* v2)
329 {
330   return v1->compInY(v2);
331 }
332 
findDiagonals(Int total_num_edges,directedLine ** sortedVertices,sweepRange ** ranges,Int & num_diagonals,directedLine ** diagonal_vertices)333 void findDiagonals(Int total_num_edges, directedLine** sortedVertices, sweepRange** ranges, Int& num_diagonals, directedLine** diagonal_vertices)
334 {
335   Int i,j,k;
336 
337   k=0;
338 
339   for(i=0; i<total_num_edges; i++)
340     {
341       directedLine* vert =sortedVertices[i];
342       directedLine* thisEdge = vert;
343       directedLine* prevEdge = vert->getPrev();
344 /*
345 printf("find i=%i\n", i);
346 printf("the vertex is\n");
347 vert->printSingle();
348 */
349       if(isBelow(vert, thisEdge) && isBelow(vert, prevEdge) && compEdges(prevEdge, thisEdge)<0)
350 	{
351 	  /*this is an upward interior cusp*/
352 	  diagonal_vertices[k++] = vert;
353 
354 	  for(j=i+1; j<total_num_edges; j++)
355 	    if(sweepRangeEqual(ranges[i], ranges[j]))
356 	      {
357 		diagonal_vertices[k++] = sortedVertices[j];
358 		break;
359 	      }
360 	  assert(j<total_num_edges);
361 
362 
363 	}
364       else if(isAbove(vert, thisEdge) && isAbove(vert, prevEdge) && compEdges(prevEdge, thisEdge)>0)
365 	{
366 	  /*this is an downward interior cusp*/
367 	  diagonal_vertices[k++] = vert;
368 	  for(j=i-1; j>=0; j--)
369 	    if(sweepRangeEqual(ranges[i], ranges[j]))
370 	      {
371 		diagonal_vertices[k++] = sortedVertices[j];
372 		break;
373 	      }
374 /*	  printf("j=%i\n", j);*/
375 	  assert(j>=0);
376 
377 
378 
379 	}
380     }
381   num_diagonals = k/2;
382 }
383 
384 /*get rid of repeated diagonlas so that each diagonal appears only once in the array
385  */
deleteRepeatDiagonals(Int num_diagonals,directedLine ** diagonal_vertices,directedLine ** new_vertices)386 Int deleteRepeatDiagonals(Int num_diagonals, directedLine** diagonal_vertices, directedLine** new_vertices)
387 {
388   Int i,k;
389   Int j,l;
390   Int index;
391   index=0;
392   for(i=0,k=0; i<num_diagonals; i++, k+=2)
393     {
394       Int isRepeated=0;
395       /*check the diagonla (diagonal_vertice[k], diagonal_vertices[k+1])
396        *is repeated or not
397        */
398       for(j=0,l=0; j<index; j++, l+=2)
399 	{
400 	  if(
401 	     (diagonal_vertices[k] == new_vertices[l] &&
402 	      diagonal_vertices[k+1] == new_vertices[l+1]
403 	      )
404 	     ||
405 	     (
406 	      diagonal_vertices[k] == new_vertices[l+1] &&
407 	      diagonal_vertices[k+1] == new_vertices[l]
408 	      )
409 	     )
410 	    {
411 	      isRepeated=1;
412 	      break;
413 	    }
414 	}
415       if(! isRepeated)
416 	{
417 	  new_vertices[index+index] = diagonal_vertices[k];
418 	  new_vertices[index+index+1] = diagonal_vertices[k+1];
419 	  index++;
420 	}
421     }
422   return index;
423 }
424 
425 /*for debug only*/
DBGfindDiagonals(directedLine * polygons,Int & num_diagonals)426 directedLine** DBGfindDiagonals(directedLine *polygons, Int& num_diagonals)
427 {
428   Int total_num_edges = 0;
429   directedLine** array = polygons->toArrayAllPolygons(total_num_edges);
430   quicksort( (void**)array, 0, total_num_edges-1, (Int (*)(void*, void*)) compInY);
431   sweepRange** ranges = (sweepRange**) malloc(sizeof(sweepRange*) * total_num_edges);
432   assert(ranges);
433 
434   sweepY(total_num_edges, array, ranges);
435 
436  directedLine** diagonal_vertices = (directedLine**) malloc(sizeof(directedLine*) * total_num_edges);
437   assert(diagonal_vertices);
438   findDiagonals(total_num_edges, array, ranges, num_diagonals, diagonal_vertices);
439 
440   num_diagonals=deleteRepeatDiagonals(num_diagonals, diagonal_vertices, diagonal_vertices);
441   return diagonal_vertices;
442 
443 }
444 
445 
446 /*partition into Y-monotone polygons*/
partitionY(directedLine * polygons,sampledLine ** retSampledLines)447 directedLine* partitionY(directedLine *polygons, sampledLine **retSampledLines)
448 {
449   Int total_num_edges = 0;
450   directedLine** array = polygons->toArrayAllPolygons(total_num_edges);
451 
452   quicksort( (void**)array, 0, total_num_edges-1, (Int (*)(void*, void*)) compInY);
453 
454   sweepRange** ranges = (sweepRange**) malloc(sizeof(sweepRange*) * (total_num_edges));
455   assert(ranges);
456 
457 
458 
459   sweepY(total_num_edges, array, ranges);
460 
461 
462 
463   /*the diagonal vertices are stored as:
464    *v0-v1: 1st diagonal
465    *v2-v3: 2nd diagonal
466    *v5-v5: 3rd diagonal
467    *...
468    */
469 
470 
471   Int num_diagonals;
472   /*number diagonals is < total_num_edges*total_num_edges*/
473   directedLine** diagonal_vertices = (directedLine**) malloc(sizeof(directedLine*) * total_num_edges*2/*total_num_edges*/);
474   assert(diagonal_vertices);
475 
476 
477 
478   findDiagonals(total_num_edges, array, ranges, num_diagonals, diagonal_vertices);
479 
480 
481 
482   directedLine* ret_polygons = polygons;
483   sampledLine* newSampledLines = NULL;
484   Int i,k;
485 
486 num_diagonals=deleteRepeatDiagonals(num_diagonals, diagonal_vertices, diagonal_vertices);
487 
488 
489 
490   Int *removedDiagonals=(Int*)malloc(sizeof(Int) * num_diagonals);
491   for(i=0; i<num_diagonals; i++)
492     removedDiagonals[i] = 0;
493 
494 
495 
496 
497 
498   for(i=0,k=0; i<num_diagonals; i++,k+=2)
499     {
500 
501 
502       directedLine* v1=diagonal_vertices[k];
503       directedLine* v2=diagonal_vertices[k+1];
504       directedLine* ret_p1;
505       directedLine* ret_p2;
506 
507       /*we ahve to determine whether v1 and v2 belong to the same polygon before
508        *their structure are modified by connectDiagonal().
509        */
510 /*
511       directedLine *root1 = v1->findRoot();
512       directedLine *root2 = v2->findRoot();
513       assert(root1);
514       assert(root2);
515 */
516 
517 directedLine*  root1 = v1->rootLinkFindRoot();
518 directedLine*  root2 = v2->rootLinkFindRoot();
519 
520       if(root1 != root2)
521 	{
522 
523 	  removedDiagonals[i] = 1;
524 	  sampledLine* generatedLine;
525 
526 
527 
528 	  v1->connectDiagonal(v1,v2, &ret_p1, &ret_p2, &generatedLine, ret_polygons);
529 
530 
531 
532 	  newSampledLines = generatedLine->insert(newSampledLines);
533 /*
534 	  ret_polygons = ret_polygons->cutoffPolygon(root1);
535 
536 	  ret_polygons = ret_polygons->cutoffPolygon(root2);
537 	  ret_polygons = ret_p1->insertPolygon(ret_polygons);
538 root1->rootLinkSet(ret_p1);
539 root2->rootLinkSet(ret_p1);
540 ret_p1->rootLinkSet(NULL);
541 ret_p2->rootLinkSet(ret_p1);
542 */
543 	  ret_polygons = ret_polygons->cutoffPolygon(root2);
544 
545 
546 
547 root2->rootLinkSet(root1);
548 ret_p1->rootLinkSet(root1);
549 ret_p2->rootLinkSet(root1);
550 
551        /*now that we have connected the diagonal v1 and v2,
552         *we have to check those unprocessed diagonals which
553         *have v1 or v2 as an end point. Notice that the head of v1
554         *has the same coodinates as the head of v2->prev, and the head of
555         *v2 has the same coordinate as the head of v1->prev.
556         *Suppose these is a diagonal (v1, x). If (v1,x) is still a valid
557         *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
558         *replaced by (v2->prev, x), that is, x is on the left of
559         * v2->prev->prev->head, v2->prev->head, v2->prev->tail.
560         */
561         Int ii, kk;
562         for(ii=0, kk=0; ii<num_diagonals; ii++, kk+=2)
563 	  if( removedDiagonals[ii]==0)
564 	    {
565 	      directedLine* d1=diagonal_vertices[kk];
566 	      directedLine* d2=diagonal_vertices[kk+1];
567 	      /*check d1, and replace diagonal_vertices[kk] if necessary*/
568 	      if(d1 == v1) {
569 		/*check if d2 is to left of v1->prev->head:v1->head:v1->tail*/
570 		if(! pointLeft2Lines(v1->getPrev()->head(),
571 				     v1->head(), v1->tail(), d2->head()))
572 		  {
573 /*
574 		    assert(pointLeft2Lines(v2->getPrev()->getPrev()->head(),
575 					   v2->getPrev()->head(),
576 					   v2->getPrev()->tail(), d2->head()));
577 */
578 		    diagonal_vertices[kk] = v2->getPrev();
579 		  }
580 	      }
581 	      if(d1 == v2) {
582 		/*check if d2 is to left of v2->prev->head:v2->head:v2->tail*/
583 		if(! pointLeft2Lines(v2->getPrev()->head(),
584 				     v2->head(), v2->tail(), d2->head()))
585 		  {
586 /*
587 		    assert(pointLeft2Lines(v1->getPrev()->getPrev()->head(),
588 					   v1->getPrev()->head(),
589 					   v1->getPrev()->tail(), d2->head()));
590 */
591 		    diagonal_vertices[kk] = v1->getPrev();
592 		  }
593 	      }
594 	      /*check d2 and replace diagonal_vertices[k+1] if necessary*/
595 	      if(d2 == v1) {
596 		/*check if d1 is to left of v1->prev->head:v1->head:v1->tail*/
597 		if(! pointLeft2Lines(v1->getPrev()->head(),
598 				     v1->head(), v1->tail(), d1->head()))
599 		  {
600 /*		    assert(pointLeft2Lines(v2->getPrev()->getPrev()->head(),
601 					   v2->getPrev()->head(),
602 					   v2->getPrev()->tail(), d1->head()));
603 */
604 		    diagonal_vertices[kk+1] = v2->getPrev();
605 		  }
606 	      }
607 	      if(d2 == v2) {
608 		/*check if d1 is to left of v2->prev->head:v2->head:v2->tail*/
609 		if(! pointLeft2Lines(v2->getPrev()->head(),
610 				     v2->head(), v2->tail(), d1->head()))
611 		  {
612 /*		    assert(pointLeft2Lines(v1->getPrev()->getPrev()->head(),
613 					   v1->getPrev()->head(),
614 					   v1->getPrev()->tail(), d1->head()));
615 */
616 		    diagonal_vertices[kk+1] = v1->getPrev();
617 		  }
618 	      }
619 	    }
620 }/*end if (root1 not equal to root 2)*/
621 }
622 
623   /*second pass,  now all diagoals should belong to the same polygon*/
624 
625 
626 
627   for(i=0,k=0; i<num_diagonals; i++, k += 2)
628     if(removedDiagonals[i] == 0)
629       {
630 
631 
632 	directedLine* v1=diagonal_vertices[k];
633 	directedLine* v2=diagonal_vertices[k+1];
634 
635 
636 
637 	directedLine* ret_p1;
638 	directedLine* ret_p2;
639 
640 	/*we ahve to determine whether v1 and v2 belong to the same polygon before
641 	 *their structure are modified by connectDiagonal().
642 	 */
643 	directedLine *root1 = v1->findRoot();
644 /*
645 	directedLine *root2 = v2->findRoot();
646 
647 
648 
649 	assert(root1);
650 	assert(root2);
651 	assert(root1 == root2);
652   */
653 	sampledLine* generatedLine;
654 
655 
656 
657 	v1->connectDiagonal(v1,v2, &ret_p1, &ret_p2, &generatedLine, ret_polygons);
658 	newSampledLines = generatedLine->insert(newSampledLines);
659 
660 	ret_polygons = ret_polygons->cutoffPolygon(root1);
661 
662 	ret_polygons = ret_p1->insertPolygon(ret_polygons);
663 
664 	ret_polygons = ret_p2->insertPolygon(ret_polygons);
665 
666 
667 
668 	for(Int j=i+1; j<num_diagonals; j++)
669 	  {
670 	    if(removedDiagonals[j] ==0)
671 	      {
672 
673 		directedLine* temp1=diagonal_vertices[2*j];
674 		directedLine* temp2=diagonal_vertices[2*j+1];
675                if(temp1==v1 || temp1==v2 || temp2==v1 || temp2==v2)
676 		if(! temp1->samePolygon(temp1, temp2))
677 		  {
678 		    /*if temp1 and temp2 are in different polygons,
679 		     *then one of them must be v1 or v2.
680 		     */
681 
682 
683 
684 		    assert(temp1==v1 || temp1 == v2 || temp2==v1 || temp2 ==v2);
685 		    if(temp1==v1)
686 		      {
687 			diagonal_vertices[2*j] = v2->getPrev();
688 		      }
689 		    if(temp2==v1)
690 		      {
691 			diagonal_vertices[2*j+1] = v2->getPrev();
692 		      }
693 		    if(temp1==v2)
694 		      {
695 			diagonal_vertices[2*j] = v1->getPrev();
696 		      }
697 		    if(temp2==v2)
698 		      {
699 			diagonal_vertices[2*j+1] = v1->getPrev();
700 		      }
701 		  }
702 	      }
703 	  }
704 
705       }
706 
707   /*clean up spaces*/
708   free(array);
709   free(ranges);
710   free(diagonal_vertices);
711   free(removedDiagonals);
712 
713   *retSampledLines = newSampledLines;
714   return ret_polygons;
715 }
716 
717 /*given a set of simple polygons where the interior
718  *is decided by left-hand principle,
719  *return a range (sight) for each vertex. This is called
720  *Trapezoidalization.
721  */
sweepY(Int nVertices,directedLine ** sortedVertices,sweepRange ** ret_ranges)722 void sweepY(Int nVertices, directedLine** sortedVertices, sweepRange** ret_ranges)
723 {
724   Int i;
725   /*for each vertex in the sorted list, update the binary search tree.
726    *and store the range information for each vertex.
727    */
728   treeNode* searchTree = NULL;
729   for(i=0; i<nVertices;i++)
730     {
731 
732       directedLine* vert = sortedVertices[i];
733 
734       directedLine* thisEdge = vert;
735       directedLine* prevEdge = vert->getPrev();
736 
737       if(isBelow(vert, thisEdge) && isAbove(vert, prevEdge))
738 	{
739 
740 	  /*case 1: this < v < prev
741 	   *the polygon is going down at v, the interior is to
742 	   *the right hand side.
743 	   * find the edge to the right of thisEdge for right range.
744            * delete thisEdge
745            * insert prevEdge
746 	   */
747 	  treeNode* thisNode = TreeNodeFind(searchTree, thisEdge, ( Int (*) (void *, void *))compEdges);
748 	  assert(thisNode);
749 
750 	  treeNode* succ = TreeNodeSuccessor(thisNode);
751 	  assert(succ);
752 	  searchTree = TreeNodeDeleteSingleNode(searchTree, thisNode);
753 	  searchTree = TreeNodeInsert(searchTree, TreeNodeMake(prevEdge), ( Int (*) (void *, void *))compEdges);
754 
755 
756 	  ret_ranges[i] = sweepRangeMake(vert, 0, (directedLine*) (succ->key), 1);
757 
758 	}
759       else if(isAbove(vert, thisEdge) && isBelow(vert, prevEdge))
760 	{
761 
762 	  /*case 2: this > v > prev
763 	   *the polygon is going up at v, the interior is to
764 	   *the left hand side.
765 	   * find the edge to the left of thisEdge for left range.
766            * delete prevEdge
767            * insert thisEdge
768 	   */
769 	  treeNode* prevNode = TreeNodeFind(searchTree, prevEdge, ( Int (*) (void *, void *))compEdges);
770 	  assert(prevNode);
771 	  treeNode* pred = TreeNodePredecessor(prevNode);
772 	  searchTree = TreeNodeDeleteSingleNode(searchTree, prevNode);
773 	  searchTree = TreeNodeInsert(searchTree, TreeNodeMake(thisEdge), ( Int (*) (void *, void *))compEdges);
774 	  ret_ranges[i] = sweepRangeMake((directedLine*)(pred->key), 1, vert, 0);
775 	}
776       else if(isAbove(vert, thisEdge) && isAbove(vert, prevEdge))
777 	{
778 
779 	  /*case 3: insert both edges*/
780 	  treeNode* thisNode = TreeNodeMake(thisEdge);
781 	  treeNode* prevNode = TreeNodeMake(prevEdge);
782 	  searchTree = TreeNodeInsert(searchTree, thisNode, ( Int (*) (void *, void *))compEdges);
783 	  searchTree = TreeNodeInsert(searchTree, prevNode, ( Int (*) (void *, void *))compEdges);
784 	  if(compEdges(thisEdge, prevEdge)<0)  /*interior cusp*/
785 	    {
786 
787 	      treeNode* leftEdge = TreeNodePredecessor(thisNode);
788 	      treeNode* rightEdge = TreeNodeSuccessor(prevNode);
789 	      ret_ranges[i] = sweepRangeMake( (directedLine*) leftEdge->key, 1,
790 					     (directedLine*) rightEdge->key, 1
791 					     );
792 	    }
793 	  else /*exterior cusp*/
794 	    {
795 
796 	      ret_ranges[i] = sweepRangeMake( prevEdge, 1, thisEdge, 1);
797 	    }
798 	}
799       else if(isBelow(vert, thisEdge) && isBelow(vert, prevEdge))
800 	{
801 
802 	  /*case 4: delete both edges*/
803 	  treeNode* thisNode = TreeNodeFind(searchTree, thisEdge, ( Int (*) (void *, void *))compEdges);
804 	  treeNode* prevNode = TreeNodeFind(searchTree, prevEdge, ( Int (*) (void *, void *))compEdges);
805 	  if(compEdges(thisEdge, prevEdge)>0) /*interior cusp*/
806 	    {
807 	      treeNode* leftEdge = TreeNodePredecessor(prevNode);
808 	      treeNode* rightEdge = TreeNodeSuccessor(thisNode);
809 	      ret_ranges[i] = sweepRangeMake( (directedLine*) leftEdge->key, 1,
810 					     (directedLine*) rightEdge->key, 1
811 					     );
812 	    }
813 	  else /*exterior cusp*/
814 	    {
815 	      ret_ranges[i] = sweepRangeMake( thisEdge, 1, prevEdge, 1);
816 	    }
817 	  searchTree = TreeNodeDeleteSingleNode(searchTree, thisNode);
818 	  searchTree = TreeNodeDeleteSingleNode(searchTree, prevNode);
819 	}
820       else
821 	{
822 	  fprintf(stderr,"error in partitionY.C, invalid case\n");
823 	  printf("vert is\n");
824 	  vert->printSingle();
825 	  printf("thisEdge is\n");
826 	  thisEdge->printSingle();
827 	  printf("prevEdge is\n");
828 	  prevEdge->printSingle();
829 
830 	  exit(1);
831 	}
832     }
833 
834   /*finaly clean up space: delete the search tree*/
835   TreeNodeDeleteWholeTree(searchTree);
836 }
837