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
IndexOutputHDF5Stream.cpp
Go to the documentation of this file.
1 /**
2  * @file IndexOutputHDF5Stream.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 data based on index senor mask into
10  * the output HDF5 file.
11  *
12  * @version kspaceFirstOrder3D 3.4
13  *
14  * @date 29 August 2014, 10:10 (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 
35 
36 #include <Parameters/Parameters.h>
37 #include <Logger/Logger.h>
38 
39 //------------------------------------------------------------------------------------------------//
40 //------------------------------------------ Constants -------------------------------------------//
41 //------------------------------------------------------------------------------------------------//
42 
43 //------------------------------------------------------------------------------------------------//
44 //--------------------------------------- Public methods -----------------------------------------//
45 //------------------------------------------------------------------------------------------------//
46 
47 
48 /**
49  * Constructor - links the HDF5 dataset, source (sampled matrix), sensor mask and the reduce
50  * operator together. The constructor DOES NOT allocate memory because the size of the sensor mask
51  * is not known at the time the instance of the class is being created.
52  *
53  * @param [in] file - Handle to the HDF5 (output) file
54  * @param [in] datasetName - The dataset's name (index based sensor data stored in a single dataset)
55  * @param [in] sourceMatrix - The source matrix (only real matrices are supported)
56  * @param [in] sensorMask - Index based sensor mask
57  * @param [in] reduceOp - Reduce operator
58  */
60  TMatrixName& datasetName,
61  const TRealMatrix& sourceMatrix,
62  const TIndexMatrix& sensorMask,
63  const TReduceOperator reduceOp)
64  : TBaseOutputHDF5Stream(file, datasetName, sourceMatrix, reduceOp),
65  sensorMask(sensorMask),
66  dataset(H5I_BADID),
67  sampledTimeStep(0),
68  eventSamplingFinished()
69 {
70  // Create event for sampling
71  checkCudaErrors(cudaEventCreate(&eventSamplingFinished));
72 }// end of TIndexOutputHDF5Stream
73 //--------------------------------------------------------------------------------------------------
74 
75 /**
76  * Destructor.
77  * If the file is still opened, it applies the post processing and flush the data.
78  * Then, the object memory is freed and the object destroyed.
79  */
81 {
82  // Destroy sampling event
83  checkCudaErrors(cudaEventDestroy(eventSamplingFinished));
84 
85  Close();
86  // free memory
87  FreeMemory();
88 }// end of Destructor
89 //--------------------------------------------------------------------------------------------------
90 
91 /**
92  * Create a HDF5 stream, create a dataset, and allocate data for it.
93  */
95 {
96 
97  size_t nSampledElementsPerStep = sensorMask.GetElementCount();
98 
99  const TParameters& params = TParameters::GetInstance();
100 
101  // Derive dataset dimension sizes
102  TDimensionSizes datasetSize(nSampledElementsPerStep,
103  (reduceOp == NONE) ? params.Get_nt() - params.GetStartTimeIndex() : 1,
104  1);
105 
106  // Set HDF5 chunk size
107  TDimensionSizes chunkSize(nSampledElementsPerStep, 1, 1);
108  // for chunks bigger than 32 MB
109  if (nSampledElementsPerStep > (CHUNK_SIZE_4MB * 8))
110  {
111  chunkSize.nx = CHUNK_SIZE_4MB; // set chunk size to MB
112  }
113 
114  // Create a dataset under the root group
117  datasetSize,
118  chunkSize,
119  params.GetCompressionLevel());
120 
121  // Write dataset parameters
123  file.WriteMatrixDataType (file.GetRootGroup(), rootObjectName, THDF5_File::FLOAT);
124 
125  // Sampled time step
126  sampledTimeStep = 0;
127 
128  // Set buffer size
129  bufferSize = nSampledElementsPerStep;
130 
131  // Allocate memory
132  AllocateMemory();
133 }// end of Create
134 //--------------------------------------------------------------------------------------------------
135 
136 /**
137  * Reopen the output stream after restart.
138  *
139  */
141 {
142  // Get parameters
143  const TParameters& params = TParameters::GetInstance();
144 
145  // Set buffer size
147 
148  // Allocate memory
149  AllocateMemory();
150 
151  // Reopen the dataset
153 
154 
155  if (reduceOp == NONE)
156  { // raw time series - just seek to the right place in the dataset
157  sampledTimeStep = (params.Get_t_index() < params.GetStartTimeIndex()) ?
158  0 : (params.Get_t_index() - params.GetStartTimeIndex());
159 
160  }
161  else
162  { // aggregated quantities - reload data
163  sampledTimeStep = 0;
164 
165  // Read data from disk only if there were anything stored there (t_index >= start_index)
167  {
168  // Since there is only a single timestep in the dataset, I can read the whole dataset
172  hostBuffer);
173 
174  // Send data to device
176  }
177  }
178 }// end of Reopen
179 //--------------------------------------------------------------------------------------------------
180 
181 
182 /**
183  * Sample grid points, line them up in the buffer, if necessary a reduce operator is applied.
184  *
185  * @warning data is not flushed, there is no sync.
186  */
188 {
189  switch (reduceOp)
190  {
191  case NONE :
192  {
193  OutputStreamsCUDAKernels::SampleIndex<NONE>
194  (deviceBuffer,
198 
199  // Record an event when the data has been copied over.
200  checkCudaErrors(cudaEventRecord(eventSamplingFinished));
201 
202  break;
203  }// case NONE
204 
205  case RMS :
206  {
207  OutputStreamsCUDAKernels::SampleIndex<RMS>
208  (deviceBuffer,
212 
213  break;
214  }// case RMS
215 
216  case MAX :
217  {
218  OutputStreamsCUDAKernels::SampleIndex<MAX>
219  (deviceBuffer,
223  break;
224  }// case MAX
225 
226  case MIN :
227  {
228  OutputStreamsCUDAKernels::SampleIndex<MIN>
229  (deviceBuffer,
233  break;
234  } //case MIN
235  }// switch
236 }// end of Sample
237 //--------------------------------------------------------------------------------------------------
238 
239 
240 /**
241  * Flush data for the timestep. Only applicable on RAW data series.
242  */
244 {
245  if (reduceOp == NONE)
246  {
247  // make sure the data has been copied from the GPU
248  cudaEventSynchronize(eventSamplingFinished);
249 
250  // only raw time series are flushed down to the disk every time step
252  }
253 }// end of FlushRaw
254 //--------------------------------------------------------------------------------------------------
255 
256 
257 /**
258  * Apply post-processing on the buffer and flush it to the file.
259  */
261 {
262  // run inherited method
264 
265  // When no reduction operator is applied, the data is flushed after every time step
266  // which means it has been done before
267  if (reduceOp != NONE)
268  {
269  // Copy data from GPU matrix
271  // flush to disk
273  }
274 }// end of PostProcessing
275 //-------------------------------------------------------------------------------------------------
276 
277 
278 /**
279  * Checkpoint the stream and close.
280  *
281  */
283 {
284  // raw data has already been flushed, others has to be flushed here
285  if (reduceOp != NONE)
286  {
287  // copy data from the device
289  // flush to disk
291  }
292 }// end of Checkpoint
293 //--------------------------------------------------------------------------------------------------
294 
295 /**
296  * Close stream (apply post-processing if necessary, flush data and close).
297  */
299 {
300  // the dataset is still opened
301  if (dataset != H5I_BADID)
302  {
304  }
305 
306  dataset = H5I_BADID;
307 }// end of Close
308 //--------------------------------------------------------------------------------------------------
309 
310 //------------------------------------------------------------------------------------------------//
311 //-------------------------------------- Protected methods ---------------------------------------//
312 //------------------------------------------------------------------------------------------------//
313 
314 
315 /**
316  * Flush the buffer down to the file at the actual position
317  */
319 {
323  hostBuffer);
324  sampledTimeStep++;
325 }// end of FlushToFile
326 //--------------------------------------------------------------------------------------------------
327 
328 //------------------------------------------------------------------------------------------------//
329 //--------------------------------------- Private methods ----------------------------------------//
330 //------------------------------------------------------------------------------------------------//
331 
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 size_t * GetDeviceData()
Get raw GPU data out of the class (for direct GPU kernel access).
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
const TReduceOperator reduceOp
Reduce operator.
virtual void Reopen()
Reopen the output stream after restart and reload data.
virtual void Sample()
Sample data into buffer, apply reduction or copy to the CPU side.
The header file containing the parameters of the simulation.
The header file of cuda kernels used for data sampling (output streams).
static const size_t CHUNK_SIZE_4MB
chunk size of 4MB in number of float elements
const std::string TMatrixName
Datatype for matrix names.
Definition: MatrixNames.h:45
#define checkCudaErrors(val)
Macro checking cuda errors and printing the file name and line. Inspired by CUDA common checking rout...
Definition: Logger.h:209
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.
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
TIndexOutputHDF5Stream(THDF5_File &file, TMatrixName &datasetName, const TRealMatrix &sourceMatrix, const TIndexMatrix &sensorMask, const TReduceOperator reduceOp)
Constructor.
virtual void Checkpoint()
Checkpoint the stream.
virtual void FlushBufferToFile()
Flush the buffer to the file.
The header file containing a class responsible for printing out info and error messages (stdout...
Class storing all parameters of the simulation.
Definition: Parameters.h:49
hid_t GetRootGroup() const
Get handle to the root group.
Definition: HDF5_File.h:573
TReduceOperator
How to aggregate data.
float * deviceBuffer
Temporary buffer on the GPU side - only for aggregated quantities.
store actual data (time series)
cudaEvent_t eventSamplingFinished
Has the sampling finished?
Abstract base class for output data streams (sampled data).
virtual void Close()
Close stream (apply post-processing if necessary, flush data and close).
virtual void FreeMemory()
A generic function to free memory - not used in the base class.
virtual size_t GetElementCount() const
Get total element count of the matrix.
hid_t dataset
Handle to a HDF5 dataset.
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.
The class for 64b unsigned integers (indices). It is used for sensor_mask_index or sensor_corners_mas...
Definition: IndexMatrix.h:46
The header file of the class saving data based on the index senor mask into the output HDF5 file...
virtual void FlushRaw()
Flush data to disk (from raw streams only).
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.
size_t sampledTimeStep
Time step to store (N/A for aggregated)
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
size_t Get_nt() const
Get Nt value.
Definition: Parameters.h:100
virtual void PostProcess()
Apply post-processing on the buffer and flush it to the file.
void CloseDataset(const hid_t dataset)
Close the HDF5 dataset.
Definition: HDF5_File.cpp:457
virtual ~TIndexOutputHDF5Stream()
Destructor.
Class wrapping the HDF5 routines.
Definition: HDF5_File.h:500
virtual void Create()
Create a HDF5 stream and allocate data for it.
const TIndexMatrix & sensorMask
Sensor mask to sample data.
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