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 
36 #ifndef vtkNetCDFCFReader_h
37 #define vtkNetCDFCFReader_h
38 
39 #include "vtkIONetCDFModule.h" // For export macro
40 #include "vtkNetCDFReader.h"
41 
42 #include "vtkStdString.h" // Used for ivars.
43 
44 class vtkImageData;
45 class vtkPoints;
46 class vtkRectilinearGrid;
47 class vtkStructuredGrid;
49 
50 class VTKIONETCDF_EXPORT vtkNetCDFCFReader : public vtkNetCDFReader
51 {
52 public:
54  static vtkNetCDFCFReader *New();
55  void PrintSelf(ostream &os, vtkIndent indent) override;
56 
58 
63  vtkGetMacro(SphericalCoordinates, vtkTypeBool);
64  vtkSetMacro(SphericalCoordinates, vtkTypeBool);
65  vtkBooleanMacro(SphericalCoordinates, vtkTypeBool);
67 
69 
80  vtkGetMacro(VerticalScale, double);
81  vtkSetMacro(VerticalScale, double);
82  vtkGetMacro(VerticalBias, double);
83  vtkSetMacro(VerticalBias, double);
85 
87 
94  vtkGetMacro(OutputType, int);
95  virtual void SetOutputType(int type);
96  void SetOutputTypeToAutomatic() { this->SetOutputType(-1); }
97  void SetOutputTypeToImage() { this->SetOutputType(VTK_IMAGE_DATA); }
98  void SetOutputTypeToRectilinear() {this->SetOutputType(VTK_RECTILINEAR_GRID);}
99  void SetOutputTypeToStructured() { this->SetOutputType(VTK_STRUCTURED_GRID); }
101  this->SetOutputType(VTK_UNSTRUCTURED_GRID);
102  }
104 
108  static int CanReadFile(const char *filename);
109 
110 protected:
112  ~vtkNetCDFCFReader() override;
113 
115 
117  double VerticalBias;
118 
120 
121  int RequestDataObject(vtkInformation *request,
122  vtkInformationVector **inputVector,
123  vtkInformationVector *outputVector) override;
124 
125  int RequestInformation(vtkInformation *request,
126  vtkInformationVector **inputVector,
127  vtkInformationVector *outputVector) override;
128 
129  int RequestData(vtkInformation *request,
130  vtkInformationVector **inputVector,
131  vtkInformationVector *outputVector) override;
132 
134 
137  int ReadMetaData(int ncFD) override;
138  int IsTimeDimension(int ncFD, int dimId) override;
139  vtkSmartPointer<vtkDoubleArray> GetTimeValues(int ncFD, int dimId) override;
141 
143  public:
145  vtkDimensionInfo(int ncFD, int id);
146  const char *GetName() const { return this->Name.c_str(); }
147  enum UnitsEnum {
152  VERTICAL_UNITS
153  };
154  UnitsEnum GetUnits() const { return this->Units; }
155  vtkSmartPointer<vtkDoubleArray> GetCoordinates() {return this->Coordinates;}
156  vtkSmartPointer<vtkDoubleArray> GetBounds() { return this->Bounds; }
157  bool GetHasRegularSpacing() const { return this->HasRegularSpacing; }
158  double GetOrigin() const { return this->Origin; }
159  double GetSpacing() const { return this->Spacing; }
161  return this->SpecialVariables;
162  }
163  protected:
165  int DimId;
170  double Origin, Spacing;
172  int LoadMetaData(int ncFD);
173  };
174  class vtkDimensionInfoVector;
175  friend class vtkDimensionInfoVector;
176  vtkDimensionInfoVector *DimensionInfo;
177  vtkDimensionInfo *GetDimensionInfo(int dimension);
178 
180  public:
181  vtkDependentDimensionInfo() : Valid(false) { };
182  vtkDependentDimensionInfo(int ncFD, int varId, vtkNetCDFCFReader *parent);
183  bool GetValid() const { return this->Valid; }
184  bool GetHasBounds() const { return this->HasBounds; }
185  bool GetCellsUnstructured() const { return this->CellsUnstructured; }
187  return this->GridDimensions;
188  }
190  return this->LongitudeCoordinates;
191  }
193  return this->LatitudeCoordinates;
194  }
196  return this->SpecialVariables;
197  }
198  protected:
199  bool Valid;
200  bool HasBounds;
206  int LoadMetaData(int ncFD, int varId, vtkNetCDFCFReader *parent);
207  int LoadCoordinateVariable(int ncFD, int varId, vtkDoubleArray *coords);
208  int LoadBoundsVariable(int ncFD, int varId, vtkDoubleArray *coords);
209  int LoadUnstructuredBoundsVariable(int ncFD, int varId,
210  vtkDoubleArray *coords);
211  };
213  class vtkDependentDimensionInfoVector;
214  friend class vtkDependentDimensionInfoVector;
215  vtkDependentDimensionInfoVector *DependentDimensionInfo;
216 
217  // Finds the dependent dimension information for the given set of dimensions.
218  // Returns nullptr if no information has been recorded.
219  vtkDependentDimensionInfo *FindDependentDimensionInfo(vtkIntArray *dims);
220 
226  virtual void IdentifySphericalCoordinates(vtkIntArray *dimensions,
227  int &longitudeDim,
228  int &latitudeDim,
229  int &verticalDim);
230 
240  COORDS_SPHERICAL_PSIDED_CELLS
241  };
242 
248  CoordinateTypesEnum CoordinateType(vtkIntArray *dimensions);
249 
253  bool DimensionsAreForPointData(vtkIntArray *dimensions) override;
254 
260  void ExtentForDimensionsAndPiece(int pieceNumber,
261  int numberOfPieces,
262  int ghostLevels,
263  int extent[6]);
264 
268  void GetUpdateExtentForOutput(vtkDataSet *output, int extent[6]) override;
269 
271 
274  void AddRectilinearCoordinates(vtkImageData *imageOutput);
275  void AddRectilinearCoordinates(vtkRectilinearGrid *rectilinearOutput);
276  void FakeRectilinearCoordinates(vtkRectilinearGrid *rectilinearOutput);
277  void Add1DRectilinearCoordinates(vtkPoints *points, const int extent[6]);
278  void Add2DRectilinearCoordinates(vtkPoints *points, const int extent[6]);
279  void Add1DRectilinearCoordinates(vtkStructuredGrid *structuredOutput);
280  void Add2DRectilinearCoordinates(vtkStructuredGrid *structuredOutput);
281  void FakeStructuredCoordinates(vtkStructuredGrid *structuredOutput);
282  void Add1DRectilinearCoordinates(vtkUnstructuredGrid *unstructuredOutput,
283  const int extent[6]);
284  void Add2DRectilinearCoordinates(vtkUnstructuredGrid *unstructuredOutput,
285  const int extent[6]);
287 
289 
292  void Add1DSphericalCoordinates(vtkPoints *points, const int extent[6]);
293  void Add2DSphericalCoordinates(vtkPoints *points, const int extent[6]);
294  void Add1DSphericalCoordinates(vtkStructuredGrid *structuredOutput);
295  void Add2DSphericalCoordinates(vtkStructuredGrid *structuredOutput);
296  void Add1DSphericalCoordinates(vtkUnstructuredGrid *unstructuredOutput,
297  const int extent[6]);
298  void Add2DSphericalCoordinates(vtkUnstructuredGrid *unstructuredOutput,
299  const int extent[6]);
301 
305  void AddStructuredCells(vtkUnstructuredGrid *unstructuredOutput,
306  const int extent[6]);
307 
309 
312  void AddUnstructuredRectilinearCoordinates(
313  vtkUnstructuredGrid *unstructuredOutput,
314  const int extent[6]);
315  void AddUnstructuredSphericalCoordinates(
316  vtkUnstructuredGrid *unstructuredOutput,
317  const int extent[6]);
319 
320 
321 private:
322  vtkNetCDFCFReader(const vtkNetCDFCFReader &) = delete;
323  void operator=(const vtkNetCDFCFReader &) = delete;
324 };
325 
326 #endif //vtkNetCDFCFReader_h
327 
Wrapper around std::string to keep symbols short.
Definition: vtkStdString.h:40
vtkSmartPointer< vtkDoubleArray > Coordinates
vtkSmartPointer< vtkDoubleArray > GetCoordinates()
a dataset that is topologically regular with variable spacing in the three coordinate directions ...
#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:62
vtkDimensionInfoVector * DimensionInfo
vtkSmartPointer< vtkStringArray > SpecialVariables
void SetOutputTypeToUnstructured()
Set/get the data type of the output.
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:45
void SetOutputTypeToAutomatic()
Set/get the data type of the output.
static vtkNetCDFReader * New()
a simple class to control print indentation
Definition: vtkIndent.h:39
vtkSmartPointer< vtkDoubleArray > GetLatitudeCoordinates() const
topologically and geometrically regular array of data
Definition: vtkImageData.h:45
dataset represents arbitrary combinations of all possible cell types
vtkSmartPointer< vtkStringArray > GetSpecialVariables() 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< vtkDoubleArray > GetLongitudeCoordinates() const
vtkSmartPointer< vtkDoubleArray > Bounds
vtkSmartPointer< vtkIntArray > GridDimensions
#define VTK_STRUCTURED_GRID
Definition: vtkType.h:93
represent and manipulate 3D points
Definition: vtkPoints.h:39
#define VTK_UNSTRUCTURED_GRID
Definition: vtkType.h:95
vtkSmartPointer< vtkStringArray > GetSpecialVariables() const
vtkSmartPointer< vtkIntArray > GetGridDimensions() const