0026106: BRepMesh - revision of data model
[occt.git] / src / IVtkOCC / IVtkOCC_ShapeMesher.cxx
1 // Created on: 2011-10-14 
2 // Created by: Roman KOZLOV
3 // Copyright (c) 2011-2014 OPEN CASCADE SAS 
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 #include <Adaptor3d_IsoCurve.hxx>
17 #include <Bnd_Box.hxx>
18 #include <BRep_Tool.hxx>
19 #include <BRepBndLib.hxx>
20 #include <BRepMesh_DiscretFactory.hxx>
21 #include <BRepMesh_DiscretRoot.hxx>
22 #include <BRepTools.hxx>
23 #include <Hatch_Hatcher.hxx>
24 #include <GCPnts_QuasiUniformDeflection.hxx>
25 #include <GCPnts_TangentialDeflection.hxx>
26 #include <Geom_BezierSurface.hxx>
27 #include <Geom_BSplineSurface.hxx>
28 #include <Geom2dAdaptor_Curve.hxx>
29 #include <GeomAdaptor_Curve.hxx>
30 #include <gp_Dir2d.hxx>
31 #include <gp_Pnt2d.hxx>
32 #include <IVtkOCC_ShapeMesher.hxx>
33 #include <NCollection_Array1.hxx>
34 #include <Poly_Polygon3D.hxx>
35 #include <Poly_PolygonOnTriangulation.hxx>
36 #include <Poly_Triangulation.hxx>
37 #include <Precision.hxx>
38 #include <Prs3d.hxx>
39 #include <Prs3d_Drawer.hxx>
40 #include <Standard_ErrorHandler.hxx>
41 #include <TColgp_SequenceOfPnt2d.hxx>
42 #include <TColStd_Array1OfReal.hxx>
43 #include <TopExp.hxx>
44 #include <TopExp_Explorer.hxx>
45
46 IMPLEMENT_STANDARD_RTTIEXT(IVtkOCC_ShapeMesher,IVtk_IShapeMesher)
47
48 // Handle implementation
49
50
51 //================================================================
52 // Function : internalBuild
53 // Purpose  : 
54 //================================================================
55 void IVtkOCC_ShapeMesher::internalBuild()
56 {
57   // TODO: do we need any protection here so as not to triangualte
58   // the shape twice??? This can be done e.g. by checking if
59   // triangulation exists for TopoDS_Shape..
60   meshShape();
61
62   // Free vertices and free edges should always be shown.
63   // Shared edges are needed in WF representation only.
64   // TODO: how to filter free edges at visualization level????
65   addFreeVertices();
66   addEdges();
67
68   // Build wireframe points and cells (lines for isolines)
69   addWireFrameFaces();
70
71   // Build shaded representation (based on Poly_Triangulation)
72   addShadedFaces();
73 }
74
75 //================================================================
76 // Function : GetShapeObj
77 // Purpose  : 
78 //================================================================
79 const IVtkOCC_Shape::Handle IVtkOCC_ShapeMesher::GetShapeObj() const
80 {
81   return (IVtkOCC_Shape::Handle::DownCast(myShapeObj));
82 }
83
84 //================================================================
85 // Function : GetDeflection
86 // Purpose  : Returns absolute deflection used by this algorithm.
87 //================================================================
88 Standard_Real IVtkOCC_ShapeMesher::GetDeflection() const
89 {
90   if (myDeflection < Precision::Confusion()) // if not yet initialized
91   {
92     Handle(Prs3d_Drawer) aDefDrawer = new Prs3d_Drawer();
93     aDefDrawer->SetTypeOfDeflection (Aspect_TOD_RELATIVE);
94     aDefDrawer->SetDeviationCoefficient (GetDeviationCoeff());
95     myDeflection = Prs3d::GetDeflection (GetShapeObj()->GetShape(), aDefDrawer);
96   }
97
98   return myDeflection;
99 }
100
101 //================================================================
102 // Function : meshShape
103 // Purpose  : 
104 //================================================================
105 void IVtkOCC_ShapeMesher::meshShape()
106 {
107   const TopoDS_Shape& anOcctShape = GetShapeObj()->GetShape();
108   if (anOcctShape.IsNull())
109   {
110     return;
111   }
112
113   //Clean triangulation before compute incremental mesh
114   BRepTools::Clean (anOcctShape);
115
116   //Compute triangulation
117   Standard_Real aDeflection = GetDeflection();
118   if (aDeflection < Precision::Confusion())
119   {
120     return;
121   }
122
123   try
124   {
125     OCC_CATCH_SIGNALS
126
127     Handle(BRepMesh_DiscretRoot) anAlgo;
128     anAlgo = BRepMesh_DiscretFactory::Get().Discret (anOcctShape,
129                                                      aDeflection,
130                                                      GetDeviationAngle());
131     if (!anAlgo.IsNull())
132     {
133       anAlgo->Perform();
134     }
135   }
136   catch (Standard_Failure)
137   { }
138 }
139
140 //================================================================
141 // Function : addFreeVertices
142 // Purpose  : 
143 //================================================================
144 void IVtkOCC_ShapeMesher::addFreeVertices()
145 {
146   TopTools_IndexedDataMapOfShapeListOfShape aVertexMap;
147   TopExp::MapShapesAndAncestors (GetShapeObj()->GetShape(),
148                                  TopAbs_VERTEX,
149                                  TopAbs_EDGE,
150                                  aVertexMap);
151
152   Standard_Integer aVertNum = aVertexMap.Extent();
153   IVtk_MeshType aType;
154   for (Standard_Integer anIt = 1; anIt <= aVertNum; anIt++)
155   {
156     if (aVertexMap.FindFromIndex(anIt).IsEmpty())
157     {
158       aType = MT_FreeVertex;
159     }
160     else
161     {
162       aType = MT_SharedVertex;
163     }
164     const TopoDS_Vertex& aVertex = TopoDS::Vertex (aVertexMap.FindKey (anIt));
165     addVertex (aVertex, GetShapeObj()->GetSubShapeId (aVertex), aType);
166   }
167 }
168
169 //================================================================
170 // Function : addEdges
171 // Purpose  : 
172 //================================================================
173 void IVtkOCC_ShapeMesher::addEdges()
174 {
175   TopTools_IndexedDataMapOfShapeListOfShape anEdgesMap;
176   TopExp::MapShapesAndAncestors (GetShapeObj()->GetShape(), 
177                                  TopAbs_EDGE, 
178                                  TopAbs_FACE,
179                                  anEdgesMap);
180
181   int aNbFaces;
182   IVtk_MeshType aType;
183   myEdgesTypes.Clear();
184
185   TopExp_Explorer anEdgeIter (GetShapeObj()->GetShape(), TopAbs_EDGE);
186   for (; anEdgeIter.More(); anEdgeIter.Next())
187   {
188     const TopoDS_Edge& anOcctEdge = TopoDS::Edge (anEdgeIter.Current());
189     aNbFaces = anEdgesMap.FindFromKey (anOcctEdge).Extent();
190     if (aNbFaces == 0)
191     {
192       aType = MT_FreeEdge;
193     }
194     else if (aNbFaces == 1)
195     {
196       aType = MT_BoundaryEdge;
197     }
198     else
199     {
200       aType = MT_SharedEdge;
201     }
202     addEdge (anOcctEdge, GetShapeObj()->GetSubShapeId (anOcctEdge), aType);
203     myEdgesTypes.Bind (anOcctEdge, aType);
204   }
205 }
206
207 //================================================================
208 // Function : addWireFrameFaces
209 // Purpose  : 
210 //================================================================
211 void IVtkOCC_ShapeMesher::addWireFrameFaces()
212 {
213   // Check the deflection value once for all faces
214   if (GetDeflection() < Precision::Confusion())
215   {
216     return;
217   }
218
219   TopExp_Explorer aFaceIter (GetShapeObj()->GetShape(), TopAbs_FACE);
220   for (; aFaceIter.More(); aFaceIter.Next())
221   {
222     const TopoDS_Face& anOcctFace = TopoDS::Face (aFaceIter.Current());
223     try
224     {
225       OCC_CATCH_SIGNALS
226       addWFFace (anOcctFace, 
227                  GetShapeObj()->GetSubShapeId (anOcctFace));
228     }
229     catch (Standard_Failure)
230     { }
231   }
232 }
233
234 //================================================================
235 // Function : addShadedFaces
236 // Purpose  : 
237 //================================================================
238 void IVtkOCC_ShapeMesher::addShadedFaces()
239 {
240   TopExp_Explorer aFaceIter (GetShapeObj()->GetShape(), TopAbs_FACE);
241   for (; aFaceIter.More(); aFaceIter.Next())
242   {
243     const TopoDS_Face& anOcctFace = TopoDS::Face (aFaceIter.Current());
244     addShadedFace (anOcctFace,
245                    GetShapeObj()->GetSubShapeId (anOcctFace));
246   }
247 }
248
249 //================================================================
250 // Function : addVertex
251 // Purpose  : 
252 //================================================================
253 void IVtkOCC_ShapeMesher::addVertex (const TopoDS_Vertex& theVertex,
254                                      const IVtk_IdType    theShapeId,
255                                      const IVtk_MeshType  theMeshType)
256 {
257   if (theVertex.IsNull())
258   {
259     return;
260   }
261
262   gp_Pnt aPnt3d = BRep_Tool::Pnt (theVertex);
263
264   IVtk_PointId anId = 
265     myShapeData->InsertCoordinate (aPnt3d.X(), aPnt3d.Y(), aPnt3d.Z());
266   myShapeData->InsertVertex (theShapeId, anId, theMeshType);
267
268 }
269
270 //================================================================
271 // Function : processPolyline
272 // Purpose  : 
273 //================================================================
274 void IVtkOCC_ShapeMesher::processPolyline (Standard_Integer          theNbNodes,
275                                       const TColgp_Array1OfPnt&      thePoints,
276                                       const TColStd_Array1OfInteger& thePointIds,
277                                       const IVtk_IdType              theOcctId,
278                                       bool                           theNoTransform,
279                                       gp_Trsf                        theTransformation,
280                                       const IVtk_MeshType            theMeshType)
281 {
282   if (theNbNodes < 2)
283   {
284     return;
285   }
286
287   IVtk_PointIdList aPolyPointIds;
288
289   IVtk_PointId anId;
290   for (Standard_Integer aJ = 0; aJ < theNbNodes; aJ++)
291   {
292     Standard_Integer aPntId = thePointIds (aJ + 1);
293     gp_Pnt point = thePoints (aPntId);
294
295     if (!theNoTransform)
296     {
297       // Apply the transformation to points
298       point.Transform (theTransformation);
299     }
300
301     anId = myShapeData->InsertCoordinate (point.X(), point.Y(), point.Z());
302     aPolyPointIds.Append (anId);
303   }
304
305   myShapeData->InsertLine (theOcctId, &aPolyPointIds, theMeshType);
306 }
307
308 //================================================================
309 // Function : addEdge
310 // Purpose  : 
311 //================================================================
312 void IVtkOCC_ShapeMesher::addEdge (const TopoDS_Edge&  theEdge,
313                                    const IVtk_IdType   theShapeId,
314                                    const IVtk_MeshType theMeshType)
315 {
316   if (theEdge.IsNull() || BRep_Tool::Degenerated (theEdge))
317   {
318     return;
319   }
320
321   // Two discrete representations of an OCCT edge are possible:
322   // 1. Polygon on trinagulation - holds Ids of points
323   // contained in Poly_Triangulation object
324   Handle(Poly_PolygonOnTriangulation) aPolyOnTriangulation;
325   Handle(Poly_Triangulation) aTriangulation;
326   TopLoc_Location aLocation;
327   BRep_Tool::PolygonOnTriangulation (theEdge,
328                                      aPolyOnTriangulation,
329                                      aTriangulation,
330                                      aLocation,
331                                      1);
332
333   // 2. 3D polygon - holds 3D points
334   Handle(Poly_Polygon3D) aPoly3d;
335   if (aPolyOnTriangulation.IsNull())
336   {
337     aPoly3d = BRep_Tool::Polygon3D (theEdge, aLocation);
338   }
339
340   if (aPoly3d.IsNull() && aPolyOnTriangulation.IsNull())
341   {
342     return;
343   }
344
345   // Handle a non-identity transofmation applied to the edge
346   gp_Trsf anEdgeTransf;
347   bool noTransform = true;
348   if (!aLocation.IsIdentity())
349   {
350     noTransform = false;
351     anEdgeTransf = aLocation.Transformation();
352   }
353
354   if (!aPoly3d.IsNull())
355   {
356     Standard_Integer aNbNodes = aPoly3d->NbNodes();
357     const TColgp_Array1OfPnt& aPoints = aPoly3d->Nodes();
358     TColStd_Array1OfInteger aPointIds (1, aNbNodes);
359
360     for (Standard_Integer anI = 1; anI <= aNbNodes; anI++)
361     {
362       aPointIds.SetValue (anI, anI);
363     }
364
365     processPolyline (aNbNodes,
366                      aPoints,
367                      aPointIds,
368                      theShapeId,
369                      noTransform,
370                      anEdgeTransf,
371                      theMeshType);
372   }
373   else
374   {
375     Standard_Integer aNbNodes = aPolyOnTriangulation->NbNodes();
376     const TColStd_Array1OfInteger& aPointIds = aPolyOnTriangulation->Nodes();
377     const TColgp_Array1OfPnt& aPoints = aTriangulation->Nodes();
378
379     processPolyline (aNbNodes,
380                      aPoints,
381                      aPointIds,
382                      theShapeId,
383                      noTransform,
384                      anEdgeTransf,
385                      theMeshType);
386   }
387 }
388
389
390 //================================================================
391 // Function : FindLimits
392 // Purpose  : Static internal function, finds parametrical limits of the curve.
393 //! @param [in] theCurve 3D curve adaptor used to retrieve the curve geometry
394 //! @param [in] theLimit maximum allowed absolute parameter value
395 //! @param [out] theFirst minimum parameter value for the curve
396 //! @param [out] theLast maximum parameter value for the curve
397 //================================================================
398 static void FindLimits (const Adaptor3d_Curve& theCurve,
399                         const Standard_Real&   theLimit,
400                         Standard_Real&         theFirst,
401                         Standard_Real&         theLast)
402 {
403   theFirst = Max(theCurve.FirstParameter(), theFirst);
404   theLast  = Min(theCurve.LastParameter(), theLast);
405   Standard_Boolean isFirstInf = Precision::IsNegativeInfinite (theFirst);
406   Standard_Boolean isLastInf  = Precision::IsPositiveInfinite (theLast);
407
408   if (isFirstInf || isLastInf)
409   {
410     gp_Pnt aP1, aP2;
411     Standard_Real aDelta = 1;
412     if (isFirstInf && isLastInf)
413     {
414       do
415       {
416         aDelta *= 2;
417         theFirst = - aDelta;
418         theLast  =   aDelta;
419         theCurve.D0 (theFirst, aP1);
420         theCurve.D0 (theLast, aP2);
421       } while (aP1.Distance (aP2) < theLimit);
422     }
423     else if (isFirstInf)
424     {
425       theCurve.D0 (theLast, aP2);
426       do {
427         aDelta *= 2;
428         theFirst = theLast - aDelta;
429         theCurve.D0 (theFirst, aP1);
430       } while (aP1.Distance(aP2) < theLimit);
431     }
432     else if (isLastInf)
433     {
434       theCurve.D0 (theFirst, aP1);
435       do
436       {
437         aDelta *= 2;
438         theLast = theFirst + aDelta;
439         theCurve.D0 (theLast, aP2);
440       } while (aP1.Distance (aP2) < theLimit);
441     }
442   }
443 }
444
445 //================================================================
446 // Function : FindLimits
447 // Purpose  :Static helper function, builds a discrete representation
448 //! (sequence of points) for the given curve.
449 //!
450 //! @param [in] theCurve 3D curve adaptor used to retrieve the curve geometry
451 //! @param [in] theDeflection absolute deflection value
452 //! @param [in] theAngle deviation angle value
453 //! @param [in] theU1 minimal curve parameter value
454 //! @param [in] theU2 maximal curve parameter value
455 //! @param [out] thePoints the container for generated polyline
456 //================================================================
457 static void DrawCurve (Adaptor3d_Curve&    theCurve,
458                        const Standard_Real theDeflection,
459                        const Standard_Real theAngle,
460                        const Standard_Real theU1,
461                        const Standard_Real theU2,
462                        IVtk_Polyline&      thePoints)
463 {
464   switch (theCurve.GetType())
465   {
466   case GeomAbs_Line:
467     {
468       gp_Pnt aPnt = theCurve.Value(theU1);
469       thePoints.Append (aPnt);
470
471       aPnt = theCurve.Value(0.5 * (theU1 + theU2));
472       thePoints.Append (aPnt);
473
474       aPnt = theCurve.Value (theU2);
475       thePoints.Append(aPnt);
476    }
477     break;
478   default:
479     {
480       Standard_Integer aNbInter = theCurve.NbIntervals(GeomAbs_C1);
481       Standard_Integer anI, aJ;
482       TColStd_Array1OfReal aParams(1, aNbInter+1);
483       theCurve.Intervals(aParams, GeomAbs_C1);
484       Standard_Real anU1, anU2;
485       Standard_Integer NumberOfPoints;
486
487       for (aJ = 1; aJ <= aNbInter; aJ++)
488       {
489         anU1 = aParams (aJ); anU2 = aParams (aJ + 1);
490         if (anU2 > anU1 && anU1 < anU2)
491         {
492           anU1 = Max(anU1, anU1);
493           anU2 = Min(anU2, anU2);
494
495           GCPnts_TangentialDeflection anAlgo (theCurve, anU1, anU2, theAngle, theDeflection);
496           NumberOfPoints = anAlgo.NbPoints();
497
498           if (NumberOfPoints > 0)
499           {
500             for (anI = 1; anI < NumberOfPoints; anI++)
501             {
502               thePoints.Append(anAlgo.Value (anI));
503             }
504             if (aJ == aNbInter)
505             {
506               thePoints.Append (anAlgo.Value (NumberOfPoints));
507             }
508           }
509         }
510       }
511     }
512   }
513 }
514
515 //================================================================
516 // Function : buildIsoLines
517 // Purpose  : 
518 //================================================================
519 void IVtkOCC_ShapeMesher::buildIsoLines (const Handle(BRepAdaptor_HSurface)& theFace,
520                                          const Standard_Boolean          theIsDrawUIso,
521                                          const Standard_Boolean          theIsDrawVIso,
522                                          const Standard_Integer          theNBUiso,
523                                          const Standard_Integer          theNBViso,
524                                          IVtk_PolylineList&              thePolylines)
525 {
526   Standard_Real anUF, anUL, aVF, aVL;
527   anUF = theFace->FirstUParameter();
528   anUL = theFace->LastUParameter();
529   aVF = theFace->FirstVParameter();
530   aVL = theFace->LastVParameter();
531
532   // Restrict maximal parameter value
533   // in OCCT it's 5e+5 by default
534   const Standard_Real aLimit = 5e+5;
535
536   // compute bounds of the restriction
537   Standard_Real anUMin, anUMax, aVMin, aVMax;
538   Standard_Integer anI;
539
540   anUMin = Max (anUF, -aLimit);
541   anUMax = Min (anUL, aLimit);
542   aVMin = Max (aVF, -aLimit);
543   aVMax = Min (aVL, aLimit);
544
545   // update min max for the hatcher.
546   gp_Pnt2d aP1,aP2;
547   gp_Pnt aDummyPnt;
548
549   Standard_Real aDdefle = Max (anUMax - anUMin, aVMax - aVMin) * GetDeviationCoeff();
550   TColgp_SequenceOfPnt2d aTabPoints;
551
552   anUMin = aVMin = 1.e100;
553   anUMax = aVMax = -1.e100;
554
555   // Process the edges
556   TopExp_Explorer aToolRst;
557   TopoDS_Face aTopoFace (((BRepAdaptor_Surface*)&(theFace->Surface()))->Face());
558   for (aToolRst.Init (aTopoFace, TopAbs_EDGE); aToolRst.More(); aToolRst.Next())
559   {
560     TopAbs_Orientation anOrient = aToolRst.Current().Orientation();
561     // Skip INTERNAL and EXTERNAL edges
562     if (anOrient == TopAbs_FORWARD || anOrient == TopAbs_REVERSED)
563     {
564       Standard_Real anU1, anU2;
565       const Handle(Geom2d_Curve)& aCurve = 
566         BRep_Tool::CurveOnSurface (TopoDS::Edge (aToolRst.Current()),
567                                    aTopoFace,
568                                    anU1, anU2);
569       if (aCurve.IsNull())
570       {
571         continue;
572       }
573
574       Geom2dAdaptor_Curve aRCurve;
575       aRCurve.Load (aCurve, anU1, anU2);
576       if (aRCurve.GetType() != GeomAbs_Line)
577       {
578         GCPnts_QuasiUniformDeflection aUDP(aRCurve, aDdefle);
579         if (aUDP.IsDone())
580         {
581           Standard_Integer NumberOfPoints = aUDP.NbPoints();
582           if ( NumberOfPoints >= 2 )
583           {
584             aDummyPnt = aUDP.Value (1);
585             aP2.SetCoord (aDummyPnt.X(), aDummyPnt.Y());
586             anUMin = Min (aP2.X(), anUMin);
587             anUMax = Max (aP2.X(), anUMax);
588             aVMin = Min (aP2.Y(), aVMin);
589             aVMax = Max (aP2.Y(), aVMax);
590             for (anI = 2; anI <= NumberOfPoints; anI++)
591             {
592               aP1 = aP2;
593               aDummyPnt = aUDP.Value (anI);
594               aP2.SetCoord (aDummyPnt.X(), aDummyPnt.Y());
595               anUMin = Min(aP2.X(), anUMin);
596               anUMax = Max(aP2.X(), anUMax);
597               aVMin = Min(aP2.Y(), aVMin);
598               aVMax = Max(aP2.Y(), aVMax);
599
600               if(anOrient == TopAbs_FORWARD )
601               {
602                 //isobuild.Trim(P1,P2);
603                 aTabPoints.Append (aP1);
604                 aTabPoints.Append (aP2);
605               }
606               else
607               {
608                 //isobuild.Trim(P2,P1);
609                 aTabPoints.Append (aP2);
610                 aTabPoints.Append (aP1);
611               }
612             }
613           }
614         }
615         else
616         {
617           cout << "Cannot evaluate curve on surface"<<endl;
618         }
619       }
620       else
621       {
622         anU1 = aRCurve.FirstParameter();
623         anU2 = aRCurve.LastParameter();
624         // MSV 17.08.06 OCC13144: U2 occured less than U1, to overcome it
625         // ensure that distance U2-U1 is not greater than aLimit*2,
626         // if greater then choose an origin and use aLimit to define
627         // U1 and U2 anew
628         Standard_Real aOrigin = 0.;
629         if (!Precision::IsNegativeInfinite(anU1) || !Precision::IsPositiveInfinite (anU2))
630         {
631           if (Precision::IsNegativeInfinite (anU1))
632           {
633             aOrigin = anU2 - aLimit;
634           }
635           else if (Precision::IsPositiveInfinite (anU2))
636           {
637             aOrigin = anU1 + aLimit;
638           }
639           else
640           {
641             aOrigin = (anU1 + anU2) * 0.5;
642           }
643         }
644
645         anU1 = Max (aOrigin - aLimit, anU1);
646         anU2 = Min (aOrigin + aLimit, anU2);
647         aP1 = aRCurve.Value (anU1);
648         aP2 = aRCurve.Value (anU2);
649         anUMin = Min(aP1.X(), anUMin);
650         anUMax = Max(aP1.X(), anUMax);
651         aVMin = Min(aP1.Y(), aVMin);
652         aVMax = Max(aP1.Y(), aVMax);
653         anUMin = Min(aP2.X(), anUMin);
654         anUMax = Max(aP2.X(), anUMax);
655         aVMin = Min(aP2.Y(), aVMin);
656         aVMax = Max(aP2.Y(), aVMax);
657         if(anOrient == TopAbs_FORWARD )
658         {
659          // isobuild.Trim(P1,P2);
660           aTabPoints.Append (aP1);
661           aTabPoints.Append (aP2);
662         }
663         else
664         {
665           //isobuild.Trim(P2,P1);
666           aTabPoints.Append (aP2);
667           aTabPoints.Append (aP1);
668         }
669       }
670     }
671   }
672
673   // load the isos
674   const Standard_Real anIntersectionTolerance = 1.e-5;
675   Hatch_Hatcher anIsoBuild (anIntersectionTolerance, Standard_True );
676
677   Standard_Boolean isUClosed = theFace->IsUClosed();
678   Standard_Boolean isVClosed = theFace->IsVClosed();
679
680   if (!isUClosed)
681   {
682     anUMin = anUMin + (anUMax - anUMin) / 1000.0;
683     anUMax = anUMax - (anUMax - anUMin) /1000.0;
684   }
685
686   if (!isVClosed)
687   {
688     aVMin = aVMin + (aVMax - aVMin) /1000.0;
689     aVMax = aVMax - (aVMax - aVMin) /1000.0;
690   }
691
692   if (theIsDrawUIso)
693   {
694     if (theNBUiso > 0)
695     {
696       isUClosed = Standard_False;
697       Standard_Real aDu= isUClosed ? (anUMax - anUMin) / theNBUiso : (anUMax - anUMin) / (1 + theNBUiso);
698       for (anI = 1; anI <= theNBUiso; anI++)
699       {
700         anIsoBuild.AddXLine (anUMin + aDu*anI);
701       }
702     }
703   }
704   if (theIsDrawVIso)
705   {
706     if (theNBViso > 0)
707     {
708       isVClosed = Standard_False;
709       Standard_Real aDv= isVClosed ? (aVMax - aVMin) / theNBViso : (aVMax - aVMin) / (1 + theNBViso);
710       for (anI = 1; anI <= theNBViso; anI++)
711       {
712         anIsoBuild.AddYLine (aVMin + aDv*anI);
713       }
714     }
715   }
716
717   Standard_Integer aLength = aTabPoints.Length();
718   for (anI = 1; anI <= aLength; anI += 2)
719   {
720     anIsoBuild.Trim (aTabPoints (anI),aTabPoints (anI + 1));
721   }
722
723   // Create the polylines for isos
724   Adaptor3d_IsoCurve anIso;
725   anIso.Load(theFace);
726   Handle(Geom_Curve) aBCurve;
727   const BRepAdaptor_Surface& aBSurf = *(BRepAdaptor_Surface*)&(theFace->Surface());
728   GeomAbs_SurfaceType aType = theFace->GetType();
729
730   Standard_Integer aNumberOfLines = anIsoBuild.NbLines();
731   Handle(Geom_Surface) aGeomSurf;
732   if (aType == GeomAbs_BezierSurface)
733   {
734     aGeomSurf = aBSurf.Bezier();
735   }
736   else if (aType == GeomAbs_BSplineSurface)
737   {
738     aGeomSurf = aBSurf.BSpline();
739   }
740
741   Standard_Real aDeflection = GetDeflection();
742   Standard_Real anAngle     = GetDeviationAngle();
743   for (anI = 1; anI <= aNumberOfLines; anI++)
744   {
745     Standard_Integer aNumberOfIntervals = anIsoBuild.NbIntervals(anI);
746     Standard_Real aCoord = anIsoBuild.Coordinate(anI);
747     for (Standard_Integer aJ = 1; aJ <= aNumberOfIntervals; aJ++)
748     {
749       Standard_Real aB1 = anIsoBuild.Start (anI, aJ);
750       Standard_Real aB2 = anIsoBuild.End(anI, aJ);
751
752       if (!aGeomSurf.IsNull())
753       {
754         if (anIsoBuild.IsXLine (anI))
755         {
756           aBCurve = aGeomSurf->UIso (aCoord);
757         }
758         else
759         {
760           aBCurve = aGeomSurf->VIso (aCoord);
761         }
762
763         GeomAdaptor_Curve aGeomCurve (aBCurve);
764         FindLimits (aGeomCurve, aLimit, aB1, aB2);
765         if (aB2 - aB1 > Precision::Confusion())
766         {
767           IVtk_Polyline aPoints;
768           DrawCurve (aGeomCurve, aDeflection, anAngle, aB1, aB2, aPoints);
769           thePolylines.Append (aPoints);
770         }
771       }
772       else
773       {
774         if (anIsoBuild.IsXLine (anI))
775         {
776           anIso.Load (GeomAbs_IsoU, aCoord, aB1, aB2);
777         }
778         else
779         {
780           anIso.Load (GeomAbs_IsoV, aCoord, aB1, aB2);
781         }
782         FindLimits (anIso, aLimit, aB1, aB2);
783         if (aB2 - aB1>Precision::Confusion())
784         {
785           IVtk_Polyline aPoints;
786           DrawCurve (anIso, aDeflection, anAngle, aB1, aB2, aPoints);
787           thePolylines.Append (aPoints);
788         }
789       }
790     }
791   }
792 }
793
794 //================================================================
795 // Function : addWFFace
796 // Purpose  : 
797 //================================================================
798 void IVtkOCC_ShapeMesher::addWFFace (const TopoDS_Face& theFace,
799                                      const IVtk_IdType  theShapeId)
800 {
801   if (theFace.IsNull())
802   {
803     return;
804   }
805
806   TopoDS_Face aFaceToMesh = theFace;
807   aFaceToMesh.Orientation (TopAbs_FORWARD);
808
809   // The code that builds wireframe representation for a TopoDS_Face
810   // has been adapted from some OCCT 6.5.1 methods:
811   // - Prs3d_WFShape::Add()
812   // - StdPrs_WFDeflectionRestrictedFace::Add()
813   // - StdPrs_DeflectionCurve::Add()
814
815   // Add face's edges here but with the face ID
816   TopExp_Explorer anEdgeIter (aFaceToMesh, TopAbs_EDGE );
817   for (; anEdgeIter.More(); anEdgeIter.Next())
818   {
819     const TopoDS_Edge& anOcctEdge = TopoDS::Edge (anEdgeIter.Current());
820     addEdge (anOcctEdge, theShapeId, myEdgesTypes (anOcctEdge));
821   }
822
823   TopLoc_Location aLoc;
824   const Handle(Geom_Surface)& aGeomSurf = BRep_Tool::Surface (aFaceToMesh, aLoc);
825   if (aGeomSurf.IsNull())
826   {
827     return;
828   }
829
830   BRepAdaptor_Surface aSurf;
831   aSurf.Initialize (aFaceToMesh);
832   Handle(BRepAdaptor_HSurface) aSurfAdaptor = new BRepAdaptor_HSurface (aSurf);
833
834   IVtk_PolylineList aPolylines;
835   gp_Trsf aDummyTrsf;
836
837   // Building U isolines
838   // Introducing a local scope here to simplify variable naming
839   {
840     buildIsoLines (aSurfAdaptor,
841                    myNbIsos[0] != 0,
842                    Standard_False,
843                    myNbIsos[0],
844                    0,
845                    aPolylines);
846
847     IVtk_PolylineList::Iterator anIt (aPolylines);
848     for (; anIt.More(); anIt.Next())
849     {
850       const IVtk_Polyline& aPntSeq = anIt.Value();
851       Standard_Integer aNbNodes = aPntSeq.Length();
852       TColgp_Array1OfPnt aPoints (1, aNbNodes);
853       for (Standard_Integer aJ = 1; aJ <= aNbNodes; aJ++)
854       {
855         aPoints.SetValue (aJ, aPntSeq.Value(aJ));
856       }
857
858       TColStd_Array1OfInteger aPointIds (1, aNbNodes);
859       for (Standard_Integer anI = 1; anI <= aNbNodes; anI++)
860       {
861         aPointIds.SetValue (anI, anI);
862       }
863
864       processPolyline (aNbNodes,
865                        aPoints,
866                        aPointIds,
867                        theShapeId,
868                        true,
869                        aDummyTrsf, 
870                        MT_IsoLine);
871     }
872   }
873
874   // Building V isolines
875   {
876     aPolylines.Clear();
877     buildIsoLines (aSurfAdaptor,
878                    Standard_False,
879                    myNbIsos[1] != 0,
880                    0,
881                    myNbIsos[1],
882                    aPolylines);
883
884     IVtk_PolylineList::Iterator anIt (aPolylines);
885     for (; anIt.More(); anIt.Next())
886     {
887       const IVtk_Polyline& aPntSeq = anIt.Value();
888       Standard_Integer aNbNodes = aPntSeq.Length();
889       TColgp_Array1OfPnt aPoints (1, aNbNodes);
890       for (int aJ = 1; aJ <= aNbNodes; aJ++)
891       {
892         aPoints.SetValue (aJ, aPntSeq.Value (aJ));
893       }
894
895       TColStd_Array1OfInteger aPointIds (1, aNbNodes);
896       for (Standard_Integer anI = 1; anI <= aNbNodes; anI++)
897       {
898         aPointIds.SetValue (anI, anI);
899       }
900
901       processPolyline (aNbNodes,
902                        aPoints,
903                        aPointIds,
904                        theShapeId,
905                        true,
906                        aDummyTrsf, 
907                        MT_IsoLine);
908     }
909   }
910 }
911
912 //================================================================
913 // Function : addShadedFace
914 // Purpose  : 
915 //================================================================
916 void IVtkOCC_ShapeMesher::addShadedFace (const TopoDS_Face& theFace,
917                                          const IVtk_IdType  theShapeId)
918 {
919   if (theFace.IsNull())
920   {
921     return;
922   }
923
924   // Build triangulation of the face.
925   TopLoc_Location aLoc;
926   Handle(Poly_Triangulation) anOcctTriangulation = BRep_Tool::Triangulation (theFace, aLoc);
927   if (anOcctTriangulation.IsNull())
928   {
929     return;
930   }
931
932   gp_Trsf aPntTransform;
933   Standard_Boolean noTransform = Standard_True;
934   if (!aLoc.IsIdentity())
935   {
936     noTransform = Standard_False;
937     aPntTransform = aLoc.Transformation();
938   }
939
940   // Get triangulation points.
941   const TColgp_Array1OfPnt& aPoints = anOcctTriangulation->Nodes();
942   Standard_Integer aNbPoints = anOcctTriangulation->NbNodes();
943
944   // Keep inserted points id's of triangulation in an array.
945   NCollection_Array1<IVtk_PointId> aPointIds (1, aNbPoints);
946   IVtk_PointId anId;
947
948   Standard_Integer anI;
949   for (anI = 1; anI <= aNbPoints; anI++)
950   {
951     gp_Pnt aPoint = aPoints (anI);
952
953     if (!noTransform)
954     {
955       aPoint.Transform (aPntTransform);
956     }
957
958     // Add a point into output shape data and keep its id in the array.
959     anId = myShapeData->InsertCoordinate (aPoint.X(), aPoint.Y(), aPoint.Z());
960     aPointIds.SetValue (anI, anId);
961   }
962
963   // Create triangles on the created triangulation points.
964   const Poly_Array1OfTriangle& aTriangles = anOcctTriangulation->Triangles();
965   Standard_Integer aNbTriangles = anOcctTriangulation->NbTriangles();
966   Standard_Integer aN1, aN2, aN3;
967   for (anI = 1; anI <= aNbTriangles; anI++)
968   {
969     aTriangles(anI).Get (aN1, aN2, aN3); // get indexes of triangle's points
970     // Insert new triangle on these points into output shape data.
971     myShapeData->InsertTriangle (
972       theShapeId, aPointIds(aN1), aPointIds(aN2), aPointIds(aN3), MT_ShadedFace);
973   }
974 }