NEMoSys  0.63.0
A modular, extensible resource with robust automated mesh generation, mesh quality analysis, adaptive mesh refinement, and data transfer between arbitrary meshes.
hdf5Reader Class Reference

Class for reading data from HDF5 files. More...

Detailed Description

Definition at line 60 of file hdf5Reader.H.

Public Member Functions

 hdf5Reader (std::string fname, std::string _outputFname, int _verb=0)
 Construct HDF5 reader object. More...
 
virtual ~hdf5Reader ()
 
int openFile ()
 Open HDF5 file. More...
 
int closeFile ()
 Close HDF5 file. More...
 
int readTopDataSet (std::string dsetName, std::string &buffer)
 Read string data from top level HDF5 DataSet. More...
 
int readDataSet (DataSet dset, std::string &buffer)
 Read string data from generic HDF5 DataSet. More...
 
template<class T >
int readTopDataSet (std::string dsetName, std::vector< T > &buffer)
 Read numeric data from top level HDF5 DataSet. More...
 
template<class T >
int readDataSet (DataSet dset, std::vector< T > &buffer)
 Read numeric data from generic HDF5 DataSet. More...
 
int readGroup (std::string groupName, Group &group)
 Read existing HDF5 Group. More...
 
template<class T >
int readGroupDataSet (std::string dsetName, Group &group, std::vector< T > &buffer)
 Read numeric data from HDF5 group. More...
 
std::string getFileName ()
 Get HDF5 filename. More...
 
template<class T >
void flattenBuffer2d (std::vector< std::vector< T >> originalBuffer, std::vector< T > &flatBuffer)
 Flatten 2d data buffer. More...
 
template<class T >
void unflattenBuffer2d (std::vector< T > flatBuffer, std::vector< std::vector< T >> &originalBuffer)
 Unflatten 2d data buffer. More...
 
void exportToVTKMesh ()
 Export to VTK format. More...
 
vtkSmartPointer< vtkDataSet > getVTKMesh ()
 Get VTK mesh object. More...
 
void exportToMeshBase ()
 Export HDF5 mesh to meshBase object. More...
 
meshBasegetMeshBase ()
 Get meshBase object. More...
 
void writeVTK ()
 Write mesh to VTK file. More...
 
void setFields (meshBase *myMeshBase, std::vector< std::string > dataNames, std::vector< std::vector< double >> data, int pointOrCell)
 Set meshBase field data. More...
 
void setNumberOfVertices (int numVertices)
 Set number of vertices in meshBase. More...
 
void setNumberOfElements (int numElements)
 Set number of elements in meshBase. More...
 
void setNumberOfDimensions (int numDimensions)
 Set number of spatial dimensions in meshBase. More...
 
void setCoordinates (std::vector< double > xCrd)
 Set meshBase 1d coordinate data. More...
 
void setCoordinates (std::vector< double > xCrd, std::vector< double > yCrd)
 Set meshBase 2d coordinate data. More...
 
void setCoordinates (std::vector< double > xCrd, std::vector< double > yCrd, std::vector< double > zCrd)
 Set meshBase 3d coordinate data. More...
 
void setConnectivities (std::vector< std::map< int, std::vector< int >>> elementConnectivity)
 Set meshBase element connectivity table. More...
 
void setElementTypes (std::vector< int > vtkElementTypes)
 Set VTK element types for each element. More...
 
void setElementTypesList (std::vector< int > vtkElementTypesList)
 Set list of each unique type of VTK element in mesh. More...
 
std::vector< int > getElementConnectivity (int elementId)
 Get element connectivity for global element ID. More...
 

Protected Attributes

std::string h5FileName
 HDF5 filename. More...
 
std::string baseH5FileName
 base of HDF5 filename More...
 
std::string outputFname
 VTK output filename. More...
 
int verb
 verbosity flag More...
 
H5File file
 HDF5 file object. More...
 
vtkSmartPointer< vtkDataSet > vtkMesh
 VTK DataSet object. More...
 
meshBasemyMeshBase
 meshBase object More...
 
std::vector< double > xCrd
 
std::vector< double > yCrd
 
std::vector< double > zCrd
 vector of coordinates More...
 
std::vector< std::map< int, std::vector< int > > > elementConnectivity
 element connectivity table More...
 
int numVertices
 number of vertices in mesh More...
 
int numElements
 number of elements in mesh More...
 
int numDimensions
 number of spatial dimensions More...
 
std::vector< int > vtkElementTypes
 vector of VTK element type for each element More...
 
std::vector< int > vtkElementTypesList
 vector storing list of each unique VTK element type More...
 

Inherited by proteusHdf5.

Constructor & Destructor Documentation

◆ hdf5Reader()

hdf5Reader::hdf5Reader ( std::string  fname,
std::string  _outputFname,
int  _verb = 0 
)
inline
Parameters
fnameHDF5 filename
_outputFNameFilename for VTK output file
_verbBoolean verbosity flag

Definition at line 67 of file hdf5Reader.H.

68  : // Initialization of class variables
69  h5FileName(fname),
70  outputFname(_outputFname),
71  verb(_verb),
72  vtkMesh(0) {
73  // Parse filename
74  std::size_t _loc = h5FileName.find_last_of(".");
75  baseH5FileName = h5FileName.substr(0, _loc);
76  };
std::string baseH5FileName
base of HDF5 filename
Definition: hdf5Reader.H:373
vtkSmartPointer< vtkDataSet > vtkMesh
VTK DataSet object.
Definition: hdf5Reader.H:377
std::string h5FileName
HDF5 filename.
Definition: hdf5Reader.H:372
int verb
verbosity flag
Definition: hdf5Reader.H:375
std::string outputFname
VTK output filename.
Definition: hdf5Reader.H:374

◆ ~hdf5Reader()

virtual hdf5Reader::~hdf5Reader ( )
inlinevirtual

Definition at line 79 of file hdf5Reader.H.

79 {};

Member Function Documentation

◆ closeFile()

int hdf5Reader::closeFile ( )
Returns
Integer flag: 0 (failure) or 1 (success)

Definition at line 63 of file hdf5Reader.C.

References file, and h5FileName.

Referenced by proteusHdf5::proteusHdf5().

63  {
64  std::cout << "Closing HDF5 file " << h5FileName << std::endl;
65 
66  try {
67  // Close HDF5 file
68  file.close();
69  } catch (FileIException error) {
70 #if H5_VERS_MAJOR >= 1 && H5_VERS_MINOR >= 8 && H5_VERS_RELEASE >= 20
71  error.printErrorStack();
72 #elif H5_VERS_MAJOR >= 1 && H5_VERS_MINOR >= 10 && H5_VERS_RELEASE >= 2
73  error.printErrorStack();
74 #elif H5_VERS_MAJOR >= 1 && H5_VERS_MINOR >= 11 && H5_VERS_RELEASE >= 0
75  error.printErrorStack();
76 #else
77  error.printError();
78 #endif
79  return -1;
80  }
81  return 0;
82 }
H5File file
HDF5 file object.
Definition: hdf5Reader.H:376
std::string h5FileName
HDF5 filename.
Definition: hdf5Reader.H:372

◆ exportToMeshBase()

void hdf5Reader::exportToMeshBase ( )

Definition at line 313 of file hdf5Reader.C.

References meshBase::Create(), exportToVTKMesh(), myMeshBase, and outputFname.

Referenced by proteusHdf5::proteusHdf5(), and writeVTK().

313  {
314  std::cout << "Exporting hdf5Reader object to meshBase format..." << std::endl;
315  if (!vtkMesh) {
316  exportToVTKMesh();
317  }
319 }
void exportToVTKMesh()
Export to VTK format.
Definition: hdf5Reader.C:241
static meshBase * Create(const std::string &fname)
Construct vtkMesh from filename.
Definition: meshBase.C:78
meshBase * myMeshBase
meshBase object
Definition: hdf5Reader.H:378
std::string outputFname
VTK output filename.
Definition: hdf5Reader.H:374

◆ exportToVTKMesh()

void hdf5Reader::exportToVTKMesh ( )

Definition at line 241 of file hdf5Reader.C.

References getElementConnectivity(), NEM::MSH::New(), numDimensions, numElements, numVertices, points, vtkElementTypes, xCrd, yCrd, and zCrd.

Referenced by exportToMeshBase(), getVTKMesh(), and proteusHdf5::proteusHdf5().

241  {
242  std::cout << "Exporting hdf5Reader object to VTK format..." << std::endl;
243 
244  if (!vtkMesh) {
245  // declare vtk dataset
246  auto dataSet_tmp = vtkSmartPointer<vtkUnstructuredGrid>::New();
247 
248  // points to be pushed into dataSet
250  // set double precision points
251  points->SetDataTypeToDouble();
252 
253  // allocate size for vtk point container
254  points->SetNumberOfPoints(numVertices);
255 
256  if (numDimensions == 2) {
257  for (int i = 0; i < numVertices; ++i) {
258  points->SetPoint(i, xCrd[i], yCrd[i], 0.0);
259  }
260  } else {
261  for (int i = 0; i < numVertices; ++i) {
262  points->SetPoint(i, xCrd[i], yCrd[i], zCrd[i]);
263  }
264  }
265  // add points to vtk mesh data structure
266  std::cout << " setting VTK point data..." << std::endl;
267  dataSet_tmp->SetPoints(points);
268  // allocate space for elements
269  std::cout << " allocating memory for elements..." << std::endl;
270  dataSet_tmp->Allocate(numElements);
271  // add the elements
272  std::cout << " adding VTK elements..." << std::endl;
273  for (int i = 0; i < numElements; ++i) {
274  auto vtkElmIds = vtkSmartPointer<vtkIdList>::New();
275  std::vector<int> elmIds(getElementConnectivity(i));
276  vtkElmIds->SetNumberOfIds(elmIds.size());
277  for (int j = 0; j < elmIds.size(); ++j) {
278  vtkElmIds->SetId(j, elmIds[j] - 0);
279  }
280  switch (vtkElementTypes[i]) {
281 #if VTK_MAJOR_VERSION > 8 || (VTK_MAJOR_VERSION == 8 && VTK_MINOR_VERSION >= 1)
282  case VTK_LAGRANGE_QUADRILATERAL:
283  dataSet_tmp->InsertNextCell(VTK_LAGRANGE_QUADRILATERAL, vtkElmIds);
284  break;
285 #endif
286  case VTK_QUAD:
287  dataSet_tmp->InsertNextCell(VTK_QUAD, vtkElmIds);
288  break;
289  case VTK_TRIANGLE:
290  dataSet_tmp->InsertNextCell(VTK_TRIANGLE, vtkElmIds);
291  break;
292  default:
293  std::cerr << "Unknown element type " << vtkElementTypes[i]
294  << std::endl;
295  break;
296  }
297  }
298  // build links - element incidences
299  dataSet_tmp->BuildLinks();
300  vtkMesh = dataSet_tmp;
301  }
302 }
geoMeshBase * New(MeshType meshType)
Create a new mesh object.
std::vector< double > yCrd
Definition: hdf5Reader.H:381
std::vector< double > xCrd
Definition: hdf5Reader.H:381
int numVertices
number of vertices in mesh
Definition: hdf5Reader.H:384
std::vector< vtkIdType > points
points given by id in .inp file
Definition: inpGeoMesh.C:133
int numDimensions
number of spatial dimensions
Definition: hdf5Reader.H:386
std::vector< double > zCrd
vector of coordinates
Definition: hdf5Reader.H:381
std::vector< int > getElementConnectivity(int elementId)
Get element connectivity for global element ID.
Definition: hdf5Reader.C:195
int numElements
number of elements in mesh
Definition: hdf5Reader.H:385
std::vector< int > vtkElementTypes
vector of VTK element type for each element
Definition: hdf5Reader.H:388

◆ flattenBuffer2d()

template<class T >
void hdf5Reader::flattenBuffer2d ( std::vector< std::vector< T >>  originalBuffer,
std::vector< T > &  flatBuffer 
)
inline
Parameters
originalBufferOriginal unflattened data
flatBufferOutput flattened data

Definition at line 241 of file hdf5Reader.H.

242  {
243  // Populate flat buffer
244  for (auto bufItr1 = originalBuffer.begin(); bufItr1 != originalBuffer.end();
245  bufItr1++) {
246  for (auto bufItr2 = (*bufItr1).begin(); bufItr2 != (*bufItr1).end();
247  bufItr2++) {
248  flatBuffer.push_back(*bufItr2);
249  }
250  }
251  }

◆ getElementConnectivity()

std::vector< int > hdf5Reader::getElementConnectivity ( int  elementId)
Parameters
elementIdElement ID
Returns
Element connectivity table

Definition at line 195 of file hdf5Reader.C.

References elementConnectivity, vtkElementTypes, and vtkElementTypesList.

Referenced by exportToVTKMesh().

195  {
196  int nVerts;
197  int connId = -1;
198 #if VTK_MAJOR_VERSION > 8 || (VTK_MAJOR_VERSION == 8 && VTK_MINOR_VERSION >= 1)
199  if (vtkElementTypes[elementId] == VTK_LAGRANGE_QUADRILATERAL) {
200  nVerts = 16;
201  auto itr = std::find(vtkElementTypesList.begin(), vtkElementTypesList.end(),
202  VTK_LAGRANGE_QUADRILATERAL);
203  if (itr != vtkElementTypesList.end()) {
204  connId = distance(vtkElementTypesList.begin(), itr);
205  }
206  } else if (vtkElementTypes[elementId] == VTK_QUAD)
207 #else
208  if (vtkElementTypes[elementId] == VTK_QUAD)
209 #endif
210  {
211  nVerts = 4;
212  auto itr = std::find(vtkElementTypesList.begin(), vtkElementTypesList.end(),
213  VTK_QUAD);
214  if (itr != vtkElementTypesList.end()) {
215  connId = distance(vtkElementTypesList.begin(), itr);
216  }
217  } else if (vtkElementTypes[elementId] == VTK_TRIANGLE) {
218  nVerts = 3;
219  auto itr = std::find(vtkElementTypesList.begin(), vtkElementTypesList.end(),
220  VTK_TRIANGLE);
221  if (itr != vtkElementTypesList.end()) {
222  connId = distance(vtkElementTypesList.begin(), itr);
223  }
224  } else {
225  std::cerr << "Element type not supported: " << vtkElementTypes[elementId]
226  << std::endl;
227  exit(-1);
228  }
229  if (connId != -1) {
230  std::vector<int> subConnectivity(
231  elementConnectivity[connId][elementId].begin(),
232  elementConnectivity[connId][elementId].begin() + nVerts);
233  return subConnectivity;
234  } else {
235  std::cerr << "Element type " << vtkElementTypes[elementId]
236  << " not found in connectivity table" << std::endl;
237  exit(-1);
238  }
239 }
std::vector< int > vtkElementTypesList
vector storing list of each unique VTK element type
Definition: hdf5Reader.H:389
std::vector< std::map< int, std::vector< int > > > elementConnectivity
element connectivity table
Definition: hdf5Reader.H:383
std::vector< int > vtkElementTypes
vector of VTK element type for each element
Definition: hdf5Reader.H:388

◆ getFileName()

std::string hdf5Reader::getFileName ( )
Returns
HDF5 filename

Definition at line 39 of file hdf5Reader.C.

References h5FileName.

39 { return h5FileName; }
std::string h5FileName
HDF5 filename.
Definition: hdf5Reader.H:372

◆ getMeshBase()

meshBase * hdf5Reader::getMeshBase ( )
Returns
meshBase object

Definition at line 321 of file hdf5Reader.C.

References myMeshBase.

Referenced by proteusHdf5::proteusHdf5().

321 { return myMeshBase; }
meshBase * myMeshBase
meshBase object
Definition: hdf5Reader.H:378

◆ getVTKMesh()

vtkSmartPointer< vtkDataSet > hdf5Reader::getVTKMesh ( )
Returns
Pointer to vtkDataSet

Definition at line 304 of file hdf5Reader.C.

References exportToVTKMesh(), and vtkMesh.

304  {
305  if (vtkMesh)
306  return vtkMesh;
307  else {
308  exportToVTKMesh();
309  return vtkMesh;
310  }
311 }
void exportToVTKMesh()
Export to VTK format.
Definition: hdf5Reader.C:241
vtkSmartPointer< vtkDataSet > vtkMesh
VTK DataSet object.
Definition: hdf5Reader.H:377

◆ openFile()

int hdf5Reader::openFile ( )
Returns
Integer flag: 0 (failure) or 1 (success)

Definition at line 41 of file hdf5Reader.C.

References file, and h5FileName.

Referenced by proteusHdf5::proteusHdf5().

41  {
42  std::cout << "Opening HDF5 file " << h5FileName << std::endl;
43 
44  try {
45  // Open HDF5 file
46  H5std_string h5FileName_(h5FileName);
47  file.openFile(h5FileName_, H5F_ACC_RDWR);
48  } catch (FileIException error) {
49 #if H5_VERS_MAJOR >= 1 && H5_VERS_MINOR >= 8 && H5_VERS_RELEASE >= 20
50  error.printErrorStack();
51 #elif H5_VERS_MAJOR >= 1 && H5_VERS_MINOR >= 10 && H5_VERS_RELEASE >= 2
52  error.printErrorStack();
53 #elif H5_VERS_MAJOR >= 1 && H5_VERS_MINOR >= 11 && H5_VERS_RELEASE >= 0
54  error.printErrorStack();
55 #else
56  error.printError();
57 #endif
58  return -1;
59  }
60  return 0;
61 }
H5File file
HDF5 file object.
Definition: hdf5Reader.H:376
std::string h5FileName
HDF5 filename.
Definition: hdf5Reader.H:372

◆ readDataSet() [1/2]

int hdf5Reader::readDataSet ( DataSet  dset,
std::string &  buffer 
)
Parameters
dsetHDF5 dataset
bufferBuffer into which data is read
Returns
Integer flag: 0 (failure) or 1 (success)

Definition at line 109 of file hdf5Reader.C.

Referenced by readTopDataSet().

109  {
110  try {
111  // Get type of data
112  DataType dtype = dset.getDataType();
113  // Read data
114  dset.read(buffer, dtype);
115  // Close dataset
116  dset.close();
117  }
118  // catch failure caused by the H5File operations
119  catch (DataSetIException error) {
120 #if H5_VERS_MAJOR >= 1 && H5_VERS_MINOR >= 8 && H5_VERS_RELEASE >= 20
121  error.printErrorStack();
122 #elif H5_VERS_MAJOR >= 1 && H5_VERS_MINOR >= 10 && H5_VERS_RELEASE >= 2
123  error.printErrorStack();
124 #elif H5_VERS_MAJOR >= 1 && H5_VERS_MINOR >= 11 && H5_VERS_RELEASE >= 0
125  error.printErrorStack();
126 #else
127  error.printError();
128 #endif
129  return -1;
130  }
131  return 0;
132 }

◆ readDataSet() [2/2]

template<class T >
int hdf5Reader::readDataSet ( DataSet  dset,
std::vector< T > &  buffer 
)
inline
Parameters
dsetHDF5 dataset
bufferBuffer into which numeric data is read
Returns
Integer flag: 0 (failure) or 1 (success)

Definition at line 148 of file hdf5Reader.H.

148  {
149  try {
150  // Get type of data
151  DataType dtype = dset.getDataType();
152  // Read data
153  dset.read(&buffer[0], dtype);
154  // Close dataset
155  dset.close();
156  }
157  // catch failure caused by the H5File operations
158  catch (DataSetIException error) {
159 #if H5_VERS_MAJOR >= 1 && H5_VERS_MINOR >= 8 && H5_VERS_RELEASE >= 20
160  error.printErrorStack();
161 #elif H5_VERS_MAJOR >= 1 && H5_VERS_MINOR >= 10 && H5_VERS_RELEASE >= 2
162  error.printErrorStack();
163 #elif H5_VERS_MAJOR >= 1 && H5_VERS_MINOR >= 11 && H5_VERS_RELEASE >= 0
164  error.printErrorStack();
165 #else
166  error.printError();
167 #endif
168  return -1;
169  }
170 
171  return 0;
172  }

◆ readGroup()

int hdf5Reader::readGroup ( std::string  groupName,
Group &  group 
)
Parameters
groupNameHDF5 group object name
groupHDF5 group object
Returns
Integer flag: 0 (failure) or 1 (success)

Definition at line 135 of file hdf5Reader.C.

References file.

Referenced by proteusHdf5::getBlocks().

135  {
136  try {
137  // Open group
138  group = Group(file.openGroup(groupName));
139  } catch (FileIException error) {
140 #if H5_VERS_MAJOR >= 1 && H5_VERS_MINOR >= 8 && H5_VERS_RELEASE >= 20
141  error.printErrorStack();
142 #elif H5_VERS_MAJOR >= 1 && H5_VERS_MINOR >= 10 && H5_VERS_RELEASE >= 2
143  error.printErrorStack();
144 #elif H5_VERS_MAJOR >= 1 && H5_VERS_MINOR >= 11 && H5_VERS_RELEASE >= 0
145  error.printErrorStack();
146 #else
147  error.printError();
148 #endif
149  return -1;
150  }
151  return 0;
152 }
H5File file
HDF5 file object.
Definition: hdf5Reader.H:376

◆ readGroupDataSet()

template<class T >
int hdf5Reader::readGroupDataSet ( std::string  dsetName,
Group &  group,
std::vector< T > &  buffer 
)
inline
Parameters
dsetNameDataSet name
groupHDF5 group object
bufferBuffer into which numeric data is read
Returns
Integer flag: 0 (failure) or 1 (success)

Definition at line 190 of file hdf5Reader.H.

Referenced by proteusHdf5::getBlockElementData(), proteusHdf5::getBlockGlobalID(), proteusHdf5::getBlockInfo(), proteusHdf5::getBlockVertexData(), and proteusHdf5::getBlockXYZ().

191  {
192  try {
193  // Open dataset
194  DataSet dset = group.openDataSet(dsetName);
195  // Get type of data
196  DataType dtype = dset.getDataType();
197  // Read data
198  dset.read(&buffer[0], dtype);
199  // Close dataset
200  dset.close();
201  }
202  // catch failure caused by the H5File operations
203  catch (GroupIException error) {
204 #if H5_VERS_MAJOR >= 1 && H5_VERS_MINOR >= 8 && H5_VERS_RELEASE >= 20
205  error.printErrorStack();
206 #elif H5_VERS_MAJOR >= 1 && H5_VERS_MINOR >= 10 && H5_VERS_RELEASE >= 2
207  error.printErrorStack();
208 #elif H5_VERS_MAJOR >= 1 && H5_VERS_MINOR >= 11 && H5_VERS_RELEASE >= 0
209  error.printErrorStack();
210 #else
211  error.printError();
212 #endif
213  return -1;
214  } catch (DataSetIException error) {
215 #if H5_VERS_MAJOR >= 1 && H5_VERS_MINOR >= 8 && H5_VERS_RELEASE >= 20
216  error.printErrorStack();
217 #elif H5_VERS_MAJOR >= 1 && H5_VERS_MINOR >= 10 && H5_VERS_RELEASE >= 2
218  error.printErrorStack();
219 #elif H5_VERS_MAJOR >= 1 && H5_VERS_MINOR >= 11 && H5_VERS_RELEASE >= 0
220  error.printErrorStack();
221 #else
222  error.printError();
223 #endif
224  return -1;
225  }
226 
227  return 0;
228  }

◆ readTopDataSet() [1/2]

int hdf5Reader::readTopDataSet ( std::string  dsetName,
std::string &  buffer 
)
Parameters
dsetNameDataset name
bufferBuffer into which string data is read
Returns
Integer flag: 0 (failure) or 1 (success)

Definition at line 85 of file hdf5Reader.C.

References file, and readDataSet().

Referenced by proteusHdf5::getControlInformation(), and proteusHdf5::getVectorNames().

85  {
86  try {
87  // Open dataset
88  DataSet dset = file.openDataSet(dsetName);
89  // Read dataset
90  readDataSet(dset, buffer);
91  }
92  // catch failure caused by the H5File operations
93  catch (FileIException error) {
94 #if H5_VERS_MAJOR >= 1 && H5_VERS_MINOR >= 8 && H5_VERS_RELEASE >= 20
95  error.printErrorStack();
96 #elif H5_VERS_MAJOR >= 1 && H5_VERS_MINOR >= 10 && H5_VERS_RELEASE >= 2
97  error.printErrorStack();
98 #elif H5_VERS_MAJOR >= 1 && H5_VERS_MINOR >= 11 && H5_VERS_RELEASE >= 0
99  error.printErrorStack();
100 #else
101  error.printError();
102 #endif
103  return -1;
104  }
105  return 0;
106 }
H5File file
HDF5 file object.
Definition: hdf5Reader.H:376
int readDataSet(DataSet dset, std::string &buffer)
Read string data from generic HDF5 DataSet.
Definition: hdf5Reader.C:109

◆ readTopDataSet() [2/2]

template<class T >
int hdf5Reader::readTopDataSet ( std::string  dsetName,
std::vector< T > &  buffer 
)
inline
Parameters
dsetNameDataset name
bufferBuffer into which numeric data is read
Returns
Integer flag: 0 (failure) or 1 (success)

Definition at line 115 of file hdf5Reader.H.

115  {
116  try {
117  std::cout << "H5_VERS_MAJOR = " << H5_VERS_MAJOR << std::endl;
118  std::cout << "H5_VERS_MINOR = " << H5_VERS_MINOR << std::endl;
119  std::cout << "H5_VERS_RELEASE = " << H5_VERS_RELEASE << std::endl;
120  std::cout << "H5_VERS_SUBRELEASE = " << H5_VERS_SUBRELEASE << std::endl;
121  // Open dataset
122  DataSet dset = file.openDataSet(dsetName);
123  readDataSet(dset, buffer);
124  }
125  // catch failure caused by the H5File operations
126  catch (FileIException error) {
127 #if H5_VERS_MAJOR >= 1 && H5_VERS_MINOR >= 8 && H5_VERS_RELEASE >= 20
128  error.printErrorStack();
129 #elif H5_VERS_MAJOR >= 1 && H5_VERS_MINOR >= 10 && H5_VERS_RELEASE >= 2
130  error.printErrorStack();
131 #elif H5_VERS_MAJOR >= 1 && H5_VERS_MINOR >= 11 && H5_VERS_RELEASE >= 0
132  error.printErrorStack();
133 #else
134  error.printError();
135 #endif
136  return -1;
137  }
138 
139  return 0;
140  }
H5File file
HDF5 file object.
Definition: hdf5Reader.H:376
int readDataSet(DataSet dset, std::string &buffer)
Read string data from generic HDF5 DataSet.
Definition: hdf5Reader.C:109

◆ setConnectivities()

void hdf5Reader::setConnectivities ( std::vector< std::map< int, std::vector< int >>>  elementConnectivity)
Parameters
elementConnectivityElement connectivity table

Definition at line 190 of file hdf5Reader.C.

References elementConnectivity.

Referenced by proteusHdf5::proteusHdf5().

191  {
192  elementConnectivity = _elementConnectivity;
193 }

◆ setCoordinates() [1/3]

void hdf5Reader::setCoordinates ( std::vector< double >  xCrd)
Parameters
xCrdvector of x-coordinates

Definition at line 166 of file hdf5Reader.C.

References xCrd.

Referenced by proteusHdf5::proteusHdf5().

166 { xCrd = _xCrd; }

◆ setCoordinates() [2/3]

void hdf5Reader::setCoordinates ( std::vector< double >  xCrd,
std::vector< double >  yCrd 
)
Parameters
xCrdvector of x-coordinates
yCrdvector of y-coordinates

Definition at line 168 of file hdf5Reader.C.

References xCrd, and yCrd.

169  {
170  xCrd = _xCrd;
171  yCrd = _yCrd;
172 }

◆ setCoordinates() [3/3]

void hdf5Reader::setCoordinates ( std::vector< double >  xCrd,
std::vector< double >  yCrd,
std::vector< double >  zCrd 
)
Parameters
xCrdvector of x-coordinates
yCrdvector of y-coordinates
zCrdvector of z-coordinates

Definition at line 174 of file hdf5Reader.C.

References xCrd, yCrd, and zCrd.

176  {
177  xCrd = _xCrd;
178  yCrd = _yCrd;
179  zCrd = _zCrd;
180 }

◆ setElementTypes()

void hdf5Reader::setElementTypes ( std::vector< int >  vtkElementTypes)
Parameters
vtkElementTypesVector of VTK element types for each element

Definition at line 182 of file hdf5Reader.C.

References vtkElementTypes.

Referenced by proteusHdf5::proteusHdf5().

182  {
183  vtkElementTypes = _vtkElementTypes;
184 }

◆ setElementTypesList()

void hdf5Reader::setElementTypesList ( std::vector< int >  vtkElementTypesList)
Parameters
vtkElementTypesListVector of each unique VTK element type in mesh

Definition at line 186 of file hdf5Reader.C.

References vtkElementTypesList.

Referenced by proteusHdf5::proteusHdf5().

186  {
187  vtkElementTypesList = _vtkElementTypesList;
188 }

◆ setFields()

void hdf5Reader::setFields ( meshBase myMeshBase,
std::vector< std::string >  dataNames,
std::vector< std::vector< double >>  data,
int  pointOrCell 
)
Parameters
myMeshBasemeshBase object
dataNamesVector of names of data fields
dataVector storing data
pointOrCellBoolean indicating point or cell data

Definition at line 323 of file hdf5Reader.C.

References data, meshBase::setCellDataArray(), and meshBase::setPointDataArray().

Referenced by proteusHdf5::proteusHdf5().

326  {
327  if (pointOrCell == 0) {
328  std::cout << "Setting hdf5Reader object point fields..." << std::endl;
329  } else if (pointOrCell == 1) {
330  std::cout << "Setting hdf5Reader object cell fields..." << std::endl;
331  }
332 
333  // Set element data
334  auto dataNameItr = dataNames.begin();
335  auto dataItr = data.begin();
336  while (dataNameItr != dataNames.end() || dataItr != data.end()) {
337  std::cout << " " << *dataNameItr << std::endl;
338 
339  if (pointOrCell == 0) {
340  myMeshBase->setPointDataArray((*dataNameItr).c_str(), *dataItr);
341  } else if (pointOrCell == 1) {
342  myMeshBase->setCellDataArray((*dataNameItr).c_str(), *dataItr);
343  }
344 
345  if (dataNameItr != dataNames.end()) {
346  dataNameItr++;
347  }
348  if (dataItr != data.end()) {
349  dataItr++;
350  }
351  }
352 }
virtual void setPointDataArray(const std::string &name, const std::vector< std::vector< double >> &data)
register data to dataSet&#39;s point data
Definition: meshBase.H:319
virtual void setCellDataArray(const std::string &name, const std::vector< std::vector< double >> &data)
register data to dataSet&#39;s cell data
Definition: meshBase.H:334

◆ setNumberOfDimensions()

void hdf5Reader::setNumberOfDimensions ( int  numDimensions)
Parameters
numDimensionsNumber of spatial dimensions

Definition at line 162 of file hdf5Reader.C.

References numDimensions.

Referenced by proteusHdf5::proteusHdf5().

162  {
163  numDimensions = _numDimensions;
164 }
int numDimensions
number of spatial dimensions
Definition: hdf5Reader.H:386

◆ setNumberOfElements()

void hdf5Reader::setNumberOfElements ( int  numElements)
Parameters
numElementsNumber of elements

Definition at line 158 of file hdf5Reader.C.

References numElements.

Referenced by proteusHdf5::proteusHdf5().

158  {
159  numElements = _numElements;
160 }
int numElements
number of elements in mesh
Definition: hdf5Reader.H:385

◆ setNumberOfVertices()

void hdf5Reader::setNumberOfVertices ( int  numVertices)
Parameters
numVerticesNumber of vertices

Definition at line 154 of file hdf5Reader.C.

References numVertices.

Referenced by proteusHdf5::proteusHdf5().

154  {
155  numVertices = _numVertices;
156 }
int numVertices
number of vertices in mesh
Definition: hdf5Reader.H:384

◆ unflattenBuffer2d()

template<class T >
void hdf5Reader::unflattenBuffer2d ( std::vector< T >  flatBuffer,
std::vector< std::vector< T >> &  originalBuffer 
)
inline
Parameters
flatBufferOriginal flattened data
originalBufferOutput unflattened data

Definition at line 258 of file hdf5Reader.H.

References data.

Referenced by proteusHdf5::getBlockElementData(), and proteusHdf5::getBlockVertexData().

259  {
260  // Populate original buffer
261  auto flatItr = flatBuffer.begin();
262  for (auto bufItr1 = originalBuffer.begin(); bufItr1 != originalBuffer.end();
263  bufItr1++) {
264  for (auto bufItr2 = (*bufItr1).begin(); bufItr2 != (*bufItr1).end();
265  bufItr2++) {
266  *bufItr2 = *flatItr;
267  flatItr++;
268  }
269  }
270  }

◆ writeVTK()

void hdf5Reader::writeVTK ( )

Definition at line 354 of file hdf5Reader.C.

References exportToMeshBase(), myMeshBase, and meshBase::write().

354  {
355  std::cout << "Writing hdf5Reader object to VTK file" << std::endl;
356  if (!myMeshBase) {
358  }
359  myMeshBase->write();
360 }
virtual void write() const
write the mesh to file named after the private var &#39;filename&#39;.
Definition: meshBase.H:598
void exportToMeshBase()
Export HDF5 mesh to meshBase object.
Definition: hdf5Reader.C:313
meshBase * myMeshBase
meshBase object
Definition: hdf5Reader.H:378

Member Data Documentation

◆ baseH5FileName

std::string hdf5Reader::baseH5FileName
protected

Definition at line 373 of file hdf5Reader.H.

◆ elementConnectivity

std::vector<std::map<int, std::vector<int> > > hdf5Reader::elementConnectivity
protected

◆ file

H5File hdf5Reader::file
protected

Definition at line 376 of file hdf5Reader.H.

Referenced by closeFile(), openFile(), readGroup(), and readTopDataSet().

◆ h5FileName

std::string hdf5Reader::h5FileName
protected

Definition at line 372 of file hdf5Reader.H.

Referenced by closeFile(), getFileName(), and openFile().

◆ myMeshBase

meshBase* hdf5Reader::myMeshBase
protected

Definition at line 378 of file hdf5Reader.H.

Referenced by exportToMeshBase(), getMeshBase(), proteusHdf5::proteusHdf5(), and writeVTK().

◆ numDimensions

int hdf5Reader::numDimensions
protected

Definition at line 386 of file hdf5Reader.H.

Referenced by exportToVTKMesh(), and setNumberOfDimensions().

◆ numElements

int hdf5Reader::numElements
protected

Definition at line 385 of file hdf5Reader.H.

Referenced by exportToVTKMesh(), and setNumberOfElements().

◆ numVertices

int hdf5Reader::numVertices
protected

Definition at line 384 of file hdf5Reader.H.

Referenced by exportToVTKMesh(), and setNumberOfVertices().

◆ outputFname

std::string hdf5Reader::outputFname
protected

Definition at line 374 of file hdf5Reader.H.

Referenced by exportToMeshBase().

◆ verb

int hdf5Reader::verb
protected

Definition at line 375 of file hdf5Reader.H.

◆ vtkElementTypes

std::vector<int> hdf5Reader::vtkElementTypes
protected

Definition at line 388 of file hdf5Reader.H.

Referenced by exportToVTKMesh(), getElementConnectivity(), and setElementTypes().

◆ vtkElementTypesList

std::vector<int> hdf5Reader::vtkElementTypesList
protected

Definition at line 389 of file hdf5Reader.H.

Referenced by getElementConnectivity(), and setElementTypesList().

◆ vtkMesh

vtkSmartPointer<vtkDataSet> hdf5Reader::vtkMesh
protected

Definition at line 377 of file hdf5Reader.H.

Referenced by getVTKMesh().

◆ xCrd

std::vector<double> hdf5Reader::xCrd
protected

Definition at line 381 of file hdf5Reader.H.

Referenced by exportToVTKMesh(), and setCoordinates().

◆ yCrd

std::vector<double> hdf5Reader::yCrd
protected

Definition at line 381 of file hdf5Reader.H.

Referenced by exportToVTKMesh(), and setCoordinates().

◆ zCrd

std::vector<double> hdf5Reader::zCrd
protected

Definition at line 381 of file hdf5Reader.H.

Referenced by exportToVTKMesh(), and setCoordinates().


The documentation for this class was generated from the following files: