]> SALOME platform Git repositories - modules/smesh.git/blob - src/OBJECT/SMESH_DeviceActor.cxx
Salome HOME
0022100: EDF 2413 SMESH: Take into account TRIA7
[modules/smesh.git] / src / OBJECT / SMESH_DeviceActor.cxx
1 // Copyright (C) 2007-2013  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.
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 //  SMESH OBJECT : interactive object for SMESH visualization
24 //  File   : SMESH_DeviceActor.cxx
25 //  Author : 
26 //  Module : SMESH
27 //
28 #include "SMESH_DeviceActor.h"
29 #include "SMESH_ScalarBarActor.h"
30 #include "SMESH_ExtractGeometry.h"
31 #include "SMESH_ControlsDef.hxx"
32 #include "SMESH_ActorUtils.h"
33 #include "SMESH_FaceOrientationFilter.h"
34 #include "VTKViewer_CellLocationsArray.h"
35 #include "VTKViewer_PolyDataMapper.h"
36
37 #include <VTKViewer_Transform.h>
38 #include <VTKViewer_TransformFilter.h>
39 #include <VTKViewer_ExtractUnstructuredGrid.h>
40
41 // VTK Includes
42 #include <vtkObjectFactory.h>
43 #include <vtkShrinkFilter.h>
44 #include <vtkShrinkPolyData.h>
45
46 #include <vtkProperty.h>
47 #include <vtkPolyData.h>
48 #include <vtkMergeFilter.h>
49 #include <vtkPolyDataMapper.h>
50 #include <vtkUnstructuredGrid.h>
51
52 #include <vtkLookupTable.h>
53 #include <vtkDoubleArray.h>
54 #include <vtkCellData.h>
55
56 #include <vtkCell.h>
57 #include <vtkIdList.h>
58 #include <vtkCellArray.h>
59 #include <vtkUnsignedCharArray.h>
60
61 #include <vtkImplicitBoolean.h>
62 #include <vtkPassThroughFilter.h>
63
64 #include <vtkRenderer.h>
65
66 #include "utilities.h"
67
68 #ifdef _DEBUG_
69 static int MYDEBUG = 0;
70 #else
71 static int MYDEBUG = 0;
72 #endif
73
74 using namespace std;
75
76
77 vtkStandardNewMacro(SMESH_DeviceActor);
78
79
80 SMESH_DeviceActor
81 ::SMESH_DeviceActor()
82 {
83   if(MYDEBUG) MESSAGE("SMESH_DeviceActor - "<<this);
84
85   myIsShrinkable = false;
86   myIsShrunk = false;
87   myIsHighlited = false;
88
89   myRepresentation = eSurface;
90
91   myProperty = vtkProperty::New();
92   myMapper = VTKViewer_PolyDataMapper::New();
93
94   vtkMapper::GetResolveCoincidentTopologyPolygonOffsetParameters(myPolygonOffsetFactor,
95                                                                  myPolygonOffsetUnits);
96
97   myMapper->UseLookupTableScalarRangeOn();
98   myMapper->SetColorModeToMapScalars();
99
100   myShrinkFilter = vtkShrinkFilter::New();
101
102   myStoreClippingMapping = false;
103
104   myExtractGeometry = SMESH_ExtractGeometry::New();
105   myExtractGeometry->SetReleaseDataFlag(true);
106   myIsImplicitFunctionUsed = false;
107
108   myExtractUnstructuredGrid = VTKViewer_ExtractUnstructuredGrid::New();
109     
110   myMergeFilter = vtkMergeFilter::New();
111
112   myGeomFilter = VTKViewer_GeometryFilter::New();
113
114   myTransformFilter = VTKViewer_TransformFilter::New();
115
116   for(int i = 0; i < 6; i++)
117     myPassFilter.push_back(vtkPassThroughFilter::New());
118
119   // Orientation of faces
120   myIsFacesOriented = false;
121
122   double anRGB[3] = { 1, 1, 1 };
123   SMESH::GetColor( "SMESH", "orientation_color", anRGB[0], anRGB[1], anRGB[2], QColor( 255, 255, 255 ) );
124
125   myFaceOrientationFilter = SMESH_FaceOrientationFilter::New();
126
127   myFaceOrientationDataMapper = vtkPolyDataMapper::New();
128   myFaceOrientationDataMapper->SetInputConnection(myFaceOrientationFilter->GetOutputPort());
129
130   myFaceOrientation = vtkActor::New();
131   myFaceOrientation->SetMapper(myFaceOrientationDataMapper);
132   myFaceOrientation->GetProperty()->SetColor(anRGB[0], anRGB[1], anRGB[2]);
133 }
134
135
136 SMESH_DeviceActor
137 ::~SMESH_DeviceActor()
138 {
139   if(MYDEBUG) MESSAGE("~SMESH_DeviceActor - "<<this);
140
141   myMapper->Delete();
142
143   myProperty->Delete();
144
145   myExtractGeometry->Delete();
146
147   myMergeFilter->Delete();
148   myExtractUnstructuredGrid->Delete();
149
150   // Orientation of faces
151   myFaceOrientationFilter->Delete();
152   myFaceOrientationDataMapper->RemoveAllInputs();
153   myFaceOrientationDataMapper->Delete();
154   myFaceOrientation->Delete();
155
156   myGeomFilter->Delete();
157
158   myTransformFilter->Delete();
159
160   for(int i = 0, iEnd = myPassFilter.size(); i < iEnd; i++)
161     myPassFilter[i]->Delete();
162
163   myShrinkFilter->Delete();
164 }
165
166
167 void
168 SMESH_DeviceActor
169 ::SetStoreGemetryMapping(bool theStoreMapping)
170 {
171   myGeomFilter->SetStoreMapping(theStoreMapping);
172   SetStoreClippingMapping(theStoreMapping);
173 }
174
175
176 void
177 SMESH_DeviceActor
178 ::SetStoreClippingMapping(bool theStoreMapping)
179 {
180   myStoreClippingMapping = theStoreMapping;
181   myExtractGeometry->SetStoreMapping(theStoreMapping && myIsImplicitFunctionUsed);
182   SetStoreIDMapping(theStoreMapping);
183 }
184
185
186 void
187 SMESH_DeviceActor
188 ::SetStoreIDMapping(bool theStoreMapping)
189 {
190   myExtractUnstructuredGrid->SetStoreMapping(theStoreMapping);
191 }
192
193
194 void 
195 SMESH_DeviceActor
196 ::Init(TVisualObjPtr theVisualObj, 
197        vtkImplicitBoolean* theImplicitBoolean)
198 {
199   myVisualObj = theVisualObj;
200   myExtractGeometry->SetImplicitFunction(theImplicitBoolean);
201   SetUnstructuredGrid(myVisualObj->GetUnstructuredGrid());
202 }
203
204
205 void
206 SMESH_DeviceActor
207 ::SetImplicitFunctionUsed(bool theIsImplicitFunctionUsed)
208 {
209   int anId = 0;
210   if(theIsImplicitFunctionUsed)
211     myPassFilter[ anId ]->SetInputConnection( myExtractGeometry->GetOutputPort() );
212   else
213     myPassFilter[ anId ]->SetInputConnection( myMergeFilter->GetOutputPort() );
214     
215   myIsImplicitFunctionUsed = theIsImplicitFunctionUsed;
216   SetStoreClippingMapping(myStoreClippingMapping);
217 }
218
219
220 void
221 SMESH_DeviceActor
222 ::SetUnstructuredGrid(vtkUnstructuredGrid* theGrid)
223 {
224   if(theGrid){
225     //myIsShrinkable = theGrid->GetNumberOfCells() > 10;
226     myIsShrinkable = true;
227
228     myExtractUnstructuredGrid->SetInputData(theGrid);
229
230     myMergeFilter->SetGeometryConnection(myExtractUnstructuredGrid->GetOutputPort());
231
232     myExtractGeometry->SetInputConnection(myMergeFilter->GetOutputPort());
233
234     int anId = 0;
235     SetImplicitFunctionUsed(myIsImplicitFunctionUsed);
236     myPassFilter[ anId + 1]->SetInputConnection( myPassFilter[ anId ]->GetOutputPort() );
237     
238     anId++; // 1
239     myTransformFilter->SetInputConnection( myPassFilter[ anId ]->GetOutputPort() );
240
241     anId++; // 2
242     myPassFilter[ anId ]->SetInputConnection( myTransformFilter->GetOutputPort() );
243     myPassFilter[ anId + 1 ]->SetInputConnection( myPassFilter[ anId ]->GetOutputPort() );
244
245     anId++; // 3
246     myGeomFilter->SetInputConnection( myPassFilter[ anId ]->GetOutputPort() );
247
248     anId++; // 4
249     myPassFilter[ anId ]->SetInputConnection( myGeomFilter->GetOutputPort() ); 
250     myPassFilter[ anId + 1 ]->SetInputConnection( myPassFilter[ anId ]->GetOutputPort() );
251
252     anId++; // 5
253     myMapper->SetInputConnection( myPassFilter[ anId ]->GetOutputPort() );
254
255     vtkLODActor::SetMapper( myMapper );
256     Modified();
257   }
258 }
259
260
261 VTKViewer_ExtractUnstructuredGrid* 
262 SMESH_DeviceActor
263 ::GetExtractUnstructuredGrid()
264 {
265   return myExtractUnstructuredGrid;
266 }
267
268
269 vtkUnstructuredGrid* 
270 SMESH_DeviceActor
271 ::GetUnstructuredGrid()
272 {
273   myExtractUnstructuredGrid->Update();
274   return myExtractUnstructuredGrid->GetOutput();
275 }
276
277
278 void
279 SMESH_DeviceActor
280 ::SetControlMode(SMESH::Controls::FunctorPtr theFunctor,
281                  SMESH_ScalarBarActor* theScalarBarActor,
282                  vtkLookupTable* theLookupTable)
283 {
284   bool anIsInitialized = theFunctor;
285   if(anIsInitialized){
286     vtkUnstructuredGrid* aDataSet = vtkUnstructuredGrid::New();
287
288     SetStoreIDMapping(true);
289     myExtractUnstructuredGrid->Update();
290     vtkUnstructuredGrid* aGrid = myExtractUnstructuredGrid->GetOutput();
291
292     aDataSet->ShallowCopy(aGrid);
293     
294     vtkDoubleArray *aScalars = vtkDoubleArray::New();
295     vtkIdType aNbCells = aGrid->GetNumberOfCells();
296     aScalars->SetNumberOfComponents(1);
297     aScalars->SetNumberOfTuples(aNbCells);
298     double* range;// = aScalars->GetRange();
299     
300     myVisualObj->UpdateFunctor(theFunctor);
301
302     using namespace SMESH::Controls;
303     if(NumericalFunctor* aNumericalFunctor = dynamic_cast<NumericalFunctor*>(theFunctor.get()))
304     {
305       for(vtkIdType i = 0; i < aNbCells; i++){
306         vtkIdType anId = myExtractUnstructuredGrid->GetInputId(i);
307         vtkIdType anObjId = myVisualObj->GetElemObjId(anId);
308         double aValue = aNumericalFunctor->GetValue(anObjId);
309         aScalars->SetValue(i,aValue);
310       }
311       range = aScalars->GetRange();
312       if ( range[1] - range[0] < ( qMax(qAbs(range[0]),qAbs(range[1])) + 1e-100 ) * 1e-6 )
313       {
314         range[1] = range[0];
315         for(vtkIdType i = 0; i < aNbCells; i++)
316           aScalars->SetValue(i,range[0]);
317       }
318     }
319     else if(Predicate* aPredicate = dynamic_cast<Predicate*>(theFunctor.get()))
320     {
321       for(vtkIdType i = 0; i < aNbCells; i++){
322         vtkIdType anId = myExtractUnstructuredGrid->GetInputId(i);
323         vtkIdType anObjId = myVisualObj->GetElemObjId(anId);
324         bool aValue = aPredicate->IsSatisfy(anObjId);
325         aScalars->SetValue(i,aValue);
326       }
327       range = aScalars->GetRange();
328     }
329
330     aDataSet->GetCellData()->SetScalars(aScalars);
331     aScalars->Delete();
332
333     theLookupTable->SetRange( range );
334     theLookupTable->SetNumberOfTableValues(theScalarBarActor->GetMaximumNumberOfColors());
335     theLookupTable->Build();
336     
337     myMergeFilter->SetScalarsData(aDataSet);
338     aDataSet->Delete();
339   }
340   GetMapper()->SetScalarVisibility(anIsInitialized);
341   theScalarBarActor->SetVisibility(anIsInitialized);
342 }
343
344 void
345 SMESH_DeviceActor
346 ::SetExtControlMode(SMESH::Controls::FunctorPtr theFunctor,
347                     SMESH_ScalarBarActor* theScalarBarActor,
348                     vtkLookupTable* theLookupTable)
349 {
350   bool anIsInitialized = theFunctor;
351   myExtractUnstructuredGrid->ClearRegisteredCells();
352   myExtractUnstructuredGrid->ClearRegisteredCellsWithType();
353   myExtractUnstructuredGrid->SetModeOfChanging(VTKViewer_ExtractUnstructuredGrid::ePassAll);
354   myVisualObj->UpdateFunctor(theFunctor);
355
356   using namespace SMESH::Controls;
357   if (anIsInitialized){
358     if (Length2D* aLength2D = dynamic_cast<Length2D*>(theFunctor.get())){
359       SMESH::Controls::Length2D::TValues aValues;
360
361       aLength2D->GetValues(aValues);
362       vtkUnstructuredGrid* aDataSet = vtkUnstructuredGrid::New();
363       vtkUnstructuredGrid* aGrid = myVisualObj->GetUnstructuredGrid();
364
365       aDataSet->SetPoints(aGrid->GetPoints());
366       
367       vtkIdType aNbCells = aValues.size();
368       
369       vtkDoubleArray *aScalars = vtkDoubleArray::New();
370       aScalars->SetNumberOfComponents(1);
371       aScalars->SetNumberOfTuples(aNbCells);
372
373       vtkIdType aCellsSize = 3*aNbCells;
374       vtkCellArray* aConnectivity = vtkCellArray::New();
375       aConnectivity->Allocate( aCellsSize, 0 );
376       
377       vtkUnsignedCharArray* aCellTypesArray = vtkUnsignedCharArray::New();
378       aCellTypesArray->SetNumberOfComponents( 1 );
379       aCellTypesArray->Allocate( aNbCells * aCellTypesArray->GetNumberOfComponents() );
380       
381       vtkIdList *anIdList = vtkIdList::New();
382       anIdList->SetNumberOfIds(2);
383       
384       Length2D::TValues::const_iterator anIter = aValues.begin();
385       for(vtkIdType aVtkId = 0; anIter != aValues.end(); anIter++,aVtkId++){
386         const Length2D::Value& aValue = *anIter;
387         int aNode[2] = {
388           myVisualObj->GetNodeVTKId(aValue.myPntId[0]),
389           myVisualObj->GetNodeVTKId(aValue.myPntId[1])
390         };
391         if(aNode[0] >= 0 && aNode[1] >= 0){
392           anIdList->SetId( 0, aNode[0] );
393           anIdList->SetId( 1, aNode[1] );
394           aConnectivity->InsertNextCell( anIdList );
395           aCellTypesArray->InsertNextValue( VTK_LINE );
396           aScalars->SetValue(aVtkId,aValue.myLength);
397         }
398       }
399       
400       VTKViewer_CellLocationsArray* aCellLocationsArray = VTKViewer_CellLocationsArray::New();
401       aCellLocationsArray->SetNumberOfComponents( 1 );
402       aCellLocationsArray->SetNumberOfTuples( aNbCells );
403       
404       aConnectivity->InitTraversal();
405       for( vtkIdType idType = 0, *pts, npts; aConnectivity->GetNextCell( npts, pts ); idType++ )
406         aCellLocationsArray->SetValue( idType, aConnectivity->GetTraversalLocation( npts ) );
407       
408       aDataSet->SetCells( aCellTypesArray, aCellLocationsArray,aConnectivity );
409       SetUnstructuredGrid(aDataSet);
410
411       aDataSet->GetCellData()->SetScalars(aScalars);
412       aScalars->Delete();
413       
414       theLookupTable->SetRange(aScalars->GetRange());
415       theLookupTable->Build();
416       
417       myMergeFilter->SetScalarsData(aDataSet);
418       aDataSet->Delete();
419     }
420     else if (MultiConnection2D* aMultiConnection2D = dynamic_cast<MultiConnection2D*>(theFunctor.get())){
421       SMESH::Controls::MultiConnection2D::MValues aValues;
422
423       aMultiConnection2D->GetValues(aValues);
424       vtkUnstructuredGrid* aDataSet = vtkUnstructuredGrid::New();
425       vtkUnstructuredGrid* aGrid = myVisualObj->GetUnstructuredGrid();
426       aDataSet->SetPoints(aGrid->GetPoints());
427       
428       vtkIdType aNbCells = aValues.size();
429       vtkDoubleArray *aScalars = vtkDoubleArray::New();
430       aScalars->SetNumberOfComponents(1);
431       aScalars->SetNumberOfTuples(aNbCells);
432
433       vtkIdType aCellsSize = 3*aNbCells;
434       vtkCellArray* aConnectivity = vtkCellArray::New();
435       aConnectivity->Allocate( aCellsSize, 0 );
436       
437       vtkUnsignedCharArray* aCellTypesArray = vtkUnsignedCharArray::New();
438       aCellTypesArray->SetNumberOfComponents( 1 );
439       aCellTypesArray->Allocate( aNbCells * aCellTypesArray->GetNumberOfComponents() );
440       
441       vtkIdList *anIdList = vtkIdList::New();
442       anIdList->SetNumberOfIds(2);
443       
444       MultiConnection2D::MValues::const_iterator anIter = aValues.begin();
445       for(vtkIdType aVtkId = 0; anIter != aValues.end(); anIter++,aVtkId++){
446         const MultiConnection2D::Value& aValue = (*anIter).first;
447         int aNode[2] = {
448           myVisualObj->GetNodeVTKId(aValue.myPntId[0]),
449           myVisualObj->GetNodeVTKId(aValue.myPntId[1])
450         };
451         if(aNode[0] >= 0 && aNode[1] >= 0){
452           anIdList->SetId( 0, aNode[0] );
453           anIdList->SetId( 1, aNode[1] );
454           aConnectivity->InsertNextCell( anIdList );
455           aCellTypesArray->InsertNextValue( VTK_LINE );
456           aScalars->SetValue(aVtkId,(*anIter).second);
457         }
458       }
459       
460       VTKViewer_CellLocationsArray* aCellLocationsArray = VTKViewer_CellLocationsArray::New();
461       aCellLocationsArray->SetNumberOfComponents( 1 );
462       aCellLocationsArray->SetNumberOfTuples( aNbCells );
463       
464       aConnectivity->InitTraversal();
465       for( vtkIdType idType = 0, *pts, npts; aConnectivity->GetNextCell( npts, pts ); idType++ )
466         aCellLocationsArray->SetValue( idType, aConnectivity->GetTraversalLocation( npts ) );
467       
468       aDataSet->SetCells( aCellTypesArray, aCellLocationsArray,aConnectivity );
469       SetUnstructuredGrid(aDataSet);
470
471       aDataSet->GetCellData()->SetScalars(aScalars);
472       aScalars->Delete();
473       
474       theLookupTable->SetRange(aScalars->GetRange());
475       theLookupTable->Build();
476       
477       myMergeFilter->SetScalarsData(aDataSet);
478       aDataSet->Delete();
479     }
480   }
481   GetMapper()->SetScalarVisibility(anIsInitialized);
482   theScalarBarActor->SetVisibility(anIsInitialized);
483 }
484
485 void
486 SMESH_DeviceActor
487 ::SetExtControlMode(SMESH::Controls::FunctorPtr theFunctor)
488 {
489   myExtractUnstructuredGrid->ClearRegisteredCells();
490   myExtractUnstructuredGrid->ClearRegisteredCellsWithType();
491   myExtractUnstructuredGrid->SetModeOfChanging(VTKViewer_ExtractUnstructuredGrid::ePassAll);
492   myVisualObj->UpdateFunctor(theFunctor);
493
494   using namespace SMESH::Controls;
495   if ( dynamic_cast<FreeBorders          *>(theFunctor.get()) ||
496        dynamic_cast<FreeFaces            *>(theFunctor.get()) ||
497        dynamic_cast<BareBorderVolume     *>(theFunctor.get()) ||
498        dynamic_cast<BareBorderFace       *>(theFunctor.get()) ||
499        dynamic_cast<OverConstrainedVolume*>(theFunctor.get()) ||
500        dynamic_cast<CoincidentElements1D *>(theFunctor.get()) ||
501        dynamic_cast<CoincidentElements2D *>(theFunctor.get()) ||
502        dynamic_cast<CoincidentElements3D *>(theFunctor.get()) ||
503        dynamic_cast<OverConstrainedFace  *>(theFunctor.get()))
504   {
505     Predicate* aPredicate = dynamic_cast<Predicate*>(theFunctor.get());
506     myExtractUnstructuredGrid->SetModeOfChanging(VTKViewer_ExtractUnstructuredGrid::eAdding);
507     vtkUnstructuredGrid* aGrid = myVisualObj->GetUnstructuredGrid();
508     vtkIdType aNbCells = aGrid->GetNumberOfCells();
509     for( vtkIdType i = 0; i < aNbCells; i++ ){
510       vtkIdType anObjId = myVisualObj->GetElemObjId(i);
511       if(aPredicate->IsSatisfy(anObjId))
512         myExtractUnstructuredGrid->RegisterCell(i);
513     }
514     if(!myExtractUnstructuredGrid->IsCellsRegistered())
515       myExtractUnstructuredGrid->RegisterCell(-1);
516     SetUnstructuredGrid(myVisualObj->GetUnstructuredGrid());
517   }
518   else if(FreeEdges* aFreeEdges = dynamic_cast<FreeEdges*>(theFunctor.get()))
519   {
520     SMESH::Controls::FreeEdges::TBorders aBorders;
521     aFreeEdges->GetBoreders(aBorders);
522     vtkUnstructuredGrid* aDataSet = vtkUnstructuredGrid::New();
523     vtkUnstructuredGrid* aGrid = myVisualObj->GetUnstructuredGrid();
524     aDataSet->SetPoints(aGrid->GetPoints());
525
526     vtkIdType aNbCells = aBorders.size();
527     vtkIdType aCellsSize = 3*aNbCells;
528     vtkCellArray* aConnectivity = vtkCellArray::New();
529     aConnectivity->Allocate( aCellsSize, 0 );
530     
531     vtkUnsignedCharArray* aCellTypesArray = vtkUnsignedCharArray::New();
532     aCellTypesArray->SetNumberOfComponents( 1 );
533     aCellTypesArray->Allocate( aNbCells * aCellTypesArray->GetNumberOfComponents() );
534     
535     vtkIdList *anIdList = vtkIdList::New();
536     anIdList->SetNumberOfIds(2);
537     
538     FreeEdges::TBorders::const_iterator anIter = aBorders.begin();
539     for(; anIter != aBorders.end(); anIter++){
540       const FreeEdges::Border& aBorder = *anIter;
541       int aNode[2] = {
542         myVisualObj->GetNodeVTKId(aBorder.myPntId[0]),
543         myVisualObj->GetNodeVTKId(aBorder.myPntId[1])
544       };
545       //cout<<"aNode = "<<aBorder.myPntId[0]<<"; "<<aBorder.myPntId[1]<<endl;
546       if(aNode[0] >= 0 && aNode[1] >= 0){
547         anIdList->SetId( 0, aNode[0] );
548         anIdList->SetId( 1, aNode[1] );
549         aConnectivity->InsertNextCell( anIdList );
550         aCellTypesArray->InsertNextValue( VTK_LINE );
551       }
552     }
553     
554     VTKViewer_CellLocationsArray* aCellLocationsArray = VTKViewer_CellLocationsArray::New();
555     aCellLocationsArray->SetNumberOfComponents( 1 );
556     aCellLocationsArray->SetNumberOfTuples( aNbCells );
557     
558     aConnectivity->InitTraversal();
559     for( vtkIdType idType = 0, *pts, npts; aConnectivity->GetNextCell( npts, pts ); idType++ )
560       aCellLocationsArray->SetValue( idType, aConnectivity->GetTraversalLocation( npts ) );
561     
562     aDataSet->SetCells( aCellTypesArray, aCellLocationsArray,aConnectivity );
563
564     SetUnstructuredGrid(aDataSet);
565     aDataSet->Delete();
566   }
567   else if(dynamic_cast<FreeNodes      *>(theFunctor.get()) ||
568           dynamic_cast<CoincidentNodes*>(theFunctor.get()))
569   {
570     Predicate* aPredicate = dynamic_cast<Predicate*>(theFunctor.get());
571     myExtractUnstructuredGrid->SetModeOfChanging(VTKViewer_ExtractUnstructuredGrid::eAdding);
572     vtkIdType aNbNodes = myVisualObj->GetNbEntities(SMDSAbs_Node);
573     for( vtkIdType i = 0; i < aNbNodes; i++ ){
574       vtkIdType anObjId = myVisualObj->GetNodeObjId(i);
575       if(aPredicate->IsSatisfy(anObjId))
576         myExtractUnstructuredGrid->RegisterCell(i);
577     }
578     if(!myExtractUnstructuredGrid->IsCellsRegistered())
579       myExtractUnstructuredGrid->RegisterCell(-1);
580     SetUnstructuredGrid(myVisualObj->GetUnstructuredGrid());
581   }
582 }
583
584
585
586
587 unsigned long int 
588 SMESH_DeviceActor
589 ::GetMTime()
590 {
591   unsigned long mTime = this->Superclass::GetMTime();
592   mTime = max(mTime,myExtractGeometry->GetMTime());
593   mTime = max(mTime,myExtractUnstructuredGrid->GetMTime());
594   mTime = max(mTime,myMergeFilter->GetMTime());
595   mTime = max(mTime,myGeomFilter->GetMTime());
596   mTime = max(mTime,myTransformFilter->GetMTime());
597   mTime = max(mTime,myFaceOrientationFilter->GetMTime());
598   return mTime;
599 }
600
601
602 void
603 SMESH_DeviceActor
604 ::SetTransform(VTKViewer_Transform* theTransform)
605 {
606   myTransformFilter->SetTransform(theTransform);
607 }
608
609
610 void
611 SMESH_DeviceActor
612 ::SetShrink() 
613 {
614   if ( !myIsShrinkable ) return;
615   if ( vtkAlgorithmOutput* aDataSet = myPassFilter[ 0 ]->GetOutputPort() )
616   {
617     myShrinkFilter->SetInputConnection( aDataSet );
618     myPassFilter[ 1 ]->SetInputConnection( myShrinkFilter->GetOutputPort() );
619     myIsShrunk = true;
620   }
621 }
622
623 void
624 SMESH_DeviceActor
625 ::UnShrink() 
626 {
627   if ( !myIsShrunk ) return;
628   if ( vtkAlgorithmOutput* aDataSet = myPassFilter[ 0 ]->GetOutputPort() )
629   {    
630     myPassFilter[ 1 ]->SetInputConnection( aDataSet );
631     myPassFilter[ 1 ]->Modified();
632     myIsShrunk = false;
633     Modified();
634   }
635 }
636
637
638 void
639 SMESH_DeviceActor
640 ::SetFacesOriented(bool theIsFacesOriented) 
641 {
642   if ( vtkAlgorithmOutput* aDataSet = myTransformFilter->GetOutputPort() )
643   {
644     myIsFacesOriented = theIsFacesOriented;
645     if( theIsFacesOriented )
646       myFaceOrientationFilter->SetInputConnection( aDataSet );
647     UpdateFaceOrientation();
648   }
649 }
650
651 void
652 SMESH_DeviceActor
653 ::SetFacesOrientationColor(double r,double g,double b)
654 {
655   myFaceOrientation->GetProperty()->SetColor( r, g, b );
656 }
657
658 void
659 SMESH_DeviceActor
660 ::GetFacesOrientationColor(double& r,double& g,double& b)
661 {
662   myFaceOrientation->GetProperty()->GetColor( r, g, b );
663 }
664
665 void
666 SMESH_DeviceActor
667 ::SetFacesOrientationScale(double theScale)
668 {
669   myFaceOrientationFilter->SetOrientationScale( theScale );
670 }
671
672 double
673 SMESH_DeviceActor
674 ::GetFacesOrientationScale()
675 {
676   return myFaceOrientationFilter->GetOrientationScale();
677 }
678
679 void
680 SMESH_DeviceActor
681 ::SetFacesOrientation3DVectors(bool theState)
682 {
683   myFaceOrientationFilter->Set3dVectors( theState );
684 }
685
686 bool
687 SMESH_DeviceActor
688 ::GetFacesOrientation3DVectors()
689 {
690   return myFaceOrientationFilter->Get3dVectors();
691 }
692
693 void
694 SMESH_DeviceActor
695 ::UpdateFaceOrientation()
696 {
697   bool aShowFaceOrientation = myIsFacesOriented;
698   aShowFaceOrientation &= GetVisibility();
699   aShowFaceOrientation &= myRepresentation == eSurface;
700   myFaceOrientation->SetVisibility(aShowFaceOrientation);
701 }
702
703
704 void
705 SMESH_DeviceActor
706 ::SetRepresentation(EReperesent theMode)
707 {
708   switch(theMode){
709   case ePoint:
710     myGeomFilter->SetInside(true);
711     myGeomFilter->SetWireframeMode(false);
712     GetProperty()->SetRepresentation(0);
713     break;
714   case eWireframe:
715     myGeomFilter->SetInside(false);
716     myGeomFilter->SetWireframeMode(true);
717     GetProperty()->SetRepresentation(theMode);
718     break;
719   case eInsideframe:
720     myGeomFilter->SetInside(true);
721     myGeomFilter->SetWireframeMode(true);
722     GetProperty()->SetRepresentation(1);
723     break;
724   case eSurface:
725     myGeomFilter->SetInside(false);
726     myGeomFilter->SetWireframeMode(false);
727     GetProperty()->SetRepresentation(theMode);
728   }
729   SetMarkerEnabled(theMode == ePoint);
730   myRepresentation = theMode;
731   UpdateFaceOrientation();
732   GetProperty()->Modified();
733   myMapper->Modified();
734   Modified();
735 }
736
737
738 void
739 SMESH_DeviceActor
740 ::SetVisibility(int theMode)
741 {
742   if(!myExtractUnstructuredGrid->GetInput() || 
743      GetUnstructuredGrid()->GetNumberOfCells())
744   {
745     vtkLODActor::SetVisibility(theMode);
746   }else{
747     vtkLODActor::SetVisibility(false);
748   }
749   UpdateFaceOrientation();
750 }
751
752
753 int
754 SMESH_DeviceActor
755 ::GetVisibility()
756 {
757   if(!GetUnstructuredGrid()->GetNumberOfCells()){
758     vtkLODActor::SetVisibility(false);
759   }
760   return vtkLODActor::GetVisibility();
761 }
762
763
764 void
765 SMESH_DeviceActor
766 ::AddToRender(vtkRenderer* theRenderer)
767 {
768   theRenderer->AddActor(this);
769   theRenderer->AddActor(myFaceOrientation);
770 }
771
772 void
773 SMESH_DeviceActor
774 ::RemoveFromRender(vtkRenderer* theRenderer)
775 {
776   theRenderer->RemoveActor(this);
777   theRenderer->RemoveActor(myFaceOrientation);
778 }
779
780
781 int
782 SMESH_DeviceActor
783 ::GetNodeObjId(int theVtkID)
784 {
785   vtkIdType anID = theVtkID;
786
787   if(IsImplicitFunctionUsed())
788     anID = myExtractGeometry->GetNodeObjId(theVtkID);
789
790   vtkIdType aRetID = myVisualObj->GetNodeObjId(anID);
791   if(MYDEBUG) MESSAGE("GetNodeObjId - theVtkID = "<<theVtkID<<"; anID = "<<anID<<"; aRetID = "<<aRetID);
792   return aRetID;
793 }
794
795 double* 
796 SMESH_DeviceActor
797 ::GetNodeCoord(int theObjID)
798 {
799   vtkDataSet* aDataSet = myMergeFilter->GetOutput();
800   vtkIdType anID = myVisualObj->GetNodeVTKId(theObjID);
801   double* aCoord = (anID >=0) ? aDataSet->GetPoint(anID) : NULL;
802   if(MYDEBUG) MESSAGE("GetNodeCoord - theObjID = "<<theObjID<<"; anID = "<<anID);
803   return aCoord;
804 }
805
806
807 int
808 SMESH_DeviceActor
809 ::GetElemObjId(int theVtkID)
810 {
811   vtkIdType anId = myGeomFilter->GetElemObjId(theVtkID);
812   if(anId < 0) 
813     return -1;
814
815   vtkIdType anId2 = anId;
816   if(IsImplicitFunctionUsed())
817     anId2 = myExtractGeometry->GetElemObjId(anId);
818   if(anId2 < 0) 
819     return -1;
820
821   vtkIdType anId3 = myExtractUnstructuredGrid->GetInputId(anId2);
822   if(anId3 < 0) 
823     return -1;
824
825   vtkIdType aRetID = myVisualObj->GetElemObjId(anId3);
826   if(MYDEBUG) 
827      MESSAGE("GetElemObjId - theVtkID = "<<theVtkID<<"; anId2 = "<<anId2<<"; anId3 = "<<anId3<<"; aRetID = "<<aRetID);
828   return aRetID;
829 }
830
831 vtkCell* 
832 SMESH_DeviceActor
833 ::GetElemCell(int theObjID)
834 {
835   vtkDataSet* aDataSet = myVisualObj->GetUnstructuredGrid();
836   vtkIdType aGridID = myVisualObj->GetElemVTKId(theObjID);
837   vtkCell* aCell = (aGridID >= 0 ) ? aDataSet->GetCell(aGridID) : NULL;
838   if(MYDEBUG) 
839     MESSAGE("GetElemCell - theObjID = "<<theObjID<<"; aGridID = "<<aGridID);
840   return aCell;
841 }
842
843
844 double 
845 SMESH_DeviceActor
846 ::GetShrinkFactor()
847 {
848   return myShrinkFilter->GetShrinkFactor();
849 }
850
851 void
852 SMESH_DeviceActor
853 ::SetShrinkFactor(double theValue)
854 {
855   theValue = theValue > 0.1? theValue: 0.8;
856   myShrinkFilter->SetShrinkFactor(theValue);
857   Modified();
858 }
859
860
861 void
862 SMESH_DeviceActor
863 ::SetHighlited(bool theIsHighlited)
864 {
865   if ( myIsHighlited == theIsHighlited )
866     return;
867   myIsHighlited = theIsHighlited;
868   Modified();
869 }
870
871 void
872 SMESH_DeviceActor
873 ::Render(vtkRenderer *ren, vtkMapper* m)
874 {
875   int aResolveCoincidentTopology = vtkMapper::GetResolveCoincidentTopology();
876   double aStoredFactor, aStoredUnit; 
877   vtkMapper::GetResolveCoincidentTopologyPolygonOffsetParameters(aStoredFactor,aStoredUnit);
878
879   vtkMapper::SetResolveCoincidentTopologyToPolygonOffset();
880   double aFactor = myPolygonOffsetFactor, aUnits = myPolygonOffsetUnits;
881   if(myIsHighlited){
882     static double EPS = .01;
883     aUnits *= (1.0-EPS);
884   }
885   vtkMapper::SetResolveCoincidentTopologyPolygonOffsetParameters(aFactor,aUnits);
886   vtkLODActor::Render(ren,m);
887
888   vtkMapper::SetResolveCoincidentTopologyPolygonOffsetParameters(aStoredFactor,aStoredUnit);
889   vtkMapper::SetResolveCoincidentTopology(aResolveCoincidentTopology);
890 }
891
892
893 void
894 SMESH_DeviceActor
895 ::SetPolygonOffsetParameters(double factor, 
896                              double units)
897 {
898   myPolygonOffsetFactor = factor;
899   myPolygonOffsetUnits = units;
900 }
901
902 /*!
903  * On/Off representation 2D quadratic element as arked polygon
904  */
905 void SMESH_DeviceActor::SetQuadraticArcMode(bool theFlag){
906   myGeomFilter->SetQuadraticArcMode(theFlag);
907 }
908
909 /*!
910  * Return true if 2D quadratic element displayed as arked polygon
911  */
912 bool SMESH_DeviceActor::GetQuadraticArcMode(){
913   return myGeomFilter->GetQuadraticArcMode();
914 }
915 /*!
916  * Set Max angle for representation 2D quadratic element as arked polygon
917  */
918 void SMESH_DeviceActor::SetQuadraticArcAngle(double theMaxAngle){
919   myGeomFilter->SetQuadraticArcAngle(theMaxAngle);
920 }
921
922 /*!
923  * Return Max angle of the representation 2D quadratic element as arked polygon
924  */
925 double SMESH_DeviceActor::GetQuadraticArcAngle(){
926   return myGeomFilter->GetQuadraticArcAngle();
927 }
928
929 /*!
930  * Set point marker enabled
931  * \param theMarkerEnabled flag to enable/disable point marker
932  */
933 void SMESH_DeviceActor::SetMarkerEnabled( bool theMarkerEnabled )
934 {
935   myMapper->SetMarkerEnabled( theMarkerEnabled );
936 }
937
938 /*!
939  * Set standard point marker
940  * \param theMarkerType type of the marker
941  */
942 void SMESH_DeviceActor::SetMarkerStd( VTK::MarkerType theMarkerType, VTK::MarkerScale theMarkerScale )
943 {
944   myMapper->SetMarkerStd( theMarkerType, theMarkerScale );
945 }
946
947 /*!
948  * Set custom point marker
949  * \param theMarkerId id of the marker texture
950  * \param theMarkerTexture marker texture
951  */
952 void SMESH_DeviceActor::SetMarkerTexture( int theMarkerId, VTK::MarkerTexture theMarkerTexture )
953 {
954   myMapper->SetMarkerTexture( theMarkerId, theMarkerTexture );
955 }
956
957 /*!
958  * Get type of the point marker
959  * \return type of the point marker
960  */
961 VTK::MarkerType SMESH_DeviceActor::GetMarkerType()
962 {
963   return myMapper->GetMarkerType();
964 }
965
966 /*!
967   Get scale of the point marker
968   \return scale of the point marker
969 */
970 VTK::MarkerScale SMESH_DeviceActor::GetMarkerScale()
971 {
972   return myMapper->GetMarkerScale();
973 }
974
975 /*!
976  * Get texture identifier of the point marker
977  * \return texture identifier of the point marker
978  */
979 int SMESH_DeviceActor::GetMarkerTexture()
980 {
981   return myMapper->GetMarkerTexture();
982 }
983
984 void SMESH_DeviceActor::SetCoincident3DAllowed(bool theFlag) {
985   myGeomFilter->SetAppendCoincident3D(theFlag);
986 }
987
988 bool SMESH_DeviceActor::IsCoincident3DAllowed() const {
989   return myGeomFilter->GetAppendCoincident3D();
990 }