0031035: Coding - uninitialized class fields reported by Visual Studio Code Analysis
[occt.git] / src / BRepMesh / BRepMesh_DelaunayDeflectionControlMeshAlgo.hxx
1 // Created on: 2016-07-07
2 // Copyright (c) 2016 OPEN CASCADE SAS
3 // Created by: Oleg AGASHIN
4 //
5 // This file is part of Open CASCADE Technology software library.
6 //
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.
12 //
13 // Alternatively, this file may be used under the terms of Open CASCADE
14 // commercial license or contractual agreement.
15
16 #ifndef _BRepMeshTools_DelaunayDeflectionControlMeshAlgo_HeaderFile
17 #define _BRepMeshTools_DelaunayDeflectionControlMeshAlgo_HeaderFile
18
19 #include <BRepMesh_DelaunayNodeInsertionMeshAlgo.hxx>
20 #include <BRepMesh_GeomTool.hxx>
21 #include <GeomLib.hxx>
22
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>
27 {
28 private:
29   // Typedef for OCCT RTTI
30   typedef BRepMesh_DelaunayNodeInsertionMeshAlgo<RangeSplitter, BaseAlgo> DelaunayInsertionBaseClass;
31
32 public:
33
34   //! Constructor.
35   BRepMesh_DelaunayDeflectionControlMeshAlgo()
36     : myMaxSqDeflection(-1.),
37       myIsAllDegenerated(Standard_False),
38       myCircles(NULL)
39   {
40   }
41
42   //! Destructor.
43   virtual ~BRepMesh_DelaunayDeflectionControlMeshAlgo()
44   {
45   }
46
47 protected:
48
49   //! Perfroms processing of generated mesh. Generates surface nodes and inserts them into structure.
50   virtual void postProcessMesh(BRepMesh_Delaun& theMesher) Standard_OVERRIDE
51   {
52     // Insert surface nodes.
53     DelaunayInsertionBaseClass::postProcessMesh(theMesher);
54
55     if (this->getParameters().ControlSurfaceDeflection &&
56         this->getStructure()->ElementsOfDomain().Extent() > 0)
57     {
58       optimizeMesh(theMesher);
59     }
60   }
61
62   //! Checks deviation of a mesh from geometrical surface.
63   //! Inserts additional nodes in case of huge deviation.
64   virtual void optimizeMesh(BRepMesh_Delaun& theMesher)
65   {
66     Handle(NCollection_IncAllocator) aTmpAlloc =
67       new NCollection_IncAllocator(IMeshData::MEMORY_BLOCK_SIZE_HUGE);
68
69     myCouplesMap   = new IMeshData::MapOfOrientedEdges(3 * this->getStructure()->ElementsOfDomain().Extent(), aTmpAlloc);
70     myControlNodes = new IMeshData::ListOfPnt2d(aTmpAlloc);
71     myCircles      = &theMesher.Circles();
72
73     const Standard_Integer aIterationsNb = 11;
74     Standard_Boolean isInserted = Standard_True;
75     for (Standard_Integer aPass = 1; aPass <= aIterationsNb && isInserted && !myIsAllDegenerated; ++aPass)
76     {
77       // Reset stop condition
78       myMaxSqDeflection = -1.;
79       myIsAllDegenerated = Standard_True;
80       myControlNodes->Clear();
81
82       if (this->getStructure()->ElementsOfDomain().Extent() < 1)
83       {
84         break;
85       }
86
87       // Iterate on current triangles
88       IMeshData::IteratorOfMapOfInteger aTriangleIt(this->getStructure()->ElementsOfDomain());
89       for (; aTriangleIt.More(); aTriangleIt.Next())
90       {
91         const BRepMesh_Triangle& aTriangle = this->getStructure()->GetElement(aTriangleIt.Key());
92         splitTriangleGeometry(aTriangle);
93       }
94
95       isInserted = this->insertNodes(myControlNodes, theMesher);
96     }
97
98     myCouplesMap.Nullify();
99     myControlNodes.Nullify();
100
101     if (!(myMaxSqDeflection < 0.))
102     {
103       this->getDFace()->SetDeflection(Sqrt(myMaxSqDeflection));
104     }
105   }
106
107 private:
108   //! Contains geometrical data related to node of triangle.
109   struct TriangleNodeInfo
110   {
111     TriangleNodeInfo()
112     : isFrontierLink(Standard_False)
113     {
114     }
115
116     gp_XY            Point2d;
117     gp_XYZ           Point;
118     Standard_Boolean isFrontierLink;
119   };
120
121   //! Functor computing deflection of a point from surface.
122   class NormalDeviation
123   {
124   public:
125     NormalDeviation(
126       const gp_Pnt& theRefPnt,
127       const gp_Vec& theNormal)
128       : myRefPnt(theRefPnt),
129         myNormal(theNormal)
130     {
131     }
132
133     Standard_Real SquareDeviation(const gp_Pnt& thePoint) const
134     {
135       const Standard_Real aDeflection = Abs(myNormal.Dot(gp_Vec(myRefPnt, thePoint)));
136       return aDeflection * aDeflection;
137     }
138
139   private:
140
141     NormalDeviation (const NormalDeviation& theOther);
142
143     void operator= (const NormalDeviation& theOther);
144
145   private:
146
147     const gp_Pnt& myRefPnt;
148     const gp_Vec& myNormal;
149   };
150
151   //! Functor computing deflection of a point on triangle link from surface.
152   class LineDeviation
153   {
154   public:
155
156     LineDeviation(
157       const gp_Pnt& thePnt1,
158       const gp_Pnt& thePnt2)
159       : myPnt1(thePnt1),
160         myPnt2(thePnt2)
161     {
162     }
163
164     Standard_Real SquareDeviation(const gp_Pnt& thePoint) const
165     {
166       return BRepMesh_GeomTool::SquareDeflectionOfSegment(myPnt1, myPnt2, thePoint);
167     }
168
169   private:
170
171     LineDeviation (const LineDeviation& theOther);
172
173     void operator= (const LineDeviation& theOther);
174
175   private:
176     const gp_Pnt& myPnt1;
177     const gp_Pnt& myPnt2;
178   };
179
180   //! Returns nodes info of the given triangle.
181   inline void getTriangleInfo(
182     const BRepMesh_Triangle& theTriangle,
183     const Standard_Integer (&theNodesIndices)[3],
184     TriangleNodeInfo       (&theInfo)[3]) const
185   {
186     const Standard_Integer(&e)[3] = theTriangle.myEdges;
187     for (Standard_Integer i = 0; i < 3; ++i)
188     {
189       const BRepMesh_Vertex& aVertex = this->getStructure()->GetNode(theNodesIndices[i]);
190       theInfo[i].Point2d        = this->getRangeSplitter().Scale(aVertex.Coord(), Standard_False).XY();
191       theInfo[i].Point          = this->getNodesMap()->Value(aVertex.Location3d()).XYZ();
192       theInfo[i].isFrontierLink = (this->getStructure()->GetLink(e[i]).Movability() == BRepMesh_Frontier);
193     }
194   }
195
196   // Check geometry of the given triangle. If triangle does not suit specified deflection, inserts new point.
197   void splitTriangleGeometry(const BRepMesh_Triangle& theTriangle)
198   {
199     if (theTriangle.Movability() != BRepMesh_Deleted)
200     {
201       Standard_Integer aNodexIndices[3];
202       this->getStructure()->ElementNodes(theTriangle, aNodexIndices);
203
204       TriangleNodeInfo aNodesInfo[3];
205       getTriangleInfo(theTriangle, aNodexIndices, aNodesInfo);
206
207       gp_Vec aNormal;
208       gp_Vec aLinkVec[3];
209       if (computeTriangleGeometry(aNodesInfo, aLinkVec, aNormal))
210       {
211         myIsAllDegenerated = Standard_False;
212
213         const gp_XY aCenter2d = (aNodesInfo[0].Point2d +
214                                  aNodesInfo[1].Point2d +
215                                  aNodesInfo[2].Point2d) / 3.;
216
217         usePoint(aCenter2d, NormalDeviation(aNodesInfo[0].Point, aNormal));
218         splitLinks(aNodesInfo, aNodexIndices);
219       }
220     }
221   }
222
223   //! Updates array of links vectors.
224   //! @return False on degenerative triangle.
225   inline Standard_Boolean computeTriangleGeometry(
226     const TriangleNodeInfo(&theNodesInfo)[3],
227     gp_Vec                (&theLinks)[3],
228     gp_Vec                 &theNormal)
229   {
230     if (checkTriangleForDegenerativityAndGetLinks(theNodesInfo, theLinks))
231     {
232       if (checkTriangleArea2d(theNodesInfo))
233       {
234         if (computeNormal(theLinks[0], theLinks[1], theNormal))
235         {
236           return Standard_True;
237         }
238       }
239     }
240
241     return Standard_False;
242   }
243
244   //! Updates array of links vectors.
245   //! @return False on degenerative triangle.
246   inline Standard_Boolean checkTriangleForDegenerativityAndGetLinks(
247     const TriangleNodeInfo (&theNodesInfo)[3],
248     gp_Vec                 (&theLinks)[3])
249   {
250     const Standard_Real MinimalSqLength3d = 1.e-12;
251     for (Standard_Integer i = 0; i < 3; ++i)
252     {
253       theLinks[i] = theNodesInfo[(i + 1) % 3].Point - theNodesInfo[i].Point;
254       if (theLinks[i].SquareMagnitude() < MinimalSqLength3d)
255       {
256         return Standard_False;
257       }
258     }
259
260     return Standard_True;
261   }
262
263   //! Checks area of triangle in parametric space for degenerativity.
264   //! @return False on degenerative triangle.
265   inline Standard_Boolean checkTriangleArea2d(
266     const TriangleNodeInfo (&theNodesInfo)[3])
267   {
268     const gp_Vec2d aLink2d1(theNodesInfo[0].Point2d, theNodesInfo[1].Point2d);
269     const gp_Vec2d aLink2d2(theNodesInfo[1].Point2d, theNodesInfo[2].Point2d);
270
271     const Standard_Real MinimalArea2d = 1.e-9;
272     return (Abs(aLink2d1 ^ aLink2d2) > MinimalArea2d);
273   }
274
275   //! Computes normal using two link vectors.
276   //! @return True on success, False in case of normal of null magnitude.
277   inline Standard_Boolean computeNormal(const gp_Vec& theLink1,
278                                         const gp_Vec& theLink2,
279                                         gp_Vec&       theNormal)
280   {
281     const gp_Vec aNormal(theLink1 ^ theLink2);
282     if (aNormal.SquareMagnitude() > gp::Resolution())
283     {
284       theNormal = aNormal.Normalized();
285       return Standard_True;
286     }
287
288     return Standard_False;
289   }
290
291   //! Computes deflection of midpoints of triangles links.
292   //! @return True if point fits specified deflection.
293   inline void splitLinks(
294     const TriangleNodeInfo (&theNodesInfo)[3],
295     const Standard_Integer (&theNodesIndices)[3])
296   {
297     // Check deflection at triangle links
298     for (Standard_Integer i = 0; i < 3; ++i)
299     {
300       if (theNodesInfo[i].isFrontierLink)
301       {
302         continue;
303       }
304
305       const Standard_Integer j = (i + 1) % 3;
306       // Check if this link was already processed
307       Standard_Integer aFirstVertex, aLastVertex;
308       if (theNodesIndices[i] < theNodesIndices[j])
309       {
310         aFirstVertex = theNodesIndices[i];
311         aLastVertex  = theNodesIndices[j];
312       }
313       else
314       {
315         aFirstVertex = theNodesIndices[j];
316         aLastVertex  = theNodesIndices[i];
317       }
318
319       if (myCouplesMap->Add(BRepMesh_OrientedEdge(aFirstVertex, aLastVertex)))
320       {
321         const gp_XY aMidPnt2d = (theNodesInfo[i].Point2d +
322                                  theNodesInfo[j].Point2d) / 2.;
323
324         if (!usePoint (aMidPnt2d, LineDeviation (theNodesInfo[i].Point, 
325                                                  theNodesInfo[j].Point)))
326         {
327           if (!checkLinkEndsForAngularDeviation(theNodesInfo[i], 
328                                                 theNodesInfo[j],
329                                                 aMidPnt2d))
330           {
331             myControlNodes->Append(aMidPnt2d);
332           }
333         }
334       }
335     }
336   }
337
338   //! Checks the given point (located between the given nodes)
339   //! for specified angular deviation.
340   Standard_Boolean checkLinkEndsForAngularDeviation(const TriangleNodeInfo& theNodeInfo1,
341                                                     const TriangleNodeInfo& theNodeInfo2,
342                                                     const gp_XY& /*theMidPoint*/)
343   {
344     gp_Dir aNorm1, aNorm2;
345     const Handle(Geom_Surface)& aSurf = 
346       this->getDFace()->GetSurface()->ChangeSurface().Surface().Surface();
347     
348     if ((GeomLib::NormEstim(aSurf, theNodeInfo1.Point2d, Precision::Confusion(), aNorm1) == 0) &&
349         (GeomLib::NormEstim(aSurf, theNodeInfo2.Point2d, Precision::Confusion(), aNorm2) == 0))
350     {
351       Standard_Real anAngle = aNorm1.Angle(aNorm2);
352       if (anAngle > this->getParameters().AngleInterior)
353         return Standard_False;
354     }
355 #if 0
356     else if (GeomLib::NormEstim(aSurf, theMidPoint, Precision::Confusion(), aNorm1) != 0)
357     {
358       // It is better to consider the singular point as a node of triangulation.
359       // However, it leads to hangs up meshing some faces (including faces with
360       // degenerated edges). E.g. tests "mesh standard_incmesh Q6".
361       // So, this code fragment is better to implement in the future.
362       return Standard_False;
363     }
364 #endif
365
366     return Standard_True;
367   }
368
369   //! Computes deflection of the given point and caches it for
370   //! insertion in case if it overflows deflection.
371   //! @return True if point has been cached for insertion.
372   template<class DeflectionFunctor>
373   inline Standard_Boolean usePoint(
374     const gp_XY&             thePnt2d,
375     const DeflectionFunctor& theDeflectionFunctor)
376   {
377     gp_Pnt aPnt;
378     this->getDFace()->GetSurface()->D0(thePnt2d.X(), thePnt2d.Y(), aPnt);
379     if (!checkDeflectionOfPointAndUpdateCache(thePnt2d, aPnt, theDeflectionFunctor.SquareDeviation(aPnt)))
380     {
381       myControlNodes->Append(thePnt2d);
382       return Standard_True;
383     }
384
385     return Standard_False;
386   }
387
388   //! Checks the given point for specified linear deflection.
389   //! Updates value of total mesh defleciton.
390   Standard_Boolean checkDeflectionOfPointAndUpdateCache(
391     const gp_XY&        thePnt2d,
392     const gp_Pnt&       thePnt3d,
393     const Standard_Real theSqDeflection)
394   {
395     if (theSqDeflection > myMaxSqDeflection)
396     {
397       myMaxSqDeflection = theSqDeflection;
398     }
399
400     const Standard_Real aSqDeflection = 
401       this->getDFace()->GetDeflection() * this->getDFace()->GetDeflection();
402     if (theSqDeflection < aSqDeflection)
403     {
404       return Standard_True;
405     }
406
407     return rejectByMinSize(thePnt2d, thePnt3d);
408   }
409
410   //! Checks the given node for 
411   Standard_Boolean rejectByMinSize(
412     const gp_XY&  thePnt2d,
413     const gp_Pnt& thePnt3d)
414   {
415     const Standard_Real aSqMinSize = 
416       this->getParameters().MinSize * this->getParameters().MinSize;
417
418     IMeshData::MapOfInteger aUsedNodes;
419     IMeshData::ListOfInteger& aCirclesList =
420       const_cast<BRepMesh_CircleTool&>(*myCircles).Select(
421         this->getRangeSplitter().Scale(thePnt2d, Standard_True).XY());
422
423     IMeshData::ListOfInteger::Iterator aCircleIt(aCirclesList);
424     for (; aCircleIt.More(); aCircleIt.Next())
425     {
426       const BRepMesh_Triangle& aTriangle = this->getStructure()->GetElement(aCircleIt.Value());
427
428       Standard_Integer aNodes[3];
429       this->getStructure()->ElementNodes(aTriangle, aNodes);
430
431       for (Standard_Integer i = 0; i < 3; ++i)
432       {
433         if (!aUsedNodes.Contains(aNodes[i]))
434         {
435           aUsedNodes.Add(aNodes[i]);
436           const BRepMesh_Vertex& aVertex = this->getStructure()->GetNode(aNodes[i]);
437           const gp_Pnt& aPoint = this->getNodesMap()->Value(aVertex.Location3d());
438
439           if (thePnt3d.SquareDistance(aPoint) < aSqMinSize)
440           {
441             return Standard_True;
442           }
443         }
444       }
445     }
446
447     return Standard_False;
448   }
449
450 private:
451   Standard_Real                         myMaxSqDeflection;
452   Standard_Boolean                      myIsAllDegenerated;
453   Handle(IMeshData::MapOfOrientedEdges) myCouplesMap;
454   Handle(IMeshData::ListOfPnt2d)        myControlNodes;
455   const BRepMesh_CircleTool*            myCircles;
456 };
457
458 #endif