VTK  9.2.6
vtkMultiThreshold.h
Go to the documentation of this file.
1/*=========================================================================
2
3 Program: Visualization Toolkit
4 Module: vtkMultiThreshold.h
5
6 Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
7 All rights reserved.
8 See Copyright.txt or http://www.kitware.com/Copyright.htm for details.
9
10 This software is distributed WITHOUT ANY WARRANTY; without even
11 the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
12 PURPOSE. See the above copyright notice for more information.
13
14=========================================================================*/
15/*----------------------------------------------------------------------------
16 Copyright (c) Sandia Corporation
17 See Copyright.txt or http://www.paraview.org/HTML/Copyright.html for details.
18----------------------------------------------------------------------------*/
19
106#ifndef vtkMultiThreshold_h
107#define vtkMultiThreshold_h
108
109#include "vtkFiltersGeneralModule.h" // For export macro
110#include "vtkMath.h" // for Inf() and NegInf()
112
113#include <map> // for IntervalRules map
114#include <set> // for UpdateDependents()
115#include <string> // for holding array names in NormKey
116#include <vector> // for lists of threshold rules
117
118class vtkCell;
119class vtkCellData;
120class vtkDataArray;
121class vtkGenericCell;
122class vtkPointSet;
124
125class VTKFILTERSGENERAL_EXPORT vtkMultiThreshold : public vtkMultiBlockDataSetAlgorithm
126{
127public:
130 void PrintSelf(ostream& os, vtkIndent indent) override;
131
134 {
135 OPEN = 0,
136 CLOSED = 1
137 };
139 enum Norm
140 {
141 LINFINITY_NORM = -3,
142 L2_NORM = -2,
143 L1_NORM = -1
144 };
148 {
153 NAND
154 };
155
157
209 int AddIntervalSet(double xmin, double xmax, int omin, int omax, int assoc, const char* arrayName,
210 int component, int allScalars);
211 int AddIntervalSet(double xmin, double xmax, int omin, int omax, int assoc, int attribType,
212 int component, int allScalars);
214
216
223 int AddLowpassIntervalSet(
224 double xmax, int assoc, const char* arrayName, int component, int allScalars);
225 int AddHighpassIntervalSet(
226 double xmin, int assoc, const char* arrayName, int component, int allScalars);
227 int AddBandpassIntervalSet(
228 double xmin, double xmax, int assoc, const char* arrayName, int component, int allScalars);
229 int AddNotchIntervalSet(
230 double xlo, double xhi, int assoc, const char* arrayName, int component, int allScalars);
232
236 int AddBooleanSet(int operation, int numInputs, int* inputs);
237
241 int OutputSet(int setId);
242
246 void Reset();
247
250 typedef double (*TupleNorm)(vtkDataArray* arr, vtkIdType tuple, int component);
251
252 // NormKey must be able to use TupleNorm typedef:
253 class NormKey;
254
255 // Interval must be able to use NormKey typedef:
256 class Interval;
257
258 // Set needs to refer to boolean set pointers
259 class BooleanSet;
260
263 {
264 public:
265 int Association; // vtkDataObject::FIELD_ASSOCIATION_POINTS or
266 // vtkDataObject::FIELD_ASSOCIATION_CELLS
267 int Type; // -1 => use Name, otherwise: vtkDataSetAttributes::{SCALARS, VECTORS, TENSORS,
268 // NORMALS, TCOORDS, GLOBALIDS}
269 std::string Name; // Either empty or (when ArrayType == -1) an input array name
270 int Component; // LINFINITY_NORM, L1_NORM, L2_NORM or an integer component number
271 int AllScalars; // For Association == vtkDataObject::FIELD_ASSOCIATION_POINTS, must all points
272 // be in the interval?
273 int InputArrayIndex; // The number passed to vtkAlgorithm::SetInputArrayToProcess()
274 TupleNorm NormFunction; // A function pointer to compute the norm (or fetcht the correct
275 // component) of a tuple.
276
280 vtkIdType cellId, vtkCell* cell, vtkDataArray* array, double cellNorm[2]) const;
281
284 bool operator<(const NormKey& other) const
285 {
286 if (this->Association < other.Association)
287 return true;
288 else if (this->Association > other.Association)
289 return false;
290
291 if (this->Component < other.Component)
292 return true;
293 else if (this->Component > other.Component)
294 return false;
295
296 if ((!this->AllScalars) && other.AllScalars)
297 return true;
298 else if (this->AllScalars && (!other.AllScalars))
299 return false;
300
301 if (this->Type == -1)
302 {
303 if (other.Type == -1)
304 return this->Name < other.Name;
305 return true;
306 }
307 else
308 {
309 return this->Type < other.Type;
310 }
311 }
312 };
313
318 class Set
319 {
320 public:
321 int Id;
324
327 Set() { this->OutputId = -1; }
329 virtual ~Set() = default;
331 virtual void PrintNodeName(ostream& os);
333 virtual void PrintNode(ostream& os) = 0;
335 virtual BooleanSet* GetBooleanSetPointer();
336 virtual Interval* GetIntervalPointer();
337 };
338
340 class Interval : public Set
341 {
342 public:
344 double EndpointValues[2];
346 int EndpointClosures[2];
349
355 int Match(double cellNorm[2]);
356
357 ~Interval() override = default;
358 void PrintNode(ostream& os) override;
359 Interval* GetIntervalPointer() override;
360 };
361
363 class BooleanSet : public Set
364 {
365 public:
369 std::vector<int> Inputs;
370
372 BooleanSet(int sId, int op, int* inBegin, int* inEnd)
373 : Inputs(inBegin, inEnd)
374 {
375 this->Id = sId;
376 this->Operator = op;
377 }
378 ~BooleanSet() override = default;
379 void PrintNode(ostream& os) override;
380 BooleanSet* GetBooleanSetPointer() override;
381 };
382
383protected:
386
402 {
403 INCONCLUSIVE = -1,
404 INCLUDE = -2,
405 EXCLUDE = -3
406 };
407
412
418 int FillInputPortInformation(int port, vtkInformation* info) override;
419
426
431
433 typedef std::vector<Interval*> IntervalList;
435 typedef std::map<NormKey, IntervalList> RuleMap;
436
437 typedef std::vector<int> TruthTreeValues;
438 typedef std::vector<TruthTreeValues> TruthTree;
439
444
450 std::vector<Set*> Sets;
451
459
463 void UpdateDependents(int id, std::set<int>& unresolvedOutputs, TruthTreeValues& setStates,
464 vtkCellData* inCellData, vtkIdType inCell, vtkGenericCell* cell,
465 std::vector<vtkUnstructuredGrid*>& outv);
466
470 int AddIntervalSet(NormKey& nk, double xmin, double xmax, int omin, int omax);
471
475 void PrintGraph(ostream& os);
476
478 void operator=(const vtkMultiThreshold&) = delete;
479};
480
482 double xmax, int assoc, const char* arrayName, int component, int allScalars)
483{
484 return this->AddIntervalSet(
485 vtkMath::NegInf(), xmax, CLOSED, CLOSED, assoc, arrayName, component, allScalars);
486}
487
489 double xmin, int assoc, const char* arrayName, int component, int allScalars)
490{
491 return this->AddIntervalSet(
492 xmin, vtkMath::Inf(), CLOSED, CLOSED, assoc, arrayName, component, allScalars);
493}
494
496 double xmin, double xmax, int assoc, const char* arrayName, int component, int allScalars)
497{
498 return this->AddIntervalSet(xmin, xmax, CLOSED, CLOSED, assoc, arrayName, component, allScalars);
499}
500
502 double xlo, double xhi, int assoc, const char* arrayName, int component, int allScalars)
503{
504 int band =
505 this->AddIntervalSet(xlo, xhi, CLOSED, CLOSED, assoc, arrayName, component, allScalars);
506 if (band < 0)
507 {
508 return -1;
509 }
510 return this->AddBooleanSet(NAND, 1, &band);
511}
512
517
522
527
532
533#endif // vtkMultiThreshold_h
represent and manipulate cell attribute data
Definition vtkCellData.h:42
abstract class to specify cell behavior
Definition vtkCell.h:61
abstract superclass for arrays of numeric data
provides thread-safe access to cells
a simple class to control print indentation
Definition vtkIndent.h:40
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
static vtkMultiThreshold * New()
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::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.
std::vector< TruthTreeValues > TruthTree
std::vector< int > TruthTreeValues
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.
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::map< NormKey, IntervalList > RuleMap
A map describing the IntervalSets that share a common attribute and norm.
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:70
dataset represents arbitrary combinations of all possible cell types
int vtkIdType
Definition vtkType.h:332