VTK  9.1.0
vtkStreamTracer.h
Go to the documentation of this file.
1/*=========================================================================
2
3 Program: Visualization Toolkit
4 Module: vtkStreamTracer.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=========================================================================*/
86#ifndef vtkStreamTracer_h
87#define vtkStreamTracer_h
88
89#include "vtkFiltersFlowPathsModule.h" // For export macro
91
92#include "vtkInitialValueProblemSolver.h" // Needed for constants
93
96class vtkDataArray;
98class vtkDoubleArray;
99class vtkExecutive;
100class vtkGenericCell;
101class vtkIdList;
102class vtkIntArray;
103class vtkPoints;
104
105#include <vector> // for std::vector
106
107class VTKFILTERSFLOWPATHS_EXPORT vtkStreamTracer : public vtkPolyDataAlgorithm
108{
109public:
111 void PrintSelf(ostream& os, vtkIndent indent) override;
112
121
123
128 vtkSetVector3Macro(StartPosition, double);
129 vtkGetVector3Macro(StartPosition, double);
131
133
142
148
149 // The previously-supported TIME_UNIT is excluded in this current
150 // enumeration definition because the underlying step size is ALWAYS in
151 // arc length unit (LENGTH_UNIT) while the 'real' time interval (virtual
152 // for steady flows) that a particle actually takes to trave in a single
153 // step is obtained by dividing the arc length by the LOCAL speed. The
154 // overall elapsed time (i.e., the life span) of the particle is the sum
155 // of those individual step-wise time intervals. The arc-length-to-time
156 // conversion only occurs for vorticity computation and for generating a
157 // point data array named 'IntegrationTime'.
158 enum Units
159 {
160 LENGTH_UNIT = 1,
161 CELL_LENGTH_UNIT = 2
162 };
163
165 {
170 UNKNOWN
171 };
172
174 {
178 OUT_OF_LENGTH = 4,
179 OUT_OF_STEPS = 5,
180 STAGNATION = 6,
181 FIXED_REASONS_FOR_TERMINATION_COUNT
182 };
183
185
196 vtkGetObjectMacro(Integrator, vtkInitialValueProblemSolver);
199 void SetIntegratorTypeToRungeKutta2() { this->SetIntegratorType(RUNGE_KUTTA2); }
200 void SetIntegratorTypeToRungeKutta4() { this->SetIntegratorType(RUNGE_KUTTA4); }
201 void SetIntegratorTypeToRungeKutta45() { this->SetIntegratorType(RUNGE_KUTTA45); }
203
209
215
217
220 vtkSetMacro(MaximumPropagation, double);
221 vtkGetMacro(MaximumPropagation, double);
223
231 int GetIntegrationStepUnit() { return this->IntegrationStepUnit; }
232
234
241 vtkSetMacro(InitialIntegrationStep, double);
242 vtkGetMacro(InitialIntegrationStep, double);
244
246
252 vtkSetMacro(MinimumIntegrationStep, double);
253 vtkGetMacro(MinimumIntegrationStep, double);
255
257
263 vtkSetMacro(MaximumIntegrationStep, double);
264 vtkGetMacro(MaximumIntegrationStep, double);
266
268
271 vtkSetMacro(MaximumError, double);
272 vtkGetMacro(MaximumError, double);
274
276
279 vtkSetMacro(MaximumNumberOfSteps, vtkIdType);
280 vtkGetMacro(MaximumNumberOfSteps, vtkIdType);
282
284
287 vtkSetMacro(TerminalSpeed, double);
288 vtkGetMacro(TerminalSpeed, double);
290
292
295 vtkGetMacro(SurfaceStreamlines, bool);
296 vtkSetMacro(SurfaceStreamlines, bool);
297 vtkBooleanMacro(SurfaceStreamlines, bool);
299
300 enum
301 {
304 BOTH
305 };
306
307 enum
308 {
310 INTERPOLATOR_WITH_CELL_LOCATOR
311 };
312
314
318 vtkSetClampMacro(IntegrationDirection, int, FORWARD, BOTH);
319 vtkGetMacro(IntegrationDirection, int);
320 void SetIntegrationDirectionToForward() { this->SetIntegrationDirection(FORWARD); }
321 void SetIntegrationDirectionToBackward() { this->SetIntegrationDirection(BACKWARD); }
322 void SetIntegrationDirectionToBoth() { this->SetIntegrationDirection(BOTH); }
324
326
331 vtkSetMacro(ComputeVorticity, bool);
332 vtkGetMacro(ComputeVorticity, bool);
334
336
340 vtkSetMacro(RotationScale, double);
341 vtkGetMacro(RotationScale, double);
343
345
353 vtkSetMacro(UseLocalSeedSource, bool);
354 vtkGetMacro(UseLocalSeedSource, bool);
355 vtkBooleanMacro(UseLocalSeedSource, bool);
357
363
373 void SetInterpolatorType(int interpType);
374
384 typedef bool (*CustomTerminationCallbackType)(
385 void* clientdata, vtkPoints* points, vtkDataArray* velocity, int integrationDirection);
395 CustomTerminationCallbackType callback, void* clientdata, int reasonForTermination);
396
397protected:
400
401 // Create a default executive.
403
404 // hide the superclass' AddInput() from the user and the compiler
406 {
407 vtkErrorMacro(<< "AddInput() must be called with a vtkDataSet not a vtkDataObject.");
408 }
409
412
414 vtkGenericCell* cell, double pcoords[3], vtkDoubleArray* cellVectors, double vorticity[3]);
415 void Integrate(vtkPointData* inputData, vtkPolyData* output, vtkDataArray* seedSource,
416 vtkIdList* seedIds, vtkIntArray* integrationDirections, double lastPoint[3],
417 vtkAbstractInterpolatedVelocityField* func, int maxCellSize, int vecType,
418 const char* vecFieldName, double& propagation, vtkIdType& numSteps, double& integrationTime);
419 double SimpleIntegrate(double seed[3], double lastPoint[3], double stepSize,
422 void GenerateNormals(vtkPolyData* output, double* firstNormal, const char* vecName);
423
425
426 // starting from global x-y-z position
427 double StartPosition[3];
428
429 static const double EPSILON;
431
433
435 {
436 double Interval;
437 int Unit;
438 };
439
444
446 double& step, double& minStep, double& maxStep, int direction, double cellLength);
447 static double ConvertToLength(double interval, int unit, double cellLength);
448 static double ConvertToLength(IntervalInformation& interval, double cellLength);
449
451 void InitializeSeeds(vtkDataArray*& seeds, vtkIdList*& seedIds,
452 vtkIntArray*& integrationDirections, vtkDataSet* source);
453
456
457 // Prototype showing the integrator type to be set by the user.
459
462
465
466 // Compute streamlines only on surface.
468
469 // Only relevant for the parallel version of this filter (see vtkPStreamTracer)
470 bool UseLocalSeedSource = true;
471
473
475 bool
476 HasMatchingPointAttributes; // does the point data in the multiblocks have the same attributes?
477 std::vector<CustomTerminationCallbackType> CustomTerminationCallback;
478 std::vector<void*> CustomTerminationClientData;
480
481 friend class PStreamTracerUtils;
482
483private:
484 vtkStreamTracer(const vtkStreamTracer&) = delete;
485 void operator=(const vtkStreamTracer&) = delete;
486};
487
488#endif
An abstract class for obtaining the interpolated velocity values at a point.
Proxy object to connect input/output ports.
abstract superclass for composite (multi-block or AMR) datasets
abstract superclass for arrays of numeric data
Definition: vtkDataArray.h:50
general representation of visualization data
Definition: vtkDataObject.h:60
represent and manipulate attribute data in a dataset
abstract class to specify dataset behavior
Definition: vtkDataSet.h:57
dynamic, self-adjusting array of double
Superclass for all pipeline executives in VTK.
Definition: vtkExecutive.h:47
provides thread-safe access to cells
list of point or cell ids
Definition: vtkIdList.h:31
a simple class to control print indentation
Definition: vtkIndent.h:34
Store zero or more vtkInformation instances.
Store vtkAlgorithm input/output information.
Integrate a set of ordinary differential equations (initial value problem) in time.
dynamic, self-adjusting array of int
Definition: vtkIntArray.h:40
represent and manipulate point attribute data
Definition: vtkPointData.h:33
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:86
Streamline generator.
void SetIntegratorTypeToRungeKutta45()
Set/get the integrator type to be used for streamline generation.
int FillInputPortInformation(int, vtkInformation *) override
Fill the input port information objects for this algorithm.
int SetupOutput(vtkInformation *inInfo, vtkInformation *outInfo)
std::vector< void * > CustomTerminationClientData
void Integrate(vtkPointData *inputData, vtkPolyData *output, vtkDataArray *seedSource, vtkIdList *seedIds, vtkIntArray *integrationDirections, double lastPoint[3], vtkAbstractInterpolatedVelocityField *func, int maxCellSize, int vecType, const char *vecFieldName, double &propagation, vtkIdType &numSteps, double &integrationTime)
@ INTERPOLATOR_WITH_DATASET_POINT_LOCATOR
void SetSourceData(vtkDataSet *source)
Specify the source object used to generate starting points (seeds).
vtkDataSet * GetSource()
Specify the source object used to generate starting points (seeds).
double InitialIntegrationStep
vtkAbstractInterpolatedVelocityField * InterpolatorPrototype
void SetInterpolatorTypeToCellLocator()
Set the velocity field interpolator type to the one involving a cell locator.
void PrintSelf(ostream &os, vtkIndent indent) override
Methods invoked by print to print information about the object including superclasses.
void CalculateVorticity(vtkGenericCell *cell, double pcoords[3], vtkDoubleArray *cellVectors, double vorticity[3])
double MinimumIntegrationStep
void SetIntegratorTypeToRungeKutta4()
Set/get the integrator type to be used for streamline generation.
static double ConvertToLength(double interval, int unit, double cellLength)
void SetIntegrator(vtkInitialValueProblemSolver *)
Set/get the integrator type to be used for streamline generation.
void SetSourceConnection(vtkAlgorithmOutput *algOutput)
Specify the source object used to generate starting points (seeds).
std::vector< int > CustomReasonForTermination
static double ConvertToLength(IntervalInformation &interval, double cellLength)
int CheckInputs(vtkAbstractInterpolatedVelocityField *&func, int *maxCellSize)
void ConvertIntervals(double &step, double &minStep, double &maxStep, int direction, double cellLength)
void GenerateNormals(vtkPolyData *output, double *firstNormal, const char *vecName)
static const double EPSILON
vtkIdType MaximumNumberOfSteps
void SetIntegrationDirectionToForward()
Specify whether the streamline is integrated in the upstream or downstream direction.
std::vector< CustomTerminationCallbackType > CustomTerminationCallback
static vtkStreamTracer * New()
Construct object to start from position (0,0,0), with forward integration, terminal speed 1....
bool HasMatchingPointAttributes
vtkCompositeDataSet * InputData
void SetInterpolatorType(int interpType)
Set the type of the velocity field interpolator to determine whether vtkInterpolatedVelocityField (IN...
double MaximumIntegrationStep
int RequestData(vtkInformation *, vtkInformationVector **, vtkInformationVector *) override
This is called by the superclass.
bool GenerateNormalsInIntegrate
vtkExecutive * CreateDefaultExecutive() override
Create a default executive.
void SetIntegrationDirectionToBackward()
Specify whether the streamline is integrated in the upstream or downstream direction.
void SetInterpolatorTypeToDataSetPointLocator()
Set the velocity field interpolator type to the one involving a dataset point locator.
int GetIntegratorType()
Set/get the integrator type to be used for streamline generation.
void AddCustomTerminationCallback(CustomTerminationCallbackType callback, void *clientdata, int reasonForTermination)
Adds a custom termination callback.
void InitializeSeeds(vtkDataArray *&seeds, vtkIdList *&seedIds, vtkIntArray *&integrationDirections, vtkDataSet *source)
void SetIntegratorTypeToRungeKutta2()
Set/get the integrator type to be used for streamline generation.
void SetIntegrationDirectionToBoth()
Specify whether the streamline is integrated in the upstream or downstream direction.
double SimpleIntegrate(double seed[3], double lastPoint[3], double stepSize, vtkAbstractInterpolatedVelocityField *func)
~vtkStreamTracer() override
void AddInput(vtkDataObject *)
vtkInitialValueProblemSolver * Integrator
void SetInterpolatorPrototype(vtkAbstractInterpolatedVelocityField *ivf)
The object used to interpolate the velocity field during integration is of the same class as this pro...
void SetIntegrationStepUnit(int unit)
Specify a uniform integration step unit for MinimumIntegrationStep, InitialIntegrationStep,...
void SetIntegratorType(int type)
Set/get the integrator type to be used for streamline generation.
int GetIntegrationStepUnit()
@ points
Definition: vtkX3D.h:452
@ direction
Definition: vtkX3D.h:266
@ type
Definition: vtkX3D.h:522
boost::graph_traits< vtkGraph * >::vertex_descriptor source(boost::graph_traits< vtkGraph * >::edge_descriptor e, vtkGraph *)
int vtkIdType
Definition: vtkType.h:332