0023106: BRepMesh_IncrementalMesh returns wrong status
[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_EdgeParameterProvider.hxx>
27 #include <BRepMesh_IEdgeTool.hxx>
28 #include <BRepMesh_EdgeTessellator.hxx>
29 #include <BRepMesh_EdgeTessellationExtractor.hxx>
30
31 #include <BRepAdaptor_Curve.hxx>
32 #include <BRepAdaptor_Surface.hxx>
33 #include <BRepAdaptor_HSurface.hxx>
34
35 #include <Bnd_Box.hxx>
36 #include <BRepTools.hxx>
37 #include <BRepBndLib.hxx>
38 #include <BndLib_Add3dCurve.hxx>
39 #include <Poly_Triangulation.hxx>
40 #include <Poly_PolygonOnTriangulation.hxx>
41
42 #include <Precision.hxx>
43 #include <Geom2d_Curve.hxx>
44 #include <Geom_Surface.hxx>
45 #include <Geom_Plane.hxx>
46 #include <GeomAbs_SurfaceType.hxx>
47 #include <Extrema_LocateExtPC.hxx>
48
49 #include <TColStd_Array1OfInteger.hxx>
50 #include <TColStd_HArray1OfReal.hxx>
51 #include <TColgp_Array1OfPnt2d.hxx>
52 #include <TColGeom2d_SequenceOfCurve.hxx>
53 #include <SortTools_ShellSortOfReal.hxx>
54 #include <TCollection_CompareOfReal.hxx>
55
56 #include <TopTools_SequenceOfShape.hxx>
57 #include <TopTools_ListIteratorOfListOfShape.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(
89   const Standard_Real                              theDefle,
90   const Standard_Real                              theAngl,
91   const Bnd_Box&                                   theBox,
92   const Standard_Boolean                           theWithShare,
93   const Standard_Boolean                           theInshape,
94   const Standard_Boolean                           theRelative,
95   const Standard_Boolean                           theShapetrigu,
96   const Standard_Boolean                           isInParallel)
97 : myAngle (theAngl),
98   myDeflection (theDefle),
99   myWithShare (theWithShare),
100   myInParallel (isInParallel),
101   myRelative (theRelative),
102   myShapetrigu (theShapetrigu), 
103   myInshape (theInshape)
104 {
105   if ( myRelative )
106     BRepMesh_ShapeTool::BoxMaxDimension(theBox, myDtotale);
107 }
108
109 //=======================================================================
110 //function : BRepMesh_FastDiscret
111 //purpose  : 
112 //=======================================================================
113 BRepMesh_FastDiscret::BRepMesh_FastDiscret(const TopoDS_Shape&    theShape,
114                                            const Standard_Real    theDefle,
115                                            const Standard_Real    theAngl,
116                                            const Bnd_Box&         theBox,
117                                            const Standard_Boolean theWithShare,
118                                            const Standard_Boolean theInshape,
119                                            const Standard_Boolean theRelative,
120                                            const Standard_Boolean theShapetrigu,
121                                            const Standard_Boolean isInParallel)
122 : myAngle (theAngl),
123   myDeflection (theDefle),
124   myWithShare (theWithShare),
125   myInParallel (isInParallel),
126   myRelative (theRelative),
127   myShapetrigu (theShapetrigu),
128   myInshape (theInshape)
129 {
130   if ( myRelative )
131     BRepMesh_ShapeTool::BoxMaxDimension(theBox, myDtotale);
132
133   Perform(theShape);
134 }
135
136 //=======================================================================
137 //function : InitSharedFaces
138 //purpose  : 
139 //=======================================================================
140 void BRepMesh_FastDiscret::InitSharedFaces(const TopoDS_Shape& theShape)
141 {
142   TopExp::MapShapesAndAncestors(theShape, TopAbs_EDGE, TopAbs_FACE, mySharedFaces);
143 }
144
145 //=======================================================================
146 //function : Perform(shape)
147 //purpose  : 
148 //=======================================================================
149 void BRepMesh_FastDiscret::Perform(const TopoDS_Shape& theShape)
150 {
151   InitSharedFaces(theShape);
152
153   std::vector<TopoDS_Face> aFaces;
154   TopExp_Explorer anExplorer(theShape, TopAbs_FACE);
155   for (; anExplorer.More(); anExplorer.Next())
156   {
157     const TopoDS_Face& aFace = TopoDS::Face(anExplorer.Current());
158     Add(aFace);
159     aFaces.push_back(aFace);
160   }
161
162 #ifdef HAVE_TBB
163   if ( myInParallel )
164   {
165     tbb::parallel_for_each(aFaces.begin(), aFaces.end(), *this);
166   }
167   else
168   {
169 #endif
170     std::vector<TopoDS_Face>::const_iterator anIt(aFaces.begin());
171     for (; anIt != aFaces.end(); anIt++)
172       Process(*anIt);
173 #ifdef HAVE_TBB
174   }
175 #endif
176 }
177
178
179 //=======================================================================
180 //function : Process
181 //purpose  : 
182 //=======================================================================
183 void BRepMesh_FastDiscret::Process(const TopoDS_Face& theFace) const
184 {
185   Handle(BRepMesh_FaceAttribute) anAttribute;
186   if (GetFaceAttribute(theFace, anAttribute))
187   {
188     try
189     {
190       OCC_CATCH_SIGNALS
191
192       BRepMesh_FastDiscretFace aTool(GetAngle(), WithShare());
193       aTool.Add(anAttribute);
194     }
195     catch (Standard_Failure)
196     {
197       anAttribute->SetStatus(BRepMesh_Failure);
198     }
199   }
200 }
201
202 //=======================================================================
203 //function : Add(face)
204 //purpose  : 
205 //=======================================================================
206 Standard_Integer BRepMesh_FastDiscret::Add(const TopoDS_Face& theFace)
207 {
208   try
209   {
210     OCC_CATCH_SIGNALS
211
212     // Initialize face attributes
213     myAttribute.Nullify();
214     GetFaceAttribute(theFace, myAttribute);
215     if (myAttribute.IsNull())
216     {
217       myAttribute = new BRepMesh_FaceAttribute(theFace,
218         myBoundaryVertices, myBoundaryPoints);
219
220       myAttributes.Bind(theFace, myAttribute);
221     }
222
223     myStructure     = myAttribute->ResetStructure();
224     myVertexEdgeMap = myAttribute->ChangeVertexEdgeMap();
225     myInternalEdges = myAttribute->ChangeInternalEdges();
226
227     if (!myWithShare)
228     {
229       myEdges.Clear();
230       myBoundaryVertices.Clear();
231       myBoundaryPoints.Clear();
232     }
233
234     Standard_Real defedge;
235     Standard_Integer nbEdge = 0;
236     Standard_Real savangle = myAngle;
237     Standard_Real cdef;
238     Standard_Real maxdef = 2.* BRepMesh_ShapeTool::MaxFaceTolerance(theFace);
239
240     Standard_Real defface = 0.;
241     if (!myRelative)
242       defface = Max(myDeflection, maxdef);
243
244     N_SEQUENCE<EdgePCurve>  aPCurves;
245     N_SEQUENCE<TopoDS_Edge> aFaceEdges;
246
247     const TopoDS_Face&                  aFace = myAttribute->Face();
248     const Handle(BRepAdaptor_HSurface)& gFace = myAttribute->Surface();
249     TopExp_Explorer aWireIt(aFace, TopAbs_WIRE);
250     for (; aWireIt.More(); aWireIt.Next())
251     {
252       TopExp_Explorer aEdgeIt(aWireIt.Current(), TopAbs_EDGE);
253       for (; aEdgeIt.More(); aEdgeIt.Next(), ++nbEdge)
254       {
255         const TopoDS_Edge& aEdge = TopoDS::Edge(aEdgeIt.Current());
256         if (!myMapdefle.IsBound(aEdge))
257         {
258           if (myRelative)
259           {
260             if (myEdges.IsBound(aEdge))
261             {
262               const BRepMesh_PairOfPolygon& aPair = myEdges.Find(aEdge);
263               const Handle(Poly_PolygonOnTriangulation)& aPolygon = aPair.First();
264               defedge = aPolygon->Deflection();
265             }
266             else
267             {
268               defedge = BRepMesh_ShapeTool::RelativeEdgeDeflection(
269                 aEdge, myDeflection, myDtotale, cdef);
270
271               myAngle = savangle * cdef;
272             }
273
274             defface += defedge;
275             defface = Max(maxdef, defface);
276           }
277           else
278           {
279             defedge = myDeflection;
280           }
281
282           defedge = Max(maxdef, defedge);
283           defedge = Max(UVDEFLECTION, defedge);
284           myMapdefle.Bind(aEdge, defedge);
285         }
286         else
287         {
288           defedge = myMapdefle(aEdge);
289           if ( myRelative )
290           {
291             defface += defedge;
292             defface = Max(maxdef, defface);
293           }
294         }
295
296         Standard_Real aFirstParam, aLastParam;
297         Handle(Geom2d_Curve) aCurve2d = 
298           BRep_Tool::CurveOnSurface(aEdge, aFace, aFirstParam, aLastParam);
299
300         if (aCurve2d.IsNull())
301           continue;
302
303         EdgePCurve aPCurve = { aCurve2d, aFirstParam, aLastParam };
304         aPCurves.Append(aPCurve);
305         aFaceEdges.Append(aEdge);
306
307         add(aEdge, aPCurve, defedge);
308         myAngle = savangle;
309       }
310     }
311
312     if ( nbEdge == 0 || myVertexEdgeMap->Extent() < 3 )
313     {
314       myAttribute->SetStatus(BRepMesh_Failure);
315       return myAttribute->GetStatus();
316     }
317
318     if ( myRelative )
319     {
320       defface = defface / nbEdge;
321     }
322     else
323     {
324       defface = myDeflection;
325     }
326
327     if ( myWithShare )
328       defface = Max(maxdef, defface);
329
330     TopLoc_Location aLoc;
331     Handle(Poly_Triangulation) aTriangulation = BRep_Tool::Triangulation(aFace, aLoc);
332
333     if (!myShapetrigu || aTriangulation.IsNull())
334     {
335       Standard_Real xCur, yCur;
336       Standard_Real maxX, minX, maxY, minY;
337
338       minX = minY = 1.e100;
339       maxX = maxY =-1.e100;
340
341       Standard_Integer ipn = 0;
342       Standard_Integer i1 = 1;
343       for ( i1 = 1; i1 <= myVertexEdgeMap->Extent(); ++i1 )
344       {
345         const BRepMesh_Vertex& aVertex = myStructure->GetNode(myVertexEdgeMap->FindKey(i1));
346
347         ++ipn;
348
349         xCur = aVertex.Coord().X();
350         yCur = aVertex.Coord().Y();
351
352         minX = Min(xCur, minX);
353         maxX = Max(xCur, maxX);
354         minY = Min(yCur, minY);
355         maxY = Max(yCur, maxY);
356       }
357
358       Standard_Real myumin = minX;
359       Standard_Real myumax = maxX;
360       Standard_Real myvmin = minY;
361       Standard_Real myvmax = maxY;
362
363       const Standard_Real umin = gFace->FirstUParameter();
364       const Standard_Real umax = gFace->LastUParameter();
365       const Standard_Real vmin = gFace->FirstVParameter();
366       const Standard_Real vmax = gFace->LastVParameter();
367
368       if (myumin < umin || myumax > umax)
369       {
370         if (gFace->IsUPeriodic())
371         {
372           if ((myumax - myumin) > (umax - umin))
373             myumax = myumin + (umax - umin);
374         }
375         else
376         {
377           if (umin > myumin)
378             myumin = umin;
379
380           if (umax < myumax)
381             myumax = umax;
382         }
383       }
384
385       if (myvmin < vmin || myvmax > vmax)
386       {
387         if (gFace->IsVPeriodic())
388         {
389           if ((myvmax - myvmin) > (vmax - vmin))
390             myvmax = myvmin + (vmax - vmin);
391         }
392         else
393         {
394           if ( vmin > myvmin )
395             myvmin = vmin;
396
397           if (vmax < myvmax)
398             myvmax = vmax;
399         }
400       }
401
402       GeomAbs_SurfaceType aSurfType = gFace->GetType();
403       // Fast verification of the validity of calculated limits.
404       // If wrong, sure a problem of pcurve.
405       if (aSurfType == GeomAbs_BezierSurface &&
406          (myumin < -0.5 || myumax > 1.5 || myvmin < -0.5 || myvmax > 1.5) )
407       {
408         myAttribute->SetStatus(BRepMesh_Failure);
409         return myAttribute->GetStatus();
410       }
411
412       //define parameters for correct parametrics
413       Standard_Real deltaX = 1.0;
414       Standard_Real deltaY = 1.0;
415
416       {
417         BRepMeshCol::HClassifier& aClassifier = myAttribute->ChangeClassifier();
418         BRepMesh_WireChecker aDFaceChecker(aFace, Precision::PConfusion(),
419           myInternalEdges, myVertexEdgeMap, myStructure,
420           myumin, myumax, myvmin, myvmax, myInParallel );
421
422         aDFaceChecker.ReCompute(aClassifier);
423         BRepMesh_Status aCheckStatus = aDFaceChecker.Status();
424
425         if (aCheckStatus == BRepMesh_SelfIntersectingWire)
426         {
427           Standard_Integer nbmaill = 0;
428           Standard_Real eps = Precision::Confusion();
429           while (nbmaill < 5 && aCheckStatus != BRepMesh_ReMesh)
430           {
431             ++nbmaill;
432
433             // Clear the structure of links
434             myStructure = myAttribute->ResetStructure();
435
436             
437             for (Standard_Integer j = 1; j <= aFaceEdges.Length(); ++j)
438             {
439               const TopoDS_Edge& anEdge = aFaceEdges(j);
440               if (myEdges.IsBound(anEdge))
441                 myEdges.UnBind(anEdge);
442
443               defedge = Max(myMapdefle(anEdge) / 3.0, eps);
444               myMapdefle.Bind(anEdge, defedge);
445
446               add(anEdge, aPCurves(j), defedge);
447             }
448
449             aDFaceChecker.ReCompute(aClassifier);
450             if (aDFaceChecker.Status() == BRepMesh_NoError)
451               aCheckStatus = BRepMesh_ReMesh;
452           }
453         }
454
455         myAttribute->SetStatus(aCheckStatus);
456         if (!myAttribute->IsValid())
457           //RemoveFaceAttribute(theFace);
458           return myAttribute->GetStatus();
459       }
460
461       // try to find the real length:
462       // akm (bug OCC16) : We must calculate these measures in non-singular
463       //     parts of face. Let's try to compute average value of three
464       //     (umin, (umin+umax)/2, umax), and respectively for v.
465       //                 vvvvv
466       Standard_Real longu = 0.0, longv = 0.0; //, last , first;
467       gp_Pnt P11, P12, P21, P22, P31, P32;
468
469       Standard_Real du = 0.05 * ( myumax - myumin );
470       Standard_Real dv = 0.05 * ( myvmax - myvmin );
471       Standard_Real dfuave = 0.5 * ( myumin + myumax );
472       Standard_Real dfvave = 0.5 * ( myvmin + myvmax );
473       Standard_Real dfucur, dfvcur;
474
475       // U loop
476       gFace->D0(myumin, myvmin, P11);
477       gFace->D0(myumin, dfvave, P21);
478       gFace->D0(myumin, myvmax, P31);
479       for (i1=1, dfucur=myumin+du; i1 <= 20; i1++, dfucur+=du)
480       {
481         gFace->D0(dfucur, myvmin, P12);
482         gFace->D0(dfucur, dfvave, P22);
483         gFace->D0(dfucur, myvmax, P32);
484         longu += ( P11.Distance(P12) + P21.Distance(P22) + P31.Distance(P32) );
485         P11 = P12;
486         P21 = P22;
487         P31 = P32;
488       }
489
490       // V loop
491       gFace->D0(myumin, myvmin, P11);
492       gFace->D0(dfuave, myvmin, P21);
493       gFace->D0(myumax, myvmin, P31);
494       for (i1=1, dfvcur=myvmin+dv; i1 <= 20; i1++, dfvcur+=dv)
495       {
496         gFace->D0(myumin, dfvcur, P12);
497         gFace->D0(dfuave, dfvcur, P22);
498         gFace->D0(myumax, dfvcur, P32);
499         longv += ( P11.Distance(P12) + P21.Distance(P22) + P31.Distance(P32) );
500         P11 = P12;
501         P21 = P22;
502         P31 = P32;
503       }
504
505       longu /= 3.;
506       longv /= 3.;
507       // akm (bug OCC16) ^^^^^
508
509       if (longu <= 1.e-16 || longv <= 1.e-16)
510       {
511         //yes, it is seen!!
512         myAttribute->SetStatus(BRepMesh_Failure);
513         return myAttribute->GetStatus();
514       }
515
516
517       if (aSurfType == GeomAbs_Torus)
518       {
519         gp_Torus Tor = gFace->Torus();
520         Standard_Real r = Tor.MinorRadius(), R = Tor.MajorRadius();
521         Standard_Real Du, Dv;//, pasu, pasv;
522
523         Dv = Max(1.0e0 - (defface/r),0.0e0) ;
524         Standard_Real oldDv = 2.0 * ACos (Dv);
525         oldDv = Min(oldDv, myAngle);
526         Dv  =  0.9*oldDv; //TWOTHIRD * oldDv;
527         Dv = oldDv;
528
529         Standard_Integer nbV = Max((Standard_Integer)((myvmax-myvmin)/Dv), 2);
530         Dv = (myvmax-myvmin)/(nbV+1);
531
532         Standard_Real ru = R + r;
533         if ( ru > 1.e-16 )
534         {
535           Du = Max(1.0e0 - (defface/ru),0.0e0);
536           Du  = (2.0 * ACos (Du));
537           Du = Min(Du, myAngle);
538           Standard_Real aa = sqrt(Du*Du + oldDv*oldDv);
539
540           if (aa < gp::Resolution())
541             return myAttribute->GetStatus();
542
543           Du = Du * Min(oldDv, Du) / aa;
544         }
545         else
546         {
547           Du = Dv;
548         }
549
550         Standard_Integer nbU = Max((Standard_Integer)((myumax-myumin)/Du), 2);
551         nbU = Max(nbU, (Standard_Integer)(nbV*(myumax-myumin)*R/((myvmax-myvmin)*r)/5.));
552       
553         Du = (myumax-myumin)/(nbU+1);
554         //-- DeltaX and DeltaY are chosen so that to avoid "jumping" 
555         //-- of points on the grid
556         deltaX = Du;
557         deltaY = Dv;
558       }
559       else if (aSurfType == GeomAbs_Cylinder)
560       {
561         gp_Cylinder Cyl = gFace->Cylinder();
562         Standard_Real R = Cyl.Radius();
563
564         // Calculate parameters for iteration in U direction
565         Standard_Real Du = 1.0 - (defface/R);
566         if (Du < 0.0)
567           Du = 0.0;
568
569         Du = 2.0 * ACos (Du);
570         if (Du > GetAngle())
571           Du = GetAngle();
572
573         deltaX = Du / longv;
574         deltaY = 1.;
575       }
576       else
577       {
578         deltaX = (myumax-myumin)/longu;
579         deltaY = (myvmax-myvmin)/longv;
580       }
581
582       // Restore face attribute
583       myAttribute->SetDefFace(defface);
584       myAttribute->SetUMax(myumax);
585       myAttribute->SetVMax(myvmax);
586       myAttribute->SetUMin(myumin);
587       myAttribute->SetVMin(myvmin);
588       myAttribute->SetDeltaX(deltaX);
589       myAttribute->SetDeltaY(deltaY);
590     }
591   }
592   catch(Standard_Failure)
593   {
594     myAttribute->SetStatus(BRepMesh_Failure);
595   }
596
597   return myAttribute->GetStatus();
598 }
599
600 //=======================================================================
601 //function : addLinkToMesh
602 //purpose  :
603 //=======================================================================
604 void BRepMesh_FastDiscret::addLinkToMesh(
605   const Standard_Integer   theFirstNodeId,
606   const Standard_Integer   theLastNodeId,
607   const TopAbs_Orientation theOrientation)
608 {
609   if (theOrientation == TopAbs_FORWARD)
610     myStructure->AddLink(BRepMesh_Edge(theFirstNodeId, theLastNodeId, BRepMesh_Frontier));
611   else if (theOrientation == TopAbs_REVERSED)
612     myStructure->AddLink(BRepMesh_Edge(theLastNodeId, theFirstNodeId, BRepMesh_Frontier));
613   else if (theOrientation == TopAbs_INTERNAL)
614     myStructure->AddLink(BRepMesh_Edge(theFirstNodeId, theLastNodeId, BRepMesh_Fixed));
615 }
616
617 //=======================================================================
618 //function : getEdgeAttributes
619 //purpose  :
620 //=======================================================================
621 Standard_Boolean BRepMesh_FastDiscret::getEdgeAttributes(
622   const TopoDS_Edge&                      theEdge,
623   const BRepMesh_FastDiscret::EdgePCurve& thePCurve,
624   const Standard_Real                     theDefEdge,
625   BRepMesh_FastDiscret::EdgeAttributes&   theAttributes) const
626 {
627   EdgeAttributes& aEAttr = theAttributes;
628
629   // Get vertices
630   TopExp::Vertices(theEdge, aEAttr.FirstVertex, aEAttr.LastVertex);
631   if (aEAttr.FirstVertex.IsNull() || aEAttr.LastVertex.IsNull())
632     return Standard_False;
633
634   // Get range on 2d curve
635   const TopoDS_Face& aFace = myAttribute->Face();
636   BRep_Tool::Range(theEdge, aFace, aEAttr.FirstParam, aEAttr.LastParam);
637
638   // Get end points on 2d curve
639   BRep_Tool::UVPoints(theEdge, aFace, aEAttr.FirstUV, aEAttr.LastUV);
640
641   aEAttr.IsSameUV =
642     aEAttr.FirstUV.IsEqual(aEAttr.LastUV, Precision::PConfusion());
643
644   //Control tolerance of vertices
645   const Handle(BRepAdaptor_HSurface)& gFace = myAttribute->Surface();
646   gp_Pnt pFirst = gFace->Value(aEAttr.FirstUV.X(), aEAttr.FirstUV.Y());
647   gp_Pnt pLast  = gFace->Value(aEAttr.LastUV.X(),  aEAttr.LastUV.Y());
648
649   aEAttr.MinDist = 10. * Max(pFirst.Distance(BRep_Tool::Pnt(aEAttr.FirstVertex)),
650                              pLast .Distance(BRep_Tool::Pnt(aEAttr.LastVertex)));
651
652   if (aEAttr.MinDist < BRep_Tool::Tolerance(aEAttr.FirstVertex) ||
653       aEAttr.MinDist < BRep_Tool::Tolerance(aEAttr.LastVertex))
654   {
655     aEAttr.MinDist = theDefEdge;
656   }
657
658   if (aEAttr.IsSameUV)
659   {
660     // 1. is it really sameUV without being degenerated
661     gp_Pnt2d uvF, uvL;
662     thePCurve.Curve2d->D0(thePCurve.FirstParam, uvF);
663     thePCurve.Curve2d->D0(thePCurve.LastParam,  uvL);
664
665     if (!aEAttr.FirstUV.IsEqual(uvF, Precision::PConfusion()))
666       aEAttr.FirstUV = uvF;
667
668     if (!aEAttr.LastUV.IsEqual(uvL, Precision::PConfusion()))
669       aEAttr.LastUV = uvL;
670   }
671
672   return Standard_True;
673 }
674
675 //=======================================================================
676 //function : registerEdgeVertices
677 //purpose  : 
678 //=======================================================================
679 void BRepMesh_FastDiscret::registerEdgeVertices(
680   BRepMesh_FastDiscret::EdgeAttributes& theAttributes,
681   Standard_Integer&                     ipf,
682   Standard_Integer&                     ivf,
683   Standard_Integer&                     isvf,
684   Standard_Integer&                     ipl,
685   Standard_Integer&                     ivl,
686   Standard_Integer&                     isvl)
687 {
688   EdgeAttributes& aEAttr = theAttributes;
689   if (aEAttr.FirstVExtractor.IsNull())
690   {
691     // Use edge geometry to produce tesselation.
692     aEAttr.FirstVExtractor = 
693       new TopoDSVExplorer(aEAttr.FirstVertex, aEAttr.IsSameUV, aEAttr.LastVertex);
694   }
695
696   if (aEAttr.LastVExtractor.IsNull())
697   {
698     // Use edge geometry to produce tesselation.
699     aEAttr.LastVExtractor = 
700       new TopoDSVExplorer(aEAttr.LastVertex, aEAttr.IsSameUV, aEAttr.FirstVertex);
701   }
702
703   gp_XY aTmpUV;
704   // Process first vertex
705   ipf = myAttribute->GetVertexIndex(aEAttr.FirstVExtractor, Standard_True);
706   aTmpUV = BRepMesh_ShapeTool::FindUV(ipf, aEAttr.FirstUV, aEAttr.FirstVertex, 
707     aEAttr.MinDist, myAttribute);
708
709   myAttribute->AddNode(ipf, aTmpUV, BRepMesh_Frontier, ivf, isvf);
710
711   // Process last vertex
712   ipl = aEAttr.LastVertex.IsSame(aEAttr.FirstVertex) ? ipf :
713     myAttribute->GetVertexIndex(aEAttr.LastVExtractor, Standard_True);
714   aTmpUV = BRepMesh_ShapeTool::FindUV(ipl, aEAttr.LastUV, aEAttr.LastVertex, 
715     aEAttr.MinDist, myAttribute);
716
717   myAttribute->AddNode(ipl, aTmpUV, BRepMesh_Frontier, ivl, isvl);
718 }
719
720 //=======================================================================
721 //function : add
722 //purpose  : 
723 //=======================================================================
724 void BRepMesh_FastDiscret::add(
725   const TopoDS_Edge&                      theEdge,
726   const BRepMesh_FastDiscret::EdgePCurve& thePCurve,
727   const Standard_Real                     theDefEdge)
728 {
729   const TopAbs_Orientation orEdge = theEdge.Orientation();
730   if (orEdge == TopAbs_EXTERNAL)
731     return;
732
733   EdgeAttributes aEAttr;
734   if (!getEdgeAttributes(theEdge, thePCurve, theDefEdge, aEAttr))
735     return;
736
737   if (!myEdges.IsBound(theEdge))
738   {
739     update(theEdge, thePCurve.Curve2d, theDefEdge, aEAttr);
740     return;
741   }
742
743   Standard_Integer ipf, ivf, isvf, ipl, ivl, isvl;
744   registerEdgeVertices(aEAttr, ipf, ivf, isvf, ipl, ivl, isvl);
745
746   // If this Edge has been already checked and it is not degenerated, 
747   // the points of the polygon calculated at the first check are retrieved :
748
749   // retrieve the polygone:
750   const BRepMesh_PairOfPolygon&              aPair    = myEdges.Find(theEdge);
751   const Handle(Poly_PolygonOnTriangulation)& aPolygon = aPair.First();
752   const TColStd_Array1OfInteger&             aNodes   = aPolygon->Nodes();
753   Handle(TColStd_HArray1OfReal)              aParams  = aPolygon->Parameters();
754
755   // creation anew:
756   const Standard_Integer  aNodesNb = aNodes.Length();
757   TColStd_Array1OfInteger aNewNodes (1, aNodesNb);
758   TColStd_Array1OfReal    aNewParams(1, aNodesNb);
759
760   aNewNodes (1) = isvf;
761   aNewParams(1) = aEAttr.FirstParam;
762
763   aNewNodes (aNodesNb) = isvl;
764   aNewParams(aNodesNb) = aEAttr.LastParam;
765
766   const TopoDS_Face& aFace = myAttribute->Face();
767   if (!BRepMesh_ShapeTool::IsDegenerated(theEdge, aFace))
768   {
769     BRepMesh_EdgeParameterProvider aProvider(theEdge, aFace, aParams);
770     for (Standard_Integer i = 2; i < aNodesNb; ++i)
771     {
772       const Standard_Integer aPointId = aNodes(i);
773       gp_Pnt& aPnt = myBoundaryPoints(aPointId);
774
775       const Standard_Real aParam = aProvider.Parameter(i, aPnt);
776       gp_Pnt2d aUV = thePCurve.Curve2d->Value(aParam);
777
778       Standard_Integer iv2, isv;
779       myAttribute->AddNode(aPointId, aUV.Coord(), BRepMesh_OnCurve, iv2, isv);
780
781       aNewNodes (i) = isv;
782       aNewParams(i) = aParam;
783
784       addLinkToMesh(ivf, iv2, orEdge);
785       ivf = iv2;
786     }
787
788     // last point
789     if (ivf != ivl)
790       addLinkToMesh(ivf, ivl, orEdge);
791   }
792
793   Handle(Poly_PolygonOnTriangulation) P1 = 
794     new Poly_PolygonOnTriangulation(aNewNodes, aNewParams);
795
796   storePolygon(theEdge, P1, theDefEdge);
797 }
798
799 //=======================================================================
800 //function : update(edge)
801 //purpose  :
802 //=======================================================================
803 void BRepMesh_FastDiscret::update(
804   const TopoDS_Edge&                                theEdge,
805   const Handle(Geom2d_Curve)&                       theC2d,
806   const Standard_Real                               theDefEdge,
807   BRepMesh_FastDiscret::EdgeAttributes&             theAttributes)
808 {
809   EdgeAttributes& aEAttr = theAttributes;
810
811   const TopoDS_Face& aFace = myAttribute->Face();
812   Handle(BRepMesh_IEdgeTool) aEdgeTool;
813   // Try to find existing tessellation.
814   for (Standard_Integer i = 1; aEdgeTool.IsNull(); ++i)
815   {
816     TopLoc_Location aLoc;
817     Handle(Poly_Triangulation) aTriangulation;
818     Handle(Poly_PolygonOnTriangulation) aPolygon;
819     BRep_Tool::PolygonOnTriangulation(theEdge, aPolygon, aTriangulation, aLoc, i);
820
821     if (aPolygon.IsNull())
822       break;
823
824     if (aTriangulation.IsNull() || !aPolygon->HasParameters())
825       continue;
826
827     if (aPolygon->Deflection() > 1.1 * theDefEdge)
828       continue;
829
830     const TColgp_Array1OfPnt&      aNodes   = aTriangulation->Nodes();
831     const TColStd_Array1OfInteger& aIndices = aPolygon->Nodes();
832     Handle(TColStd_HArray1OfReal)  aParams  = aPolygon->Parameters();
833
834     aEAttr.FirstVExtractor = new PolyVExplorer(aEAttr.FirstVertex, 
835       aEAttr.IsSameUV, aEAttr.LastVertex, aIndices(1), aNodes, aLoc);
836
837     aEAttr.LastVExtractor = new PolyVExplorer(aEAttr.LastVertex, 
838       aEAttr.IsSameUV, aEAttr.FirstVertex, aIndices(aIndices.Length()), aNodes, aLoc);
839
840     aEdgeTool = new BRepMesh_EdgeTessellationExtractor(theEdge, theC2d, 
841       aFace, aTriangulation, aPolygon, aLoc);
842   }
843
844   if (aEdgeTool.IsNull())
845   {
846     aEdgeTool = new BRepMesh_EdgeTessellator(theEdge, myAttribute, 
847       mySharedFaces, theDefEdge, myAngle);
848   }
849
850   Standard_Integer ipf, ivf, isvf, ipl, ivl, isvl;
851   registerEdgeVertices(aEAttr, ipf, ivf, isvf, ipl, ivl, isvl);
852
853   TopAbs_Orientation orEdge = theEdge.Orientation();
854   Handle(Poly_PolygonOnTriangulation) P1, P2;
855   if (BRepMesh_ShapeTool::IsDegenerated(theEdge, aFace))
856   {
857     const Standard_Integer  aNodesNb = 2;
858     TColStd_Array1OfInteger aNewNodes      (1, aNodesNb);
859     TColStd_Array1OfInteger aNewNodInStruct(1, aNodesNb);
860     TColStd_Array1OfReal    aNewParams     (1, aNodesNb);
861
862     aNewNodInStruct(1) = ipf;
863     aNewNodes      (1) = isvf;
864     aNewParams     (1) = aEAttr.FirstParam;
865
866     aNewNodInStruct(aNodesNb) = ipl;
867     aNewNodes      (aNodesNb) = isvl;
868     aNewParams     (aNodesNb) = aEAttr.LastParam;
869
870     P1 = new Poly_PolygonOnTriangulation(aNewNodes,       aNewParams);
871     P2 = new Poly_PolygonOnTriangulation(aNewNodInStruct, aNewParams);
872   }
873   else
874   {
875     const Standard_Integer  aNodesNb = aEdgeTool->NbPoints();
876     TColStd_Array1OfInteger aNewNodes      (1, aNodesNb);
877     TColStd_Array1OfInteger aNewNodInStruct(1, aNodesNb);
878     TColStd_Array1OfReal    aNewParams     (1, aNodesNb);
879
880     aNewNodInStruct(1) = ipf;
881     aNewNodes      (1) = isvf;
882     aNewParams     (1) = aEAttr.FirstParam;
883
884     aNewNodInStruct(aNodesNb) = ipl;
885     aNewNodes      (aNodesNb) = isvl;
886     aNewParams     (aNodesNb) = aEAttr.LastParam;
887
888     Standard_Integer aLastPointId = myAttribute->LastPointId();
889     for (Standard_Integer i = 2; i < aNodesNb; ++i)
890     {
891       gp_Pnt        aPnt;
892       gp_Pnt2d      aUV;
893       Standard_Real aParam;
894       aEdgeTool->Value(i, aParam, aPnt, aUV);
895       myBoundaryPoints.Bind(++aLastPointId, aPnt);
896
897       Standard_Integer iv2, isv;
898       myAttribute->AddNode(aLastPointId, aUV.Coord(), BRepMesh_Frontier, iv2, isv);
899
900       aNewNodInStruct(i) = aLastPointId;
901       aNewNodes      (i) = isv;
902       aNewParams     (i) = aParam;
903
904       addLinkToMesh(ivf, iv2, orEdge);
905       ivf = iv2;
906     }
907
908     P1 = new Poly_PolygonOnTriangulation(aNewNodes,       aNewParams);
909     P2 = new Poly_PolygonOnTriangulation(aNewNodInStruct, aNewParams);
910   }
911
912   // last point
913   if (ivf != ivl)
914     addLinkToMesh(ivf, ivl, orEdge);
915
916   storePolygon(theEdge, P1, theDefEdge);
917   storePolygonSharedData(theEdge, P2, theDefEdge);
918 }
919
920 //=======================================================================
921 //function : storeInternalPolygon
922 //purpose  : 
923 //=======================================================================
924 void BRepMesh_FastDiscret::storePolygon(
925   const TopoDS_Edge&                         theEdge,
926   Handle(Poly_PolygonOnTriangulation)&       thePolygon,
927   const Standard_Real                        theDeflection)
928 {
929   thePolygon->Deflection(theDeflection);
930   if (myInternalEdges->IsBound(theEdge))
931   {
932     BRepMesh_PairOfPolygon& aPair = myInternalEdges->ChangeFind(theEdge);
933     if (theEdge.Orientation() == TopAbs_REVERSED)
934       aPair.Append(thePolygon);
935     else
936       aPair.Prepend(thePolygon);
937
938     return;
939   }
940
941   BRepMesh_PairOfPolygon aPair;
942   aPair.Append(thePolygon);
943   myInternalEdges->Bind(theEdge, aPair);
944 }
945
946 //=======================================================================
947 //function : storePolygonSharedData
948 //purpose  : 
949 //=======================================================================
950 void BRepMesh_FastDiscret::storePolygonSharedData(
951   const TopoDS_Edge&                         theEdge,
952   Handle(Poly_PolygonOnTriangulation)&       thePolygon,
953   const Standard_Real                        theDeflection)
954 {
955   thePolygon->Deflection(theDeflection);
956   BRepMesh_PairOfPolygon aPair;
957   aPair.Append(thePolygon);
958   myEdges.Bind(theEdge, aPair);
959 }
960
961 //=======================================================================
962 //function : GetFaceAttribute
963 //purpose  : 
964 //=======================================================================
965 Standard_Boolean BRepMesh_FastDiscret::GetFaceAttribute(
966   const TopoDS_Face&              theFace,
967   Handle(BRepMesh_FaceAttribute)& theAttribute ) const
968 {
969   if (myAttributes.IsBound(theFace))
970   {
971     theAttribute = myAttributes(theFace);
972     return Standard_True;
973   }
974
975   return Standard_False;
976 }
977
978 //=======================================================================
979 //function : RemoveFaceAttribute
980 //purpose  : 
981 //=======================================================================
982 void BRepMesh_FastDiscret::RemoveFaceAttribute(const TopoDS_Face& theFace)
983 {
984   if (myAttributes.IsBound(theFace))
985     myAttributes.UnBind(theFace);
986 }