VTK  9.3.20240328
vtkSLACReader.h
Go to the documentation of this file.
1 // SPDX-FileCopyrightText: Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
2 // SPDX-FileCopyrightText: Copyright 2008 Sandia Corporation
3 // SPDX-License-Identifier: LicenseRef-BSD-3-Clause-Sandia-LANL-California-USGov
4 
20 #ifndef vtkSLACReader_h
21 #define vtkSLACReader_h
22 
23 #include "vtkIONetCDFModule.h" // For export macro
25 
26 #include "vtkSmartPointer.h" // For internal method.
27 
28 VTK_ABI_NAMESPACE_BEGIN
30 class vtkDoubleArray;
31 class vtkIdTypeArray;
34 
35 class VTKIONETCDF_EXPORT vtkSLACReader : public vtkMultiBlockDataSetAlgorithm
36 {
37 public:
39  static vtkSLACReader* New();
40  void PrintSelf(ostream& os, vtkIndent indent) override;
41 
42  vtkGetFilePathMacro(MeshFileName);
43  vtkSetFilePathMacro(MeshFileName);
44 
46 
51  virtual void AddModeFileName(VTK_FILEPATH const char* fname);
52  virtual void RemoveAllModeFileNames();
53  virtual unsigned int GetNumberOfModeFileNames();
54  virtual VTK_FILEPATH const char* GetModeFileName(unsigned int idx);
56 
58 
61  vtkGetMacro(ReadInternalVolume, vtkTypeBool);
62  vtkSetMacro(ReadInternalVolume, vtkTypeBool);
63  vtkBooleanMacro(ReadInternalVolume, vtkTypeBool);
65 
67 
70  vtkGetMacro(ReadExternalSurface, vtkTypeBool);
71  vtkSetMacro(ReadExternalSurface, vtkTypeBool);
72  vtkBooleanMacro(ReadExternalSurface, vtkTypeBool);
74 
76 
80  vtkGetMacro(ReadMidpoints, vtkTypeBool);
81  vtkSetMacro(ReadMidpoints, vtkTypeBool);
82  vtkBooleanMacro(ReadMidpoints, vtkTypeBool);
84 
86 
90  virtual const char* GetVariableArrayName(int index);
91  virtual int GetVariableArrayStatus(const char* name);
92  virtual void SetVariableArrayStatus(const char* name, int status);
94 
96 
99  virtual void ResetFrequencyScales();
100  virtual void SetFrequencyScale(int index, double scale);
102 
104 
107  virtual void ResetPhaseShifts();
108  virtual void SetPhaseShift(int index, double shift);
110 
112 
118 
122  static int CanReadFile(VTK_FILEPATH const char* filename);
123 
129 
135 
137 
146 
148 
152  class VTKIONETCDF_EXPORT EdgeEndpoints
153  {
154  public:
156  : MinEndPoint(-1)
157  , MaxEndPoint(-1)
158  {
159  }
160  EdgeEndpoints(vtkIdType endpointA, vtkIdType endpointB)
161  {
162  if (endpointA < endpointB)
163  {
164  this->MinEndPoint = endpointA;
165  this->MaxEndPoint = endpointB;
166  }
167  else
168  {
169  this->MinEndPoint = endpointB;
170  this->MaxEndPoint = endpointA;
171  }
172  }
173  inline vtkIdType GetMinEndPoint() const { return this->MinEndPoint; }
174  inline vtkIdType GetMaxEndPoint() const { return this->MaxEndPoint; }
175  inline bool operator==(const EdgeEndpoints& other) const
176  {
177  return ((this->GetMinEndPoint() == other.GetMinEndPoint()) &&
178  (this->GetMaxEndPoint() == other.GetMaxEndPoint()));
179  }
180 
181  protected:
184  };
186 
188 
191  class VTKIONETCDF_EXPORT MidpointCoordinates
192  {
193  public:
194  MidpointCoordinates() = default;
195  MidpointCoordinates(const double coord[3], vtkIdType id)
196  {
197  this->Coordinate[0] = coord[0];
198  this->Coordinate[1] = coord[1];
199  this->Coordinate[2] = coord[2];
200  this->ID = id;
201  }
202  double Coordinate[3];
204  };
206 
207  enum
208  {
209  SURFACE_OUTPUT = 0,
210  VOLUME_OUTPUT = 1,
211  NUM_OUTPUTS = 2
212  };
213 
214 protected:
216  ~vtkSLACReader() override;
217 
218  class vtkInternal;
219  vtkInternal* Internal;
220 
221  // Friend so vtkInternal can access MidpointIdMap
222  // (so Sun CC compiler doesn't complain).
223  friend class vtkInternal;
224 
226 
230 
235 
240 
245 
247  vtkInformationVector* outputVector) override;
248 
249  int RequestData(vtkInformation* request, vtkInformationVector** inputVector,
250  vtkInformationVector* outputVector) override;
251 
256  vtkObject* caller, unsigned long eid, void* clientdata, void* calldata);
257 
265  virtual vtkIdType GetNumTuplesInVariable(int ncFD, int varId, int expectedNumComponents);
266 
271  virtual int CheckTetrahedraWinding(int meshFD);
272 
277  virtual int ReadConnectivity(
278  int meshFD, vtkMultiBlockDataSet* surfaceOutput, vtkMultiBlockDataSet* volumeOutput);
279 
281 
284  virtual int ReadTetrahedronInteriorArray(int meshFD, vtkIdTypeArray* connectivity);
285  virtual int ReadTetrahedronExteriorArray(int meshFD, vtkIdTypeArray* connectivity);
287 
291  virtual vtkSmartPointer<vtkDataArray> ReadPointDataArray(int ncFD, int varId);
292 
296  enum
297  {
298  NumPerTetInt = 5,
299  NumPerTetExt = 9
300  };
301 
303 
306  class VTKIONETCDF_EXPORT MidpointCoordinateMap
307  {
308  public:
312 
313  void AddMidpoint(const EdgeEndpoints& edge, const MidpointCoordinates& midpoint);
314  void RemoveMidpoint(const EdgeEndpoints& edge);
317 
323 
324  protected:
325  class vtkInternal;
326  vtkInternal* Internal;
327 
328  private:
329  // Too lazy to implement these.
331  void operator=(const MidpointCoordinateMap&) = delete;
332  };
333 
335 
338  class VTKIONETCDF_EXPORT MidpointIdMap
339  {
340  public:
344 
345  void AddMidpoint(const EdgeEndpoints& edge, vtkIdType midpoint);
346  void RemoveMidpoint(const EdgeEndpoints& edge);
349 
354 
362  bool GetNextMidpoint(EdgeEndpoints& edge, vtkIdType& midpoint);
363 
364  protected:
365  class vtkInternal;
366  vtkInternal* Internal;
367 
368  private:
369  // Too lazy to implement these.
370  MidpointIdMap(const MidpointIdMap&) = delete;
371  void operator=(const MidpointIdMap&) = delete;
372  };
373 
378  virtual int ReadCoordinates(int meshFD, vtkMultiBlockDataSet* output);
379 
386  int meshFD, vtkMultiBlockDataSet* output, MidpointCoordinateMap& map);
387 
393  virtual int ReadMidpointData(
394  int meshFD, vtkMultiBlockDataSet* output, MidpointIdMap& midpointIds);
395 
400  virtual int RestoreMeshCache(vtkMultiBlockDataSet* surfaceOutput,
401  vtkMultiBlockDataSet* volumeOutput, vtkMultiBlockDataSet* compositeOutput);
402 
407  virtual int ReadFieldData(const int* modeFDArray, int numModeFDs, vtkMultiBlockDataSet* output);
408 
414 
421 
426  virtual int MeshUpToDate();
427 
428 private:
429  vtkSLACReader(const vtkSLACReader&) = delete;
430  void operator=(const vtkSLACReader&) = delete;
431 };
432 
433 VTK_ABI_NAMESPACE_END
434 #endif // vtkSLACReader_h
Store on/off settings for data arrays, etc.
dynamic, self-adjusting array of double
dynamic, self-adjusting array of vtkIdType
a simple class to control print indentation
Definition: vtkIndent.h:108
Key for integer values in vtkInformation.
Key for vtkObjectBase values.
Store zero or more vtkInformation instances.
Store vtkAlgorithm input/output information.
Superclass for algorithms that produce only vtkMultiBlockDataSet as output.
Composite dataset that organizes datasets into blocks.
abstract base class for most VTK objects
Definition: vtkObject.h:161
Simple class used internally to define an edge based on the endpoints.
vtkIdType GetMinEndPoint() const
vtkIdType GetMaxEndPoint() const
bool operator==(const EdgeEndpoints &other) const
EdgeEndpoints(vtkIdType endpointA, vtkIdType endpointB)
Manages a map from edges to midpoint coordinates.
MidpointCoordinates * FindMidpoint(const EdgeEndpoints &edge)
Finds the coordinates for the given edge or returns nullptr if it does not exist.
void RemoveMidpoint(const EdgeEndpoints &edge)
void AddMidpoint(const EdgeEndpoints &edge, const MidpointCoordinates &midpoint)
Simple class used internally for holding midpoint information.
MidpointCoordinates(const double coord[3], vtkIdType id)
Manages a map from edges to the point id of the midpoint.
vtkIdType * FindMidpoint(const EdgeEndpoints &edge)
Finds the id for the given edge or returns nullptr if it does not exist.
vtkIdType GetNumberOfMidpoints() const
void AddMidpoint(const EdgeEndpoints &edge, vtkIdType midpoint)
void InitTraversal()
Initialize iteration.
bool GetNextMidpoint(EdgeEndpoints &edge, vtkIdType &midpoint)
Get the next midpoint in the iteration.
void RemoveMidpoint(const EdgeEndpoints &edge)
A reader for a data format used by Omega3p, Tau3p, and several other tools used at the Standford Line...
Definition: vtkSLACReader.h:36
virtual void ResetFrequencyScales()
Sets the scale factor for each mode.
int RequestData(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector) override
This is called by the superclass.
bool TimeStepModes
True if "mode" files are a sequence of time steps.
virtual vtkDoubleArray * GetPhaseShifts()
NOTE: This is not thread-safe.
static vtkInformationIntegerKey * IS_EXTERNAL_SURFACE()
This key is attached to the metadata information of all data sets in the output that are part of the ...
vtkInternal * Internal
virtual void SetPhaseShift(int index, double shift)
Sets the phase offset for each mode.
virtual void SetFrequencyScale(int index, double scale)
Sets the scale factor for each mode.
virtual int ReadMidpointCoordinates(int meshFD, vtkMultiBlockDataSet *output, MidpointCoordinateMap &map)
Reads in the midpoint coordinate data from the mesh file and returns a map from edges to midpoints.
static void SelectionModifiedCallback(vtkObject *caller, unsigned long eid, void *clientdata, void *calldata)
Callback registered with the VariableArraySelection.
static int CanReadFile(VTK_FILEPATH const char *filename)
Returns true if the given file can be read by this reader.
vtkTypeBool ReadMidpoints
virtual vtkDoubleArray * GetFrequencyScales()
NOTE: This is not thread-safe.
~vtkSLACReader() override
vtkGetFilePathMacro(MeshFileName)
static vtkInformationObjectBaseKey * POINTS()
All the data sets stored in the multiblock output share the same point data.
virtual int RestoreMeshCache(vtkMultiBlockDataSet *surfaceOutput, vtkMultiBlockDataSet *volumeOutput, vtkMultiBlockDataSet *compositeOutput)
Instead of reading data from the mesh file, restore the data from the previous mesh file read.
virtual int ReadTetrahedronExteriorArray(int meshFD, vtkIdTypeArray *connectivity)
Reads tetrahedron connectivity arrays.
virtual int ReadFieldData(const int *modeFDArray, int numModeFDs, vtkMultiBlockDataSet *output)
Read in the field data from the mode file.
virtual void ResetPhaseShifts()
Sets the phase offset for each mode.
static vtkInformationIntegerKey * IS_INTERNAL_VOLUME()
This key is attached to the metadata information of all data sets in the output that are part of the ...
virtual int MeshUpToDate()
Returns 1 if the mesh is up to date, 0 if the mesh needs to be read from disk.
char * MeshFileName
virtual int ReadMidpointData(int meshFD, vtkMultiBlockDataSet *output, MidpointIdMap &midpointIds)
Read in the midpoint data from the mesh file.
virtual unsigned int GetNumberOfModeFileNames()
There may be one mode file (usually for actual modes) or multiple mode files (which usually actually ...
bool FrequencyModes
True if mode files describe vibrating fields.
virtual VTK_FILEPATH const char * GetModeFileName(unsigned int idx)
There may be one mode file (usually for actual modes) or multiple mode files (which usually actually ...
static vtkInformationObjectBaseKey * POINT_DATA()
All the data sets stored in the multiblock output share the same point data.
virtual int ReadCoordinates(int meshFD, vtkMultiBlockDataSet *output)
Read in the point coordinate data from the mesh file.
virtual vtkSmartPointer< vtkDataArray > ReadPointDataArray(int ncFD, int varId)
Reads point data arrays.
virtual int ReadTetrahedronInteriorArray(int meshFD, vtkIdTypeArray *connectivity)
Reads tetrahedron connectivity arrays.
static vtkSLACReader * New()
vtkTimeStamp MeshReadTime
A time stamp for the last time the mesh file was read.
virtual int ReadConnectivity(int meshFD, vtkMultiBlockDataSet *surfaceOutput, vtkMultiBlockDataSet *volumeOutput)
Read the connectivity information from the mesh file.
bool ReadModeData
True if reading from a proper mode file.
virtual int GetVariableArrayStatus(const char *name)
Variable array selection.
int RequestInformation(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector) override
This is called by the superclass.
virtual vtkIdType GetNumTuplesInVariable(int ncFD, int varId, int expectedNumComponents)
Convenience function that checks the dimensions of a 2D netCDF array that is supposed to be a set of ...
vtkTypeBool ReadInternalVolume
virtual const char * GetVariableArrayName(int index)
Variable array selection.
vtkTypeBool ReadExternalSurface
virtual int InterpolateMidpointData(vtkMultiBlockDataSet *output, MidpointIdMap &map)
Takes the data read on the fields and interpolates data for the midpoints.
virtual void AddModeFileName(VTK_FILEPATH const char *fname)
There may be one mode file (usually for actual modes) or multiple mode files (which usually actually ...
virtual int GetNumberOfVariableArrays()
Variable array selection.
virtual void RemoveAllModeFileNames()
There may be one mode file (usually for actual modes) or multiple mode files (which usually actually ...
void PrintSelf(ostream &os, vtkIndent indent) override
Methods invoked by print to print information about the object including superclasses.
vtkSetFilePathMacro(MeshFileName)
virtual int CheckTetrahedraWinding(int meshFD)
Checks the winding of the tetrahedra in the mesh file.
virtual void SetVariableArrayStatus(const char *name, int status)
Variable array selection.
record modification and/or execution time
Definition: vtkTimeStamp.h:44
@ Coordinate
Definition: vtkX3D.h:44
@ scale
Definition: vtkX3D.h:229
@ name
Definition: vtkX3D.h:219
@ index
Definition: vtkX3D.h:246
int vtkTypeBool
Definition: vtkABI.h:64
int vtkIdType
Definition: vtkType.h:315
#define VTK_FILEPATH