1 // Created on: 2016-07-07
2 // Copyright (c) 2016 OPEN CASCADE SAS
3 // Created by: Oleg AGASHIN
5 // This file is part of Open CASCADE Technology software library.
7 // This library is free software; you can redistribute it and/or modify it under
8 // the terms of the GNU Lesser General Public License version 2.1 as published
9 // by the Free Software Foundation, with special exception defined in the file
10 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
11 // distribution for complete text of the license and disclaimer of any warranty.
13 // Alternatively, this file may be used under the terms of Open CASCADE
14 // commercial license or contractual agreement.
16 #ifndef _BRepMeshTools_DelaunayDeflectionControlMeshAlgo_HeaderFile
17 #define _BRepMeshTools_DelaunayDeflectionControlMeshAlgo_HeaderFile
19 #include <BRepMesh_DelaunayNodeInsertionMeshAlgo.hxx>
20 #include <BRepMesh_GeomTool.hxx>
21 #include <GeomLib.hxx>
23 //! Extends node insertion Delaunay meshing algo in order to control
24 //! deflection of generated trianges. Splits triangles failing the check.
25 template<class RangeSplitter, class BaseAlgo>
26 class BRepMesh_DelaunayDeflectionControlMeshAlgo : public BRepMesh_DelaunayNodeInsertionMeshAlgo<RangeSplitter, BaseAlgo>
29 // Typedef for OCCT RTTI
30 typedef BRepMesh_DelaunayNodeInsertionMeshAlgo<RangeSplitter, BaseAlgo> DelaunayInsertionBaseClass;
35 BRepMesh_DelaunayDeflectionControlMeshAlgo()
36 : myMaxSqDeflection(-1.),
37 myIsAllDegenerated(Standard_False),
43 virtual ~BRepMesh_DelaunayDeflectionControlMeshAlgo()
49 //! Performs processing of generated mesh. Generates surface nodes and inserts them into structure.
50 virtual void postProcessMesh (BRepMesh_Delaun& theMesher,
51 const Message_ProgressRange& theRange) Standard_OVERRIDE
53 Message_ProgressScope aPS(theRange, "Post process mesh", 2);
54 // Insert surface nodes.
55 DelaunayInsertionBaseClass::postProcessMesh (theMesher, aPS.Next());
61 if (this->getParameters().ControlSurfaceDeflection &&
62 this->getStructure()->ElementsOfDomain().Extent() > 0)
64 optimizeMesh(theMesher, aPS.Next());
72 //! Checks deviation of a mesh from geometrical surface.
73 //! Inserts additional nodes in case of huge deviation.
74 virtual void optimizeMesh (BRepMesh_Delaun& theMesher,
75 const Message_ProgressRange& theRange)
77 Handle(NCollection_IncAllocator) aTmpAlloc =
78 new NCollection_IncAllocator(IMeshData::MEMORY_BLOCK_SIZE_HUGE);
80 myCouplesMap = new IMeshData::MapOfOrientedEdges(3 * this->getStructure()->ElementsOfDomain().Extent(), aTmpAlloc);
81 myControlNodes = new IMeshData::ListOfPnt2d(aTmpAlloc);
82 myCircles = &theMesher.Circles();
84 const Standard_Integer aIterationsNb = 11;
85 Standard_Boolean isInserted = Standard_True;
86 Message_ProgressScope aPS(theRange, "Iteration", aIterationsNb);
87 for (Standard_Integer aPass = 1; aPass <= aIterationsNb && isInserted && !myIsAllDegenerated; ++aPass)
93 // Reset stop condition
94 myMaxSqDeflection = -1.;
95 myIsAllDegenerated = Standard_True;
96 myControlNodes->Clear();
98 if (this->getStructure()->ElementsOfDomain().Extent() < 1)
102 // Iterate on current triangles
103 IMeshData::IteratorOfMapOfInteger aTriangleIt(this->getStructure()->ElementsOfDomain());
104 for (; aTriangleIt.More(); aTriangleIt.Next())
106 const BRepMesh_Triangle& aTriangle = this->getStructure()->GetElement(aTriangleIt.Key());
107 splitTriangleGeometry(aTriangle);
110 isInserted = this->insertNodes(myControlNodes, theMesher, aPS.Next());
113 myCouplesMap.Nullify();
114 myControlNodes.Nullify();
116 if (!(myMaxSqDeflection < 0.))
118 this->getDFace()->SetDeflection(Sqrt(myMaxSqDeflection));
123 //! Contains geometrical data related to node of triangle.
124 struct TriangleNodeInfo
127 : isFrontierLink(Standard_False)
133 Standard_Boolean isFrontierLink;
136 //! Functor computing deflection of a point from surface.
137 class NormalDeviation
141 const gp_Pnt& theRefPnt,
142 const gp_Vec& theNormal)
143 : myRefPnt(theRefPnt),
148 Standard_Real SquareDeviation(const gp_Pnt& thePoint) const
150 const Standard_Real aDeflection = Abs(myNormal.Dot(gp_Vec(myRefPnt, thePoint)));
151 return aDeflection * aDeflection;
156 NormalDeviation (const NormalDeviation& theOther);
158 void operator= (const NormalDeviation& theOther);
162 const gp_Pnt& myRefPnt;
163 const gp_Vec& myNormal;
166 //! Functor computing deflection of a point on triangle link from surface.
172 const gp_Pnt& thePnt1,
173 const gp_Pnt& thePnt2)
179 Standard_Real SquareDeviation(const gp_Pnt& thePoint) const
181 return BRepMesh_GeomTool::SquareDeflectionOfSegment(myPnt1, myPnt2, thePoint);
186 LineDeviation (const LineDeviation& theOther);
188 void operator= (const LineDeviation& theOther);
191 const gp_Pnt& myPnt1;
192 const gp_Pnt& myPnt2;
195 //! Returns nodes info of the given triangle.
196 void getTriangleInfo(
197 const BRepMesh_Triangle& theTriangle,
198 const Standard_Integer (&theNodesIndices)[3],
199 TriangleNodeInfo (&theInfo)[3]) const
201 const Standard_Integer(&e)[3] = theTriangle.myEdges;
202 for (Standard_Integer i = 0; i < 3; ++i)
204 const BRepMesh_Vertex& aVertex = this->getStructure()->GetNode(theNodesIndices[i]);
205 theInfo[i].Point2d = this->getRangeSplitter().Scale(aVertex.Coord(), Standard_False).XY();
206 theInfo[i].Point = this->getNodesMap()->Value(aVertex.Location3d()).XYZ();
207 theInfo[i].isFrontierLink = (this->getStructure()->GetLink(e[i]).Movability() == BRepMesh_Frontier);
211 // Check geometry of the given triangle. If triangle does not suit specified deflection, inserts new point.
212 void splitTriangleGeometry(const BRepMesh_Triangle& theTriangle)
214 if (theTriangle.Movability() != BRepMesh_Deleted)
216 Standard_Integer aNodexIndices[3];
217 this->getStructure()->ElementNodes(theTriangle, aNodexIndices);
219 TriangleNodeInfo aNodesInfo[3];
220 getTriangleInfo(theTriangle, aNodexIndices, aNodesInfo);
224 if (computeTriangleGeometry(aNodesInfo, aLinkVec, aNormal))
226 myIsAllDegenerated = Standard_False;
228 const gp_XY aCenter2d = (aNodesInfo[0].Point2d +
229 aNodesInfo[1].Point2d +
230 aNodesInfo[2].Point2d) / 3.;
232 usePoint(aCenter2d, NormalDeviation(aNodesInfo[0].Point, aNormal));
233 splitLinks(aNodesInfo, aNodexIndices);
238 //! Updates array of links vectors.
239 //! @return False on degenerative triangle.
240 Standard_Boolean computeTriangleGeometry(
241 const TriangleNodeInfo(&theNodesInfo)[3],
242 gp_Vec (&theLinks)[3],
245 if (checkTriangleForDegenerativityAndGetLinks(theNodesInfo, theLinks))
247 if (checkTriangleArea2d(theNodesInfo))
249 if (computeNormal(theLinks[0], theLinks[1], theNormal))
251 return Standard_True;
256 return Standard_False;
259 //! Updates array of links vectors.
260 //! @return False on degenerative triangle.
261 Standard_Boolean checkTriangleForDegenerativityAndGetLinks(
262 const TriangleNodeInfo (&theNodesInfo)[3],
263 gp_Vec (&theLinks)[3])
265 const Standard_Real MinimalSqLength3d = 1.e-12;
266 for (Standard_Integer i = 0; i < 3; ++i)
268 theLinks[i] = theNodesInfo[(i + 1) % 3].Point - theNodesInfo[i].Point;
269 if (theLinks[i].SquareMagnitude() < MinimalSqLength3d)
271 return Standard_False;
275 return Standard_True;
278 //! Checks area of triangle in parametric space for degenerativity.
279 //! @return False on degenerative triangle.
280 Standard_Boolean checkTriangleArea2d(
281 const TriangleNodeInfo (&theNodesInfo)[3])
283 const gp_Vec2d aLink2d1(theNodesInfo[0].Point2d, theNodesInfo[1].Point2d);
284 const gp_Vec2d aLink2d2(theNodesInfo[1].Point2d, theNodesInfo[2].Point2d);
286 const Standard_Real MinimalArea2d = 1.e-9;
287 return (Abs(aLink2d1 ^ aLink2d2) > MinimalArea2d);
290 //! Computes normal using two link vectors.
291 //! @return True on success, False in case of normal of null magnitude.
292 Standard_Boolean computeNormal(const gp_Vec& theLink1,
293 const gp_Vec& theLink2,
296 const gp_Vec aNormal(theLink1 ^ theLink2);
297 if (aNormal.SquareMagnitude() > gp::Resolution())
299 theNormal = aNormal.Normalized();
300 return Standard_True;
303 return Standard_False;
306 //! Computes deflection of midpoints of triangles links.
307 //! @return True if point fits specified deflection.
309 const TriangleNodeInfo (&theNodesInfo)[3],
310 const Standard_Integer (&theNodesIndices)[3])
312 // Check deflection at triangle links
313 for (Standard_Integer i = 0; i < 3; ++i)
315 if (theNodesInfo[i].isFrontierLink)
320 const Standard_Integer j = (i + 1) % 3;
321 // Check if this link was already processed
322 Standard_Integer aFirstVertex, aLastVertex;
323 if (theNodesIndices[i] < theNodesIndices[j])
325 aFirstVertex = theNodesIndices[i];
326 aLastVertex = theNodesIndices[j];
330 aFirstVertex = theNodesIndices[j];
331 aLastVertex = theNodesIndices[i];
334 if (myCouplesMap->Add(BRepMesh_OrientedEdge(aFirstVertex, aLastVertex)))
336 const gp_XY aMidPnt2d = (theNodesInfo[i].Point2d +
337 theNodesInfo[j].Point2d) / 2.;
339 if (!usePoint (aMidPnt2d, LineDeviation (theNodesInfo[i].Point,
340 theNodesInfo[j].Point)))
342 if (!checkLinkEndsForAngularDeviation(theNodesInfo[i],
346 myControlNodes->Append(aMidPnt2d);
353 //! Checks the given point (located between the given nodes)
354 //! for specified angular deviation.
355 Standard_Boolean checkLinkEndsForAngularDeviation(const TriangleNodeInfo& theNodeInfo1,
356 const TriangleNodeInfo& theNodeInfo2,
357 const gp_XY& /*theMidPoint*/)
359 gp_Dir aNorm1, aNorm2;
360 const Handle(Geom_Surface)& aSurf = this->getDFace()->GetSurface()->Surface().Surface();
362 if ((GeomLib::NormEstim(aSurf, theNodeInfo1.Point2d, Precision::Confusion(), aNorm1) == 0) &&
363 (GeomLib::NormEstim(aSurf, theNodeInfo2.Point2d, Precision::Confusion(), aNorm2) == 0))
365 Standard_Real anAngle = aNorm1.Angle(aNorm2);
366 if (anAngle > this->getParameters().AngleInterior)
367 return Standard_False;
370 else if (GeomLib::NormEstim(aSurf, theMidPoint, Precision::Confusion(), aNorm1) != 0)
372 // It is better to consider the singular point as a node of triangulation.
373 // However, it leads to hangs up meshing some faces (including faces with
374 // degenerated edges). E.g. tests "mesh standard_incmesh Q6".
375 // So, this code fragment is better to implement in the future.
376 return Standard_False;
380 return Standard_True;
383 //! Computes deflection of the given point and caches it for
384 //! insertion in case if it overflows deflection.
385 //! @return True if point has been cached for insertion.
386 template<class DeflectionFunctor>
387 Standard_Boolean usePoint(
388 const gp_XY& thePnt2d,
389 const DeflectionFunctor& theDeflectionFunctor)
392 this->getDFace()->GetSurface()->D0(thePnt2d.X(), thePnt2d.Y(), aPnt);
393 if (!checkDeflectionOfPointAndUpdateCache(thePnt2d, aPnt, theDeflectionFunctor.SquareDeviation(aPnt)))
395 myControlNodes->Append(thePnt2d);
396 return Standard_True;
399 return Standard_False;
402 //! Checks the given point for specified linear deflection.
403 //! Updates value of total mesh defleciton.
404 Standard_Boolean checkDeflectionOfPointAndUpdateCache(
405 const gp_XY& thePnt2d,
406 const gp_Pnt& thePnt3d,
407 const Standard_Real theSqDeflection)
409 if (theSqDeflection > myMaxSqDeflection)
411 myMaxSqDeflection = theSqDeflection;
414 const Standard_Real aSqDeflection =
415 this->getDFace()->GetDeflection() * this->getDFace()->GetDeflection();
416 if (theSqDeflection < aSqDeflection)
418 return Standard_True;
421 return rejectByMinSize(thePnt2d, thePnt3d);
424 //! Checks the given node for
425 Standard_Boolean rejectByMinSize(
426 const gp_XY& thePnt2d,
427 const gp_Pnt& thePnt3d)
429 const Standard_Real aSqMinSize =
430 this->getParameters().MinSize * this->getParameters().MinSize;
432 IMeshData::MapOfInteger aUsedNodes;
433 IMeshData::ListOfInteger& aCirclesList =
434 const_cast<BRepMesh_CircleTool&>(*myCircles).Select(
435 this->getRangeSplitter().Scale(thePnt2d, Standard_True).XY());
437 IMeshData::ListOfInteger::Iterator aCircleIt(aCirclesList);
438 for (; aCircleIt.More(); aCircleIt.Next())
440 const BRepMesh_Triangle& aTriangle = this->getStructure()->GetElement(aCircleIt.Value());
442 Standard_Integer aNodes[3];
443 this->getStructure()->ElementNodes(aTriangle, aNodes);
445 for (Standard_Integer i = 0; i < 3; ++i)
447 if (!aUsedNodes.Contains(aNodes[i]))
449 aUsedNodes.Add(aNodes[i]);
450 const BRepMesh_Vertex& aVertex = this->getStructure()->GetNode(aNodes[i]);
451 const gp_Pnt& aPoint = this->getNodesMap()->Value(aVertex.Location3d());
453 if (thePnt3d.SquareDistance(aPoint) < aSqMinSize)
455 return Standard_True;
461 return Standard_False;
465 Standard_Real myMaxSqDeflection;
466 Standard_Boolean myIsAllDegenerated;
467 Handle(IMeshData::MapOfOrientedEdges) myCouplesMap;
468 Handle(IMeshData::ListOfPnt2d) myControlNodes;
469 const BRepMesh_CircleTool* myCircles;