src/BRepMesh/BRepMesh_Delaun.cxx

   1 // Created on: 1993-05-12
   2 // Created by: Didier PIFFAULT
   3 // Copyright (c) 1993-1999 Matra Datavision
   4 // Copyright (c) 1999-2014 OPEN CASCADE SAS
   5 //
   6 // This file is part of Open CASCADE Technology software library.
   7 //
   8 // This library is free software; you can redistribute it and/or modify it under
   9 // the terms of the GNU Lesser General Public License version 2.1 as published
  10 // by the Free Software Foundation, with special exception defined in the file
  11 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
  12 // distribution for complete text of the license and disclaimer of any warranty.
  13 //
  14 // Alternatively, this file may be used under the terms of Open CASCADE
  15 // commercial license or contractual agreement.
  16
  17 #include <BRepMesh_Delaun.hxx>
  18
  19 #include <gp.hxx>
  20 #include <gp_XY.hxx>
  21 #include <gp_Pnt2d.hxx>
  22 #include <gp_Vec2d.hxx>
  23
  24 #include <Precision.hxx>
  25 #include <Bnd_Box2d.hxx>
  26 #include <Bnd_B2d.hxx>
  27
  28 #include <BRepMesh_SelectorOfDataStructureOfDelaun.hxx>
  29
  30 #include <BRepMesh_Edge.hxx>
  31 #include <BRepMesh_Vertex.hxx>
  32 #include <BRepMesh_Triangle.hxx>
  33
  34 #include <NCollection_Vector.hxx>
  35
  36 #include <algorithm>
  37 #include <stack>
  38
  39 const Standard_Real AngDeviation1Deg  = M_PI/180.;
  40 const Standard_Real AngDeviation90Deg = 90 * AngDeviation1Deg;
  41 const Standard_Real Angle2PI          = 2 * M_PI;
  42
  43 const Standard_Real Precision    = Precision::PConfusion();
  44 const Standard_Real Precision2   = Precision * Precision;
  45
  46 namespace {
  47   //! Sort two points in projection on vector (1,1)
  48   struct ComparatorOfVertexOfDelaun
  49   {
  50     bool operator() (const BRepMesh_Vertex& theLeft, const BRepMesh_Vertex& theRight)
  51     {
  52       return theLeft.Coord().X() + theLeft.Coord().Y() < theRight.Coord().X() + theRight.Coord().Y();
  53     }
  54   };
  55
  56   //! Sort two points in projection on vector (1,1)
  57   struct ComparatorOfIndexedVertexOfDelaun
  58   {
  59   public:
  60     ComparatorOfIndexedVertexOfDelaun (const Handle(BRepMesh_DataStructureOfDelaun)& theDS)
  61       : myStructure(theDS) {}
  62
  63     bool operator() (Standard_Integer theLeft, Standard_Integer theRight)
  64     {
  65       const BRepMesh_Vertex& aLeft  = myStructure->GetNode(theLeft);
  66       const BRepMesh_Vertex& aRight = myStructure->GetNode(theRight);
  67       return ComparatorOfVertexOfDelaun() (aLeft, aRight);
  68     }
  69
  70   private:
  71     Handle(BRepMesh_DataStructureOfDelaun) myStructure;
  72   };
  73
  74   inline void UpdateBndBox(const gp_XY& thePnt1, const gp_XY& thePnt2, Bnd_B2d& theBox)
  75   {
  76     theBox.Add( thePnt1 );
  77     theBox.Add( thePnt2 );
  78     theBox.Enlarge(Precision);
  79   }
  80 } // anonymous namespace
  81
  82 //=======================================================================
  83 //function : BRepMesh_Delaun
  84 //purpose  :
  85 //=======================================================================
  86 BRepMesh_Delaun::BRepMesh_Delaun (
  87   const Handle(BRepMesh_DataStructureOfDelaun)& theOldMesh,
  88   const Standard_Integer                        theCellsCountU,
  89   const Standard_Integer                        theCellsCountV,
  90   const Standard_Boolean                        isFillCircles)
  91 : myMeshData ( theOldMesh ),
  92   myCircles (new NCollection_IncAllocator(
  93              IMeshData::MEMORY_BLOCK_SIZE_HUGE))
  94 {
  95   memset (mySupVert, 0, sizeof (mySupVert));
  96   if (isFillCircles)
  97   {
  98     InitCirclesTool (theCellsCountU, theCellsCountV);
  99   }
 100 }
 101
 102 //=======================================================================
 103 //function : BRepMesh_Delaun
 104 //purpose  : Creates the triangulation with an empty Mesh data structure
 105 //=======================================================================
 106 BRepMesh_Delaun::BRepMesh_Delaun(IMeshData::Array1OfVertexOfDelaun& theVertices)
 107 : myCircles (theVertices.Length(), new NCollection_IncAllocator(
 108              IMeshData::MEMORY_BLOCK_SIZE_HUGE))
 109 {
 110   memset (mySupVert, 0, sizeof (mySupVert));
 111   if ( theVertices.Length() > 2 )
 112   {
 113     myMeshData = new BRepMesh_DataStructureOfDelaun(
 114       new NCollection_IncAllocator(IMeshData::MEMORY_BLOCK_SIZE_HUGE),
 115       theVertices.Length() );
 116     Init( theVertices );
 117   }
 118 }
 119
 120 //=======================================================================
 121 //function : BRepMesh_Delaun
 122 //purpose  : Creates the triangulation with and existent Mesh data structure
 123 //=======================================================================
 124 BRepMesh_Delaun::BRepMesh_Delaun(
 125   const Handle(BRepMesh_DataStructureOfDelaun)& theOldMesh,
 126   IMeshData::Array1OfVertexOfDelaun&            theVertices)
 127 : myMeshData( theOldMesh ),
 128   myCircles ( theVertices.Length(), new NCollection_IncAllocator(
 129              IMeshData::MEMORY_BLOCK_SIZE_HUGE))
 130 {
 131   memset (mySupVert, 0, sizeof (mySupVert));
 132   if ( theVertices.Length() > 2 )
 133   {
 134     Init( theVertices );
 135   }
 136 }
 137
 138 //=======================================================================
 139 //function : BRepMesh_Delaun
 140 //purpose  : Creates the triangulation with and existent Mesh data structure
 141 //=======================================================================
 142 BRepMesh_Delaun::BRepMesh_Delaun(
 143   const Handle(BRepMesh_DataStructureOfDelaun)& theOldMesh,
 144   IMeshData::VectorOfInteger&                   theVertexIndices)
 145 : myMeshData( theOldMesh ),
 146   myCircles ( theVertexIndices.Length(), new NCollection_IncAllocator(
 147              IMeshData::MEMORY_BLOCK_SIZE_HUGE))
 148 {
 149   memset (mySupVert, 0, sizeof (mySupVert));
 150   perform(theVertexIndices);
 151 }
 152
 153 //=======================================================================
 154 //function : BRepMesh_Delaun
 155 //purpose  : Creates the triangulation with and existent Mesh data structure
 156 //=======================================================================
 157 BRepMesh_Delaun::BRepMesh_Delaun (const Handle (BRepMesh_DataStructureOfDelaun)& theOldMesh,
 158                                   IMeshData::VectorOfInteger&                    theVertexIndices,
 159                                   const Standard_Integer                         theCellsCountU,
 160                                   const Standard_Integer                         theCellsCountV)
 161 : myMeshData (theOldMesh),
 162   myCircles (theVertexIndices.Length (), new NCollection_IncAllocator(
 163              IMeshData::MEMORY_BLOCK_SIZE_HUGE))
 164 {
 165   memset (mySupVert, 0, sizeof (mySupVert));
 166   perform (theVertexIndices, theCellsCountU, theCellsCountV);
 167 }
 168
 169 //=======================================================================
 170 //function : Init
 171 //purpose  : Initializes the triangulation with an Array of Vertex
 172 //=======================================================================
 173 void BRepMesh_Delaun::Init(IMeshData::Array1OfVertexOfDelaun& theVertices)
 174 {
 175   Standard_Integer aLowerIdx  = theVertices.Lower();
 176   Standard_Integer anUpperIdx = theVertices.Upper();
 177   IMeshData::VectorOfInteger aVertexIndexes(theVertices.Size());
 178
 179   Standard_Integer anIndex = aLowerIdx;
 180   for ( ; anIndex <= anUpperIdx; ++anIndex )
 181   {
 182     aVertexIndexes.Append(myMeshData->AddNode( theVertices( anIndex ) ));
 183   }
 184
 185   perform( aVertexIndexes );
 186 }
 187
 188 //=======================================================================
 189 //function : InitCirclesTool
 190 //purpose  :
 191 //=======================================================================
 192 void BRepMesh_Delaun::InitCirclesTool (const Standard_Integer theCellsCountU,
 193                                        const Standard_Integer theCellsCountV)
 194 {
 195   Bnd_Box2d aBox;
 196   for (Standard_Integer aNodeIt = 1; aNodeIt <= myMeshData->NbNodes(); ++aNodeIt)
 197   {
 198     aBox.Add (gp_Pnt2d (GetVertex (aNodeIt).Coord ()));
 199   }
 200   aBox.Enlarge (Precision);
 201
 202   initCirclesTool (aBox, theCellsCountU, theCellsCountV);
 203
 204   IMeshData::IteratorOfMapOfInteger aTriangleIt (myMeshData->ElementsOfDomain());
 205   for (; aTriangleIt.More(); aTriangleIt.Next())
 206   {
 207     Standard_Integer aNodesIndices[3];
 208     const BRepMesh_Triangle& aTriangle = myMeshData->GetElement (aTriangleIt.Key());
 209     myMeshData->ElementNodes (aTriangle, aNodesIndices);
 210     myCircles.Bind (aTriangleIt.Key(),
 211                     GetVertex( aNodesIndices[0] ).Coord(),
 212                     GetVertex( aNodesIndices[1] ).Coord(),
 213                     GetVertex( aNodesIndices[2] ).Coord());
 214   }
 215 }
 216
 217 //=======================================================================
 218 //function : initCirclesTool
 219 //purpose  :
 220 //=======================================================================
 221 void BRepMesh_Delaun::initCirclesTool (const Bnd_Box2d&       theBox,
 222                                        const Standard_Integer theCellsCountU,
 223                                        const Standard_Integer theCellsCountV)
 224 {
 225   Standard_Real aMinX, aMinY, aMaxX, aMaxY;
 226   theBox.Get  ( aMinX, aMinY, aMaxX, aMaxY );
 227   const Standard_Real aDeltaX = aMaxX - aMinX;
 228   const Standard_Real aDeltaY = aMaxY - aMinY;
 229
 230   Standard_Integer aScaler = 2;
 231   if ( myMeshData->NbNodes() > 100 )
 232   {
 233     aScaler = 5;
 234   }
 235   else if( myMeshData->NbNodes() > 1000 )
 236   {
 237     aScaler = 7;
 238   }
 239
 240   myCircles.SetMinMaxSize( gp_XY( aMinX, aMinY ), gp_XY( aMaxX, aMaxY ) );
 241   myCircles.SetCellSize  ( aDeltaX / Max (theCellsCountU, aScaler),
 242                            aDeltaY / Max (theCellsCountV, aScaler));
 243 }
 244
 245 //=======================================================================
 246 //function : perform
 247 //purpose  : Create super mesh and run triangulation procedure
 248 //=======================================================================
 249 void BRepMesh_Delaun::perform(IMeshData::VectorOfInteger& theVertexIndices,
 250                               const Standard_Integer      theCellsCountU /* = -1 */,
 251                               const Standard_Integer      theCellsCountV /* = -1 */)
 252 {
 253   if (theVertexIndices.Length () <= 2)
 254   {
 255     return;
 256   }
 257
 258   Bnd_Box2d aBox;
 259   Standard_Integer anIndex = theVertexIndices.Lower ();
 260   Standard_Integer anUpper = theVertexIndices.Upper ();
 261   for (; anIndex <= anUpper; ++anIndex)
 262   {
 263     aBox.Add (gp_Pnt2d (GetVertex (theVertexIndices (anIndex)).Coord ()));
 264   }
 265
 266   aBox.Enlarge (Precision);
 267
 268   initCirclesTool (aBox, theCellsCountU, theCellsCountV);
 269   superMesh       (aBox);
 270
 271   ComparatorOfIndexedVertexOfDelaun aCmp(myMeshData);
 272   std::make_heap(theVertexIndices.begin(), theVertexIndices.end(), aCmp);
 273   std::sort_heap(theVertexIndices.begin(), theVertexIndices.end(), aCmp);
 274
 275   compute( theVertexIndices );
 276 }
 277
 278 //=======================================================================
 279 //function : superMesh
 280 //purpose  : Build the super mesh
 281 //=======================================================================
 282 void BRepMesh_Delaun::superMesh(const Bnd_Box2d& theBox)
 283 {
 284   Standard_Real aMinX, aMinY, aMaxX, aMaxY;
 285   theBox.Get  ( aMinX, aMinY, aMaxX, aMaxY );
 286   Standard_Real aDeltaX = aMaxX - aMinX;
 287   Standard_Real aDeltaY = aMaxY - aMinY;
 288
 289   Standard_Real aDeltaMin = Min( aDeltaX, aDeltaY );
 290   Standard_Real aDeltaMax = Max( aDeltaX, aDeltaY );
 291   Standard_Real aDelta    = aDeltaX + aDeltaY;
 292
 293   mySupVert[0] = myMeshData->AddNode(
 294     BRepMesh_Vertex( ( aMinX + aMaxX ) / 2, aMaxY + aDeltaMax, BRepMesh_Free ) );
 295
 296   mySupVert[1] = myMeshData->AddNode(
 297     BRepMesh_Vertex( aMinX - aDelta, aMinY - aDeltaMin, BRepMesh_Free ) );
 298
 299   mySupVert[2] = myMeshData->AddNode(
 300     BRepMesh_Vertex( aMaxX + aDelta, aMinY - aDeltaMin, BRepMesh_Free ) );
 301
 302   Standard_Integer e[3];
 303   Standard_Boolean o[3];
 304   for (Standard_Integer aNodeId = 0; aNodeId < 3; ++aNodeId)
 305   {
 306     Standard_Integer aFirstNode = aNodeId;
 307     Standard_Integer aLastNode  = (aNodeId + 1) % 3;
 308     Standard_Integer aLinkIndex = myMeshData->AddLink( BRepMesh_Edge(
 309       mySupVert[aFirstNode], mySupVert[aLastNode], BRepMesh_Free ) );
 310
 311     e[aNodeId] = Abs(aLinkIndex);
 312     o[aNodeId] = (aLinkIndex > 0);
 313   }
 314
 315   mySupTrian = BRepMesh_Triangle(e, o, BRepMesh_Free);
 316 }
 317
 318 //=======================================================================
 319 //function : deleteTriangle
 320 //purpose  : Deletes the triangle with the given index and adds the free
 321 //           edges into the map.
 322 //           When an edge is suppressed more than one time it is destroyed.
 323 //=======================================================================
 324 void BRepMesh_Delaun::deleteTriangle(const Standard_Integer          theIndex,
 325                                      IMeshData::MapOfIntegerInteger& theLoopEdges )
 326 {
 327   if (!myCircles.IsEmpty())
 328   {
 329     myCircles.Delete (theIndex);
 330   }
 331
 332   const BRepMesh_Triangle& aElement = GetTriangle(theIndex);
 333   const Standard_Integer(&e)[3] = aElement.myEdges;
 334   const Standard_Boolean(&o)[3] = aElement.myOrientations;
 335
 336   myMeshData->RemoveElement( theIndex );
 337
 338   for ( Standard_Integer i = 0; i < 3; ++i )
 339   {
 340     if ( !theLoopEdges.Bind( e[i], o[i] ) )
 341     {
 342       theLoopEdges.UnBind( e[i] );
 343       myMeshData->RemoveLink( e[i] );
 344     }
 345   }
 346 }
 347
 348 //=======================================================================
 349 //function : compute
 350 //purpose  : Computes the triangulation and add the vertices edges and
 351 //           triangles to the Mesh data structure
 352 //=======================================================================
 353 void BRepMesh_Delaun::compute(IMeshData::VectorOfInteger& theVertexIndexes)
 354 {
 355   // Insertion of edges of super triangles in the list of free edges:
 356   Handle(NCollection_IncAllocator) aAllocator = new NCollection_IncAllocator(
 357     IMeshData::MEMORY_BLOCK_SIZE_HUGE);
 358
 359   IMeshData::MapOfIntegerInteger aLoopEdges(10, aAllocator);
 360   const Standard_Integer(&e)[3] = mySupTrian.myEdges;
 361
 362   aLoopEdges.Bind( e[0], Standard_True );
 363   aLoopEdges.Bind( e[1], Standard_True );
 364   aLoopEdges.Bind( e[2], Standard_True );
 365
 366   if ( theVertexIndexes.Length() > 0 )
 367   {
 368     // Creation of 3 trianglers with the first node and the edges of the super triangle:
 369     Standard_Integer anVertexIdx = theVertexIndexes.Lower();
 370     createTriangles( theVertexIndexes( anVertexIdx ), aLoopEdges );
 371
 372     // Add other nodes to the mesh
 373     createTrianglesOnNewVertices( theVertexIndexes );
 374   }
 375
 376   // Destruction of triangles containing a top of the super triangle
 377   BRepMesh_SelectorOfDataStructureOfDelaun aSelector( myMeshData );
 378   for (Standard_Integer aSupVertId = 0; aSupVertId < 3; ++aSupVertId)
 379     aSelector.NeighboursOfNode( mySupVert[aSupVertId] );
 380
 381   aLoopEdges.Clear();
 382   IMeshData::IteratorOfMapOfInteger aFreeTriangles( aSelector.Elements() );
 383   for ( ; aFreeTriangles.More(); aFreeTriangles.Next() )
 384     deleteTriangle( aFreeTriangles.Key(), aLoopEdges );
 385
 386   // All edges that remain free are removed from aLoopEdges;
 387   // only the boundary edges of the triangulation remain there
 388   IMeshData::MapOfIntegerInteger::Iterator aFreeEdges( aLoopEdges );
 389   for ( ; aFreeEdges.More(); aFreeEdges.Next() )
 390   {
 391     if ( myMeshData->ElementsConnectedTo( aFreeEdges.Key() ).IsEmpty() )
 392       myMeshData->RemoveLink( aFreeEdges.Key() );
 393   }
 394
 395   // The tops of the super triangle are destroyed
 396   for (Standard_Integer aSupVertId = 0; aSupVertId < 3; ++aSupVertId)
 397     myMeshData->RemoveNode( mySupVert[aSupVertId] );
 398 }
 399
 400 //=======================================================================
 401 //function : createTriangles
 402 //purpose  : Creates the triangles beetween the node and the polyline.
 403 //=======================================================================
 404 void BRepMesh_Delaun::createTriangles(const Standard_Integer          theVertexIndex,
 405                                       IMeshData::MapOfIntegerInteger& thePoly)
 406 {
 407   IMeshData::ListOfInteger aLoopEdges, anExternalEdges;
 408   const gp_XY& aVertexCoord = myMeshData->GetNode( theVertexIndex ).Coord();
 409
 410   IMeshData::MapOfIntegerInteger::Iterator anEdges( thePoly );
 411   for ( ; anEdges.More(); anEdges.Next() )
 412   {
 413     Standard_Integer     anEdgeId = anEdges.Key();
 414     const BRepMesh_Edge& anEdge   = GetEdge( anEdgeId );
 415
 416     Standard_Boolean isPositive = thePoly( anEdgeId ) != 0;
 417
 418     Standard_Integer aNodes[3];
 419     if ( isPositive )
 420     {
 421       aNodes[0] = anEdge.FirstNode();
 422       aNodes[2] = anEdge.LastNode();
 423     }
 424     else
 425     {
 426       aNodes[0] = anEdge.LastNode();
 427       aNodes[2] = anEdge.FirstNode();
 428     }
 429     aNodes[1] = theVertexIndex;
 430
 431     const BRepMesh_Vertex& aFirstVertex = GetVertex( aNodes[0] );
 432     const BRepMesh_Vertex& aLastVertex  = GetVertex( aNodes[2] );
 433
 434     gp_XY anEdgeDir( aLastVertex.Coord() - aFirstVertex.Coord() );
 435     Standard_Real anEdgeLen = anEdgeDir.Modulus();
 436     if ( anEdgeLen < Precision )
 437       continue;
 438
 439     anEdgeDir.SetCoord( anEdgeDir.X() / anEdgeLen,
 440                         anEdgeDir.Y() / anEdgeLen );
 441
 442     gp_XY aFirstLinkDir( aFirstVertex.Coord() - aVertexCoord );
 443     gp_XY aLastLinkDir ( aVertexCoord         - aLastVertex.Coord() );
 444
 445     Standard_Real aDist12 = aFirstLinkDir ^ anEdgeDir;
 446     Standard_Real aDist23 = anEdgeDir     ^ aLastLinkDir;
 447     if ( Abs( aDist12 ) < Precision ||
 448          Abs( aDist23 ) < Precision )
 449     {
 450       continue;
 451     }
 452
 453     BRepMesh_Edge aFirstLink( aNodes[1], aNodes[0], BRepMesh_Free );
 454     BRepMesh_Edge aLastLink ( aNodes[2], aNodes[1], BRepMesh_Free );
 455
 456     Standard_Integer anEdgesInfo[3] = {
 457       myMeshData->AddLink( aFirstLink ),
 458       isPositive ? anEdgeId : -anEdgeId,
 459       myMeshData->AddLink( aLastLink ) };
 460
 461     Standard_Boolean isSensOK = (aDist12 > 0. && aDist23 > 0.);
 462     if (isSensOK)
 463     {
 464       Standard_Integer anEdgeIds[3];
 465       Standard_Boolean anEdgesOri[3];
 466       for ( Standard_Integer aTriLinkIt = 0; aTriLinkIt < 3; ++aTriLinkIt )
 467       {
 468         const Standard_Integer& anEdgeInfo = anEdgesInfo[aTriLinkIt];
 469         anEdgeIds[aTriLinkIt]    = Abs( anEdgeInfo );
 470         anEdgesOri[aTriLinkIt] = anEdgeInfo > 0;
 471       }
 472
 473       addTriangle(anEdgeIds, anEdgesOri, aNodes );
 474     }
 475     else
 476     {
 477       if ( isPositive )
 478         aLoopEdges.Append(  anEdges.Key() );
 479       else
 480         aLoopEdges.Append( -anEdges.Key() );
 481
 482       if ( aFirstLinkDir.SquareModulus() > aLastLinkDir.SquareModulus() )
 483         anExternalEdges.Append( Abs( anEdgesInfo[0] ) );
 484       else
 485         anExternalEdges.Append( Abs( anEdgesInfo[2] ) );
 486     }
 487   }
 488
 489   thePoly.Clear();
 490   while ( !anExternalEdges.IsEmpty() )
 491   {
 492     const BRepMesh_PairOfIndex& aPair =
 493       myMeshData->ElementsConnectedTo( Abs( anExternalEdges.First() ) );
 494
 495
 496     if ( !aPair.IsEmpty() )
 497       deleteTriangle( aPair.FirstIndex(), thePoly );
 498
 499     anExternalEdges.RemoveFirst();
 500   }
 501
 502   for ( anEdges.Initialize( thePoly ); anEdges.More(); anEdges.Next() )
 503   {
 504     if ( myMeshData->ElementsConnectedTo( anEdges.Key() ).IsEmpty() )
 505       myMeshData->RemoveLink( anEdges.Key() );
 506   }
 507
 508   while ( !aLoopEdges.IsEmpty() )
 509   {
 510     const BRepMesh_Edge& anEdge = GetEdge( Abs( aLoopEdges.First() ) );
 511     if ( anEdge.Movability() != BRepMesh_Deleted )
 512     {
 513       Standard_Integer anEdgeIdx = aLoopEdges.First();
 514       meshLeftPolygonOf( Abs( anEdgeIdx ), ( anEdgeIdx > 0 ) );
 515     }
 516
 517     aLoopEdges.RemoveFirst();
 518   }
 519 }
 520
 521 //=======================================================================
 522 //function : createTrianglesOnNewVertices
 523 //purpose  : Creation of triangles from the new nodes
 524 //=======================================================================
 525 void BRepMesh_Delaun::createTrianglesOnNewVertices(
 526   IMeshData::VectorOfInteger& theVertexIndexes)
 527 {
 528   Handle(NCollection_IncAllocator) aAllocator =
 529     new NCollection_IncAllocator(IMeshData::MEMORY_BLOCK_SIZE_HUGE);
 530
 531   Standard_Real aTolU, aTolV;
 532   myMeshData->Data()->GetTolerance(aTolU, aTolV);
 533   const Standard_Real aSqTol = aTolU * aTolU + aTolV * aTolV;
 534
 535   // Insertion of nodes :
 536   Standard_Boolean isModify = Standard_True;
 537
 538   Standard_Integer anIndex = theVertexIndexes.Lower();
 539   Standard_Integer anUpper = theVertexIndexes.Upper();
 540   for( ; anIndex <= anUpper; ++anIndex )
 541   {
 542     aAllocator->Reset(Standard_False);
 543     IMeshData::MapOfIntegerInteger aLoopEdges(10, aAllocator);
 544
 545     Standard_Integer aVertexIdx = theVertexIndexes( anIndex );
 546     const BRepMesh_Vertex& aVertex = GetVertex( aVertexIdx );
 547
 548     // Iterator in the list of indexes of circles containing the node
 549     IMeshData::ListOfInteger& aCirclesList = myCircles.Select( aVertex.Coord() );
 550
 551     Standard_Integer onEgdeId = 0, aTriangleId = 0;
 552     IMeshData::ListOfInteger::Iterator aCircleIt( aCirclesList );
 553     for ( ; aCircleIt.More(); aCircleIt.Next() )
 554     {
 555       // To add a node in the mesh it is necessary to check conditions:
 556       // - the node should be within the boundaries of the mesh and so in an existing triangle
 557       // - all adjacent triangles should belong to a component connected with this triangle
 558       if ( Contains( aCircleIt.Value(), aVertex, aSqTol, onEgdeId ) )
 559       {
 560         if (onEgdeId != 0 && GetEdge(onEgdeId).Movability() != BRepMesh_Free)
 561         {
 562           // We can skip free vertex too close to the frontier edge.
 563           if (aVertex.Movability() == BRepMesh_Free)
 564             continue;
 565
 566           // However, we should add vertex that have neighboring frontier edges.
 567         }
 568
 569         // Remove triangle even if it contains frontier edge in order
 570         // to prevent appearance of incorrect configurations like free
 571         // edge glued with frontier #26407
 572         aTriangleId = aCircleIt.Value();
 573         aCirclesList.Remove( aCircleIt );
 574         break;
 575       }
 576     }
 577
 578     if ( aTriangleId > 0 )
 579     {
 580       deleteTriangle( aTriangleId, aLoopEdges );
 581
 582       isModify = Standard_True;
 583       while ( isModify && !aCirclesList.IsEmpty() )
 584       {
 585         isModify = Standard_False;
 586         IMeshData::ListOfInteger::Iterator aCircleIt1( aCirclesList );
 587         for ( ; aCircleIt1.More(); aCircleIt1.Next() )
 588         {
 589           const BRepMesh_Triangle& aElement = GetTriangle(aCircleIt1.Value());
 590           const Standard_Integer(&e)[3] = aElement.myEdges;
 591
 592           if ( aLoopEdges.IsBound( e[0] ) ||
 593                aLoopEdges.IsBound( e[1] ) ||
 594                aLoopEdges.IsBound( e[2] ) )
 595           {
 596             isModify = Standard_True;
 597             deleteTriangle( aCircleIt1.Value(), aLoopEdges );
 598             aCirclesList.Remove( aCircleIt1 );
 599             break;
 600           }
 601         }
 602       }
 603
 604       // Creation of triangles with the current node and free edges
 605       // and removal of these edges from the list of free edges
 606       createTriangles( aVertexIdx, aLoopEdges );
 607     }
 608   }
 609
 610   ProcessConstraints();
 611 }
 612
 613 //=======================================================================
 614 //function : insertInternalEdges
 615 //purpose  :
 616 //=======================================================================
 617 void BRepMesh_Delaun::insertInternalEdges()
 618 {
 619   Handle(IMeshData::MapOfInteger) anInternalEdges = InternalEdges();
 620
 621   // Destruction of triancles intersecting internal edges
 622   // and their replacement by makeshift triangles
 623   IMeshData::IteratorOfMapOfInteger anInernalEdgesIt( *anInternalEdges );
 624   for ( ; anInernalEdgesIt.More(); anInernalEdgesIt.Next() )
 625   {
 626     const Standard_Integer aLinkIndex = anInernalEdgesIt.Key();
 627     const BRepMesh_PairOfIndex& aPair = myMeshData->ElementsConnectedTo(aLinkIndex);
 628
 629     // Check both sides of link for adjusted triangle.
 630     Standard_Boolean isGo[2] = { Standard_True, Standard_True };
 631     for (Standard_Integer aTriangleIt = 1; aTriangleIt <= aPair.Extent(); ++aTriangleIt)
 632     {
 633       const BRepMesh_Triangle& aElement = GetTriangle(aPair.Index(aTriangleIt));
 634       const Standard_Integer(&e)[3] = aElement.myEdges;
 635       const Standard_Boolean(&o)[3] = aElement.myOrientations;
 636
 637       for (Standard_Integer i = 0; i < 3; ++i)
 638       {
 639         if (e[i] == aLinkIndex)
 640         {
 641           isGo[o[i] ? 0 : 1] = Standard_False;
 642           break;
 643         }
 644       }
 645     }
 646
 647     if (isGo[0])
 648     {
 649       meshLeftPolygonOf(aLinkIndex, Standard_True);
 650     }
 651
 652     if (isGo[1])
 653     {
 654       meshLeftPolygonOf(aLinkIndex, Standard_False);
 655     }
 656   }
 657 }
 658
 659 //=======================================================================
 660 //function : isBoundToFrontier
 661 //purpose  :
 662 //=======================================================================
 663 Standard_Boolean BRepMesh_Delaun::isBoundToFrontier(
 664   const Standard_Integer theRefNodeId,
 665   const Standard_Integer theRefLinkId)
 666 {
 667   std::stack<Standard_Integer> aLinkStack;
 668   TColStd_PackedMapOfInteger   aVisitedLinks;
 669
 670   aLinkStack.push (theRefLinkId);
 671   while (!aLinkStack.empty ())
 672   {
 673     const Standard_Integer aCurrentLinkId = aLinkStack.top ();
 674     aLinkStack.pop ();
 675
 676     const BRepMesh_PairOfIndex& aPair = myMeshData->ElementsConnectedTo (aCurrentLinkId);
 677     if (aPair.IsEmpty ())
 678       return Standard_False;
 679
 680     const Standard_Integer aNbElements = aPair.Extent ();
 681     for (Standard_Integer anElemIt = 1; anElemIt <= aNbElements; ++anElemIt)
 682     {
 683       const Standard_Integer aTriId = aPair.Index (anElemIt);
 684       if (aTriId < 0)
 685         continue;
 686
 687       const BRepMesh_Triangle& aElement = GetTriangle (aTriId);
 688       const Standard_Integer (&anEdges)[3] = aElement.myEdges;
 689
 690       for (Standard_Integer anEdgeIt = 0; anEdgeIt < 3; ++anEdgeIt)
 691       {
 692         const Standard_Integer anEdgeId = anEdges[anEdgeIt];
 693         if (anEdgeId == aCurrentLinkId)
 694           continue;
 695
 696         const BRepMesh_Edge& anEdge = GetEdge (anEdgeId);
 697         if (anEdge.FirstNode () != theRefNodeId &&
 698             anEdge.LastNode  () != theRefNodeId)
 699         {
 700           continue;
 701         }
 702
 703         if (anEdge.Movability () != BRepMesh_Free)
 704           return Standard_True;
 705
 706         if (aVisitedLinks.Add (anEdgeId))
 707         {
 708           aLinkStack.push (anEdgeId);
 709         }
 710       }
 711     }
 712   }
 713
 714   return Standard_False;
 715 }
 716
 717 //=======================================================================
 718 //function : cleanupMesh
 719 //purpose  : Cleanup mesh from the free triangles
 720 //=======================================================================
 721 void BRepMesh_Delaun::cleanupMesh()
 722 {
 723   Handle(NCollection_IncAllocator) aAllocator =
 724     new NCollection_IncAllocator(IMeshData::MEMORY_BLOCK_SIZE_HUGE);
 725
 726   for(;;)
 727   {
 728     aAllocator->Reset(Standard_False);
 729     IMeshData::MapOfIntegerInteger aLoopEdges(10, aAllocator);
 730     IMeshData::MapOfInteger aDelTriangles;
 731
 732     Handle(IMeshData::MapOfInteger) aFreeEdges = FreeEdges();
 733     IMeshData::IteratorOfMapOfInteger aFreeEdgesIt( *aFreeEdges );
 734     for ( ; aFreeEdgesIt.More(); aFreeEdgesIt.Next() )
 735     {
 736       const Standard_Integer& aFreeEdgeId = aFreeEdgesIt.Key();
 737       const BRepMesh_Edge&    anEdge      = GetEdge( aFreeEdgeId );
 738       if ( anEdge.Movability() == BRepMesh_Frontier )
 739         continue;
 740
 741       const BRepMesh_PairOfIndex& aPair =
 742         myMeshData->ElementsConnectedTo( aFreeEdgeId );
 743       if ( aPair.IsEmpty() )
 744       {
 745         aLoopEdges.Bind( aFreeEdgeId, Standard_True );
 746         continue;
 747       }
 748
 749       Standard_Integer aTriId = aPair.FirstIndex();
 750
 751       // Check that the connected triangle is not surrounded by another triangles
 752       const BRepMesh_Triangle& aElement = GetTriangle(aTriId);
 753       const Standard_Integer(&anEdges)[3] = aElement.myEdges;
 754
 755       Standard_Boolean isCanNotBeRemoved = Standard_True;
 756       for ( Standard_Integer aCurEdgeIdx = 0; aCurEdgeIdx < 3; ++aCurEdgeIdx )
 757       {
 758         if ( anEdges[aCurEdgeIdx] != aFreeEdgeId )
 759           continue;
 760
 761         for ( Standard_Integer anOtherEdgeIt = 1; anOtherEdgeIt <= 2 && isCanNotBeRemoved; ++anOtherEdgeIt )
 762         {
 763           Standard_Integer anOtherEdgeId = ( aCurEdgeIdx + anOtherEdgeIt ) % 3;
 764           const BRepMesh_PairOfIndex& anOtherEdgePair =
 765             myMeshData->ElementsConnectedTo( anEdges[anOtherEdgeId] );
 766
 767           if ( anOtherEdgePair.Extent() < 2 )
 768           {
 769             isCanNotBeRemoved = Standard_False;
 770           }
 771           else
 772           {
 773             for (int aTriIdx = 1; aTriIdx <= anOtherEdgePair.Extent () && isCanNotBeRemoved; ++aTriIdx)
 774             {
 775               if (anOtherEdgePair.Index (aTriIdx) == aTriId)
 776                 continue;
 777
 778               Standard_Integer v[3];
 779               const BRepMesh_Triangle& aCurTriangle = GetTriangle (anOtherEdgePair.Index (aTriIdx));
 780               myMeshData->ElementNodes (aCurTriangle, v);
 781               for (int aNodeIdx = 0; aNodeIdx < 3 && isCanNotBeRemoved; ++aNodeIdx)
 782               {
 783                 if (v[aNodeIdx] == mySupVert[0] ||
 784                     v[aNodeIdx] == mySupVert[1] ||
 785                     v[aNodeIdx] == mySupVert[2])
 786                 {
 787                   isCanNotBeRemoved = Standard_False;
 788                 }
 789               }
 790             }
 791           }
 792         }
 793
 794         break;
 795       }
 796
 797       if ( isCanNotBeRemoved )
 798         continue;
 799
 800       Standard_Boolean isConnected[2] = { Standard_False, Standard_False };
 801       for ( Standard_Integer aLinkNodeIt = 0; aLinkNodeIt < 2; ++aLinkNodeIt )
 802       {
 803         isConnected[aLinkNodeIt] = isBoundToFrontier( ( aLinkNodeIt == 0 ) ?
 804           anEdge.FirstNode() : anEdge.LastNode(), aFreeEdgeId);
 805       }
 806
 807       if ( !isConnected[0] || !isConnected[1] )
 808         aDelTriangles.Add( aTriId );
 809     }
 810
 811     // Destruction of triangles :
 812     Standard_Integer aDeletedTrianglesNb = 0;
 813     IMeshData::IteratorOfMapOfInteger aDelTrianglesIt( aDelTriangles );
 814     for ( ; aDelTrianglesIt.More(); aDelTrianglesIt.Next() )
 815     {
 816       deleteTriangle( aDelTrianglesIt.Key(), aLoopEdges );
 817       aDeletedTrianglesNb++;
 818     }
 819
 820     // Destruction of remaining hanging edges
 821     IMeshData::MapOfIntegerInteger::Iterator aLoopEdgesIt( aLoopEdges );
 822     for ( ; aLoopEdgesIt.More(); aLoopEdgesIt.Next() )
 823     {
 824       if ( myMeshData->ElementsConnectedTo( aLoopEdgesIt.Key() ).IsEmpty() )
 825         myMeshData->RemoveLink( aLoopEdgesIt.Key() );
 826     }
 827
 828     if ( aDeletedTrianglesNb == 0 )
 829       break;
 830   }
 831 }
 832
 833 //=======================================================================
 834 //function : frontierAdjust
 835 //purpose  : Adjust the mesh on the frontier
 836 //=======================================================================
 837 void BRepMesh_Delaun::frontierAdjust()
 838 {
 839   Handle(IMeshData::MapOfInteger)  aFrontier = Frontier();
 840
 841   Handle(NCollection_IncAllocator) aAllocator =
 842     new NCollection_IncAllocator(IMeshData::MEMORY_BLOCK_SIZE_HUGE);
 843
 844   IMeshData::VectorOfInteger      aFailedFrontiers(256, aAllocator);
 845   IMeshData::MapOfIntegerInteger  aLoopEdges(10, aAllocator);
 846   Handle(IMeshData::MapOfInteger) aIntFrontierEdges = new IMeshData::MapOfInteger;
 847
 848   for ( Standard_Integer aPass = 1; aPass <= 2; ++aPass )
 849   {
 850     // 1 pass): find external triangles on boundary edges;
 851     // 2 pass): find external triangles on boundary edges appeared
 852     //          during triangles replacement.
 853
 854     IMeshData::IteratorOfMapOfInteger aFrontierIt( *aFrontier );
 855     for ( ; aFrontierIt.More(); aFrontierIt.Next() )
 856     {
 857       Standard_Integer aFrontierId = aFrontierIt.Key();
 858       const BRepMesh_PairOfIndex& aPair = myMeshData->ElementsConnectedTo( aFrontierId );
 859       Standard_Integer aNbElem = aPair.Extent();
 860       for( Standard_Integer aElemIt = 1; aElemIt <= aNbElem; ++aElemIt )
 861       {
 862         const Standard_Integer aPriorElemId = aPair.Index( aElemIt );
 863         if( aPriorElemId < 0 )
 864           continue;
 865
 866         const BRepMesh_Triangle& aElement = GetTriangle(aPriorElemId);
 867         const Standard_Integer(&e)[3] = aElement.myEdges;
 868         const Standard_Boolean(&o)[3] = aElement.myOrientations;
 869
 870         Standard_Boolean isTriangleFound = Standard_False;
 871         for ( Standard_Integer n = 0; n < 3; ++n )
 872         {
 873           if ( aFrontierId == e[n] && !o[n] )
 874           {
 875             // Destruction  of external triangles on boundary edges
 876             isTriangleFound = Standard_True;
 877             deleteTriangle( aPriorElemId, aLoopEdges );
 878             break;
 879           }
 880         }
 881
 882         if ( isTriangleFound )
 883           break;
 884       }
 885     }
 886
 887     // destrucrion of remaining hanging edges :
 888     IMeshData::MapOfIntegerInteger::Iterator aLoopEdgesIt( aLoopEdges );
 889     for ( ; aLoopEdgesIt.More(); aLoopEdgesIt.Next() )
 890     {
 891       Standard_Integer aLoopEdgeId = aLoopEdgesIt.Key();
 892       if (myMeshData->ElementsConnectedTo( aLoopEdgeId ).IsEmpty() )
 893         myMeshData->RemoveLink( aLoopEdgeId );
 894     }
 895
 896     // destruction of triangles crossing the boundary edges and
 897     // their replacement by makeshift triangles
 898     for ( aFrontierIt.Reset(); aFrontierIt.More(); aFrontierIt.Next() )
 899     {
 900       Standard_Integer aFrontierId = aFrontierIt.Key();
 901       if ( !myMeshData->ElementsConnectedTo( aFrontierId ).IsEmpty() )
 902         continue;
 903
 904       Standard_Boolean isSuccess =
 905         meshLeftPolygonOf( aFrontierId, Standard_True, aIntFrontierEdges );
 906
 907       if ( aPass == 2 && !isSuccess )
 908         aFailedFrontiers.Append( aFrontierId );
 909     }
 910   }
 911
 912   cleanupMesh();
 913
 914   // When the mesh has been cleaned up, try to process frontier edges
 915   // once again to fill the possible gaps that might be occured in case of "saw" -
 916   // situation when frontier edge has a triangle at a right side, but its free
 917   // links cross another frontieres  and meshLeftPolygonOf itself can't collect
 918   // a closed polygon.
 919   IMeshData::VectorOfInteger::Iterator aFailedFrontiersIt( aFailedFrontiers );
 920   for ( ; aFailedFrontiersIt.More(); aFailedFrontiersIt.Next() )
 921   {
 922     Standard_Integer aFrontierId = aFailedFrontiersIt.Value();
 923     if ( !myMeshData->ElementsConnectedTo( aFrontierId ).IsEmpty() )
 924       continue;
 925
 926     meshLeftPolygonOf( aFrontierId, Standard_True, aIntFrontierEdges );
 927   }
 928 }
 929
 930 //=======================================================================
 931 //function : fillBndBox
 932 //purpose  : Add boundig box for edge defined by start & end point to
 933 //           the given vector of bounding boxes for triangulation edges
 934 //=======================================================================
 935 void BRepMesh_Delaun::fillBndBox(IMeshData::SequenceOfBndB2d& theBoxes,
 936                                  const BRepMesh_Vertex&       theV1,
 937                                  const BRepMesh_Vertex&       theV2)
 938 {
 939   Bnd_B2d aBox;
 940   UpdateBndBox(theV1.Coord(), theV2.Coord(), aBox);
 941   theBoxes.Append( aBox );
 942 }
 943
 944 //=======================================================================
 945 //function : meshLeftPolygonOf
 946 //purpose  : Collect the polygon at the left of the given edge (material side)
 947 //=======================================================================
 948 Standard_Boolean BRepMesh_Delaun::meshLeftPolygonOf(
 949   const Standard_Integer          theStartEdgeId,
 950   const Standard_Boolean          isForward,
 951   Handle(IMeshData::MapOfInteger) theSkipped)
 952 {
 953   if ( !theSkipped.IsNull() && theSkipped->Contains( theStartEdgeId ) )
 954     return Standard_True;
 955
 956   const BRepMesh_Edge& aRefEdge = GetEdge( theStartEdgeId );
 957
 958   IMeshData::SequenceOfInteger aPolygon;
 959   Standard_Integer aStartNode, aPivotNode;
 960   if ( isForward )
 961   {
 962     aPolygon.Append( theStartEdgeId );
 963     aStartNode = aRefEdge.FirstNode();
 964     aPivotNode = aRefEdge.LastNode();
 965   }
 966   else
 967   {
 968     aPolygon.Append( -theStartEdgeId );
 969     aStartNode = aRefEdge.LastNode();
 970     aPivotNode = aRefEdge.FirstNode();
 971   }
 972
 973
 974   const BRepMesh_Vertex& aStartEdgeVertexS = GetVertex( aStartNode );
 975   BRepMesh_Vertex        aPivotVertex      = GetVertex( aPivotNode );
 976
 977   gp_Vec2d aRefLinkDir( aPivotVertex.Coord() -
 978     aStartEdgeVertexS.Coord() );
 979
 980   if ( aRefLinkDir.SquareMagnitude() < Precision2 )
 981     return Standard_True;
 982
 983   // Auxilary structures.
 984   // Bounding boxes of polygon links to be used for preliminary
 985   // analysis of intersections
 986   IMeshData::SequenceOfBndB2d aBoxes;
 987   fillBndBox( aBoxes, aStartEdgeVertexS, aPivotVertex );
 988
 989   // Hanging ends
 990   IMeshData::MapOfInteger aDeadLinks;
 991
 992   // Links are temporarily excluded from consideration
 993   IMeshData::MapOfInteger aLeprousLinks;
 994   aLeprousLinks.Add( theStartEdgeId );
 995
 996   Standard_Boolean isSkipLeprous = Standard_True;
 997   Standard_Integer aFirstNode    = aStartNode;
 998   while ( aPivotNode != aFirstNode )
 999   {
1000     Bnd_B2d          aNextLinkBndBox;
1001     gp_Vec2d         aNextLinkDir;
1002     Standard_Integer aNextPivotNode = 0;
1003
1004     Standard_Integer aNextLinkId = findNextPolygonLink(
1005       aFirstNode,
1006       aPivotNode,     aPivotVertex,  aRefLinkDir,
1007       aBoxes,         aPolygon,      theSkipped,
1008       isSkipLeprous,  aLeprousLinks, aDeadLinks,
1009       aNextPivotNode, aNextLinkDir,  aNextLinkBndBox );
1010
1011     if ( aNextLinkId != 0 )
1012     {
1013       aStartNode   = aPivotNode;
1014       aRefLinkDir  = aNextLinkDir;
1015
1016       aPivotNode   = aNextPivotNode;
1017       aPivotVertex = GetVertex( aNextPivotNode );
1018
1019       aBoxes.Append  ( aNextLinkBndBox );
1020       aPolygon.Append( aNextLinkId );
1021
1022       isSkipLeprous = Standard_True;
1023     }
1024     else
1025     {
1026       // Nothing to do
1027       if ( aPolygon.Length() == 1 )
1028         return Standard_False;
1029
1030       // Return to the previous point
1031       Standard_Integer aDeadLinkId = Abs( aPolygon.Last() );
1032       aDeadLinks.Add      ( aDeadLinkId );
1033
1034       aLeprousLinks.Remove( aDeadLinkId );
1035       aPolygon.Remove     ( aPolygon.Length() );
1036       aBoxes.Remove       ( aBoxes.Length() );
1037
1038       Standard_Integer aPrevLinkInfo = aPolygon.Last();
1039       const BRepMesh_Edge& aPrevLink = GetEdge( Abs( aPrevLinkInfo ) );
1040
1041       if( aPrevLinkInfo > 0 )
1042       {
1043         aStartNode = aPrevLink.FirstNode();
1044         aPivotNode = aPrevLink.LastNode();
1045       }
1046       else
1047       {
1048         aStartNode = aPrevLink.LastNode();
1049         aPivotNode = aPrevLink.FirstNode();
1050       }
1051
1052       aPivotVertex = GetVertex( aPivotNode );
1053       aRefLinkDir =
1054         aPivotVertex.Coord() - GetVertex( aStartNode ).Coord();
1055
1056       isSkipLeprous = Standard_False;
1057     }
1058   }
1059
1060   if ( aPolygon.Length() < 3 )
1061     return Standard_False;
1062
1063   cleanupPolygon( aPolygon, aBoxes );
1064   meshPolygon   ( aPolygon, aBoxes, theSkipped );
1065
1066   return Standard_True;
1067 }
1068
1069 //=======================================================================
1070 //function : findNextPolygonLink
1071 //purpose  : Find next link starting from the given node and has maximum
1072 //           angle respect the given reference link.
1073 //           Each time the next link is found other neighbor links at the
1074 //           pivot node are marked as leprous and will be excluded from
1075 //           consideration next time until a hanging end is occured.
1076 //=======================================================================
1077 Standard_Integer BRepMesh_Delaun::findNextPolygonLink(
1078   const Standard_Integer&               theFirstNode,
1079   const Standard_Integer&               thePivotNode,
1080   const BRepMesh_Vertex&                thePivotVertex,
1081   const gp_Vec2d&                       theRefLinkDir,
1082   const IMeshData::SequenceOfBndB2d&    theBoxes,
1083   const IMeshData::SequenceOfInteger&   thePolygon,
1084   const Handle(IMeshData::MapOfInteger) theSkipped,
1085   const Standard_Boolean&               isSkipLeprous,
1086   IMeshData::MapOfInteger&              theLeprousLinks,
1087   IMeshData::MapOfInteger&              theDeadLinks,
1088   Standard_Integer&                     theNextPivotNode,
1089   gp_Vec2d&                             theNextLinkDir,
1090   Bnd_B2d&                              theNextLinkBndBox )
1091 {
1092   // Find the next link having the greatest angle
1093   // respect to a direction of a reference one
1094   Standard_Real aMaxAngle = RealFirst();
1095
1096   Standard_Integer aNextLinkId = 0;
1097   IMeshData::ListOfInteger::Iterator aLinkIt( myMeshData->LinksConnectedTo( thePivotNode ) );
1098   for ( ; aLinkIt.More(); aLinkIt.Next() )
1099   {
1100     const Standard_Integer& aNeighbourLinkInfo = aLinkIt.Value();
1101     Standard_Integer        aNeighbourLinkId   = Abs( aNeighbourLinkInfo );
1102
1103     if ( theDeadLinks.Contains( aNeighbourLinkId ) ||
1104        ( !theSkipped.IsNull() && theSkipped->Contains( aNeighbourLinkId ) ) )
1105     {
1106       continue;
1107     }
1108
1109     Standard_Boolean isLeprous = theLeprousLinks.Contains( aNeighbourLinkId );
1110     if ( isSkipLeprous && isLeprous )
1111       continue;
1112
1113     const BRepMesh_Edge& aNeighbourLink = GetEdge( aNeighbourLinkId );
1114
1115     // Determine whether the link belongs to the mesh
1116     if ( aNeighbourLink.Movability() == BRepMesh_Free &&
1117          myMeshData->ElementsConnectedTo( aNeighbourLinkInfo ).IsEmpty() )
1118     {
1119       theDeadLinks.Add( aNeighbourLinkId );
1120       continue;
1121     }
1122
1123     Standard_Integer anOtherNode = aNeighbourLink.FirstNode();
1124     if ( anOtherNode == thePivotNode )
1125       anOtherNode = aNeighbourLink.LastNode();
1126
1127     gp_Vec2d aCurLinkDir( GetVertex( anOtherNode ).Coord() -
1128       thePivotVertex.Coord() );
1129
1130     if ( aCurLinkDir.SquareMagnitude() < Precision2 )
1131     {
1132       theDeadLinks.Add( aNeighbourLinkId );
1133       continue;
1134     }
1135
1136     if ( !isLeprous )
1137       theLeprousLinks.Add( aNeighbourLinkId );
1138
1139     Standard_Real    anAngle    = theRefLinkDir.Angle( aCurLinkDir );
1140     Standard_Boolean isFrontier =
1141       ( aNeighbourLink.Movability() == BRepMesh_Frontier );
1142
1143     Standard_Boolean isCheckPointOnEdge = Standard_True;
1144     if ( isFrontier )
1145     {
1146       if ( Abs( Abs(anAngle) - M_PI ) < Precision::Angular() )
1147       {
1148         // Glued constrains - don't check intersection
1149         isCheckPointOnEdge = Standard_False;
1150         anAngle            = Abs( anAngle );
1151       }
1152     }
1153
1154     if (anAngle <= aMaxAngle)
1155       continue;
1156
1157     Standard_Boolean isCheckEndPoints = ( anOtherNode != theFirstNode );
1158
1159     Bnd_B2d aBox;
1160     Standard_Boolean isNotIntersect =
1161       checkIntersection( aNeighbourLink, thePolygon, theBoxes,
1162       isCheckEndPoints, isCheckPointOnEdge, Standard_True, aBox );
1163
1164     if( isNotIntersect )
1165     {
1166       aMaxAngle         = anAngle;
1167
1168       theNextLinkDir    = aCurLinkDir;
1169       theNextPivotNode  = anOtherNode;
1170       theNextLinkBndBox = aBox;
1171
1172       aNextLinkId       = ( aNeighbourLink.FirstNode() == thePivotNode ) ?
1173         aNeighbourLinkId : -aNeighbourLinkId;
1174     }
1175   }
1176
1177   return aNextLinkId;
1178 }
1179
1180 //=======================================================================
1181 //function : checkIntersection
1182 //purpose  : Check is the given link intersects the polygon boundaries.
1183 //           Returns bounding box for the given link trough the
1184 //           <theLinkBndBox> parameter.
1185 //=======================================================================
1186 Standard_Boolean BRepMesh_Delaun::checkIntersection(
1187   const BRepMesh_Edge&                theLink,
1188   const IMeshData::SequenceOfInteger& thePolygon,
1189   const IMeshData::SequenceOfBndB2d&  thePolyBoxes,
1190   const Standard_Boolean              isConsiderEndPointTouch,
1191   const Standard_Boolean              isConsiderPointOnEdge,
1192   const Standard_Boolean              isSkipLastEdge,
1193   Bnd_B2d&                            theLinkBndBox ) const
1194 {
1195   UpdateBndBox(GetVertex(theLink.FirstNode()).Coord(),
1196     GetVertex(theLink.LastNode()).Coord(), theLinkBndBox);
1197
1198   Standard_Integer aPolyLen = thePolygon.Length();
1199   // Don't check intersection with the last link
1200   if ( isSkipLastEdge )
1201     --aPolyLen;
1202
1203   Standard_Boolean isFrontier =
1204     ( theLink.Movability() == BRepMesh_Frontier );
1205
1206   for ( Standard_Integer aPolyIt = 1; aPolyIt <= aPolyLen; ++aPolyIt )
1207   {
1208     if ( !theLinkBndBox.IsOut( thePolyBoxes.Value( aPolyIt ) ) )
1209     {
1210       // intersection is possible...
1211       Standard_Integer aPolyLinkId   = Abs( thePolygon( aPolyIt ) );
1212       const BRepMesh_Edge& aPolyLink = GetEdge( aPolyLinkId );
1213
1214       // skip intersections between frontier edges
1215       if ( aPolyLink.Movability() == BRepMesh_Frontier && isFrontier )
1216         continue;
1217
1218       gp_Pnt2d anIntPnt;
1219       BRepMesh_GeomTool::IntFlag aIntFlag = intSegSeg( theLink, aPolyLink,
1220         isConsiderEndPointTouch, isConsiderPointOnEdge, anIntPnt );
1221
1222       if ( aIntFlag != BRepMesh_GeomTool::NoIntersection )
1223         return Standard_False;
1224     }
1225   }
1226
1227   // Ok, no intersection
1228   return Standard_True;
1229 }
1230
1231 //=======================================================================
1232 //function : addTriangle
1233 //purpose  : Add a triangle based on the given oriented edges into mesh
1234 //=======================================================================
1235 void BRepMesh_Delaun::addTriangle( const Standard_Integer (&theEdgesId)[3],
1236                                    const Standard_Boolean (&theEdgesOri)[3],
1237                                    const Standard_Integer (&theNodesId)[3] )
1238 {
1239   Standard_Integer aNewTriangleId =
1240     myMeshData->AddElement(BRepMesh_Triangle(theEdgesId,
1241       theEdgesOri, BRepMesh_Free));
1242
1243   Standard_Boolean isAdded = myCircles.Bind(
1244     aNewTriangleId,
1245     GetVertex( theNodesId[0] ).Coord(),
1246     GetVertex( theNodesId[1] ).Coord(),
1247     GetVertex( theNodesId[2] ).Coord() );
1248
1249   if ( !isAdded )
1250     myMeshData->RemoveElement( aNewTriangleId );
1251 }
1252
1253 //=======================================================================
1254 //function : cleanupPolygon
1255 //purpose  : Remove internal triangles from the given polygon
1256 //=======================================================================
1257 void BRepMesh_Delaun::cleanupPolygon(const IMeshData::SequenceOfInteger& thePolygon,
1258                                      const IMeshData::SequenceOfBndB2d&  thePolyBoxes )
1259 {
1260   Standard_Integer aPolyLen = thePolygon.Length();
1261   if ( aPolyLen < 3 )
1262     return;
1263
1264   Handle(NCollection_IncAllocator) aAllocator =
1265     new NCollection_IncAllocator(IMeshData::MEMORY_BLOCK_SIZE_HUGE);
1266
1267   IMeshData::MapOfIntegerInteger aLoopEdges(10, aAllocator);
1268   IMeshData::MapOfInteger    anIgnoredEdges;
1269   IMeshData::MapOfInteger    aPolyVerticesFindMap;
1270   IMeshData::VectorOfInteger aPolyVertices(256, aAllocator);
1271   // Collect boundary vertices of the polygon
1272   for ( Standard_Integer aPolyIt = 1; aPolyIt <= aPolyLen; ++aPolyIt )
1273   {
1274     Standard_Integer aPolyEdgeInfo = thePolygon( aPolyIt );
1275     Standard_Integer aPolyEdgeId   = Abs( aPolyEdgeInfo );
1276     anIgnoredEdges.Add( aPolyEdgeId );
1277
1278     Standard_Boolean isForward = ( aPolyEdgeInfo > 0 );
1279     const BRepMesh_PairOfIndex& aPair =
1280       myMeshData->ElementsConnectedTo( aPolyEdgeId );
1281
1282     Standard_Integer anElemIt = 1;
1283     for ( ; anElemIt <= aPair.Extent(); ++anElemIt )
1284     {
1285       Standard_Integer anElemId = aPair.Index( anElemIt );
1286       if ( anElemId < 0 )
1287         continue;
1288
1289       const BRepMesh_Triangle& aElement = GetTriangle(anElemId);
1290       const Standard_Integer(&anEdges)[3] = aElement.myEdges;
1291       const Standard_Boolean(&anEdgesOri)[3] = aElement.myOrientations;
1292
1293       Standard_Integer isTriangleFound = Standard_False;
1294       for ( Standard_Integer anEdgeIt = 0; anEdgeIt < 3; ++anEdgeIt )
1295       {
1296         if ( anEdges[anEdgeIt]    == aPolyEdgeId &&
1297              anEdgesOri[anEdgeIt] == isForward )
1298         {
1299           isTriangleFound = Standard_True;
1300           deleteTriangle( anElemId, aLoopEdges );
1301           break;
1302         }
1303       }
1304
1305       if ( isTriangleFound )
1306         break;
1307     }
1308
1309     // Skip a neighbor link to extract unique vertices each time
1310     if ( aPolyIt % 2 )
1311     {
1312       const BRepMesh_Edge& aPolyEdge = GetEdge( aPolyEdgeId );
1313       Standard_Integer aFirstVertex  = aPolyEdge.FirstNode();
1314       Standard_Integer aLastVertex   = aPolyEdge.LastNode();
1315
1316       aPolyVerticesFindMap.Add( aFirstVertex );
1317       aPolyVerticesFindMap.Add( aLastVertex );
1318
1319       if ( aPolyEdgeInfo > 0 )
1320       {
1321         aPolyVertices.Append( aFirstVertex );
1322         aPolyVertices.Append( aLastVertex );
1323       }
1324       else
1325       {
1326         aPolyVertices.Append( aLastVertex );
1327         aPolyVertices.Append( aFirstVertex );
1328       }
1329     }
1330   }
1331
1332   // Make closed sequence
1333   if ( aPolyVertices.First() != aPolyVertices.Last() )
1334     aPolyVertices.Append( aPolyVertices.First() );
1335
1336   IMeshData::MapOfInteger aSurvivedLinks( anIgnoredEdges );
1337
1338   Standard_Integer aPolyVertIt          = 0;
1339   Standard_Integer anUniqueVerticesNum  = aPolyVertices.Length() - 1;
1340   for ( ; aPolyVertIt < anUniqueVerticesNum; ++aPolyVertIt )
1341   {
1342     killTrianglesAroundVertex( aPolyVertices( aPolyVertIt ),
1343       aPolyVertices, aPolyVerticesFindMap, thePolygon,
1344       thePolyBoxes, aSurvivedLinks, aLoopEdges );
1345   }
1346
1347   IMeshData::MapOfIntegerInteger::Iterator aLoopEdgesIt( aLoopEdges );
1348   for ( ; aLoopEdgesIt.More(); aLoopEdgesIt.Next() )
1349   {
1350     const Standard_Integer& aLoopEdgeId = aLoopEdgesIt.Key();
1351     if ( anIgnoredEdges.Contains( aLoopEdgeId ) )
1352       continue;
1353
1354     if ( myMeshData->ElementsConnectedTo( aLoopEdgeId ).IsEmpty() )
1355       myMeshData->RemoveLink( aLoopEdgesIt.Key() );
1356   }
1357 }
1358
1359 //=======================================================================
1360 //function : killTrianglesAroundVertex
1361 //purpose  : Remove all triangles and edges that are placed
1362 //           inside the polygon or crossed it.
1363 //=======================================================================
1364 void BRepMesh_Delaun::killTrianglesAroundVertex(
1365   const Standard_Integer              theZombieNodeId,
1366   const IMeshData::VectorOfInteger&   thePolyVertices,
1367   const IMeshData::MapOfInteger&      thePolyVerticesFindMap,
1368   const IMeshData::SequenceOfInteger& thePolygon,
1369   const IMeshData::SequenceOfBndB2d&  thePolyBoxes,
1370   IMeshData::MapOfInteger&            theSurvivedLinks,
1371   IMeshData::MapOfIntegerInteger&     theLoopEdges )
1372 {
1373   IMeshData::ListOfInteger::Iterator aNeighborsIt =
1374     myMeshData->LinksConnectedTo( theZombieNodeId );
1375
1376   // Try to infect neighbor nodes
1377   IMeshData::VectorOfInteger aVictimNodes;
1378   for ( ; aNeighborsIt.More(); aNeighborsIt.Next() )
1379   {
1380     const Standard_Integer& aNeighborLinkId = aNeighborsIt.Value();
1381     if ( theSurvivedLinks.Contains( aNeighborLinkId ) )
1382       continue;
1383
1384     const BRepMesh_Edge& aNeighborLink = GetEdge( aNeighborLinkId );
1385     if ( aNeighborLink.Movability() == BRepMesh_Frontier )
1386     {
1387       // Though, if it lies onto the polygon boundary -
1388       // take its triangles
1389       Bnd_B2d aBox;
1390       Standard_Boolean isNotIntersect =
1391         checkIntersection( aNeighborLink, thePolygon,
1392         thePolyBoxes, Standard_False, Standard_True,
1393         Standard_False, aBox );
1394
1395       if ( isNotIntersect )
1396       {
1397         // Don't touch it
1398         theSurvivedLinks.Add( aNeighborLinkId );
1399         continue;
1400       }
1401     }
1402     else
1403     {
1404       Standard_Integer anOtherNode = aNeighborLink.FirstNode();
1405       if ( anOtherNode == theZombieNodeId )
1406         anOtherNode = aNeighborLink.LastNode();
1407
1408       // Possible sacrifice
1409       if ( !thePolyVerticesFindMap.Contains( anOtherNode ) )
1410       {
1411         if ( isVertexInsidePolygon( anOtherNode, thePolyVertices ) )
1412         {
1413           // Got you!
1414           aVictimNodes.Append( anOtherNode );
1415         }
1416         else
1417         {
1418           // Lucky. But it may intersect the polygon boundary.
1419           // Let's check it!
1420           killTrianglesOnIntersectingLinks( aNeighborLinkId,
1421             aNeighborLink, anOtherNode, thePolygon,
1422             thePolyBoxes, theSurvivedLinks, theLoopEdges );
1423
1424           continue;
1425         }
1426       }
1427     }
1428
1429     // Add link to the survivers to avoid cycling
1430     theSurvivedLinks.Add( aNeighborLinkId );
1431     killLinkTriangles( aNeighborLinkId, theLoopEdges );
1432   }
1433
1434   // Go and do your job!
1435   IMeshData::VectorOfInteger::Iterator aVictimIt( aVictimNodes );
1436   for ( ; aVictimIt.More(); aVictimIt.Next() )
1437   {
1438     killTrianglesAroundVertex( aVictimIt.Value(), thePolyVertices,
1439       thePolyVerticesFindMap, thePolygon, thePolyBoxes,
1440       theSurvivedLinks, theLoopEdges );
1441   }
1442 }
1443
1444 //=======================================================================
1445 //function : isVertexInsidePolygon
1446 //purpose  : Checks is the given vertex lies inside the polygon
1447 //=======================================================================
1448 Standard_Boolean BRepMesh_Delaun::isVertexInsidePolygon(
1449   const Standard_Integer&           theVertexId,
1450   const IMeshData::VectorOfInteger& thePolygonVertices ) const
1451 {
1452   Standard_Integer aPolyLen = thePolygonVertices.Length();
1453   if ( aPolyLen < 3 )
1454     return Standard_False;
1455
1456
1457   const gp_XY aCenterPointXY = GetVertex( theVertexId ).Coord();
1458
1459   const BRepMesh_Vertex& aFirstVertex = GetVertex( thePolygonVertices( 0 ) );
1460   gp_Vec2d aPrevVertexDir( aFirstVertex.Coord() - aCenterPointXY );
1461   if ( aPrevVertexDir.SquareMagnitude() < Precision2 )
1462     return Standard_True;
1463
1464   Standard_Real aTotalAng = 0.0;
1465   for ( Standard_Integer aPolyIt = 1; aPolyIt < aPolyLen; ++aPolyIt )
1466   {
1467     const BRepMesh_Vertex& aPolyVertex = GetVertex( thePolygonVertices( aPolyIt ) );
1468
1469     gp_Vec2d aCurVertexDir( aPolyVertex.Coord() - aCenterPointXY );
1470     if ( aCurVertexDir.SquareMagnitude() < Precision2 )
1471       return Standard_True;
1472
1473     aTotalAng     += aCurVertexDir.Angle( aPrevVertexDir );
1474     aPrevVertexDir = aCurVertexDir;
1475   }
1476
1477   if ( Abs( Angle2PI - aTotalAng ) > Precision::Angular() )
1478     return Standard_False;
1479
1480   return Standard_True;
1481 }
1482
1483 //=======================================================================
1484 //function : killTrianglesOnIntersectingLinks
1485 //purpose  : Checks is the given link crosses the polygon boundary.
1486 //           If yes, kills its triangles and checks neighbor links on
1487 //           boundary intersection. Does nothing elsewhere.
1488 //=======================================================================
1489 void BRepMesh_Delaun::killTrianglesOnIntersectingLinks(
1490   const Standard_Integer&             theLinkToCheckId,
1491   const BRepMesh_Edge&                theLinkToCheck,
1492   const Standard_Integer&             theEndPoint,
1493   const IMeshData::SequenceOfInteger& thePolygon,
1494   const IMeshData::SequenceOfBndB2d&  thePolyBoxes,
1495   IMeshData::MapOfInteger&            theSurvivedLinks,
1496   IMeshData::MapOfIntegerInteger&     theLoopEdges )
1497 {
1498   if ( theSurvivedLinks.Contains( theLinkToCheckId ) )
1499     return;
1500
1501   Bnd_B2d aBox;
1502   Standard_Boolean isNotIntersect =
1503     checkIntersection( theLinkToCheck, thePolygon,
1504       thePolyBoxes, Standard_False, Standard_False,
1505       Standard_False, aBox );
1506
1507   theSurvivedLinks.Add( theLinkToCheckId );
1508
1509   if ( isNotIntersect )
1510     return;
1511
1512   killLinkTriangles( theLinkToCheckId, theLoopEdges );
1513
1514   IMeshData::ListOfInteger::Iterator aNeighborsIt(
1515     myMeshData->LinksConnectedTo(theEndPoint));
1516
1517   for ( ; aNeighborsIt.More(); aNeighborsIt.Next() )
1518   {
1519     const Standard_Integer& aNeighborLinkId = aNeighborsIt.Value();
1520     const BRepMesh_Edge&    aNeighborLink   = GetEdge( aNeighborLinkId );
1521     Standard_Integer anOtherNode = aNeighborLink.FirstNode();
1522     if ( anOtherNode == theEndPoint )
1523       anOtherNode = aNeighborLink.LastNode();
1524
1525     killTrianglesOnIntersectingLinks( aNeighborLinkId,
1526       aNeighborLink, anOtherNode, thePolygon,
1527       thePolyBoxes, theSurvivedLinks, theLoopEdges );
1528   }
1529 }
1530
1531 //=======================================================================
1532 //function : killLinkTriangles
1533 //purpose  : Kill triangles bound to the given link.
1534 //=======================================================================
1535 void BRepMesh_Delaun::killLinkTriangles(
1536   const Standard_Integer&         theLinkId,
1537   IMeshData::MapOfIntegerInteger& theLoopEdges )
1538 {
1539   const BRepMesh_PairOfIndex& aPair =
1540     myMeshData->ElementsConnectedTo( theLinkId );
1541
1542   Standard_Integer anElemNb = aPair.Extent();
1543   for ( Standard_Integer aPairIt = 1; aPairIt <= anElemNb; ++aPairIt )
1544   {
1545     Standard_Integer anElemId = aPair.FirstIndex();
1546     if ( anElemId < 0 )
1547       continue;
1548
1549     deleteTriangle( anElemId, theLoopEdges );
1550   }
1551 }
1552
1553 //=======================================================================
1554 //function : getOrientedNodes
1555 //purpose  : Returns start and end nodes of the given edge in respect to
1556 //           its orientation.
1557 //=======================================================================
1558 void BRepMesh_Delaun::getOrientedNodes(const BRepMesh_Edge&   theEdge,
1559                                        const Standard_Boolean isForward,
1560                                        Standard_Integer       *theNodes) const
1561 {
1562   if ( isForward )
1563   {
1564     theNodes[0] = theEdge.FirstNode();
1565     theNodes[1] = theEdge.LastNode();
1566   }
1567   else
1568   {
1569     theNodes[0] = theEdge.LastNode();
1570     theNodes[1] = theEdge.FirstNode();
1571   }
1572 }
1573
1574 //=======================================================================
1575 //function : processLoop
1576 //purpose  : Processes loop within the given polygon formed by range of
1577 //           its links specified by start and end link indices.
1578 //=======================================================================
1579 void BRepMesh_Delaun::processLoop(const Standard_Integer              theLinkFrom,
1580                                   const Standard_Integer              theLinkTo,
1581                                   const IMeshData::SequenceOfInteger& thePolygon,
1582                                   const IMeshData::SequenceOfBndB2d&  thePolyBoxes)
1583 {
1584   Standard_Integer aNbOfLinksInLoop = theLinkTo - theLinkFrom - 1;
1585   if ( aNbOfLinksInLoop < 3 )
1586     return;
1587
1588   IMeshData::SequenceOfInteger aPolygon;
1589   IMeshData::SequenceOfBndB2d  aPolyBoxes;
1590   for ( ; aNbOfLinksInLoop > 0; --aNbOfLinksInLoop )
1591   {
1592     Standard_Integer aLoopLinkIndex = theLinkFrom + aNbOfLinksInLoop;
1593     aPolygon  .Prepend( thePolygon  ( aLoopLinkIndex ) );
1594     aPolyBoxes.Prepend( thePolyBoxes( aLoopLinkIndex ) );
1595   }
1596   meshPolygon( aPolygon, aPolyBoxes );
1597 }
1598
1599 //=======================================================================
1600 //function : createAndReplacePolygonLink
1601 //purpose  : Creates new link based on the given nodes and updates the
1602 //           given polygon.
1603 //=======================================================================
1604 Standard_Integer BRepMesh_Delaun::createAndReplacePolygonLink(
1605   const Standard_Integer       *theNodes,
1606   const gp_Pnt2d               *thePnts,
1607   const Standard_Integer        theRootIndex,
1608   const ReplaceFlag             theReplaceFlag,
1609   IMeshData::SequenceOfInteger& thePolygon,
1610   IMeshData::SequenceOfBndB2d&  thePolyBoxes )
1611 {
1612   Standard_Integer aNewEdgeId =
1613     myMeshData->AddLink( BRepMesh_Edge(
1614     theNodes[0], theNodes[1], BRepMesh_Free ) );
1615
1616   Bnd_B2d aNewBox;
1617   UpdateBndBox(thePnts[0].Coord(), thePnts[1].Coord(), aNewBox);
1618
1619   switch ( theReplaceFlag )
1620   {
1621   case BRepMesh_Delaun::Replace:
1622     thePolygon  .SetValue( theRootIndex, aNewEdgeId );
1623     thePolyBoxes.SetValue( theRootIndex, aNewBox );
1624     break;
1625
1626   case BRepMesh_Delaun::InsertAfter:
1627     thePolygon  .InsertAfter( theRootIndex, aNewEdgeId );
1628     thePolyBoxes.InsertAfter( theRootIndex, aNewBox );
1629     break;
1630
1631   case BRepMesh_Delaun::InsertBefore:
1632     thePolygon  .InsertBefore( theRootIndex, aNewEdgeId );
1633     thePolyBoxes.InsertBefore( theRootIndex, aNewBox );
1634     break;
1635   }
1636
1637   return aNewEdgeId;
1638 }
1639
1640 //=======================================================================
1641 //function : meshPolygon
1642 //purpose  :
1643 //=======================================================================
1644 void BRepMesh_Delaun::meshPolygon(IMeshData::SequenceOfInteger&   thePolygon,
1645                                   IMeshData::SequenceOfBndB2d&    thePolyBoxes,
1646                                   Handle(IMeshData::MapOfInteger) theSkipped)
1647 {
1648   // Check is the source polygon elementary
1649   if ( meshElementaryPolygon( thePolygon ) )
1650     return;
1651
1652   // Check and correct boundary edges
1653   Standard_Integer aPolyLen = thePolygon.Length();
1654   const Standard_Real aPolyArea      = Abs( polyArea( thePolygon, 1, aPolyLen ) );
1655   const Standard_Real aSmallLoopArea = 0.001 * aPolyArea;
1656   for ( Standard_Integer aPolyIt = 1; aPolyIt < aPolyLen; ++aPolyIt )
1657   {
1658     Standard_Integer aCurEdgeInfo = thePolygon( aPolyIt );
1659     Standard_Integer aCurEdgeId   = Abs( aCurEdgeInfo );
1660     const BRepMesh_Edge* aCurEdge = &GetEdge( aCurEdgeId );
1661     if ( aCurEdge->Movability() != BRepMesh_Frontier )
1662       continue;
1663
1664     Standard_Integer aCurNodes[2];
1665     getOrientedNodes( *aCurEdge, aCurEdgeInfo > 0, aCurNodes );
1666
1667     gp_Pnt2d aCurPnts[2] = {
1668       GetVertex(aCurNodes[0]).Coord(),
1669       GetVertex(aCurNodes[1]).Coord()
1670     };
1671
1672     gp_Vec2d aCurVec( aCurPnts[0], aCurPnts[1] );
1673
1674     // check further links
1675     Standard_Integer aNextPolyIt = aPolyIt + 1;
1676     for ( ; aNextPolyIt <= aPolyLen; ++aNextPolyIt )
1677     {
1678       Standard_Integer aNextEdgeInfo = thePolygon( aNextPolyIt );
1679       Standard_Integer aNextEdgeId   = Abs( aNextEdgeInfo );
1680       const BRepMesh_Edge* aNextEdge = &GetEdge( aNextEdgeId );
1681       if ( aNextEdge->Movability() != BRepMesh_Frontier )
1682         continue;
1683
1684       Standard_Integer aNextNodes[2];
1685       getOrientedNodes( *aNextEdge, aNextEdgeInfo > 0, aNextNodes );
1686
1687       gp_Pnt2d aNextPnts[2] = {
1688         GetVertex(aNextNodes[0]).Coord(),
1689         GetVertex(aNextNodes[1]).Coord()
1690       };
1691
1692       gp_Pnt2d anIntPnt;
1693       BRepMesh_GeomTool::IntFlag aIntFlag = intSegSeg( *aCurEdge, *aNextEdge,
1694         Standard_False, Standard_True, anIntPnt );
1695
1696       if ( aIntFlag == BRepMesh_GeomTool::NoIntersection )
1697         continue;
1698
1699       Standard_Boolean isRemoveFromFirst  = Standard_False;
1700       Standard_Boolean isAddReplacingEdge = Standard_True;
1701       Standard_Integer aIndexToRemoveTo   = aNextPolyIt;
1702       if ( aIntFlag == BRepMesh_GeomTool::Cross )
1703       {
1704         Standard_Real aLoopArea = polyArea( thePolygon, aPolyIt + 1, aNextPolyIt );
1705         gp_Vec2d aVec1( anIntPnt, aCurPnts [1] );
1706         gp_Vec2d aVec2( anIntPnt, aNextPnts[0] );
1707
1708         aLoopArea += ( aVec1 ^ aVec2 ) / 2.;
1709         if ( Abs( aLoopArea ) > aSmallLoopArea )
1710         {
1711           aNextNodes[1] = aCurNodes[0];
1712           aNextPnts [1] = aCurPnts [0];
1713
1714           aNextEdgeId = Abs( createAndReplacePolygonLink( aNextNodes, aNextPnts,
1715             aNextPolyIt, BRepMesh_Delaun::Replace, thePolygon, thePolyBoxes ) );
1716
1717           processLoop( aPolyIt, aNextPolyIt, thePolygon, thePolyBoxes );
1718           return;
1719         }
1720
1721         Standard_Real aDist1 = anIntPnt.SquareDistance(aNextPnts[0]);
1722         Standard_Real aDist2 = anIntPnt.SquareDistance(aNextPnts[1]);
1723
1724         // Choose node with lower distance
1725         const Standard_Boolean isCloseToStart = ( aDist1 < aDist2 );
1726         const Standard_Integer aEndPointIndex = isCloseToStart ? 0 : 1;
1727         aCurNodes[1] = aNextNodes[aEndPointIndex];
1728         aCurPnts [1] = aNextPnts [aEndPointIndex];
1729
1730         if ( isCloseToStart )
1731           --aIndexToRemoveTo;
1732
1733         // In this context only intersections between frontier edges
1734         // are possible. If intersection between edges of different
1735         // types occured - treat this case as invalid (i.e. result
1736         // might not reflect the expectations).
1737         if ( !theSkipped.IsNull() )
1738         {
1739           Standard_Integer aSkippedLinkIt = aPolyIt;
1740           for ( ; aSkippedLinkIt <= aIndexToRemoveTo; ++aSkippedLinkIt )
1741             theSkipped->Add( Abs( thePolygon( aSkippedLinkIt ) ) );
1742         }
1743       }
1744       else if ( aIntFlag == BRepMesh_GeomTool::PointOnSegment )
1745       {
1746         // Indentify chopping link
1747         Standard_Boolean isFirstChopping = Standard_False;
1748         Standard_Integer aCheckPointIt = 0;
1749         for ( ; aCheckPointIt < 2; ++aCheckPointIt )
1750         {
1751           gp_Pnt2d& aRefPoint = aCurPnts[aCheckPointIt];
1752           // Check is second link touches the first one
1753           gp_Vec2d aVec1( aRefPoint, aNextPnts[0] );
1754           gp_Vec2d aVec2( aRefPoint, aNextPnts[1] );
1755           if ( Abs( aVec1 ^ aVec2 ) < Precision )
1756           {
1757             isFirstChopping = Standard_True;
1758             break;
1759           }
1760         }
1761
1762         if ( isFirstChopping )
1763         {
1764           // Split second link
1765           isAddReplacingEdge = Standard_False;
1766           isRemoveFromFirst  = ( aCheckPointIt == 0 );
1767
1768           Standard_Integer aSplitLink[3] = {
1769             aNextNodes[0],
1770             aCurNodes [aCheckPointIt],
1771             aNextNodes[1]
1772           };
1773
1774           gp_Pnt2d aSplitPnts[3] = {
1775             aNextPnts[0],
1776             aCurPnts [aCheckPointIt],
1777             aNextPnts[1]
1778           };
1779
1780           Standard_Integer aSplitLinkIt = 0;
1781           for ( ; aSplitLinkIt < 2; ++aSplitLinkIt )
1782           {
1783             createAndReplacePolygonLink( &aSplitLink[aSplitLinkIt],
1784               &aSplitPnts[aSplitLinkIt], aNextPolyIt, ( aSplitLinkIt == 0 ) ?
1785               BRepMesh_Delaun::Replace : BRepMesh_Delaun::InsertAfter,
1786               thePolygon, thePolyBoxes );
1787           }
1788
1789           processLoop( aPolyIt + aCheckPointIt, aIndexToRemoveTo,
1790             thePolygon, thePolyBoxes );
1791         }
1792         else
1793         {
1794           // Split first link
1795           Standard_Integer aSplitLinkNodes[2] = {
1796             aNextNodes[1],
1797             aCurNodes [1]
1798           };
1799
1800           gp_Pnt2d aSplitLinkPnts[2] = {
1801             aNextPnts[1],
1802             aCurPnts [1]
1803           };
1804           createAndReplacePolygonLink( aSplitLinkNodes, aSplitLinkPnts,
1805             aPolyIt, BRepMesh_Delaun::InsertAfter, thePolygon, thePolyBoxes );
1806
1807           aCurNodes[1] = aNextNodes[1];
1808           aCurPnts [1] = aNextPnts [1];
1809           ++aIndexToRemoveTo;
1810
1811           processLoop( aPolyIt + 1, aIndexToRemoveTo,
1812             thePolygon, thePolyBoxes );
1813         }
1814       }
1815       else if ( aIntFlag == BRepMesh_GeomTool::Glued )
1816       {
1817         if ( aCurNodes[1] == aNextNodes[0] )
1818         {
1819           aCurNodes[1] = aNextNodes[1];
1820           aCurPnts [1] = aNextPnts [1];
1821         }
1822         // TODO: Non-adjacent glued links within the polygon
1823       }
1824       else if ( aIntFlag == BRepMesh_GeomTool::Same )
1825       {
1826         processLoop( aPolyIt, aNextPolyIt, thePolygon, thePolyBoxes );
1827
1828         isRemoveFromFirst  = Standard_True;
1829         isAddReplacingEdge = Standard_False;
1830       }
1831       else
1832         continue; // Not supported type
1833
1834       if ( isAddReplacingEdge )
1835       {
1836         aCurEdgeId = Abs( createAndReplacePolygonLink( aCurNodes, aCurPnts,
1837           aPolyIt, BRepMesh_Delaun::Replace, thePolygon, thePolyBoxes ) );
1838
1839         aCurEdge = &GetEdge( aCurEdgeId );
1840         aCurVec  = gp_Vec2d( aCurPnts[0], aCurPnts[1] );
1841       }
1842
1843       Standard_Integer aIndexToRemoveFrom =
1844         isRemoveFromFirst ? aPolyIt : aPolyIt + 1;
1845
1846       thePolygon  .Remove( aIndexToRemoveFrom, aIndexToRemoveTo );
1847       thePolyBoxes.Remove( aIndexToRemoveFrom, aIndexToRemoveTo );
1848
1849       aPolyLen = thePolygon.Length();
1850       if ( isRemoveFromFirst )
1851       {
1852         --aPolyIt;
1853         break;
1854       }
1855
1856       aNextPolyIt = aPolyIt;
1857     }
1858   }
1859
1860   IMeshData::SequenceOfInteger* aPolygon1   = &thePolygon;
1861   IMeshData::SequenceOfBndB2d*  aPolyBoxes1 = &thePolyBoxes;
1862
1863   Handle(IMeshData::SequenceOfInteger) aPolygon2   = new IMeshData::SequenceOfInteger;
1864   Handle(IMeshData::SequenceOfBndB2d)  aPolyBoxes2 = new IMeshData::SequenceOfBndB2d;
1865
1866   NCollection_Sequence<Handle(IMeshData::SequenceOfInteger)> aPolyStack;
1867   NCollection_Sequence<Handle(IMeshData::SequenceOfBndB2d)>  aPolyBoxStack;
1868   for (;;)
1869   {
1870     decomposeSimplePolygon(*aPolygon1, *aPolyBoxes1, *aPolygon2, *aPolyBoxes2);
1871     if (!aPolygon2->IsEmpty())
1872     {
1873       aPolyStack.Append(aPolygon2);
1874       aPolyBoxStack.Append(aPolyBoxes2);
1875
1876       aPolygon2   = new IMeshData::SequenceOfInteger;
1877       aPolyBoxes2 = new IMeshData::SequenceOfBndB2d;
1878     }
1879
1880     if (aPolygon1->IsEmpty())
1881     {
1882       if (!aPolyStack.IsEmpty() && aPolygon1 == &(*aPolyStack.First()))
1883       {
1884         aPolyStack.Remove(1);
1885         aPolyBoxStack.Remove(1);
1886       }
1887
1888       if (aPolyStack.IsEmpty())
1889         return;
1890
1891       aPolygon1   = &(*aPolyStack.ChangeFirst());
1892       aPolyBoxes1 = &(*aPolyBoxStack.ChangeFirst());
1893     }
1894   }
1895 }
1896
1897 //=======================================================================
1898 //function : meshElementaryPolygon
1899 //purpose  : Triangulation of closed polygon containing only three edges.
1900 //=======================================================================
1901 Standard_Boolean BRepMesh_Delaun::meshElementaryPolygon(
1902   const IMeshData::SequenceOfInteger& thePolygon)
1903 {
1904   Standard_Integer aPolyLen = thePolygon.Length();
1905   if ( aPolyLen < 3 )
1906     return Standard_True;
1907   else if ( aPolyLen > 3 )
1908     return Standard_False;
1909
1910   // Just create a triangle
1911   Standard_Integer anEdges[3];
1912   Standard_Boolean anEdgesOri[3];
1913
1914   for ( Standard_Integer anEdgeIt = 0; anEdgeIt < 3; ++anEdgeIt )
1915   {
1916     Standard_Integer anEdgeInfo = thePolygon( anEdgeIt + 1 );
1917     anEdges[anEdgeIt]           = Abs( anEdgeInfo );
1918     anEdgesOri[anEdgeIt]        = ( anEdgeInfo > 0 );
1919   }
1920
1921   const BRepMesh_Edge& anEdge1 = GetEdge( anEdges[0] );
1922   const BRepMesh_Edge& anEdge2 = GetEdge( anEdges[1] );
1923
1924   Standard_Integer aNodes[3] = { anEdge1.FirstNode(),
1925                                  anEdge1.LastNode(),
1926                                  anEdge2.FirstNode() };
1927   if ( aNodes[2] == aNodes[0] ||
1928        aNodes[2] == aNodes[1] )
1929   {
1930     aNodes[2] = anEdge2.LastNode();
1931   }
1932
1933   addTriangle( anEdges, anEdgesOri, aNodes );
1934   return Standard_True;
1935 }
1936
1937 //=======================================================================
1938 //function : meshSimplePolygon
1939 //purpose  :
1940 //=======================================================================
1941 void BRepMesh_Delaun::decomposeSimplePolygon(
1942   IMeshData::SequenceOfInteger& thePolygon,
1943   IMeshData::SequenceOfBndB2d&  thePolyBoxes,
1944   IMeshData::SequenceOfInteger& thePolygonCut,
1945   IMeshData::SequenceOfBndB2d&  thePolyBoxesCut)
1946 {
1947   // Check is the given polygon elementary
1948   if ( meshElementaryPolygon( thePolygon ) )
1949   {
1950     thePolygon.Clear();
1951     thePolyBoxes.Clear();
1952     return;
1953   }
1954
1955   // Polygon contains more than 3 links
1956   Standard_Integer aFirstEdgeInfo = thePolygon(1);
1957   const BRepMesh_Edge& aFirstEdge = GetEdge( Abs( aFirstEdgeInfo ) );
1958
1959   Standard_Integer aNodes[3];
1960   getOrientedNodes( aFirstEdge, aFirstEdgeInfo > 0, aNodes );
1961
1962   gp_Pnt2d aRefVertices[3];
1963   aRefVertices[0] = GetVertex( aNodes[0] ).Coord();
1964   aRefVertices[1] = GetVertex( aNodes[1] ).Coord();
1965
1966   gp_Vec2d aRefEdgeDir( aRefVertices[0], aRefVertices[1] );
1967
1968   Standard_Real aRefEdgeLen = aRefEdgeDir.Magnitude();
1969   if ( aRefEdgeLen < Precision )
1970   {
1971     thePolygon.Clear();
1972     thePolyBoxes.Clear();
1973     return;
1974   }
1975
1976   aRefEdgeDir /= aRefEdgeLen;
1977
1978   // Find a point with minimum distance respect
1979   // the end of reference link
1980   Standard_Integer aUsedLinkId = 0;
1981   Standard_Real    aOptAngle   = 0.0;
1982   Standard_Real    aMinDist    = RealLast();
1983   Standard_Integer aPivotNode  = aNodes[1];
1984   Standard_Integer aPolyLen    = thePolygon.Length();
1985   for ( Standard_Integer aLinkIt = 3; aLinkIt <= aPolyLen; ++aLinkIt )
1986   {
1987     Standard_Integer aLinkInfo = thePolygon( aLinkIt );
1988     const BRepMesh_Edge& aNextEdge = GetEdge( Abs( aLinkInfo ) );
1989
1990     aPivotNode = aLinkInfo > 0 ?
1991       aNextEdge.FirstNode() :
1992       aNextEdge.LastNode();
1993
1994     gp_Pnt2d aPivotVertex = GetVertex( aPivotNode ).Coord();
1995     gp_Vec2d aDistanceDir( aRefVertices[1], aPivotVertex );
1996
1997     Standard_Real aDist     = aRefEdgeDir ^ aDistanceDir;
1998     Standard_Real aAngle    = Abs( aRefEdgeDir.Angle(aDistanceDir) );
1999     Standard_Real anAbsDist = Abs( aDist );
2000     if (anAbsDist < Precision || aDist < 0.)
2001       continue;
2002
2003     if ( ( anAbsDist >= aMinDist                             ) &&
2004          ( aAngle <= aOptAngle || aAngle > AngDeviation90Deg ) )
2005     {
2006       continue;
2007     }
2008
2009     // Check is the test link crosses the polygon boudaries
2010     Standard_Boolean isIntersect = Standard_False;
2011     for ( Standard_Integer aRefLinkNodeIt = 0; aRefLinkNodeIt < 2; ++aRefLinkNodeIt )
2012     {
2013       const Standard_Integer& aLinkFirstNode   = aNodes[aRefLinkNodeIt];
2014       const gp_Pnt2d&         aLinkFirstVertex = aRefVertices[aRefLinkNodeIt];
2015
2016       Bnd_B2d aBox;
2017       UpdateBndBox(aLinkFirstVertex.Coord(), aPivotVertex.Coord(), aBox);
2018
2019       BRepMesh_Edge aCheckLink( aLinkFirstNode, aPivotNode, BRepMesh_Free );
2020
2021       Standard_Integer aCheckLinkIt = 2;
2022       for ( ; aCheckLinkIt <= aPolyLen; ++aCheckLinkIt )
2023       {
2024         if( aCheckLinkIt == aLinkIt )
2025           continue;
2026
2027         if ( !aBox.IsOut( thePolyBoxes.Value( aCheckLinkIt ) ) )
2028         {
2029           const BRepMesh_Edge& aPolyLink =
2030             GetEdge( Abs( thePolygon( aCheckLinkIt ) ) );
2031
2032           if ( aCheckLink.IsEqual( aPolyLink ) )
2033             continue;
2034
2035           // intersection is possible...
2036           gp_Pnt2d anIntPnt;
2037           BRepMesh_GeomTool::IntFlag aIntFlag = intSegSeg( aCheckLink, aPolyLink,
2038             Standard_False, Standard_False, anIntPnt );
2039
2040           if( aIntFlag != BRepMesh_GeomTool::NoIntersection )
2041           {
2042             isIntersect = Standard_True;
2043             break;
2044           }
2045         }
2046       }
2047
2048       if ( isIntersect )
2049         break;
2050     }
2051
2052     if( isIntersect )
2053       continue;
2054
2055
2056     aOptAngle       = aAngle;
2057     aMinDist        = anAbsDist;
2058     aNodes[2]       = aPivotNode;
2059     aRefVertices[2] = aPivotVertex;
2060     aUsedLinkId     = aLinkIt;
2061   }
2062
2063   if ( aUsedLinkId == 0 )
2064   {
2065     thePolygon.Clear();
2066     thePolyBoxes.Clear();
2067     return;
2068   }
2069
2070
2071   BRepMesh_Edge aNewEdges[2] = {
2072     BRepMesh_Edge( aNodes[1], aNodes[2], BRepMesh_Free ),
2073     BRepMesh_Edge( aNodes[2], aNodes[0], BRepMesh_Free ) };
2074
2075   Standard_Integer aNewEdgesInfo[3] = {
2076     aFirstEdgeInfo,
2077     myMeshData->AddLink( aNewEdges[0] ),
2078     myMeshData->AddLink( aNewEdges[1] ) };
2079
2080
2081   Standard_Integer anEdges[3];
2082   Standard_Boolean anEdgesOri[3];
2083   for ( Standard_Integer aTriEdgeIt = 0; aTriEdgeIt < 3; ++aTriEdgeIt )
2084   {
2085     const Standard_Integer& anEdgeInfo = aNewEdgesInfo[aTriEdgeIt];
2086     anEdges[aTriEdgeIt]    = Abs( anEdgeInfo );
2087     anEdgesOri[aTriEdgeIt] = anEdgeInfo > 0;
2088   }
2089   addTriangle( anEdges, anEdgesOri, aNodes );
2090
2091   // Create triangle and split the source polygon on two
2092   // parts (if possible) and mesh each part as independent
2093   // polygon.
2094   if ( aUsedLinkId < aPolyLen )
2095   {
2096     thePolygon.Split(aUsedLinkId, thePolygonCut);
2097     thePolygonCut.Prepend( -aNewEdgesInfo[2] );
2098     thePolyBoxes.Split(aUsedLinkId, thePolyBoxesCut);
2099
2100     Bnd_B2d aBox;
2101     UpdateBndBox(aRefVertices[0].Coord(), aRefVertices[2].Coord(), aBox);
2102     thePolyBoxesCut.Prepend( aBox );
2103   }
2104   else
2105   {
2106     thePolygon.Remove  ( aPolyLen );
2107     thePolyBoxes.Remove( aPolyLen );
2108   }
2109
2110   if ( aUsedLinkId > 3 )
2111   {
2112     thePolygon.SetValue( 1, -aNewEdgesInfo[1] );
2113
2114     Bnd_B2d aBox;
2115     UpdateBndBox(aRefVertices[1].Coord(), aRefVertices[2].Coord(), aBox);
2116     thePolyBoxes.SetValue( 1, aBox );
2117   }
2118 }
2119
2120 //=======================================================================
2121 //function : RemoveVertex
2122 //purpose  : Removes a vertex from the triangulation
2123 //=======================================================================
2124 void BRepMesh_Delaun::RemoveVertex( const BRepMesh_Vertex& theVertex )
2125 {
2126   BRepMesh_SelectorOfDataStructureOfDelaun aSelector( myMeshData );
2127   aSelector.NeighboursOf( theVertex );
2128
2129   IMeshData::MapOfIntegerInteger aLoopEdges;//( 10, myMeshData->Allocator() );
2130
2131   // Loop on triangles to be destroyed :
2132   IMeshData::IteratorOfMapOfInteger aTriangleIt( aSelector.Elements() );
2133   for ( ; aTriangleIt.More(); aTriangleIt.Next() )
2134     deleteTriangle( aTriangleIt.Key(), aLoopEdges );
2135
2136   IMeshData::SequenceOfBndB2d  aBoxes;
2137   IMeshData::SequenceOfInteger aPolygon;
2138   Standard_Integer aLoopEdgesCount = aLoopEdges.Extent();
2139   IMeshData::MapOfIntegerInteger::Iterator aLoopEdgesIt( aLoopEdges );
2140
2141   if ( aLoopEdgesIt.More() )
2142   {
2143     const BRepMesh_Edge& anEdge = GetEdge( aLoopEdgesIt.Key() );
2144     Standard_Integer aFirstNode = anEdge.FirstNode();
2145     Standard_Integer aLastNode;
2146     Standard_Integer aPivotNode = anEdge.LastNode();
2147     Standard_Integer anEdgeId   = aLoopEdgesIt.Key();
2148
2149     Standard_Boolean isPositive = aLoopEdges( anEdgeId ) != 0;
2150     if ( !isPositive )
2151     {
2152       Standard_Integer aTmp;
2153       aTmp       = aFirstNode;
2154       aFirstNode = aPivotNode;
2155       aPivotNode = aTmp;
2156
2157       aPolygon.Append( -anEdgeId );
2158     }
2159     else
2160       aPolygon.Append( anEdgeId );
2161
2162     fillBndBox( aBoxes, GetVertex( aFirstNode ), GetVertex( aPivotNode ) );
2163
2164     aLoopEdges.UnBind( anEdgeId );
2165
2166     aLastNode = aFirstNode;
2167     while ( aPivotNode != aLastNode )
2168     {
2169       IMeshData::ListOfInteger::Iterator aLinkIt( myMeshData->LinksConnectedTo( aPivotNode ) );
2170       for ( ; aLinkIt.More(); aLinkIt.Next() )
2171       {
2172         if ( aLinkIt.Value() != anEdgeId &&
2173              aLoopEdges.IsBound( aLinkIt.Value() ) )
2174         {
2175           Standard_Integer aCurrentNode;
2176           anEdgeId = aLinkIt.Value();
2177           const BRepMesh_Edge& anEdge1 = GetEdge( anEdgeId );
2178
2179           aCurrentNode = anEdge1.LastNode();
2180           if ( aCurrentNode != aPivotNode )
2181           {
2182             aCurrentNode = anEdge1.FirstNode();
2183             aPolygon.Append( -anEdgeId );
2184           }
2185           else
2186             aPolygon.Append( anEdgeId );
2187
2188           fillBndBox( aBoxes, GetVertex( aCurrentNode ), GetVertex( aPivotNode ) );
2189
2190           aPivotNode = aCurrentNode;
2191           aLoopEdges.UnBind( anEdgeId );
2192           break;
2193         }
2194       }
2195
2196       if ( aLoopEdgesCount <= 0 )
2197         break;
2198       --aLoopEdgesCount;
2199     }
2200
2201     meshPolygon( aPolygon, aBoxes );
2202   }
2203 }
2204
2205
2206 //=======================================================================
2207 //function : AddVertices
2208 //purpose  : Adds some vertices in the triangulation.
2209 //=======================================================================
2210 void BRepMesh_Delaun::AddVertices(IMeshData::VectorOfInteger& theVertices)
2211 {
2212   ComparatorOfIndexedVertexOfDelaun aCmp(myMeshData);
2213   std::make_heap(theVertices.begin(), theVertices.end(), aCmp);
2214   std::sort_heap(theVertices.begin(), theVertices.end(), aCmp);
2215
2216   createTrianglesOnNewVertices(theVertices);
2217 }
2218
2219 //=======================================================================
2220 //function : UseEdge
2221 //purpose  : Modify mesh to use the edge. Return True if done
2222 //=======================================================================
2223 Standard_Boolean BRepMesh_Delaun::UseEdge( const Standard_Integer /*theIndex*/ )
2224 {
2225   /*
2226   const BRepMesh_PairOfIndex& aPair = myMeshData->ElemConnectedTo( theIndex );
2227   if ( aPair.Extent() == 0 )
2228   {
2229     const BRepMesh_Edge& anEdge = GetEdge( theIndex );
2230
2231     Standard_Integer aStartNode, aPivotNode, anOtherNode;
2232     aStartNode = anEdge.FirstNode();
2233     aPivotNode = anEdge.LastNode();
2234
2235     const BRepMesh_ListOfInteger& aStartNodeNeighbors = myMeshData->LinkNeighboursOf( aStartNode );
2236     const BRepMesh_ListOfInteger& aPivotNodeNeighbors = myMeshData->LinkNeighboursOf( aPivotNode );
2237
2238     if ( aStartNodeNeighbors.Extent() > 0 &&
2239          aPivotNodeNeighbors.Extent() > 0 )
2240     {
2241       const BRepMesh_Vertex& aStartVertex = GetVertex( aStartNode );
2242       const BRepMesh_Vertex& aPivotVertex = GetVertex( aPivotNode );
2243
2244       gp_XY aVEdge   ( aPivotVertex.Coord() );
2245       aVEdge.Subtract( aStartVertex.Coord() );
2246
2247       Standard_Real    anAngle    = 0.;
2248       Standard_Real    anAngleMin = RealLast();
2249       Standard_Real    anAngleMax = RealFirst();
2250       Standard_Integer aLeftEdge  = 0, aRightEdge = 0;
2251
2252       BRepMesh_ListOfInteger::Iterator aNeighborIt( aPivotNodeNeighbors );
2253       for ( ; aNeighborIt.More(); aNeighborIt.Next() )
2254       {
2255         Standard_Integer anEdgeId = aNeighborIt.Value();
2256         if ( anEdgeId != theIndex )
2257         {
2258           const BRepMesh_Edge& aNextEdge = GetEdge( anEdgeId );
2259
2260           Standard_Boolean isInMesh = Standard_True;
2261           if ( aNextEdge.Movability() == BRepMesh_Free )
2262           {
2263             if ( myMeshData->ElemConnectedTo( anEdgeId ).IsEmpty() )
2264               isInMesh = Standard_False;
2265           }
2266
2267           if ( isInMesh )
2268           {
2269             anOtherNode = aNextEdge.FirstNode();
2270             if ( anOtherNode == aPivotNode )
2271               anOtherNode = aNextEdge.LastNode();
2272
2273             gp_XY aVEdgeCur = GetVertex( anOtherNode ).Coord();
2274             aVEdgeCur.Subtract( aPivotVertex.Coord() );
2275
2276             anAngle = gp_Vec2d( aVEdge ).Angle( gp_Vec2d( aVEdgeCur ) );
2277           }
2278
2279           if ( anAngle > anAngleMax )
2280           {
2281             anAngleMax = anAngle;
2282             aLeftEdge  = anEdgeId;
2283           }
2284           if ( anAngle < anAngleMin )
2285           {
2286             anAngleMin = anAngle;
2287             aRightEdge = anEdgeId;
2288           }
2289         }
2290       }
2291
2292       if ( aLeftEdge > 0 )
2293       {
2294         if (aLeftEdge==aRightEdge)
2295         {
2296         }
2297         else
2298         {
2299         }
2300       }
2301     }
2302   }
2303   */
2304   return Standard_False;
2305 }
2306
2307 //=======================================================================
2308 //function : getEdgesByType
2309 //purpose  : Gives the list of edges with type defined by input parameter
2310 //=======================================================================
2311 Handle(IMeshData::MapOfInteger) BRepMesh_Delaun::getEdgesByType(
2312   const BRepMesh_DegreeOfFreedom theEdgeType ) const
2313 {
2314   Handle(NCollection_IncAllocator) anAlloc = new NCollection_IncAllocator;
2315   Handle(IMeshData::MapOfInteger) aResult = new IMeshData::MapOfInteger;
2316   IMeshData::IteratorOfMapOfInteger anEdgeIt( myMeshData->LinksOfDomain() );
2317
2318   for ( ; anEdgeIt.More(); anEdgeIt.Next() )
2319   {
2320     Standard_Integer anEdge = anEdgeIt.Key();
2321     Standard_Boolean isToAdd = (theEdgeType == BRepMesh_Free) ?
2322       (myMeshData->ElementsConnectedTo( anEdge ).Extent() <= 1) :
2323       (GetEdge( anEdge ).Movability() == theEdgeType);
2324
2325     if (isToAdd)
2326       aResult->Add( anEdge );
2327   }
2328
2329   return aResult;
2330 }
2331
2332 //=======================================================================
2333 //function : calculateDist
2334 //purpose  : Calculates distances between the given point and edges of
2335 //           triangle
2336 //=======================================================================
2337 Standard_Real BRepMesh_Delaun::calculateDist( const gp_XY            theVEdges[3],
2338                                               const gp_XY            thePoints[3],
2339                                               const BRepMesh_Vertex& theVertex,
2340                                               Standard_Real          theDistance[3],
2341                                               Standard_Real          theSqModulus[3],
2342                                               Standard_Integer&      theEdgeOn ) const
2343 {
2344   Standard_Real aMinDist = RealLast();
2345   for( Standard_Integer i = 0; i < 3; ++i )
2346   {
2347     theSqModulus[i] = theVEdges[i].SquareModulus();
2348     if ( theSqModulus[i] <= Precision2 )
2349       return -1;
2350
2351     theDistance[i] = theVEdges[i] ^ ( theVertex.Coord() - thePoints[i] );
2352
2353     Standard_Real aDist = theDistance[i] * theDistance[i];
2354     aDist /= theSqModulus[i];
2355
2356     if ( aDist < aMinDist )
2357     {
2358       theEdgeOn = i;
2359       aMinDist  = aDist;
2360     }
2361   }
2362
2363   return aMinDist;
2364 }
2365
2366 //=======================================================================
2367 //function : Contains
2368 //purpose  : Test if triangle of index <TrianIndex> contains geometricaly
2369 //           <theVertex>. If <theEdgeOn> is != 0 then theVertex is on Edge
2370 //           of  index <theEdgeOn>
2371 //=======================================================================
2372 Standard_Boolean BRepMesh_Delaun::Contains( const Standard_Integer theTriangleId,
2373                                             const BRepMesh_Vertex& theVertex,
2374                                             const Standard_Real    theSqTolerance,
2375                                             Standard_Integer&      theEdgeOn) const
2376 {
2377   theEdgeOn = 0;
2378
2379   Standard_Integer p[3];
2380
2381   const BRepMesh_Triangle& aElement = GetTriangle( theTriangleId );
2382   const Standard_Integer(&e)[3] = aElement.myEdges;
2383
2384   const BRepMesh_Edge* anEdges[3] = { &GetEdge( e[0] ),
2385                                       &GetEdge( e[1] ),
2386                                       &GetEdge( e[2] ) };
2387
2388   myMeshData->ElementNodes(aElement, p);
2389
2390   gp_XY aPoints[3];
2391   aPoints[0] = GetVertex( p[0] ).Coord();
2392   aPoints[1] = GetVertex( p[1] ).Coord();
2393   aPoints[2] = GetVertex( p[2] ).Coord();
2394
2395   gp_XY aVEdges[3];
2396   aVEdges[0] = aPoints[1];
2397   aVEdges[0].Subtract( aPoints[0] );
2398
2399   aVEdges[1] = aPoints[2];
2400   aVEdges[1].Subtract( aPoints[1] );
2401
2402   aVEdges[2] = aPoints[0];
2403   aVEdges[2].Subtract( aPoints[2] );
2404
2405   Standard_Real aDistance[3];
2406   Standard_Real aSqModulus[3];
2407
2408   Standard_Real aSqMinDist;
2409   Standard_Integer aEdgeOnId;
2410   aSqMinDist = calculateDist( aVEdges, aPoints, theVertex, aDistance, aSqModulus, aEdgeOnId );
2411   if ( aSqMinDist < 0 )
2412     return Standard_False;
2413
2414   const Standard_Boolean isNotFree = (anEdges[aEdgeOnId]->Movability() != BRepMesh_Free);
2415   if ( aSqMinDist > theSqTolerance )
2416   {
2417     if (isNotFree && aDistance[aEdgeOnId] < ( aSqModulus[aEdgeOnId] / 5. ))
2418       theEdgeOn = e[aEdgeOnId];
2419   }
2420   else if (isNotFree)
2421     return Standard_False;
2422   else
2423     theEdgeOn = e[aEdgeOnId];
2424
2425   return (aDistance[0] >= 0. && aDistance[1] >= 0. && aDistance[2] >= 0.);
2426 }
2427
2428 //=============================================================================
2429 //function : intSegSeg
2430 //purpose  : Checks intersection between the two segments.
2431 //=============================================================================
2432 BRepMesh_GeomTool::IntFlag BRepMesh_Delaun::intSegSeg(
2433   const BRepMesh_Edge&   theEdg1,
2434   const BRepMesh_Edge&   theEdg2,
2435   const Standard_Boolean isConsiderEndPointTouch,
2436   const Standard_Boolean isConsiderPointOnEdge,
2437   gp_Pnt2d&              theIntPnt) const
2438 {
2439   gp_XY p1, p2, p3, p4;
2440   p1 = GetVertex( theEdg1.FirstNode() ).Coord();
2441   p2 = GetVertex( theEdg1.LastNode()  ).Coord();
2442   p3 = GetVertex( theEdg2.FirstNode() ).Coord();
2443   p4 = GetVertex( theEdg2.LastNode()  ).Coord();
2444
2445   return BRepMesh_GeomTool::IntSegSeg(p1, p2, p3, p4,
2446     isConsiderEndPointTouch, isConsiderPointOnEdge, theIntPnt);
2447 }
2448
2449 //=============================================================================
2450 //function : polyArea
2451 //purpose  : Returns area of the loop of the given polygon defined by indices
2452 //           of its start and end links.
2453 //=============================================================================
2454 Standard_Real BRepMesh_Delaun::polyArea(const IMeshData::SequenceOfInteger& thePolygon,
2455                                         const Standard_Integer              theStartIndex,
2456                                         const Standard_Integer              theEndIndex) const
2457 {
2458   Standard_Real aArea = 0.0;
2459   Standard_Integer aPolyLen = thePolygon.Length();
2460   if ( theStartIndex >= theEndIndex ||
2461        theStartIndex > aPolyLen )
2462   {
2463     return aArea;
2464   }
2465   Standard_Integer aCurEdgeInfo = thePolygon( theStartIndex );
2466   Standard_Integer aCurEdgeId   = Abs( aCurEdgeInfo );
2467   const BRepMesh_Edge* aCurEdge = &GetEdge( aCurEdgeId );
2468
2469   Standard_Integer aNodes[2];
2470   getOrientedNodes( *aCurEdge, aCurEdgeInfo > 0, aNodes );
2471
2472   gp_Pnt2d aRefPnt = GetVertex( aNodes[0] ).Coord();
2473   Standard_Integer aPolyIt = theStartIndex + 1;
2474   for ( ; aPolyIt <= theEndIndex; ++aPolyIt )
2475   {
2476     aCurEdgeInfo = thePolygon( aPolyIt );
2477     aCurEdgeId   = Abs( aCurEdgeInfo );
2478     aCurEdge     = &GetEdge( aCurEdgeId );
2479
2480     getOrientedNodes( *aCurEdge, aCurEdgeInfo > 0, aNodes );
2481     gp_Vec2d aVec1( aRefPnt, GetVertex( aNodes[0] ).Coord() );
2482     gp_Vec2d aVec2( aRefPnt, GetVertex( aNodes[1] ).Coord() );
2483
2484     aArea += aVec1 ^ aVec2;
2485   }
2486
2487   return aArea / 2.;
2488 }
2489
2490 #ifdef OCCT_DEBUG
2491 //=======================================================================
2492 //function : BRepMesh_DumpPoly
2493 //purpose  :
2494 //=======================================================================
2495 #include <TopoDS_Compound.hxx>
2496 #include <BRep_Builder.hxx>
2497 #include <Standard_ErrorHandler.hxx>
2498 #include <BRepBuilderAPI_MakeEdge.hxx>
2499 #include <BRepTools.hxx>
2500 Standard_CString BRepMesh_DumpPoly(void*            thePolygon,
2501                                    void*            theMeshHandlePtr,
2502                                    Standard_CString theFileNameStr)
2503 {
2504   if (thePolygon == 0 || theFileNameStr == 0)
2505   {
2506     return "Error: file name or polygon data is null";
2507   }
2508
2509   IMeshData::SequenceOfInteger& aPolygon = *(IMeshData::SequenceOfInteger*)thePolygon;
2510
2511   Handle(BRepMesh_DataStructureOfDelaun) aMeshData =
2512     *(Handle(BRepMesh_DataStructureOfDelaun)*)theMeshHandlePtr;
2513
2514   if (aMeshData.IsNull())
2515     return "Error: mesh data is empty";
2516
2517   TopoDS_Compound aMesh;
2518   BRep_Builder aBuilder;
2519   aBuilder.MakeCompound(aMesh);
2520
2521   try
2522   {
2523     OCC_CATCH_SIGNALS
2524
2525     IMeshData::SequenceOfInteger::Iterator aLinksIt(aPolygon);
2526     for (; aLinksIt.More(); aLinksIt.Next())
2527     {
2528       const BRepMesh_Edge& aLink = aMeshData->GetLink(Abs(aLinksIt.Value()));
2529
2530       gp_Pnt aPnt[2];
2531       for (Standard_Integer i = 0; i < 2; ++i)
2532       {
2533         const Standard_Integer aNodeId =
2534           (i == 0) ? aLink.FirstNode() : aLink.LastNode();
2535
2536         const gp_XY& aNode = aMeshData->GetNode(aNodeId).Coord();
2537         aPnt[i] = gp_Pnt(aNode.X(), aNode.Y(), 0.);
2538       }
2539
2540       if (aPnt[0].SquareDistance(aPnt[1]) < Precision::SquareConfusion())
2541         continue;
2542
2543       aBuilder.Add(aMesh, BRepBuilderAPI_MakeEdge(aPnt[0], aPnt[1]));
2544     }
2545
2546     if (!BRepTools::Write(aMesh, theFileNameStr))
2547       return "Error: write failed";
2548   }
2549   catch (Standard_Failure const& anException)
2550   {
2551     return anException.GetMessageString();
2552   }
2553
2554   return theFileNameStr;
2555 }
2556 #endif