source/zbufinternals.cc

   1 /* This file is part of Shapes.
   2  *
   3  * Shapes is free software: you can redistribute it and/or modify
   4  * it under the terms of the GNU General Public License as published by
   5  * the Free Software Foundation, either version 3 of the License, or
   6  * any later version.
   7  *
   8  * Shapes is distributed in the hope that it will be useful,
   9  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  10  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  11  * GNU General Public License for more details.
  12  *
  13  * You should have received a copy of the GNU General Public License
  14  * along with Shapes.  If not, see <http://www.gnu.org/licenses/>.
  15  *
  16  * Copyright 2008 Henrik Tidefelt
  17  */
  18
  19 #include <cmath>
  20
  21 #include "shapestypes.h"
  22 #include "shapesexceptions.h"
  23 #include "astexpr.h"
  24 #include "consts.h"
  25 #include "angleselect.h"
  26 #include "astvar.h"
  27 #include "astclass.h"
  28 #include "statetypes.h"
  29 #include "lighttypes.h"
  30 #include "shadingtypes.h"
  31 #include "globals.h"
  32 #include "trianglefunctions.h"
  33 #include "constructorrepresentation.h"
  34
  35 #include <ctype.h>
  36 #include <list>
  37 #include <algorithm>
  38
  39 using namespace Shapes;
  40
  41 namespace Shapes
  42 {
  43         namespace Lang
  44         {
  45
  46         class ZBuf::IndexRepresentation
  47         {
  48                 char end_;
  49                 bool pushed_;
  50         public:
  51                 size_t i0;
  52                 size_t i1;
  53                 size_t i2;
  54
  55                 IndexRepresentation( );
  56                 void push_back( size_t i );
  57                 Lang::ZBuf::PolyIndices toPolyIndices( const std::vector< Concrete::Coords2D > & cornerMem ) const;
  58                 size_t getDifferent( size_t notThis1, size_t notThis2 ) const;
  59                 void addEdges( Computation::UndirectedEdgeMatrix< std::pair< bool, std::list< Lang::ZBuf::IndexRepresentation * > > > * edgeMatrix, std::list< Computation::UndirectedEdge > * jobQueue, const std::vector< Concrete::Coords2D > & cornerMem );
  60                 void removeEdges( Computation::UndirectedEdgeMatrix< std::pair< bool, std::list< Lang::ZBuf::IndexRepresentation * > > > * edgeMatrix );
  61                 void pushTriangle( std::list< Computation::ZBufTriangle > * triangles, const std::vector< Concrete::Coords2D > & cornerMem, const Computation::PaintedPolygon3D * painter, const RefCountPtr< const Computation::ZBufTriangle::ZMap > zMap );
  62                 void show( std::ostream & os ) const;
  63                 void show( std::ostream & os, const std::vector< Concrete::Coords2D > & cornerMem ) const;
  64                 Computation::UndirectedEdge getEdge( char pos ) const;
  65                 size_t getDifferent( char pos ) const;
  66                 bool isMember( const size_t i ) const;
  67         };
  68         }
  69 }
  70
  71 Lang::ZBuf::IndexRepresentation::IndexRepresentation( )
  72         : end_( 0 ), pushed_( false )
  73 { }
  74
  75 void
  76 Lang::ZBuf::IndexRepresentation::push_back( size_t i )
  77 {
  78         switch( end_ )
  79                 {
  80                 case 0:
  81                         i0 = i;
  82                         ++end_;
  83                         break;
  84                 case 1:
  85                         i1 = i;
  86                         ++end_;
  87                         break;
  88                 case 2:
  89                         i2 = i;
  90                         ++end_;
  91                         break;
  92                 default:
  93                         throw Exceptions::InternalError( "IndexRepresentation overflow" );
  94                 }
  95 }
  96
  97 Lang::ZBuf::PolyIndices
  98 Lang::ZBuf::IndexRepresentation::toPolyIndices( const std::vector< Concrete::Coords2D > & cornerMem ) const
  99 {
 100         PolyIndices res;
 101         const Concrete::Coords2D d1 = cornerMem[ i1 ] - cornerMem[ i0 ];
 102         const Concrete::Coords2D d2 = cornerMem[ i2 ] - cornerMem[ i0 ];
 103         if( d1.x_ * d2.y_ - d1.y_ * d2.x_ > 0 )
 104                 {
 105                         // we have counter-clockwise ordering
 106                         res.push_back( i0 );
 107                         res.push_back( i1 );
 108                         res.push_back( i2 );
 109                 }
 110         else
 111                 {
 112                         res.push_back( i0 );
 113                         res.push_back( i2 );
 114                         res.push_back( i1 );
 115                 }
 116         return res;
 117 }
 118
 119 size_t
 120 Lang::ZBuf::IndexRepresentation::getDifferent( size_t notThis1, size_t notThis2 ) const
 121 {
 122         if( i0 != notThis1 && i0 != notThis2 )
 123                 {
 124                         return i0;
 125                 }
 126         if( i1 != notThis1 && i1 != notThis2 )
 127                 {
 128                         return i1;
 129                 }
 130         return i2;
 131 }
 132
 133 bool
 134 Lang::ZBuf::IndexRepresentation::isMember( const size_t i ) const
 135 {
 136         switch( end_ )
 137                 {
 138                 case 0:
 139                         return false;
 140                 case 1:
 141                         return i == i0;
 142                 case 2:
 143                         return i == i0 || i == i1;
 144                 case 3:
 145                         return i == i0 || i == i1 || i == i2;
 146                 default:
 147                         throw Exceptions::InternalError( "IndexRepresentation end out of range." );
 148                 }
 149 }
 150
 151 void
 152 Lang::ZBuf::IndexRepresentation::addEdges( Computation::UndirectedEdgeMatrix< std::pair< bool, std::list< Lang::ZBuf::IndexRepresentation * > > > * edgeMatrix, std::list< Computation::UndirectedEdge > * jobQueue, const std::vector< Concrete::Coords2D > & cornerMem )
 153 {
 154         // When adding edges, we must be prepared to accept that the triangle may be rejected at any edge because it
 155         // competes with a bigger triangle on the same side of the edge.
 156
 157         // We should remember that the triangle may first knock out other triangles when being added to the first edges,
 158         // and then be knocked out itself at a later edge.      Then it seems that the triangles it knocked out first shouldn't
 159         // have been knocked out.       However, they wouldn't have been knocked out if they were not very small, so it's
 160         // probably not a big loss anyway.
 161
 162         if( ! addEdge( edgeMatrix, Computation::UndirectedEdge( i0, i1 ), this, jobQueue, cornerMem ) )
 163                 {
 164                         return;
 165                 }
 166         if( ! addEdge( edgeMatrix, Computation::UndirectedEdge( i1, i2 ), this, jobQueue, cornerMem ) )
 167                 {
 168                         // We know that this was push_back'ed, so we don't have to search using removeEdge.
 169                         {
 170                                 std::pair< bool, std::list< Lang::ZBuf::IndexRepresentation * > > & p = (*edgeMatrix)[ Computation::UndirectedEdge( i0, i1 ) ];
 171                                 p.first = false;
 172                                 p.second.pop_back( );
 173                         }
 174                         return;
 175                 }
 176         if( ! addEdge( edgeMatrix, Computation::UndirectedEdge( i2, i0 ), this, jobQueue, cornerMem ) )
 177                 {
 178                         // We know that this was push_back'ed, so we don't have to search using removeEdge.
 179                         {
 180                                 std::pair< bool, std::list< Lang::ZBuf::IndexRepresentation * > > & p = (*edgeMatrix)[ Computation::UndirectedEdge( i0, i1 ) ];
 181                                 p.first = false;
 182                                 p.second.pop_back( );
 183                         }
 184                         {
 185                                 std::pair< bool, std::list< Lang::ZBuf::IndexRepresentation * > > & p = (*edgeMatrix)[ Computation::UndirectedEdge( i1, i2 ) ];
 186                                 p.first = false;
 187                                 p.second.pop_back( );
 188                         }
 189                         return;
 190                 }
 191 }
 192
 193 void
 194 Lang::ZBuf::IndexRepresentation::removeEdges( Computation::UndirectedEdgeMatrix< std::pair< bool, std::list< Lang::ZBuf::IndexRepresentation * > > > * edgeMatrix )
 195 {
 196         removeEdge( edgeMatrix, Computation::UndirectedEdge( i0, i1 ), this );
 197         removeEdge( edgeMatrix, Computation::UndirectedEdge( i1, i2 ), this );
 198         removeEdge( edgeMatrix, Computation::UndirectedEdge( i2, i0 ), this );
 199 }
 200
 201 void
 202 Lang::ZBuf::IndexRepresentation::pushTriangle( std::list< Computation::ZBufTriangle > * triangles, const std::vector< Concrete::Coords2D > & cornerMem, const Computation::PaintedPolygon3D * painter, const RefCountPtr< const Computation::ZBufTriangle::ZMap > zMap )
 203 {
 204         if( pushed_ )
 205                 {
 206                         return;
 207                 }
 208         pushed_ = true;
 209
 210         triangles->push_back( Computation::ZBufTriangle( painter,
 211                                                                                                                                                                                          zMap,
 212                                                                                                                                                                                          cornerMem[ i0 ], cornerMem[ i1 ], cornerMem[ i2 ] ) );
 213 }
 214
 215 void
 216 Lang::ZBuf::IndexRepresentation::show( std::ostream & os ) const
 217 {
 218         switch( end_ )
 219                 {
 220                 case 0:
 221                         os << "No points" ;
 222                         break;
 223                 case 1:
 224                         os << "1 point: " << i0 ;
 225                         break;
 226                 case 2:
 227                         os << "2 points: " << i0 << " -- " << i1 ;
 228                         break;
 229                 case 3:
 230                         os << "3 points: " << i0 << " -- " << i1 << " -- " << i2 ;
 231                         break;
 232                 default:
 233                         os << "end_ out of range!" ;
 234                 }
 235 }
 236
 237 void
 238 Lang::ZBuf::IndexRepresentation::show( std::ostream & os, const std::vector< Concrete::Coords2D > & cornerMem ) const
 239 {
 240         switch( end_ )
 241                 {
 242                 case 0:
 243                         os << "|** No points" ;
 244                         break;
 245                 case 1:
 246                         os << "|** 1 point:     " << cornerMem[ i0 ] ;
 247                         break;
 248                 case 2:
 249                         os << "|** 2 points: " << cornerMem[ i0 ] << " -- " << cornerMem[ i1 ] ;
 250                         break;
 251                 case 3:
 252                         os << "@<< @width:0.3bp | stroke [] " << Helpers::shapesFormat( cornerMem[ i0 ] ) << "--" << Helpers::shapesFormat( cornerMem[ i1 ] ) << "--" << Helpers::shapesFormat( cornerMem[ i2 ] ) << "--cycle" << std::endl ;
 253                         break;
 254                 default:
 255                         os << "end_ out of range!" ;
 256                 }
 257 }
 258
 259 Computation::UndirectedEdge
 260 Lang::ZBuf::IndexRepresentation::getEdge( char pos ) const
 261 {
 262         switch( pos )
 263                 {
 264                 case 0:
 265                         return Computation::UndirectedEdge( i0, i1 );
 266                 case 1:
 267                         return Computation::UndirectedEdge( i1, i2 );
 268                 case 2:
 269                         return Computation::UndirectedEdge( i2, i0 );
 270                 default:
 271                         throw Exceptions::InternalError( "getEdge: pos out of range." );
 272                 }
 273 }
 274
 275 size_t
 276 Lang::ZBuf::IndexRepresentation::getDifferent( char pos ) const
 277 {
 278         switch( pos )
 279                 {
 280                 case 0:
 281                         return i2;
 282                 case 1:
 283                         return i0;
 284                 case 2:
 285                         return i1;
 286                 default:
 287                         throw Exceptions::InternalError( "getEdge: pos out of range." );
 288                 }
 289 }
 290
 291
 292 void
 293 Lang::ZBuf::mergeTriangles( std::list< Computation::ZBufTriangle > * triangles )
 294 {
 295         if( triangles->size( ) <= 1 )
 296                 {
 297                         return;
 298                 }
 299
 300         const Computation::PaintedPolygon3D * painter = triangles->front( ).painter_;
 301         const RefCountPtr< const Computation::ZBufTriangle::ZMap > zMap = triangles->front( ).zMap_;
 302
 303         std::vector< Lang::ZBuf::IndexRepresentation > indexRepMem;
 304         indexRepMem.reserve( 2 * triangles->size( ) ); // This is crucial!
 305
 306         // One needs to be careful here.        The reuse of points based on distances being less than the tolerance
 307         // may cause two points on the same triangle to be approximated by the very _same_ point.       This will
 308         // lead to a degenerate triangle, and must be avoided in order to prevent confusion later on
 309         // in the algorithm.
 310
 311         // Even worse, the triangles being "disjoint" does not mean they are perfectly disjoint (perhaps due to some
 312         // point approximations somewhere), but only overlaping by at most some tolerance value.        Hence, we may
 313         // even find triangles included in other triangles and similar phenomena.
 314         // This would be interesting to design better with more powerful invariants so that one could guarantee
 315         // stronger properties. However, as the current state of affairs when writing this desperate note is
 316         // that we may find ourselves with more than two triangles sharing the same edge!
 317         // The solution to this problem is to discard, on each side of the edge (a zero tolerance test is used to determine the side
 318         // of a triangle) discard all but the triangle with the biggest area.   Luckily, this can be done incrementaly by
 319         // always keeping track of the biggest triangle on each side.
 320         // Need I say it feels like this algorithm is running out of control?
 321
 322         std::vector< Concrete::Coords2D > cornerMem;
 323         while( ! triangles->empty( ) )
 324                 {
 325                         const std::vector< Concrete::Coords2D > & points = triangles->front( ).points_;
 326                         indexRepMem.push_back( IndexRepresentation( ) );
 327                         Lang::ZBuf::IndexRepresentation & current = indexRepMem.back( );
 328                         for( std::vector< Concrete::Coords2D >::const_iterator i = points.begin( ); i != points.end( ); ++i )
 329                                 {
 330                                         Concrete::Coords2D tmp = *i;
 331                                         bool reuse = false;
 332                                         size_t memIdx = 0;
 333                                         typedef typeof cornerMem MemType;
 334                                         for( MemType::const_iterator j = cornerMem.begin( ); j != cornerMem.end( ); ++j, ++memIdx )
 335                                                 {
 336                                                         if( ( tmp - *j ).norm( ) < Computation::theTrixelizeSplicingTol )
 337                                                                 {
 338                                                                         if( current.isMember( memIdx ) )
 339                                                                                 {
 340                                                                                         // This is the error situation when we are just about to reuse a point which is already part of this triangle.
 341                                                                                         // The solution is to abort the creation of the current triangle.
 342                                                                                         indexRepMem.pop_back( );
 343                                                                                         goto abortCurrent;
 344                                                                                 }
 345                                                                         reuse = true;
 346                                                                         current.push_back( memIdx );
 347                                                                         break;
 348                                                                 }
 349                                                 }
 350                                         if( ! reuse )
 351                                                 {
 352                                                         current.push_back( cornerMem.size( ) );
 353                                                         cornerMem.push_back( tmp );
 354                                                 }
 355                                 }
 356                 abortCurrent:
 357                         triangles->pop_front( );
 358                 }
 359
 360         const size_t pointCount = cornerMem.size( );
 361
 362         Computation::UndirectedEdgeMatrix< std::pair< bool, std::list< Lang::ZBuf::IndexRepresentation * > > >
 363                 edgeMatrix( pointCount,
 364                                                                 std::pair< bool, std::list< Lang::ZBuf::IndexRepresentation * > >( false,
 365                                                                                                                                                                                                                                                                                                                                          std::list< Lang::ZBuf::IndexRepresentation * >( ) ) );
 366
 367         std::list< Computation::UndirectedEdge > jobQueue;
 368
 369         {
 370                 typedef typeof indexRepMem ListType;
 371                 for( ListType::iterator i = indexRepMem.begin( ); i != indexRepMem.end( ); ++i )
 372                         {
 373                                 i->addEdges( & edgeMatrix, & jobQueue, cornerMem );
 374                         }
 375         }
 376
 377         while( ! jobQueue.empty( ) )
 378                 {
 379                         Computation::UndirectedEdge edge = jobQueue.front( );
 380                         jobQueue.pop_front( );
 381
 382                         // Elements can be removed from the jobQueue by unflagging them; only flagged elements that appear in the jobQueue
 383                         // are really waiting to be processed.
 384                         if( ! edgeMatrix[ edge ].first )
 385                                 {
 386                                         continue;
 387                                 }
 388                         edgeMatrix[ edge ].first = false;       // Flag that this is no longer waiting to be processed.
 389
 390                         typedef std::list< Lang::ZBuf::IndexRepresentation * > ListType;
 391                         ListType & meetingTriangles = edgeMatrix[ edge ].second;
 392
 393                         if( meetingTriangles.size( ) < 2 )
 394                                 {
 395                                         // One situation I can think of when this could happen is when three triangles have been
 396                                         // found to be on the same side of the edge during addEdge.     Then two of them are removed,
 397                                         // but the edge is still in the jobQueue since it was queued when the second triange was added
 398                                         // to the edge.
 399                                         continue;
 400                                 }
 401
 402                         if( meetingTriangles.size( ) > 2 )
 403                                 {
 404                                         std::ostringstream oss;
 405                                         oss << "Only two (not " << meetingTriangles.size( ) << ") edges can share a pair of points!" ;
 406                                         oss << std::endl << "|** Triangles in queue: " ;
 407                                         for( ListType::const_iterator i = meetingTriangles.begin( ); i != meetingTriangles.end( ); ++i )
 408                                                 {
 409                                                         oss << std::endl ;
 410                                                         (*i)->show( oss, cornerMem );
 411                                                 }
 412                                         throw Exceptions::InternalError( oss );
 413                                 }
 414
 415                         Concrete::Coords2D cornerLow = cornerMem[ edge.low( ) ];
 416                         Concrete::Coords2D cornerHigh = cornerMem[ edge.high( ) ];
 417
 418                         ListType::iterator src = meetingTriangles.begin( );
 419                         Lang::ZBuf::IndexRepresentation * triangle1 = *src;
 420                         Concrete::Coords2D differentPoint1 = cornerMem[ triangle1->getDifferent( edge.low( ), edge.high( ) ) ];
 421
 422                         ++src;
 423                         Lang::ZBuf::IndexRepresentation * triangle2 = *src;
 424                         Concrete::Coords2D differentPoint2 = cornerMem[ triangle2->getDifferent( edge.low( ), edge.high( ) ) ];
 425                         bool aligned = false;
 426                         size_t keepPoint = 0;   /* Initialized with dummy value only because it is hard for the compiler to see that it won't be used uninitialized. */
 427                         if( areAligned( cornerLow, differentPoint1, differentPoint2 ) )
 428                                 {
 429                                         aligned = true;
 430                                         keepPoint = edge.high( );
 431                                 }
 432                         else if( areAligned( cornerHigh, differentPoint1, differentPoint2 ) )
 433                                 {
 434                                         aligned = true;
 435                                         keepPoint = edge.low( );
 436                                 }
 437
 438                         if( aligned )
 439                                 {
 440                                         triangle1->removeEdges( & edgeMatrix );
 441                                         triangle2->removeEdges( & edgeMatrix );
 442                                         indexRepMem.push_back( IndexRepresentation( ) );
 443                                         Lang::ZBuf::IndexRepresentation & current = indexRepMem.back( );
 444                                         current.push_back( triangle1->getDifferent( edge.low( ), edge.high( ) ) );
 445                                         current.push_back( triangle2->getDifferent( edge.low( ), edge.high( ) ) );
 446                                         current.push_back( keepPoint );
 447                                         current.addEdges( & edgeMatrix, & jobQueue, cornerMem );
 448                                 }
 449                 }
 450
 451         for( Computation::UndirectedEdge edge = Computation::UndirectedEdge::begin( ); edge != Computation::UndirectedEdge::end( ); edge.increment( pointCount ) )
 452                 {
 453                         std::list< Lang::ZBuf::IndexRepresentation * > & lst = edgeMatrix[ edge ].second;
 454                         for( std::list< Lang::ZBuf::IndexRepresentation * >::iterator i = lst.begin( ); i != lst.end( ); ++i )
 455                                 {
 456                                         // The IndexRepresentation keeps a state to ensure that is is really pushed only once.
 457                                         (*i)->pushTriangle( triangles, cornerMem, painter, zMap );
 458                                 }
 459                 }
 460 }
 461
 462 // recombineTriangles generalizes mergeTriangles, but is more expensive and should be called less often.
 463 void
 464 Lang::ZBuf::recombineTriangles( std::list< Computation::ZBufTriangle > * mergedTriangles )
 465 {
 466         if( mergedTriangles->size( ) <= 1 )
 467                 {
 468                         return;
 469                 }
 470
 471         const Computation::PaintedPolygon3D * painter = mergedTriangles->front( ).painter_;
 472         const RefCountPtr< const Computation::ZBufTriangle::ZMap > zMap = mergedTriangles->front( ).zMap_;
 473
 474         std::list< Lang::ZBuf::IndexRepresentation > indexRepMem;
 475
 476         std::vector< Concrete::Coords2D > cornerMem;
 477         while( ! mergedTriangles->empty( ) )
 478                 {
 479                         const std::vector< Concrete::Coords2D > & points = mergedTriangles->front( ).points_;
 480                         indexRepMem.push_back( IndexRepresentation( ) );
 481                         Lang::ZBuf::IndexRepresentation & current = indexRepMem.back( );
 482                         for( std::vector< Concrete::Coords2D >::const_iterator i = points.begin( ); i != points.end( ); ++i )
 483                                 {
 484                                         Concrete::Coords2D tmp = *i;
 485                                         bool reuse = false;
 486                                         size_t memIdx = 0;
 487                                         typedef typeof cornerMem MemType;
 488                                         for( MemType::const_iterator j = cornerMem.begin( ); j != cornerMem.end( ); ++j, ++memIdx )
 489                                                 {
 490                                                         if( ( tmp - *j ).norm( ) < Computation::theTrixelizeSplicingTol )
 491                                                                 {
 492                                                                         if( current.isMember( memIdx ) )
 493                                                                                 {
 494                                                                                         // This is the error situation when we are just about to reuse a point which is already part of this triangle.
 495                                                                                         // The solution is to abort the creation of the current triangle.
 496                                                                                         indexRepMem.pop_back( );
 497                                                                                         goto abortCurrent;
 498                                                                                 }
 499                                                                         reuse = true;
 500                                                                         current.push_back( memIdx );
 501                                                                         break;
 502                                                                 }
 503                                                 }
 504                                         if( ! reuse )
 505                                                 {
 506                                                         current.push_back( cornerMem.size( ) );
 507                                                         cornerMem.push_back( tmp );
 508                                                 }
 509                                 }
 510                 abortCurrent:
 511                         mergedTriangles->pop_front( );
 512                 }
 513
 514         const size_t pointCount = cornerMem.size( );
 515
 516         std::vector< PolyIndices > polyMem;
 517         polyMem.reserve( indexRepMem.size( ) ); // this is crucial to avoid reallocation!
 518
 519         Computation::UndirectedEdgeMatrix< std::list< PolyIndices * > > edgeMatrix( pointCount );
 520
 521         while( ! indexRepMem.empty( ) )
 522         {
 523                 bool foundExtension = false;
 524                 typedef typeof indexRepMem TrianglesType;
 525                 for( TrianglesType::iterator i = indexRepMem.begin( ); i != indexRepMem.end( ) && ! foundExtension; )
 526                         {
 527                                 bool erase = false;
 528                                 for( char e = 0; e < 3; ++e )
 529                                         {
 530                                                 Computation::UndirectedEdge edge = i->getEdge( e );
 531                                                 std::list< PolyIndices * > & polys = edgeMatrix[ edge ];
 532                                                 if( polys.size( ) == 1 )        // By the way, it is an error if the size is greater than 1
 533                                                         {
 534                                                                 if( ! extendPoly( polys.front( ), edge, i->getDifferent( e ), & edgeMatrix, cornerMem, true ) ) // true for convex
 535                                                                         {
 536                                                                                 continue;
 537                                                                         }
 538                                                                 foundExtension = true;
 539                                                                 erase = true;
 540                                                                 break;
 541                                                         }
 542                                         }
 543                                 if( erase ){
 544                                         TrianglesType::iterator iErase = i;
 545                                         ++i;
 546                                         indexRepMem.erase( iErase );
 547                                 }else{
 548                                         ++i;
 549                                 }
 550                         }
 551                 if( ! foundExtension )
 552                         {
 553                                 polyMem.push_back( indexRepMem.front( ).toPolyIndices( cornerMem ) );
 554                                 indexRepMem.pop_front( );
 555                                 PolyIndices * newPoly = & polyMem.back( );
 556                                 addEdges( newPoly, & edgeMatrix );
 557                         }
 558         }
 559
 560         while( ! polyMem.empty( ) )
 561                 {
 562                         PolyIndices & poly = polyMem.back( );
 563                         // First we remove unnecessary points on straight segments
 564                         if( poly.size( ) > 3 )
 565                                 {
 566                                         for( PolyIndices::iterator i0 = poly.begin( ); i0 != poly.end( ); ++i0 )
 567                                                 {
 568                                                         Concrete::Coords2D p0 = cornerMem[ *i0 ];
 569                                                         PolyIndices::iterator i1 = i0;
 570                                                         ++i1;
 571                                                         if( i1 == poly.end( ) )
 572                                                                 {
 573                                                                         i1 = poly.begin( );
 574                                                                 }
 575                                                         Concrete::Coords2D p1 = cornerMem[ *i1 ];
 576                                                         PolyIndices::iterator i2 = i1;
 577                                                         ++i2;
 578                                                         if( i2 == poly.end( ) )
 579                                                                 {
 580                                                                         i2 = poly.begin( );
 581                                                                 }
 582                                                         Concrete::Coords2D p2 = cornerMem[ *i2 ];
 583
 584                                                         // We do the usual inscribed circle test to test alignment
 585                                                         while( Computation::triangleArea( p0, p1, p2 ) < Computation::theTrixelizeOverlapTol * Computation::triangleSemiPerimeter( p0, p1, p2 ) )
 586                                                                 {
 587                                                                         PolyIndices::iterator iRemove = i1;
 588                                                                         i1 = i2;
 589                                                                         p1 = p2;
 590                                                                         poly.erase( iRemove );
 591                                                                         if( poly.size( ) < 3 )
 592                                                                                 {
 593                                                                                         goto done;
 594                                                                                 }
 595                                                                         ++i2;
 596                                                                         if( i2 == poly.end( ) )
 597                                                                                 {
 598                                                                                         i2 = poly.begin( );
 599                                                                                 }
 600                                                                         p2 = cornerMem[ *i2 ];
 601                                                                 }
 602                                                 }
 603                                 done:
 604                                         if( poly.size( ) < 3 )
 605                                                 {
 606                                                         polyMem.pop_back( );
 607                                                         continue;
 608                                                 }
 609                                 }
 610
 611
 612                         PolyIndices::const_iterator src = poly.begin( );
 613                         size_t i0 = *src;
 614                         ++src;
 615                         size_t i1 = *src;
 616                         ++src;
 617                         for( ; src != poly.end( ); ++src )
 618                                 {
 619                                         size_t i2 = *src;
 620                                         mergedTriangles->push_back( Computation::ZBufTriangle( painter,
 621                                                                                                                                                                                                                                                  zMap,
 622                                                                                                                                                                                                                                                  cornerMem[ i0 ], cornerMem[ i1 ], cornerMem[ i2 ] ) );
 623                                         i1 = i2;
 624                                 }
 625                         polyMem.pop_back( );
 626                 }
 627
 628 }
 629
 630 void
 631 Lang::ZBuf::trianglesToPolys( std::list< Computation::ZBufTriangle > * triangles, Lang::Clipped2D * dst )
 632 {
 633         if( triangles->empty( ) )
 634                 {
 635                         return;
 636                 }
 637
 638         const bool DEBUGME = false;
 639         if( DEBUGME )
 640         {
 641                 while( ! triangles->empty( ) )
 642                         {
 643                                 const std::vector< Concrete::Coords2D > & points = triangles->front( ).points_;
 644
 645                                 Lang::ElementaryPath2D * newPath = new Lang::ElementaryPath2D;
 646
 647                                 for( std::vector< Concrete::Coords2D >::const_iterator src = points.begin( ); src != points.end( ); ++src )
 648                                 {
 649                                         newPath->push_back( new Concrete::PathPoint2D( new Concrete::Coords2D( *src ) ) );
 650                                 }
 651                                 newPath->close( );
 652                                 dst->addSubPath( RefCountPtr< const Lang::ElementaryPath2D >( newPath ) );
 653
 654                                 triangles->pop_front( );
 655                         }
 656                 return;
 657         }
 658
 659         std::list< Lang::ZBuf::IndexRepresentation > indexRepMem;
 660         std::vector< Concrete::Coords2D > cornerMem;
 661
 662         {
 663                 std::list< Lang::ZBuf::IndexRepresentation > indexRepQueue;
 664                 while( ! triangles->empty( ) )
 665                         {
 666                                 const std::vector< Concrete::Coords2D > & points = triangles->front( ).points_;
 667                                 indexRepQueue.push_back( IndexRepresentation( ) );
 668                                 Lang::ZBuf::IndexRepresentation & current = indexRepQueue.back( );
 669                                 for( std::vector< Concrete::Coords2D >::const_iterator i = points.begin( ); i != points.end( ); ++i )
 670                                         {
 671                                                 Concrete::Coords2D tmp = *i;
 672                                                 bool reuse = false;
 673                                                 size_t memIdx = 0;
 674                                                 typedef typeof cornerMem MemType;
 675                                                 for( MemType::const_iterator j = cornerMem.begin( ); j != cornerMem.end( ); ++j, ++memIdx )
 676                                                         {
 677                                                                 if( ( tmp - *j ).norm( ) < Computation::theTrixelizeSplicingTol )
 678                                                                         {
 679                                                                                 if( current.isMember( memIdx ) )
 680                                                                                         {
 681                                                                                                 // This is the error situation when we are just about to reuse a point which is already part of this triangle.
 682                                                                                                 // The solution is to abort the creation of the current triangle.
 683                                                                                                 indexRepMem.pop_back( );
 684                                                                                                 goto abortCurrent;
 685                                                                                         }
 686                                                                                 reuse = true;
 687                                                                                 current.push_back( memIdx );
 688                                                                                 break;
 689                                                                         }
 690                                                         }
 691                                                 if( ! reuse )
 692                                                         {
 693                                                                 current.push_back( cornerMem.size( ) );
 694                                                                 cornerMem.push_back( tmp );
 695                                                         }
 696                                         }
 697                         abortCurrent:
 698                                 triangles->pop_front( );
 699                         }
 700
 701                 while( ! indexRepQueue.empty( ) )
 702                         {
 703                                 Lang::ZBuf::IndexRepresentation & current = indexRepQueue.front( );
 704
 705                                 bool wasSplit = false;
 706                                 for( char e = 0; e < 3 && ! wasSplit; ++e )
 707                                         {
 708                                                 Computation::UndirectedEdge edge = current.getEdge( e );
 709                                                 size_t iLow = edge.low( );
 710                                                 size_t iHigh = edge.high( );
 711                                                 size_t iOther = current.getDifferent( e );
 712                                                 Concrete::Coords2D pLow = cornerMem[ iLow ];
 713                                                 Concrete::Coords2D pHigh = cornerMem[ iHigh ];
 714                                                 for( size_t i = 0; i < cornerMem.size( ) && ! wasSplit; ++i )
 715                                                         {
 716                                                                 if( i == iLow || i == iHigh || i == iOther )
 717                                                                         {
 718                                                                                 continue;
 719                                                                         }
 720                                                                 if( areAlignedAndOrdered( pLow, cornerMem[ i ], pHigh ) )
 721                                                                         {
 722                                                                                 // To avoid infinite loops, we make sure that the triangles get smaller each time
 723                                                                                 // they are splitted.   Otherwise, four points in a row can cause an infinite loop.
 724
 725                                                                                 Concrete::Coords2D pOther = cornerMem[ iOther ];
 726                                                                                 Concrete::Coords2D pi = cornerMem[ i ];
 727
 728                                                                                 Concrete::Area area0 = Computation::triangleArea( pLow, pHigh, pOther );
 729
 730                                                                                 const double AREA_TOL = 1e-3;
 731
 732                                                                                 {
 733                                                                                         const double relArea = Computation::triangleArea( pOther, pLow, pi ) / area0;
 734                                                                                         if( AREA_TOL < relArea && relArea < 1 - AREA_TOL )
 735                                                                                                 {
 736                                                                                                         indexRepQueue.push_back( IndexRepresentation( ) );
 737                                                                                                         Lang::ZBuf::IndexRepresentation & newTriangle = indexRepQueue.back( );
 738                                                                                                         newTriangle.push_back( iOther );
 739                                                                                                         newTriangle.push_back( iLow );
 740                                                                                                         newTriangle.push_back( i );
 741                                                                                                         wasSplit = true;
 742                                                                                                 }
 743                                                                                 }
 744                                                                                 {
 745                                                                                         const double relArea = Computation::triangleArea( pOther, pi, pHigh ) / area0;
 746                                                                                         if( AREA_TOL < relArea && relArea < 1 - AREA_TOL )
 747                                                                                                 {
 748                                                                                                         indexRepQueue.push_back( IndexRepresentation( ) );
 749                                                                                                         Lang::ZBuf::IndexRepresentation & newTriangle = indexRepQueue.back( );
 750                                                                                                         newTriangle.push_back( iOther );
 751                                                                                                         newTriangle.push_back( i );
 752                                                                                                         newTriangle.push_back( iHigh );
 753                                                                                                         wasSplit = true;
 754                                                                                                 }
 755                                                                                 }
 756                                                                         }
 757                                                         }
 758                                         }
 759
 760                                 if( ! wasSplit )
 761                                         {
 762                                                 indexRepMem.push_back( current );
 763                                         }
 764
 765                                 indexRepQueue.pop_front( );
 766                         }
 767         }
 768
 769         const size_t pointCount = cornerMem.size( );
 770
 771         std::vector< PolyIndices > polyMem;
 772         polyMem.reserve( indexRepMem.size( ) ); // this is crucial to avoid reallocation!
 773
 774         Computation::UndirectedEdgeMatrix< std::list< PolyIndices * > > edgeMatrix( pointCount );
 775
 776         while( ! indexRepMem.empty( ) )
 777         {
 778                 bool foundExtension = false;
 779                 typedef typeof indexRepMem TrianglesType;
 780                 for( TrianglesType::iterator i = indexRepMem.begin( ); i != indexRepMem.end( ) && ! foundExtension; )
 781                         {
 782                                 bool erase = false;
 783                                 for( char e = 0; e < 3; ++e )
 784                                         {
 785                                                 Computation::UndirectedEdge edge = i->getEdge( e );
 786                                                 std::list< PolyIndices * > & polys = edgeMatrix[ edge ];
 787                                                 if( polys.size( ) == 1 )        // By the way, it is an error if the size is greater than 1
 788                                                         {
 789                                                                 if( ! extendPoly( polys.front( ), edge, i->getDifferent( e ), & edgeMatrix, cornerMem, false ) ) // false to not require convex
 790                                                                         {
 791                                                                                 continue;
 792                                                                         }
 793                                                                 foundExtension = true;
 794                                                                 erase = true;
 795                                                                 break;
 796                                                         }
 797                                         }
 798                                 if( erase ){
 799                                         TrianglesType::iterator iErase = i;
 800                                         ++i;
 801                                         indexRepMem.erase( iErase );
 802                                 }else{
 803                                         ++i;
 804                                 }
 805                         }
 806                 if( ! foundExtension )
 807                         {
 808                                 polyMem.push_back( indexRepMem.front( ).toPolyIndices( cornerMem ) );
 809                                 indexRepMem.pop_front( );
 810                                 PolyIndices * newPoly = & polyMem.back( );
 811                                 addEdges( newPoly, & edgeMatrix );
 812                         }
 813         }
 814
 815         while( ! polyMem.empty( ) )
 816                 {
 817                         PolyIndices & poly = polyMem.back( );
 818                         if( poly.size( ) > 3 )
 819                                 {
 820                                         // First we remove insections, that is sequences like 4-5-4 et c.
 821                                         {
 822                                                 PolyIndices::iterator theEnd = poly.begin( );
 823                                                 PolyIndices::iterator i0 = poly.begin( );
 824                                                 while( poly.size( ) >= 5 )      // it would be pathological to have insection in polygons with four corners or less
 825                                                         {
 826                                                                 PolyIndices::iterator i1 = i0;
 827                                                                 ++i1;
 828                                                                 if( i1 == poly.end( ) )
 829                                                                         {
 830                                                                                 i1 = poly.begin( );
 831                                                                         }
 832                                                                 PolyIndices::iterator i2 = i1;
 833                                                                 ++i2;
 834                                                                 if( i2 == poly.end( ) )
 835                                                                         {
 836                                                                                 i2 = poly.begin( );
 837                                                                         }
 838
 839                                                                 if( *i0 == *i2 )
 840                                                                         {
 841                                                                                 // We found an insection
 842                                                                                 poly.erase( i1 );
 843                                                                                 poly.erase( i2 );
 844                                                                                 theEnd = i0;
 845                                                                         }
 846
 847                                                                 ++i0;
 848                                                                 if( i0 == poly.end( ) )
 849                                                                         {
 850                                                                                 i0 = poly.begin( );
 851                                                                         }
 852                                                                 if( i0 == theEnd )
 853                                                                         {
 854                                                                                 break;
 855                                                                         }
 856                                                         }
 857                                         }
 858                                         // Then we remove unnecessary points on straight segments
 859                                         for( PolyIndices::iterator i0 = poly.begin( ); i0 != poly.end( ); ++i0 )
 860                                                 {
 861                                                         Concrete::Coords2D p0 = cornerMem[ *i0 ];
 862                                                         PolyIndices::iterator i1 = i0;
 863                                                         ++i1;
 864                                                         if( i1 == poly.end( ) )
 865                                                                 {
 866                                                                         i1 = poly.begin( );
 867                                                                 }
 868                                                         Concrete::Coords2D p1 = cornerMem[ *i1 ];
 869                                                         PolyIndices::iterator i2 = i1;
 870                                                         ++i2;
 871                                                         if( i2 == poly.end( ) )
 872                                                                 {
 873                                                                         i2 = poly.begin( );
 874                                                                 }
 875                                                         Concrete::Coords2D p2 = cornerMem[ *i2 ];
 876
 877                                                         // We do the usual inscribed circle test to test alignment
 878                                                         while( Computation::triangleArea( p0, p1, p2 ) < Computation::theTrixelizeOverlapTol * Computation::triangleSemiPerimeter( p0, p1, p2 ) )
 879                                                                 {
 880                                                                         PolyIndices::iterator iRemove = i1;
 881                                                                         i1 = i2;
 882                                                                         p1 = p2;
 883                                                                         poly.erase( iRemove );
 884                                                                         if( poly.size( ) < 3 )
 885                                                                                 {
 886                                                                                         goto done;
 887                                                                                 }
 888                                                                         ++i2;
 889                                                                         if( i2 == poly.end( ) )
 890                                                                                 {
 891                                                                                         i2 = poly.begin( );
 892                                                                                 }
 893                                                                         p2 = cornerMem[ *i2 ];
 894                                                                 }
 895                                                 }
 896                                 done:
 897                                         if( poly.size( ) < 3 )
 898                                                 {
 899                                                         polyMem.pop_back( );
 900                                                         continue;
 901                                                 }
 902                                 }
 903
 904
 905                         {
 906                                 Lang::ElementaryPath2D * newPath = new Lang::ElementaryPath2D;
 907
 908                                 for( PolyIndices::const_iterator src = poly.begin( ); src != poly.end( ); ++src )
 909                                 {
 910                                         newPath->push_back( new Concrete::PathPoint2D( new Concrete::Coords2D( cornerMem[ *src ] ) ) );
 911                                 }
 912                                 newPath->close( );
 913                                 dst->addSubPath( RefCountPtr< const Lang::ElementaryPath2D >( newPath ) );
 914                         }
 915                         polyMem.pop_back( );
 916                 }
 917
 918 }
 919
 920 bool
 921 Lang::ZBuf::addEdge( Computation::UndirectedEdgeMatrix< std::pair< bool, std::list< Lang::ZBuf::IndexRepresentation * > > > * edgeMatrix, const Computation::UndirectedEdge & edge, Lang::ZBuf::IndexRepresentation * indRep, std::list< Computation::UndirectedEdge > * jobQueue, const std::vector< Concrete::Coords2D > & cornerMem )
 922 {
 923         typedef typeof (*edgeMatrix)[ edge ].second ListType;
 924         ListType & trianglesAtEdge = (*edgeMatrix)[ edge ].second;
 925         trianglesAtEdge.push_back( indRep );
 926
 927         bool alreadyInQueue = false;
 928
 929         if( trianglesAtEdge.size( ) > 2 )
 930                 {
 931                         // This is a problematic case.  I've written about it elsewhere.        The decision I took was to discard
 932                         // the smaller if the triangles being on the same side of the edge.
 933                         // Since we only do this when there are already three triangles sharing the edge, there's a risk
 934                         // that all three are on the same side. Then we may even drop two of them.
 935
 936                         alreadyInQueue = true;
 937
 938                         Concrete::UnitFloatPair edgeNormalNotNormalized = cornerMem[ edge.low( ) ].unNormalizedOrthogonal( cornerMem[ edge.high( ) ] );
 939                         double m = Concrete::innerScalar( cornerMem[ edge.low( ) ], edgeNormalNotNormalized );
 940
 941                         Lang::ZBuf::IndexRepresentation * posMax = 0;
 942                         Lang::ZBuf::IndexRepresentation * negMax = 0;
 943                         Concrete::Area posArea = -1;
 944                         Concrete::Area negArea = -1;
 945
 946                         for( ListType::iterator i = trianglesAtEdge.begin( ); i != trianglesAtEdge.end( ); ++i )
 947                                 {
 948                                         Concrete::Area tmpArea = Computation::triangleArea( cornerMem[ (*i)->i0 ], cornerMem[ (*i)->i1 ], cornerMem[ (*i)->i2 ] );
 949                                         if( Concrete::innerScalar( cornerMem[ (*i)->getDifferent( edge.low( ), edge.high( ) ) ], edgeNormalNotNormalized ) > m )
 950                                                 {
 951                                                         if( tmpArea > posArea )
 952                                                                 {
 953                                                                         posMax = *i;
 954                                                                         posArea = tmpArea;
 955                                                                 }
 956                                                 }
 957                                         else
 958                                                 {
 959                                                         if( tmpArea > negArea )
 960                                                                 {
 961                                                                         negMax = *i;
 962                                                                         negArea = tmpArea;
 963                                                                 }
 964                                                 }
 965                                 }
 966
 967                         // *** It is required that the new triangle appears at the end of trianglesAtEdge in case we add it. ***
 968
 969                         // We temporarily remove the newly added triangle...
 970                         trianglesAtEdge.pop_back( );
 971
 972                         // Now we remove those of the old triangles that were knocked out.
 973                         // Note that removing triangles from edges clears the job flags at their edges, so we don't have to do that here.
 974                         {
 975                                 // We must be careful since trianglesAtEdge might change when triangles are removed!
 976                                 ListType::iterator i = trianglesAtEdge.begin( );
 977                                 Lang::ZBuf::IndexRepresentation * t1 = *i;
 978                                 ++i;
 979                                 Lang::ZBuf::IndexRepresentation * t2 = *i;
 980
 981                                 if( t1 != posMax && t1 != negMax )
 982                                         {
 983                                                 t1->removeEdges( edgeMatrix );
 984                                         }
 985                                 if( t2 != posMax && t2 != negMax )
 986                                         {
 987                                                 t2->removeEdges( edgeMatrix );
 988                                         }
 989                         }
 990
 991                         // Next, we consider the new triangle.
 992                         if( indRep != posMax && indRep != negMax )
 993                                 {
 994                                         // The new triangle couldn't compete with the old ones, so we shall return false.
 995                                         return false;
 996                                 }
 997
 998                         // ... and now we put it back if it passed the test.
 999                         trianglesAtEdge.push_back( indRep );
1000                 }
1001
1002         if( jobQueue != 0 )
1003                 {
1004                         if( trianglesAtEdge.size( ) == 2 && ! (*edgeMatrix)[ edge ].first )
1005                                 {
1006                                         if( ! alreadyInQueue )
1007                                                 {
1008                                                         jobQueue->push_back( edge );
1009                                                 }
1010                                         (*edgeMatrix)[ edge ].first = true; // Flag that this pair is in the queue.
1011                                 }
1012                 }
1013         return true;
1014 }
1015
1016 void
1017 Lang::ZBuf::removeEdge( Computation::UndirectedEdgeMatrix< std::pair< bool, std::list< Lang::ZBuf::IndexRepresentation * > > > * edgeMatrix, const Computation::UndirectedEdge & edge, Lang::ZBuf::IndexRepresentation * indRep )
1018 {
1019         std::pair< bool, std::list< Lang::ZBuf::IndexRepresentation * > > & p = (*edgeMatrix)[ edge ];
1020
1021         // An edge can have at most two triangles, so when one is removed, the edge cannot be queued for work.
1022         // If this was the only triangle at the edge, this will be a no-op.
1023         p.first = false;
1024
1025         std::list< Lang::ZBuf::IndexRepresentation * > & lst = p.second;
1026         std::list< Lang::ZBuf::IndexRepresentation * >::iterator i = find( lst.begin( ), lst.end( ), indRep );
1027         if( i == lst.end( ) )
1028                 {
1029                         throw Exceptions::InternalError( "Attempt to remove nonexistent edge." );
1030                 }
1031         lst.erase( i );
1032 }
1033
1034 bool
1035 Lang::ZBuf::areAligned( const Concrete::Coords2D p0, const Concrete::Coords2D p1, const Concrete::Coords2D p2 )
1036 {
1037         // We test alignment by looking at the radius of the inscribed circle.
1038         return Computation::triangleArea( p0, p1, p2 ) < Computation::theTrixelizeSplicingTol * Computation::triangleSemiPerimeter( p0, p1, p2 );
1039 }
1040
1041 bool
1042 Lang::ZBuf::areAlignedAndOrdered( const Concrete::Coords2D p0, const Concrete::Coords2D p1, const Concrete::Coords2D p2 )
1043 {
1044         // We test alignment by looking at the radius of the inscribed circle.
1045         if( Computation::triangleArea( p0, p1, p2 ) > Computation::theTrixelizeSplicingTol * Computation::triangleSemiPerimeter( p0, p1, p2 ) )
1046                 {
1047                         return false;
1048                 }
1049
1050         // We then test that p1 is between p0 and p2 by computing a projection onto an unnormalized vector...
1051         // With dHat = d.normalized( ) = d * ( 1 / d.norm( ) ), the condition
1052         //       dHat.inner( p1 - p0 ) < d.norm( )
1053         // is equivalent to
1054         //       d.inner( p1 - p0 ) < d.norm( ) * d.norm( )
1055         Concrete::UnitFloatPair d( ( p2.x_ - p0.x_ ).offtype< 1, 0 >( ), ( p2.y_ - p0.y_ ).offtype< 1, 0 >( ), bool( ) );
1056         Concrete::Length c = Concrete::inner( d, p1 - p0 );
1057         if( c < Concrete::ZERO_LENGTH )
1058                 {
1059                         return false;
1060                 }
1061         if( c > d.normSquaredThatOughtToBeOne( ) )
1062                 {
1063                         return false;
1064                 }
1065         return true;
1066 }
1067
1068 void
1069 Lang::ZBuf::addEdges( PolyIndices * poly, Computation::UndirectedEdgeMatrix< std::list< PolyIndices * > > * edgeMatrix )
1070 {
1071         PolyIndices::const_iterator src1 = poly->begin( );
1072         PolyIndices::const_iterator src2 = src1;
1073         ++src2;
1074         for( ; src2 != poly->end( ); ++src1, ++src2 )
1075                 {
1076                         addEdge( edgeMatrix, *src1, *src2, poly );
1077                 }
1078         src2 = poly->begin( );
1079         addEdge( edgeMatrix, *src1, *src2, poly );
1080 }
1081
1082 void
1083 Lang::ZBuf::addEdge( Computation::UndirectedEdgeMatrix< std::list< PolyIndices * > > * edgeMatrix, size_t ia, size_t ib, PolyIndices * indRep )
1084 {
1085         (*edgeMatrix)[ Computation::UndirectedEdge( ia, ib ) ].push_back( indRep );
1086 }
1087
1088 // Return true when extension was really performed
1089 bool
1090 Lang::ZBuf::extendPoly( PolyIndices * poly, const Computation::UndirectedEdge & edge, size_t iNew, Computation::UndirectedEdgeMatrix< std::list< PolyIndices * > > * edgeMatrix, const std::vector< Concrete::Coords2D > & cornerMem, bool convex )
1091 {
1092         const PolyIndices::iterator src0 = std::find( poly->begin( ), poly->end( ), edge.low( ) );
1093         if( src0 == poly->end( ) )
1094                 {
1095                         return false;
1096                 }
1097
1098         PolyIndices::iterator src1;
1099         PolyIndices::iterator src2;
1100         // First we try the previous point
1101         if( src0 == poly->begin( ) )
1102                 {
1103                         if( poly->back( ) == edge.high( ) )
1104                                 {
1105                                         src2 = src0;
1106                                         src1 = poly->end( );
1107                                         --src1;
1108                                         goto doEdges;
1109                                 }
1110                 }
1111         else
1112                 {
1113                         src1 = src0;
1114                         --src1;
1115                         src2 = src0;
1116                         if( *src1 == edge.high( ) )
1117                                 {
1118                                         goto doEdges;
1119                                 }
1120                 }
1121
1122         // Then we try the next point
1123         {
1124                 src1 = src0;
1125                 src2 = src0;
1126                 ++src2;
1127                 if( src2 == poly->end( ) )
1128                         {
1129                                 src2 = poly->begin( );
1130                         }
1131                 if( *src2 == edge.high( ) )
1132                         {
1133                                 goto doEdges;
1134                         }
1135         }
1136         return false;
1137
1138  doEdges:
1139         // Before we replace the old edge by the new ones we check convexity, if that's what the caller wants
1140         if( convex )
1141                 {
1142                         PolyIndices::iterator srcBack = src1;
1143                         if( srcBack == poly->begin( ) )
1144                                 {
1145                                         srcBack = poly->end( );
1146                                 }
1147                         --srcBack;
1148                         {
1149                                 Concrete::UnitFloatPair inHat = cornerMem[ *srcBack ].normalizedOrthogonal( cornerMem[ *src1 ] );
1150                                 if( Concrete::inner( inHat, cornerMem[ iNew ] - cornerMem[ *src1 ] ) < - Computation::theTrixelizeSplicingTol )
1151                                         {
1152                                                 return false;
1153                                         }
1154                         }
1155
1156                         PolyIndices::iterator srcForward = src2;
1157                         ++srcForward;
1158                         if( srcForward == poly->end( ) )
1159                                 {
1160                                         srcForward = poly->begin( );
1161                                 }
1162                         {
1163                                 Concrete::UnitFloatPair inHat = cornerMem[ *src2 ].normalizedOrthogonal( cornerMem[ *srcForward ] );
1164                                 if( Concrete::inner( inHat, cornerMem[ iNew ] - cornerMem[ *src2 ] ) < - Computation::theTrixelizeSplicingTol )
1165                                         {
1166                                                 return false;
1167                                         }
1168                         }
1169                 }
1170
1171         poly->insert( src2, iNew );
1172         (*edgeMatrix)[ edge ].remove( poly );
1173         (*edgeMatrix)[ Computation::UndirectedEdge( edge.low( ), iNew ) ].push_back( poly );
1174         (*edgeMatrix)[ Computation::UndirectedEdge( edge.high( ), iNew ) ].push_back( poly );
1175         return true;
1176 }