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