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