VTK
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
95class vtkDataArray;
96class vtkDoubleArray;
97class vtkExecutive;
98class vtkGenericCell;
99class vtkIdList;
100class vtkIntArray;
102
103class VTKFILTERSFLOWPATHS_EXPORT vtkStreamTracer : public vtkPolyDataAlgorithm
104{
105public:
107 void PrintSelf(ostream& os, vtkIndent indent);
108
117
119
124 vtkSetVector3Macro(StartPosition, double);
125 vtkGetVector3Macro(StartPosition, double);
127
129
138
144
145 // The previously-supported TIME_UNIT is excluded in this current
146 // enumeration definition because the underlying step size is ALWAYS in
147 // arc length unit (LENGTH_UNIT) while the 'real' time interval (virtual
148 // for steady flows) that a particle actually takes to trave in a single
149 // step is obtained by dividing the arc length by the LOCAL speed. The
150 // overall elapsed time (i.e., the life span) of the particle is the sum
151 // of those individual step-wise time intervals. The arc-length-to-time
152 // conversion only occurs for vorticity computation and for generating a
153 // point data array named 'IntegrationTime'.
154 enum Units
155 {
156 LENGTH_UNIT = 1,
157 CELL_LENGTH_UNIT = 2
158 };
159
161 {
166 UNKNOWN
167 };
168
170 {
174 OUT_OF_LENGTH = 4,
175 OUT_OF_STEPS = 5,
176 STAGNATION = 6
177 };
178
180
191 vtkGetObjectMacro ( Integrator, vtkInitialValueProblemSolver );
195 {this->SetIntegratorType(RUNGE_KUTTA2);};
197 {this->SetIntegratorType(RUNGE_KUTTA4);};
199 {this->SetIntegratorType(RUNGE_KUTTA45);};
201
207
213
215
218 vtkSetMacro(MaximumPropagation, double);
219 vtkGetMacro(MaximumPropagation, double);
221
228 void SetIntegrationStepUnit( int unit );
229 int GetIntegrationStepUnit() { return this->IntegrationStepUnit; }
230
232
239 vtkSetMacro(InitialIntegrationStep, double);
240 vtkGetMacro(InitialIntegrationStep, double);
242
244
250 vtkSetMacro(MinimumIntegrationStep, double);
251 vtkGetMacro(MinimumIntegrationStep, double);
253
255
261 vtkSetMacro(MaximumIntegrationStep, double);
262 vtkGetMacro(MaximumIntegrationStep, double);
264
266
269 vtkSetMacro(MaximumError, double);
270 vtkGetMacro(MaximumError, double);
272
274
277 vtkSetMacro(MaximumNumberOfSteps, vtkIdType);
278 vtkGetMacro(MaximumNumberOfSteps, vtkIdType);
280
282
285 vtkSetMacro(TerminalSpeed, double);
286 vtkGetMacro(TerminalSpeed, double);
288
290
293 vtkGetMacro(SurfaceStreamlines, bool);
294 vtkSetMacro(SurfaceStreamlines, bool);
295 vtkBooleanMacro(SurfaceStreamlines, bool);
297
298 enum
299 {
302 BOTH
303 };
304
305 enum
306 {
308 INTERPOLATOR_WITH_CELL_LOCATOR
309 };
310
312
316 vtkSetClampMacro(IntegrationDirection, int, FORWARD, BOTH);
317 vtkGetMacro(IntegrationDirection, int);
319 {this->SetIntegrationDirection(FORWARD);};
321 {this->SetIntegrationDirection(BACKWARD);};
323 {this->SetIntegrationDirection(BOTH);};
325
327
332 vtkSetMacro(ComputeVorticity, bool);
333 vtkGetMacro(ComputeVorticity, bool);
335
337
341 vtkSetMacro(RotationScale, double);
342 vtkGetMacro(RotationScale, double);
344
350
360 void SetInterpolatorType( int interpType );
361
362protected:
363
366
367 // Create a default executive.
369
370 // hide the superclass' AddInput() from the user and the compiler
372 { vtkErrorMacro( << "AddInput() must be called with a vtkDataSet not a vtkDataObject."); };
373
376
377 void CalculateVorticity( vtkGenericCell* cell, double pcoords[3],
378 vtkDoubleArray* cellVectors, double vorticity[3] );
379 void Integrate(vtkPointData *inputData,
380 vtkPolyData* output,
381 vtkDataArray* seedSource,
382 vtkIdList* seedIds,
383 vtkIntArray* integrationDirections,
384 double lastPoint[3],
386 int maxCellSize,
387 int vecType,
388 const char *vecFieldName,
389 double& propagation,
390 vtkIdType& numSteps,
391 double& integrationTime);
392 double SimpleIntegrate(double seed[3],
393 double lastPoint[3],
394 double stepSize,
397 int* maxCellSize);
398 void GenerateNormals(vtkPolyData* output, double* firstNormal, const char *vecName);
399
401
402 // starting from global x-y-z position
403 double StartPosition[3];
404
405 static const double EPSILON;
407
409
411 {
412 double Interval;
413 int Unit;
414 };
415
420
421 void ConvertIntervals( double& step, double& minStep, double& maxStep,
422 int direction, double cellLength );
423 static double ConvertToLength( double interval, int unit, double cellLength );
424 static double ConvertToLength( IntervalInformation& interval, double cellLength );
425
427 vtkInformation* outInfo);
429 vtkIdList*& seedIds,
430 vtkIntArray*& integrationDirections,
432
435
436 // Prototype showing the integrator type to be set by the user.
438
441
444
445 // Compute streamlines only on surface.
447
449
451 bool HasMatchingPointAttributes; //does the point data in the multiblocks have the same attributes?
452
453 friend class PStreamTracerUtils;
454
455private:
456 vtkStreamTracer(const vtkStreamTracer&) VTK_DELETE_FUNCTION;
457 void operator=(const vtkStreamTracer&) VTK_DELETE_FUNCTION;
458};
459
460
461#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:55
general representation of visualization data
Definition: vtkDataObject.h:65
abstract class to specify dataset behavior
Definition: vtkDataSet.h:63
dynamic, self-adjusting array of double
Superclass for all pipeline executives in VTK.
Definition: vtkExecutive.h:50
provides thread-safe access to cells
list of point or cell ids
Definition: vtkIdList.h:37
a simple class to control print indentation
Definition: vtkIndent.h:40
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:46
represent and manipulate point attribute data
Definition: vtkPointData.h:38
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()
int SetupOutput(vtkInformation *inInfo, vtkInformation *outInfo)
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)
void SetSourceData(vtkDataSet *source)
Specify the source object used to generate starting points (seeds).
vtkDataSet * GetSource()
double InitialIntegrationStep
vtkAbstractInterpolatedVelocityField * InterpolatorPrototype
void SetInterpolatorTypeToCellLocator()
Set the velocity field interpolator type to the one involving a cell locator.
void CalculateVorticity(vtkGenericCell *cell, double pcoords[3], vtkDoubleArray *cellVectors, double vorticity[3])
double MinimumIntegrationStep
void SetIntegratorTypeToRungeKutta4()
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).
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
@ INTERPOLATOR_WITH_DATASET_POINT_LOCATOR
vtkIdType MaximumNumberOfSteps
void SetIntegrationDirectionToForward()
static vtkStreamTracer * New()
Construct object to start from position (0,0,0), with forward integration, terminal speed 1....
bool HasMatchingPointAttributes
vtkCompositeDataSet * InputData
virtual int FillInputPortInformation(int, vtkInformation *)
Fill the input port information objects for this algorithm.
virtual vtkExecutive * CreateDefaultExecutive()
Create a default executive.
void SetInterpolatorType(int interpType)
Set the type of the velocity field interpolator to determine whether vtkInterpolatedVelocityField (IN...
double MaximumIntegrationStep
bool GenerateNormalsInIntegrate
void SetIntegrationDirectionToBackward()
void SetInterpolatorTypeToDataSetPointLocator()
Set the velocity field interpolator type to the one involving a dataset point locator.
int GetIntegratorType()
void PrintSelf(ostream &os, vtkIndent indent)
Methods invoked by print to print information about the object including superclasses.
virtual int RequestData(vtkInformation *, vtkInformationVector **, vtkInformationVector *)
This is called by the superclass.
void InitializeSeeds(vtkDataArray *&seeds, vtkIdList *&seedIds, vtkIntArray *&integrationDirections, vtkDataSet *source)
void SetIntegratorTypeToRungeKutta2()
void SetIntegrationDirectionToBoth()
double SimpleIntegrate(double seed[3], double lastPoint[3], double stepSize, vtkAbstractInterpolatedVelocityField *func)
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)
int GetIntegrationStepUnit()
@ direction
Definition: vtkX3D.h:260
@ type
Definition: vtkX3D.h:516
boost::graph_traits< vtkGraph * >::vertex_descriptor source(boost::graph_traits< vtkGraph * >::edge_descriptor e, vtkGraph *)
vtkSetMacro(IgnoreDriverBugs, bool)
Updates the extensions string.
vtkBooleanMacro(IgnoreDriverBugs, bool)
Updates the extensions string.
int vtkIdType
Definition: vtkType.h:287