VTK
vtkNetCDFCFReader.h
Go to the documentation of this file.
1 // -*- c++ -*-
2 /*=========================================================================
3 
4  Program: Visualization Toolkit
5  Module: vtkNetCDFCFReader.h
6 
7  Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
8  All rights reserved.
9  See Copyright.txt or http://www.kitware.com/Copyright.htm for details.
10 
11  This software is distributed WITHOUT ANY WARRANTY; without even
12  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
13  PURPOSE. See the above copyright notice for more information.
14 
15 =========================================================================*/
16 
17 /*-------------------------------------------------------------------------
18  Copyright 2008 Sandia Corporation.
19  Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
20  the U.S. Government retains certain rights in this software.
21 -------------------------------------------------------------------------*/
22 
33 #ifndef vtkNetCDFCFReader_h
34 #define vtkNetCDFCFReader_h
35 
36 #include "vtkIONetCDFModule.h" // For export macro
37 #include "vtkNetCDFReader.h"
38 
39 #include "vtkStdString.h" // Used for ivars.
40 
41 class vtkImageData;
42 class vtkPoints;
43 class vtkRectilinearGrid;
44 class vtkStructuredGrid;
46 
47 class VTKIONETCDF_EXPORT vtkNetCDFCFReader : public vtkNetCDFReader
48 {
49 public:
51  static vtkNetCDFCFReader *New();
52  void PrintSelf(ostream &os, vtkIndent indent) override;
53 
55 
60  vtkGetMacro(SphericalCoordinates, vtkTypeBool);
61  vtkSetMacro(SphericalCoordinates, vtkTypeBool);
62  vtkBooleanMacro(SphericalCoordinates, vtkTypeBool);
64 
66 
77  vtkGetMacro(VerticalScale, double);
78  vtkSetMacro(VerticalScale, double);
79  vtkGetMacro(VerticalBias, double);
80  vtkSetMacro(VerticalBias, double);
82 
84 
91  vtkGetMacro(OutputType, int);
92  virtual void SetOutputType(int type);
93  void SetOutputTypeToAutomatic() { this->SetOutputType(-1); }
94  void SetOutputTypeToImage() { this->SetOutputType(VTK_IMAGE_DATA); }
95  void SetOutputTypeToRectilinear() {this->SetOutputType(VTK_RECTILINEAR_GRID);}
96  void SetOutputTypeToStructured() { this->SetOutputType(VTK_STRUCTURED_GRID); }
98  this->SetOutputType(VTK_UNSTRUCTURED_GRID);
99  }
101 
105  static int CanReadFile(const char *filename);
106 
107 protected:
109  ~vtkNetCDFCFReader() override;
110 
112 
114  double VerticalBias;
115 
117 
118  int RequestDataObject(vtkInformation *request,
119  vtkInformationVector **inputVector,
120  vtkInformationVector *outputVector) override;
121 
122  int RequestInformation(vtkInformation *request,
123  vtkInformationVector **inputVector,
124  vtkInformationVector *outputVector) override;
125 
126  int RequestData(vtkInformation *request,
127  vtkInformationVector **inputVector,
128  vtkInformationVector *outputVector) override;
129 
131 
134  int ReadMetaData(int ncFD) override;
135  int IsTimeDimension(int ncFD, int dimId) override;
136  vtkSmartPointer<vtkDoubleArray> GetTimeValues(int ncFD, int dimId) override;
138 
140  public:
142  vtkDimensionInfo(int ncFD, int id);
143  const char *GetName() const { return this->Name.c_str(); }
144  enum UnitsEnum {
149  VERTICAL_UNITS
150  };
151  UnitsEnum GetUnits() const { return this->Units; }
152  vtkSmartPointer<vtkDoubleArray> GetCoordinates() {return this->Coordinates;}
153  vtkSmartPointer<vtkDoubleArray> GetBounds() { return this->Bounds; }
154  bool GetHasRegularSpacing() const { return this->HasRegularSpacing; }
155  double GetOrigin() const { return this->Origin; }
156  double GetSpacing() const { return this->Spacing; }
158  return this->SpecialVariables;
159  }
160  protected:
162  int DimId;
167  double Origin, Spacing;
169  int LoadMetaData(int ncFD);
170  };
171  class vtkDimensionInfoVector;
172  friend class vtkDimensionInfoVector;
173  vtkDimensionInfoVector *DimensionInfo;
174  vtkDimensionInfo *GetDimensionInfo(int dimension);
175 
177  public:
178  vtkDependentDimensionInfo() : Valid(false) { };
179  vtkDependentDimensionInfo(int ncFD, int varId, vtkNetCDFCFReader *parent);
180  bool GetValid() const { return this->Valid; }
181  bool GetHasBounds() const { return this->HasBounds; }
182  bool GetCellsUnstructured() const { return this->CellsUnstructured; }
184  return this->GridDimensions;
185  }
187  return this->LongitudeCoordinates;
188  }
190  return this->LatitudeCoordinates;
191  }
193  return this->SpecialVariables;
194  }
195  protected:
196  bool Valid;
197  bool HasBounds;
203  int LoadMetaData(int ncFD, int varId, vtkNetCDFCFReader *parent);
204  int LoadCoordinateVariable(int ncFD, int varId, vtkDoubleArray *coords);
205  int LoadBoundsVariable(int ncFD, int varId, vtkDoubleArray *coords);
206  int LoadUnstructuredBoundsVariable(int ncFD, int varId,
207  vtkDoubleArray *coords);
208  };
210  class vtkDependentDimensionInfoVector;
211  friend class vtkDependentDimensionInfoVector;
212  vtkDependentDimensionInfoVector *DependentDimensionInfo;
213 
214  // Finds the dependent dimension information for the given set of dimensions.
215  // Returns nullptr if no information has been recorded.
216  vtkDependentDimensionInfo *FindDependentDimensionInfo(vtkIntArray *dims);
217 
223  virtual void IdentifySphericalCoordinates(vtkIntArray *dimensions,
224  int &longitudeDim,
225  int &latitudeDim,
226  int &verticalDim);
227 
237  COORDS_SPHERICAL_PSIDED_CELLS
238  };
239 
245  CoordinateTypesEnum CoordinateType(vtkIntArray *dimensions);
246 
250  bool DimensionsAreForPointData(vtkIntArray *dimensions) override;
251 
257  void ExtentForDimensionsAndPiece(int pieceNumber,
258  int numberOfPieces,
259  int ghostLevels,
260  int extent[6]);
261 
265  void GetUpdateExtentForOutput(vtkDataSet *output, int extent[6]) override;
266 
268 
271  void AddRectilinearCoordinates(vtkImageData *imageOutput);
272  void AddRectilinearCoordinates(vtkRectilinearGrid *rectilinearOutput);
273  void FakeRectilinearCoordinates(vtkRectilinearGrid *rectilinearOutput);
274  void Add1DRectilinearCoordinates(vtkPoints *points, const int extent[6]);
275  void Add2DRectilinearCoordinates(vtkPoints *points, const int extent[6]);
276  void Add1DRectilinearCoordinates(vtkStructuredGrid *structuredOutput);
277  void Add2DRectilinearCoordinates(vtkStructuredGrid *structuredOutput);
278  void FakeStructuredCoordinates(vtkStructuredGrid *structuredOutput);
279  void Add1DRectilinearCoordinates(vtkUnstructuredGrid *unstructuredOutput,
280  const int extent[6]);
281  void Add2DRectilinearCoordinates(vtkUnstructuredGrid *unstructuredOutput,
282  const int extent[6]);
284 
286 
289  void Add1DSphericalCoordinates(vtkPoints *points, const int extent[6]);
290  void Add2DSphericalCoordinates(vtkPoints *points, const int extent[6]);
291  void Add1DSphericalCoordinates(vtkStructuredGrid *structuredOutput);
292  void Add2DSphericalCoordinates(vtkStructuredGrid *structuredOutput);
293  void Add1DSphericalCoordinates(vtkUnstructuredGrid *unstructuredOutput,
294  const int extent[6]);
295  void Add2DSphericalCoordinates(vtkUnstructuredGrid *unstructuredOutput,
296  const int extent[6]);
298 
302  void AddStructuredCells(vtkUnstructuredGrid *unstructuredOutput,
303  const int extent[6]);
304 
306 
309  void AddUnstructuredRectilinearCoordinates(
310  vtkUnstructuredGrid *unstructuredOutput,
311  const int extent[6]);
312  void AddUnstructuredSphericalCoordinates(
313  vtkUnstructuredGrid *unstructuredOutput,
314  const int extent[6]);
316 
317 
318 private:
319  vtkNetCDFCFReader(const vtkNetCDFCFReader &) = delete;
320  void operator=(const vtkNetCDFCFReader &) = delete;
321 };
322 
323 #endif //vtkNetCDFCFReader_h
324 
Wrapper around std::string to keep symbols short.
Definition: vtkStdString.h:34
vtkSmartPointer< vtkDoubleArray > Coordinates
vtkSmartPointer< vtkDoubleArray > GetCoordinates()
a dataset that is topologically regular with variable spacing in the three coordinate directions
vtkSmartPointer< vtkStringArray > GetSpecialVariables() const
#define VTK_IMAGE_DATA
Definition: vtkType.h:97
int RequestData(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector) override
vtkSmartPointer< vtkDoubleArray > LatitudeCoordinates
void SetOutputTypeToStructured()
Set/get the data type of the output.
#define VTK_RECTILINEAR_GRID
Definition: vtkType.h:94
Store vtkAlgorithm input/output information.
Reads netCDF files that follow the CF convention.
abstract class to specify dataset behavior
Definition: vtkDataSet.h:56
virtual bool DimensionsAreForPointData(vtkIntArray *vtkNotUsed(dimensions))
Called internally to determine whether a variable with the given set of dimensions should be loaded a...
vtkDimensionInfoVector * DimensionInfo
vtkSmartPointer< vtkStringArray > SpecialVariables
void SetOutputTypeToUnstructured()
Set/get the data type of the output.
vtkSmartPointer< vtkDoubleArray > GetLongitudeCoordinates() const
void SetOutputTypeToImage()
Set/get the data type of the output.
virtual int ReadMetaData(int ncFD)
Reads meta data and populates ivars.
A superclass for reading netCDF files.
dynamic, self-adjusting array of double
int RequestInformation(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector) override
int RequestDataObject(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector) override
This is called by the superclass.
int vtkTypeBool
Definition: vtkABI.h:69
dynamic, self-adjusting array of int
Definition: vtkIntArray.h:39
void SetOutputTypeToAutomatic()
Set/get the data type of the output.
static vtkNetCDFReader * New()
virtual void GetUpdateExtentForOutput(vtkDataSet *output, int extent[6])
Retrieves the update extent for the output object.
a simple class to control print indentation
Definition: vtkIndent.h:33
topologically and geometrically regular array of data
Definition: vtkImageData.h:39
dataset represents arbitrary combinations of all possible cell types
vtkSmartPointer< vtkDoubleArray > GetLatitudeCoordinates() const
vtkSmartPointer< vtkIntArray > GetGridDimensions() const
virtual int IsTimeDimension(int ncFD, int dimId)
Determines whether the given variable is a time dimension.
vtkSmartPointer< vtkStringArray > SpecialVariables
vtkTypeBool SphericalCoordinates
vtkSmartPointer< vtkDoubleArray > GetBounds()
virtual vtkSmartPointer< vtkDoubleArray > GetTimeValues(int ncFD, int dimId)
Given a dimension already determined to be a time dimension (via a call to IsTimeDimension) returns a...
void PrintSelf(ostream &os, vtkIndent indent) override
Methods invoked by print to print information about the object including superclasses.
vtkSmartPointer< vtkDoubleArray > LongitudeCoordinates
vtkDependentDimensionInfoVector * DependentDimensionInfo
topologically regular array of data
void SetOutputTypeToRectilinear()
Set/get the data type of the output.
Store zero or more vtkInformation instances.
vtkSmartPointer< vtkStringArray > GetSpecialVariables() const
vtkSmartPointer< vtkDoubleArray > Bounds
vtkSmartPointer< vtkIntArray > GridDimensions
#define VTK_STRUCTURED_GRID
Definition: vtkType.h:93
represent and manipulate 3D points
Definition: vtkPoints.h:33
#define VTK_UNSTRUCTURED_GRID
Definition: vtkType.h:95