VTK  9.2.6
vtkVectorFieldTopology.h
Go to the documentation of this file.
1 /*=========================================================================
2 
3  Program: Visualization Toolkit
4  Module: vtkVectorFieldTopology.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 =========================================================================*/
30 #ifndef vtkVectorFieldTopology_h
31 #define vtkVectorFieldTopology_h
32 
33 #include "vtkFiltersFlowPathsModule.h" // For export macro
34 #include "vtkPolyDataAlgorithm.h"
35 #include "vtkStreamTracer.h" // for vtkStreamSurface::CELL_LENGTH_UNIT
36 
37 class vtkGradientFilter;
38 class vtkImageData;
39 class vtkPolyData;
40 class vtkStreamSurface;
42 
43 class VTKFILTERSFLOWPATHS_EXPORT vtkVectorFieldTopology : public vtkPolyDataAlgorithm
44 {
45 public:
46  static vtkVectorFieldTopology* New();
48  void PrintSelf(ostream& os, vtkIndent indent) override;
49 
51 
57  vtkSetMacro(IntegrationStepUnit, int);
58  vtkGetMacro(IntegrationStepUnit, int);
60 
62 
65  vtkSetMacro(MaxNumSteps, int);
66  vtkGetMacro(MaxNumSteps, int);
68 
70 
74  vtkSetMacro(IntegrationStepSize, double);
75  vtkGetMacro(IntegrationStepSize, double);
77 
79 
83  vtkSetMacro(SeparatrixDistance, double);
84  vtkGetMacro(SeparatrixDistance, double);
86 
88 
91  vtkSetMacro(UseIterativeSeeding, bool);
92  vtkGetMacro(UseIterativeSeeding, bool);
94 
96 
99  vtkSetMacro(ComputeSurfaces, bool);
100  vtkGetMacro(ComputeSurfaces, bool);
102 
104 
107  vtkSetMacro(ExcludeBoundary, bool);
108  vtkGetMacro(ExcludeBoundary, bool);
110 
112 
115  vtkSetMacro(UseBoundarySwitchPoints, bool);
116  vtkGetMacro(UseBoundarySwitchPoints, bool);
118 
120 
129  vtkSetMacro(VectorAngleThreshold, double);
130  vtkGetMacro(VectorAngleThreshold, double);
132 
134 
137  vtkSetMacro(OffsetAwayFromBoundary, double);
138  vtkGetMacro(OffsetAwayFromBoundary, double);
140 
142 
145  vtkSetMacro(EpsilonCriticalPoint, double);
146  vtkGetMacro(EpsilonCriticalPoint, double);
148 
154  void SetInterpolatorType(int interpType);
155 
159  void SetInterpolatorTypeToCellLocator();
160 
164  void SetInterpolatorTypeToDataSetPointLocator();
165 
166 protected:
168  ~vtkVectorFieldTopology() override;
169 
170  int FillInputPortInformation(int port, vtkInformation* info) override;
173 
174 private:
176  void operator=(const vtkVectorFieldTopology&) = delete;
177 
181  int Validate();
182 
190  int ImageDataPrepare(vtkDataSet* dataSetInput, vtkUnstructuredGrid* tridataset);
191 
199  int UnstructuredGridPrepare(vtkDataSet* dataSetInput, vtkUnstructuredGrid* tridataset);
200 
206  int RemoveBoundary(vtkSmartPointer<vtkUnstructuredGrid> tridataset);
207 
215  int ComputeCriticalPoints2D(
217 
225  int ComputeCriticalPoints3D(
227 
238  static void InterpolateVector(
239  double x0, double x1, double x, const double v0[3], const double v1[3], double v[3]);
240 
247  int ComputeBoundarySwitchPoints(
248  vtkPolyData* boundarySwitchPoints, vtkUnstructuredGrid* tridataset);
249 
265  int ComputeSeparatricesBoundarySwitchPoints(vtkPolyData* boundarySwitchPoints,
266  vtkPolyData* separatrices, vtkDataSet* dataset, vtkPoints* interestPoints,
267  int integrationStepUnit, double dist, int maxNumSteps);
268 
286  int ComputeSeparatricesBoundarySwitchLines(vtkPolyData* boundarySwitchLines,
287  vtkPolyData* separatrices, vtkDataSet* dataset, int integrationStepUnit, double dist,
288  int maxNumSteps, bool computeSurfaces, bool useIterativeSeeding);
289 
308  int ComputeSeparatrices(vtkPolyData* criticalPoints, vtkPolyData* separatrices,
309  vtkPolyData* surfaces, vtkDataSet* dataset, vtkPoints* interestPoints, int integrationStepUnit,
310  double dist, double stepSize, int maxNumSteps, bool computeSurfaces, bool useIterativeSeeding);
311 
328  int ComputeSurface(int numberOfSeparatingSurfaces, bool isBackward, double normal[3],
329  double zeroPos[3], vtkPolyData* streamSurfaces, vtkDataSet* dataset, int integrationStepUnit,
330  double dist, double stepSize, int maxNumSteps, bool useIterativeSeeding);
331 
336  enum CriticalType2D
337  {
338  DEGENERATE_2D = -1,
339  SINK_2D = 0,
340  SADDLE_2D = 1,
341  SOURCE_2D = 2,
342  CENTER_2D = 3
343  };
344 
349  enum CriticalTypeDetailed2D
350  {
351  // DEGENERATE2D = -1,
352  ATTRACTING_NODE_2D = 0,
353  ATTRACTING_FOCUS_2D = 1,
354  NODE_SADDLE_2D = 2,
355  REPELLING_NODE_2D = 3,
356  REPELLING_FOCUS_2D = 4,
357  CENTER_DETAILED_2D = 5
358  };
359 
364  enum CriticalType3D
365  {
366  DEGENERATE_3D = -1,
367  SINK_3D = 0,
368  SADDLE_1_3D = 1,
369  SADDLE_2_3D = 2,
370  SOURCE_3D = 3,
371  CENTER_3D = 4
372  };
373 
378  enum CriticalTypeDetailed3D
379  {
380  ATTRACTING_NODE_3D = 0,
381  ATTRACTING_FOCUS_3D = 1,
382  NODE_SADDLE_1_3D = 2,
383  FOCUS_SADDLE_1_3D = 3,
384  NODE_SADDLE_2_3D = 4,
385  FOCUS_SADDLE_2_3D = 5,
386  REPELLING_NODE_3D = 6,
387  REPELLING_FOCUS_3D = 7,
388  CENTER_DETAILED_3D = 8
389  };
390 
398  static int Classify2D(int countComplex, int countPos, int countNeg);
399 
408  static int ClassifyDetailed2D(int countComplex, int countPos, int countNeg);
409 
418  static int Classify3D(int countComplex, int countPos, int countNeg);
419 
429  static int ClassifyDetailed3D(int countComplex, int countPos, int countNeg);
430 
434  int MaxNumSteps = 100;
435 
439  double IntegrationStepSize = 1;
440 
444  double SeparatrixDistance = 1;
445 
449  bool UseIterativeSeeding = false;
450 
454  bool ComputeSurfaces = false;
455 
459  const char* NameOfVectorArray;
460 
465  bool ExcludeBoundary = false;
466 
470  int Dimension = 2;
471 
479  int IntegrationStepUnit = vtkStreamTracer::CELL_LENGTH_UNIT;
480 
487  bool UseBoundarySwitchPoints = false;
488 
494 
503  double VectorAngleThreshold = 1;
504 
512  double OffsetAwayFromBoundary = 1e-3;
513 
517  double EpsilonCriticalPoint = 1e-10;
518 
519  vtkNew<vtkStreamSurface> StreamSurface;
520 };
521 #endif
Store vtkAlgorithm input/output information.
abstract class to specify dataset behavior
Definition: vtkDataSet.h:62
virtual int RequestData(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector)
This is called by the superclass.
Extract the topological skeleton as output datasets.
Hold a reference to a vtkObjectBase instance.
Definition: vtkMeta.h:32
concrete dataset represents vertices, lines, polygons, and triangle strips
Definition: vtkPolyData.h:90
static vtkPolyDataAlgorithm * New()
Superclass for algorithms that produce only polydata as output.
a simple class to control print indentation
Definition: vtkIndent.h:39
topologically and geometrically regular array of data
Definition: vtkImageData.h:53
dataset represents arbitrary combinations of all possible cell types
A general filter for gradient estimation.
int FillOutputPortInformation(int port, vtkInformation *info) override
Fill the output port information objects for this algorithm.
int FillInputPortInformation(int port, vtkInformation *info) override
Fill the input port information objects for this algorithm.
Store zero or more vtkInformation instances.
Advect a stream surface in a vector field.
represent and manipulate 3D points
Definition: vtkPoints.h:39
void PrintSelf(ostream &os, vtkIndent indent) override
Methods invoked by print to print information about the object including superclasses.