]> SALOME platform Git repositories - modules/geom.git/blob - src/OBJECT/GEOM_OCCReader.cxx
Salome HOME
bos #29484 Merge branch 'vsr/29484'
[modules/geom.git] / src / OBJECT / GEOM_OCCReader.cxx
1 // Copyright (C) 2007-2021  CEA/DEN, EDF R&D, OPEN CASCADE
2 //
3 // Copyright (C) 2003-2007  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
4 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
5 //
6 // This library is free software; you can redistribute it and/or
7 // modify it under the terms of the GNU Lesser General Public
8 // License as published by the Free Software Foundation; either
9 // version 2.1 of the License, or (at your option) any later version.
10 //
11 // This library is distributed in the hope that it will be useful,
12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
14 // Lesser General Public License for more details.
15 //
16 // You should have received a copy of the GNU Lesser General Public
17 // License along with this library; if not, write to the Free Software
18 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
19 //
20 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
21 //
22
23 //  GEOM OBJECT : interactive object for Geometry entities visualization
24 //  File   : GEOM_OCCReader.h
25 //  Author : Christophe ATTANASIO
26 //  Module : GEOM
27
28 #include "GEOM_OCCReader.h"
29
30 #include <GEOMUtils_Hatcher.hxx>
31
32 // VTK Includes
33 #include <vtkPoints.h>
34 #include <vtkCellArray.h>
35
36 #include <vtkObjectFactory.h>
37 #include <vtkPolyData.h>
38 #include <vtkInformation.h>
39 #include <vtkInformationVector.h>
40
41 // OpenCASCADE Includes
42 #include <Adaptor3d_HCurve.hxx>
43 #include <Poly_Triangulation.hxx>
44 #include <Poly_Polygon3D.hxx>
45 #include <Poly_PolygonOnTriangulation.hxx>
46 #include <TopoDS.hxx>
47 #include <TopoDS_Face.hxx>
48 #include <TopoDS_Edge.hxx>
49 #include <Precision.hxx>
50 #include <BRep_Tool.hxx>
51 #include <TColStd_Array1OfInteger.hxx>
52 #include <TColStd_Array1OfReal.hxx>
53
54 #include "utilities.h"
55
56
57 static Standard_Integer lastVTKpoint = 0;
58 static Standard_Integer PlotCount = 0;
59 static Standard_Real IsoRatio = 1.001;
60 static Standard_Integer MaxPlotCount = 5; 
61
62 //=======================================================================
63 // Function : New
64 // Purpose  : 
65 //=======================================================================
66
67 GEOM_OCCReader* GEOM_OCCReader::New()
68 {
69   vtkObject* ret = vtkObjectFactory::CreateInstance("GEOM_OCCReader");
70   if(ret) {
71     return (GEOM_OCCReader*)ret;
72   }
73   return new GEOM_OCCReader;
74 }
75
76 //=======================================================================
77 // Function : GEOM_OCCReader
78 // Purpose  : 
79 //=======================================================================
80
81 GEOM_OCCReader::GEOM_OCCReader()
82 {
83   //this->myShape = NULL;
84   this->amode = 0;
85   this->forced = Standard_False;
86   this->discretiso = 15;
87   this->nbisos = 1;
88 }
89 //=======================================================================
90 // Function : ~GEOM_OCCReader
91 // Purpose  : 
92 //=======================================================================
93
94 GEOM_OCCReader::~GEOM_OCCReader()
95 {
96 }
97
98
99 //=======================================================================
100 // Function : RequestData
101 // Purpose  : 
102 //=======================================================================
103
104 int GEOM_OCCReader::RequestData(vtkInformation *vtkNotUsed(request),
105                                 vtkInformationVector **/*inputVector*/,
106                                 vtkInformationVector *outputVector)
107 {
108   vtkInformation *outInfo = outputVector->GetInformationObject(0);
109   vtkPolyData *output = vtkPolyData::SafeDownCast(
110     outInfo->Get(vtkDataObject::DATA_OBJECT()));
111
112   vtkPoints* Pts = NULL;
113   vtkCellArray* Cells = NULL;
114   TopLoc_Location aLoc;
115
116   // Allocation
117   Pts = vtkPoints::New();
118   Cells = vtkCellArray::New();
119         
120   //Compute number of triangles and points
121   Standard_Integer nbpoly=0,nbpts=0;
122
123   if(amode==1) {
124     //for shading
125     
126     if(myShape.ShapeType() == TopAbs_FACE) {
127       // whole FACE 
128       const TopoDS_Face& aFace = TopoDS::Face(myShape);
129       Handle(Poly_Triangulation) aPoly = BRep_Tool::Triangulation(aFace,aLoc);
130       if(aPoly.IsNull()) {
131         Pts->Delete();
132         Cells->Delete();
133         return 0;
134       }
135
136       nbpts = aPoly->NbNodes();
137       nbpoly = aPoly->NbTriangles();
138
139       Pts->SetNumberOfPoints(nbpts);
140       Cells->Allocate(Cells->EstimateSize(nbpoly,3));
141     }
142     else { 
143         Cells->Delete();
144         Pts->Delete();
145         return 0; 
146     }
147   }
148
149   // Start computation
150   if(amode == 0) {
151     ComputeWireframe(Pts,Cells);
152     output->SetPoints(Pts);
153     output->SetLines(Cells);
154     output->Squeeze();
155   }
156   else {
157     if(myShape.ShapeType() == TopAbs_FACE) {
158       ComputeShading(Pts,Cells);
159
160       output->SetPoints(Pts);
161       output->SetPolys(Cells);
162       output->Squeeze();
163     }
164   }
165   Pts->Delete();
166   Cells->Delete();
167   return 1;
168 }
169
170 //=======================================================================
171 // Function : ComputeWireframe
172 // Purpose  : Compute the shape in CAD wireframe mode
173 //=======================================================================
174
175 void GEOM_OCCReader::ComputeWireframe(vtkPoints* Pts,vtkCellArray* Cells){
176
177   // Check the type of the shape:
178   if(myShape.ShapeType() == TopAbs_FACE) {
179     // Face
180     TransferFaceWData(TopoDS::Face(myShape),Pts,Cells);
181   } else if(myShape.ShapeType() == TopAbs_EDGE) {
182     // Edge
183     TransferEdgeWData(TopoDS::Edge(myShape),Pts,Cells);
184   } else {
185     if(myShape.ShapeType() == TopAbs_VERTEX) {
186       // Vertex
187       TransferVertexWData(TopoDS::Vertex(myShape),Pts,Cells);
188     }
189   }
190 }
191
192 //=======================================================================
193 // Function : TransferFaceWData
194 // Purpose  : Transfer wireframe data for FACE
195 //=======================================================================
196
197 void GEOM_OCCReader::TransferFaceWData(const TopoDS_Face& aFace,
198                                          vtkPoints* Pts,
199                                          vtkCellArray* Cells) 
200 {
201   TopoDS_Face aCopyFace = aFace; 
202   aCopyFace.Orientation (TopAbs_FORWARD);
203   createISO(aCopyFace,1,Pts,Cells);
204 }
205
206 //=======================================================================
207 // Function : createISO
208 // Purpose  : Create ISO for Face Wireframe representation 
209 //=======================================================================
210
211 void GEOM_OCCReader::createISO (const TopoDS_Face      &TopologicalFace,
212                                 const Standard_Integer  NbIsos,
213                                       vtkPoints        *Pts,
214                                       vtkCellArray     *Cell)
215 {
216   GEOMUtils::Hatcher aHatcher(TopologicalFace);
217
218   aHatcher.Init(NbIsos);
219   aHatcher.Perform();
220
221   if (aHatcher.IsDone()) {
222     // Push iso lines in vtk kernel
223     Standard_Integer  pt_start_idx = 0;
224
225     createIsos(aHatcher, Standard_True, pt_start_idx, Pts, Cell);
226     createIsos(aHatcher, Standard_False, pt_start_idx, Pts, Cell);
227   }
228 }
229
230 //=======================================================================
231 // Function : createIsos
232 // Purpose  : Create isolines obtained from hatcher.
233 //=======================================================================
234 void GEOM_OCCReader::createIsos(const GEOMUtils::Hatcher &theHatcher,
235                                 const Standard_Boolean   IsUIso,
236                                       Standard_Integer  &pt_start_idx,
237                                       vtkPoints         *Pts,
238                                       vtkCellArray      *Cell)
239 {
240   // Push iso lines in vtk kernel
241   Handle(TColStd_HArray1OfInteger) anIndices;
242   Handle(TColStd_HArray1OfReal)    aParams;
243
244   if (IsUIso) {
245     // U-isolines
246     anIndices = theHatcher.GetUIndices();
247     aParams   = theHatcher.GetUParams();
248   } else {
249     // V-isolines
250     anIndices = theHatcher.GetVIndices();
251     aParams   = theHatcher.GetVParams();
252   }
253
254   if (anIndices.IsNull() || aParams.IsNull()) {
255     if (IsUIso) {
256       MESSAGE("GEOMUtils_Hatcher: null U-isoline indices")
257     } else {
258       MESSAGE("GEOMUtils_Hatcher: null V-isoline indices")
259     }
260   } else {
261     const GeomAbs_IsoType aType    = (IsUIso ? GeomAbs_IsoU : GeomAbs_IsoV);
262     Standard_Integer      anIsoInd = anIndices->Lower();
263
264     for (; anIsoInd <= anIndices->Upper(); anIsoInd++) {
265       const Standard_Integer aHatchingIndex = anIndices->Value(anIsoInd);
266
267       if (aHatchingIndex != 0) {
268         const Standard_Real    aParam     = aParams->Value(anIsoInd);
269         const Standard_Integer aNbDomains =
270           theHatcher.GetNbDomains(aHatchingIndex);
271
272         if (aNbDomains < 0) {
273           if (IsUIso) {
274             MESSAGE("GEOMUtils_Hatcher: U iso of parameter: "<<aParam)
275           } else {
276             MESSAGE("GEOMUtils_Hatcher: V iso of parameter: "<<aParam)
277           }
278
279           switch (theHatcher.GetHatcher().Status (aHatchingIndex)) {
280           case HatchGen_NoProblem          :
281             MESSAGE("No Problem")          ; break ;
282           case HatchGen_TrimFailure        :
283             MESSAGE("Trim Failure")        ; break ;
284           case HatchGen_TransitionFailure  :
285             MESSAGE("Transition Failure")  ; break ;
286           case HatchGen_IncoherentParity   :
287             MESSAGE("Incoherent Parity")   ; break ;
288           case HatchGen_IncompatibleStates :
289             MESSAGE("Incompatible States") ; break ;
290           }
291         } else {
292           Standard_Integer anIDom = 1;
293           Standard_Real    aV1;
294           Standard_Real    aV2;
295
296           for (; anIDom <= aNbDomains; anIDom++) {
297             if (theHatcher.GetDomain(aHatchingIndex, anIDom, aV1, aV2)) {
298               DrawIso(aType, aParam, aV1, aV2, Pts, Cell,pt_start_idx);
299             }
300           }
301         }
302       }
303     }
304   }
305 }
306
307 //=======================================================================
308 // Function : MoveTo
309 // Purpose  : Init VTK ISO PLOT
310 //=======================================================================
311 void GEOM_OCCReader::MoveTo(gp_Pnt P,
312                               vtkPoints* Pts)
313 {    
314   float coord[3];
315
316   coord[0] = P.X(); coord[1] = P.Y(); coord[2] = P.Z();
317   lastVTKpoint = Pts->InsertNextPoint(coord);
318     
319
320
321 //=======================================================================
322 // Function : DrawTo
323 // Purpose  : Plot point in VTK
324 //=======================================================================
325 void GEOM_OCCReader::DrawTo(gp_Pnt P,
326                               vtkPoints* Pts,
327                               vtkCellArray* Cells)
328 {
329   float coord[3];
330   coord[0] = P.X(); coord[1] = P.Y(); coord[2] = P.Z();
331   Standard_Integer NewVTKpoint =  Pts->InsertNextPoint(coord);
332
333   vtkIdType pts[2];
334   pts[0] = lastVTKpoint;
335   pts[1] = NewVTKpoint;
336
337   Cells->InsertNextCell(2,pts);
338     
339   lastVTKpoint = NewVTKpoint;
340 }
341
342
343 //=======================================================================
344 // Function : DrawIso
345 // Purpose  : Draw an iso on vtk
346 //=======================================================================
347 void GEOM_OCCReader::DrawIso(GeomAbs_IsoType T, 
348                                Standard_Real Par, 
349                                Standard_Real T1,
350                                Standard_Real T2,
351                                vtkPoints* Pts,
352                                vtkCellArray* Cells,
353                                Standard_Integer& /*startidx*/)
354 {
355
356   Standard_Boolean halt = Standard_False;
357   Standard_Integer j,myDiscret = discretiso;
358   Standard_Real U1,U2,V1,V2,stepU=0.,stepV=0.;
359   gp_Pnt P;
360   TopLoc_Location l;
361
362   const Handle(Geom_Surface)& S = BRep_Tool::Surface(TopoDS::Face(myShape),l);
363   if (!S.IsNull()) {
364     BRepAdaptor_Surface S(TopoDS::Face(myShape),Standard_False);
365       
366     GeomAbs_SurfaceType SurfType = S.GetType();
367
368     GeomAbs_CurveType CurvType = GeomAbs_OtherCurve;
369
370     Standard_Integer Intrv, nbIntv;
371     Standard_Integer nbUIntv = S.NbUIntervals(GeomAbs_CN);
372     Standard_Integer nbVIntv = S.NbVIntervals(GeomAbs_CN);
373     TColStd_Array1OfReal TI(1,Max(nbUIntv, nbVIntv)+1);
374
375
376     if (T == GeomAbs_IsoU) {
377       S.VIntervals(TI, GeomAbs_CN);
378       V1 = Max(T1, TI(1));
379       V2 = Min(T2, TI(2));
380       U1 = Par;
381       U2 = Par;
382       stepU = 0;
383       nbIntv = nbVIntv;
384     }
385     else {
386       S.UIntervals(TI, GeomAbs_CN);
387       U1 = Max(T1, TI(1));
388       U2 = Min(T2, TI(2));
389       V1 = Par;
390       V2 = Par;
391       stepV = 0;
392       nbIntv = nbUIntv;
393     }   
394         
395     S.D0(U1,V1,P);
396     MoveTo(P,Pts);
397
398     for (Intrv = 1; Intrv <= nbIntv; Intrv++) {
399
400       if (TI(Intrv) <= T1 && TI(Intrv + 1) <= T1)
401         continue;
402       if (TI(Intrv) >= T2 && TI(Intrv + 1) >= T2)
403         continue;
404       if (T == GeomAbs_IsoU) {
405         V1 = Max(T1, TI(Intrv));
406         V2 = Min(T2, TI(Intrv + 1));
407         stepV = (V2 - V1) / myDiscret;
408       }
409       else {
410         U1 = Max(T1, TI(Intrv));
411         U2 = Min(T2, TI(Intrv + 1));
412         stepU = (U2 - U1) / myDiscret;
413       }
414
415       switch (SurfType) {
416         //-------------GeomAbs_Plane---------------
417       case GeomAbs_Plane :
418         break;
419         //----GeomAbs_Cylinder   GeomAbs_Cone------
420       case GeomAbs_Cylinder :
421       case GeomAbs_Cone :
422         if (T == GeomAbs_IsoV) {
423           for (j = 1; j < myDiscret; j++) {
424             U1 += stepU;
425             V1 += stepV;
426             S.D0(U1,V1,P);
427             DrawTo(P,Pts,Cells);
428           }
429         }
430         break;
431         //---GeomAbs_Sphere   GeomAbs_Torus--------
432         //GeomAbs_BezierSurface GeomAbs_BezierSurface
433       case GeomAbs_Sphere :
434       case GeomAbs_Torus :
435       case GeomAbs_OffsetSurface :
436       case GeomAbs_OtherSurface :
437         for (j = 1; j < myDiscret; j++) {
438           U1 += stepU;
439           V1 += stepV;
440           S.D0(U1,V1,P);
441           DrawTo(P,Pts,Cells);
442         }
443         break;
444         //-------------GeomAbs_BSplineSurface------
445       case GeomAbs_BezierSurface :
446       case GeomAbs_BSplineSurface :
447         for (j = 1; j <= myDiscret/2; j++) {
448
449           PlotCount = 0;
450
451           PlotIso ( S, T, U1, V1, (T == GeomAbs_IsoV) ? stepU*2. : stepV*2., halt, Pts, Cells);
452           U1 += stepU*2.;
453           V1 += stepV*2.;
454         }
455         break;
456         //-------------GeomAbs_SurfaceOfExtrusion--
457         //-------------GeomAbs_SurfaceOfRevolution-
458       case GeomAbs_SurfaceOfExtrusion :
459       case GeomAbs_SurfaceOfRevolution :
460         if ((T == GeomAbs_IsoV && SurfType == GeomAbs_SurfaceOfRevolution) ||
461             (T == GeomAbs_IsoU && SurfType == GeomAbs_SurfaceOfExtrusion)) {
462           if (SurfType == GeomAbs_SurfaceOfExtrusion) break;
463           for (j = 1; j < myDiscret; j++) {
464             U1 += stepU;
465             V1 += stepV;
466             S.D0(U1,V1,P);
467             DrawTo(P,Pts,Cells);
468           }
469         } else {
470           CurvType = (S.BasisCurve())->GetType();
471           switch (CurvType) {
472           case GeomAbs_Line :
473             break;
474           case GeomAbs_Circle :
475           case GeomAbs_Ellipse :
476             for (j = 1; j < myDiscret; j++) {
477               U1 += stepU;
478               V1 += stepV;
479               S.D0(U1,V1,P);
480               DrawTo(P,Pts,Cells);
481             }
482             break;
483           case GeomAbs_Parabola :
484           case GeomAbs_Hyperbola :
485           case GeomAbs_BezierCurve :
486           case GeomAbs_BSplineCurve :
487           case GeomAbs_OffsetCurve:
488           case GeomAbs_OtherCurve :
489             for (j = 1; j <= myDiscret/2; j++) {
490
491               PlotCount = 0;
492
493               PlotIso ( S, T, U1, V1,(T == GeomAbs_IsoV) ? stepU*2. : stepV*2., halt, Pts, Cells);
494               U1 += stepU*2.;
495               V1 += stepV*2.;
496             }
497             break;
498           }
499         }
500       }
501     }
502     S.D0(U2,V2,P);
503     DrawTo(P,Pts,Cells);
504   }  
505 }
506
507 //=======================================================================
508 // Function : PlotIso
509 // Purpose  : Plot iso for other surface
510 //=======================================================================
511
512 void GEOM_OCCReader::PlotIso (BRepAdaptor_Surface& S, 
513                                 GeomAbs_IsoType T,
514                                 Standard_Real& U, 
515                                 Standard_Real& V, 
516                                 Standard_Real Step, 
517                                 Standard_Boolean& halt,
518                                 vtkPoints* Pts,
519                                 vtkCellArray* Cells)
520 {
521
522   ++PlotCount; 
523
524   gp_Pnt Pl, Pr, Pm;
525
526   if (T == GeomAbs_IsoU) {
527     S.D0(U, V, Pl);
528     S.D0(U, V + Step/2., Pm);
529     S.D0(U, V + Step, Pr);
530   } else {
531     S.D0(U, V, Pl);
532     S.D0(U + Step/2., V, Pm);
533     S.D0(U + Step, V, Pr);
534   }
535
536   if (PlotCount > MaxPlotCount) {
537     DrawTo(Pr,Pts,Cells);
538     return;
539   }
540
541   if (Pm.Distance(Pl) + Pm.Distance(Pr) <= IsoRatio*Pl.Distance(Pr)) {
542     DrawTo(Pr,Pts,Cells);
543   } else 
544     if (T == GeomAbs_IsoU) {
545       PlotIso ( S, T, U, V, Step/2, halt, Pts, Cells);
546       Standard_Real aLocalV = V + Step/2 ;
547       PlotIso ( S, T, U, aLocalV , Step/2, halt, Pts, Cells);
548     } else {
549       PlotIso ( S, T, U, V, Step/2, halt, Pts, Cells);
550       Standard_Real aLocalU = U + Step/2 ;
551       PlotIso ( S, T, aLocalU , V, Step/2, halt, Pts, Cells);
552     }
553 }
554
555 //=======================================================================
556 // Function : TransferEdgeWData
557 // Purpose  : Transfer wireframe data for EDGE
558 //=======================================================================
559
560 void GEOM_OCCReader::TransferEdgeWData(const TopoDS_Edge& aEdge,
561                                          vtkPoints* Pts,
562                                          vtkCellArray* Cells)
563 {
564   Handle(Poly_PolygonOnTriangulation) aEdgePoly;
565   Standard_Integer i = 1;
566   Handle(Poly_Triangulation) T;
567   TopLoc_Location aEdgeLoc;
568   BRep_Tool::PolygonOnTriangulation(aEdge, aEdgePoly, T, aEdgeLoc, i);
569   
570   Handle(Poly_Polygon3D) P;
571   if(aEdgePoly.IsNull()) {
572     P = BRep_Tool::Polygon3D(aEdge, aEdgeLoc);
573   }
574
575   if(P.IsNull() && aEdgePoly.IsNull())
576     return;
577   
578   // Location edges
579   //---------------
580   
581   gp_Trsf edgeTransf;
582   Standard_Boolean isidtrsf = true;
583   if(!aEdgeLoc.IsIdentity())  {
584     isidtrsf = false;
585     edgeTransf = aEdgeLoc.Transformation();
586   }
587
588   gp_Pnt aP1, aP2;
589
590   Standard_Integer nbnodes;
591   if (aEdgePoly.IsNull()) {
592     nbnodes = P->NbNodes();
593     const TColgp_Array1OfPnt& theNodesP = P->Nodes();
594
595     aP1 = theNodesP(1);
596     aP2 = theNodesP(nbnodes);
597
598     float coord[3];
599     vtkIdType pts[2];
600
601     for(int j=1;j<nbnodes;j++) {
602       gp_Pnt pt1 = theNodesP(j);
603       gp_Pnt pt2 = theNodesP(j+1);
604     
605       if(!isidtrsf) {
606         // apply edge transformation
607         pt1.Transform(edgeTransf);
608         pt2.Transform(edgeTransf);
609       }
610       
611       // insert pt1
612       coord[0] = pt1.X(); coord[1] = pt1.Y(); coord[2] = pt1.Z();
613       pts[0] = Pts->InsertNextPoint(coord);
614       
615       // insert pt2
616       coord[0] = pt2.X(); coord[1] = pt2.Y(); coord[2] = pt2.Z();
617       pts[1] = Pts->InsertNextPoint(coord);
618       
619       // insert line (pt1,pt2)
620       Cells->InsertNextCell(2,pts);
621     }
622   } else {
623     nbnodes = aEdgePoly->NbNodes();
624     const TColStd_Array1OfInteger& Nodesidx = aEdgePoly->Nodes();
625     const TColgp_Array1OfPnt& theNodesPoly = T->Nodes();
626
627     aP1 = theNodesPoly(1);
628     aP2 = theNodesPoly(nbnodes);
629
630     float coord[3];
631     vtkIdType pts[2];
632     
633     for(int j=1;j<nbnodes;j++) {
634       Standard_Integer id1 = Nodesidx(j);
635       Standard_Integer id2 = Nodesidx(j+1);
636       
637       gp_Pnt pt1 = theNodesPoly(id1);
638       gp_Pnt pt2 = theNodesPoly(id2);
639           
640       if(!isidtrsf) {
641         // apply edge transformation
642         pt1.Transform(edgeTransf);
643         pt2.Transform(edgeTransf);
644       }
645       
646       // insert pt1
647       coord[0] = pt1.X(); coord[1] = pt1.Y(); coord[2] = pt1.Z();
648       pts[0] = Pts->InsertNextPoint(coord);
649       
650       // insert pt2
651       coord[0] = pt2.X(); coord[1] = pt2.Y(); coord[2] = pt2.Z();
652       pts[1] = Pts->InsertNextPoint(coord);
653       
654       // insert line (pt1,pt2)
655       Cells->InsertNextCell(2,pts);
656     }
657   }
658
659   // vector representation has an arrow on its end
660   if (myIsVector)
661   {
662     if (!isidtrsf) {
663       // apply edge transformation
664       aP1.Transform(edgeTransf);
665       aP2.Transform(edgeTransf);
666     }
667
668     // draw an arrow
669     gp_Vec aDirVec (aP1, aP2);
670     Standard_Real aDist = aDirVec.Magnitude();
671     if (aDist < gp::Resolution()) return;
672     gp_Dir aDirection (aDirVec);
673
674     Standard_Real anAngle = M_PI/180. * 5.;
675     Standard_Real aLength = aDist/10.;
676
677     Standard_Real dx,dy,dz;
678     aDirection.Coord(dx,dy,dz);
679
680     // Pointe de la fleche
681     Standard_Real xo,yo,zo;
682     aP2.Coord(xo,yo,zo);
683
684     // Centre du cercle base de la fleche
685     gp_XYZ aPc = aP2.XYZ() - aDirection.XYZ() * aLength;
686
687     // Construction d'un repere i,j pour le cercle
688     gp_Dir aDirN;
689     if      (Abs(dx) <= Abs(dy) && Abs(dx) <= Abs(dz)) aDirN = gp::DX();
690     else if (Abs(dy) <= Abs(dz) && Abs(dy) <= Abs(dx)) aDirN = gp::DY();
691     else aDirN = gp::DZ();
692
693     gp_Dir aDirI = aDirection ^ aDirN;
694     gp_Dir aDirJ = aDirection ^ aDirI;
695
696     // Add points and segments, composing the arrow
697     Standard_Real cosinus, sinus, Tg = tan(anAngle);
698
699     float coord[3];
700     coord[0] = xo; coord[1] = yo; coord[2] = zo;
701
702     int ptLoc = Pts->InsertNextPoint(coord);
703     int ptFirst = 0;
704     int ptPrev = 0;
705     int ptCur = 0;
706
707     vtkIdType pts[2];
708
709     int NbPoints = 15;
710     for (int i = 1; i <= NbPoints; i++, ptPrev = ptCur)
711     {
712       cosinus = cos(2. * M_PI / NbPoints * (i-1));   
713       sinus   = sin(2. * M_PI / NbPoints * (i-1));
714
715       gp_XYZ aP = aPc + (aDirI.XYZ() * cosinus + aDirJ.XYZ() * sinus) * aLength * Tg;
716       coord[0] = aP.X();
717       coord[1] = aP.Y();
718       coord[2] = aP.Z();
719
720       // insert pts
721       ptCur = Pts->InsertNextPoint(coord);
722       pts[0] = ptCur;
723
724       if (i == 1) {
725         ptFirst = ptCur;
726       }
727       else {
728         // insert line (ptCur,ptPrev)
729         pts[1] = ptPrev;
730         Cells->InsertNextCell(2,pts);
731       }
732
733       // insert line (ptCur,ptLoc)
734       pts[1] = ptLoc;
735       Cells->InsertNextCell(2,pts);
736     }
737
738     // insert line (ptCur,ptFirst)
739     pts[0] = ptCur;
740     pts[1] = ptFirst;
741     Cells->InsertNextCell(2,pts);
742   }
743 }
744
745 /*  Standard_Integer nbnodes = aEdgePoly->NbNodes();
746   const TColStd_Array1OfInteger& Nodesidx = aEdgePoly->Nodes();
747   const TColgp_Array1OfPnt& theNodes = T->Nodes();
748     
749   float coord[3];
750   vtkIdType pts[2];
751     
752
753   // PUSH NODES
754   for(i=1;i<=nbnodes;i++) {
755     Standard_Integer id = Nodesidx(i);
756     gp_Pnt pt = theNodes(id);
757  
758     float coord[3];
759     if(!isidtrsf) pt.Transform(edgeTransf);
760
761     coord[0] = pt.X(); coord[1] = pt.Y(); coord[2] = pt.Z();
762
763     Pts->SetPoint(id-1,coord);
764
765   }
766
767   // PUSH EDGES
768   for(i=1;i<nbnodes;i++) {
769       
770     Standard_Integer id1 = Nodesidx(i);
771     Standard_Integer id2 = Nodesidx(i+1);
772     
773     vtkIdType pts[2];
774     pts[0] = id1-1; pts[1] = id2-1;
775
776     // insert line (pt1,pt2)
777     Cells->InsertNextCell(2,pts);
778   }
779  
780   
781   }*/
782
783 //=======================================================================
784 // Function : TransferVertexWData
785 // Purpose  : Transfer wireframe data for VERTEX
786 //=======================================================================
787
788 void GEOM_OCCReader::TransferVertexWData(const TopoDS_Vertex& /*aVertex*/,
789                                          vtkPoints*           Pts,
790                                          vtkCellArray*        Cells)
791 {
792 #define ZERO_COORD coord[0] = coord[1] = coord[2] = 0.0
793   
794   // gp_Pnt P = BRep_Tool::Pnt( aVertex ); ??????????????????????????
795   float delta = 1, coord[3];
796   vtkIdType pts[2];
797   // insert pt
798   ZERO_COORD; coord[0] = +delta;
799   pts[0] = Pts->InsertNextPoint(coord);
800   coord[0] = -delta;
801   pts[1] = Pts->InsertNextPoint(coord);
802   // insert line (pt1,pt2)
803   Cells->InsertNextCell(2,pts);
804
805   ZERO_COORD; coord[1] = +delta;
806   pts[0] = Pts->InsertNextPoint(coord);
807   coord[1] = -delta;
808   pts[1] = Pts->InsertNextPoint(coord);
809   // insert line (pt1,pt2)
810   Cells->InsertNextCell(2,pts);
811
812   ZERO_COORD; coord[2] = +delta;
813   pts[0] = Pts->InsertNextPoint(coord);
814   coord[2] = -delta;
815   pts[1] = Pts->InsertNextPoint(coord);
816   // insert line (pt1,pt2)
817   Cells->InsertNextCell(2,pts);
818
819 #undef ZERO_COORD
820 }       
821
822 //=======================================================================
823 // Function : TransferEdgeSData(
824 // Purpose  : Transfer shading data for EDGE
825 //=======================================================================
826
827 void GEOM_OCCReader::TransferEdgeSData(const TopoDS_Edge& /*aFace*/,
828                                        vtkPoints* /*Pts*/,
829                                        vtkCellArray* /*Cells*/) 
830 {
831 }
832
833
834 //=======================================================================
835 // Function : TransferFaceSData
836 // Purpose  : Transfer shading data for FACE
837 //=======================================================================
838 void GEOM_OCCReader::TransferFaceSData(const TopoDS_Face& aFace,
839                                          vtkPoints* Pts,
840                                          vtkCellArray* Cells) {
841
842   TopLoc_Location aLoc;
843   Handle(Poly_Triangulation) aPoly = BRep_Tool::Triangulation(aFace,aLoc);
844   if(aPoly.IsNull()) return;
845   else {
846     
847     gp_Trsf myTransf;
848     Standard_Boolean identity = true;
849     if(!aLoc.IsIdentity())  {
850       identity = false;
851       myTransf = aLoc.Transformation();
852     }
853
854     Standard_Integer nbNodesInFace = aPoly->NbNodes();
855     Standard_Integer nbTriInFace = aPoly->NbTriangles();
856                 
857     const Poly_Array1OfTriangle& Triangles = aPoly->Triangles();
858     const TColgp_Array1OfPnt& Nodes = aPoly->Nodes();
859       
860     Standard_Integer i;
861     for(i=1;i<=nbNodesInFace;i++) {
862       gp_Pnt P = Nodes(i);
863       float coord[3];
864       if(!identity) P.Transform(myTransf);
865       coord[0] = P.X(); coord[1] = P.Y(); coord[2] = P.Z();
866       Pts->SetPoint(i-1,coord);
867     }
868
869     for(i=1;i<=nbTriInFace;i++) {
870       // Get the triangle
871         
872       Standard_Integer N1,N2,N3;
873       Triangles(i).Get(N1,N2,N3);
874         
875       vtkIdType pts[3];
876       pts[0] = N1-1; pts[1] = N2-1; pts[2] = N3-1;
877       Cells->InsertNextCell(3,pts);
878
879     }
880   } 
881 }
882
883 //=======================================================================
884 // Function : ComputeShading
885 // Purpose  : Compute the shape in shading mode
886 //=======================================================================
887 void GEOM_OCCReader::ComputeShading(vtkPoints* Pts,vtkCellArray* Cells){
888
889   // Check the type of the shape:
890   if(myShape.ShapeType() == TopAbs_FACE) {
891     // Face
892     TransferFaceSData(TopoDS::Face(myShape),Pts,Cells);
893   }
894   else {
895     if(myShape.ShapeType() == TopAbs_EDGE) {
896       // Edge
897       TransferEdgeSData(TopoDS::Edge(myShape),Pts,Cells);
898     }
899     else {
900     }
901
902   } 
903
904 }
905
906 //=======================================================================
907 // Function : 
908 // Purpose  : Set parameters
909 //=======================================================================
910 void GEOM_OCCReader::setDisplayMode(int thenewmode) {
911   amode = thenewmode;
912 }
913
914 void GEOM_OCCReader::setTopo(const TopoDS_Shape& aShape, bool isVector) {
915   myShape = aShape;
916   myIsVector = isVector;
917 }
918
919 void GEOM_OCCReader::setForceUpdate(Standard_Boolean bol) {
920   forced = bol;
921 }
922
923 //=======================================================================
924 // Function : 
925 // Purpose  : Get parameters
926 //=======================================================================
927 const TopoDS_Shape& GEOM_OCCReader::getTopo() {
928   return myShape;
929 }
930
931 int GEOM_OCCReader::getDisplayMode() {
932   return amode;
933 }