kspaceFirstOrder3D-CUDA  1.1
The CUDA/C++ implementation of the k-wave toolbox for the time-domain simulation of acoustic wave fields in 3D
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
WholeDomainOutputHDF5Stream.cpp
Go to the documentation of this file.
1 /**
2  * @file WholeDomainOutputHDF5Stream.cpp
3  *
4  * @author Jiri Jaros \n
5  * Faculty of Information Technology \n
6  * Brno University of Technology \n
7  * jarosjir@fit.vutbr.cz
8  *
9  * @brief The implementation file of the class saving RealMatrix data into the output
10  * HDF5 file, e.g. p_max_all.
11  *
12  * @version kspaceFirstOrder3D 3.4
13  *
14  * @date 28 August 2014, 11:15 (created)
15  * 29 July 2016, 16:56 (revised)
16  *
17  * @section License
18  * This file is part of the C++ extension of the k-Wave Toolbox
19  * (http://www.k-wave.org).\n Copyright (C) 2016 Jiri Jaros and Bradley Treeby.
20  *
21  * This file is part of the k-Wave. k-Wave is free software: you can redistribute it and/or modify
22  * it under the terms of the GNU Lesser General Public License as published by the Free Software
23  * Foundation, either version 3 of the License, or (at your option) any later version.
24  *
25  * k-Wave is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even
26  * the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser
27  * General Public License for more details.
28  *
29  * You should have received a copy of the GNU Lesser General Public License along with k-Wave.
30  * If not, see http://www.gnu.org/licenses/.
31  */
32 
34 #include <Parameters/Parameters.h>
35 
37 
38 //------------------------------------------------------------------------------------------------//
39 //------------------------------------------ Constants -------------------------------------------//
40 //------------------------------------------------------------------------------------------------//
41 
42 //------------------------------------------------------------------------------------------------//
43 //--------------------------------------- Public methods -----------------------------------------//
44 //------------------------------------------------------------------------------------------------//
45 
46 /**
47  * Constructor - links the HDF5 dataset and SourceMatrix
48  * @param [in] file - HDF5 file to write the output to
49  * @param [in] datasetName - The name of the HDF5 group containing datasets for particular cuboids
50  * @param [in] sourceMatrix - Source matrix to be sampled
51  * @param [in] reduceOp - Reduce operator
52  */
54  TMatrixName& datasetName,
55  TRealMatrix& sourceMatrix,
56  const TReduceOperator reduceOp)
57  : TBaseOutputHDF5Stream(file, datasetName, sourceMatrix, reduceOp),
58  dataset(H5I_BADID),
59  sampledTimeStep(0)
60 {
61 
62 }// end of TWholeDomainOutputHDF5Stream
63 //--------------------------------------------------------------------------------------------------
64 
65 
66 /**
67  * Destructor.
68  * If the file is still opened, it applies the post processing and flush the data.
69  * Then, the object memory is freed and the object destroyed.
70  */
72 {
73  Close();
74  // free memory
75  FreeMemory();
76 }// end of Destructor
77 //--------------------------------------------------------------------------------------------------
78 
79 
80 
81 /**
82  * Create a HDF5 stream for the whole domain and allocate data for it.
83  */
85 {
88  1);
89 
90  // Create a dataset under the root group
94  chunkSize,
96 
97  // Write dataset parameters
100 
101  // Set buffer size
103 
104  // Allocate memory
105  AllocateMemory();
106 }//end of Create
107 //--------------------------------------------------------------------------------------------------
108 
109 
110 /**
111  * Reopen the output stream after restart and reload data
112  */
114 {
115  const TParameters& params = TParameters::GetInstance();
116 
117  // Set buffer size
119 
120  // Allocate memory
121  AllocateMemory();
122 
123  // Open the dataset under the root group
125 
126  sampledTimeStep = 0;
127  if (reduceOp == NONE)
128  { // seek in the dataset
129  sampledTimeStep = (params.Get_t_index() < params.GetStartTimeIndex()) ?
130  0 : (params.Get_t_index() - params.GetStartTimeIndex());
131  }
132  else
133  { // reload data
134  // Read data from disk only if there were anything stored there (t_index > start_index)
135  //(one step ahead)
136  if (params.Get_t_index() > params.GetStartTimeIndex())
137  {
141  hostBuffer);
142  // Send data to device
144  }
145  }
146 }// end of Reopen
147 //--------------------------------------------------------------------------------------------------
148 
149 
150 /**
151  * Sample all grid points, line them up in the buffer an flush to the disk unless a reduce operator
152  * is applied
153  */
155 {
156  switch (reduceOp)
157  {
158  case NONE :
159  {
160  // Copy all data from GPU to CPU (no need to use a kernel)
161  // this violates the const prerequisite, however this routine is still NOT used in the code
162  const_cast<TRealMatrix&> (sourceMatrix).CopyFromDevice();
163 
164  // We use here direct HDF5 offload using MEMSPACE - seems to be faster for bigger datasets
165  const TDimensionSizes datasetPosition(0, 0, 0, sampledTimeStep); //4D position in the dataset
166 
167  TDimensionSizes cuboidSize(sourceMatrix.GetDimensionSizes());// Size of the cuboid
168  cuboidSize.nt = 1;
169 
170  // iterate over all cuboid to be sampled
172  datasetPosition,
173  TDimensionSizes(0,0,0,0), // position in the SourceMatrix
174  cuboidSize,
177 
178  sampledTimeStep++; // Move forward in time
179 
180  break;
181  }// case NONE
182 
183  case RMS :
184  {
185  OutputStreamsCUDAKernels::SampleAll<RMS>
186  (deviceBuffer,
189  break;
190  }// case RMS
191 
192  case MAX :
193  {
194  OutputStreamsCUDAKernels::SampleAll<MAX>
195  (deviceBuffer,
198  break;
199  }//case MAX
200 
201  case MIN :
202  {
203  OutputStreamsCUDAKernels::SampleAll<MIN>
204  (deviceBuffer,
207  break;
208  } //case MIN
209  }// switch
210 }// end of Sample
211 //--------------------------------------------------------------------------------------------------
212 
213 
214 /**
215  * Apply post-processing on the buffer and flush it to the file.
216  */
218 {
219  // run inherited method
221 
222  // When no reduce operator is applied, the data is flushed after every time step
223  // which means it has been done before
224  if (reduceOp != NONE)
225  {
226  // Copy data from GPU matrix
228 
230  }
231 }// end of PostProcessing
232 //--------------------------------------------------------------------------------------------------
233 
234 /**
235  * Checkpoint the stream.
236  */
238 {
239  // copy data on the CPU
241 
242  // raw data has already been flushed, others has to be flushed here
243  if (reduceOp != NONE) FlushBufferToFile();
244 }// end of Checkpoint
245 //-------------------------------------------------------------------------------------------------
246 
247 /**
248  * Close stream (apply post-processing if necessary, flush data and close).
249  */
251 {
252  // the dataset is still opened
253  if (dataset != H5I_BADID)
254  {
256  }
257 
258  dataset = H5I_BADID;
259 }// end of Close
260 //--------------------------------------------------------------------------------------------------
261 
262 
263 //------------------------------------------------------------------------------------------------//
264 //-------------------------------------- Protected methods ---------------------------------------//
265 //------------------------------------------------------------------------------------------------//
266 
267 
268 /**
269  * Flush the buffer down to the file at the actual position.
270  */
272 {
274  TDimensionSizes position(0,0,0);
275 
276  // Not used for NONE now!
277  if (reduceOp == NONE)
278  {
279  position.nt = sampledTimeStep;
280  size.nt = sampledTimeStep;
281  }
282 
283  file.WriteHyperSlab(dataset, position, size, hostBuffer);
284  sampledTimeStep++;
285 }// end of FlushToFile
286 //--------------------------------------------------------------------------------------------------
287 
288 //------------------------------------------------------------------------------------------------//
289 //--------------------------------------- Private methods ----------------------------------------//
290 //------------------------------------------------------------------------------------------------//
291 
size_t nx
number of elements in the x direction
virtual void FlushBufferToFile()
Flush the buffer to the file.
THDF5_File & file
HDF5 file handle.
hid_t CreateFloatDataset(const hid_t parentGroup, TMatrixName &datasetName, const TDimensionSizes &dimensionSizes, const TDimensionSizes &chunkSizes, const size_t compressionLevel)
Create a float HDF5 dataset at a specified place in the file tree (3D/4D).
Definition: HDF5_File.cpp:296
size_t GetStartTimeIndex() const
Get start time index for sensor recording.
Definition: Parameters.h:225
void WriteMatrixDomainType(const hid_t parentGroup, TMatrixName &datasetName, const TMatrixDomainType &matrixDomainType)
Write matrix domain type into the dataset under the root group.
Definition: HDF5_File.cpp:1048
virtual void Reopen()
Reopen the output stream after restart and reload data.
size_t nt
Number of time steps (for time series datasets).
virtual void CopyDataFromDevice()
Copy data deviceBuffer -> hostBuffer.
calculate root mean square
std::string rootObjectName
Dataset name.
static TParameters & GetInstance()
Get instance of the singleton class.
Definition: Parameters.cpp:70
virtual void Close()
Close stream (apply post-processing if necessary, flush data and close).
virtual size_t GetElementCount() const
Get element count of the matrix.
const TReduceOperator reduceOp
Reduce operator.
void WriteCuboidToHyperSlab(const hid_t dataset, const TDimensionSizes &hyperslabPosition, const TDimensionSizes &cuboidPosition, const TDimensionSizes &cuboidSize, const TDimensionSizes &matrixDimensions, const float *mMatrixData)
Write a cuboid selected within the matrixData into a hyperslab.
Definition: HDF5_File.cpp:618
virtual void Sample()
Sample data (copy from GPU memory and then flush - no overlapping implemented!)
The header file containing the parameters of the simulation.
virtual float * GetHostData()
Get raw CPU data out of the class (for direct CPU kernel access).
The header file of cuda kernels used for data sampling (output streams).
const std::string TMatrixName
Datatype for matrix names.
Definition: MatrixNames.h:45
virtual void PostProcess()
Apply post-processing on the buffer and flush it to the file.
virtual void AllocateMemory()
A generic function to allocate memory - not used in the base class.
virtual TDimensionSizes GetDimensionSizes() const
Get dimension sizes of the matrix.
void WriteMatrixDataType(const hid_t parentGroup, TMatrixName &datasetName, const TMatrixDataType &matrixDataType)
Write matrix data type into the dataset under a specified group.
Definition: HDF5_File.cpp:1029
const TRealMatrix & sourceMatrix
Source matrix to be sampled.
size_t Get_t_index() const
Get simulation time step.
Definition: Parameters.h:102
virtual ~TWholeDomainOutputHDF5Stream()
Destructor.
Class storing all parameters of the simulation.
Definition: Parameters.h:49
TWholeDomainOutputHDF5Stream(THDF5_File &file, TMatrixName &datasetName, TRealMatrix &sourceMatrix, const TReduceOperator reduceOp)
Constructor.
size_t sampledTimeStep
Time step to store (N/A for aggregated).
hid_t GetRootGroup() const
Get handle to the root group.
Definition: HDF5_File.h:573
The header file of the class saving whole RealMatrix into the output HDF5 file, e.g. p_max_all.
TReduceOperator
How to aggregate data.
float * deviceBuffer
Temporary buffer on the GPU side - only for aggregated quantities.
store actual data (time series)
size_t ny
number of elements in the y direction
Abstract base class for output data streams (sampled data).
virtual void FreeMemory()
A generic function to free memory - not used in the base class.
The class for real matrices.
Definition: RealMatrix.h:45
float * hostBuffer
Temporary buffer for store on the GPU side.
virtual void CopyDataToDevice()
Copy data hostBuffer-> deviceBuffer.
virtual void PostProcess()
Apply post-processing on the buffer and flush it to the file.
TDimensionSizes position
Position in the dataset.
hid_t dataset
Handle to a HDF5 dataset.
size_t GetCompressionLevel() const
Get compression level.
Definition: Parameters.h:218
void WriteHyperSlab(const hid_t dataset, const TDimensionSizes &position, const TDimensionSizes &size, const float *data)
Write a hyper-slab into the dataset - float dataset.
Definition: HDF5_File.cpp:473
virtual float * GetDeviceData()
Get raw GPU data out of the class (for direct GPU kernel access).
size_t bufferSize
Buffer size.
virtual void Create()
Create a HDF5 stream and allocate data for it.
void ReadCompleteDataset(const hid_t parentGroup, TMatrixName &datasetName, const TDimensionSizes &dimensionSizes, float *data)
Read data from the dataset under a specified group, float dataset.
Definition: HDF5_File.cpp:887
void CloseDataset(const hid_t dataset)
Close the HDF5 dataset.
Definition: HDF5_File.cpp:457
Class wrapping the HDF5 routines.
Definition: HDF5_File.h:500
Structure with 4D dimension sizes (3 in space and 1 in time).
hid_t OpenDataset(const hid_t parentGroup, TMatrixName &datasetName)
Open the HDF5 dataset at a specified place in the file tree.
Definition: HDF5_File.cpp:268