1 // Created by: Ekaterina SMIRNOVA
2 // Copyright (c) 2008-2014 OPEN CASCADE SAS
4 // This file is part of Open CASCADE Technology software library.
6 // This library is free software; you can redistribute it and/or modify it under
7 // the terms of the GNU Lesser General Public License version 2.1 as published
8 // by the Free Software Foundation, with special exception defined in the file
9 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
10 // distribution for complete text of the license and disclaimer of any warranty.
12 // Alternatively, this file may be used under the terms of Open CASCADE
13 // commercial license or contractual agreement.
15 #include <BRepMesh_FastDiscretFace.hxx>
17 #include <BRepMesh_PairOfPolygon.hxx>
18 #include <BRepMesh_ShapeTool.hxx>
19 #include <Poly_PolygonOnTriangulation.hxx>
20 #include <Poly_Triangulation.hxx>
22 #include <BRepAdaptor_Surface.hxx>
23 #include <BRepAdaptor_HSurface.hxx>
24 #include <BRepAdaptor_Curve.hxx>
25 #include <Adaptor3d_IsoCurve.hxx>
27 #include <BRep_ListIteratorOfListOfPointRepresentation.hxx>
28 #include <BRep_PointRepresentation.hxx>
29 #include <BRep_TVertex.hxx>
30 #include <BRep_Tool.hxx>
32 #include <GeomLib.hxx>
33 #include <Geom_Surface.hxx>
34 #include <Geom_BSplineSurface.hxx>
35 #include <GCPnts_TangentialDeflection.hxx>
37 #include <Standard_ErrorHandler.hxx>
38 #include <Standard_Failure.hxx>
39 #include <SortTools_ShellSortOfReal.hxx>
40 #include <TCollection_CompareOfReal.hxx>
41 #include <TColStd_Array1OfReal.hxx>
42 #include <TColStd_ListOfInteger.hxx>
43 #include <TColStd_SequenceOfReal.hxx>
44 #include <TColStd_Array1OfInteger.hxx>
45 #include <TColStd_HArray1OfReal.hxx>
46 #include <TColgp_Array1OfPnt2d.hxx>
47 #include <TopTools_DataMapOfShapeReal.hxx>
49 #include <TopExp_Explorer.hxx>
51 #include <TopoDS_Vertex.hxx>
54 #include <NCollection_Map.hxx>
55 #include <Bnd_Box2d.hxx>
59 IMPLEMENT_STANDARD_HANDLE (BRepMesh_FastDiscretFace, Standard_Transient)
60 IMPLEMENT_STANDARD_RTTIEXT(BRepMesh_FastDiscretFace, Standard_Transient)
62 static Standard_Real FUN_CalcAverageDUV(TColStd_Array1OfReal& P, const Standard_Integer PLen)
64 Standard_Integer i, j, n = 0;
65 Standard_Real p, result = 0.;
67 for(i = 1; i <= PLen; i++)
70 for(j = i + 1; j <= PLen; j++)
82 p = Abs(P(i) - P(i-1));
90 return (n? (result / (Standard_Real) n) : -1.);
95 //! Auxiliary class used to extract geometrical parameters of fixed TopoDS_Vertex.
100 DEFINE_STANDARD_ALLOC
102 FixedVExplorer(const TopoDS_Vertex& theVertex)
103 : myVertex(theVertex)
107 const TopoDS_Vertex& Vertex() const
112 Standard_Boolean IsSameUV() const
114 return Standard_False;
117 TopoDS_Vertex SameVertex() const
119 return TopoDS_Vertex();
122 virtual gp_Pnt Point() const
124 return BRep_Tool::Pnt(myVertex);
129 void operator =(const FixedVExplorer& /*theOther*/)
134 const TopoDS_Vertex& myVertex;
139 //=======================================================================
140 //function : BRepMesh_FastDiscretFace
142 //=======================================================================
143 BRepMesh_FastDiscretFace::BRepMesh_FastDiscretFace
144 (const Standard_Real theAngle,
145 const Standard_Boolean theWithShare) :
146 myAngle(theAngle), myWithShare(theWithShare),
147 myInternalVerticesMode(Standard_True)
149 myAllocator = new NCollection_IncAllocator(64000);
152 //=======================================================================
153 //function : Add(face)
155 //=======================================================================
156 void BRepMesh_FastDiscretFace::Add(const Handle(BRepMesh_FaceAttribute)& theAttribute)
158 if (!theAttribute->IsValid())
161 myAttribute = theAttribute;
162 myStructure = myAttribute->ChangeStructure();
163 myVertexEdgeMap = myAttribute->ChangeVertexEdgeMap();
164 myClassifier = myAttribute->ChangeClassifier();
165 mySurfacePoints = myAttribute->ChangeSurfacePoints();
167 Standard_Real aTolU = myAttribute->ToleranceU();
168 Standard_Real aTolV = myAttribute->ToleranceV();
169 Standard_Real uCellSize = 14.0 * aTolU;
170 Standard_Real vCellSize = 14.0 * aTolV;
172 const Handle(BRepAdaptor_HSurface)& gFace = myAttribute->Surface();
173 GeomAbs_SurfaceType thetype = gFace->GetType();
175 Standard_Integer i = 1;
177 ////////////////////////////////////////////////////////////
178 //add internal vertices after self-intersection check
179 Standard_Integer nbVertices = 0;
180 if ( myInternalVerticesMode )
182 TopExp_Explorer anExplorer(myAttribute->Face(), TopAbs_VERTEX, TopAbs_EDGE);
183 for ( ; anExplorer.More(); anExplorer.Next() )
184 add(TopoDS::Vertex(anExplorer.Current()));
185 nbVertices = myVertexEdgeMap->Extent();
188 // essai de determination de la longueur vraie:
189 // akm (bug OCC16) : We must calculate these measures in non-singular
190 // parts of face. Let`s try to compute average value of three
191 // (umin, (umin+umax)/2, umax), and respectively for v.
193 //Standard_Real longu = 0.0, longv = 0.0; //, last , first;
194 //gp_Pnt P11, P12, P21, P22, P31, P32;
196 const Standard_Real deltaX = myAttribute->GetDeltaX();
197 const Standard_Real deltaY = myAttribute->GetDeltaY();
199 BRepMeshCol::Array1OfInteger tabvert_corr(1, nbVertices);
201 // Check the necessity to fill the map of parameters
202 const Standard_Boolean useUVParam = (thetype == GeomAbs_Torus ||
203 thetype == GeomAbs_BezierSurface ||
204 thetype == GeomAbs_BSplineSurface);
208 BRepMesh_VertexTool aMoveNodes(nbVertices, myAllocator);
210 aMoveNodes.SetCellSize ( uCellSize / deltaX, vCellSize / deltaY);
211 aMoveNodes.SetTolerance( aTolU / deltaX, aTolV / deltaY);
213 for ( i = 1; i <= myStructure->NbNodes(); ++i )
215 const BRepMesh_Vertex& v = myStructure->GetNode(i);
216 gp_XY aPnt2d = v.Coord();
220 myUParam.Add(aPnt2d.X());
221 myVParam.Add(aPnt2d.Y());
224 aPnt2d = myAttribute->Scale(aPnt2d, Standard_True);
225 BRepMesh_Vertex v_new(aPnt2d, v.Location3d(), v.Movability());
226 const BRepMeshCol::ListOfInteger& alist = myStructure->LinksConnectedTo(i);
227 aMoveNodes.Add(v_new, alist);
230 myStructure->ReplaceNodes(aMoveNodes);
232 Standard_Boolean rajout =
233 (thetype == GeomAbs_Sphere || thetype == GeomAbs_Torus);
235 BRepMesh_Delaun trigu(myStructure, tabvert_corr);
237 //removed all free edges from triangulation
238 const Standard_Integer nbLinks = myStructure->NbLinks();
239 for( i = 1; i <= nbLinks; i++ )
241 if( myStructure->ElementsConnectedTo(i).Extent() < 1 )
243 BRepMesh_Edge& anEdge = (BRepMesh_Edge&)trigu.GetEdge(i);
244 if ( anEdge.Movability() == BRepMesh_Deleted )
247 anEdge.SetMovability(BRepMesh_Free);
248 myStructure->RemoveLink(i);
252 const Standard_Real umax = myAttribute->GetUMax();
253 const Standard_Real umin = myAttribute->GetUMin();
254 const Standard_Real vmax = myAttribute->GetVMax();
255 const Standard_Real vmin = myAttribute->GetVMin();
257 Standard_Boolean isaline =
258 ((umax - umin) < Precision::PConfusion() ||
259 (vmax - vmin) < Precision::PConfusion());
261 Standard_Real aDef = -1;
262 if ( !isaline && myStructure->ElementsOfDomain().Extent() > 0 )
264 BRepMeshCol::ListOfVertex aNewVertices;
267 aDef = control(aNewVertices, trigu, Standard_True);
269 if( aDef > myAttribute->GetDefFace() || aDef < 0.)
270 rajout = Standard_True;
273 if ( !rajout && useUVParam )
275 if ( myVParam.Extent() > 2 &&
276 (gFace->IsUClosed() ||
279 rajout = Standard_True;
285 insertInternalVertices(aNewVertices, trigu);
287 //control internal points
288 aDef = control(aNewVertices, trigu, Standard_False);
292 //modify myStructure back
293 aMoveNodes.SetCellSize ( uCellSize, vCellSize );
294 aMoveNodes.SetTolerance( aTolU , aTolV );
296 for ( i = 1; i <= myStructure->NbNodes(); ++i )
298 const BRepMesh_Vertex& v = myStructure->GetNode(i);
299 gp_XY aPnt2d = myAttribute->Scale(v.Coord(), Standard_False);
300 BRepMesh_Vertex v_new(aPnt2d, v.Location3d(), v.Movability());
301 const BRepMeshCol::ListOfInteger& alist = myStructure->LinksConnectedTo(i);
302 aMoveNodes.Add(v_new, alist);
304 // Register internal nodes used in triangulation
305 if (!alist.IsEmpty() && myVertexEdgeMap->FindIndex(i) == 0)
306 myVertexEdgeMap->Add(i);
308 myStructure->ReplaceNodes(aMoveNodes);
311 myAttribute->SetDefFace(aDef);
314 //=======================================================================
315 //function : addVerticesToMesh
317 //=======================================================================
318 Standard_Boolean BRepMesh_FastDiscretFace::addVerticesToMesh(
319 const BRepMeshCol::ListOfVertex& theVertices,
320 BRepMesh_Delaun& theMeshBuilder)
322 if (theVertices.IsEmpty())
323 return Standard_False;
325 BRepMeshCol::Array1OfVertexOfDelaun aArrayOfNewVertices(1, theVertices.Extent());
326 BRepMeshCol::ListOfVertex::Iterator aVertexIt(theVertices);
327 for (Standard_Integer aVertexId = 0; aVertexIt.More(); aVertexIt.Next())
328 aArrayOfNewVertices(++aVertexId) = aVertexIt.Value();
330 theMeshBuilder.AddVertices(aArrayOfNewVertices);
331 return Standard_True;
334 //=======================================================================
335 //function : insertInternalVertices
337 //=======================================================================
338 static void filterParameters(const BRepMeshCol::IMapOfReal& theParams,
339 const Standard_Real theMinDist,
340 const Standard_Real theFilterDist,
341 BRepMeshCol::SequenceOfReal& theResult)
343 // Sort sequence of parameters
344 BRepMeshCol::SequenceOfReal aParamTmp;
345 Standard_Integer aParamLength = 1;
346 const Standard_Integer anInitLen = theParams.Extent();
348 TColStd_Array1OfReal aParamArray(1, anInitLen);
350 for (j = 1; j <= anInitLen; j++)
351 aParamArray(j) = theParams(j);
353 TCollection_CompareOfReal aCompare;
354 SortTools_ShellSortOfReal::Sort(aParamArray, aCompare);
356 // mandadory pre-filtering using the first (minimal) filter value
357 Standard_Real aP1, aP2;
358 aP1 = aParamArray(1);
359 aParamTmp.Append(aP1);
360 for (j = 2; j <= anInitLen; j++)
362 aP2 = aParamArray(j);
363 if ((aP2-aP1) > theMinDist)
365 aParamTmp.Append(aP2);
371 //add last point if required
372 if(aParamArray(anInitLen)-theParams(aParamLength) > theMinDist)
374 aParamTmp.Append(aParamArray(anInitLen));
378 //perform filtering on series
379 Standard_Real aLastAdded, aLastCandidate;
380 Standard_Boolean isCandidateDefined = Standard_False;
381 aLastAdded = aParamTmp.First();
382 aLastCandidate = aLastAdded;
383 theResult.Append(aParamTmp.First());
385 for(j=2;j<aParamTmp.Length();j++)
387 Standard_Real aVal = aParamTmp.Value(j);
388 if(aVal-aLastAdded > theFilterDist)
391 if(isCandidateDefined) {
392 aLastAdded = aLastCandidate;
393 isCandidateDefined = Standard_False;
400 theResult.Append(aLastAdded);
404 aLastCandidate = aVal;
405 isCandidateDefined = Standard_True;
407 theResult.Append(aParamTmp.Last());
409 if( theResult.Length() == 2 )
411 Standard_Real dist = theResult.Last() - theResult.First();
412 Standard_Integer nbint = (Standard_Integer)((dist / theFilterDist) + 0.5);
416 //Five points more is maximum
423 Standard_Real dU = dist / nbint;
424 for( i = 1; i < nbint; i++ )
426 theResult.InsertAfter(i, theResult.First()+i*dU);
432 void BRepMesh_FastDiscretFace::insertInternalVertices(
433 BRepMeshCol::ListOfVertex& theNewVertices,
434 BRepMesh_Delaun& theMeshBuilder)
436 const Handle(BRepAdaptor_HSurface)& gFace = myAttribute->Surface();
437 switch (gFace->GetType())
440 insertInternalVerticesSphere(theNewVertices);
443 case GeomAbs_Cylinder:
444 insertInternalVerticesCylinder(theNewVertices);
448 insertInternalVerticesCone(theNewVertices);
452 insertInternalVerticesTorus(theNewVertices);
455 case GeomAbs_BezierSurface:
456 case GeomAbs_BSplineSurface:
457 insertInternalVerticesBSpline(theNewVertices);
461 insertInternalVerticesOther(theNewVertices);
465 addVerticesToMesh(theNewVertices, theMeshBuilder);
468 //=======================================================================
469 //function : insertInternalVerticesSphere
471 //=======================================================================
472 void BRepMesh_FastDiscretFace::insertInternalVerticesSphere(
473 BRepMeshCol::ListOfVertex& theNewVertices)
475 const Standard_Real umax = myAttribute->GetUMax();
476 const Standard_Real umin = myAttribute->GetUMin();
477 const Standard_Real vmax = myAttribute->GetVMax();
478 const Standard_Real vmin = myAttribute->GetVMin();
480 gp_Sphere S = myAttribute->Surface()->Sphere();
481 const Standard_Real R = S.Radius();
483 // Calculate parameters for iteration in V direction
484 Standard_Real Dv = 1.0 - (myAttribute->GetDefFace() / R);
485 if (Dv < 0.0) Dv = 0.0;
486 Standard_Real oldDv = 2.0 * ACos(Dv);
487 Dv = .7 * oldDv; //.7 ~= sqrt(2.) - Dv is hypotenuse of triangle when oldDv is legs
488 const Standard_Real sv = vmax - vmin;
489 Dv = sv / ((Standard_Integer)(sv / Dv) + 1);
490 const Standard_Real pasvmax = vmax - Dv*0.5;
492 //Du can be defined from relation: 2*r*Sin(Du/2) = 2*R*Sin(Dv/2), r = R*Cos(v)
493 //here approximate relation r*Du = R*Dv is used
495 Standard_Real Du, pasu, pasv; //, ru;
496 const Standard_Real su = umax - umin;
497 Standard_Boolean Shift = Standard_False;
498 for (pasv = vmin + Dv; pasv < pasvmax; pasv += Dv)
500 // Calculate parameters for iteration in U direction
501 // 1.-.365*pasv*pasv is simple approximation of Cos(pasv)
502 // with condition that it gives ~.1 when pasv = pi/2
503 Du = Dv / (1. - .365*pasv*pasv);
504 Du = su / ((Standard_Integer)(su / Du) + 1);
506 const Standard_Real d = (Shift) ? Du*.5 : 0.;
507 const Standard_Real pasumax = umax - Du*0.5 + d;
508 for (pasu = umin + Du - d; pasu < pasumax; pasu += Du)
510 tryToInsertAnalyticVertex(gp_Pnt2d(pasu, pasv), S, theNewVertices);
515 //=======================================================================
516 //function : insertInternalVerticesCylinder
518 //=======================================================================
519 void BRepMesh_FastDiscretFace::insertInternalVerticesCylinder(
520 BRepMeshCol::ListOfVertex& theNewVertices)
522 const Standard_Real umax = myAttribute->GetUMax();
523 const Standard_Real umin = myAttribute->GetUMin();
524 const Standard_Real vmax = myAttribute->GetVMax();
525 const Standard_Real vmin = myAttribute->GetVMin();
527 gp_Cylinder S = myAttribute->Surface()->Cylinder();
528 const Standard_Real R = S.Radius();
530 // Calculate parameters for iteration in U direction
531 Standard_Real Du = 1.0 - (myAttribute->GetDefFace() / R);
532 if (Du < 0.0) Du = 0.0;
534 if (Du > myAngle) Du = myAngle;
535 const Standard_Real su = umax - umin;
536 const Standard_Integer nbU = (Standard_Integer)(su / Du);
539 // Calculate parameters for iteration in V direction
540 const Standard_Real sv = vmax - vmin;
541 Standard_Integer nbV = (Standard_Integer)(nbU*sv / (su*R));
542 nbV = Min(nbV, 100 * nbU);
543 Standard_Real Dv = sv / (nbV + 1);
545 Standard_Real pasu, pasv, pasvmax = vmax - Dv*0.5, pasumax = umax - Du*0.5;
546 for (pasv = vmin + Dv; pasv < pasvmax; pasv += Dv)
548 for (pasu = umin + Du; pasu < pasumax; pasu += Du)
550 tryToInsertAnalyticVertex(gp_Pnt2d(pasu, pasv), S, theNewVertices);
555 //=======================================================================
556 //function : insertInternalVerticesCone
558 //=======================================================================
559 void BRepMesh_FastDiscretFace::insertInternalVerticesCone(
560 BRepMeshCol::ListOfVertex& theNewVertices)
562 const Standard_Real umax = myAttribute->GetUMax();
563 const Standard_Real umin = myAttribute->GetUMin();
564 const Standard_Real vmax = myAttribute->GetVMax();
565 const Standard_Real vmin = myAttribute->GetVMin();
567 Standard_Real R, RefR, SAng;
568 gp_Cone C = myAttribute->Surface()->Cone();
569 RefR = C.RefRadius();
570 SAng = C.SemiAngle();
571 R = Max(Abs(RefR + vmin*Sin(SAng)), Abs(RefR + vmax*Sin(SAng)));
572 Standard_Real Du, Dv, pasu, pasv;
573 Du = Max(1.0e0 - (myAttribute->GetDefFace() / R), 0.0e0);
574 Du = (2.0 * ACos(Du));
575 Standard_Integer nbU = (Standard_Integer)((umax - umin) / Du);
576 Standard_Integer nbV = (Standard_Integer)(nbU*(vmax - vmin) / ((umax - umin)*R));
577 Du = (umax - umin) / (nbU + 1);
578 Dv = (vmax - vmin) / (nbV + 1);
580 Standard_Real pasvmax = vmax - Dv*0.5, pasumax = umax - Du*0.5;
581 for (pasv = vmin + Dv; pasv < pasvmax; pasv += Dv)
583 for (pasu = umin + Du; pasu < pasumax; pasu += Du)
585 tryToInsertAnalyticVertex(gp_Pnt2d(pasu, pasv), C, theNewVertices);
590 //=======================================================================
591 //function : insertInternalVerticesTorus
593 //=======================================================================
594 void BRepMesh_FastDiscretFace::insertInternalVerticesTorus(
595 BRepMeshCol::ListOfVertex& theNewVertices)
597 const Standard_Real umax = myAttribute->GetUMax();
598 const Standard_Real umin = myAttribute->GetUMin();
599 const Standard_Real vmax = myAttribute->GetVMax();
600 const Standard_Real vmin = myAttribute->GetVMin();
601 const Standard_Real deltaX = myAttribute->GetDeltaX();
602 const Standard_Real deltaY = myAttribute->GetDeltaY();
603 const Standard_Real aDefFace = myAttribute->GetDefFace();
605 gp_Torus T = myAttribute->Surface()->Torus();
607 Standard_Boolean insert;
608 Standard_Integer i, j, ParamULength, ParamVLength;
609 Standard_Real pp, pasu, pasv;
610 Standard_Real r = T.MinorRadius(), R = T.MajorRadius();
612 BRepMeshCol::SequenceOfReal ParamU, ParamV;
614 Standard_Real Du, Dv;//, pasu, pasv;
615 Dv = Max(1.0e0 - (aDefFace / r), 0.0e0);
616 Standard_Real oldDv = 2.0 * ACos(Dv);
617 oldDv = Min(oldDv, myAngle);
618 Dv = 0.9*oldDv; //TWOTHIRD * oldDv;
621 Standard_Integer nbV = Max((Standard_Integer)((vmax - vmin) / Dv), 2);
622 Dv = (vmax - vmin) / (nbV + 1);
623 Standard_Real ru = R + r;
626 Du = 2.0 * ACos(Max(1.0 - (aDefFace / ru), 0.0));
627 if (myAngle < Du) Du = myAngle;
628 Standard_Real aa = sqrt(Du*Du + oldDv*oldDv);
629 if (aa < gp::Resolution())
631 Du *= Min(oldDv, Du) / aa;
635 Standard_Integer nbU = Max((Standard_Integer)((umax - umin) / Du), 2);
636 nbU = Max(nbU, (int)(nbV*(umax - umin)*R / ((vmax - vmin)*r) / 5.));
637 Du = (umax - umin) / (nbU + 1);
641 // As the points of edges are returned.
642 // in this case, the points are not representative.
644 //-- Choose DeltaX and DeltaY so that to avoid skipping points on the grid
645 for (i = 0; i <= nbU; i++) ParamU.Append(umin + i* Du);
650 // Number of mapped U parameters
651 const Standard_Integer LenU = myUParam.Extent();
652 // Fill array of U parameters
653 TColStd_Array1OfReal Up(1, LenU);
654 for (j = 1; j <= LenU; j++) Up(j) = myUParam(j);
656 // Calculate DU, leave array of parameters
657 Standard_Real aDU = FUN_CalcAverageDUV(Up, LenU);
658 aDU = Max(aDU, Abs(umax - umin) / (Standard_Real)nbU / 2.);
659 Standard_Real dUstd = Abs(umax - umin) / (Standard_Real)LenU;
660 if (aDU > dUstd) dUstd = aDU;
662 for (j = 1; j <= LenU; j++)
665 insert = Standard_True;
666 ParamULength = ParamU.Length();
667 for (i = 1; i <= ParamULength && insert; i++)
669 insert = (Abs(ParamU.Value(i) - pp) > (0.5*dUstd));
671 if (insert) ParamU.Append(pp);
676 // Number of mapped V parameters
677 const Standard_Integer LenV = myVParam.Extent();
678 // Fill array of V parameters
679 TColStd_Array1OfReal Vp(1, LenV);
680 for (j = 1; j <= LenV; j++) Vp(j) = myVParam(j);
681 // Calculate DV, sort array of parameters
682 Standard_Real aDV = FUN_CalcAverageDUV(Vp, LenV);
683 aDV = Max(aDV, Abs(vmax - vmin) / (Standard_Real)nbV / 2.);
685 Standard_Real dVstd = Abs(vmax - vmin) / (Standard_Real)LenV;
686 if (aDV > dVstd) dVstd = aDV;
688 for (j = 1; j <= LenV; j++)
692 insert = Standard_True;
693 ParamVLength = ParamV.Length();
694 for (i = 1; i <= ParamVLength && insert; i++)
696 insert = (Abs(ParamV.Value(i) - pp) > (dVstd*2. / 3.));
698 if (insert) ParamV.Append(pp);
701 Standard_Integer Lu = ParamU.Length(), Lv = ParamV.Length();
702 Standard_Real uminnew = umin + deltaY*0.1;
703 Standard_Real vminnew = vmin + deltaX*0.1;
704 Standard_Real umaxnew = umax - deltaY*0.1;
705 Standard_Real vmaxnew = vmax - deltaX*0.1;
707 for (i = 1; i <= Lu; i++)
709 pasu = ParamU.Value(i);
710 if (pasu >= uminnew && pasu < umaxnew)
712 for (j = 1; j <= Lv; j++)
714 pasv = ParamV.Value(j);
715 if (pasv >= vminnew && pasv < vmaxnew)
717 tryToInsertAnalyticVertex(gp_Pnt2d(pasu, pasv), T, theNewVertices);
724 //=======================================================================
725 //function : insertInternalVerticesBSpline
727 //=======================================================================
728 void BRepMesh_FastDiscretFace::insertInternalVerticesBSpline(
729 BRepMeshCol::ListOfVertex& theNewVertices)
731 const Standard_Real aRange[2][2] = {
732 { myAttribute->GetUMax(), myAttribute->GetUMin() },
733 { myAttribute->GetVMax(), myAttribute->GetVMin() }
736 const Standard_Real aDelta[2] = {
737 myAttribute->GetDeltaX(),
738 myAttribute->GetDeltaY()
741 const Standard_Real aDefFace = myAttribute->GetDefFace();
742 const Handle(BRepAdaptor_HSurface)& gFace = myAttribute->Surface();
744 BRepMeshCol::SequenceOfReal aParams[2];
745 for (Standard_Integer i = 0; i < 2; ++i)
747 Standard_Boolean isU = (i == 0);
748 Standard_Real aRes = isU ?
749 gFace->UResolution(aDefFace) :
750 gFace->VResolution(aDefFace);
752 // Sort and filter sequence of parameters
753 Standard_Real aMinDiff = Precision::PConfusion();
755 aMinDiff /= aDelta[i];
757 Standard_Real aRangeDiff = aRange[i][0] - aRange[i][1];
758 Standard_Real aDiffMaxLim = 0.1 * aRangeDiff;
759 Standard_Real aDiff = Min(aDiffMaxLim, Max(0.005 * aRangeDiff, 2. * aRes));
760 filterParameters(isU ? myUParam : myVParam, aMinDiff, aDiff, aParams[i]);
763 // check intermediate isolines
764 Handle(Geom_Surface) aBSpline;
765 GeomAbs_SurfaceType aSurfType = gFace->GetType();
766 if (aSurfType == GeomAbs_BezierSurface)
767 aBSpline = gFace->Bezier();
768 else if (aSurfType == GeomAbs_BSplineSurface)
769 aBSpline = gFace->BSpline();
771 // precision for compare square distances
772 const Standard_Real aPrecision = Precision::Confusion();
773 const Standard_Real aSqPrecision = aPrecision * aPrecision;
774 for (Standard_Integer k = 0; k < 2; ++k)
776 BRepMeshCol::SequenceOfReal& aParams1 = aParams[k];
777 BRepMeshCol::SequenceOfReal& aParams2 = aParams[(k + 1) % 2];
778 const Standard_Boolean isU = (k == 0);
779 Standard_Integer aStartIndex, aEndIndex;
783 aEndIndex = aParams1.Length();
788 aEndIndex = aParams1.Length() - 1;
791 for (Standard_Integer i = aStartIndex; i <= aEndIndex; ++i)
793 const Standard_Real aParam1 = aParams1(i);
794 Handle(Geom_Curve) aIso = isU ?
795 aBSpline->UIso(aParam1) : aBSpline->VIso(aParam1);
797 Standard_Real aPrevParam2 = aParams2(1);
798 gp_Pnt aPrevPnt2 = aIso->Value(aPrevParam2);
799 for (Standard_Integer j = 2; j <= aParams2.Length();)
801 Standard_Real aParam2 = aParams2(j);
802 gp_Pnt aPnt2 = aIso->Value(aParam2);
803 Standard_Real aMidParam = 0.5 * (aPrevParam2 + aParam2);
804 gp_Pnt aMidPnt = aIso->Value(aMidParam);
806 // 23.03.2010 skl for OCC21645 - change precision for comparison
808 if (aPrevPnt2.SquareDistance(aPnt2) > aSqPrecision)
810 gp_Lin aLin(aPrevPnt2, gp_Dir(gp_Vec(aPrevPnt2, aPnt2)));
811 aDist = aLin.Distance(aMidPnt);
814 aDist = aPrevPnt2.Distance(aMidPnt);
816 if (aDist > aDefFace)
819 aParams2.InsertBefore(j, aMidParam);
823 //put regular grig for normals
824 gp_Pnt2d aStPnt1, aStPnt2;
827 aStPnt1 = gp_Pnt2d(aParam1, aPrevParam2);
828 aStPnt2 = gp_Pnt2d(aParam1, aMidParam);
832 aStPnt1 = gp_Pnt2d(aPrevParam2, aParam1);
833 aStPnt2 = gp_Pnt2d(aMidParam, aParam1);
836 gp_Dir N1(0, 0, 1), N2(0, 0, 1);
837 Standard_Boolean aSt1 = GeomLib::NormEstim(aBSpline, aStPnt1, aPrecision, N1);
838 Standard_Boolean aSt2 = GeomLib::NormEstim(aBSpline, aStPnt2, aPrecision, N2);
840 Standard_Real aAngle = N2.Angle(N1);
841 if (aSt1 < 1 && aSt2 < 1 && aAngle > myAngle)
844 aParams2.InsertBefore(j, aMidParam);
848 aPrevParam2 = aParam2;
851 if (!isU && j < aParams2.Length())
853 // Update point parameter.
854 aStPnt1.SetX(aPrevParam2);
856 // Classify intersection point
857 if (myClassifier->Perform(aStPnt1) == TopAbs_IN)
859 insertVertex(aPrevPnt2, aStPnt1.Coord(), theNewVertices);
871 //=======================================================================
872 //function : insertInternalVerticesOther
874 //=======================================================================
875 void BRepMesh_FastDiscretFace::insertInternalVerticesOther(
876 BRepMeshCol::ListOfVertex& theNewVertices)
878 const Standard_Real aAngle = 0.35;
879 const Standard_Real aRange[2][2] = {
880 { myAttribute->GetUMax(), myAttribute->GetUMin() },
881 { myAttribute->GetVMax(), myAttribute->GetVMin() }
884 const Standard_Real aDefFace = myAttribute->GetDefFace();
885 const Handle(BRepAdaptor_HSurface)& gFace = myAttribute->Surface();
887 BRepMeshCol::SequenceOfReal aParams[2];
888 const Standard_Integer aIsoPointsNb = 11;
889 for (Standard_Integer k = 0; k < 2; ++k)
891 const Standard_Boolean isU = (k == 0);
892 const GeomAbs_IsoType aIsoType = isU ? GeomAbs_IsoU : GeomAbs_IsoV;
893 const Standard_Integer aOtherParamsIndex = (k + 1) % 2;
894 const Standard_Real (&aRange1)[2] = aRange[k];
895 const Standard_Real (&aRange2)[2] = aRange[aOtherParamsIndex];
897 GCPnts_TangentialDeflection aDiscretIso[aIsoPointsNb];
898 const Standard_Real aStepWidth = (aRange1[0] - aRange1[1]) / aIsoPointsNb;
900 // Find the most curved Iso.
901 Standard_Integer aMaxIndex = 1, aMaxPointsNb = 0;
902 for (Standard_Integer aIsoIt = 0; aIsoIt < aIsoPointsNb; ++aIsoIt)
904 Standard_Real aParam = aRange1[1] + aIsoIt * aStepWidth;
905 Adaptor3d_IsoCurve aIso(gFace, aIsoType, aParam);
907 Standard_Real aFirstParam = Max(aRange2[1], aIso.FirstParameter());
908 Standard_Real aLastParam = Min(aRange2[0], aIso.LastParameter());
910 aDiscretIso[aIsoIt].Initialize(aIso, aFirstParam, aLastParam,
911 aAngle, 0.7 * aDefFace, 2);
913 const Standard_Integer aPointsNb = aDiscretIso[aIsoIt].NbPoints();
914 if (aPointsNb > aMaxPointsNb)
916 aMaxPointsNb = aPointsNb;
921 BRepMeshCol::SequenceOfReal& aParams2 = aParams[aOtherParamsIndex];
922 GCPnts_TangentialDeflection& aDIso = aDiscretIso[aMaxIndex];
923 for (Standard_Integer i = 1; i <= aMaxPointsNb; ++i)
924 aParams2.Append(aDIso.Parameter(i));
926 if (aParams2.Length() == 2)
927 aParams2.InsertAfter(1, 0.5 * (aRange2[1] + aRange2[0]));
930 Adaptor3d_IsoCurve aIsoV;
933 const Standard_Integer aUPointsNb = aParams[0].Length();
934 const Standard_Integer aVPointsNb = aParams[1].Length();
935 for (Standard_Integer i = 2; i < aVPointsNb; ++i)
937 const Standard_Real aV = aParams[1].Value(i);
938 aIsoV.Load(GeomAbs_IsoV, aV);
940 for (Standard_Integer j = 2; j < aUPointsNb; ++j)
942 const Standard_Real aU = aParams[0].Value(j);
944 const gp_Pnt2d aNewPoint(aU, aV);
945 if (myClassifier->Perform(aNewPoint) == TopAbs_IN)
946 insertVertex(aIsoV.Value(aU), aNewPoint.Coord(), theNewVertices);
951 //=======================================================================
954 //=======================================================================
955 Standard_Real BRepMesh_FastDiscretFace::control(
956 BRepMeshCol::ListOfVertex& theNewVertices,
957 BRepMesh_Delaun& theTrigu,
958 const Standard_Boolean theIsFirst)
960 #define CHECK_DEF_AND_INSERT_CURRENT(isSkipped) \
961 if (aSqDef > aMaxSqDef) \
962 aMaxSqDef = aSqDef; \
964 (isSkipped) = Standard_False; \
965 if (aSqDef > aSqDefFace) \
967 (isSkipped) = theIsFirst; \
969 insertVertex(pDef, mi2d, theNewVertices); \
973 Standard_Integer aTrianglesNb = myStructure->ElementsOfDomain().Extent();
974 if (aTrianglesNb < 1)
977 //IMPORTANT: Constants used in calculations
978 const Standard_Real MinimalArea2d = 1.e-9;
979 const Standard_Real MinimalSqLength3d = 1.e-12;
980 const Standard_Real aSqDefFace = myAttribute->GetDefFace() * myAttribute->GetDefFace();
982 const Handle(BRepAdaptor_HSurface)& gFace = myAttribute->Surface();
984 Handle(Geom_Surface) aBSpline;
985 GeomAbs_SurfaceType aSurfType = gFace->GetType();
986 if (aSurfType == GeomAbs_BezierSurface)
987 aBSpline = gFace->Bezier();
988 else if (aSurfType == GeomAbs_BSplineSurface)
989 aBSpline = gFace->BSpline();
991 NCollection_DataMap<Standard_Integer, gp_Dir> aNorMap;
992 BRepMeshCol::MapOfIntegerInteger aStatMap;
993 NCollection_Map<BRepMesh_OrientedEdge> aCouples(3 * aTrianglesNb);
995 // Perform refinement passes
996 // Define the number of iterations
997 Standard_Integer aIterationsNb = 11;
998 const Standard_Integer aPassesNb = (theIsFirst ? 1 : aIterationsNb);
999 // Initialize stop condition
1000 Standard_Real aMaxSqDef = -1.;
1001 Standard_Integer aPass = 1, aInsertedNb = 1;
1002 Standard_Boolean isAllDegenerated = Standard_False;
1003 for (; aPass <= aPassesNb && aInsertedNb && !isAllDegenerated; ++aPass)
1005 theNewVertices.Clear();
1007 // Reset stop condition
1010 isAllDegenerated = Standard_True;
1012 aTrianglesNb = myStructure->ElementsOfDomain().Extent();
1013 if (aTrianglesNb < 1)
1016 // Iterate on current triangles
1017 const BRepMeshCol::MapOfInteger& aTriangles = myStructure->ElementsOfDomain();
1018 BRepMeshCol::MapOfInteger::Iterator aTriangleIt(aTriangles);
1019 for (; aTriangleIt.More(); aTriangleIt.Next())
1021 const Standard_Integer aTriangleId = aTriangleIt.Key();
1022 const BRepMesh_Triangle& aCurrentTriangle = myStructure->GetElement(aTriangleId);
1024 if (aCurrentTriangle.Movability() == BRepMesh_Deleted)
1027 Standard_Integer v[3];
1028 myStructure->ElementNodes(aCurrentTriangle, v);
1030 Standard_Integer e[3];
1031 Standard_Boolean o[3];
1032 aCurrentTriangle.Edges(e, o);
1036 Standard_Boolean m[3];
1037 for (Standard_Integer i = 0; i < 3; ++i)
1039 m[i] = (myStructure->GetLink(e[i]).Movability() == BRepMesh_Frontier);
1041 const BRepMesh_Vertex& aVertex = myStructure->GetNode(v[i]);
1042 xy[i] = myAttribute->Scale(aVertex.Coord(), Standard_False);
1043 p [i] = myAttribute->GetPoint(aVertex).Coord();
1047 Standard_Boolean isDegeneratedTri = Standard_False;
1048 for (Standard_Integer i = 0; i < 3 && !isDegeneratedTri; ++i)
1050 aLinkVec[i] = p[(i + 1) % 3] - p[i];
1051 isDegeneratedTri = (aLinkVec[i].SquareModulus() < MinimalSqLength3d);
1054 if (isDegeneratedTri)
1057 isAllDegenerated = Standard_False;
1059 // Check triangle area in 2d
1060 if (Abs((xy[1]-xy[0])^(xy[2]-xy[1])) < MinimalArea2d)
1063 // Check triangle normal
1065 Standard_Real aSqDef = -1.;
1066 Standard_Boolean isSkipped = Standard_False;
1067 gp_XYZ normal(aLinkVec[0] ^ aLinkVec[1]);
1074 // Check deflection on triangle
1075 gp_XY mi2d = (xy[0] + xy[1] + xy[2]) / 3.0;
1076 gFace->D0(mi2d.X(), mi2d.Y(), pDef);
1077 aSqDef = Abs(normal * (pDef.XYZ() - p[0]));
1080 CHECK_DEF_AND_INSERT_CURRENT(isSkipped);
1084 catch (Standard_Failure)
1089 for (Standard_Integer i = 0; i < 3 && !isSkipped; ++i)
1091 if (m[i]) // is a boundary
1094 Standard_Integer j = (i + 1) % 3;
1095 // Check if this link was already processed
1096 Standard_Integer aFirstVertex, aLastVertex;
1099 aFirstVertex = v[i];
1104 aFirstVertex = v[j];
1108 if (aCouples.Add(BRepMesh_OrientedEdge(aFirstVertex, aLastVertex)))
1110 // Check deflection on edge 1
1111 gp_XY mi2d = (xy[i] + xy[j]) * 0.5;
1112 gFace->D0(mi2d.X(), mi2d.Y(), pDef);
1113 gp_Lin aLin(p[i], gp_Vec(p[i], p[j]));
1114 aSqDef = aLin.SquareDistance(pDef);
1116 CHECK_DEF_AND_INSERT_CURRENT(isSkipped);
1123 //check normal on bsplines
1124 if (theIsFirst && !aBSpline.IsNull())
1126 gp_Dir N[3] = { gp::DZ(), gp::DZ(), gp::DZ() };
1127 Standard_Integer aSt[3];
1129 for (Standard_Integer i = 0; i < 3; ++i)
1131 if (aNorMap.IsBound(v[i]))
1133 aSt[i] = aStatMap.Find(v[i]);
1134 N[i] = aNorMap.Find(v[i]);
1138 aSt[i] = GeomLib::NormEstim(aBSpline, gp_Pnt2d(xy[i]), Precision::Confusion(), N[i]);
1139 aStatMap.Bind(v[i], aSt[i]);
1140 aNorMap.Bind(v[i], N[i]);
1144 Standard_Real aAngle[3];
1145 for (Standard_Integer i = 0; i < 3; ++i)
1146 aAngle[i] = N[(i + 1) % 3].Angle(N[i]);
1148 if (aSt[0] < 1 && aSt[1] < 1 && aSt[2] < 1)
1150 if (aAngle[0] > myAngle || aAngle[1] > myAngle || aAngle[2] > myAngle)
1162 if (addVerticesToMesh(theNewVertices, theTrigu))
1166 return (aMaxSqDef < 0) ? aMaxSqDef : Sqrt(aMaxSqDef);
1169 //=======================================================================
1172 //=======================================================================
1173 void BRepMesh_FastDiscretFace::add(const TopoDS_Vertex& theVertex)
1175 if (theVertex.Orientation() != TopAbs_INTERNAL)
1182 gp_Pnt2d aPnt2d = BRep_Tool::Parameters(theVertex, myAttribute->Face());
1184 N_HANDLE<FixedVExplorer> aFixedVExplorer = new FixedVExplorer(theVertex);
1185 Standard_Integer aIndex = myAttribute->GetVertexIndex(aFixedVExplorer);
1186 gp_XY anUV = BRepMesh_ShapeTool::FindUV(aIndex, aPnt2d,
1187 theVertex, BRep_Tool::Tolerance(theVertex), myAttribute);
1189 Standard_Integer aTmpId1, aTmpId2;
1190 myAttribute->AddNode(aIndex, anUV, BRepMesh_Fixed, aTmpId1, aTmpId2);
1192 catch (Standard_Failure)
1197 //=======================================================================
1198 //function : insertVertex
1200 //=======================================================================
1201 void BRepMesh_FastDiscretFace::insertVertex(
1202 const gp_Pnt& thePnt3d,
1204 BRepMeshCol::ListOfVertex& theVertices)
1206 Standard_Integer aNbLocat = myAttribute->LastPointId();
1207 mySurfacePoints->Bind(++aNbLocat, thePnt3d);
1209 gp_XY aPnt2d = myAttribute->Scale(theUV, Standard_True);
1210 BRepMesh_Vertex aVertex(aPnt2d, aNbLocat, BRepMesh_Free);
1211 theVertices.Append(aVertex);