VTK  9.3.0
vtkMultiThreshold.h
Go to the documentation of this file.
1 // SPDX-FileCopyrightText: Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
2 // SPDX-FileCopyrightText: Copyright (c) Sandia Corporation
3 // SPDX-License-Identifier: BSD-3-Clause
4 
100 #ifndef vtkMultiThreshold_h
101 #define vtkMultiThreshold_h
102 
103 #include "vtkFiltersGeneralModule.h" // For export macro
104 #include "vtkMath.h" // for Inf() and NegInf()
106 
107 #include <map> // for IntervalRules map
108 #include <set> // for UpdateDependents()
109 #include <string> // for holding array names in NormKey
110 #include <vector> // for lists of threshold rules
111 
112 VTK_ABI_NAMESPACE_BEGIN
113 class vtkCell;
114 class vtkCellData;
115 class vtkDataArray;
116 class vtkGenericCell;
117 class vtkPointSet;
118 class vtkUnstructuredGrid;
119 
120 class VTKFILTERSGENERAL_EXPORT vtkMultiThreshold : public vtkMultiBlockDataSetAlgorithm
121 {
122 public:
125  void PrintSelf(ostream& os, vtkIndent indent) override;
126 
128  enum Closure
129  {
130  OPEN = 0,
131  CLOSED = 1
132  };
134  enum Norm
135  {
136  LINFINITY_NORM = -3,
137  L2_NORM = -2,
138  L1_NORM = -1
139  };
143  {
144  AND,
145  OR,
146  XOR,
147  WOR,
148  NAND
149  };
150 
152 
204  int AddIntervalSet(double xmin, double xmax, int omin, int omax, int assoc, const char* arrayName,
205  int component, int allScalars);
206  int AddIntervalSet(double xmin, double xmax, int omin, int omax, int assoc, int attribType,
207  int component, int allScalars);
209 
211 
218  int AddLowpassIntervalSet(
219  double xmax, int assoc, const char* arrayName, int component, int allScalars);
220  int AddHighpassIntervalSet(
221  double xmin, int assoc, const char* arrayName, int component, int allScalars);
222  int AddBandpassIntervalSet(
223  double xmin, double xmax, int assoc, const char* arrayName, int component, int allScalars);
224  int AddNotchIntervalSet(
225  double xlo, double xhi, int assoc, const char* arrayName, int component, int allScalars);
227 
231  int AddBooleanSet(int operation, int numInputs, int* inputs);
232 
236  int OutputSet(int setId);
237 
241  void Reset();
242 
245  typedef double (*TupleNorm)(vtkDataArray* arr, vtkIdType tuple, int component);
246 
247  // NormKey must be able to use TupleNorm typedef:
248  class NormKey;
249 
250  // Interval must be able to use NormKey typedef:
251  class Interval;
252 
253  // Set needs to refer to boolean set pointers
254  class BooleanSet;
255 
257  class NormKey
258  {
259  public:
260  int Association; // vtkDataObject::FIELD_ASSOCIATION_POINTS or
261  // vtkDataObject::FIELD_ASSOCIATION_CELLS
262  int Type; // -1 => use Name, otherwise: vtkDataSetAttributes::{SCALARS, VECTORS, TENSORS,
263  // NORMALS, TCOORDS, GLOBALIDS}
264  std::string Name; // Either empty or (when ArrayType == -1) an input array name
265  int Component; // LINFINITY_NORM, L1_NORM, L2_NORM or an integer component number
266  int AllScalars; // For Association == vtkDataObject::FIELD_ASSOCIATION_POINTS, must all points
267  // be in the interval?
268  int InputArrayIndex; // The number passed to vtkAlgorithm::SetInputArrayToProcess()
269  TupleNorm NormFunction; // A function pointer to compute the norm (or fetcht the correct
270  // component) of a tuple.
271 
275  vtkIdType cellId, vtkCell* cell, vtkDataArray* array, double cellNorm[2]) const;
276 
279  bool operator<(const NormKey& other) const
280  {
281  if (this->Association < other.Association)
282  return true;
283  else if (this->Association > other.Association)
284  return false;
285 
286  if (this->Component < other.Component)
287  return true;
288  else if (this->Component > other.Component)
289  return false;
290 
291  if ((!this->AllScalars) && other.AllScalars)
292  return true;
293  else if (this->AllScalars && (!other.AllScalars))
294  return false;
295 
296  if (this->Type == -1)
297  {
298  if (other.Type == -1)
299  return this->Name < other.Name;
300  return true;
301  }
302  else
303  {
304  return this->Type < other.Type;
305  }
306  }
307  };
308 
313  class Set
314  {
315  public:
316  int Id;
317  int OutputId;
319 
322  Set() { this->OutputId = -1; }
324  virtual ~Set() = default;
326  virtual void PrintNodeName(ostream& os);
328  virtual void PrintNode(ostream& os) = 0;
330  virtual BooleanSet* GetBooleanSetPointer();
331  virtual Interval* GetIntervalPointer();
332  };
333 
335  class Interval : public Set
336  {
337  public:
339  double EndpointValues[2];
341  int EndpointClosures[2];
344 
350  int Match(double cellNorm[2]);
351 
352  ~Interval() override = default;
353  void PrintNode(ostream& os) override;
354  Interval* GetIntervalPointer() override;
355  };
356 
358  class BooleanSet : public Set
359  {
360  public:
362  int Operator;
364  std::vector<int> Inputs;
365 
367  BooleanSet(int sId, int op, int* inBegin, int* inEnd)
368  : Inputs(inBegin, inEnd)
369  {
370  this->Id = sId;
371  this->Operator = op;
372  }
373  ~BooleanSet() override = default;
374  void PrintNode(ostream& os) override;
375  BooleanSet* GetBooleanSetPointer() override;
376  };
377 
378 protected:
380  ~vtkMultiThreshold() override;
381 
396  enum Ruling
397  {
398  INCONCLUSIVE = -1,
399  INCLUDE = -2,
400  EXCLUDE = -3
401  };
402 
407 
414 
421 
426 
428  typedef std::vector<Interval*> IntervalList;
430  typedef std::map<NormKey, IntervalList> RuleMap;
431 
432  typedef std::vector<int> TruthTreeValues;
433  typedef std::vector<TruthTreeValues> TruthTree;
434 
439 
445  std::vector<Set*> Sets;
446 
454 
458  void UpdateDependents(int id, std::set<int>& unresolvedOutputs, TruthTreeValues& setStates,
459  vtkCellData* inCellData, vtkIdType inCell, vtkGenericCell* cell,
460  std::vector<vtkUnstructuredGrid*>& outv);
461 
465  int AddIntervalSet(NormKey& nk, double xmin, double xmax, int omin, int omax);
466 
470  void PrintGraph(ostream& os);
471 
473  void operator=(const vtkMultiThreshold&) = delete;
474 };
475 
477  double xmax, int assoc, const char* arrayName, int component, int allScalars)
478 {
479  return this->AddIntervalSet(
480  vtkMath::NegInf(), xmax, CLOSED, CLOSED, assoc, arrayName, component, allScalars);
481 }
482 
484  double xmin, int assoc, const char* arrayName, int component, int allScalars)
485 {
486  return this->AddIntervalSet(
487  xmin, vtkMath::Inf(), CLOSED, CLOSED, assoc, arrayName, component, allScalars);
488 }
489 
491  double xmin, double xmax, int assoc, const char* arrayName, int component, int allScalars)
492 {
493  return this->AddIntervalSet(xmin, xmax, CLOSED, CLOSED, assoc, arrayName, component, allScalars);
494 }
495 
497  double xlo, double xhi, int assoc, const char* arrayName, int component, int allScalars)
498 {
499  int band =
500  this->AddIntervalSet(xlo, xhi, CLOSED, CLOSED, assoc, arrayName, component, allScalars);
501  if (band < 0)
502  {
503  return -1;
504  }
505  return this->AddBooleanSet(NAND, 1, &band);
506 }
507 
509 {
510  return nullptr;
511 }
512 
514 {
515  return nullptr;
516 }
517 
519 {
520  return this;
521 }
522 
524 {
525  return this;
526 }
527 
528 VTK_ABI_NAMESPACE_END
529 #endif // vtkMultiThreshold_h
represent and manipulate cell attribute data
Definition: vtkCellData.h:40
abstract class to specify cell behavior
Definition: vtkCell.h:59
abstract superclass for arrays of numeric data
Definition: vtkDataArray.h:54
provides thread-safe access to cells
a simple class to control print indentation
Definition: vtkIndent.h:38
Store zero or more vtkInformation instances.
Store vtkAlgorithm input/output information.
static double NegInf()
Special IEEE-754 number used to represent negative infinity.
static double Inf()
Special IEEE-754 number used to represent positive infinity.
Superclass for algorithms that produce only vtkMultiBlockDataSet as output.
A subset of a mesh represented as a boolean set operation.
int Operator
The boolean operation that will be performed on the inputs to obtain the output.
BooleanSet(int sId, int op, int *inBegin, int *inEnd)
Construct a new set with the given ID, operator, and inputs.
void PrintNode(ostream &os) override
Print a graphviz node name for use in an edge statement.
std::vector< int > Inputs
A list of input sets. These may be IntervalSets or BooleanSets.
~BooleanSet() override=default
BooleanSet * GetBooleanSetPointer() override
Avoid dynamic_casts. Subclasses must override.
A subset of a mesh represented by a range of acceptable attribute values.
void PrintNode(ostream &os) override
Print a graphviz node name for use in an edge statement.
~Interval() override=default
int Match(double cellNorm[2])
Does the specified range fall inside the interval? For cell-centered attributes, only cellNorm[0] is ...
NormKey Norm
This contains information about the attribute over which the interval is defined.
Interval * GetIntervalPointer() override
A class with comparison operator used to index input array norms used in threshold rules.
void ComputeNorm(vtkIdType cellId, vtkCell *cell, vtkDataArray *array, double cellNorm[2]) const
Compute the norm of a cell by calling NormFunction for all its points or for its single cell-centered...
bool operator<(const NormKey &other) const
A partial ordering of NormKey objects is required for them to serve as keys in the vtkMultiThreshold:...
A base class for representing threshold sets.
int OutputId
A unique identifier for this set.
virtual Interval * GetIntervalPointer()
virtual void PrintNodeName(ostream &os)
Print a graphviz node label statement (with fancy node name and shape).
virtual void PrintNode(ostream &os)=0
Print a graphviz node name for use in an edge statement.
virtual BooleanSet * GetBooleanSetPointer()
Avoid dynamic_casts. Subclasses must override.
Set()
The index of the output mesh that will hold this set or -1 if the set is not output.
virtual ~Set()=default
Virtual destructor since we have virtual members.
Threshold cells within multiple intervals.
~vtkMultiThreshold() override
TruthTree DependentSets
A list of boolean sets whose values depend on the given set.
int NumberOfOutputs
The number of output datasets.
vtkMultiThreshold(const vtkMultiThreshold &)=delete
void operator=(const vtkMultiThreshold &)=delete
Ruling
When an interval is evaluated, its value is used to update a truth table.
Norm
Norms that can be used to threshold vector attributes.
std::map< NormKey, IntervalList > RuleMap
A map describing the IntervalSets that share a common attribute and norm.
std::vector< Set * > Sets
A list of rules keyed by their unique integer ID.
std::vector< Interval * > IntervalList
A list of pointers to IntervalSets.
int NextArrayIndex
A variable used to store the next index to use when calling SetInputArrayToProcess.
RuleMap IntervalRules
A set of threshold rules sorted by the attribute+norm to which they are applied.
int AddBandpassIntervalSet(double xmin, double xmax, int assoc, const char *arrayName, int component, int allScalars)
These convenience members make it easy to insert closed intervals.
static vtkMultiThreshold * New()
std::vector< int > TruthTreeValues
void Reset()
Remove all the intervals currently defined.
int AddLowpassIntervalSet(double xmax, int assoc, const char *arrayName, int component, int allScalars)
These convenience members make it easy to insert closed intervals.
void PrintSelf(ostream &os, vtkIndent indent) override
Methods invoked by print to print information about the object including superclasses.
void PrintGraph(ostream &os)
Print out a graphviz-formatted text description of all the sets.
int RequestData(vtkInformation *, vtkInformationVector **, vtkInformationVector *) override
This function performs the actual thresholding.
int AddBooleanSet(int operation, int numInputs, int *inputs)
Create a new mesh subset using boolean operations on pre-existing sets.
int OutputSet(int setId)
Create an output mesh containing a boolean or interval subset of the input mesh.
void UpdateDependents(int id, std::set< int > &unresolvedOutputs, TruthTreeValues &setStates, vtkCellData *inCellData, vtkIdType inCell, vtkGenericCell *cell, std::vector< vtkUnstructuredGrid * > &outv)
Recursively update the setStates and unresolvedOutputs vectors based on this->DependentSets.
int AddIntervalSet(double xmin, double xmax, int omin, int omax, int assoc, int attribType, int component, int allScalars)
Add a mesh subset to be computed by thresholding an attribute of the input mesh.
int AddIntervalSet(double xmin, double xmax, int omin, int omax, int assoc, const char *arrayName, int component, int allScalars)
Add a mesh subset to be computed by thresholding an attribute of the input mesh.
int AddNotchIntervalSet(double xlo, double xhi, int assoc, const char *arrayName, int component, int allScalars)
These convenience members make it easy to insert closed intervals.
SetOperation
Operations that can be performed on sets to generate another set.
@ WOR
Include elements that belong to an odd number of input sets (a kind of "winding XOR")
@ XOR
Include an element if it belongs to exactly one input set.
@ NAND
Only include elements that don't belong to any input set.
@ AND
Only include an element if it belongs to all the input sets.
@ OR
Include an element if it belongs to any input set.
int AddHighpassIntervalSet(double xmin, int assoc, const char *arrayName, int component, int allScalars)
These convenience members make it easy to insert closed intervals.
Closure
Whether the endpoint value of an interval should be included or excluded.
@ CLOSED
Specify a closed interval.
int FillInputPortInformation(int port, vtkInformation *info) override
We accept any mesh that is descended from vtkPointSet.
std::vector< TruthTreeValues > TruthTree
int AddIntervalSet(NormKey &nk, double xmin, double xmax, int omin, int omax)
A utility method called by the public AddInterval members.
concrete class for storing a set of points
Definition: vtkPointSet.h:68
dataset represents arbitrary combinations of all possible cell types
@ component
Definition: vtkX3D.h:175
@ info
Definition: vtkX3D.h:376
@ port
Definition: vtkX3D.h:447
@ string
Definition: vtkX3D.h:490
int vtkIdType
Definition: vtkType.h:315