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::cerr << "|** Triangles in queue: " << std::endl ;
 405                                         for( ListType::const_iterator i = meetingTriangles.begin( ); i != meetingTriangles.end( ); ++i )
 406                                                 {
 407                                                         (*i)->show( std::cerr, cornerMem );
 408                                                         std::cerr << std::endl ;
 409                                                 }
 410                                         std::ostringstream oss;
 411                                         oss << "Only two (not " << meetingTriangles.size( ) << ") edges can share a pair of points!" ;
 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; ++i )
 526                         {
 527                                 for( char e = 0; e < 3; ++e )
 528                                         {
 529                                                 Computation::UndirectedEdge edge = i->getEdge( e );
 530                                                 std::list< PolyIndices * > & polys = edgeMatrix[ edge ];
 531                                                 if( polys.size( ) == 1 )        // By the way, it is an error if the size is greater than 1
 532                                                         {
 533                                                                 if( ! extendPoly( polys.front( ), edge, i->getDifferent( e ), & edgeMatrix, cornerMem, true ) ) // true for convex
 534                                                                         {
 535                                                                                 continue;
 536                                                                         }
 537                                                                 foundExtension = true;
 538                                                                 indexRepMem.erase( i );
 539                                                                 break;
 540                                                         }
 541                                         }
 542                         }
 543                 if( ! foundExtension )
 544                         {
 545                                 polyMem.push_back( indexRepMem.front( ).toPolyIndices( cornerMem ) );
 546                                 indexRepMem.pop_front( );
 547                                 PolyIndices * newPoly = & polyMem.back( );
 548                                 addEdges( newPoly, & edgeMatrix );
 549                         }
 550         }
 551
 552         while( ! polyMem.empty( ) )
 553                 {
 554                         PolyIndices & poly = polyMem.back( );
 555                         // First we remove unnecessary points on straight segments
 556                         if( poly.size( ) > 3 )
 557                                 {
 558                                         for( PolyIndices::iterator i0 = poly.begin( ); i0 != poly.end( ); ++i0 )
 559                                                 {
 560                                                         Concrete::Coords2D p0 = cornerMem[ *i0 ];
 561                                                         PolyIndices::iterator i1 = i0;
 562                                                         ++i1;
 563                                                         if( i1 == poly.end( ) )
 564                                                                 {
 565                                                                         i1 = poly.begin( );
 566                                                                 }
 567                                                         Concrete::Coords2D p1 = cornerMem[ *i1 ];
 568                                                         PolyIndices::iterator i2 = i1;
 569                                                         ++i2;
 570                                                         if( i2 == poly.end( ) )
 571                                                                 {
 572                                                                         i2 = poly.begin( );
 573                                                                 }
 574                                                         Concrete::Coords2D p2 = cornerMem[ *i2 ];
 575
 576                                                         // We do the usual inscribed circle test to test alignment
 577                                                         while( Computation::triangleArea( p0, p1, p2 ) < Computation::theTrixelizeOverlapTol * Computation::triangleSemiPerimeter( p0, p1, p2 ) )
 578                                                                 {
 579                                                                         PolyIndices::iterator iRemove = i1;
 580                                                                         i1 = i2;
 581                                                                         p1 = p2;
 582                                                                         poly.erase( iRemove );
 583                                                                         if( poly.size( ) < 3 )
 584                                                                                 {
 585                                                                                         goto done;
 586                                                                                 }
 587                                                                         ++i2;
 588                                                                         if( i2 == poly.end( ) )
 589                                                                                 {
 590                                                                                         i2 = poly.begin( );
 591                                                                                 }
 592                                                                         p2 = cornerMem[ *i2 ];
 593                                                                 }
 594                                                 }
 595                                 done:
 596                                         if( poly.size( ) < 3 )
 597                                                 {
 598                                                         polyMem.pop_back( );
 599                                                         continue;
 600                                                 }
 601                                 }
 602
 603
 604                         PolyIndices::const_iterator src = poly.begin( );
 605                         size_t i0 = *src;
 606                         ++src;
 607                         size_t i1 = *src;
 608                         ++src;
 609                         for( ; src != poly.end( ); ++src )
 610                                 {
 611                                         size_t i2 = *src;
 612                                         mergedTriangles->push_back( Computation::ZBufTriangle( painter,
 613                                                                                                                                                                                                                                                  zMap,
 614                                                                                                                                                                                                                                                  cornerMem[ i0 ], cornerMem[ i1 ], cornerMem[ i2 ] ) );
 615                                         i1 = i2;
 616                                 }
 617                         polyMem.pop_back( );
 618                 }
 619
 620 }
 621
 622 void
 623 Lang::ZBuf::trianglesToPolys( std::list< Computation::ZBufTriangle > * triangles, Lang::Clipped2D * dst )
 624 {
 625         if( triangles->empty( ) )
 626                 {
 627                         return;
 628                 }
 629
 630         const bool DEBUGME = false;
 631         if( DEBUGME )
 632         {
 633                 while( ! triangles->empty( ) )
 634                         {
 635                                 const std::vector< Concrete::Coords2D > & points = triangles->front( ).points_;
 636
 637                                 Lang::ElementaryPath2D * newPath = new Lang::ElementaryPath2D;
 638
 639                                 for( std::vector< Concrete::Coords2D >::const_iterator src = points.begin( ); src != points.end( ); ++src )
 640                                 {
 641                                         newPath->push_back( new Concrete::PathPoint2D( new Concrete::Coords2D( *src ) ) );
 642                                 }
 643                                 newPath->close( );
 644                                 dst->addSubPath( RefCountPtr< const Lang::ElementaryPath2D >( newPath ) );
 645
 646                                 triangles->pop_front( );
 647                         }
 648                 return;
 649         }
 650
 651         std::list< Lang::ZBuf::IndexRepresentation > indexRepMem;
 652         std::vector< Concrete::Coords2D > cornerMem;
 653
 654         {
 655                 std::list< Lang::ZBuf::IndexRepresentation > indexRepQueue;
 656                 while( ! triangles->empty( ) )
 657                         {
 658                                 const std::vector< Concrete::Coords2D > & points = triangles->front( ).points_;
 659                                 indexRepQueue.push_back( IndexRepresentation( ) );
 660                                 Lang::ZBuf::IndexRepresentation & current = indexRepQueue.back( );
 661                                 for( std::vector< Concrete::Coords2D >::const_iterator i = points.begin( ); i != points.end( ); ++i )
 662                                         {
 663                                                 Concrete::Coords2D tmp = *i;
 664                                                 bool reuse = false;
 665                                                 size_t memIdx = 0;
 666                                                 typedef typeof cornerMem MemType;
 667                                                 for( MemType::const_iterator j = cornerMem.begin( ); j != cornerMem.end( ); ++j, ++memIdx )
 668                                                         {
 669                                                                 if( ( tmp - *j ).norm( ) < Computation::theTrixelizeSplicingTol )
 670                                                                         {
 671                                                                                 if( current.isMember( memIdx ) )
 672                                                                                         {
 673                                                                                                 // This is the error situation when we are just about to reuse a point which is already part of this triangle.
 674                                                                                                 // The solution is to abort the creation of the current triangle.
 675                                                                                                 indexRepMem.pop_back( );
 676                                                                                                 goto abortCurrent;
 677                                                                                         }
 678                                                                                 reuse = true;
 679                                                                                 current.push_back( memIdx );
 680                                                                                 break;
 681                                                                         }
 682                                                         }
 683                                                 if( ! reuse )
 684                                                         {
 685                                                                 current.push_back( cornerMem.size( ) );
 686                                                                 cornerMem.push_back( tmp );
 687                                                         }
 688                                         }
 689                         abortCurrent:
 690                                 triangles->pop_front( );
 691                         }
 692
 693                 while( ! indexRepQueue.empty( ) )
 694                         {
 695                                 Lang::ZBuf::IndexRepresentation & current = indexRepQueue.front( );
 696
 697                                 bool wasSplit = false;
 698                                 for( char e = 0; e < 3 && ! wasSplit; ++e )
 699                                         {
 700                                                 Computation::UndirectedEdge edge = current.getEdge( e );
 701                                                 size_t iLow = edge.low( );
 702                                                 size_t iHigh = edge.high( );
 703                                                 size_t iOther = current.getDifferent( e );
 704                                                 Concrete::Coords2D pLow = cornerMem[ iLow ];
 705                                                 Concrete::Coords2D pHigh = cornerMem[ iHigh ];
 706                                                 for( size_t i = 0; i < cornerMem.size( ) && ! wasSplit; ++i )
 707                                                         {
 708                                                                 if( i == iLow || i == iHigh || i == iOther )
 709                                                                         {
 710                                                                                 continue;
 711                                                                         }
 712                                                                 if( areAlignedAndOrdered( pLow, cornerMem[ i ], pHigh ) )
 713                                                                         {
 714                                                                                 // To avoid infinite loops, we make sure that the triangles get smaller each time
 715                                                                                 // they are splitted.   Otherwise, four points in a row can cause an infinite loop.
 716
 717                                                                                 Concrete::Coords2D pOther = cornerMem[ iOther ];
 718                                                                                 Concrete::Coords2D pi = cornerMem[ i ];
 719
 720                                                                                 Concrete::Area area0 = Computation::triangleArea( pLow, pHigh, pOther );
 721
 722                                                                                 const double AREA_TOL = 1e-3;
 723
 724                                                                                 {
 725                                                                                         const double relArea = Computation::triangleArea( pOther, pLow, pi ) / area0;
 726                                                                                         if( AREA_TOL < relArea && relArea < 1 - AREA_TOL )
 727                                                                                                 {
 728                                                                                                         indexRepQueue.push_back( IndexRepresentation( ) );
 729                                                                                                         Lang::ZBuf::IndexRepresentation & newTriangle = indexRepQueue.back( );
 730                                                                                                         newTriangle.push_back( iOther );
 731                                                                                                         newTriangle.push_back( iLow );
 732                                                                                                         newTriangle.push_back( i );
 733                                                                                                         wasSplit = true;
 734                                                                                                 }
 735                                                                                 }
 736                                                                                 {
 737                                                                                         const double relArea = Computation::triangleArea( pOther, pi, pHigh ) / area0;
 738                                                                                         if( AREA_TOL < relArea && relArea < 1 - AREA_TOL )
 739                                                                                                 {
 740                                                                                                         indexRepQueue.push_back( IndexRepresentation( ) );
 741                                                                                                         Lang::ZBuf::IndexRepresentation & newTriangle = indexRepQueue.back( );
 742                                                                                                         newTriangle.push_back( iOther );
 743                                                                                                         newTriangle.push_back( i );
 744                                                                                                         newTriangle.push_back( iHigh );
 745                                                                                                         wasSplit = true;
 746                                                                                                 }
 747                                                                                 }
 748                                                                         }
 749                                                         }
 750                                         }
 751
 752                                 if( ! wasSplit )
 753                                         {
 754                                                 indexRepMem.push_back( current );
 755                                         }
 756
 757                                 indexRepQueue.pop_front( );
 758                         }
 759         }
 760
 761         const size_t pointCount = cornerMem.size( );
 762
 763         std::vector< PolyIndices > polyMem;
 764         polyMem.reserve( indexRepMem.size( ) ); // this is crucial to avoid reallocation!
 765
 766         Computation::UndirectedEdgeMatrix< std::list< PolyIndices * > > edgeMatrix( pointCount );
 767
 768         while( ! indexRepMem.empty( ) )
 769         {
 770                 bool foundExtension = false;
 771                 typedef typeof indexRepMem TrianglesType;
 772                 for( TrianglesType::iterator i = indexRepMem.begin( ); i != indexRepMem.end( ) && ! foundExtension; ++i )
 773                         {
 774                                 for( char e = 0; e < 3; ++e )
 775                                         {
 776                                                 Computation::UndirectedEdge edge = i->getEdge( e );
 777                                                 std::list< PolyIndices * > & polys = edgeMatrix[ edge ];
 778                                                 if( polys.size( ) == 1 )        // By the way, it is an error if the size is greater than 1
 779                                                         {
 780                                                                 if( ! extendPoly( polys.front( ), edge, i->getDifferent( e ), & edgeMatrix, cornerMem, false ) ) // false to not require convex
 781                                                                         {
 782                                                                                 continue;
 783                                                                         }
 784                                                                 foundExtension = true;
 785                                                                 indexRepMem.erase( i );
 786                                                                 break;
 787                                                         }
 788                                         }
 789                         }
 790                 if( ! foundExtension )
 791                         {
 792                                 polyMem.push_back( indexRepMem.front( ).toPolyIndices( cornerMem ) );
 793                                 indexRepMem.pop_front( );
 794                                 PolyIndices * newPoly = & polyMem.back( );
 795                                 addEdges( newPoly, & edgeMatrix );
 796                         }
 797         }
 798
 799         while( ! polyMem.empty( ) )
 800                 {
 801                         PolyIndices & poly = polyMem.back( );
 802                         if( poly.size( ) > 3 )
 803                                 {
 804                                         // First we remove insections, that is sequences like 4-5-4 et c.
 805                                         {
 806                                                 PolyIndices::iterator theEnd = poly.begin( );
 807                                                 PolyIndices::iterator i0 = poly.begin( );
 808                                                 while( poly.size( ) >= 5 )      // it would be pathological to have insection in polygons with four corners or less
 809                                                         {
 810                                                                 PolyIndices::iterator i1 = i0;
 811                                                                 ++i1;
 812                                                                 if( i1 == poly.end( ) )
 813                                                                         {
 814                                                                                 i1 = poly.begin( );
 815                                                                         }
 816                                                                 PolyIndices::iterator i2 = i1;
 817                                                                 ++i2;
 818                                                                 if( i2 == poly.end( ) )
 819                                                                         {
 820                                                                                 i2 = poly.begin( );
 821                                                                         }
 822
 823                                                                 if( *i0 == *i2 )
 824                                                                         {
 825                                                                                 // We found an insection
 826                                                                                 poly.erase( i1 );
 827                                                                                 poly.erase( i2 );
 828                                                                                 theEnd = i0;
 829                                                                         }
 830
 831                                                                 ++i0;
 832                                                                 if( i0 == poly.end( ) )
 833                                                                         {
 834                                                                                 i0 = poly.begin( );
 835                                                                         }
 836                                                                 if( i0 == theEnd )
 837                                                                         {
 838                                                                                 break;
 839                                                                         }
 840                                                         }
 841                                         }
 842                                         // Then we remove unnecessary points on straight segments
 843                                         for( PolyIndices::iterator i0 = poly.begin( ); i0 != poly.end( ); ++i0 )
 844                                                 {
 845                                                         Concrete::Coords2D p0 = cornerMem[ *i0 ];
 846                                                         PolyIndices::iterator i1 = i0;
 847                                                         ++i1;
 848                                                         if( i1 == poly.end( ) )
 849                                                                 {
 850                                                                         i1 = poly.begin( );
 851                                                                 }
 852                                                         Concrete::Coords2D p1 = cornerMem[ *i1 ];
 853                                                         PolyIndices::iterator i2 = i1;
 854                                                         ++i2;
 855                                                         if( i2 == poly.end( ) )
 856                                                                 {
 857                                                                         i2 = poly.begin( );
 858                                                                 }
 859                                                         Concrete::Coords2D p2 = cornerMem[ *i2 ];
 860
 861                                                         // We do the usual inscribed circle test to test alignment
 862                                                         while( Computation::triangleArea( p0, p1, p2 ) < Computation::theTrixelizeOverlapTol * Computation::triangleSemiPerimeter( p0, p1, p2 ) )
 863                                                                 {
 864                                                                         PolyIndices::iterator iRemove = i1;
 865                                                                         i1 = i2;
 866                                                                         p1 = p2;
 867                                                                         poly.erase( iRemove );
 868                                                                         if( poly.size( ) < 3 )
 869                                                                                 {
 870                                                                                         goto done;
 871                                                                                 }
 872                                                                         ++i2;
 873                                                                         if( i2 == poly.end( ) )
 874                                                                                 {
 875                                                                                         i2 = poly.begin( );
 876                                                                                 }
 877                                                                         p2 = cornerMem[ *i2 ];
 878                                                                 }
 879                                                 }
 880                                 done:
 881                                         if( poly.size( ) < 3 )
 882                                                 {
 883                                                         polyMem.pop_back( );
 884                                                         continue;
 885                                                 }
 886                                 }
 887
 888
 889                         {
 890                                 Lang::ElementaryPath2D * newPath = new Lang::ElementaryPath2D;
 891
 892                                 for( PolyIndices::const_iterator src = poly.begin( ); src != poly.end( ); ++src )
 893                                 {
 894                                         newPath->push_back( new Concrete::PathPoint2D( new Concrete::Coords2D( cornerMem[ *src ] ) ) );
 895                                 }
 896                                 newPath->close( );
 897                                 dst->addSubPath( RefCountPtr< const Lang::ElementaryPath2D >( newPath ) );
 898                         }
 899                         polyMem.pop_back( );
 900                 }
 901
 902 }
 903
 904 bool
 905 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 )
 906 {
 907         typedef typeof (*edgeMatrix)[ edge ].second ListType;
 908         ListType & trianglesAtEdge = (*edgeMatrix)[ edge ].second;
 909         trianglesAtEdge.push_back( indRep );
 910
 911         bool alreadyInQueue = false;
 912
 913         if( trianglesAtEdge.size( ) > 2 )
 914                 {
 915                         // This is a problematic case.  I've written about it elsewhere.        The decision I took was to discard
 916                         // the smaller if the triangles being on the same side of the edge.
 917                         // Since we only do this when there are already three triangles sharing the edge, there's a risk
 918                         // that all three are on the same side. Then we may even drop two of them.
 919
 920                         alreadyInQueue = true;
 921
 922                         Concrete::UnitFloatPair edgeNormalNotNormalized = cornerMem[ edge.low( ) ].unNormalizedOrthogonal( cornerMem[ edge.high( ) ] );
 923                         double m = Concrete::innerScalar( cornerMem[ edge.low( ) ], edgeNormalNotNormalized );
 924
 925                         Lang::ZBuf::IndexRepresentation * posMax = 0;
 926                         Lang::ZBuf::IndexRepresentation * negMax = 0;
 927                         Concrete::Area posArea = -1;
 928                         Concrete::Area negArea = -1;
 929
 930                         for( ListType::iterator i = trianglesAtEdge.begin( ); i != trianglesAtEdge.end( ); ++i )
 931                                 {
 932                                         Concrete::Area tmpArea = Computation::triangleArea( cornerMem[ (*i)->i0 ], cornerMem[ (*i)->i1 ], cornerMem[ (*i)->i2 ] );
 933                                         if( Concrete::innerScalar( cornerMem[ (*i)->getDifferent( edge.low( ), edge.high( ) ) ], edgeNormalNotNormalized ) > m )
 934                                                 {
 935                                                         if( tmpArea > posArea )
 936                                                                 {
 937                                                                         posMax = *i;
 938                                                                         posArea = tmpArea;
 939                                                                 }
 940                                                 }
 941                                         else
 942                                                 {
 943                                                         if( tmpArea > negArea )
 944                                                                 {
 945                                                                         negMax = *i;
 946                                                                         negArea = tmpArea;
 947                                                                 }
 948                                                 }
 949                                 }
 950
 951                         // *** It is required that the new triangle appears at the end of trianglesAtEdge in case we add it. ***
 952
 953                         // We temporarily remove the newly added triangle...
 954                         trianglesAtEdge.pop_back( );
 955
 956                         // Now we remove those of the old triangles that were knocked out.
 957                         // Note that removing triangles from edges clears the job flags at their edges, so we don't have to do that here.
 958                         {
 959                                 // We must be careful since trianglesAtEdge might change when triangles are removed!
 960                                 ListType::iterator i = trianglesAtEdge.begin( );
 961                                 Lang::ZBuf::IndexRepresentation * t1 = *i;
 962                                 ++i;
 963                                 Lang::ZBuf::IndexRepresentation * t2 = *i;
 964
 965                                 if( t1 != posMax && t1 != negMax )
 966                                         {
 967                                                 t1->removeEdges( edgeMatrix );
 968                                         }
 969                                 if( t2 != posMax && t2 != negMax )
 970                                         {
 971                                                 t2->removeEdges( edgeMatrix );
 972                                         }
 973                         }
 974
 975                         // Next, we consider the new triangle.
 976                         if( indRep != posMax && indRep != negMax )
 977                                 {
 978                                         // The new triangle couldn't compete with the old ones, so we shall return false.
 979                                         return false;
 980                                 }
 981
 982                         // ... and now we put it back if it passed the test.
 983                         trianglesAtEdge.push_back( indRep );
 984                 }
 985
 986         if( jobQueue != 0 )
 987                 {
 988                         if( trianglesAtEdge.size( ) == 2 && ! (*edgeMatrix)[ edge ].first )
 989                                 {
 990                                         if( ! alreadyInQueue )
 991                                                 {
 992                                                         jobQueue->push_back( edge );
 993                                                 }
 994                                         (*edgeMatrix)[ edge ].first = true; // Flag that this pair is in the queue.
 995                                 }
 996                 }
 997         return true;
 998 }
 999
1000 void
1001 Lang::ZBuf::removeEdge( Computation::UndirectedEdgeMatrix< std::pair< bool, std::list< Lang::ZBuf::IndexRepresentation * > > > * edgeMatrix, const Computation::UndirectedEdge & edge, Lang::ZBuf::IndexRepresentation * indRep )
1002 {
1003         std::pair< bool, std::list< Lang::ZBuf::IndexRepresentation * > > & p = (*edgeMatrix)[ edge ];
1004
1005         // An edge can have at most two triangles, so when one is removed, the edge cannot be queued for work.
1006         // If this was the only triangle at the edge, this will be a no-op.
1007         p.first = false;
1008
1009         std::list< Lang::ZBuf::IndexRepresentation * > & lst = p.second;
1010         std::list< Lang::ZBuf::IndexRepresentation * >::iterator i = find( lst.begin( ), lst.end( ), indRep );
1011         if( i == lst.end( ) )
1012                 {
1013                         throw Exceptions::InternalError( "Attempt to remove nonexistent edge." );
1014                 }
1015         lst.erase( i );
1016 }
1017
1018 bool
1019 Lang::ZBuf::areAligned( const Concrete::Coords2D p0, const Concrete::Coords2D p1, const Concrete::Coords2D p2 )
1020 {
1021         // We test alignment by looking at the radius of the inscribed circle.
1022         return Computation::triangleArea( p0, p1, p2 ) < Computation::theTrixelizeSplicingTol * Computation::triangleSemiPerimeter( p0, p1, p2 );
1023 }
1024
1025 bool
1026 Lang::ZBuf::areAlignedAndOrdered( const Concrete::Coords2D p0, const Concrete::Coords2D p1, const Concrete::Coords2D p2 )
1027 {
1028         // We test alignment by looking at the radius of the inscribed circle.
1029         if( Computation::triangleArea( p0, p1, p2 ) > Computation::theTrixelizeSplicingTol * Computation::triangleSemiPerimeter( p0, p1, p2 ) )
1030                 {
1031                         return false;
1032                 }
1033
1034         // We then test that p1 is between p0 and p2 by computing a projection onto an unnormalized vector...
1035         // With dHat = d.normalized( ) = d * ( 1 / d.norm( ) ), the condition
1036         //       dHat.inner( p1 - p0 ) < d.norm( )
1037         // is equivalent to
1038         //       d.inner( p1 - p0 ) < d.norm( ) * d.norm( )
1039         Concrete::UnitFloatPair d( ( p2.x_ - p0.x_ ).offtype< 1, 0 >( ), ( p2.y_ - p0.y_ ).offtype< 1, 0 >( ), bool( ) );
1040         Concrete::Length c = Concrete::inner( d, p1 - p0 );
1041         if( c < Concrete::ZERO_LENGTH )
1042                 {
1043                         return false;
1044                 }
1045         if( c > d.normSquaredThatOughtToBeOne( ) )
1046                 {
1047                         return false;
1048                 }
1049         return true;
1050 }
1051
1052 void
1053 Lang::ZBuf::addEdges( PolyIndices * poly, Computation::UndirectedEdgeMatrix< std::list< PolyIndices * > > * edgeMatrix )
1054 {
1055         PolyIndices::const_iterator src1 = poly->begin( );
1056         PolyIndices::const_iterator src2 = src1;
1057         ++src2;
1058         for( ; src2 != poly->end( ); ++src1, ++src2 )
1059                 {
1060                         addEdge( edgeMatrix, *src1, *src2, poly );
1061                 }
1062         src2 = poly->begin( );
1063         addEdge( edgeMatrix, *src1, *src2, poly );
1064 }
1065
1066 void
1067 Lang::ZBuf::addEdge( Computation::UndirectedEdgeMatrix< std::list< PolyIndices * > > * edgeMatrix, size_t ia, size_t ib, PolyIndices * indRep )
1068 {
1069         (*edgeMatrix)[ Computation::UndirectedEdge( ia, ib ) ].push_back( indRep );
1070 }
1071
1072 // Return true when extension was really performed
1073 bool
1074 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 )
1075 {
1076         const PolyIndices::iterator src0 = std::find( poly->begin( ), poly->end( ), edge.low( ) );
1077         if( src0 == poly->end( ) )
1078                 {
1079                         return false;
1080                 }
1081
1082         PolyIndices::iterator src1;
1083         PolyIndices::iterator src2;
1084         // First we try the previous point
1085         if( src0 == poly->begin( ) )
1086                 {
1087                         if( poly->back( ) == edge.high( ) )
1088                                 {
1089                                         src2 = src0;
1090                                         src1 = poly->end( );
1091                                         --src1;
1092                                         goto doEdges;
1093                                 }
1094                 }
1095         else
1096                 {
1097                         src1 = src0;
1098                         --src1;
1099                         src2 = src0;
1100                         if( *src1 == edge.high( ) )
1101                                 {
1102                                         goto doEdges;
1103                                 }
1104                 }
1105
1106         // Then we try the next point
1107         {
1108                 src1 = src0;
1109                 src2 = src0;
1110                 ++src2;
1111                 if( src2 == poly->end( ) )
1112                         {
1113                                 src2 = poly->begin( );
1114                         }
1115                 if( *src2 == edge.high( ) )
1116                         {
1117                                 goto doEdges;
1118                         }
1119         }
1120         return false;
1121
1122  doEdges:
1123         // Before we replace the old edge by the new ones we check convexity, if that's what the caller wants
1124         if( convex )
1125                 {
1126                         PolyIndices::iterator srcBack = src1;
1127                         if( srcBack == poly->begin( ) )
1128                                 {
1129                                         srcBack = poly->end( );
1130                                 }
1131                         --srcBack;
1132                         {
1133                                 Concrete::UnitFloatPair inHat = cornerMem[ *srcBack ].normalizedOrthogonal( cornerMem[ *src1 ] );
1134                                 if( Concrete::inner( inHat, cornerMem[ iNew ] - cornerMem[ *src1 ] ) < - Computation::theTrixelizeSplicingTol )
1135                                         {
1136                                                 return false;
1137                                         }
1138                         }
1139
1140                         PolyIndices::iterator srcForward = src2;
1141                         ++srcForward;
1142                         if( srcForward == poly->end( ) )
1143                                 {
1144                                         srcForward = poly->begin( );
1145                                 }
1146                         {
1147                                 Concrete::UnitFloatPair inHat = cornerMem[ *src2 ].normalizedOrthogonal( cornerMem[ *srcForward ] );
1148                                 if( Concrete::inner( inHat, cornerMem[ iNew ] - cornerMem[ *src2 ] ) < - Computation::theTrixelizeSplicingTol )
1149                                         {
1150                                                 return false;
1151                                         }
1152                         }
1153                 }
1154
1155         poly->insert( src2, iNew );
1156         (*edgeMatrix)[ edge ].remove( poly );
1157         (*edgeMatrix)[ Computation::UndirectedEdge( edge.low( ), iNew ) ].push_back( poly );
1158         (*edgeMatrix)[ Computation::UndirectedEdge( edge.high( ), iNew ) ].push_back( poly );
1159         return true;
1160 }