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