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