VTK  9.3.1
vtkMPASReader.h
Go to the documentation of this file.
1 // SPDX-FileCopyrightText: Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
2 // SPDX-FileCopyrightText: Copyright (c) 2002-2005 Los Alamos National Laboratory
3 // SPDX-License-Identifier: BSD-3-Clause-Sandia-LANL-California-USGov
39 #ifndef vtkMPASReader_h
40 #define vtkMPASReader_h
41 
42 #include "vtkIONetCDFModule.h" // For export macro
44 
45 #include <string> // for std::string
46 
47 VTK_ABI_NAMESPACE_BEGIN
48 class vtkCallbackCommand;
50 class vtkDoubleArray;
51 class vtkStringArray;
52 
53 class VTKIONETCDF_EXPORT vtkMPASReader : public vtkUnstructuredGridAlgorithm
54 {
55 public:
56  static vtkMPASReader* New();
58  void PrintSelf(ostream& os, vtkIndent indent) override;
59 
61 
64  vtkSetFilePathMacro(FileName);
65  vtkGetFilePathMacro(FileName);
67 
69 
72  vtkGetMacro(MaximumCells, int);
74 
76 
79  vtkGetMacro(MaximumPoints, int);
81 
83 
86  virtual int GetNumberOfCellVars();
87  virtual int GetNumberOfPointVars();
89 
91 
97 
99 
105  vtkSetMacro(UseDimensionedArrayNames, bool);
106  vtkGetMacro(UseDimensionedArrayNames, bool);
107  vtkBooleanMacro(UseDimensionedArrayNames, bool);
109 
111 
116  int GetNumberOfPointArrays();
117  const char* GetPointArrayName(int index);
118  int GetPointArrayStatus(const char* name);
119  void SetPointArrayStatus(const char* name, int status);
120  void DisableAllPointArrays();
121  void EnableAllPointArrays();
123 
124  int GetNumberOfCellArrays();
125  const char* GetCellArrayName(int index);
126  int GetCellArrayStatus(const char* name);
127  void SetCellArrayStatus(const char* name, int status);
128  void DisableAllCellArrays();
129  void EnableAllCellArrays();
130 
132 
140  vtkIdType GetNumberOfDimensions();
141  std::string GetDimensionName(int idx);
142  vtkStringArray* GetAllDimensions();
143  int GetDimensionCurrentIndex(const std::string& dim);
144  void SetDimensionCurrentIndex(const std::string& dim, int idx);
145  int GetDimensionSize(const std::string& dim);
147 
149 
153  vtkSetMacro(VerticalDimension, std::string);
154  vtkGetMacro(VerticalDimension, std::string);
156 
158 
162  void SetVerticalLevel(int level);
163  int GetVerticalLevel();
165 
166  vtkGetVector2Macro(VerticalLevelRange, int);
167 
168  vtkSetMacro(LayerThickness, int);
169  vtkGetMacro(LayerThickness, int);
170  vtkGetVector2Macro(LayerThicknessRange, int);
171 
172  void SetCenterLon(int val);
173  vtkGetVector2Macro(CenterLonRange, int);
174 
175  vtkSetMacro(ProjectLatLon, bool);
176  vtkGetMacro(ProjectLatLon, bool);
177 
178  vtkSetMacro(IsAtmosphere, bool);
179  vtkGetMacro(IsAtmosphere, bool);
180 
181  vtkSetMacro(IsZeroCentered, bool);
182  vtkGetMacro(IsZeroCentered, bool);
183 
184  vtkSetMacro(ShowMultilayerView, bool);
185  vtkGetMacro(ShowMultilayerView, bool);
186 
190  static int CanReadFile(VTK_FILEPATH const char* filename);
191 
192  vtkMTimeType GetMTime() override;
193 
194 protected:
195  vtkMPASReader();
196  ~vtkMPASReader() override;
197  void ReleaseNcData();
198  void DestroyData();
199 
200  char* FileName; // First field part file giving path
201 
202  size_t NumberOfTimeSteps; // Temporal domain
203  double DTime; // The current time
204 
205  // Observer to modify this object when array selections are modified
207 
210 
211  static void SelectionCallback(
212  vtkObject* caller, unsigned long eid, void* clientdata, void* calldata);
213 
214  // Selected field of interest
217 
222  void UpdateDimensions(bool force = false);
223 
225  int VerticalLevelRange[2];
226 
228  int LayerThicknessRange[2];
229 
231  int CenterLonRange[2];
232 
234  {
237  Planar
238  };
239 
241 
242  bool ProjectLatLon; // User option
243  bool OnASphere; // Data file attribute
247 
249  bool DoBugFix;
250  double CenterRad;
251 
253 
254  // geometry
259  size_t PointOffset;
261  size_t CurrentExtraPoint; // current extra point
262  size_t CurrentExtraCell; // current extra cell
263  double* PointX; // x coord of point
264  double* PointY; // y coord of point
265  double* PointZ; // z coord of point
266  size_t ModNumPoints;
267  size_t ModNumCells;
268  int* OrigConnections; // original connections
269  int* ModConnections; // modified connections
270  size_t* CellMap; // maps from added cell to original cell #
271  size_t* PointMap; // maps from added point to original point #
273  int MaximumCells; // max cells
274  int MaximumPoints; // max points
275 
276  void SetDefaults();
277  int GetNcDims();
278  int GetNcAtts();
279  int CheckParams();
280  int GetNcVars(const char* cellDimName, const char* pointDimName);
281  int ReadAndOutputGrid();
282  int BuildVarArrays();
283  int AllocSphericalGeometry();
284  int AllocProjectedGeometry();
285  int AllocPlanarGeometry();
286  void ShiftLonData();
287  int AddMirrorPoint(int index, double dividerX, double offset);
288  void FixPoints();
289  int EliminateXWrap();
290  void OutputPoints();
291  void OutputCells();
292  unsigned char GetCellType();
293 
294  vtkDataArray* LoadPointVarData(int variable);
295  vtkDataArray* LoadCellVarData(int variable);
296  vtkDataArray* LookupPointDataArray(int varIdx);
297  vtkDataArray* LookupCellDataArray(int varIdx);
298 
307  void LoadTimeFieldData(vtkUnstructuredGrid* dataset);
308 
309 private:
310  vtkMPASReader(const vtkMPASReader&) = delete;
311  void operator=(const vtkMPASReader&) = delete;
312 
313  class Internal;
314  Internal* Internals;
315 };
316 
317 VTK_ABI_NAMESPACE_END
318 #endif
GeometryType Geometry
virtual int RequestData(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector)
This is called by the superclass.
size_t PointsPerCell
size_t ModNumPoints
abstract base class for most VTK objects
Definition: vtkObject.h:51
Store vtkAlgorithm input/output information.
vtkTypeUInt32 vtkMTimeType
Definition: vtkType.h:270
size_t * CellMap
size_t NumberOfTimeSteps
static vtkUnstructuredGridAlgorithm * New()
size_t PointOffset
double * PointZ
a vtkAbstractArray subclass for strings
bool IncludeTopography
virtual int RequestInformation(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector)
size_t CurrentExtraPoint
int * OrigConnections
int vtkIdType
Definition: vtkType.h:315
vtkDataArraySelection * PointDataArraySelection
Read an MPAS netCDF file.
Definition: vtkMPASReader.h:53
dynamic, self-adjusting array of double
double * PointX
bool UseDimensionedArrayNames
void PrintSelf(ostream &os, vtkIndent indent) override
Methods invoked by print to print information about the object including superclasses.
supports function callbacks
double * PointY
a simple class to control print indentation
Definition: vtkIndent.h:28
Store on/off settings for data arrays, etc.
dataset represents arbitrary combinations of all possible cell types
abstract superclass for arrays of numeric data
Definition: vtkDataArray.h:44
virtual vtkMTimeType GetMTime()
Return this object's modified time.
size_t NumberOfPoints
size_t CurrentExtraCell
Superclass for algorithms that produce only unstructured grid as output.
std::string VerticalDimension
int GetCellType(const Ioss::ElementTopology *topology)
Returns VTK celltype for a Ioss topology element.
#define VTK_FILEPATH
size_t MaximumNVertLevels
vtkCallbackCommand * SelectionObserver
Store zero or more vtkInformation instances.
int * MaximumLevelPoint
vtkUnstructuredGrid * GetOutput()
Get the output data object for a port on this algorithm.
int * ModConnections
size_t ModNumCells
vtkDataArraySelection * CellDataArraySelection
size_t NumberOfCells
size_t * PointMap
bool ShowMultilayerView