src/BRepMesh/BRepMesh_FastDiscret.cxx

   1 // Created on: 1996-02-27
   2 // Created by: Ekaterina SMIRNOVA
   3 // Copyright (c) 1996-1999 Matra Datavision
   4 // Copyright (c) 1999-2012 OPEN CASCADE SAS
   5 //
   6 // The content of this file is subject to the Open CASCADE Technology Public
   7 // License Version 6.5 (the "License"). You may not use the content of this file
   8 // except in compliance with the License. Please obtain a copy of the License
   9 // at http://www.opencascade.org and read it completely before using this file.
  10 //
  11 // The Initial Developer of the Original Code is Open CASCADE S.A.S., having its
  12 // main offices at: 1, place des Freres Montgolfier, 78280 Guyancourt, France.
  13 //
  14 // The Original Code and all software distributed under the License is
  15 // distributed on an "AS IS" basis, without warranty of any kind, and the
  16 // Initial Developer hereby disclaims all such warranties, including without
  17 // limitation, any warranties of merchantability, fitness for a particular
  18 // purpose or non-infringement. Please see the License for the specific terms
  19 // and conditions governing the rights and limitations under the License.
  20
  21
  22 #include <BRepMesh_FastDiscret.ixx>
  23
  24 #include <BRepMesh_FastDiscretFace.hxx>
  25 #include <BRepMesh_FaceAttribute.hxx>
  26 #include <BRepMesh_DataStructureOfDelaun.hxx>
  27 #include <BRepMesh_ClassifierPtr.hxx>
  28 #include <BRepMesh_GeomTool.hxx>
  29 #include <BRepMesh_PairOfPolygon.hxx>
  30 #include <BRepMesh_DataMapOfShapePairOfPolygon.hxx>
  31 #include <BRepMesh_DataMapIteratorOfDataMapOfShapePairOfPolygon.hxx>
  32 #include <Geom_Plane.hxx>
  33 #include <GeomAbs_IsoType.hxx>
  34 #include <GeomAbs_SurfaceType.hxx>
  35 #include <TopAbs.hxx>
  36 #include <TColStd_HArray1OfReal.hxx>
  37 #include <Precision.hxx>
  38
  39 #include <BRep_Builder.hxx>
  40 #include <BRep_Tool.hxx>
  41 #include <Poly_Triangulation.hxx>
  42 #include <Poly_PolygonOnTriangulation.hxx>
  43 #include <Poly_Connect.hxx>
  44 #include <TColStd_SequenceOfInteger.hxx>
  45 #include <TColStd_Array1OfInteger.hxx>
  46 #include <TColStd_HArray1OfInteger.hxx>
  47
  48 #include <TColgp_Array1OfPnt.hxx>
  49 #include <TColgp_Array1OfPnt2d.hxx>
  50 #include <Precision.hxx>
  51
  52 #include <BRepAdaptor_Curve.hxx>
  53 #include <BRepAdaptor_Surface.hxx>
  54 #include <BRepAdaptor_HSurface.hxx>
  55 #include <BRepTools.hxx>
  56 #include <BndLib_Add3dCurve.hxx>
  57 #include <BRepBndLib.hxx>
  58 #include <Bnd_Box.hxx>
  59 #include <TopoDS.hxx>
  60 #include <TopExp.hxx>
  61 #include <TopExp_Explorer.hxx>
  62
  63 #include <Geom2d_Curve.hxx>
  64
  65 #include <TColStd_DataMapOfIntegerInteger.hxx>
  66 #include <BRepMesh_ShapeTool.hxx>
  67 #include <ElSLib.hxx>
  68 #include <Geom_Surface.hxx>
  69 #include <Adaptor3d_IsoCurve.hxx>
  70 #include <BRepMesh_IndexedMapOfVertex.hxx>
  71 #include <Extrema_LocateExtPC.hxx>
  72
  73 #include <BRepMesh_ListOfXY.hxx>
  74 #include <BRepMesh_ListIteratorOfListOfXY.hxx>
  75
  76 #include <TColStd_Array1OfInteger.hxx>
  77 #include <BRepMesh_IDMapOfNodeOfDataStructureOfDelaun.hxx>
  78 #include <Standard_ErrorHandler.hxx>
  79 #include <Standard_Failure.hxx>
  80 //#include <TColStd_DataMapOfInteger.hxx>
  81 #include <TColGeom2d_SequenceOfCurve.hxx>
  82 #include <TopTools_SequenceOfShape.hxx>
  83 #include <NCollection_IncAllocator.hxx>
  84
  85 #include <BRep_ListIteratorOfListOfPointRepresentation.hxx>
  86 #include <BRep_PointRepresentation.hxx>
  87 #include <BRep_TVertex.hxx>
  88 #include <TColStd_MapOfInteger.hxx>
  89 #include <SortTools_ShellSortOfReal.hxx>
  90 #include <TCollection_CompareOfReal.hxx>
  91
  92 #include <TopTools_HArray1OfShape.hxx>
  93 #include <TopTools_ListIteratorOfListOfShape.hxx>
  94
  95 #include <vector>
  96
  97 #ifdef HAVE_TBB
  98   // paralleling using Intel TBB
  99   #include <tbb/parallel_for_each.h>
 100 #endif
 101
 102 #define UVDEFLECTION 1.e-05
 103
 104 inline Standard_Real MaxFaceTol (const TopoDS_Face& theFace)
 105 {
 106   Standard_Real T, TMax = BRep_Tool::Tolerance(theFace);
 107   TopExp_Explorer Ex;
 108
 109   for (Ex.Init(theFace,TopAbs_EDGE); Ex.More(); Ex.Next())
 110   {
 111     T = BRep_Tool::Tolerance(TopoDS::Edge(Ex.Current()));
 112     if (T > TMax) TMax = T;
 113   }
 114
 115   for (Ex.Init(theFace,TopAbs_VERTEX); Ex.More(); Ex.Next())
 116   {
 117     T = BRep_Tool::Tolerance(TopoDS::Vertex(Ex.Current()));
 118     if (T > TMax) TMax = T;
 119   }
 120
 121   return TMax;
 122 }
 123
 124
 125 //=======================================================================
 126 //function : BRepMesh_FastDiscret
 127 //purpose  :
 128 //=======================================================================
 129 BRepMesh_FastDiscret::BRepMesh_FastDiscret(const Standard_Real    theDefle,
 130                                            const Standard_Real    theAngl,
 131                                            const Bnd_Box&         theBox,
 132                                            const Standard_Boolean theWithShare,
 133                                            const Standard_Boolean theInshape,
 134                                            const Standard_Boolean theRelative,
 135                                            const Standard_Boolean theShapetrigu) :
 136   myAngle (theAngl),
 137   myDeflection (theDefle),
 138   myWithShare (theWithShare),
 139   myInParallel (Standard_False),
 140   myNbLocat (0),
 141   myRelative (theRelative),
 142   myShapetrigu (theShapetrigu),
 143   myInshape (theInshape)
 144 {
 145   myAllocator = new NCollection_IncAllocator(64000);
 146   if(myRelative)
 147     BoxMaxDimension(theBox, myDtotale);
 148 }
 149
 150 //=======================================================================
 151 //function : BRepMesh_FastDiscret
 152 //purpose  :
 153 //=======================================================================
 154
 155 BRepMesh_FastDiscret::BRepMesh_FastDiscret(const Standard_Real    theDefle,
 156                                            const TopoDS_Shape&    theShape,
 157                                            const Bnd_Box&         theBox,
 158                                            const Standard_Real    theAngl,
 159                                            const Standard_Boolean theWithShare,
 160                                            const Standard_Boolean theInshape,
 161                                            const Standard_Boolean theRelative,
 162                                            const Standard_Boolean theShapetrigu):
 163   myAngle (theAngl),
 164   myDeflection (theDefle),
 165   myWithShare (theWithShare),
 166   myInParallel (Standard_False),
 167   myNbLocat (0),
 168   myRelative (theRelative),
 169   myShapetrigu (theShapetrigu),
 170   myInshape (theInshape)
 171 {
 172   myAllocator = new NCollection_IncAllocator(64000);
 173   if(myRelative)
 174     BoxMaxDimension(theBox, myDtotale);
 175   Perform(theShape);
 176 }
 177
 178 //=======================================================================
 179 //function : SetParallel
 180 //purpose  :
 181 //=======================================================================
 182 void BRepMesh_FastDiscret::SetParallel (const Standard_Boolean theInParallel)
 183 {
 184   myInParallel = theInParallel;
 185 }
 186
 187 //=======================================================================
 188 //function : IsParallel
 189 //purpose  :
 190 //=======================================================================
 191 Standard_Boolean BRepMesh_FastDiscret::IsParallel() const
 192 {
 193   return myInParallel;
 194 }
 195
 196 //=======================================================================
 197 //function : BoxMaxDimension
 198 //purpose  :
 199 //=======================================================================
 200
 201 void BRepMesh_FastDiscret::BoxMaxDimension(const Bnd_Box& theBox, Standard_Real& theMaxDim)
 202 {
 203   if(theBox.IsVoid())
 204     return;
 205   Standard_Real TXmin, TYmin, TZmin, TXmax, TYmax, TZmax;
 206   theBox.Get(TXmin, TYmin, TZmin, TXmax, TYmax, TZmax);
 207   theMaxDim = TXmax-TXmin;
 208   const Standard_Real dy = TYmax-TYmin;
 209   const Standard_Real dz = TZmax-TZmin;
 210   if (dy > theMaxDim) theMaxDim = dy;
 211   if (dz > theMaxDim) theMaxDim = dz;
 212 }
 213
 214 //=======================================================================
 215 //function : RelativeEdgeDeflection
 216 //purpose  :
 217 //=======================================================================
 218
 219 Standard_Real BRepMesh_FastDiscret::RelativeEdgeDeflection(const TopoDS_Edge& theEdge,
 220                                                                                                    const Standard_Real theDefle,
 221                                                                                                      const Standard_Real theDTotale,
 222                                                            Standard_Real& theDefCoef)
 223 {
 224   theDefCoef = 1.;
 225   Standard_Real defedge = theDefle;
 226   if(theEdge.IsNull())
 227     return defedge;
 228
 229   Bnd_Box B;
 230   BRepBndLib::Add(theEdge, B);
 231   BoxMaxDimension(B, defedge);
 232
 233   // adjusted in relation to the total size:
 234   theDefCoef = theDTotale/(2*defedge);
 235   if (theDefCoef < 0.5) theDefCoef = 0.5;
 236   if (theDefCoef > 2.) theDefCoef = 2.;
 237   defedge = theDefCoef * defedge * theDefle;
 238
 239   return defedge;
 240 }
 241
 242 //=======================================================================
 243 //function : Perform(shape)
 244 //purpose  :
 245 //=======================================================================
 246
 247 void BRepMesh_FastDiscret::Perform(const TopoDS_Shape& theShape)
 248 {
 249   TopTools_IndexedDataMapOfShapeListOfShape anAncestors;
 250   TopExp::MapShapesAndAncestors(theShape, TopAbs_EDGE, TopAbs_FACE, anAncestors);
 251   std::vector<TopoDS_Face> aFaces;
 252   for (TopExp_Explorer ex(theShape, TopAbs_FACE); ex.More(); ex.Next()) {
 253     TopoDS_Face aF = TopoDS::Face(ex.Current());
 254     Add(aF, anAncestors);
 255     aFaces.push_back(aF);
 256   }
 257
 258   if (myInParallel)
 259   {
 260   #ifdef HAVE_TBB
 261     CreateMutexesForSubShapes(theShape, TopAbs_EDGE);
 262     // mesh faces in parallel threads using TBB
 263     tbb::parallel_for_each (aFaces.begin(), aFaces.end(), *this);
 264   #else
 265     // alternative parallelization not yet available
 266     for (std::vector<TopoDS_Face>::iterator it(aFaces.begin()); it != aFaces.end(); it++)
 267       Process (*it);
 268   #endif
 269     RemoveAllMutexes();
 270   }
 271   else
 272   {
 273     for (std::vector<TopoDS_Face>::iterator it(aFaces.begin()); it != aFaces.end(); it++)
 274       Process (*it);
 275   }
 276 }
 277
 278
 279 //=======================================================================
 280 //function : Process
 281 //purpose  :
 282 //=======================================================================
 283
 284 void BRepMesh_FastDiscret::Process(const TopoDS_Face& theFace) const
 285 {
 286   //cout << "START face " << theFace.TShape().operator->() << endl << flush;
 287   Handle(BRepMesh_FaceAttribute) fattribute;
 288   if ( GetFaceAttribute (theFace, fattribute) )
 289   {
 290     BRepMesh_FastDiscretFace aTool (GetAngle(), WithShare());
 291     aTool.Add (theFace, fattribute, GetMapOfDefEdge(), myMutexProvider);
 292   }
 293   //cout << "END   face " << theFace.TShape().operator->() << endl << flush;
 294 }
 295
 296 //=======================================================================
 297 //function : Add(face)
 298 //purpose  :
 299 //=======================================================================
 300
 301 #define MESH_FAILURE(theface)        \
 302 {                                    \
 303   myFacestate = BRepMesh_Failure;    \
 304   myNottriangulated.Append(theface); \
 305   return;                            \
 306 }
 307
 308 void BRepMesh_FastDiscret::Add(const TopoDS_Face& theface,
 309                                const TopTools_IndexedDataMapOfShapeListOfShape& theAncestors)
 310 {
 311 #ifndef DEB_MESH
 312   try
 313   {
 314     OCC_CATCH_SIGNALS
 315 #endif
 316   TopoDS_Face face = theface;
 317   BRepTools::Update(face);
 318   face.Orientation(TopAbs_FORWARD);
 319   myStructure.Nullify();
 320   Handle(NCollection_IncAllocator) anAlloc = Handle(NCollection_IncAllocator)::DownCast(myAllocator);
 321   anAlloc->Reset(Standard_False);
 322   myStructure=new BRepMesh_DataStructureOfDelaun(anAlloc);
 323
 324   Standard_Real aUmin, aVmin, aUmax, aVmax;
 325   BRepTools::UVBounds (theface, aUmin, aUmax, aVmin, aVmax);
 326   Standard_Real aTolU = (aUmax - aUmin) * UVDEFLECTION;
 327   Standard_Real aTolV = (aVmax - aVmin) * UVDEFLECTION;
 328   myStructure->Data().SetCellSize ( 14 * aTolU, 14 * aTolV );
 329   myStructure->Data().SetTolerance( aTolU, aTolV );
 330
 331   BRepAdaptor_Surface  BS(face, Standard_False);
 332   Handle(BRepAdaptor_HSurface) gFace = new BRepAdaptor_HSurface(BS);
 333
 334   GeomAbs_SurfaceType thetype;
 335   thetype = BS.GetType();
 336
 337   gp_Pnt2d uvFirst, uvLast;
 338
 339   Handle(Poly_Triangulation) T;
 340   TopLoc_Location loc;
 341
 342   if (!myWithShare) {
 343     myVertices.Clear();
 344     myEdges.Clear();
 345   }
 346
 347   myVemap.Clear();
 348   myLocation2d.Clear();
 349   myInternaledges.Clear();
 350
 351   Standard_Integer i;
 352   i = 1;
 353
 354   Standard_Real defedge, defface;
 355   Standard_Integer nbEdge = 0;
 356   Standard_Real savangle = myAngle;
 357   Standard_Real cdef;
 358   Standard_Real maxdef = 2.* MaxFaceTol(theface);
 359   defface = 0.;
 360
 361   if (!myRelative) defface = Max(myDeflection, maxdef);
 362
 363   TColStd_SequenceOfReal aFSeq, aLSeq;
 364   TColGeom2d_SequenceOfCurve aCSeq;
 365   TopTools_SequenceOfShape aShSeq;
 366
 367   TopoDS_Iterator exW(face);
 368
 369   for (; exW.More(); exW.Next()) {
 370     const TopoDS_Shape& aWire = exW.Value();
 371     if (aWire.ShapeType() != TopAbs_WIRE)
 372       continue;
 373     TopoDS_Iterator ex(aWire);
 374     for(; ex.More(); ex.Next()) {
 375       const TopoDS_Edge& edge = TopoDS::Edge(ex.Value());
 376       nbEdge++;
 377       if (!myMapdefle.IsBound(edge)) {
 378         if (myRelative) {
 379           if (myEdges.IsBound(edge)) {
 380             const BRepMesh_PairOfPolygon& pair = myEdges.Find(edge);
 381             const Handle(Poly_PolygonOnTriangulation)& P = pair.First();
 382             defedge = P->Deflection();
 383           }
 384           else {
 385             defedge = RelativeEdgeDeflection(edge, myDeflection, myDtotale, cdef);
 386             myAngle = savangle * cdef;
 387           }
 388           defface = defface + defedge;
 389           defface = Max(maxdef, defface);
 390         }
 391         else defedge = myDeflection;
 392
 393         defedge = Max(maxdef, defedge);
 394         defedge = Max(UVDEFLECTION , defedge);
 395         myMapdefle.Bind(edge, defedge);
 396       }
 397       else{
 398         defedge = myMapdefle(edge);
 399         if (myRelative) {defface = defface + defedge; defface = Max(maxdef, defface);}
 400       }
 401       Standard_Real f1,l1;
 402       Handle(Geom2d_Curve) C = BRep_Tool::CurveOnSurface(edge, face, f1, l1);
 403       if (C.IsNull()) continue;
 404
 405       aFSeq.Append(f1);
 406       aLSeq.Append(l1);
 407       aCSeq.Append(C);
 408       aShSeq.Append(edge);
 409       Add(edge, face, gFace, C, theAncestors, defedge, f1, l1);
 410       myAngle = savangle;
 411     }
 412   }
 413
 414   if (nbEdge == 0 || myVemap.Extent() < 3)
 415   {
 416       MESH_FAILURE(theface);
 417   }
 418
 419   if (myRelative ) defface = defface / nbEdge;
 420   else             defface = myDeflection;
 421
 422   if (myWithShare) defface = Max(maxdef, defface);
 423
 424   T = BRep_Tool::Triangulation(face, loc);
 425
 426   if (!myShapetrigu || T.IsNull()) {
 427
 428     Standard_Real xCur, yCur;
 429     Standard_Real maxX, minX, maxY, minY;
 430     minX=minY=1.e100;
 431     maxX=maxY=-1.e100;
 432
 433     Standard_Integer ipn = 0;
 434     Standard_Integer i1 =1;
 435     for (i1 = 1; i1 <= myVemap.Extent(); i1++) {
 436       const BRepMesh_Vertex& aV = myStructure->GetNode(myVemap.FindKey(i1));
 437       ipn++;
 438       xCur=aV.Coord().X();
 439       yCur=aV.Coord().Y();
 440       minX=Min(xCur, minX);
 441       maxX=Max(xCur, maxX);
 442       minY=Min(yCur, minY);
 443       maxY=Max(yCur, maxY);
 444     }
 445     Standard_Real myumin = minX;
 446     Standard_Real myumax = maxX;
 447     Standard_Real myvmin = minY;
 448     Standard_Real myvmax = maxY;
 449
 450     const Standard_Real umin = BS.FirstUParameter();
 451     const Standard_Real umax = BS.LastUParameter();
 452     const Standard_Real vmin = BS.FirstVParameter();
 453     const Standard_Real vmax = BS.LastVParameter();
 454
 455     if (myumin < umin || myumax > umax)
 456     {
 457       if (BS.IsUPeriodic())
 458       {
 459         if ((myumax - myumin) > (umax - umin))
 460         {
 461           myumax = myumin + (umax - umin);
 462         }
 463       }
 464       else
 465       {
 466         if (umin > myumin) myumin = umin;
 467         if (umax < myumax) myumax = umax;
 468       }
 469     }
 470
 471     if (myvmin < vmin || myvmax > vmax)
 472     {
 473       if (BS.IsVPeriodic())
 474       {
 475         if ((myvmax - myvmin) > (vmax - vmin))
 476         {
 477           myvmax = myvmin + (vmax - vmin);
 478         }
 479       }
 480       else
 481       {
 482         if (vmin > myvmin) myvmin = vmin;
 483         if (vmax < myvmax) myvmax = vmax;
 484       }
 485     }
 486
 487     // fast verification of the validity of calculated limits. If wrong,
 488     // sure a problem of pcurve.
 489     if (thetype == GeomAbs_BezierSurface &&
 490         (myumin < -0.5 || myumax > 1.5 || myvmin < -0.5 || myvmax > 1.5))
 491     {
 492       MESH_FAILURE(theface);
 493     }
 494
 495     //define parameters for correct parametrics
 496
 497     Standard_Real deltaX = 1.0;
 498     Standard_Real deltaY = 1.0;
 499     Standard_Integer nbVertices = myVemap.Extent();
 500     const Standard_Real tolclass = Precision::PConfusion(); //0.03*Max(myumax-myumin, myvmax-myvmin);
 501
 502     BRepMesh_ClassifierPtr classifier (
 503       new BRepMesh_Classifier(face, tolclass,  myInternaledges, myVemap,
 504                               myStructure, myumin, myumax, myvmin, myvmax) );
 505
 506     myFacestate = classifier->State();
 507     if (myFacestate == BRepMesh_SelfIntersectingWire)
 508     {
 509       Standard_Integer nbmaill = 0;
 510       Standard_Real eps = Precision::Confusion();
 511       while (nbmaill < 5 && myFacestate != BRepMesh_ReMesh)
 512       {
 513         nbmaill++;
 514
 515         //clear the structure of links
 516         myStructure.Nullify();
 517         myStructure = new BRepMesh_DataStructureOfDelaun(anAlloc);
 518
 519         myVemap.Clear();
 520         myLocation2d.Clear();
 521         myInternaledges.Clear();
 522
 523         Standard_Integer j1;
 524         for(j1 = 1; j1 <= aShSeq.Length(); j1++)
 525         {
 526           const TopoDS_Edge& edge = TopoDS::Edge(aShSeq.Value(j1));
 527           if (myEdges.IsBound(edge))
 528           {
 529             myEdges.UnBind(edge);
 530             myInternaledges.UnBind(edge);
 531           }
 532         }
 533
 534
 535         for( j1 = 1; j1 <= aShSeq.Length(); j1++)
 536         {
 537           const TopoDS_Edge& edge = TopoDS::Edge(aShSeq.Value(j1));
 538           defedge = myMapdefle(edge) / 3.;
 539           defedge = Max(defedge, eps);
 540           myMapdefle.Bind(edge, defedge);
 541           const Handle(Geom2d_Curve)& C = aCSeq.Value(j1);
 542           Add(edge, face, gFace, C, theAncestors, defedge, aFSeq.Value(j1), aLSeq.Value(j1));
 543         }
 544
 545         classifier.Nullify();
 546
 547         classifier = new BRepMesh_Classifier(face, tolclass, myInternaledges, myVemap,
 548                                              myStructure, myumin, myumax, myvmin, myvmax);
 549
 550         if (classifier->State() == BRepMesh_NoError)
 551         {
 552           myFacestate = BRepMesh_ReMesh;
 553         }
 554         nbVertices = myVemap.Extent();
 555       }
 556     }
 557
 558     if (myFacestate != BRepMesh_NoError && myFacestate != BRepMesh_ReMesh)
 559     {
 560       myNottriangulated.Append(face);
 561       classifier.Nullify();
 562       return;
 563     }
 564
 565     // try to find the real length:
 566     // akm (bug OCC16) : We must calculate these measures in non-singular
 567     //     parts of face. Let's try to compute average value of three
 568     //     (umin, (umin+umax)/2, umax), and respectively for v.
 569     //                 vvvvv
 570     Standard_Real longu = 0.0, longv = 0.0; //, last , first;
 571     gp_Pnt P11, P12, P21, P22, P31, P32;
 572
 573     Standard_Real du = (myumax-myumin)/20;
 574     Standard_Real dv = (myvmax-myvmin)/20;
 575     Standard_Real dfuave=(myumin+myumax)/2, dfucur;
 576     Standard_Real dfvave=(myvmin+myvmax)/2, dfvcur;
 577     // U loop
 578     BS.D0 (myumin, myvmin, P11);
 579     BS.D0 (myumin, dfvave, P21);
 580     BS.D0 (myumin, myvmax, P31);
 581     for (i1=1, dfucur=myumin+du; i1 <= 20; i1++, dfucur+=du)  {
 582       BS.D0 (dfucur, myvmin, P12);
 583       BS.D0 (dfucur, dfvave, P22);
 584       BS.D0 (dfucur, myvmax, P32);
 585       longu += ( P11.Distance(P12) + P21.Distance(P22) + P31.Distance(P32) );
 586       P11 = P12;
 587       P21 = P22;
 588       P31 = P32;
 589     }
 590     // V loop
 591     BS.D0(myumin, myvmin, P11);
 592     BS.D0(dfuave, myvmin, P21);
 593     BS.D0(myumax, myvmin, P31);
 594     for (i1=1, dfvcur=myvmin+dv; i1 <= 20; i1++, dfvcur+=dv) {
 595       BS.D0 (myumin, dfvcur, P12);
 596       BS.D0 (dfuave, dfvcur, P22);
 597       BS.D0 (myumax, dfvcur, P32);
 598       longv += ( P11.Distance(P12) + P21.Distance(P22) + P31.Distance(P32) );
 599       P11 = P12;
 600       P21 = P22;
 601       P31 = P32;
 602     }
 603     longu /= 3.;
 604     longv /= 3.;
 605     // akm (bug OCC16) ^^^^^
 606
 607     if (longu <= 1.e-16 || longv <= 1.e-16) {
 608       //yes, it is seen!!
 609 #ifdef DEB_MESH_CHRONO
 610       chMaillEdges.Stop();
 611       MESH_CHRONO_TSTOP(thetype);
 612 #endif
 613       MESH_FAILURE(theface);
 614     }
 615
 616
 617     if (thetype == GeomAbs_Torus)  {
 618       gp_Torus Tor = BS.Torus();
 619       Standard_Real r = Tor.MinorRadius(), R = Tor.MajorRadius();
 620       Standard_Real Du, Dv;//, pasu, pasv;
 621
 622       Dv = Max(1.0e0 - (defface/r),0.0e0) ;
 623       Standard_Real oldDv = 2.0 * ACos (Dv);
 624       oldDv = Min(oldDv, myAngle);
 625       Dv  =  0.9*oldDv; //TWOTHIRD * oldDv;
 626       Dv = oldDv;
 627
 628       Standard_Integer nbV = Max((Standard_Integer)((myvmax-myvmin)/Dv), 2);
 629       Dv = (myvmax-myvmin)/(nbV+1);
 630
 631       Standard_Real ru = R + r;
 632       if (ru > 1.e-16) {
 633         Du = Max(1.0e0 - (defface/ru),0.0e0);
 634         Du  = (2.0 * ACos (Du));
 635         Du = Min(Du, myAngle);
 636         Standard_Real aa = sqrt(Du*Du + oldDv*oldDv);
 637         if(aa < gp::Resolution())
 638           return;
 639
 640         Du = Du * Min(oldDv, Du) / aa;
 641       }
 642       else Du = Dv;
 643
 644       Standard_Integer nbU = Max((Standard_Integer)((myumax-myumin)/Du), 2);
 645       nbU = Max(nbU, (Standard_Integer)(nbV*(myumax-myumin)*R/((myvmax-myvmin)*r)/5.));
 646
 647       Du = (myumax-myumin)/(nbU+1);
 648       //-- DeltaX and DeltaY are chosen so that to avoid "jumping"
 649       //-- of points on the grid
 650       deltaX = Du;
 651       deltaY = Dv;
 652     }
 653     else if (thetype == GeomAbs_Cylinder) {
 654       /*Standard_Real aMax = Max(myumax,myvmax);
 655       deltaX = 0.01;  //1./aMax; //0.01;
 656       deltaY = 1.0;*/
 657       gp_Cylinder Cyl = BS.Cylinder();
 658       Standard_Real R = Cyl.Radius();
 659       // Calculate parameters for iteration in U direction
 660       Standard_Real Du = 1.0 - (defface/R);
 661       if (Du < 0.0) Du = 0.0;
 662       Du = 2.0 * ACos (Du);
 663       if (Du > GetAngle()) Du = GetAngle();
 664       deltaX = Du/longv;
 665       deltaY = 1.;
 666     }
 667     else {
 668       deltaX = (myumax-myumin)/longu;
 669       deltaY = (myvmax-myvmin)/longv;
 670     }
 671
 672     if(!myMapattrib.IsBound(theface))
 673     {
 674       Handle(BRepMesh_FaceAttribute) aFA = new BRepMesh_FaceAttribute();
 675       aFA->GetDefFace() = defface;
 676       aFA->GetUMax() = myumax;
 677       aFA->GetVMax() = myvmax;
 678       aFA->GetUMin() = myumin;
 679       aFA->GetVMin() = myvmin;
 680       aFA->GetDeltaX() = deltaX;
 681       aFA->GetDeltaY() = deltaY;
 682       aFA->GetMinX() = minX;
 683       aFA->GetMinY() = minY;
 684       aFA->GetClassifier() = classifier;
 685       myMapattrib.Bind(theface, aFA);
 686     }
 687
 688     //Standard_Integer nbVertices = myVemap.Extent();
 689     T = new Poly_Triangulation(nbVertices, 1, Standard_True);
 690     TColgp_Array1OfPnt&  Nodes = T->ChangeNodes();
 691     TColgp_Array1OfPnt2d& Nodes2d = T->ChangeUVNodes();
 692
 693     Standard_Integer index;
 694     for (i = 1; i <= nbVertices; i++) {
 695       index = myVemap.FindKey(i);
 696       Nodes(i) = Pnt(index);
 697       Nodes2d(i).SetXY(Vertex(index).Coord());
 698     }
 699
 700     // storage of triangulation in the BRep.
 701     //TopLoc_Location loc = face.Location();
 702     if (!loc.IsIdentity()) {
 703       gp_Trsf tr = loc.Transformation();
 704       tr.Invert();
 705       for (i = Nodes.Lower(); i <= Nodes.Upper(); i++)
 706         Nodes(i).Transform(tr);
 707     }
 708
 709     BRep_Builder B;
 710     B.UpdateFace(face, T);
 711
 712     BRepMesh_DataMapIteratorOfDataMapOfShapePairOfPolygon It(myInternaledges);
 713     for (; It.More(); It.Next()) {
 714       const BRepMesh_PairOfPolygon& pair = It.Value();
 715       const Handle(Poly_PolygonOnTriangulation)& NOD1 = pair.First();
 716       const Handle(Poly_PolygonOnTriangulation)& NOD2 = pair.Last();
 717       if ( NOD1 == NOD2 )
 718         B.UpdateEdge(TopoDS::Edge(It.Key()), NOD1, T, loc);
 719       else
 720         B.UpdateEdge(TopoDS::Edge(It.Key()), NOD1, NOD2, T, loc);
 721     }
 722   }
 723
 724 #ifndef DEB_MESH
 725   }
 726   catch(Standard_Failure)
 727   {
 728     MESH_FAILURE(theface);
 729   }
 730 #endif // DEB_MESH
 731   myStructure.Nullify();
 732 }
 733
 734 //=======================================================================
 735 //function : splitSegment
 736 //purpose  :
 737 //=======================================================================
 738 static void splitSegment( BRepMesh_GeomTool&    theGT,
 739                           const Handle(Geom_Surface)& theSurf,
 740                           const Handle(Geom2d_Curve)& theCurve2d,
 741                           const BRepAdaptor_Curve&    theBAC,
 742                           const Standard_Real         theSquareEDef,
 743                           const Standard_Real         theFirst,
 744                           const Standard_Real         theLast,
 745                           const Standard_Integer      theNbIter)
 746 {
 747   //limit ineration depth
 748   if(theNbIter > 10)
 749     return;
 750   gp_Pnt2d uvf, uvl, uvm;
 751   gp_Pnt   P3dF, P3dL, midP3d, midP3dFromSurf;
 752   Standard_Real midpar;
 753
 754   if(Abs(theLast - theFirst) < 2*Precision::PConfusion())
 755     return;
 756
 757   theCurve2d->D0(theFirst, uvf);
 758   theCurve2d->D0(theLast, uvl);
 759
 760   P3dF = theSurf->Value(uvf.X(), uvf.Y());
 761   P3dL = theSurf->Value(uvl.X(), uvl.Y());
 762
 763   if(P3dF.SquareDistance(P3dL) < theSquareEDef)
 764     return;
 765
 766   uvm = gp_Pnt2d((uvf.XY() + uvl.XY())*0.5);
 767   midP3dFromSurf = theSurf->Value(uvm.X(), uvm.Y());
 768
 769   gp_XYZ aVec = P3dL.XYZ()-P3dF.XYZ();
 770   aVec.Normalize();
 771
 772   gp_XYZ Vec1 = midP3dFromSurf.XYZ() - P3dF.XYZ();
 773   Standard_Real aModulus = Vec1.Dot(aVec);
 774   gp_XYZ aProj = aVec*aModulus;
 775   gp_XYZ aDist = Vec1 - aProj;
 776
 777   if(aDist.SquareModulus() < theSquareEDef)
 778     return;
 779
 780   midpar = (theFirst + theLast) * 0.5;
 781   theBAC.D0(midpar, midP3d);
 782   theGT.AddPoint(midP3d, midpar, Standard_False);
 783
 784   splitSegment(theGT, theSurf, theCurve2d, theBAC, theSquareEDef, theFirst, midpar, theNbIter+1);
 785   splitSegment(theGT, theSurf, theCurve2d, theBAC, theSquareEDef, midpar, theLast,  theNbIter+1);
 786 }
 787
 788 //=======================================================================
 789 //function : Add
 790 //purpose  :
 791 //=======================================================================
 792 void BRepMesh_FastDiscret::Add( const TopoDS_Edge&                  theEdge,
 793                                 const TopoDS_Face&                  theFace,
 794                                 const Handle(BRepAdaptor_HSurface)& theGFace,
 795                                 const Handle(Geom2d_Curve)&         theC2d,
 796                                 const TopTools_IndexedDataMapOfShapeListOfShape& theAncestors,
 797                                 const Standard_Real                 theDefEdge,
 798                                 const Standard_Real                 theFirst,
 799                                 const Standard_Real                 theLast)
 800 {
 801   const TopAbs_Orientation orEdge = theEdge.Orientation();
 802   if (orEdge == TopAbs_EXTERNAL) return;
 803
 804   const Standard_Boolean isEdgeBound = myEdges.IsBound(theEdge);
 805
 806   if (!isEdgeBound)
 807   {
 808     if (Update(theEdge, theFace, theC2d, theDefEdge, theFirst, theLast))
 809     {
 810       return;
 811     }
 812   }
 813
 814   TopoDS_Vertex pBegin, pEnd;
 815   TopExp::Vertices(theEdge, pBegin, pEnd);
 816   if (pBegin.IsNull() || pEnd.IsNull())
 817   {
 818     return;
 819   }
 820
 821   Standard_Real wFirst, wLast;
 822   BRep_Tool::Range(theEdge, theFace, wFirst, wLast);
 823
 824   gp_Pnt2d uvFirst, uvLast;
 825   BRep_Tool::UVPoints(theEdge, theFace, uvFirst, uvLast);
 826
 827   const Standard_Boolean sameUV = uvFirst.IsEqual(uvLast, Precision::PConfusion());
 828
 829   //Control vertexes tolerances
 830   gp_Pnt pFirst = theGFace->Value(uvFirst.X(), uvFirst.Y());
 831   gp_Pnt pLast  = theGFace->Value(uvLast.X(), uvLast.Y());
 832
 833   Standard_Real mindist = 10. * Max(pFirst.Distance(BRep_Tool::Pnt(pBegin)),
 834                                     pLast.Distance(BRep_Tool::Pnt(pEnd)));
 835
 836   if(mindist < BRep_Tool::Tolerance(pBegin) ||
 837      mindist < BRep_Tool::Tolerance(pEnd) ) mindist = theDefEdge;
 838
 839   // control of degenerated non-coded edges.
 840
 841   Standard_Boolean degener = BRep_Tool::Degenerated(theEdge);
 842
 843   if (!degener)
 844   {
 845     if (pBegin.IsSame(pEnd))
 846     {
 847       // calculation of the length of the edge in 3D
 848       Standard_Real longueur = 0.0;
 849       Standard_Real du = (wLast-wFirst)/20;
 850       gp_Pnt P1, P2;
 851       BRepAdaptor_Curve BC(theEdge);
 852       BC.D0(wFirst, P1);
 853       Standard_Real tolV = BRep_Tool::Tolerance(pBegin);
 854       Standard_Real tolV2 = 1.2*tolV;
 855       for (Standard_Integer l = 1; l <= 20; l++) {
 856         BC.D0(wFirst + l*du, P2);
 857         longueur += P1.Distance(P2);
 858         if (longueur > tolV2) break;
 859         P1 = P2;
 860       }
 861
 862       if (longueur < tolV2)
 863       {
 864         degener = Standard_True;
 865       }
 866     }
 867   }
 868
 869   // Correct UV points
 870   if (sameUV)
 871   {
 872     // 1. is it really sameUV without being degenerated
 873     gp_Pnt2d uvF, uvL;
 874     theC2d->D0(theFirst, uvF);
 875     theC2d->D0(theLast, uvL);
 876     if (!uvFirst.IsEqual(uvF, Precision::PConfusion()))
 877     {
 878       uvFirst = uvF;
 879     }
 880     if (!uvLast.IsEqual(uvL, Precision::PConfusion()))
 881     {
 882       uvLast = uvL;
 883     }
 884   }
 885
 886   gp_XY theUV;
 887
 888   // Process first vertex
 889   Standard_Integer ipf;
 890   if (myVertices.IsBound(pBegin))
 891   {
 892     ipf = myVertices.Find(pBegin);
 893   }
 894   else
 895   {
 896     if (sameUV && myVertices.IsBound(pEnd))
 897     {
 898       ipf = myVertices.Find(pEnd);
 899     }
 900     else
 901     {
 902       myNbLocat++;
 903       ipf = myNbLocat;
 904       myLocation3d.Bind(ipf, BRep_Tool::Pnt(pBegin));
 905     }
 906     myVertices.Bind(pBegin, ipf);
 907   }
 908
 909   Handle(BRepMesh_FaceAttribute) aFaceAttribute;
 910   GetFaceAttribute ( theFace, aFaceAttribute );
 911   theUV = BRepMesh_FastDiscretFace::FindUV(pBegin, uvFirst,
 912     ipf, theGFace, mindist, aFaceAttribute, myLocation2d);
 913
 914   BRepMesh_Vertex vf(theUV, ipf, BRepMesh_Frontier);
 915   Standard_Integer ivf = myStructure->AddNode(vf);
 916
 917   // Process last vertex
 918   Standard_Integer ipl;
 919   if (pEnd.IsSame(pBegin))
 920   {
 921     ipl = ipf;
 922   }
 923   else
 924   {
 925     if (myVertices.IsBound(pEnd))
 926     {
 927       ipl = myVertices.Find(pEnd);
 928     }
 929     else
 930     {
 931       if (sameUV)
 932       {
 933         ipl = ipf;
 934       }
 935       else
 936       {
 937         myNbLocat++;
 938         ipl = myNbLocat;
 939         myLocation3d.Bind(ipl, BRep_Tool::Pnt(pEnd));
 940       }
 941       myVertices.Bind(pEnd,ipl);
 942     }
 943   }
 944
 945   theUV = BRepMesh_FastDiscretFace::FindUV(pEnd, uvLast, ipl,
 946     theGFace, mindist, aFaceAttribute, myLocation2d);
 947
 948   BRepMesh_Vertex vl(theUV, ipl, BRepMesh_Frontier);
 949   Standard_Integer ivl= myStructure->AddNode(vl);
 950
 951   Standard_Integer isvf = myVemap.FindIndex(ivf);
 952   if (isvf == 0) isvf = myVemap.Add(ivf);
 953   Standard_Integer isvl = myVemap.FindIndex(ivl);
 954   if (isvl == 0) isvl = myVemap.Add(ivl);
 955
 956   Standard_Real otherdefedge = 0.5*theDefEdge;
 957   gp_Pnt2d uv;
 958   BRepMesh_Vertex v2;
 959   gp_Pnt   P3d;
 960
 961   Handle(Poly_PolygonOnTriangulation) P1;
 962
 963   if (!isEdgeBound)
 964   {
 965     Handle(Poly_PolygonOnTriangulation) P2;
 966
 967     if (degener)
 968     {
 969       // creation anew:
 970       TColStd_Array1OfInteger Nodes(1, 2), NodInStruct(1, 2);
 971       TColStd_Array1OfReal Param(1, 2);
 972
 973       NodInStruct(1) = ipf;
 974       Nodes(1) = isvf;
 975       Param(1) = wFirst;
 976
 977       NodInStruct(2) = ipl;
 978       Nodes(2) = isvl;
 979       Param(2) = wLast;
 980
 981       P1 = new Poly_PolygonOnTriangulation(Nodes, Param);
 982       P2 = new Poly_PolygonOnTriangulation(NodInStruct, Param);
 983     }
 984     else
 985     {
 986       if (orEdge == TopAbs_INTERNAL) otherdefedge *= 0.5;
 987
 988       BRepAdaptor_Curve cons;
 989       Standard_Boolean isSameParam = BRep_Tool::SameParameter(theEdge);
 990       if (isSameParam)
 991       {
 992         cons.Initialize(theEdge);
 993       }
 994       else
 995       {
 996         cons.Initialize(theEdge, theFace);
 997       }
 998
 999       TopLoc_Location L;
1000       Standard_Integer nbpmin = 2;
1001       const GeomAbs_CurveType aCurveType = cons.GetType();
1002       if ( aCurveType == GeomAbs_Circle )
1003         nbpmin = 4; //OCC287
1004
1005       BRepMesh_GeomTool GT(cons, wFirst, wLast, 0.5 * myAngle, otherdefedge, nbpmin);
1006
1007       if ( aCurveType == GeomAbs_BSplineCurve )
1008       {
1009         const Standard_Integer aNbInt = cons.NbIntervals( GeomAbs_C1 );
1010         if ( aNbInt > 0 )
1011         {
1012           TColStd_Array1OfReal anIntervals( 1, aNbInt + 1 );
1013           cons.Intervals( anIntervals, GeomAbs_C1 );
1014           for ( Standard_Integer aIntIt = 1; aIntIt <= aNbInt; ++aIntIt )
1015           {
1016             const Standard_Real& aStartInt = anIntervals.Value( aIntIt );
1017             const Standard_Real& anEndInt  = anIntervals.Value( aIntIt + 1 );
1018
1019             BRepMesh_GeomTool aDetalizator( cons, aStartInt, anEndInt,
1020               0.5 * myAngle, otherdefedge, nbpmin );
1021
1022             Standard_Integer aNbAddNodes = aDetalizator.NbPoints();
1023             for ( Standard_Integer aNodeIt = 1; aNodeIt <= aNbAddNodes; ++aNodeIt )
1024             {
1025               Standard_Real aParam;
1026               gp_Pnt        aPoint3d;
1027               gp_Pnt2d      aPoint2d;
1028               aDetalizator.Value( cons, theGFace, aNodeIt, aParam, aPoint3d, aPoint2d );
1029               GT.AddPoint( aPoint3d, aParam, Standard_False );
1030             }
1031           }
1032         }
1033       }
1034
1035       // PTv, chl/922/G9, Take into account internal vertices
1036       // it is necessary for internal edges, which do not split other edges, by their vertex
1037       TopoDS_Iterator exV(theEdge);
1038       for ( ; exV.More(); exV.Next() )
1039       {
1040         TopoDS_Vertex aIntV = TopoDS::Vertex(exV.Value());
1041         if ( aIntV.Orientation() == TopAbs_INTERNAL )
1042          GT.AddPoint(BRep_Tool::Pnt(aIntV),
1043                      BRep_Tool::Parameter(aIntV, theEdge),
1044                      Standard_True);
1045       }
1046
1047       Standard_Integer i;
1048       Standard_Integer nbnodes = GT.NbPoints();
1049       //Check deflection in 2d space for improvement of edge tesselation.
1050       if( isSameParam && nbnodes > 1)
1051       {
1052         Standard_Real aSquareEdgeDef = otherdefedge * otherdefedge;
1053         const TopTools_ListOfShape& lf = theAncestors.FindFromKey(theEdge);
1054         TopTools_ListIteratorOfListOfShape itl(lf);
1055         for (; itl.More(); itl.Next()) {
1056           const TopoDS_Face& aFace = TopoDS::Face (itl.Value());
1057
1058           TopLoc_Location aLoc;
1059           Standard_Real aF, aL;
1060           Handle(Geom_Surface) aSurf = BRep_Tool::Surface(aFace, aLoc);
1061           const Handle(Standard_Type)& aType = aSurf->DynamicType();
1062           if(aType == STANDARD_TYPE(Geom_Plane))
1063             continue;
1064           Handle(Geom2d_Curve) aCurve2d = BRep_Tool::CurveOnSurface(theEdge, aFace, aF, aL);
1065           if(Abs(aF-wFirst)>Precision::PConfusion()||Abs(aL-wLast)>Precision::PConfusion())
1066             continue;
1067
1068           gp_Pnt2d uvf;
1069           Standard_Real parf;
1070           nbnodes = GT.NbPoints();
1071           TColStd_Array1OfReal aParamArray(1, nbnodes);
1072           for (i = 1; i <= nbnodes; i++)
1073           {
1074             GT.Value(cons, theGFace, i, parf, P3d, uvf);
1075             aParamArray.SetValue(i, parf);
1076           }
1077           for (i = 1; i < nbnodes; i++)
1078           {
1079             splitSegment(GT, aSurf, aCurve2d, cons, aSquareEdgeDef, aParamArray(i), aParamArray(i+1), 1);
1080           }
1081         }
1082       }
1083
1084       // Creation of polygons on triangulation:
1085       Standard_Real puv;
1086       nbnodes = GT.NbPoints();
1087
1088       TColStd_Array1OfInteger Nodes(1, nbnodes);
1089       TColStd_Array1OfInteger NodInStruct(1, nbnodes);
1090       TColStd_Array1OfReal Param(1, nbnodes);
1091
1092       // processing of the 1st point:
1093       Nodes(1) = isvf;
1094       NodInStruct(1) = ipf;
1095       Param(1) = wFirst;
1096
1097       // last point:
1098       Nodes(nbnodes) = isvl;
1099       NodInStruct(nbnodes) = ipl;
1100       Param(nbnodes) = wLast;
1101
1102       if(nbnodes > 2)
1103       {
1104         Standard_Integer iv2;
1105         for (i = 2; i < GT.NbPoints(); i++)
1106         {
1107           // Record 3d point
1108           GT.Value(cons, theGFace, i, puv, P3d, uv);
1109           myNbLocat++;
1110           myLocation3d.Bind(myNbLocat, P3d);
1111           NodInStruct(i) = myNbLocat;
1112           // Record 2d point
1113           v2.Initialize(uv.Coord(), myNbLocat, BRepMesh_OnCurve);
1114           iv2=myStructure->AddNode(v2);
1115
1116           Standard_Integer isv = myVemap.FindIndex(iv2);
1117           if (isv == 0) isv = myVemap.Add(iv2);
1118           Nodes(i) = isv;
1119           NodInStruct(i) = myNbLocat;
1120           Param(i) = puv;
1121
1122           if (orEdge == TopAbs_FORWARD)
1123             myStructure->AddLink(BRepMesh_Edge(ivf, iv2, BRepMesh_Frontier));
1124           else if (orEdge == TopAbs_REVERSED)
1125             myStructure->AddLink(BRepMesh_Edge(iv2, ivf, BRepMesh_Frontier));
1126           else if (orEdge == TopAbs_INTERNAL)
1127             myStructure->AddLink(BRepMesh_Edge(ivf, iv2, BRepMesh_Fixed));
1128           ivf = iv2;
1129         }
1130       }
1131
1132       P1 = new Poly_PolygonOnTriangulation(Nodes, Param);
1133       P2 = new Poly_PolygonOnTriangulation(NodInStruct, Param);
1134     }
1135
1136     P2->Deflection(otherdefedge);
1137     BRepMesh_PairOfPolygon pair;
1138     pair.Append(P2);
1139     myEdges.Bind(theEdge, pair);
1140
1141     if (ivf != ivl) {
1142       if (orEdge == TopAbs_FORWARD)
1143         myStructure->AddLink(BRepMesh_Edge(ivf, ivl, BRepMesh_Frontier));
1144       else if (orEdge == TopAbs_REVERSED)
1145         myStructure->AddLink(BRepMesh_Edge(ivl, ivf, BRepMesh_Frontier));
1146       else if (orEdge == TopAbs_INTERNAL)
1147         myStructure->AddLink(BRepMesh_Edge(ivf, ivl, BRepMesh_Fixed));
1148     }
1149
1150
1151   }
1152   // If this Edge has been already checked and it is not degenerated,
1153   // the points of the polygon calculated at the first check are retrieved :
1154   else
1155   {
1156     if (degener)
1157     {
1158       // Create 2d polygon for degenerated edge
1159       TColStd_Array1OfInteger Nodes(1, 2);
1160       TColStd_Array1OfReal PPar(1, 2);
1161
1162       Nodes(1) = isvf;
1163       PPar(1) = wFirst;
1164
1165       Nodes(2) = isvl;
1166       PPar(2) = wLast;
1167
1168       P1 = new Poly_PolygonOnTriangulation(Nodes, PPar);
1169     }
1170     else
1171     {
1172       // retrieve the polygone:
1173       const BRepMesh_PairOfPolygon& pair = myEdges.Find(theEdge);
1174       const Handle(Poly_PolygonOnTriangulation)& P = pair.First();
1175       const TColStd_Array1OfInteger& NOD = P->Nodes();
1176       Handle(TColStd_HArray1OfReal) Par = P->Parameters();
1177
1178       // creation anew:
1179       const Standard_Integer nbnodes = NOD.Length();
1180       TColStd_Array1OfInteger Nodes(1, nbnodes);
1181       TColStd_Array1OfReal PPar(1, nbnodes);
1182       Standard_Integer iv2;
1183
1184       Nodes(1) = isvf;
1185       PPar(1) = wFirst;
1186
1187       Nodes(nbnodes) = isvl;
1188       PPar(nbnodes) = wLast;
1189
1190       if (nbnodes > 2)
1191       {
1192         Standard_Integer i;
1193         if (BRep_Tool::SameParameter(theEdge))
1194         {
1195           for (i = 2; i < nbnodes; i++)
1196           {
1197             const Standard_Real puv = Par->Value(i);
1198             theC2d->D0(puv, uv);
1199             v2.Initialize(uv.Coord(), NOD(i), BRepMesh_OnCurve);
1200             iv2 = myStructure->AddNode(v2);
1201
1202             Standard_Integer isv = myVemap.FindIndex(iv2);
1203             if (isv == 0) isv = myVemap.Add(iv2);
1204             Nodes(i) = isv;
1205             PPar(i) = puv;
1206
1207             if (orEdge==TopAbs_FORWARD)
1208               myStructure->AddLink(BRepMesh_Edge(ivf, iv2, BRepMesh_Frontier));
1209             else if (orEdge == TopAbs_REVERSED)
1210               myStructure->AddLink(BRepMesh_Edge(iv2, ivf, BRepMesh_Frontier));
1211             else if (orEdge == TopAbs_INTERNAL)
1212               myStructure->AddLink(BRepMesh_Edge(ivf, iv2, BRepMesh_Fixed));
1213
1214             ivf = iv2;
1215           }
1216         }
1217         else
1218         {
1219           const Standard_Real wFold = Par->Value(Par->Lower());
1220           const Standard_Real wLold = Par->Value(Par->Upper());
1221
1222           Standard_Real wKoef = 1.;
1223           if ((wFold != wFirst || wLold != wLast) && wLold != wFold)
1224           {
1225             wKoef = (wLast - wFirst) / (wLold - wFold);
1226           }
1227
1228           BRepAdaptor_Curve cons(theEdge, theFace);
1229           Extrema_LocateExtPC pcos;
1230           pcos.Initialize(cons, cons.FirstParameter(),
1231               cons.LastParameter(), Precision::PConfusion());
1232
1233           Standard_Real wPrev;
1234           Standard_Real wCur      = wFirst;
1235           Standard_Real wCurFound = wFirst;
1236           for (i = 2; i < nbnodes; i++)
1237           {
1238             P3d = myLocation3d(NOD(i));
1239             // Record 2d point
1240             wPrev = wCur;
1241             wCur  = wFirst + wKoef*(Par->Value(i) - wFold);
1242                   wCurFound += (wCur - wPrev);
1243             pcos.Perform(P3d, wCurFound);
1244             if (pcos.IsDone()) wCurFound = pcos.Point().Parameter();
1245             theC2d->D0(wCurFound, uv);
1246             v2.Initialize(uv.Coord(), NOD(i), BRepMesh_OnCurve);
1247             iv2 = myStructure->AddNode(v2);
1248
1249             Standard_Integer isv = myVemap.FindIndex(iv2);
1250             if (isv == 0) isv = myVemap.Add(iv2);
1251                   Nodes(i) = isv;
1252             PPar(i) = wCurFound;
1253
1254             if (orEdge==TopAbs_FORWARD)
1255               myStructure->AddLink(BRepMesh_Edge(ivf, iv2, BRepMesh_Frontier));
1256             else if (orEdge == TopAbs_REVERSED)
1257               myStructure->AddLink(BRepMesh_Edge(iv2, ivf, BRepMesh_Frontier));
1258             else if (orEdge == TopAbs_INTERNAL)
1259               myStructure->AddLink(BRepMesh_Edge(ivf, iv2, BRepMesh_Fixed));
1260
1261             ivf = iv2;
1262           }
1263         }
1264       }
1265
1266     P1 = new Poly_PolygonOnTriangulation(Nodes, PPar);
1267
1268     if (ivf != ivl) {
1269       if (orEdge == TopAbs_FORWARD)
1270         myStructure->AddLink(BRepMesh_Edge(ivf, ivl, BRepMesh_Frontier));
1271       else if (orEdge == TopAbs_REVERSED)
1272         myStructure->AddLink(BRepMesh_Edge(ivl, ivf, BRepMesh_Frontier));
1273       else if (orEdge == TopAbs_INTERNAL)
1274         myStructure->AddLink(BRepMesh_Edge(ivf, ivl, BRepMesh_Fixed));
1275       }
1276     }
1277   }
1278
1279   P1->Deflection(theDefEdge);
1280   if (myInternaledges.IsBound(theEdge))
1281   {
1282     BRepMesh_PairOfPolygon& pair = myInternaledges.ChangeFind(theEdge);
1283     if (orEdge == TopAbs_REVERSED)
1284       pair.Append(P1);
1285     else
1286       pair.Prepend(P1);
1287   }
1288   else
1289   {
1290     BRepMesh_PairOfPolygon pair1;
1291     pair1.Append(P1);
1292     myInternaledges.Bind(theEdge, pair1);
1293   }
1294 }
1295
1296
1297 //=======================================================================
1298 //function : Update(edge)
1299 //purpose  :
1300 //=======================================================================
1301 Standard_Boolean BRepMesh_FastDiscret::Update(const TopoDS_Edge&          theEdge,
1302                                               const TopoDS_Face&          theFace,
1303                                               const Handle(Geom2d_Curve)& theC2d,
1304                                               const Standard_Real         theDefEdge,
1305                                               const Standard_Real         theFirst,
1306                                               const Standard_Real         theLast)
1307 {
1308   TopLoc_Location Loc;
1309   Handle(Poly_Triangulation) T, TNull;
1310   Handle(Poly_PolygonOnTriangulation) Poly, NullPoly;
1311
1312   Standard_Integer i = 1;
1313   Standard_Boolean found = Standard_False;
1314   do
1315   {
1316     BRep_Tool::PolygonOnTriangulation(theEdge,Poly,T,Loc,i);
1317     i++;
1318     if (!found && !T.IsNull() && T->HasUVNodes() &&
1319         !Poly.IsNull() && Poly->HasParameters())
1320     {
1321       if (Poly->Deflection() <= 1.1*theDefEdge)
1322       {
1323         // 2d vertex indices
1324         TopAbs_Orientation orEdge = theEdge.Orientation();
1325         Standard_Integer iv1, iv2, ivl;
1326         Standard_Integer isv1, isv, isvl;
1327
1328         // Get range on 2d curve
1329         Standard_Real wFirst, wLast;
1330         BRep_Tool::Range(theEdge, theFace, wFirst, wLast);
1331
1332         // Get end points on 2d curve
1333         gp_Pnt2d uvFirst, uvLast;
1334         BRep_Tool::UVPoints(theEdge, theFace, uvFirst, uvLast);
1335
1336         // Get vertices
1337         TopoDS_Vertex pBegin, pEnd;
1338         TopExp::Vertices(theEdge,pBegin,pEnd);
1339
1340         const Standard_Boolean sameUV =
1341           uvFirst.IsEqual(uvLast, Precision::PConfusion());
1342
1343         //Controle vertice tolerances
1344         BRepAdaptor_Surface  BS(theFace, Standard_False);
1345         Handle(BRepAdaptor_HSurface) gFace = new BRepAdaptor_HSurface(BS);
1346
1347
1348         gp_Pnt pFirst = gFace->Value(uvFirst.X(), uvFirst.Y());
1349         gp_Pnt pLast  = gFace->Value(uvLast.X(), uvLast.Y());
1350
1351         Standard_Real mindist = 10. * Max(pFirst.Distance(BRep_Tool::Pnt(pBegin)),
1352                                           pLast.Distance(BRep_Tool::Pnt(pEnd)));
1353
1354         if (mindist < BRep_Tool::Tolerance(pBegin) ||
1355             mindist < BRep_Tool::Tolerance(pEnd) ) mindist = theDefEdge;
1356
1357         if (sameUV)
1358         {
1359           // 1. is it really sameUV without being degenerated
1360           gp_Pnt2d uvF, uvL;
1361           theC2d->D0(theFirst, uvF);
1362           theC2d->D0(theLast, uvL);
1363           if (!uvFirst.IsEqual(uvF, Precision::PConfusion())) {
1364             uvFirst = uvF;
1365           }
1366           if (!uvLast.IsEqual(uvL, Precision::PConfusion())) {
1367             uvLast = uvL;
1368           }
1369         }
1370
1371         const TColgp_Array1OfPnt&      Nodes   = T->Nodes();
1372         const TColStd_Array1OfInteger& Indices = Poly->Nodes();
1373         Handle(TColStd_HArray1OfReal)  Param   = Poly->Parameters();
1374
1375         const Standard_Integer nbnodes = Indices.Length();
1376         TColStd_Array1OfInteger NewNodes(1, nbnodes);
1377         TColStd_Array1OfInteger NewNodInStruct(1, nbnodes);
1378
1379         gp_Pnt P3d;
1380         gp_XY theUV;
1381
1382         // Process first vertex
1383         Standard_Integer ipf;
1384         if (myVertices.IsBound(pBegin))
1385         {
1386           ipf = myVertices.Find(pBegin);
1387         }
1388         else
1389         {
1390           if (sameUV && myVertices.IsBound(pEnd))
1391           {
1392               ipf = myVertices.Find(pEnd);
1393           }
1394           else
1395           {
1396             P3d = Nodes(Indices(1));
1397             if (!Loc.IsIdentity())
1398               P3d.Transform(Loc.Transformation());
1399             myNbLocat++;
1400             myLocation3d.Bind(myNbLocat,P3d);
1401             ipf = myNbLocat;
1402           }
1403           myVertices.Bind(pBegin,ipf);
1404         }
1405         NewNodInStruct(1) = ipf;
1406
1407         Handle(BRepMesh_FaceAttribute) aFaceAttribute;
1408         GetFaceAttribute ( theFace, aFaceAttribute );
1409         theUV = BRepMesh_FastDiscretFace::FindUV(pBegin, uvFirst, ipf,
1410           gFace, mindist, aFaceAttribute, myLocation2d);
1411
1412         BRepMesh_Vertex vf(theUV,ipf,BRepMesh_Frontier);
1413         iv1 = myStructure->AddNode(vf);
1414         isv1 = myVemap.FindIndex(iv1);
1415         if (isv1 == 0) isv1 = myVemap.Add(iv1);
1416         NewNodes(1) = isv1;
1417
1418         // Process last vertex
1419         Standard_Integer ipl;
1420         if (pEnd.IsSame(pBegin))
1421         {
1422           ipl = ipf;
1423         }
1424         else
1425         {
1426           if (myVertices.IsBound(pEnd))
1427           {
1428             ipl = myVertices.Find(pEnd);
1429           }
1430           else
1431           {
1432             if (sameUV)
1433             {
1434               ipl = ipf;
1435               ivl = iv1;
1436               isv1 = isv1;
1437             }
1438             else
1439             {
1440               myNbLocat++;
1441               myLocation3d.Bind(myNbLocat,Nodes(Indices(nbnodes)).Transformed(Loc.Transformation()));
1442               ipl = myNbLocat;
1443             }
1444             myVertices.Bind(pEnd,ipl);
1445           }
1446         }
1447         NewNodInStruct(nbnodes) = ipl;
1448         theUV = BRepMesh_FastDiscretFace::FindUV(pEnd, uvLast, ipl,
1449           gFace, mindist, aFaceAttribute, myLocation2d);
1450
1451         BRepMesh_Vertex vl(theUV,ipl,BRepMesh_Frontier);
1452
1453         ivl = myStructure->AddNode(vl);
1454         isvl = myVemap.FindIndex(ivl);
1455         if (isvl == 0) isvl = myVemap.Add(ivl);
1456
1457         NewNodes(nbnodes) = isvl;
1458
1459         gp_Pnt2d uv;
1460         BRepMesh_Vertex v;
1461
1462         if (BRep_Tool::SameParameter(theEdge))
1463         {
1464           for (i = 2; i < Indices.Length(); i++)
1465           {
1466             // Record 3d point
1467             P3d = Nodes(Indices(i));
1468             if (!Loc.IsIdentity())
1469               P3d.Transform(Loc.Transformation());
1470             myNbLocat++;
1471             myLocation3d.Bind(myNbLocat, P3d);
1472             NewNodInStruct(i) = myNbLocat;
1473             // Record 2d point
1474             uv = theC2d->Value(Param->Value(i));
1475             v.Initialize(uv.Coord(), myNbLocat, BRepMesh_Frontier);
1476             iv2 = myStructure->AddNode(v);
1477             isv = myVemap.FindIndex(iv2);
1478             if (isv == 0) isv = myVemap.Add(iv2);
1479             NewNodes(i) = isv;
1480
1481             //add links
1482             if (orEdge == TopAbs_FORWARD)
1483               myStructure->AddLink(BRepMesh_Edge(iv1,iv2,BRepMesh_Frontier));
1484             else if (orEdge == TopAbs_REVERSED)
1485               myStructure->AddLink(BRepMesh_Edge(iv2,iv1,BRepMesh_Frontier));
1486             else if (orEdge == TopAbs_INTERNAL)
1487               myStructure->AddLink(BRepMesh_Edge(iv1,iv2,BRepMesh_Fixed));
1488             iv1 = iv2;
1489           }
1490
1491           // last point
1492           if (iv1 != ivl) {
1493             if (orEdge == TopAbs_FORWARD)
1494               myStructure->AddLink(BRepMesh_Edge(iv1,ivl,BRepMesh_Frontier));
1495             else if (orEdge == TopAbs_REVERSED)
1496               myStructure->AddLink(BRepMesh_Edge(ivl,iv1,BRepMesh_Frontier));
1497             else if (orEdge == TopAbs_INTERNAL)
1498               myStructure->AddLink(BRepMesh_Edge(iv1,ivl,BRepMesh_Fixed));
1499           }
1500
1501
1502         }
1503         else
1504         {
1505           const Standard_Real wFold = Param->Value(Param->Lower());
1506           const Standard_Real wLold = Param->Value(Param->Upper());
1507
1508           Standard_Real wKoef = 1.;
1509           if ((wFold != wFirst || wLold != wLast) && wLold != wFold)
1510           {
1511             wKoef = (wLast - wFirst) / (wLold - wFold);
1512           }
1513
1514           BRepAdaptor_Curve cons(theEdge, theFace);
1515           Extrema_LocateExtPC pcos;
1516           pcos.Initialize(cons, cons.FirstParameter(),
1517                           cons.LastParameter(), Precision::PConfusion());
1518
1519           Standard_Real wPrev;
1520           Standard_Real wCur      = wFirst;
1521           Standard_Real wCurFound = wFirst;
1522           for (i = 2; i < Indices.Length(); i++)
1523           {
1524             // Record 3d point
1525             P3d = Nodes(Indices(i));
1526             if (!Loc.IsIdentity())
1527               P3d.Transform(Loc.Transformation());
1528             myNbLocat++;
1529             myLocation3d.Bind(myNbLocat, P3d);
1530             NewNodInStruct(i) = myNbLocat;
1531             // Record 2d point
1532             wPrev = wCur;
1533             wCur  = wFirst + wKoef*(Param->Value(i) - wFold);
1534             wCurFound += (wCur - wPrev);
1535             pcos.Perform(P3d, wCurFound);
1536             if (pcos.IsDone()) wCurFound = pcos.Point().Parameter();
1537             theC2d->D0(wCurFound, uv);
1538             v.Initialize(uv.Coord(), myNbLocat, BRepMesh_Frontier);
1539             iv2 = myStructure->AddNode(v);
1540             isv = myVemap.FindIndex(iv2);
1541             if (isv == 0) isv = myVemap.Add(iv2);
1542             NewNodes(i) = isv;
1543
1544
1545             //add links
1546             if (orEdge == TopAbs_FORWARD)
1547               myStructure->AddLink(BRepMesh_Edge(iv1,iv2,BRepMesh_Frontier));
1548             else if (orEdge == TopAbs_REVERSED)
1549               myStructure->AddLink(BRepMesh_Edge(iv2,iv1,BRepMesh_Frontier));
1550             else if (orEdge == TopAbs_INTERNAL)
1551               myStructure->AddLink(BRepMesh_Edge(iv1,iv2,BRepMesh_Fixed));
1552             iv1 = iv2;
1553           }
1554
1555           // last point
1556           if (iv1 != ivl) {
1557             if (orEdge == TopAbs_FORWARD)
1558               myStructure->AddLink(BRepMesh_Edge(iv1,ivl,BRepMesh_Frontier));
1559             else if (orEdge == TopAbs_REVERSED)
1560               myStructure->AddLink(BRepMesh_Edge(ivl,iv1,BRepMesh_Frontier));
1561             else if (orEdge == TopAbs_INTERNAL)
1562               myStructure->AddLink(BRepMesh_Edge(iv1,ivl,BRepMesh_Fixed));
1563           }
1564         }
1565
1566         Handle(Poly_PolygonOnTriangulation) P1 =
1567           new Poly_PolygonOnTriangulation(NewNodes, Param->Array1());
1568         P1->Deflection(theDefEdge);
1569         if (myInternaledges.IsBound(theEdge))
1570         {
1571           BRepMesh_PairOfPolygon& pair = myInternaledges.ChangeFind(theEdge);
1572           if (theEdge.Orientation() == TopAbs_REVERSED)
1573             pair.Append(P1);
1574           else
1575             pair.Prepend(P1);
1576         }
1577         else
1578         {
1579           BRepMesh_PairOfPolygon pair1;
1580           pair1.Append(P1);
1581           myInternaledges.Bind(theEdge, pair1);
1582         }
1583
1584         Handle(Poly_PolygonOnTriangulation) P2 =
1585           new Poly_PolygonOnTriangulation(NewNodInStruct, Param->Array1());
1586         P2->Deflection(theDefEdge);
1587         BRepMesh_PairOfPolygon pair;
1588         pair.Append(P2);
1589         myEdges.Bind(theEdge, pair);
1590
1591         found = Standard_True;
1592       }
1593       else
1594       {
1595         BRep_Builder B;
1596         B.UpdateEdge(theEdge,NullPoly,T,Loc);
1597         B.UpdateFace(theFace,TNull);
1598       }
1599     }
1600     else if (!T.IsNull() && !T->HasUVNodes())
1601     {
1602       BRep_Builder B;
1603       B.UpdateEdge(theEdge,NullPoly,T,Loc);
1604       B.UpdateFace(theFace,TNull);
1605     }
1606   }
1607   while (!Poly.IsNull());
1608
1609   return found;
1610 }
1611
1612 //=======================================================================
1613 //function : output
1614 //purpose  :
1615 //=======================================================================
1616 Standard_Integer BRepMesh_FastDiscret::NbTriangles() const
1617 {
1618   return myStructure->NbElements();
1619 }
1620
1621 //=======================================================================
1622 //function : Triangle
1623 //purpose  :
1624 //=======================================================================
1625
1626 const BRepMesh_Triangle& BRepMesh_FastDiscret::Triangle
1627   (const Standard_Integer Index) const
1628 {
1629   return myStructure->GetElement(Index);
1630 }
1631
1632 //=======================================================================
1633 //function : NbEdges
1634 //purpose  :
1635 //=======================================================================
1636
1637 Standard_Integer BRepMesh_FastDiscret::NbEdges() const
1638 {
1639   return myStructure->NbLinks();
1640 }
1641
1642 //=======================================================================
1643 //function : Edge
1644 //purpose  :
1645 //=======================================================================
1646
1647 const BRepMesh_Edge& BRepMesh_FastDiscret::Edge(const Standard_Integer Index) const
1648 {
1649   return myStructure->GetLink(Index);
1650 }
1651
1652 //=======================================================================
1653 //function : NbVertices
1654 //purpose  :
1655 //=======================================================================
1656
1657 Standard_Integer BRepMesh_FastDiscret::NbVertices() const
1658 {
1659   return myStructure->NbNodes();
1660 }
1661
1662 //=======================================================================
1663 //function : Vertex
1664 //purpose  :
1665 //=======================================================================
1666
1667 const BRepMesh_Vertex& BRepMesh_FastDiscret::Vertex
1668   (const Standard_Integer Index) const
1669 {
1670   return myStructure->GetNode(Index);
1671 }
1672
1673 //=======================================================================
1674 //function : Pnt
1675 //purpose  :
1676 //=======================================================================
1677
1678 const gp_Pnt& BRepMesh_FastDiscret::Pnt(const Standard_Integer Index) const
1679 {
1680   return myLocation3d(myStructure->GetNode(Index).Location3d());
1681 }
1682
1683 //=======================================================================
1684 //function : VerticesOfDomain
1685 //purpose  :
1686 //=======================================================================
1687
1688 void BRepMesh_FastDiscret::VerticesOfDomain(BRepMesh_MapOfInteger&  Indices) const
1689 {
1690   Indices.Clear();
1691
1692   // recuperate from the map of edges.
1693   const BRepMesh_MapOfInteger& edmap = myStructure->LinkOfDomain();
1694
1695   // iterator on edges.
1696   BRepMesh_MapOfInteger::Iterator iter(edmap);
1697
1698   Standard_Integer ind_edge;
1699   for (iter.Reset(); iter.More(); iter.Next()) {
1700     ind_edge = iter.Key();
1701     const BRepMesh_Edge& Ed = Edge(ind_edge);
1702     Indices.Add(Ed.FirstNode());
1703     Indices.Add(Ed.LastNode());
1704   }
1705 }
1706
1707 BRepMesh_Status BRepMesh_FastDiscret::CurrentFaceStatus() const
1708 {
1709   return myFacestate;
1710 }
1711
1712 //=======================================================================
1713 //function : GetFaceAttribute
1714 //purpose  :
1715 //=======================================================================
1716
1717 Standard_Boolean BRepMesh_FastDiscret::GetFaceAttribute(const TopoDS_Face&              theFace,
1718                                                         Handle(BRepMesh_FaceAttribute)& theFattrib) const
1719 {
1720   if(myMapattrib.IsBound(theFace))
1721   {
1722     theFattrib = myMapattrib(theFace);
1723     return Standard_True;
1724   }
1725   return Standard_False;
1726 }
1727
1728 //=======================================================================
1729 //function : RemoveFaceAttribute
1730 //purpose  :
1731 //=======================================================================
1732
1733 void BRepMesh_FastDiscret::RemoveFaceAttribute(const TopoDS_Face& theFace)
1734 {
1735   if(myMapattrib.IsBound(theFace))
1736     myMapattrib.UnBind(theFace);
1737 }
1738
1739 //=======================================================================
1740 //function : CreateMutexesForSubShapes
1741 //purpose  :
1742 //=======================================================================
1743 void BRepMesh_FastDiscret::CreateMutexesForSubShapes(const TopoDS_Shape& theShape,
1744                                                      const TopAbs_ShapeEnum theType)
1745 {
1746   myMutexProvider.CreateMutexesForSubShapes(theShape, theType);
1747 }
1748
1749 //=======================================================================
1750 //function : RemoveAllMutexes
1751 //purpose  :
1752 //=======================================================================
1753 void BRepMesh_FastDiscret::RemoveAllMutexes()
1754 {
1755   myMutexProvider.RemoveAllMutexes();
1756 }