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