src/glu/sgi/libnurbs/nurbtess/monoTriangulation.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 ** $Date: 2005/10/28 13:09:23 $ $Revision: 1.4 $
  35 */
  36 /*
  37 ** $Header: /home/krh/git/sync/mesa-cvs-repo/Mesa/src/glu/sgi/libnurbs/nurbtess/monoTriangulation.cc,v 1.4 2005/10/28 13:09:23 brianp Exp $
  38 */
  39
  40 #include <stdlib.h>
  41 #include <stdio.h>
  42 #include "gluos.h"
  43 #include "glimports.h"
  44 #include "zlassert.h"
  45
  46 #include "monoTriangulation.h"
  47 #include "polyUtil.h" /*for area*/
  48 #include "partitionX.h"
  49 #include "monoPolyPart.h"
  50
  51
  52
  53 extern  directedLine*  polygonConvert(directedLine* polygon);
  54
  55 /*poly is NOT deleted
  56  */
  57 void monoTriangulationOpt(directedLine* poly, primStream* pStream)
  58 {
  59   Int n_cusps;
  60   Int n_edges = poly->numEdges();
  61   directedLine** cusps = (directedLine**) malloc(sizeof(directedLine*)*n_edges);
  62   assert(cusps);
  63   findInteriorCuspsX(poly, n_cusps, cusps);
  64   if(n_cusps ==0) //u monotine
  65     {
  66       monoTriangulationFun(poly, compV2InX, pStream);
  67     }
  68   else if(n_cusps == 1) // one interior cusp
  69     {
  70       directedLine* new_polygon = polygonConvert(cusps[0]);
  71       directedLine* other = findDiagonal_singleCuspX(new_polygon);
  72       //<other> should NOT be null unless there are self-intersecting
  73       //trim curves. In that case, we don't want to core dump, instead,
  74       //we triangulate anyway, and print out error message.
  75       if(other == NULL)
  76         {
  77           monoTriangulationFun(poly, compV2InX, pStream);
  78         }
  79       else
  80         {
  81           directedLine* ret_p1;
  82           directedLine* ret_p2;
  83
  84           new_polygon->connectDiagonal_2slines(new_polygon, other,
  85                                                &ret_p1,
  86                                                &ret_p2,
  87                                                new_polygon);
  88
  89           monoTriangulationFun(ret_p1, compV2InX, pStream);
  90           monoTriangulationFun(ret_p2, compV2InX, pStream);
  91
  92           ret_p1->deleteSinglePolygonWithSline();
  93           ret_p2->deleteSinglePolygonWithSline();
  94         }
  95     }
  96   else
  97     {
  98       //we need a general partitionX funtion (supposed to be in partitionX.C,
  99       //not implemented yet. XXX
 100       monoTriangulationFun(poly, compV2InY, pStream);
 101     }
 102
 103   free(cusps);
 104 }
 105
 106 void monoTriangulationRecOpt(Real* topVertex, Real* botVertex,
 107                              vertexArray* left_chain, Int left_current,
 108                              vertexArray* right_chain, Int right_current,
 109                              primStream* pStream)
 110 {
 111   Int i,j;
 112   Int n_left = left_chain->getNumElements();
 113   Int n_right = right_chain->getNumElements();
 114   if(left_current>= n_left-1 ||
 115      right_current>= n_right-1)
 116     {
 117       monoTriangulationRec(topVertex, botVertex, left_chain, left_current,
 118                            right_chain, right_current, pStream);
 119       return;
 120     }
 121   //now both left and right  have at least two vertices each.
 122   Real left_v = left_chain->getVertex(left_current)[1];
 123   Real right_v = right_chain->getVertex(right_current)[1];
 124
 125   if(left_v <= right_v) //first left vertex is below right
 126     {
 127       //find the last vertex of right which is above or equal to left
 128       for(j=right_current; j<=n_right-1; j++)
 129         {
 130           if(right_chain->getVertex(j)[1] < left_v)
 131             break;
 132         }
 133       monoTriangulationRecGen(topVertex, left_chain->getVertex(left_current),
 134                               left_chain, left_current, left_current,
 135                               right_chain, right_current, j-1,
 136                               pStream);
 137       monoTriangulationRecOpt(right_chain->getVertex(j-1),
 138                               botVertex,
 139                               left_chain, left_current,
 140                               right_chain, j,
 141                               pStream);
 142     }
 143   else //first right vertex is strictly below left
 144     {
 145       //find the last vertex of left which is strictly above right
 146       for(i=left_current; i<=n_left-1; i++)
 147         {
 148           if(left_chain->getVertex(i)[1] <= right_v)
 149             break;
 150         }
 151       monoTriangulationRecGen(topVertex, right_chain->getVertex(right_current),
 152                               left_chain, left_current, i-1,
 153                               right_chain, right_current, right_current,
 154                               pStream);
 155       monoTriangulationRecOpt(left_chain->getVertex(i-1),
 156                               botVertex,
 157                               left_chain, i,
 158                               right_chain, right_current,
 159                               pStream);
 160     }
 161 }
 162
 163
 164 void monoTriangulationRecGenTBOpt(Real* topVertex, Real* botVertex,
 165                           vertexArray* inc_chain, Int inc_current, Int inc_end,
 166                           vertexArray* dec_chain, Int dec_current, Int dec_end,
 167                           primStream* pStream)
 168 {
 169   pStream->triangle(topVertex, inc_chain->getVertex(inc_current), dec_chain->getVertex(dec_current));
 170
 171 /*printf("**(%f,%f)\n", inc_chain->getArray()[0][0],inc_chain->getArray()[0][1]);*/
 172   triangulateXYMonoTB(inc_end-inc_current+1, inc_chain->getArray()+inc_current,  dec_end-dec_current+1,  dec_chain->getArray()+dec_current, pStream);
 173
 174   pStream->triangle(botVertex, dec_chain->getVertex(dec_end), inc_chain->getVertex(inc_end));
 175 }
 176
 177
 178 /*n_left>=1
 179  *n_right>=1
 180  *the strip is going top to bottom. compared to the funtion
 181  * triangulateXYmono()
 182  */
 183 void triangulateXYMonoTB(Int n_left, Real** leftVerts,
 184                        Int n_right, Real** rightVerts,
 185                        primStream* pStream)
 186 {
 187
 188
 189   Int i,j,k,l;
 190   Real* topMostV;
 191
 192   assert(n_left>=1 && n_right>=1);
 193   if(leftVerts[0][1] >= rightVerts[0][1])
 194     {
 195       i=1;
 196       j=0;
 197       topMostV = leftVerts[0];
 198     }
 199   else
 200     {
 201       i=0;
 202       j=1;
 203       topMostV = rightVerts[0];
 204     }
 205
 206   while(1)
 207     {
 208       if(i >= n_left) /*case1: no more in left*/
 209         {
 210
 211           if(j<n_right-1) /*at least two vertices in right*/
 212             {
 213               pStream->begin();
 214               pStream->insert(topMostV);
 215               for(k=n_right-1; k>=j; k--)
 216                 pStream->insert(rightVerts[j]);
 217
 218               pStream->end(PRIMITIVE_STREAM_FAN);
 219
 220             }
 221
 222           break;
 223         }
 224       else if(j>= n_right) /*case2: no more in right*/
 225         {
 226
 227           if(i<n_left-1) /*at least two vertices in left*/
 228             {
 229               pStream->begin();
 230               pStream->insert(topMostV);
 231
 232               for(k=i; k<n_left; k++)
 233                 pStream->insert(leftVerts[k]);
 234
 235               pStream->end(PRIMITIVE_STREAM_FAN);
 236             }
 237
 238           break;
 239         }
 240       else /* case3: neither is empty, plus the topMostV, there is at least one triangle to output*/
 241         {
 242
 243           if(leftVerts[i][1] >=  rightVerts[j][1])
 244             {
 245               pStream->begin();
 246               pStream->insert(rightVerts[j]); /*the origin of this fan*/
 247
 248               pStream->insert(topMostV);
 249
 250               /*find the last k>=i such that
 251                *leftverts[k][1] >= rightverts[j][1]
 252                */
 253               k=i;
 254               while(k<n_left)
 255                 {
 256                   if(leftVerts[k][1] < rightVerts[j][1])
 257                     break;
 258                   k++;
 259                 }
 260               k--;
 261               for(l=i; l<=k; l++)
 262                 {
 263                   pStream->insert(leftVerts[l]);
 264                 }
 265
 266               pStream->end(PRIMITIVE_STREAM_FAN);
 267               //update i for next loop
 268               i = k+1;
 269               topMostV = leftVerts[k];
 270
 271             }
 272           else /*leftVerts[i][1] < rightVerts[j][1]*/
 273             {
 274               pStream->begin();
 275               pStream->insert(leftVerts[i]);/*the origion of this fan*/
 276
 277               /*find the last k>=j such that
 278                *rightverts[k][1] > leftverts[i][1]*/
 279               k=j;
 280               while(k< n_right)
 281                 {
 282                   if(rightVerts[k][1] <= leftVerts[i][1])
 283                     break;
 284                   k++;
 285                 }
 286               k--;
 287
 288               for(l=k; l>= j; l--)
 289                 pStream->insert(rightVerts[l]);
 290
 291               pStream->insert(topMostV);
 292               pStream->end(PRIMITIVE_STREAM_FAN);
 293               j=k+1;
 294               topMostV = rightVerts[j-1];
 295             }
 296         }
 297     }
 298 }
 299
 300 static int chainConvex(vertexArray* inc_chain, Int inc_current, Int inc_end)
 301 {
 302   Int i;
 303   //if there are no more than 2 vertices, return 1
 304   if(inc_current >= inc_end-1) return 1;
 305   for(i=inc_current; i<= inc_end-2; i++)
 306     {
 307       if(area(inc_chain->getVertex(i), inc_chain->getVertex(i+1), inc_chain->getVertex(i+2)) <0)
 308         return 0;
 309     }
 310   return 1;
 311 }
 312
 313 static int chainConcave(vertexArray* dec_chain, Int dec_current, Int dec_end)
 314 {
 315   Int i;
 316   //if there are no more than 2 vertices, return 1
 317   if(dec_current >= dec_end -1) return 1;
 318   for(i=dec_current; i<=dec_end-2; i++)
 319     {
 320       if(area(dec_chain->getVertex(i), dec_chain->getVertex(i+1), dec_chain->getVertex(i+2)) >0)
 321         return 0;
 322     }
 323   return 1;
 324 }
 325
 326 void monoTriangulationRecGenInU(Real* topVertex, Real* botVertex,
 327                                 vertexArray* inc_chain, Int inc_current, Int inc_end,
 328                                 vertexArray* dec_chain, Int dec_current, Int dec_end,
 329                                 primStream* pStream)
 330 {
 331
 332 }
 333
 334 void  monoTriangulationRecGenOpt(Real* topVertex, Real* botVertex,
 335                                  vertexArray* inc_chain, Int inc_current, Int inc_end,
 336                                  vertexArray* dec_chain, Int dec_current, Int dec_end,
 337                                  primStream* pStream)
 338 {
 339   Int i;
 340   //copy this to a polygon: directedLine Lioop
 341   sampledLine* sline;
 342   directedLine* dline;
 343   directedLine* poly;
 344
 345   if(inc_current <= inc_end) //at least one vertex in inc_chain
 346     {
 347       sline = new sampledLine(topVertex, inc_chain->getVertex(inc_current));
 348       poly = new directedLine(INCREASING, sline);
 349       for(i=inc_current; i<=inc_end-1; i++)
 350         {
 351           sline = new sampledLine(inc_chain->getVertex(i), inc_chain->getVertex(i+1));
 352           dline = new directedLine(INCREASING, sline);
 353           poly->insert(dline);
 354         }
 355       sline = new sampledLine(inc_chain->getVertex(inc_end), botVertex);
 356       dline = new directedLine(INCREASING, sline);
 357       poly->insert(dline);
 358     }
 359   else //inc_chian is empty
 360     {
 361       sline = new sampledLine(topVertex, botVertex);
 362       dline = new directedLine(INCREASING, sline);
 363       poly = dline;
 364     }
 365
 366   assert(poly != NULL);
 367
 368   if(dec_current <= dec_end) //at least on vertex in dec_Chain
 369     {
 370       sline = new sampledLine(botVertex, dec_chain->getVertex(dec_end));
 371       dline = new directedLine(INCREASING, sline);
 372       poly->insert(dline);
 373       for(i=dec_end; i>dec_current; i--)
 374         {
 375           sline = new sampledLine(dec_chain->getVertex(i), dec_chain->getVertex(i-1));
 376           dline = new directedLine(INCREASING, sline);
 377           poly->insert(dline);
 378         }
 379       sline = new sampledLine(dec_chain->getVertex(dec_current), topVertex);
 380       dline = new directedLine(INCREASING, sline);
 381       poly->insert(dline);
 382     }
 383   else //dec_chain  is empty
 384     {
 385       sline = new sampledLine(botVertex, topVertex);
 386       dline = new directedLine(INCREASING, sline);
 387       poly->insert(dline);
 388     }
 389
 390   {
 391     Int n_cusps;
 392     Int n_edges = poly->numEdges();
 393     directedLine** cusps = (directedLine**) malloc(sizeof(directedLine*)*n_edges);
 394     assert(cusps);
 395     findInteriorCuspsX(poly, n_cusps, cusps);
 396
 397     if(n_cusps ==0) //u monotine
 398       {
 399         monoTriangulationFun(poly, compV2InX, pStream);
 400       }
 401     else if(n_cusps == 1) // one interior cusp
 402       {
 403         directedLine* new_polygon = polygonConvert(cusps[0]);
 404         directedLine* other = findDiagonal_singleCuspX(new_polygon);
 405           //<other> should NOT be null unless there are self-intersecting
 406           //trim curves. In that case, we don't want to core dump, instead,
 407           //we triangulate anyway, and print out error message.
 408           if(other == NULL)
 409             {
 410               monoTriangulationFun(poly, compV2InX, pStream);
 411             }
 412           else
 413             {
 414               directedLine* ret_p1;
 415               directedLine* ret_p2;
 416
 417               new_polygon->connectDiagonal_2slines(new_polygon, other,
 418                                                    &ret_p1,
 419                                                    &ret_p2,
 420                                                    new_polygon);
 421
 422               monoTriangulationFun(ret_p1, compV2InX, pStream);
 423               monoTriangulationFun(ret_p2, compV2InX, pStream);
 424
 425               ret_p1->deleteSinglePolygonWithSline();
 426               ret_p2->deleteSinglePolygonWithSline();
 427             }
 428       }
 429     else
 430       {
 431         //we need a general partitionX funtion (supposed to be in partitionX.C,
 432         //not implemented yet. XXX
 433         //monoTriangulationFun(poly, compV2InY, pStream);
 434
 435         directedLine* new_polygon = polygonConvert(poly);
 436         directedLine* list = monoPolyPart(new_polygon);
 437         for(directedLine* temp = list; temp != NULL; temp = temp->getNextPolygon())
 438           {
 439             monoTriangulationFun(temp, compV2InX, pStream);
 440           }
 441         //clean up
 442         list->deletePolygonListWithSline();
 443
 444       }
 445
 446     free(cusps);
 447     /*
 448       if(numInteriorCuspsX(poly) == 0) //is u monotone
 449         monoTriangulationFun(poly, compV2InX, pStream);
 450       else //it is not u motone
 451         monoTriangulationFun(poly, compV2InY, pStream);
 452         */
 453     //clean up space
 454     poly->deleteSinglePolygonWithSline();
 455     return;
 456   }
 457
 458   //apparently the following code is not reachable,
 459   //it is for test purpose
 460   if(inc_current > inc_end || dec_current>dec_end)
 461     {
 462     monoTriangulationRecGen(topVertex, botVertex, inc_chain, inc_current, inc_end,
 463                             dec_chain, dec_current, dec_end,
 464                             pStream);
 465     return;
 466   }
 467
 468
 469   if(
 470      area(dec_chain->getVertex(dec_current),
 471           topVertex,
 472           inc_chain->getVertex(inc_current)) >=0
 473      && chainConvex(inc_chain, inc_current, inc_end)
 474      && chainConcave(dec_chain, dec_current, dec_end)
 475      && area(inc_chain->getVertex(inc_end), botVertex, dec_chain->getVertex(dec_end)) >=0
 476      )
 477     {
 478       monoTriangulationRecFunGen(topVertex, botVertex,
 479                                  inc_chain, inc_current, inc_end,
 480                                  dec_chain, dec_current, dec_end,
 481                                  compV2InX, pStream);
 482     }
 483   else
 484     {
 485       monoTriangulationRecGen(topVertex, botVertex, inc_chain, inc_current, inc_end,
 486                               dec_chain, dec_current, dec_end,
 487                               pStream);
 488     }
 489 }
 490
 491 /*if inc_current>inc_end, then inc_chain has no points to be considered
 492  *same for dec_chain
 493  */
 494 void monoTriangulationRecGen(Real* topVertex, Real* botVertex,
 495                           vertexArray* inc_chain, Int inc_current, Int inc_end,
 496                           vertexArray* dec_chain, Int dec_current, Int dec_end,
 497                           primStream* pStream)
 498 {
 499   Real** inc_array ;
 500   Real** dec_array ;
 501   Int i;
 502
 503   if(inc_current > inc_end && dec_current>dec_end)
 504     return;
 505   else if(inc_current>inc_end) /*no more vertices on inc_chain*/
 506     {
 507       dec_array = dec_chain->getArray();
 508       reflexChain rChain(100,0);
 509       /*put the top vertex into the reflex chain*/
 510       rChain.processNewVertex(topVertex, pStream);
 511       /*process all the vertices on the dec_chain*/
 512       for(i=dec_current; i<=dec_end; i++){
 513         rChain.processNewVertex(dec_array[i], pStream);
 514       }
 515       /*process the bottom vertex*/
 516       rChain.processNewVertex(botVertex, pStream);
 517     }
 518   else if(dec_current> dec_end) /*no more vertices on dec_chain*/
 519     {
 520       inc_array = inc_chain->getArray();
 521
 522       reflexChain rChain(100,1);
 523       /*put the top vertex into the reflex chain*/
 524       rChain.processNewVertex(topVertex, pStream);
 525       /*process all the vertices on the inc_chain*/
 526       for(i=inc_current; i<=inc_end; i++){
 527         rChain.processNewVertex(inc_array[i], pStream);
 528       }
 529       /*process the bottom vertex*/
 530       rChain.processNewVertex(botVertex, pStream);
 531     }
 532   else /*neither chain is empty*/
 533     {
 534       inc_array = inc_chain -> getArray();
 535       dec_array = dec_chain -> getArray();
 536
 537       /*if top of inc_chain is 'lower' than top of dec_chain, process all the
 538        *vertices on the dec_chain which are higher than top of inc_chain
 539        */
 540       if(compV2InY(inc_array[inc_current], dec_array[dec_current]) <= 0)
 541         {
 542
 543           reflexChain rChain(100, 0);
 544           rChain.processNewVertex(topVertex, pStream);
 545           for(i=dec_current; i<=dec_end; i++)
 546             {
 547               if(compV2InY(inc_array[inc_current], dec_array[i]) <= 0)
 548                 rChain.processNewVertex(dec_array[i], pStream);
 549               else
 550                 break;
 551             }
 552           rChain.outputFan(inc_array[inc_current], pStream);
 553           monoTriangulationRecGen(dec_array[i-1], botVertex,
 554                                inc_chain, inc_current, inc_end,
 555                                dec_chain, i, dec_end,
 556                                pStream);
 557         }
 558       else /*compV2InY(inc_array[inc_current], dec_array[dec_current]) > 0*/
 559         {
 560
 561           reflexChain rChain(100, 1);
 562           rChain.processNewVertex(topVertex, pStream);
 563           for(i=inc_current; i<=inc_end; i++)
 564             {
 565               if(compV2InY(inc_array[i], dec_array[dec_current]) >0)
 566                 rChain.processNewVertex(inc_array[i], pStream);
 567               else
 568                 break;
 569             }
 570           rChain.outputFan(dec_array[dec_current], pStream);
 571           monoTriangulationRecGen(inc_array[i-1], botVertex,
 572                                inc_chain, i, inc_end,
 573                                dec_chain, dec_current,dec_end,
 574                                pStream);
 575         }
 576     }/*end case neither is empty*/
 577 }
 578
 579 void monoTriangulationFun(directedLine* monoPolygon, Int (*compFun)(Real*, Real*), primStream* pStream)
 580 {
 581   Int i;
 582   /*find the top vertex, bottom vertex, inccreasing chain, and decreasing chain,
 583    *then call monoTriangulationRec
 584    */
 585   directedLine* tempV;
 586   directedLine* topV;
 587   directedLine* botV;
 588   topV = botV = monoPolygon;
 589   for(tempV = monoPolygon->getNext(); tempV != monoPolygon; tempV = tempV->getNext())
 590     {
 591       if(compFun(topV->head(), tempV->head())<0) {
 592         topV = tempV;
 593       }
 594       if(compFun(botV->head(), tempV->head())>0) {
 595         botV = tempV;
 596       }
 597     }
 598
 599   /*creat increase and decrease chains*/
 600   vertexArray inc_chain(20); /*this is a dynamic array*/
 601   for(i=1; i<=topV->get_npoints()-2; i++) { /*the first vertex is the top vertex which doesn't belong to inc_chain*/
 602     inc_chain.appendVertex(topV->getVertex(i));
 603   }
 604   for(tempV = topV->getNext(); tempV != botV; tempV = tempV->getNext())
 605     {
 606       for(i=0; i<=tempV->get_npoints()-2; i++){
 607         inc_chain.appendVertex(tempV->getVertex(i));
 608       }
 609     }
 610
 611   vertexArray dec_chain(20);
 612   for(tempV = topV->getPrev(); tempV != botV; tempV = tempV->getPrev())
 613     {
 614       for(i=tempV->get_npoints()-2; i>=0; i--){
 615         dec_chain.appendVertex(tempV->getVertex(i));
 616       }
 617     }
 618   for(i=botV->get_npoints()-2; i>=1; i--){
 619     dec_chain.appendVertex(tempV->getVertex(i));
 620   }
 621
 622   monoTriangulationRecFun(topV->head(), botV->head(), &inc_chain, 0, &dec_chain, 0, compFun, pStream);
 623
 624 }
 625
 626 void monoTriangulation(directedLine* monoPolygon, primStream* pStream)
 627 {
 628   Int i;
 629   /*find the top vertex, bottom vertex, inccreasing chain, and decreasing chain,
 630    *then call monoTriangulationRec
 631    */
 632   directedLine* tempV;
 633   directedLine* topV;
 634   directedLine* botV;
 635   topV = botV = monoPolygon;
 636   for(tempV = monoPolygon->getNext(); tempV != monoPolygon; tempV = tempV->getNext())
 637     {
 638       if(compV2InY(topV->head(), tempV->head())<0) {
 639         topV = tempV;
 640       }
 641       if(compV2InY(botV->head(), tempV->head())>0) {
 642         botV = tempV;
 643       }
 644     }
 645   /*creat increase and decrease chains*/
 646   vertexArray inc_chain(20); /*this is a dynamic array*/
 647   for(i=1; i<=topV->get_npoints()-2; i++) { /*the first vertex is the top vertex which doesn't belong to inc_chain*/
 648     inc_chain.appendVertex(topV->getVertex(i));
 649   }
 650   for(tempV = topV->getNext(); tempV != botV; tempV = tempV->getNext())
 651     {
 652       for(i=0; i<=tempV->get_npoints()-2; i++){
 653         inc_chain.appendVertex(tempV->getVertex(i));
 654       }
 655     }
 656
 657   vertexArray dec_chain(20);
 658   for(tempV = topV->getPrev(); tempV != botV; tempV = tempV->getPrev())
 659     {
 660       for(i=tempV->get_npoints()-2; i>=0; i--){
 661         dec_chain.appendVertex(tempV->getVertex(i));
 662       }
 663     }
 664   for(i=botV->get_npoints()-2; i>=1; i--){
 665     dec_chain.appendVertex(tempV->getVertex(i));
 666   }
 667
 668   monoTriangulationRec(topV->head(), botV->head(), &inc_chain, 0, &dec_chain, 0, pStream);
 669
 670 }
 671
 672 /*the chain could be increasing or decreasing, although we use the
 673  * name inc_chain.
 674  *the argument  is_increase_chain indicates whether this chain
 675  *is increasing (left chain in V-monotone case) or decreaing (right chain
 676  *in V-monotone case).
 677  */
 678 void monoTriangulation2(Real* topVertex, Real* botVertex,
 679                         vertexArray* inc_chain, Int inc_smallIndex,
 680                         Int inc_largeIndex,
 681                         Int is_increase_chain,
 682                         primStream* pStream)
 683 {
 684   assert( inc_chain != NULL);
 685   Real** inc_array ;
 686
 687   if(inc_smallIndex > inc_largeIndex)
 688     return; //no triangles
 689   if(inc_smallIndex == inc_largeIndex)
 690     {
 691       if(is_increase_chain)
 692         pStream->triangle(inc_chain->getVertex(inc_smallIndex), botVertex, topVertex);
 693       else
 694         pStream->triangle(inc_chain->getVertex(inc_smallIndex), topVertex, botVertex);
 695       return;
 696     }
 697   Int i;
 698
 699   if(is_increase_chain && botVertex[1] == inc_chain->getVertex(inc_largeIndex)[1])
 700     {
 701       pStream->triangle(botVertex, inc_chain->getVertex(inc_largeIndex-1),
 702                         inc_chain->getVertex(inc_largeIndex));
 703       monoTriangulation2(topVertex, botVertex, inc_chain, inc_smallIndex,
 704                          inc_largeIndex-1,
 705                          is_increase_chain,
 706                          pStream);
 707       return;
 708     }
 709   else if( (!is_increase_chain) && topVertex[1] == inc_chain->getVertex(inc_smallIndex)[1])
 710     {
 711       pStream->triangle(topVertex, inc_chain->getVertex(inc_smallIndex+1),
 712                         inc_chain->getVertex(inc_smallIndex));
 713       monoTriangulation2(topVertex, botVertex, inc_chain, inc_smallIndex+1,
 714                          inc_largeIndex, is_increase_chain, pStream);
 715       return ;
 716     }
 717
 718   inc_array = inc_chain->getArray();
 719
 720   reflexChain rChain(20,is_increase_chain); /*1 means the chain is increasing*/
 721
 722   rChain.processNewVertex(topVertex, pStream);
 723
 724   for(i=inc_smallIndex; i<=inc_largeIndex; i++){
 725     rChain.processNewVertex(inc_array[i], pStream);
 726   }
 727   rChain.processNewVertex(botVertex, pStream);
 728
 729 }
 730
 731 /*if compFun == compV2InY, top to bottom: V-monotone
 732  *if compFun == compV2InX, right to left: U-monotone
 733  */
 734 void monoTriangulationRecFunGen(Real* topVertex, Real* botVertex,
 735                           vertexArray* inc_chain, Int inc_current, Int inc_end,
 736                           vertexArray* dec_chain, Int dec_current, Int dec_end,
 737                           Int  (*compFun)(Real*, Real*),
 738                           primStream* pStream)
 739 {
 740   assert( inc_chain != NULL && dec_chain != NULL);
 741   assert( ! (inc_current> inc_end &&
 742              dec_current> dec_end));
 743   /*
 744   Int inc_nVertices;
 745   Int dec_nVertices;
 746   */
 747   Real** inc_array ;
 748   Real** dec_array ;
 749   Int i;
 750   assert( ! ( (inc_chain==NULL) && (dec_chain==NULL)));
 751
 752   if(inc_current> inc_end) /*no more vertices on inc_chain*/
 753     {
 754
 755       dec_array = dec_chain->getArray();
 756       reflexChain rChain(20,0);
 757       /*put the top vertex into the reflex chain*/
 758       rChain.processNewVertex(topVertex, pStream);
 759       /*process all the vertices on the dec_chain*/
 760       for(i=dec_current; i<=dec_end; i++){
 761         rChain.processNewVertex(dec_array[i], pStream);
 762       }
 763       /*process the bottom vertex*/
 764       rChain.processNewVertex(botVertex, pStream);
 765
 766     }
 767   else if(dec_current> dec_end) /*no more vertices on dec_chain*/
 768     {
 769       inc_array = inc_chain->getArray();
 770       reflexChain rChain(20,1);
 771       /*put the top vertex into the reflex chain*/
 772       rChain.processNewVertex(topVertex, pStream);
 773       /*process all the vertices on the inc_chain*/
 774       for(i=inc_current; i<=inc_end; i++){
 775         rChain.processNewVertex(inc_array[i], pStream);
 776       }
 777       /*process the bottom vertex*/
 778       rChain.processNewVertex(botVertex, pStream);
 779     }
 780   else /*neither chain is empty*/
 781     {
 782       inc_array = inc_chain -> getArray();
 783       dec_array = dec_chain -> getArray();
 784
 785       /*if top of inc_chain is 'lower' than top of dec_chain, process all the
 786        *vertices on the dec_chain which are higher than top of inc_chain
 787        */
 788       if(compFun(inc_array[inc_current], dec_array[dec_current]) <= 0)
 789         {
 790
 791           reflexChain rChain(20, 0);
 792           rChain.processNewVertex(topVertex, pStream);
 793           for(i=dec_current; i<=dec_end; i++)
 794             {
 795               if(compFun(inc_array[inc_current], dec_array[i]) <= 0)
 796                 rChain.processNewVertex(dec_array[i], pStream);
 797               else
 798                 break;
 799             }
 800           rChain.outputFan(inc_array[inc_current], pStream);
 801           monoTriangulationRecFunGen(dec_array[i-1], botVertex,
 802                                inc_chain, inc_current, inc_end,
 803                                dec_chain, i, dec_end,
 804                                compFun,
 805                                pStream);
 806         }
 807       else /*compFun(inc_array[inc_current], dec_array[dec_current]) > 0*/
 808         {
 809
 810           reflexChain rChain(20, 1);
 811           rChain.processNewVertex(topVertex, pStream);
 812           for(i=inc_current; i<=inc_end; i++)
 813             {
 814               if(compFun(inc_array[i], dec_array[dec_current]) >0)
 815                 rChain.processNewVertex(inc_array[i], pStream);
 816               else
 817                 break;
 818             }
 819           rChain.outputFan(dec_array[dec_current], pStream);
 820           monoTriangulationRecFunGen(inc_array[i-1], botVertex,
 821                                inc_chain, i,inc_end,
 822                                dec_chain, dec_current,dec_end,
 823                                compFun,
 824                                pStream);
 825         }
 826     }/*end case neither is empty*/
 827 }
 828
 829 /*if compFun == compV2InY, top to bottom: V-monotone
 830  *if compFun == compV2InX, right to left: U-monotone
 831  */
 832 void monoTriangulationRecFun(Real* topVertex, Real* botVertex,
 833                           vertexArray* inc_chain, Int inc_current,
 834                           vertexArray* dec_chain, Int dec_current,
 835                           Int  (*compFun)(Real*, Real*),
 836                           primStream* pStream)
 837 {
 838   assert( inc_chain != NULL && dec_chain != NULL);
 839   assert( ! (inc_current>=inc_chain->getNumElements() &&
 840              dec_current>=dec_chain->getNumElements()));
 841   Int inc_nVertices;
 842   Int dec_nVertices;
 843   Real** inc_array ;
 844   Real** dec_array ;
 845   Int i;
 846   assert( ! ( (inc_chain==NULL) && (dec_chain==NULL)));
 847
 848   if(inc_current>=inc_chain->getNumElements()) /*no more vertices on inc_chain*/
 849     {
 850
 851       dec_array = dec_chain->getArray();
 852       dec_nVertices = dec_chain->getNumElements();
 853       reflexChain rChain(20,0);
 854       /*put the top vertex into the reflex chain*/
 855       rChain.processNewVertex(topVertex, pStream);
 856       /*process all the vertices on the dec_chain*/
 857       for(i=dec_current; i<dec_nVertices; i++){
 858         rChain.processNewVertex(dec_array[i], pStream);
 859       }
 860       /*process the bottom vertex*/
 861       rChain.processNewVertex(botVertex, pStream);
 862
 863     }
 864   else if(dec_current>= dec_chain->getNumElements()) /*no more vertices on dec_chain*/
 865     {
 866       inc_array = inc_chain->getArray();
 867       inc_nVertices= inc_chain->getNumElements();
 868       reflexChain rChain(20,1);
 869       /*put the top vertex into the reflex chain*/
 870       rChain.processNewVertex(topVertex, pStream);
 871       /*process all the vertices on the inc_chain*/
 872       for(i=inc_current; i<inc_nVertices; i++){
 873         rChain.processNewVertex(inc_array[i], pStream);
 874       }
 875       /*process the bottom vertex*/
 876       rChain.processNewVertex(botVertex, pStream);
 877     }
 878   else /*neither chain is empty*/
 879     {
 880       inc_array = inc_chain -> getArray();
 881       dec_array = dec_chain -> getArray();
 882       inc_nVertices= inc_chain->getNumElements();
 883       dec_nVertices= dec_chain->getNumElements();
 884       /*if top of inc_chain is 'lower' than top of dec_chain, process all the
 885        *vertices on the dec_chain which are higher than top of inc_chain
 886        */
 887       if(compFun(inc_array[inc_current], dec_array[dec_current]) <= 0)
 888         {
 889
 890           reflexChain rChain(20, 0);
 891           rChain.processNewVertex(topVertex, pStream);
 892           for(i=dec_current; i<dec_nVertices; i++)
 893             {
 894               if(compFun(inc_array[inc_current], dec_array[i]) <= 0)
 895                 rChain.processNewVertex(dec_array[i], pStream);
 896               else
 897                 break;
 898             }
 899           rChain.outputFan(inc_array[inc_current], pStream);
 900           monoTriangulationRecFun(dec_array[i-1], botVertex,
 901                                inc_chain, inc_current,
 902                                dec_chain, i,
 903                                compFun,
 904                                pStream);
 905         }
 906       else /*compFun(inc_array[inc_current], dec_array[dec_current]) > 0*/
 907         {
 908
 909           reflexChain rChain(20, 1);
 910           rChain.processNewVertex(topVertex, pStream);
 911           for(i=inc_current; i<inc_nVertices; i++)
 912             {
 913               if(compFun(inc_array[i], dec_array[dec_current]) >0)
 914                 rChain.processNewVertex(inc_array[i], pStream);
 915               else
 916                 break;
 917             }
 918           rChain.outputFan(dec_array[dec_current], pStream);
 919           monoTriangulationRecFun(inc_array[i-1], botVertex,
 920                                inc_chain, i,
 921                                dec_chain, dec_current,
 922                                compFun,
 923                                pStream);
 924         }
 925     }/*end case neither is empty*/
 926 }
 927
 928
 929 void monoTriangulationRec(Real* topVertex, Real* botVertex,
 930                           vertexArray* inc_chain, Int inc_current,
 931                           vertexArray* dec_chain, Int dec_current,
 932                           primStream* pStream)
 933 {
 934   assert( inc_chain != NULL && dec_chain != NULL);
 935   assert( ! (inc_current>=inc_chain->getNumElements() &&
 936              dec_current>=dec_chain->getNumElements()));
 937   Int inc_nVertices;
 938   Int dec_nVertices;
 939   Real** inc_array ;
 940   Real** dec_array ;
 941   Int i;
 942   assert( ! ( (inc_chain==NULL) && (dec_chain==NULL)));
 943
 944   if(inc_current>=inc_chain->getNumElements()) /*no more vertices on inc_chain*/
 945     {
 946
 947       dec_array = dec_chain->getArray();
 948       dec_nVertices = dec_chain->getNumElements();
 949       reflexChain rChain(20,0);
 950       /*put the top vertex into the reflex chain*/
 951       rChain.processNewVertex(topVertex, pStream);
 952       /*process all the vertices on the dec_chain*/
 953       for(i=dec_current; i<dec_nVertices; i++){
 954         rChain.processNewVertex(dec_array[i], pStream);
 955       }
 956       /*process the bottom vertex*/
 957       rChain.processNewVertex(botVertex, pStream);
 958
 959     }
 960   else if(dec_current>= dec_chain->getNumElements()) /*no more vertices on dec_chain*/
 961     {
 962       inc_array = inc_chain->getArray();
 963       inc_nVertices= inc_chain->getNumElements();
 964       reflexChain rChain(20,1);
 965       /*put the top vertex into the reflex chain*/
 966       rChain.processNewVertex(topVertex, pStream);
 967       /*process all the vertices on the inc_chain*/
 968       for(i=inc_current; i<inc_nVertices; i++){
 969         rChain.processNewVertex(inc_array[i], pStream);
 970       }
 971       /*process the bottom vertex*/
 972       rChain.processNewVertex(botVertex, pStream);
 973     }
 974   else /*neither chain is empty*/
 975     {
 976       inc_array = inc_chain -> getArray();
 977       dec_array = dec_chain -> getArray();
 978       inc_nVertices= inc_chain->getNumElements();
 979       dec_nVertices= dec_chain->getNumElements();
 980       /*if top of inc_chain is 'lower' than top of dec_chain, process all the
 981        *vertices on the dec_chain which are higher than top of inc_chain
 982        */
 983       if(compV2InY(inc_array[inc_current], dec_array[dec_current]) <= 0)
 984         {
 985
 986           reflexChain rChain(20, 0);
 987           rChain.processNewVertex(topVertex, pStream);
 988           for(i=dec_current; i<dec_nVertices; i++)
 989             {
 990               if(compV2InY(inc_array[inc_current], dec_array[i]) <= 0)
 991                 rChain.processNewVertex(dec_array[i], pStream);
 992               else
 993                 break;
 994             }
 995           rChain.outputFan(inc_array[inc_current], pStream);
 996           monoTriangulationRec(dec_array[i-1], botVertex,
 997                                inc_chain, inc_current,
 998                                dec_chain, i,
 999                                pStream);
1000         }
1001       else /*compV2InY(inc_array[inc_current], dec_array[dec_current]) > 0*/
1002         {
1003
1004           reflexChain rChain(20, 1);
1005           rChain.processNewVertex(topVertex, pStream);
1006           for(i=inc_current; i<inc_nVertices; i++)
1007             {
1008               if(compV2InY(inc_array[i], dec_array[dec_current]) >0)
1009                 rChain.processNewVertex(inc_array[i], pStream);
1010               else
1011                 break;
1012             }
1013           rChain.outputFan(dec_array[dec_current], pStream);
1014           monoTriangulationRec(inc_array[i-1], botVertex,
1015                                inc_chain, i,
1016                                dec_chain, dec_current,
1017                                pStream);
1018         }
1019     }/*end case neither is empty*/
1020 }
1021
1022
1023
1024 /* the name here assumes that the polygon is Y-monotone, but
1025  *this function also works for X-monotone polygons.
1026  * a monotne polygon consists of two extrem verteices: topVertex and botVertex, and
1027  *two monotone chains: inc_chain, and dec_chain. The edges of the increasing chain (inc_chain)
1028  *is ordered by following pointer: next, while  the edges of the decreasing chain (dec_chain)
1029  *is ordered by following pointer: prev
1030  * inc_index index the vertex which is the toppest of the inc_chain which we are handling currently.
1031  * dec_index index the vertex which is the toppest of the dec_chain which we are handling currently.
1032  */
1033 void monoTriangulationRec(directedLine* inc_chain, Int inc_index,
1034                           directedLine* dec_chain, Int dec_index,
1035                           directedLine* topVertex, Int top_index,
1036                           directedLine* botVertex,
1037                           primStream* pStream)
1038 {
1039   Int i;
1040   directedLine *temp, *oldtemp = NULL;
1041   Int tempIndex, oldtempIndex = 0;
1042
1043   assert(inc_chain != NULL && dec_chain != NULL);
1044
1045   if(inc_chain == botVertex) {
1046     reflexChain rChain(20, 0);
1047     rChain.processNewVertex(topVertex->getVertex(top_index), pStream);
1048     for(i=dec_index; i< dec_chain->get_npoints(); i++){
1049       rChain.processNewVertex(dec_chain->getVertex(i), pStream);
1050     }
1051     for(temp = dec_chain->getPrev(); temp != botVertex; temp = temp->getPrev())
1052       {
1053         for(i=0; i<temp->get_npoints(); i++){
1054           rChain.processNewVertex(temp->getVertex(i), pStream);
1055         }
1056       }
1057   }
1058   else if(dec_chain==botVertex) {
1059     reflexChain rChain(20, 1);
1060     rChain.processNewVertex(topVertex->getVertex(top_index), pStream);
1061     for(i=inc_index; i< inc_chain->get_npoints(); i++){
1062       rChain.processNewVertex(inc_chain->getVertex(i), pStream);
1063     }
1064     for(temp = inc_chain->getPrev(); temp != botVertex; temp = temp->getNext())
1065       {
1066         for(i=0; i<temp->get_npoints(); i++){
1067           rChain.processNewVertex(temp->getVertex(i), pStream);
1068         }
1069       }
1070   }
1071   else /*neither reached the bottom*/{
1072     if(compV2InY(inc_chain->getVertex(inc_index), dec_chain->getVertex(dec_index)) <=0) {
1073       reflexChain rChain(20, 0);
1074       rChain.processNewVertex(topVertex -> getVertex(top_index), pStream);
1075       temp = dec_chain;
1076       tempIndex = dec_index;
1077       while( compV2InY(inc_chain->getVertex(inc_index), temp->getVertex(tempIndex))<=0) {
1078         oldtemp = temp;
1079         oldtempIndex = tempIndex;
1080         rChain.processNewVertex(temp->getVertex(tempIndex), pStream);
1081
1082         if(tempIndex == temp->get_npoints()-1){
1083           tempIndex = 0;
1084           temp = temp->getPrev();
1085         }
1086         else{
1087           tempIndex++;
1088         }
1089       }
1090       rChain.outputFan(inc_chain->getVertex(inc_index), pStream);
1091       monoTriangulationRec(inc_chain, inc_index, temp, tempIndex, oldtemp, oldtempIndex, botVertex, pStream);
1092     }
1093     else /* >0*/ {
1094       reflexChain rChain(20, 1);
1095       rChain.processNewVertex(topVertex -> getVertex(top_index), pStream);
1096       temp = inc_chain;
1097       tempIndex = inc_index;
1098       while( compV2InY(temp->getVertex(tempIndex), dec_chain->getVertex(dec_index))>0){
1099         oldtemp = temp;
1100         oldtempIndex = tempIndex;
1101         rChain.processNewVertex(temp->getVertex(tempIndex), pStream);
1102
1103         if(tempIndex == temp->get_npoints()-1){
1104           tempIndex = 0;
1105           temp = temp->getNext();
1106         }
1107         else{
1108           tempIndex++;
1109         }
1110       }
1111       rChain.outputFan(dec_chain->getVertex(dec_index), pStream);
1112       monoTriangulationRec(temp, tempIndex, dec_chain, dec_index, oldtemp, oldtempIndex, botVertex, pStream);
1113     }
1114   } /*end case neither reached the bottom*/
1115 }
1116
1117 /***************************vertexArray begin here**********************************/
1118 vertexArray::vertexArray(Real2* vertices, Int nVertices)
1119 {
1120   Int i;
1121   size = index = nVertices;
1122   array = (Real**) malloc(sizeof(Real*) * nVertices);
1123   assert(array);
1124   for(i=0; i<nVertices; i++)
1125     {
1126       array[i] = vertices[i];
1127       array[i] = vertices[i];
1128     }
1129 }
1130
1131 vertexArray::vertexArray(Int s)
1132 {
1133   size = s;
1134   array = (Real**) malloc(sizeof(Real*) * s);
1135   assert(array);
1136   index = 0;
1137 }
1138
1139 vertexArray::~vertexArray()
1140 {
1141   free(array);
1142 }
1143
1144 void vertexArray::appendVertex(Real* ptr)
1145 {
1146   Int i;
1147   if(index >= size){
1148     Real** temp = (Real**) malloc(sizeof(Real*) * (2*size +1));
1149     assert(temp);
1150     for(i=0; i<index; i++)
1151       temp[i] = array[i];
1152     free(array);
1153     array = temp;
1154     size = 2*size+1;
1155   }
1156   array[index++] = ptr;
1157 }
1158
1159 void vertexArray::print()
1160 {
1161   printf("vertex Array:index=%i, size=%i\n", index, size);
1162   for(Int i=0; i<index; i++)
1163     {
1164       printf("(%f,%f) ", array[i][0], array[i][1]);
1165     }
1166   printf("\n");
1167 }
1168
1169 /*find the first i such that array[i][1] >= v
1170  * and array[i+1][1] <v
1171  * if index == 0 (the array is empty, return -1.
1172  * if v is above all, return -1.
1173  * if v is below all, return index-1.
1174  */
1175 Int vertexArray::findIndexAbove(Real v)
1176 {
1177   Int i;
1178   if(index == 0)
1179     return -1;
1180   else if(array[0][1] < v)
1181     return -1;
1182   else
1183     {
1184       for(i=1; i<index; i++)
1185         {
1186           if(array[i][1] < v)
1187             break;
1188         }
1189       return i-1;
1190     }
1191 }
1192
1193 /*find the first i<=endIndex such that array[i][1] <= v
1194  * and array[i-1][1] > v
1195  *if sartIndex>endIndex, then return endIndex+1.
1196  *otherwise, startIndex<=endIndex, it is assumed that
1197  * 0<=startIndex<=endIndex<index.
1198  * if v is below all, return endIndex+1
1199  * if v is above all, return startIndex.
1200  */
1201 Int vertexArray::findIndexBelowGen(Real v, Int startIndex, Int endIndex)
1202 {
1203   Int i;
1204   if(startIndex > endIndex)
1205     return endIndex+1;
1206   else if(array[endIndex][1] >  v)
1207     return endIndex+1;
1208   else //now array[endIndex][1] <= v
1209     {
1210       for(i=endIndex-1; i>=startIndex; i--)
1211         {
1212           if(array[i][1] > v)
1213             break;
1214         }
1215       return i+1;
1216     }
1217 }
1218
1219 /*find the first i<=endIndex such that array[i-1][1] >= v
1220  * and array[i][1] < v
1221  *if sartIndex>endIndex, then return endIndex+1.
1222  *otherwise, startIndex<=endIndex, it is assumed that
1223  * 0<=startIndex<=endIndex<index.
1224  * if v is below or equal to all, return endIndex+1
1225  * if v is strictly above all, return startIndex.
1226  */
1227 Int vertexArray::findIndexStrictBelowGen(Real v, Int startIndex, Int endIndex)
1228 {
1229   Int i;
1230   if(startIndex > endIndex)
1231     return endIndex+1;
1232   else if(array[endIndex][1] >=  v)
1233     return endIndex+1;
1234   else //now array[endIndex][1] < v
1235     {
1236       for(i=endIndex-1; i>=startIndex; i--)
1237         {
1238           if(array[i][1] >= v)
1239             break;
1240         }
1241       return i+1;
1242     }
1243 }
1244
1245 /*find the first i>startIndex such that array[i-1][1] > v
1246  * and array[i][1] >=v
1247  *if sartIndex>endIndex, then return startIndex-1.
1248  *otherwise, startIndex<=endIndex, it is assumed that
1249  * 0<=startIndex<=endIndex<index.
1250  * if v is strictly above all, return startIndex-1
1251  * if v is strictly below all, return endIndex.
1252  */
1253 Int vertexArray::findIndexFirstAboveEqualGen(Real v, Int startIndex, Int endIndex)
1254 {
1255
1256   Int i;
1257   if(startIndex > endIndex)
1258     return startIndex-1;
1259   else if(array[startIndex][1] < v)
1260     return startIndex-1;
1261   else //now array[startIndex][1] >= v
1262     {
1263
1264       for(i=startIndex; i<=endIndex; i++)
1265         {
1266           if(array[i][1] <= v)
1267             break;
1268         }
1269       if(i>endIndex) // v is strictly below all
1270         return endIndex;
1271       else if(array[i][1] == v)
1272         return i;
1273       else
1274         return i-1;
1275     }
1276
1277 }
1278
1279
1280 /*find the first i>=startIndex such that array[i][1] >= v
1281  * and array[i+1][1] <v
1282  *if sartIndex>endIndex, then return startIndex-1.
1283  *otherwise, startIndex<=endIndex, it is assumed that
1284  * 0<=startIndex<=endIndex<index.
1285  * if v is above all, return startIndex-1
1286  * if v is below all, return endIndex.
1287  */
1288 Int vertexArray::findIndexAboveGen(Real v, Int startIndex, Int endIndex)
1289 {
1290   Int i;
1291   if(startIndex > endIndex)
1292     return startIndex-1;
1293   else if(array[startIndex][1] < v)
1294     return startIndex-1;
1295   else //now array[startIndex][1] >= v
1296     {
1297       for(i=startIndex+1; i<=endIndex; i++)
1298         {
1299           if(array[i][1] < v)
1300             break;
1301         }
1302       return i-1;
1303     }
1304 }
1305
1306 Int vertexArray::findDecreaseChainFromEnd(Int begin, Int end)
1307 {
1308   Int i = end;
1309   Real prevU = array[i][0];
1310   Real thisU;
1311   for(i=end-1; i>=begin; i--){
1312     thisU = array[i][0];
1313     if(thisU < prevU)
1314       prevU = thisU;
1315     else
1316       break;
1317   }
1318   return i;
1319 }
1320
1321 //if(V(start) == v, return start, other wise return the
1322 //last i so that V(i)==v
1323 Int vertexArray::skipEqualityFromStart(Real v, Int start, Int end)
1324 {
1325   Int i;
1326   if(array[start][1] != v)
1327     return start;
1328   //now array[start][1] == v
1329   for(i=start+1; i<= end; i++)
1330     if(array[i][1] != v)
1331       break;
1332   return i-1;
1333 }
1334
1335
1336 /***************************vertexArray end****************************************/
1337
1338
1339
1340 /***************************relfex chain stuff begin here*****************************/
1341
1342 reflexChain::reflexChain(Int size, Int is_increasing)
1343 {
1344   queue = (Real2*) malloc(sizeof(Real2) * size);
1345   assert(queue);
1346   index_queue = 0;
1347   size_queue = size;
1348   isIncreasing = is_increasing;
1349 }
1350
1351 reflexChain::~reflexChain()
1352 {
1353   free(queue);
1354 }
1355
1356 /*put (u,v) at the end of the queue
1357  *pay attention to space
1358  */
1359 void reflexChain::insert(Real u, Real v)
1360 {
1361   Int i;
1362   if(index_queue >= size_queue) {
1363     Real2 *temp = (Real2*) malloc(sizeof(Real2) * (2*size_queue+1));
1364     assert(temp);
1365
1366     /*copy*/
1367     for(i=0; i<index_queue; i++){
1368       temp[i][0] = queue[i][0];
1369       temp[i][1] = queue[i][1];
1370     }
1371
1372     free(queue);
1373     queue = temp;
1374     size_queue = 2*size_queue + 1;
1375   }
1376
1377   queue[index_queue][0] = u;
1378   queue[index_queue][1] = v;
1379   index_queue ++;
1380 }
1381
1382 void reflexChain::insert(Real v[2])
1383 {
1384   insert(v[0], v[1]);
1385 }
1386
1387 /*
1388 static Real area(Real A[2], Real B[2], Real C[2])
1389 {
1390   Real Bx, By, Cx, Cy;
1391   Bx = B[0] - A[0];
1392   By = B[1] - A[1];
1393   Cx = C[0] - A[0];
1394   Cy = C[1] - A[1];
1395   return Bx*Cy - Cx*By;
1396 }
1397 */
1398
1399 /*the chain is reflex, and the vertex v is
1400  *on the other side of the chain, so that
1401  *we can outout the fan with v as the
1402  *the center
1403  */
1404 void reflexChain::outputFan(Real v[2], primStream* pStream)
1405 {
1406   Int i;
1407   pStream->begin();
1408   pStream->insert(v);
1409   if(isIncreasing) {
1410     for(i=0; i<index_queue; i++)
1411       pStream->insert(queue[i]);
1412   }
1413   else {
1414     for(i=index_queue-1; i>=0; i--)
1415       pStream->insert(queue[i]);
1416   }
1417   pStream->end(PRIMITIVE_STREAM_FAN);
1418 }
1419
1420 void reflexChain::processNewVertex(Real v[2], primStream* pStream)
1421 {
1422   Int i,j,k;
1423   Int isReflex;
1424   /*if there are at most one vertex in the queue, then simply insert
1425    */
1426   if(index_queue <=1){
1427     insert(v);
1428     return;
1429   }
1430
1431   /*there are at least two vertices in the queue*/
1432   j=index_queue-1;
1433
1434   for(i=j; i>=1; i--) {
1435     if(isIncreasing) {
1436       isReflex = (area(queue[i-1], queue[i], v) <= 0.0);
1437     }
1438     else /*decreasing*/{
1439       isReflex = (area(v, queue[i], queue[i-1]) <= 0.0);
1440     }
1441     if(isReflex) {
1442       break;
1443     }
1444   }
1445
1446   /*
1447    *if i<j then vertices: i+1--j are convex
1448    * output triangle fan:
1449    *  v, and queue[i], i+1, ..., j
1450    */
1451   if(i<j)
1452     {
1453       pStream->begin();
1454       pStream->insert(v);
1455       if(isIncreasing) {
1456         for(k=i; k<=j; k++)
1457           pStream->insert(queue[k]);
1458       }
1459       else {
1460         for(k=j; k>=i; k--)
1461           pStream->insert(queue[k]);
1462       }
1463
1464       pStream->end(PRIMITIVE_STREAM_FAN);
1465     }
1466
1467   /*delete vertices i+1--j from the queue*/
1468   index_queue = i+1;
1469   /*finally insert v at the end of the queue*/
1470   insert(v);
1471
1472 }
1473
1474 void reflexChain::print()
1475 {
1476   Int i;
1477   printf("reflex chain: isIncreasing=%i\n", isIncreasing);
1478   for(i=0; i<index_queue; i++) {
1479     printf("(%f,%f) ", queue[i][0], queue[i][1]);
1480   }
1481   printf("\n");
1482 }