VTK
vtkParticleTracerBase.h
Go to the documentation of this file.
1/*=========================================================================
2
3 Program: Visualization Toolkit
4 Module: vtkParticleTracerBase.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=========================================================================*/
28#ifndef vtkParticleTracerBase_h
29#define vtkParticleTracerBase_h
30
31#include "vtkFiltersFlowPathsModule.h" // For export macro
32#include "vtkSmartPointer.h" // For protected ivars.
34
35#include <vector> // STL Header
36#include <list> // STL Header
37
40class vtkCellArray;
41class vtkCharArray;
43class vtkDataArray;
44class vtkDataSet;
45class vtkDoubleArray;
46class vtkFloatArray;
47class vtkGenericCell;
49class vtkIntArray;
52class vtkPointData;
53class vtkPoints;
54class vtkPolyData;
56
58{
59 typedef struct { double x[4]; } Position;
60 typedef struct {
61 // These are used during iteration
63 int CachedDataSetId[2];
64 vtkIdType CachedCellId[2];
66 // These are computed scalars we might display
68 int TimeStepAge; // amount of time steps the particle has advanced
70 int InjectedStepId; // time step the particle was injected
73 // These are useful to track for debugging etc
75 float age;
76 // these are needed across time steps to compute vorticity
77 float rotation;
79 float time;
80 float speed;
81 // once the partice is added, PointId is valid and is the tuple location
82 // in ProtoPD.
84 // if PointId is negative then in parallel this particle was just
85 // received and we need to get the tuple value from vtkPParticleTracerBase::Tail.
88
89 typedef std::vector<ParticleInformation> ParticleVector;
90 typedef ParticleVector::iterator ParticleIterator;
91 typedef std::list<ParticleInformation> ParticleDataList;
92 typedef ParticleDataList::iterator ParticleListIterator;
93};
94
95class VTKFILTERSFLOWPATHS_EXPORT vtkParticleTracerBase : public vtkPolyDataAlgorithm
96{
97public:
99 {
104 UNKNOWN
105 };
106
108 void PrintSelf(ostream& os, vtkIndent indent);
109 void PrintParticleHistories();
110
112
117 vtkGetMacro(ComputeVorticity, bool);
118 void SetComputeVorticity(bool);
120
122
125 vtkGetMacro(TerminalSpeed, double);
126 void SetTerminalSpeed(double);
128
130
134 vtkGetMacro(RotationScale, double);
135 void SetRotationScale(double);
137
139
143 vtkSetMacro(IgnorePipelineTime, int);
144 vtkGetMacro(IgnorePipelineTime, int);
145 vtkBooleanMacro(IgnorePipelineTime, int);
147
149
158 vtkGetMacro(ForceReinjectionEveryNSteps,int);
159 void SetForceReinjectionEveryNSteps(int);
161
163
169 void SetTerminationTime(double t);
170 vtkGetMacro(TerminationTime,double);
172
173 void SetIntegrator(vtkInitialValueProblemSolver *);
174 vtkGetObjectMacro ( Integrator, vtkInitialValueProblemSolver );
175
176 void SetIntegratorType(int type);
177 int GetIntegratorType();
178
180
184 vtkGetMacro(StartTime, double);
185 void SetStartTime(double t);
187
189
198 vtkSetMacro(StaticSeeds,int);
199 vtkGetMacro(StaticSeeds,int);
201
203
212 vtkSetMacro(StaticMesh,int);
213 vtkGetMacro(StaticMesh,int);
215
217
223 virtual void SetParticleWriter(vtkAbstractParticleWriter *pw);
224 vtkGetObjectMacro(ParticleWriter, vtkAbstractParticleWriter);
226
228
232 vtkSetStringMacro(ParticleFileName);
233 vtkGetStringMacro(ParticleFileName);
235
237
241 vtkSetMacro(EnableParticleWriting,int);
242 vtkGetMacro(EnableParticleWriting,int);
243 vtkBooleanMacro(EnableParticleWriting,int);
245
247
252 vtkSetMacro(DisableResetCache,int);
253 vtkGetMacro(DisableResetCache,int);
254 vtkBooleanMacro(DisableResetCache,int);
256
258
261 void AddSourceConnection(vtkAlgorithmOutput* input);
262 void RemoveAllSources();
264
265 protected:
266 vtkSmartPointer<vtkPolyData> Output; //managed by child classes
268
273 vtkIdType UniqueIdCounter;// global Id counter used to give particles a stamp
275 vtkSmartPointer<vtkPointData> ParticlePointData; //the current particle point data consistent
276 //with particle history
277 //Everything related to time
278 int IgnorePipelineTime; //whether to use the pipeline time for termination
279 int DisableResetCache; //whether to enable ResetCache() method
281
284
285 //
286 // Make sure the pipeline knows what type we expect as input
287 //
288 virtual int FillInputPortInformation(int port, vtkInformation* info);
289
290 //
291 // The usual suspects
292 //
293 virtual int ProcessRequest(vtkInformation* request,
294 vtkInformationVector** inputVector,
295 vtkInformationVector* outputVector);
296
297 //
298 // Store any information we need in the output and fetch what we can
299 // from the input
300 //
301 virtual int RequestInformation(vtkInformation* request,
302 vtkInformationVector** inputVector,
303 vtkInformationVector* outputVector);
304
305 //
306 // Compute input time steps given the output step
307 //
308 virtual int RequestUpdateExtent(vtkInformation* request,
309 vtkInformationVector** inputVector,
310 vtkInformationVector* outputVector);
311
312 //
313 // what the pipeline calls for each time step
314 //
315 virtual int RequestData(vtkInformation* request,
316 vtkInformationVector** inputVector,
317 vtkInformationVector* outputVector);
318
319 //
320 // these routines are internally called to actually generate the output
321 //
322 virtual int ProcessInput(vtkInformationVector** inputVector);
323
324 // This is the main part of the algorithm:
325 // * move all the particles one step
326 // * Reinject particles (by adding them to this->ParticleHistories)
327 // either at the beginning or at the end of each step (modulo this->ForceReinjectionEveryNSteps)
328 // * Output a polydata representing the moved particles
329 // Note that if the starting and the ending time coincide, the polydata is still valid.
330 virtual vtkPolyData* Execute(vtkInformationVector** inputVector);
331
332 // the RequestData will call these methods in turn
333 virtual void Initialize(){} //the first iteration
334 virtual int OutputParticles(vtkPolyData* poly)=0; //every iteration
335 virtual void Finalize(){} //the last iteration
336
341 virtual std::vector<vtkDataSet*> GetSeedSources(vtkInformationVector* inputVector, int timeStep);
342
343 //
344 // Initialization of input (vector-field) geometry
345 //
348
356 int &count);
357
359 vtkParticleTracerBaseNamespace::ParticleVector &candidates, std::vector<int> &passed);
360
368 double time, vtkDataSet *source, int sourceID, int ptId,
370 int &localAssignedCount);
371
376 virtual void AssignUniqueIds(
378
385
391 virtual bool UpdateParticleListFromOtherProcesses(){return false;}
392
398 double currenttime, double terminationtime,
399 vtkInitialValueProblemSolver* integrator);
400
401 // if the particle is added to send list, then returns value is 1,
402 // if it is kept on this process after a retry return value is 0
406 {
407 return true;
408 }
409
417 double pos[4], double p2[4], double intersection[4],
418 vtkGenericCell *cell);
419
420 //
421 // Scalar arrays that are generated as each particle is updated
422 //
424
434
435 // utility function we use to test if a point is inside any of our local datasets
436 bool InsideBounds(double point[]);
437
438 void CalculateVorticity( vtkGenericCell* cell, double pcoords[3],
439 vtkDoubleArray* cellVectors, double vorticity[3] );
440
441 //------------------------------------------------------
442
443
444 double GetCacheDataTime(int i);
446
447 virtual void ResetCache();
449
451
456 virtual bool IsPointDataValid(vtkDataObject* input);
457 bool IsPointDataValid(vtkCompositeDataSet* input, std::vector<std::string>& arrayNames);
458 void GetPointDataArrayNames(vtkDataSet* input, std::vector<std::string>& names);
460
461 vtkGetMacro(ReinjectionCounter, int);
462 vtkGetMacro(CurrentTimeValue, double);
463
468 virtual void InitializeExtraPointDataArrays(vtkPointData* vtkNotUsed(outputPD)) {}
469
471
473
478 virtual void AddRestartSeeds(vtkInformationVector** /*inputVector*/) {}
479
480 private:
484 void SetInterpolatorPrototype(vtkAbstractInterpolatedVelocityField*) {}
485
495 bool RetryWithPush(
496 vtkParticleTracerBaseNamespace::ParticleInformation &info, double* point1,double delT, int subSteps);
497
498 bool SetTerminationTimeNoModify(double t);
499
500 //Parameters of tracing
502 double IntegrationStep;
503 double MaximumError;
504 bool ComputeVorticity;
505 double RotationScale;
506 double TerminalSpeed;
507
508 // A counter to keep track of how many times we reinjected
509 int ReinjectionCounter;
510
511 // Important for Caching of Cells/Ids/Weights etc
512 int AllFixedGeometry;
513 int StaticMesh;
514 int StaticSeeds;
515
516 std::vector<double> InputTimeValues;
517 double StartTime;
518 double TerminationTime;
519 double CurrentTimeValue;
520
521 int StartTimeStep; //InputTimeValues[StartTimeStep] <= StartTime <= InputTimeValues[StartTimeStep+1]
522 int CurrentTimeStep;
523 int TerminationTimeStep; //computed from start time
524 bool FirstIteration;
525
526 //Innjection parameters
527 int ForceReinjectionEveryNSteps;
528 vtkTimeStamp ParticleInjectionTime;
529 bool HasCache;
530
531 // Particle writing to disk
532 vtkAbstractParticleWriter *ParticleWriter;
533 char *ParticleFileName;
534 int EnableParticleWriting;
535
536
537 // The main lists which are held during operation- between time step updates
539
540 // The velocity interpolator
542 vtkAbstractInterpolatedVelocityField * InterpolatorPrototype;
543
544 // Data for time step CurrentTimeStep-1 and CurrentTimeStep
546
547 // Cache bounds info for each dataset we will use repeatedly
548 typedef struct {
549 double b[6];
550 } bounds;
551 std::vector<bounds> CachedBounds[2];
552
553 // temporary variables used by Exeucte(), for convenience only
554
555 vtkSmartPointer<vtkPoints> OutputCoordinates;
558 vtkSmartPointer<vtkCharArray> ParticleSourceIds;
559 vtkSmartPointer<vtkIntArray> InjectedPointIds;
560 vtkSmartPointer<vtkIntArray> InjectedStepIds;
562 vtkSmartPointer<vtkFloatArray> ParticleVorticity;
563 vtkSmartPointer<vtkFloatArray> ParticleRotation;
564 vtkSmartPointer<vtkFloatArray> ParticleAngularVel;
566 vtkSmartPointer<vtkPointData> OutputPointData;
567 vtkSmartPointer<vtkDataSet> DataReferenceT[2];
568 vtkSmartPointer<vtkCellArray> ParticleCells;
569
570 vtkParticleTracerBase(const vtkParticleTracerBase&) VTK_DELETE_FUNCTION;
571 void operator=(const vtkParticleTracerBase&) VTK_DELETE_FUNCTION;
572 vtkTimeStamp ExecuteTime;
573
574 unsigned int NumberOfParticles();
575
578
579 static const double Epsilon;
580};
581
582#endif
An abstract class for obtaining the interpolated velocity values at a point.
abstract class to write particle data to file
Proxy object to connect input/output ports.
unsigned long ErrorCode
Definition: vtkAlgorithm.h:917
object to represent cell connectivity
Definition: vtkCellArray.h:51
dynamic, self-adjusting array of char
Definition: vtkCharArray.h:39
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
dynamic, self-adjusting array of float
Definition: vtkFloatArray.h:42
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.
Integrate a set of ordinary differential equations (initial value problem) in time.
dynamic, self-adjusting array of int
Definition: vtkIntArray.h:46
Composite dataset that organizes datasets into blocks.
Multiprocessing communication superclass.
A particle tracer for vector fields.
vtkFloatArray * GetParticleVorticity(vtkPointData *)
void IntegrateParticle(vtkParticleTracerBaseNamespace::ParticleListIterator &it, double currenttime, double terminationtime, vtkInitialValueProblemSolver *integrator)
particle between the two times supplied.
vtkIntArray * GetInjectedPointIds(vtkPointData *)
vtkFloatArray * GetParticleAge(vtkPointData *)
void TestParticles(vtkParticleTracerBaseNamespace::ParticleVector &candidates, std::vector< int > &passed)
virtual void AddRestartSeeds(vtkInformationVector **)
For restarts of particle paths, we add in the ability to add in particles from a previous computation...
bool IsPointDataValid(vtkCompositeDataSet *input, std::vector< std::string > &arrayNames)
vtkIntArray * GetInjectedStepIds(vtkPointData *)
virtual bool UpdateParticleListFromOtherProcesses()
this is used during classification of seed points and also between iterations of the main loop as par...
virtual std::vector< vtkDataSet * > GetSeedSources(vtkInformationVector *inputVector, int timeStep)
Method to get the data set seed sources.
vtkIntArray * GetErrorCodeArr(vtkPointData *)
void UpdateParticleList(vtkParticleTracerBaseNamespace::ParticleVector &candidates)
and sending between processors, into a list, which is used as the master list on this processor
virtual void AppendToExtraPointDataArrays(vtkParticleTracerBaseNamespace::ParticleInformation &)
vtkFloatArray * GetParticleAngularVel(vtkPointData *)
double GetCacheDataTime(int i)
void TestParticles(vtkParticleTracerBaseNamespace::ParticleVector &candidates, vtkParticleTracerBaseNamespace::ParticleVector &passed, int &count)
inside our data.
void CalculateVorticity(vtkGenericCell *cell, double pcoords[3], vtkDoubleArray *cellVectors, double vorticity[3])
void CreateProtoPD(vtkDataObject *input)
virtual bool SendParticleToAnotherProcess(vtkParticleTracerBaseNamespace::ParticleInformation &, vtkParticleTracerBaseNamespace::ParticleInformation &, vtkPointData *)
vtkIntArray * GetParticleIds(vtkPointData *)
int UpdateDataCache(vtkDataObject *td)
void AddParticle(vtkParticleTracerBaseNamespace::ParticleInformation &info, double *velocity)
vtkFloatArray * GetParticleRotation(vtkPointData *)
void GetPointDataArrayNames(vtkDataSet *input, std::vector< std::string > &names)
virtual void AssignSeedsToProcessors(double time, vtkDataSet *source, int sourceID, int ptId, vtkParticleTracerBaseNamespace::ParticleVector &localSeedPoints, int &localAssignedCount)
all the injection/seed points according to which processor they belong to.
bool InsideBounds(double point[])
vtkCharArray * GetParticleSourceIds(vtkPointData *)
virtual bool IsPointDataValid(vtkDataObject *input)
Methods that check that the input arrays are ordered the same on all data sets.
vtkTemporalInterpolatedVelocityField * GetInterpolator()
virtual void InitializeExtraPointDataArrays(vtkPointData *vtkNotUsed(outputPD))
Methods to append values to existing point data arrays that may only be desired on specific concrete ...
bool ComputeDomainExitLocation(double pos[4], double p2[4], double intersection[4], vtkGenericCell *cell)
This is an old routine kept for possible future use.
virtual void ResetCache()
virtual void AssignUniqueIds(vtkParticleTracerBaseNamespace::ParticleVector &localSeedPoints)
give each one a uniqu ID.
virtual int OutputParticles(vtkPolyData *poly)=0
represent and manipulate point attribute data
Definition: vtkPointData.h:38
represent and manipulate 3D points
Definition: vtkPoints.h:40
Superclass for algorithms that produce only polydata as output.
concrete dataset represents vertices, lines, polygons, and triangle strips
Definition: vtkPolyData.h:86
Hold a reference to a vtkObjectBase instance.
A helper class for interpolating between times during particle tracing.
record modification and/or execution time
Definition: vtkTimeStamp.h:36
ParticleVector::iterator ParticleIterator
std::list< ParticleInformation > ParticleDataList
ParticleDataList::iterator ParticleListIterator
std::vector< ParticleInformation > ParticleVector
@ point
Definition: vtkX3D.h:236
@ info
Definition: vtkX3D.h:376
@ port
Definition: vtkX3D.h:447
@ time
Definition: vtkX3D.h:497
@ 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.
vtkGetStringMacro(ExtensionsString)
Returns a string listing all available extensions.
int vtkIdType
Definition: vtkType.h:287