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=========================================================================*/
27#ifndef vtkVectorFieldTopology_h
28#define vtkVectorFieldTopology_h
29
30#include "vtkFiltersFlowPathsModule.h" // For export macro
32#include "vtkStreamTracer.h" // for vtkStreamSurface::CELL_LENGTH_UNIT
33
35class vtkImageData;
36class vtkPolyData;
39
40class VTKFILTERSFLOWPATHS_EXPORT vtkVectorFieldTopology : public vtkPolyDataAlgorithm
41{
42public:
45 void PrintSelf(ostream& os, vtkIndent indent) override;
46
48
54 vtkSetMacro(IntegrationStepUnit, int);
55 vtkGetMacro(IntegrationStepUnit, int);
57
59
62 vtkSetMacro(MaxNumSteps, int);
63 vtkGetMacro(MaxNumSteps, int);
65
67
71 vtkSetMacro(IntegrationStepSize, double);
72 vtkGetMacro(IntegrationStepSize, double);
74
76
80 vtkSetMacro(SeparatrixDistance, double);
81 vtkGetMacro(SeparatrixDistance, double);
83
85
88 vtkSetMacro(UseIterativeSeeding, bool);
89 vtkGetMacro(UseIterativeSeeding, bool);
91
93
96 vtkSetMacro(ComputeSurfaces, bool);
97 vtkGetMacro(ComputeSurfaces, bool);
99
101
104 vtkSetMacro(ExcludeBoundary, bool);
105 vtkGetMacro(ExcludeBoundary, bool);
107
109
112 vtkSetMacro(UseBoundarySwitchPoints, bool);
113 vtkGetMacro(UseBoundarySwitchPoints, bool);
115
117
126 vtkSetMacro(VectorAngleThreshold, double);
127 vtkGetMacro(VectorAngleThreshold, double);
129
131
134 vtkSetMacro(OffsetAwayFromBoundary, double);
135 vtkGetMacro(OffsetAwayFromBoundary, double);
137
139
142 vtkSetMacro(EpsilonCriticalPoint, double);
143 vtkGetMacro(EpsilonCriticalPoint, double);
145
151 void SetInterpolatorType(int interpType);
152
157
162
163protected:
166
167 int FillInputPortInformation(int port, vtkInformation* info) override;
168 int FillOutputPortInformation(int port, vtkInformation* info) override;
170
171private:
173 void operator=(const vtkVectorFieldTopology&) = delete;
174
178 int Validate();
179
187 int ImageDataPrepare(vtkDataSet* dataSetInput, vtkUnstructuredGrid* tridataset);
188
196 int UnstructuredGridPrepare(vtkDataSet* dataSetInput, vtkUnstructuredGrid* tridataset);
197
203 int RemoveBoundary(vtkSmartPointer<vtkUnstructuredGrid> tridataset);
204
212 int ComputeCriticalPoints2D(
214
222 int ComputeCriticalPoints3D(
224
235 static void InterpolateVector(
236 double x0, double x1, double x, const double v0[3], const double v1[3], double v[3]);
237
244 int ComputeBoundarySwitchPoints(
245 vtkPolyData* boundarySwitchPoints, vtkUnstructuredGrid* tridataset);
246
262 int ComputeSeparatricesBoundarySwitchPoints(vtkPolyData* boundarySwitchPoints,
263 vtkPolyData* separatrices, vtkDataSet* dataset, vtkPoints* interestPoints,
264 int integrationStepUnit, double dist, int maxNumSteps);
265
283 int ComputeSeparatricesBoundarySwitchLines(vtkPolyData* boundarySwitchLines,
284 vtkPolyData* separatrices, vtkDataSet* dataset, int integrationStepUnit, double dist,
285 int maxNumSteps, bool computeSurfaces, bool useIterativeSeeding);
286
305 int ComputeSeparatrices(vtkPolyData* criticalPoints, vtkPolyData* separatrices,
306 vtkPolyData* surfaces, vtkDataSet* dataset, vtkPoints* interestPoints, int integrationStepUnit,
307 double dist, double stepSize, int maxNumSteps, bool computeSurfaces, bool useIterativeSeeding);
308
325 int ComputeSurface(int numberOfSeparatingSurfaces, bool isBackward, double normal[3],
326 double zeroPos[3], vtkPolyData* streamSurfaces, vtkDataSet* dataset, int integrationStepUnit,
327 double dist, double stepSize, int maxNumSteps, bool useIterativeSeeding);
328
333 enum CriticalType2D
334 {
335 DEGENERATE_2D = -1,
336 SINK_2D = 0,
337 SADDLE_2D = 1,
338 SOURCE_2D = 2,
339 CENTER_2D = 3
340 };
341
346 enum CriticalTypeDetailed2D
347 {
348 // DEGENERATE2D = -1,
349 ATTRACTING_NODE_2D = 0,
350 ATTRACTING_FOCUS_2D = 1,
351 NODE_SADDLE_2D = 2,
352 REPELLING_NODE_2D = 3,
353 REPELLING_FOCUS_2D = 4,
354 CENTER_DETAILED_2D = 5
355 };
356
361 enum CriticalType3D
362 {
363 DEGENERATE_3D = -1,
364 SINK_3D = 0,
365 SADDLE_1_3D = 1,
366 SADDLE_2_3D = 2,
367 SOURCE_3D = 3,
368 CENTER_3D = 4
369 };
370
375 enum CriticalTypeDetailed3D
376 {
377 ATTRACTING_NODE_3D = 0,
378 ATTRACTING_FOCUS_3D = 1,
379 NODE_SADDLE_1_3D = 2,
380 FOCUS_SADDLE_1_3D = 3,
381 NODE_SADDLE_2_3D = 4,
382 FOCUS_SADDLE_2_3D = 5,
383 REPELLING_NODE_3D = 6,
384 REPELLING_FOCUS_3D = 7,
385 CENTER_DETAILED_3D = 8
386 };
387
395 static int Classify2D(int countComplex, int countPos, int countNeg);
396
405 static int ClassifyDetailed2D(int countComplex, int countPos, int countNeg);
406
415 static int Classify3D(int countComplex, int countPos, int countNeg);
416
426 static int ClassifyDetailed3D(int countComplex, int countPos, int countNeg);
427
431 int MaxNumSteps = 100;
432
436 double IntegrationStepSize = 1;
437
441 double SeparatrixDistance = 1;
442
446 bool UseIterativeSeeding = false;
447
451 bool ComputeSurfaces = false;
452
456 const char* NameOfVectorArray;
457
462 bool ExcludeBoundary = false;
463
467 int Dimension = 2;
468
476 int IntegrationStepUnit = vtkStreamTracer::CELL_LENGTH_UNIT;
477
484 bool UseBoundarySwitchPoints = false;
485
491
500 double VectorAngleThreshold = 1;
501
509 double OffsetAwayFromBoundary = 1e-3;
510
514 double EpsilonCriticalPoint = 1e-10;
515
516 vtkNew<vtkStreamSurface> StreamSurface;
517};
518#endif
abstract class to specify dataset behavior
Definition vtkDataSet.h:57
A general filter for gradient estimation.
topologically and geometrically regular array of data
a simple class to control print indentation
Definition vtkIndent.h:34
Store zero or more vtkInformation instances.
Store vtkAlgorithm input/output information.
Allocate and hold a VTK object.
Definition vtkNew.h:56
represent and manipulate 3D points
Definition vtkPoints.h:34
Superclass for algorithms that produce only polydata as output.
concrete dataset represents vertices, lines, polygons, and triangle strips
Definition vtkPolyData.h:85
Hold a reference to a vtkObjectBase instance.
Advect a stream surface in a vector field.
@ INTERPOLATOR_WITH_DATASET_POINT_LOCATOR
dataset represents arbitrary combinations of all possible cell types
Extract the topological skeleton as output datasets.
void SetInterpolatorTypeToCellLocator()
Set the velocity field interpolator type to the one involving a cell locator.
int FillInputPortInformation(int port, vtkInformation *info) override
Fill the input port information objects for this algorithm.
int RequestData(vtkInformation *, vtkInformationVector **, vtkInformationVector *) override
This is called by the superclass.
static vtkVectorFieldTopology * New()
void SetInterpolatorTypeToDataSetPointLocator()
Set the velocity field interpolator type to the one involving a dataset point locator.
~vtkVectorFieldTopology() override
void PrintSelf(ostream &os, vtkIndent indent) override
Methods invoked by print to print information about the object including superclasses.
void SetInterpolatorType(int interpType)
Set the type of the velocity field interpolator to determine whether INTERPOLATOR_WITH_DATASET_POINT_...
int FillOutputPortInformation(int port, vtkInformation *info) override
Fill the output port information objects for this algorithm.