OGS
VtkVisPointSetItem.cpp
Go to the documentation of this file.
1 
15 // ** INCLUDES **
16 #include "VtkVisPointSetItem.h"
17 
18 #include <vtkActor.h>
19 #include <vtkCellData.h>
20 #include <vtkDataSetMapper.h>
21 #include <vtkImageAlgorithm.h>
22 #include <vtkLookupTable.h>
23 #include <vtkPointData.h>
24 #include <vtkProperty.h>
25 #include <vtkRenderer.h>
26 #include <vtkSmartPointer.h>
27 #include <vtkTransform.h>
28 #include <vtkTransformFilter.h>
29 
30 #include <QObject>
31 #include <QRegExp>
32 #include <QSettings>
33 #include <QStringList>
34 #include <limits>
35 
36 #include "BaseLib/FileTools.h"
37 #include "BaseLib/Logging.h"
39 #include "QVtkDataSetMapper.h"
40 #include "VtkAlgorithmProperties.h"
42 #include "VtkCompositeFilter.h"
44 
45 // export test
46 #include <vtkDataSetAttributes.h>
47 #include <vtkPolyDataAlgorithm.h>
48 #include <vtkTriangleFilter.h>
49 #include <vtkTubeFilter.h>
50 #include <vtkUnstructuredGridAlgorithm.h>
51 #include <vtkXMLPolyDataWriter.h>
52 #include <vtkXMLUnstructuredGridWriter.h>
53 
55  vtkAlgorithm* algorithm, TreeItem* parentItem,
56  const QList<QVariant> data /*= QList<QVariant>()*/)
57  : VtkVisPipelineItem(algorithm, parentItem, data),
58  _mapper(nullptr),
59  _transformFilter(nullptr),
60  _onPointData(true),
61  _activeArrayName("")
62 {
63  auto* visParentItem = dynamic_cast<VtkVisPipelineItem*>(parentItem);
64  if (parentItem->parentItem())
65  {
66  // special case if parent is image but child is not (e.g.
67  // Image2BarChartFilter)
68  if (dynamic_cast<vtkImageAlgorithm*>(visParentItem->algorithm()))
69  {
70  _algorithm->SetInputConnection(
71  visParentItem->algorithm()->GetOutputPort());
72  }
73  else
74  {
75  auto* pointSetItem = dynamic_cast<VtkVisPointSetItem*>(parentItem);
76  if (pointSetItem)
77  {
78  _algorithm->SetInputConnection(
79  pointSetItem->transformFilter()->GetOutputPort());
80  }
81  }
82  }
83 }
84 
86  VtkCompositeFilter* compositeFilter, TreeItem* parentItem,
87  const QList<QVariant> data /*= QList<QVariant>()*/)
88  : VtkVisPipelineItem(compositeFilter, parentItem, data),
89  _mapper(nullptr),
90  _transformFilter(nullptr),
91  _onPointData(true),
92  _activeArrayName("")
93 {
94 }
95 
97 {
98  _transformFilter->Delete();
99  _mapper->Delete();
100 }
102 {
103  return _vtkProps->GetActiveAttribute();
104 }
105 
106 void VtkVisPointSetItem::Initialize(vtkRenderer* renderer)
107 {
108  // TODO vtkTransformFilter creates a new copy of the point coordinates which
109  // conflicts with VtkMappedMeshSource. Find a workaround!
110  _transformFilter = vtkTransformFilter::New();
111  vtkSmartPointer<vtkTransform> transform =
112  vtkSmartPointer<vtkTransform>::New();
113  transform->Identity();
114  _transformFilter->SetTransform(transform);
115 
116  _transformFilter->SetInputConnection(_algorithm->GetOutputPort());
117  _transformFilter->Update();
118 
119  _renderer = renderer;
121  _mapper->InterpolateScalarsBeforeMappingOff();
122  _mapper->SetColorModeToMapScalars();
123 
124  _mapper->SetInputConnection(_transformFilter->GetOutputPort());
125  _actor = vtkActor::New();
126  static_cast<vtkActor*>(_actor)->SetMapper(_mapper);
127  _renderer->AddActor(_actor);
128 
129  // Determine the right pre-set properties
130  // Order is: _algorithm, _compositeFilter, create a new one with props
131  // copied from parent
132  auto* vtkProps = dynamic_cast<VtkAlgorithmProperties*>(_algorithm);
133  if (!vtkProps)
134  {
135  vtkProps = dynamic_cast<VtkAlgorithmProperties*>(_compositeFilter);
136 
137  // Copy properties from parent or create a new VtkAlgorithmProperties
138  if (!vtkProps)
139  {
140  auto* parentItem =
141  dynamic_cast<VtkVisPipelineItem*>(this->parentItem());
142  while (parentItem)
143  {
144  VtkAlgorithmProperties* parentProps = nullptr;
145  if (dynamic_cast<VtkVisPointSetItem*>(parentItem))
146  {
147  parentProps = dynamic_cast<VtkVisPointSetItem*>(parentItem)
148  ->getVtkProperties();
149  }
150  if (parentProps)
151  {
152  vtkProps =
153  new VtkAlgorithmProperties(); // TODO memory leak?
154  vtkProps->SetScalarVisibility(
155  parentProps->GetScalarVisibility());
156  vtkProps->SetTexture(parentProps->GetTexture());
157  vtkProps->SetActiveAttribute(
158  parentProps->GetActiveAttribute());
159  parentItem = nullptr;
160  }
161  else
162  {
163  parentItem = dynamic_cast<VtkVisPipelineItem*>(
165  }
166  }
167 
168  // Has no parents
169  if (!vtkProps)
170  {
171  vtkProps = new VtkAlgorithmProperties(); // TODO memory leak?
172  }
173  }
174  }
175  _vtkProps = vtkProps;
176 
177  if (vtkProps->GetActiveAttribute().length() == 0)
178  {
179  // Get first scalar and set it to active
180  QStringList arrayNames = this->getScalarArrayNames();
181  if (arrayNames.length() > 0)
182  {
183  vtkProps->SetActiveAttribute(arrayNames[0]);
184  }
185  else
186  {
187  vtkProps->SetActiveAttribute("Solid Color");
188  }
189  }
190  this->setVtkProperties(vtkProps);
191  this->SetActiveAttribute(vtkProps->GetActiveAttribute());
192 
193  // Set global backface culling
194  QSettings settings;
195  bool backfaceCulling = settings.value("globalCullBackfaces", 0).toBool();
196  this->setBackfaceCulling(backfaceCulling);
197 
198  // Set the correct threshold range
199  if (dynamic_cast<VtkCompositeThresholdFilter*>(this->_compositeFilter))
200  {
201  double range[2];
202  this->GetRangeForActiveAttribute(range);
203  QList<QVariant> thresholdRangeList;
204  thresholdRangeList.push_back(range[0]);
205  thresholdRangeList.push_back(range[1]);
206  dynamic_cast<VtkCompositeFilter*>(this->_compositeFilter)
207  ->SetUserVectorProperty("Range", thresholdRangeList);
208  }
209 
210  // Show edges on meshes
211  if (dynamic_cast<MeshLib::VtkMappedMeshSource*>(this->_algorithm))
212  {
213  _vtkProps->GetProperties()->SetEdgeVisibility(1);
214  }
215 }
216 
218 {
220 }
221 
223 {
224  QObject::connect(vtkProps, SIGNAL(ScalarVisibilityChanged(bool)), _mapper,
225  SLOT(SetScalarVisibility(bool)));
226 
227  auto* actor = dynamic_cast<vtkActor*>(_actor);
228  if (actor)
229  {
230  if (vtkProps->GetTexture() != nullptr)
231  {
232  vtkProps->SetScalarVisibility(false);
233  actor->GetProperty()->SetColor(1, 1, 1); // don't colorise textures
234  actor->SetTexture(vtkProps->GetTexture());
235  }
236  else
237  {
238  vtkSmartPointer<vtkProperty> itemProperty =
239  vtkProps->GetProperties();
240  actor->SetProperty(itemProperty);
241  }
242 
243  if (!vtkProps->GetScalarVisibility())
244  {
245  vtkProps->SetScalarVisibility(false);
246  }
247  }
248 }
249 
250 int VtkVisPointSetItem::callVTKWriter(vtkAlgorithm* algorithm,
251  const std::string& filename) const
252 {
253  std::string file_name_cpy(filename);
254  auto* algPD = dynamic_cast<vtkPolyDataAlgorithm*>(algorithm);
255  if (algPD)
256  {
257  vtkSmartPointer<vtkXMLPolyDataWriter> pdWriter =
258  vtkSmartPointer<vtkXMLPolyDataWriter>::New();
259  pdWriter->SetInputData(algPD->GetOutputDataObject(0));
260  if (BaseLib::getFileExtension(filename) != ".vtp")
261  {
262  file_name_cpy.append(".vtp");
263  }
264  pdWriter->SetFileName(file_name_cpy.c_str());
265  return pdWriter->Write();
266  }
267 
268  auto* algUG = dynamic_cast<vtkUnstructuredGridAlgorithm*>(algorithm);
269  if (algUG)
270  {
271  vtkSmartPointer<vtkXMLUnstructuredGridWriter> ugWriter =
272  vtkSmartPointer<vtkXMLUnstructuredGridWriter>::New();
273  ugWriter->SetInputData(algUG->GetOutputDataObject(0));
274  if (BaseLib::getFileExtension(filename) != ".vtu")
275  {
276  file_name_cpy.append(".vtu");
277  }
278  ugWriter->SetFileName(file_name_cpy.c_str());
279  return ugWriter->Write();
280  }
281 
282  WARN("VtkVisPipelineItem::writeToFile(): Unknown data type.");
283  return 0;
284 }
285 
287 {
288  // Get type by identifier
289  if (name.contains(QRegExp("^P-")))
290  {
291  _onPointData = true;
292  }
293  else if (name.contains(QRegExp("^C-")))
294  {
295  _onPointData = false;
296  }
297  else if (name.contains("Solid Color"))
298  {
299  _vtkProps->SetActiveAttribute("Solid Color");
300  _mapper->ScalarVisibilityOff();
301  return;
302  }
303  else
304  {
305  return;
306  }
307 
308  // Remove type identifier
309  _activeArrayName = QString(name).remove(0, 2).toStdString();
310 
311  vtkDataSet* dataSet =
312  vtkDataSet::SafeDownCast(this->_algorithm->GetOutputDataObject(0));
313  if (!dataSet)
314  {
315  return;
316  }
317 
318  double range[2];
320  if (_onPointData)
321  {
322  vtkPointData* pointData = dataSet->GetPointData();
323  if (pointData)
324  {
325  _algorithm->SetInputArrayToProcess(
326  0, 0, 0, vtkDataObject::FIELD_ASSOCIATION_POINTS,
327  _activeArrayName.c_str());
328  _mapper->SetScalarModeToUsePointFieldData();
329  }
330  }
331  else
332  {
333  vtkCellData* cellData = dataSet->GetCellData();
334  if (cellData)
335  {
336  _algorithm->SetInputArrayToProcess(
337  0, 0, 0, vtkDataObject::FIELD_ASSOCIATION_CELLS,
338  _activeArrayName.c_str());
339  _mapper->SetScalarModeToUseCellFieldData();
340  }
341  }
342 
344 
345  _mapper->ScalarVisibilityOn();
346  _mapper->UseLookupTableScalarRangeOn();
347 
348  auto* mapper = dynamic_cast<QVtkDataSetMapper*>(_mapper);
349  if (mapper)
350  {
351  // Create a default color table when there is no lookup table for this
352  // attribute
353  vtkLookupTable* lut = _vtkProps->GetLookupTable(name);
354  if (lut == nullptr)
355  {
356  // std::cout << "Creating new lookup table for: " <<
357  // name.toStdString() << std::endl;
358  lut = vtkLookupTable::New(); // is not a memory leak, gets deleted
359  // in VtkAlgorithmProperties
360  lut->SetTableRange(range);
362  }
363  _mapper->SetLookupTable(lut);
364  }
365  _mapper->SelectColorArray(_activeArrayName.c_str());
366 }
367 
368 bool VtkVisPointSetItem::activeAttributeExists(vtkDataSetAttributes* data,
369  std::string& name)
370 {
371  bool arrayFound = false;
372  for (int i = 0; i < data->GetNumberOfArrays() && !arrayFound; i++)
373  {
374  std::string arrayName = data->GetArrayName(i);
375  if (arrayName == name)
376  {
377  arrayFound = true;
378  }
379  }
380  if (arrayFound)
381  {
382  // TODO Necessary? Currently this function is not called
383  data->SetActiveAttribute(name.c_str(), vtkDataSetAttributes::SCALARS);
384  return true;
385  }
386 
387  return false;
388 }
389 
390 void VtkVisPointSetItem::setScale(double x, double y, double z) const
391 {
392  if (this->transformFilter())
393  {
394  auto* transform =
395  static_cast<vtkTransform*>(this->_transformFilter->GetTransform());
396  double* trans = transform->GetPosition();
397  transform->Identity();
398  transform->Scale(x, y, z);
399  transform->Translate(trans[0] / x, trans[1] / y, trans[2] / z);
400  this->transformFilter()->Modified();
401  }
402 }
403 
404 void VtkVisPointSetItem::setTranslation(double x, double y, double z) const
405 {
406  if (this->transformFilter())
407  {
408  auto* transform =
409  static_cast<vtkTransform*>(this->_transformFilter->GetTransform());
410  double* scale = transform->GetScale();
411  transform->Identity();
412  transform->Scale(scale);
413  transform->Translate(x, y, z);
414  this->transformFilter()->Modified();
415  }
416 }
417 
419 {
420  return _transformFilter;
421 }
422 
424 {
425  static_cast<vtkActor*>(this->_actor)
426  ->GetProperty()
427  ->SetBackfaceCulling(static_cast<int>(enable));
428 }
429 
431 {
432  vtkDataSet* dataSet =
433  vtkDataSet::SafeDownCast(this->_algorithm->GetOutputDataObject(0));
434  if (dataSet && _activeArrayName.length() > 0)
435  {
436  if (_onPointData)
437  {
438  vtkPointData* pointData = dataSet->GetPointData();
439  if (pointData)
440  {
441  pointData->GetArray(_activeArrayName.c_str())->GetRange(range);
442  }
443  }
444  else
445  {
446  vtkCellData* cellData = dataSet->GetCellData();
447  cellData->GetArray(_activeArrayName.c_str())->GetRange(range);
448  }
449  }
450 }
Filename manipulation routines.
void WARN(char const *fmt, Args const &... args)
Definition: Logging.h:37
Definition of the QVtkDataSetMapper class.
Definition of the VtkAlgorithmProperties class.
Definition of the VtkCompositeContourFilter class.
Definition of the VtkCompositeFilter class.
Definition of the VtkCompositeThresholdFilter class.
VtkMappedMeshSource is a source class to transform OGS meshes into complete vtkUnstructuredGrids....
Definition of the VtkVisPointSetItem class.
Simply wraps vtkDataSetMapper as a Qt object to enable slot connections.
virtual void SetScalarVisibility(bool on)
Sets the scalar visibility on this mapper.
static QVtkDataSetMapper * New()
Create new objects with New() because of VTKs reference counting.
Objects nodes for the TreeModel.
Definition: TreeItem.h:28
TreeItem * parentItem() const
Definition: TreeItem.cpp:115
Contains properties for the visualization of objects as VtkVisPipelineItems.
vtkLookupTable * GetLookupTable(const QString &array_name)
Returns the colour lookup table (if one has been assigned).
void SetActiveAttribute(QString name)
Set the active attribute.
bool GetScalarVisibility() const
Returns the scalar visibility.
void SetScalarVisibility(bool on)
Sets the scalar visibility.
QString GetActiveAttribute() const
Returns the desired active attribute.
void SetLookUpTable(const QString &array_name, vtkLookupTable *lut)
Sets a colour lookup table for the given scalar array of the VtkVisPipelineItem.
vtkTexture * GetTexture()
Returns a texture (if one has been assigned).
vtkProperty * GetProperties() const
Returns the vtk properties.
Is used to combine several filter in one VtkVisPipelineItem. You can use vtk filter and custom filter...
Visualises only parts of meshes that are above/below/within given thresholds. In init() the threshold...
An item in the VtkVisPipeline containing a graphic object to be visualized.
vtkProp3D * actor() const
Returns the actor as vtkProp3D.
QStringList getScalarArrayNames() const
Returns a list of array names prefixed with P- or C- for point and cell data.
vtkAlgorithm * algorithm() const
Returns the algorithm object.
VtkAlgorithmProperties * _vtkProps
The active VtkAlgorithmProperties. From algorithm, compositeFilter, or copied from parent.
VtkAlgorithmProperties * getVtkProperties() const
Returns the VtkAlgorithmProperties.
vtkAlgorithm * _algorithm
VtkCompositeFilter * _compositeFilter
QVariant data(int column) const override
An item in the VtkVisPipeline containing a point set object to be visualized.
void setTranslation(double x, double y, double z) const override
Translates the item in visualisation-space.
void GetRangeForActiveAttribute(double range[2]) const
Get the scalar range for the active attribute.
VtkVisPointSetItem(vtkAlgorithm *algorithm, TreeItem *parentItem, const QList< QVariant > data=QList< QVariant >())
Constructor for a source/filter object.
void SetActiveAttribute(const QString &name) override
Sets the selected attribute array for the visualisation of the data set.
bool activeAttributeExists(vtkDataSetAttributes *data, std::string &name)
Checks if the selected attribute actually exists for the data set.
int callVTKWriter(vtkAlgorithm *algorithm, const std::string &filename) const override
Selects the appropriate VTK-Writer object and writes the object to a file with the given name.
void SetScalarVisibility(bool on)
std::string _activeArrayName
void setVtkProperties(VtkAlgorithmProperties *vtkProps)
Sets pre-set properties on vtkActor and on vtkMapper.
void setBackfaceCulling(bool enable) const override
Enables / disables backface culling.
QVtkDataSetMapper * _mapper
void setScale(double x, double y, double z) const override
Scales the data in visualisation-space.
QString GetActiveAttribute() const override
Gets the last selected attribute.
vtkTransformFilter * _transformFilter
vtkAlgorithm * transformFilter() const override
void Initialize(vtkRenderer *renderer) override
Initializes vtkMapper and vtkActor necessary for visualization of the item and sets the item's proper...
std::string getFileExtension(const std::string &path)
Definition: FileTools.cpp:186
void scale(PETScVector &x, double const a)
Definition: LinAlg.cpp:44