dll/win32/glu32/libnurbs/internals/intersect.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  * intersect.c++
  37  *
  38  * $Date$ $Revision: 1.1 $
  39  * $Header: /cygdrive/c/RCVS/CVS/ReactOS/reactos/lib/glu32/libnurbs/internals/intersect.cc,v 1.1 2004/02/02 16:39:11 navaraf Exp $
  40  */
  41
  42 #include "glimports.h"
  43 #include "myassert.h"
  44 #include "mystdio.h"
  45 #include "subdivider.h"
  46 #include "arc.h"
  47 #include "bin.h"
  48 #include "backend.h"
  49 #include "trimvertpool.h"
  50
  51 /*#define NOTDEF*/
  52
  53 enum i_result { INTERSECT_VERTEX, INTERSECT_EDGE };
  54
  55 /* local functions */
  56 static int              arc_classify( Arc_ptr, int, REAL );
  57 static enum i_result    pwlarc_intersect( PwlArc *, int, REAL, int, int[3] );
  58
  59
  60 void
  61 Subdivider::partition( Bin & bin, Bin & left, Bin & intersections,
  62                 Bin & right, Bin & unknown, int param, REAL value )
  63 {
  64     Bin headonleft, headonright, tailonleft, tailonright;
  65
  66     for( Arc_ptr jarc = bin.removearc(); jarc; jarc = bin.removearc() ) {
  67
  68         REAL tdiff = jarc->tail()[param] - value;
  69         REAL hdiff = jarc->head()[param] - value;
  70
  71         if( tdiff > 0.0 ) {
  72             if( hdiff > 0.0 ) {
  73                 right.addarc( jarc  );
  74             } else if( hdiff == 0.0 ) {
  75                 tailonright.addarc( jarc  );
  76             } else {
  77                 Arc_ptr jtemp;
  78                 switch( arc_split(jarc, param, value, 0) ) {
  79                     case 2:
  80                         tailonright.addarc( jarc  );
  81                         headonleft.addarc( jarc->next  );
  82                         break;
  83                     case 31:
  84                         assert( jarc->head()[param] > value );
  85                         right.addarc( jarc  );
  86                         tailonright.addarc( jtemp = jarc->next  );
  87                         headonleft.addarc( jtemp->next  );
  88                         break;
  89                     case 32:
  90                         assert( jarc->head()[param] <= value );
  91                         tailonright .addarc( jarc  );
  92                         headonleft.addarc( jtemp = jarc->next  );
  93                         left.addarc( jtemp->next  );
  94                         break;
  95                     case 4:
  96                         right.addarc( jarc  );
  97                         tailonright.addarc( jtemp = jarc->next  );
  98                         headonleft.addarc( jtemp = jtemp->next  );
  99                         left.addarc( jtemp->next  );
 100                 }
 101             }
 102         } else if( tdiff == 0.0 ) {
 103             if( hdiff > 0.0 ) {
 104                 headonright.addarc( jarc  );
 105             } else if( hdiff == 0.0 ) {
 106                 unknown.addarc( jarc  );
 107             } else {
 108                 headonleft.addarc( jarc  );
 109             }
 110         } else {
 111             if( hdiff > 0.0 ) {
 112                 Arc_ptr jtemp;
 113                 switch( arc_split(jarc, param, value, 1) ) {
 114                     case 2:
 115                         tailonleft.addarc( jarc  );
 116                         headonright.addarc( jarc->next  );
 117                         break;
 118                     case 31:
 119                         assert( jarc->head()[param] < value );
 120                         left.addarc( jarc  );
 121                         tailonleft.addarc( jtemp = jarc->next  );
 122                         headonright.addarc( jtemp->next  );
 123                         break;
 124                     case 32:
 125                         assert( jarc->head()[param] >= value );
 126                         tailonleft.addarc( jarc  );
 127                         headonright.addarc( jtemp = jarc->next  );
 128                         right.addarc( jtemp->next  );
 129                         break;
 130                     case 4:
 131                         left.addarc( jarc  );
 132                         tailonleft.addarc( jtemp = jarc->next  );
 133                         headonright.addarc( jtemp = jtemp->next  );
 134                         right.addarc( jtemp->next  );
 135                 }
 136             } else if( hdiff == 0.0 ) {
 137                 tailonleft.addarc( jarc  );
 138             } else {
 139                 left.addarc( jarc  );
 140             }
 141         }
 142     }
 143     if( param == 0 ) {
 144         classify_headonleft_s( headonleft, intersections, left, value );
 145         classify_tailonleft_s( tailonleft, intersections, left, value );
 146         classify_headonright_s( headonright, intersections, right, value );
 147         classify_tailonright_s( tailonright, intersections, right, value );
 148     } else {
 149         classify_headonleft_t( headonleft, intersections, left, value );
 150         classify_tailonleft_t( tailonleft, intersections, left, value );
 151         classify_headonright_t( headonright, intersections, right, value );
 152         classify_tailonright_t( tailonright, intersections, right, value );
 153     }
 154 }
 155
 156 inline static void
 157 vert_interp( TrimVertex *n, TrimVertex *l, TrimVertex *r, int p, REAL val )
 158 {
 159     assert( val > l->param[p]);
 160     assert( val < r->param[p]);
 161
 162     n->nuid = l->nuid;
 163
 164     n->param[p] = val;
 165     if( l->param[1-p] != r->param[1-p]  ) {
 166         REAL ratio = (val - l->param[p]) / (r->param[p] - l->param[p]);
 167         n->param[1-p] = l->param[1-p] +
 168                         ratio * (r->param[1-p] - l->param[1-p]);
 169     } else {
 170         n->param[1-p] = l->param[1-p];
 171     }
 172 }
 173
 174 int
 175 Subdivider::arc_split( Arc_ptr jarc, int param, REAL value, int dir )
 176 {
 177     int         maxvertex = jarc->pwlArc->npts;
 178     Arc_ptr     jarc1;
 179     TrimVertex* v = jarc->pwlArc->pts;
 180
 181     int         loc[3] = {0,0,0};
 182
 183     switch( pwlarc_intersect( jarc->pwlArc, param, value, dir, loc ) ) {
 184
 185                 // When the parameter value lands on a vertex, life is sweet
 186     case INTERSECT_VERTEX: {
 187             jarc1 = new(arcpool) Arc( jarc, new( pwlarcpool) PwlArc( maxvertex-loc[1], &v[loc[1]] ) );
 188             jarc->pwlArc->npts = loc[1] + 1;
 189             jarc1->next = jarc->next;
 190             jarc1->next->prev = jarc1;
 191             jarc->next = jarc1;
 192             jarc1->prev = jarc;
 193             assert(jarc->check() != 0);
 194             return 2;
 195         }
 196
 197                 // When the parameter value intersects an edge, we have to
 198                 // interpolate a new vertex.  There are special cases
 199                 // if the new vertex is adjacent to one or both of the
 200                 // endpoints of the arc.
 201     case INTERSECT_EDGE: {
 202             int i, j;
 203             if( dir == 0 ) {
 204                 i = loc[0];
 205                 j = loc[2];
 206             } else {
 207                 i = loc[2];
 208                 j = loc[0];
 209             }
 210
 211 #ifndef NOTDEF
 212             // The split is between vertices at index j and i, in that
 213             // order (j < i)
 214
 215             // JEB:  This code is my idea of how to do the split without
 216             // increasing the number of links.  I'm doing this so that
 217             // the is_rect routine can recognize rectangles created by
 218             // subdivision.  In exchange for simplifying the curve list,
 219             // however, it costs in allocated space and vertex copies.
 220
 221             TrimVertex *newjunk = trimvertexpool.get(maxvertex -i+1 /*-j*/);
 222             int k;
 223             for(k=0; k<maxvertex-i; k++)
 224               {
 225                 newjunk[k+1] = v[i+k];
 226                 newjunk[k+1].nuid = jarc->nuid;
 227               }
 228
 229             TrimVertex *vcopy = trimvertexpool.get(maxvertex);
 230             for(k=0; k<maxvertex; k++)
 231               {
 232                 vcopy[k].param[0] = v[k].param[0];
 233                 vcopy[k].param[1] = v[k].param[1];
 234               }
 235             jarc->pwlArc->pts=vcopy;
 236
 237             v[i].nuid = jarc->nuid;
 238             v[j].nuid = jarc->nuid;
 239             vert_interp( &newjunk[0], &v[loc[0]], &v[loc[2]], param, value );
 240
 241             if( showingDegenerate() )
 242                 backend.triangle( &v[i], &newjunk[0], &v[j] );
 243
 244             vcopy[j+1].param[0]=newjunk[0].param[0];
 245             vcopy[j+1].param[1]=newjunk[0].param[1];
 246
 247
 248             jarc1 = new(arcpool) Arc( jarc,
 249                         new(pwlarcpool) PwlArc(maxvertex-i+1 , newjunk ) );
 250
 251             jarc->pwlArc->npts = j+2;
 252             jarc1->next = jarc->next;
 253             jarc1->next->prev = jarc1;
 254             jarc->next = jarc1;
 255             jarc1->prev = jarc;
 256             assert(jarc->check() != 0);
 257
 258             return 2;
 259 #endif //not NOTDEF
 260                 // JEB: This is the original version:
 261 #ifdef NOTDEF
 262
 263             TrimVertex *newjunk = trimvertexpool.get(3);
 264             v[i].nuid = jarc->nuid;
 265             v[j].nuid = jarc->nuid;
 266             newjunk[0] = v[j];
 267             newjunk[2] = v[i];
 268             vert_interp( &newjunk[1], &v[loc[0]], &v[loc[2]], param, value );
 269
 270             if( showingDegenerate() )
 271                 backend.triangle( &newjunk[2], &newjunk[1], &newjunk[0] );
 272
 273                 // New vertex adjacent to both endpoints
 274             if (maxvertex == 2) {
 275                 jarc1 = new(arcpool) Arc( jarc, new(pwlarcpool) PwlArc( 2, newjunk+1 ) );
 276                 jarc->pwlArc->npts = 2;
 277                 jarc->pwlArc->pts = newjunk;
 278                 jarc1->next = jarc->next;
 279                 jarc1->next->prev = jarc1;
 280                 jarc->next = jarc1;
 281                 jarc1->prev = jarc;
 282                 assert(jarc->check() != 0);
 283
 284                 return 2;
 285
 286                 // New vertex adjacent to ending point of arc
 287             } else if (maxvertex - j == 2) {
 288                 jarc1 = new(arcpool) Arc( jarc, new(pwlarcpool) PwlArc( 2, newjunk ) );
 289                 jarc2 = new(arcpool) Arc( jarc, new(pwlarcpool) PwlArc( 2, newjunk+1 ) );
 290                 jarc->pwlArc->npts = maxvertex-1;
 291                 jarc2->next = jarc->next;
 292                 jarc2->next->prev = jarc2;
 293                 jarc->next = jarc1;
 294                 jarc1->prev = jarc;
 295                 jarc1->next = jarc2;
 296                 jarc2->prev = jarc1;
 297                 assert(jarc->check() != 0);
 298                 return 31;
 299
 300                 // New vertex adjacent to starting point of arc
 301             } else if (i == 1) {
 302                 jarc1 = new(arcpool) Arc( jarc, new(pwlarcpool) PwlArc( 2, newjunk+1 ) );
 303                 jarc2 = new(arcpool) Arc( jarc,
 304                         new(pwlarcpool) PwlArc( maxvertex-1, &jarc->pwlArc->pts[1] ) );
 305                 jarc->pwlArc->npts = 2;
 306                 jarc->pwlArc->pts = newjunk;
 307                 jarc2->next = jarc->next;
 308                 jarc2->next->prev = jarc2;
 309                 jarc->next = jarc1;
 310                 jarc1->prev = jarc;
 311                 jarc1->next = jarc2;
 312                 jarc2->prev = jarc1;
 313                 assert(jarc->check() != 0);
 314                 return 32;
 315
 316                 // It's somewhere in the middle
 317             } else {
 318                 jarc1 = new(arcpool) Arc( jarc, new(pwlarcpool) PwlArc( 2, newjunk ) );
 319                 jarc2 = new(arcpool) Arc( jarc, new(pwlarcpool) PwlArc( 2, newjunk+1 ) );
 320                 jarc3 = new(arcpool) Arc( jarc, new(pwlarcpool) PwlArc( maxvertex-i, v+i ) );
 321                 jarc->pwlArc->npts = j + 1;
 322                 jarc3->next = jarc->next;
 323                 jarc3->next->prev = jarc3;
 324                 jarc->next = jarc1;
 325                 jarc1->prev = jarc;
 326                 jarc1->next = jarc2;
 327                 jarc2->prev = jarc1;
 328                 jarc2->next = jarc3;
 329                 jarc3->prev = jarc2;
 330                 assert(jarc->check() != 0);
 331                 return 4;
 332             }
 333 #endif // NOTDEF
 334         }
 335         default:
 336         return -1; //picked -1 since it's not used
 337     }
 338 }
 339
 340 /*----------------------------------------------------------------------------
 341  * pwlarc_intersect -  find intersection of pwlArc and isoparametric line
 342  *----------------------------------------------------------------------------
 343  */
 344
 345 static enum i_result
 346 pwlarc_intersect(
 347     PwlArc *pwlArc,
 348     int param,
 349     REAL value,
 350     int dir,
 351     int loc[3] )
 352 {
 353     assert( pwlArc->npts > 0 );
 354
 355     if( dir ) {
 356         TrimVertex *v = pwlArc->pts;
 357         int imin = 0;
 358         int imax = pwlArc->npts - 1;
 359         assert( value > v[imin].param[param] );
 360         assert( value < v[imax].param[param] );
 361         while( (imax - imin) > 1 ) {
 362             int imid = (imax + imin)/2;
 363             if( v[imid].param[param] > value )
 364                 imax = imid;
 365             else if( v[imid].param[param] < value )
 366                 imin = imid;
 367             else {
 368                 loc[1] = imid;
 369                 return INTERSECT_VERTEX;
 370             }
 371         }
 372         loc[0] = imin;
 373         loc[2] = imax;
 374         return INTERSECT_EDGE;
 375     } else {
 376         TrimVertex *v = pwlArc->pts;
 377         int imax = 0;
 378         int imin = pwlArc->npts - 1;
 379         assert( value > v[imin].param[param] );
 380         assert( value < v[imax].param[param] );
 381         while( (imin - imax) > 1 ) {
 382             int imid = (imax + imin)/2;
 383             if( v[imid].param[param] > value )
 384                 imax = imid;
 385             else if( v[imid].param[param] < value )
 386                 imin = imid;
 387             else {
 388                 loc[1] = imid;
 389                 return INTERSECT_VERTEX;
 390             }
 391         }
 392         loc[0] = imin;
 393         loc[2] = imax;
 394         return INTERSECT_EDGE;
 395     }
 396 }
 397
 398 /*----------------------------------------------------------------------------
 399  * arc_classify - determine which side of a line a jarc lies
 400  *----------------------------------------------------------------------------
 401  */
 402
 403 static int
 404 arc_classify( Arc_ptr jarc, int param, REAL value )
 405 {
 406     REAL tdiff, hdiff;
 407     if( param == 0 ) {
 408         tdiff = jarc->tail()[0] - value;
 409         hdiff = jarc->head()[0] - value;
 410     } else {
 411         tdiff = jarc->tail()[1] - value;
 412         hdiff = jarc->head()[1] - value;
 413     }
 414
 415     if( tdiff > 0.0 ) {
 416         if( hdiff > 0.0 ) {
 417             return 0x11;
 418         } else if( hdiff == 0.0 ) {
 419             return 0x12;
 420         } else {
 421             return 0x10;
 422         }
 423     } else if( tdiff == 0.0 ) {
 424         if( hdiff > 0.0 ) {
 425             return 0x21;
 426         } else if( hdiff == 0.0 ) {
 427             return 0x22;
 428         } else {
 429             return 0x20;
 430         }
 431     } else {
 432         if( hdiff > 0.0 ) {
 433             return 0x01;
 434         } else if( hdiff == 0.0 ) {
 435             return 0x02;
 436         } else {
 437             return 0;
 438         }
 439     }
 440 }
 441
 442 void
 443 Subdivider::classify_tailonleft_s( Bin& bin, Bin& in, Bin& out, REAL val )
 444 {
 445     /* tail at left, head on line */
 446     Arc_ptr j;
 447
 448     while( (j = bin.removearc()) != NULL ) {
 449         assert( arc_classify( j, 0, val ) == 0x02 );
 450         j->clearitail();
 451
 452         REAL diff = j->next->head()[0] - val;
 453         if( diff > 0.0 ) {
 454             in.addarc( j );
 455         } else if( diff < 0.0 ) {
 456             if( ccwTurn_sl( j, j->next ) )
 457                 out.addarc( j );
 458             else
 459                 in.addarc( j );
 460         } else {
 461             if( j->next->tail()[1] > j->next->head()[1] )
 462                 in.addarc(j);
 463             else
 464                 out.addarc(j);
 465         }
 466     }
 467 }
 468
 469 void
 470 Subdivider::classify_tailonleft_t( Bin& bin, Bin& in, Bin& out, REAL val )
 471 {
 472     /* tail at left, head on line */
 473     Arc_ptr j;
 474
 475     while( (j = bin.removearc()) != NULL ) {
 476         assert( arc_classify( j, 1, val ) == 0x02 );
 477         j->clearitail();
 478
 479         REAL diff = j->next->head()[1] - val;
 480         if( diff > 0.0 ) {
 481             in.addarc( j );
 482         } else if( diff < 0.0 ) {
 483             if( ccwTurn_tl( j, j->next ) )
 484                 out.addarc( j );
 485             else
 486                 in.addarc( j );
 487         } else {
 488             if (j->next->tail()[0] > j->next->head()[0] )
 489                 out.addarc( j );
 490             else
 491                 in.addarc( j );
 492         }
 493     }
 494 }
 495
 496 void
 497 Subdivider::classify_headonleft_s( Bin& bin, Bin& in, Bin& out, REAL val )
 498 {
 499     /* tail on line, head at left */
 500     Arc_ptr j;
 501
 502     while( (j = bin.removearc()) != NULL ) {
 503         assert( arc_classify( j, 0, val ) == 0x20 );
 504
 505         j->setitail();
 506
 507         REAL diff = j->prev->tail()[0] - val;
 508         if( diff > 0.0 ) {
 509             out.addarc( j );
 510         } else if( diff < 0.0 ) {
 511             if( ccwTurn_sl( j->prev, j ) )
 512                 out.addarc( j );
 513             else
 514                 in.addarc( j );
 515         } else {
 516             if( j->prev->tail()[1] > j->prev->head()[1] )
 517                 in.addarc( j );
 518             else
 519                 out.addarc( j );
 520         }
 521     }
 522 }
 523
 524 void
 525 Subdivider::classify_headonleft_t( Bin& bin, Bin& in, Bin& out, REAL val )
 526 {
 527     /* tail on line, head at left */
 528     Arc_ptr j;
 529
 530     while( (j = bin.removearc()) != NULL ) {
 531         assert( arc_classify( j, 1, val ) == 0x20 );
 532         j->setitail();
 533
 534         REAL diff = j->prev->tail()[1] - val;
 535         if( diff > 0.0 ) {
 536             out.addarc( j );
 537         } else if( diff < 0.0 ) {
 538             if( ccwTurn_tl( j->prev, j ) )
 539                 out.addarc( j );
 540             else
 541                 in.addarc( j );
 542         } else {
 543             if( j->prev->tail()[0] > j->prev->head()[0] )
 544                 out.addarc( j );
 545             else
 546                 in.addarc( j );
 547         }
 548     }
 549 }
 550
 551
 552 void
 553 Subdivider::classify_tailonright_s( Bin& bin, Bin& in, Bin& out, REAL val )
 554 {
 555     /* tail at right, head on line */
 556     Arc_ptr j;
 557
 558     while( (j = bin.removearc()) != NULL ) {
 559         assert( arc_classify( j, 0, val ) == 0x12);
 560
 561         j->clearitail();
 562
 563         REAL diff = j->next->head()[0] - val;
 564         if( diff > 0.0 ) {
 565             if( ccwTurn_sr( j, j->next ) )
 566                 out.addarc( j );
 567             else
 568                 in.addarc( j );
 569         } else if( diff < 0.0 ) {
 570             in.addarc( j );
 571         } else {
 572             if( j->next->tail()[1] > j->next->head()[1] )
 573                 out.addarc( j );
 574             else
 575                 in.addarc( j );
 576         }
 577     }
 578 }
 579
 580 void
 581 Subdivider::classify_tailonright_t( Bin& bin, Bin& in, Bin& out, REAL val )
 582 {
 583     /* tail at right, head on line */
 584     Arc_ptr j;
 585
 586     while( (j = bin.removearc()) != NULL ) {
 587         assert( arc_classify( j, 1, val ) == 0x12);
 588
 589         j->clearitail();
 590
 591         REAL diff =  j->next->head()[1] - val;
 592         if( diff > 0.0 ) {
 593             if( ccwTurn_tr( j, j->next ) )
 594                 out.addarc( j );
 595             else
 596                 in.addarc( j );
 597         } else if( diff < 0.0 ) {
 598             in.addarc( j );
 599         } else {
 600             if( j->next->tail()[0] > j->next->head()[0] )
 601                 in.addarc( j );
 602             else
 603                 out.addarc( j );
 604         }
 605     }
 606 }
 607
 608 void
 609 Subdivider::classify_headonright_s( Bin& bin, Bin& in, Bin& out, REAL val )
 610 {
 611     /* tail on line, head at right */
 612     Arc_ptr j;
 613
 614     while( (j = bin.removearc()) != NULL ) {
 615         assert( arc_classify( j, 0, val ) == 0x21 );
 616
 617         j->setitail();
 618
 619         REAL diff = j->prev->tail()[0] - val;
 620         if( diff > 0.0 ) {
 621             if( ccwTurn_sr( j->prev, j ) )
 622                 out.addarc( j );
 623             else
 624                 in.addarc( j );
 625         } else if( diff < 0.0 ) {
 626             out.addarc( j );
 627         } else {
 628             if( j->prev->tail()[1] > j->prev->head()[1] )
 629                 out.addarc( j );
 630             else
 631                 in.addarc( j );
 632         }
 633     }
 634 }
 635
 636 void
 637 Subdivider::classify_headonright_t( Bin& bin, Bin& in, Bin& out, REAL val )
 638 {
 639     /* tail on line, head at right */
 640     Arc_ptr j;
 641
 642     while( (j = bin.removearc()) != NULL ) {
 643         assert( arc_classify( j, 1, val ) == 0x21 );
 644
 645         j->setitail();
 646
 647         REAL diff = j->prev->tail()[1] - val;
 648         if( diff > 0.0 ) {
 649             if( ccwTurn_tr( j->prev, j ) )
 650                 out.addarc( j );
 651             else
 652                 in.addarc( j );
 653         } else if( diff < 0.0 ) {
 654             out.addarc( j );
 655         } else {
 656             if( j->prev->tail()[0] > j->prev->head()[0] )
 657                 in.addarc( j );
 658             else
 659                 out.addarc( j );
 660         }
 661     }
 662 }
 663