Scrambler
1
|
00001 !######################################################################### 00002 ! 00003 ! Copyright (C) 2003-2012 Department of Physics and Astronomy, 00004 ! University of Rochester, 00005 ! Rochester, NY 00006 ! 00007 ! hdf5_declarations.f90 is part of AstroBEAR. 00008 ! 00009 ! AstroBEAR is free software: you can redistribute it and/or modify 00010 ! it under the terms of the GNU General Public License as published by 00011 ! the Free Software Foundation, either version 3 of the License, or 00012 ! (at your option) any later version. 00013 ! 00014 ! AstroBEAR is distributed in the hope that it will be useful, 00015 ! but WITHOUT ANY WARRANTY; without even the implied warranty of 00016 ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 00017 ! GNU General Public License for more details. 00018 ! 00019 ! You should have received a copy of the GNU General Public License 00020 ! along with AstroBEAR. If not, see <http://www.gnu.org/licenses/>. 00021 ! 00022 !######################################################################### 00025 00029 00034 MODULE HDF5Declarations 00035 00036 USE HDF5 00037 USE GlobalDeclarations 00038 USE TreeDeclarations 00039 00040 IMPLICIT NONE 00041 SAVE 00042 PUBLIC 00043 00044 ! Chombo Parameter Declarations 00045 CHARACTER(LEN = 1), PARAMETER :: S_ROOT_GROUP = '/' 00046 REAL(KIND = qPrec), PARAMETER :: DBL_TEST_REAL_DAT = 1.0 00047 INTEGER(HID_T), PARAMETER :: HID_ATTRIBUTE_RANK = 0 00048 INTEGER(HSIZE_T), DIMENSION(1), PARAMETER :: IA_SCALAR_ATTRIB_DIMS = (/1/) 00049 INTEGER, PARAMETER :: I_MAX_COMPONENTS = 100 00050 INTEGER, PARAMETER :: I_MAX_CNAME_LENGTH = 100 00051 INTEGER, PARAMETER :: I_MAX_LEVELNAME_LENGTH = 12 00052 INTEGER, PARAMETER :: I_FILENAME_LENGTH = 19 00053 INTEGER, PARAMETER :: I_DATASET_RANK = 1 00054 INTEGER, PARAMETER :: I_DEFLATE_LEVEL=6 00055 00056 ! A coordinate box requires six integers worth of space. 00057 INTEGER, PARAMETER :: I_CHOMBO_BOX_SIZE = 6 00058 00059 ! The size of dataset chunks. 00060 ! HDF5 currently can't handle chunk sizes above 4 GB. We'll do 1 GB to play it safe. 00061 INTEGER, PARAMETER :: I_MAX_CHUNK_SIZE = 1073741824 ! 1 GB 00062 00063 INTEGER(HID_T) :: hid_intvect_id 00064 INTEGER(HID_T) :: hid_floatvect_id 00065 INTEGER(HID_T) :: hid_box_id 00066 00067 LOGICAL, PARAMETER :: L_EXTENSIBLE = .TRUE. 00068 LOGICAL, PARAMETER :: L_NON_EXTENSIBLE = .FALSE. 00069 00070 00071 CONTAINS 00072 00073 00078 SUBROUTINE CatchHDF5Error(s_function_name, i_err, target_code) 00079 00080 CHARACTER(LEN=*), INTENT(IN) :: s_function_name 00081 INTEGER :: i_err 00082 INTEGER, OPTIONAL :: target_code 00083 00084 INTEGER :: success_code 00085 00086 00087 ! The default success code for most HDF5 calls is 0. If there is no designated success 00088 ! code, then just set it to 0. 00089 IF (.NOT. PRESENT(target_code)) THEN 00090 success_code = 0 00091 ELSE 00092 success_code = target_code 00093 END IF 00094 00095 IF (i_err /= success_code) THEN 00096 PRINT *, "HDF5 error: function ", s_function_name, " failed with error code ", i_err, "." 00097 STOP 00098 END IF 00099 00100 END SUBROUTINE CatchHDF5Error 00101 00103 SUBROUTINE HDF5Init() 00104 IMPLICIT NONE 00105 INTEGER :: i_err 00106 00107 ! Initializes the HDF5 library. This function is safe to call 00108 ! even if it has already been called. 00109 CALL h5open_f(i_err) 00110 CALL CatchHDF5Error("HDF5Init::h5open_f", i_err) 00111 CALL CreateHDF5IntVector(hid_intvect_id) 00112 CALL CreateHDF5FloatVector(hid_floatvect_id) 00113 CALL CreateHDF5BoxType(hid_box_id) 00114 00115 END SUBROUTINE HDF5Init 00116 00118 SUBROUTINE HDF5Finalize() 00119 00120 INTEGER :: i_err 00121 00122 CALL h5tclose_f(hid_intvect_id, i_err) 00123 CALL CatchHDF5Error("h5tclose_f(intvect_id)", i_err) 00124 CALL h5tclose_f(hid_floatvect_id, i_err) 00125 CALL CatchHDF5Error("h5tclose_f(floatvect_id)", i_err) 00126 CALL h5tclose_f(hid_box_id, i_err) 00127 CALL CatchHDF5Error("h5tclose_f(box_id)", i_err) 00128 00129 CALL h5close_f(i_err) 00130 CALL CatchHDF5Error("h5close_f", i_err) 00131 00132 END SUBROUTINE HDF5Finalize 00133 00136 LOGICAL FUNCTION HandleIsValid(handle) 00137 00138 INTEGER(HID_T) :: handle 00139 00140 INTEGER :: i_err 00141 LOGICAL :: result 00142 00143 00144 CALL h5iis_valid_f(handle, result, i_err) 00145 CALL CatchHDF5Error("HandleIsValid::h5iisvalid_f", i_err) 00146 00147 HandleIsValid = result 00148 00149 END FUNCTION HandleIsValid 00150 00151 !=============================== Begin HDF5 Type Definitions ============================== 00154 00157 SUBROUTINE CreateHDF5IntVector(hid_intvect_id) 00158 00159 USE HDF5 00160 00161 INTEGER(HID_T) :: hid_intvect_id 00162 00163 INTEGER(size_t) :: i_h5int_type_size 00164 INTEGER(size_t) :: hid_intvectid_type_size 00165 INTEGER(size_t) :: i_offset 00166 INTEGER :: i_err 00167 00168 00169 i_err = 0 00170 ! Get size of H5T_NATIVE_INTEGER data type. 00171 CALL h5tget_size_f(H5T_NATIVE_INTEGER, i_h5int_type_size, i_err) 00172 CALL CatchHDF5Error("CreateHDF5IntVector::h5tget_size_f(integer)", i_err) 00173 00174 ! Generate size of vector type 00175 hid_intvectid_type_size = 3 * i_h5int_type_size 00176 00177 ! Create the intvectid data type. 00178 CALL h5tcreate_f(H5T_COMPOUND_F, hid_intvectid_type_size, hid_intvect_id, i_err) 00179 CALL CatchHDF5Error("CreateHDF5IntVector::h5tcreate_f (intvector)", i_err) 00180 00181 ! Clear i_offset. 00182 i_offset = 0 00183 00184 ! Insert the three directional fields <i, j, k> into the integer vector data type. 00185 00186 CALL h5tinsert_f(hid_intvect_id, "intvecti", i_offset, H5T_NATIVE_INTEGER, i_err) 00187 i_offset = i_offset + i_h5int_type_size 00188 00189 CALL h5tinsert_f(hid_intvect_id, "intvectj", i_offset, H5T_NATIVE_INTEGER, i_err) 00190 i_offset = i_offset + i_h5int_type_size 00191 00192 CALL h5tinsert_f(hid_intvect_id, "intvectk", i_offset, H5T_NATIVE_INTEGER, i_err) 00193 00194 ! Die if unable to insert all three directional fields into intvectid. 00195 CALL CatchHDF5Error("CreateHDF5IntVector::h5tinsert_f(intvect_*)", i_err) 00196 00197 00198 END SUBROUTINE CreateHDF5IntVector 00199 00202 SUBROUTINE CreateHDF5FloatVector(hid_floatvect_id) 00203 00204 INTEGER(HID_T) :: hid_floatvect_id 00205 00206 INTEGER(SIZE_T) :: i_h5dbl_type_size 00207 INTEGER(SIZE_T) :: hid_floatvectid_type_size 00208 INTEGER(SIZE_T) :: i_offset 00209 INTEGER :: i_err 00210 00211 00212 ! Get size of H5T_NATIVE_DOUBLE data type. 00213 CALL h5tget_size_f(H5T_NATIVE_DOUBLE, i_h5dbl_type_size, i_err) 00214 CALL CatchHDF5Error("h5tget_size_f (double)", i_err) 00215 00216 hid_floatvectid_type_size = 3 * i_h5dbl_type_size 00217 00218 ! Create the floatvectid data type. 00219 CALL h5tcreate_f(H5T_COMPOUND_F, hid_floatvectid_type_size, hid_floatvect_id, i_err) 00220 CALL CatchHDF5Error("h5tcreate_f (floatvect)", i_err) 00221 00222 ! Clear i_offset 00223 i_offset = 0 00224 00225 ! Insert the three directional fields <i, j, k> into the floating-point vector 00226 ! data type. 00227 CALL h5tinsert_f(hid_floatvect_id, "x", i_offset, H5T_NATIVE_DOUBLE, i_err) 00228 i_offset = i_offset + i_h5dbl_type_size 00229 00230 CALL h5tinsert_f(hid_floatvect_id, "y", i_offset, H5T_NATIVE_DOUBLE, i_err) 00231 i_offset = i_offset + i_h5dbl_type_size 00232 00233 CALL h5tinsert_f(hid_floatvect_id, "z", i_offset, H5T_NATIVE_DOUBLE, i_err) 00234 00235 ! Die if unable to insert all three directional fields into floatvectid. 00236 CALL CatchHDF5Error("h5tinsert_f (floatvect_*)", i_err) 00237 00238 END SUBROUTINE CreateHDF5FloatVector 00239 00242 SUBROUTINE CreateHDF5BoxType(hid_box_id) 00243 00244 INTEGER(HID_T) :: hid_box_id 00245 00246 INTEGER(size_t) :: i_h5int_type_size 00247 INTEGER(size_t) :: hid_boxid_type_size 00248 INTEGER(size_t) :: i_offset 00249 INTEGER :: i_err 00250 00251 00252 ! Get size of H5T_NATIVE_INTEGER data type. 00253 CALL h5tget_size_f(H5T_NATIVE_INTEGER, i_h5int_type_size, i_err) 00254 CALL CatchHDF5Error("CreateHDF5BoxType::h5tget_size_f(integer)", i_err) 00255 00256 ! Set box data size 00257 hid_boxid_type_size = I_CHOMBO_BOX_SIZE * i_h5int_type_size 00258 00259 ! Get size of box type. 00260 CALL h5tcreate_f(H5T_COMPOUND_F, hid_boxid_type_size, hid_box_id, i_err) 00261 CALL CatchHDF5Error("CreateHDF5BoxType::h5tcreate_f(boxID)", i_err) 00262 00263 ! Clear i_offset 00264 i_offset = 0 00265 00266 ! Insert the six boundary-descriptors into the box data type. 00267 CALL h5tinsert_f(hid_box_id, "lo_i", i_offset, H5T_NATIVE_INTEGER, i_err) 00268 i_offset = i_offset + i_h5int_type_size 00269 00270 CALL h5tinsert_f(hid_box_id, "lo_j", i_offset, H5T_NATIVE_INTEGER, i_err) 00271 i_offset = i_offset + i_h5int_type_size 00272 00273 CALL h5tinsert_f(hid_box_id, "lo_k", i_offset, H5T_NATIVE_INTEGER, i_err) 00274 i_offset = i_offset + i_h5int_type_size 00275 00276 CALL h5tinsert_f(hid_box_id, "hi_i", i_offset, H5T_NATIVE_INTEGER, i_err) 00277 i_offset = i_offset + i_h5int_type_size 00278 00279 CALL h5tinsert_f(hid_box_id, "hi_j", i_offset, H5T_NATIVE_INTEGER, i_err) 00280 i_offset = i_offset + i_h5int_type_size 00281 00282 CALL h5tinsert_f(hid_box_id, "hi_k", i_offset, H5T_NATIVE_INTEGER, i_err) 00283 00284 ! Die if unable to insert all six directional fields into hid_box_id. 00285 CALL CatchHDF5Error("CreateHDF5BoxType::h5tinsert_f(boxID)", i_err) 00286 00287 END SUBROUTINE CreateHDF5BoxType 00288 00291 00292 !================================ Begin HDF5 Read Functions =============================== 00293 00296 00302 SUBROUTINE Read_Slab_From_Dataset_Box(s_name, hid_group_id, box_data, i_offset) 00303 00304 IMPLICIT NONE 00305 00306 ! Input parameter declarations 00307 CHARACTER(LEN=*), INTENT(IN) :: s_name ! Name of new dataset. 00308 INTEGER(HID_T), INTENT(IN) :: hid_group_id ! Group identifier. 00309 INTEGER, DIMENSION(6), INTENT(OUT) :: box_data ! Data input 00310 INTEGER(HSIZE_T), INTENT(IN) :: i_offset 00311 00312 INTEGER(HID_T) :: hid_dataset_id ! Dataset handle. 00313 INTEGER(HID_T) :: hid_dataspace_id ! Dataspace handle. 00314 INTEGER(HID_T) :: hid_memspace_id ! Dataspace handle for hyperslab reference. 00315 INTEGER(HSIZE_T) :: i_dataset_size ! Size of dataset to be retrieved. 00316 INTEGER(HSIZE_T), DIMENSION(1) :: ia_dims ! Array to hold the size of the input. 00317 INTEGER :: i_err 00318 INTEGER(HSIZE_T), DIMENSION(1) :: ia_data_offset 00319 INTEGER(HSIZE_T), DIMENSION(1) :: ia_slab_size 00320 00321 00322 ! Set the size of the hyperslab. 00323 ia_slab_size(1) = SIZE(box_data) / I_CHOMBO_BOX_SIZE 00324 00325 ! Store the size of the offset in a dimension array. 00326 ia_data_offset(1) = i_offset 00327 00328 !PRINT *, "Read_Slab_From_Dataset_Box::s_name = ", s_name 00329 !PRINT *, "Read_Slab_From_Dataset_Box::slab size = ", ia_slab_size(1) 00330 !PRINT *, "Read_Slab_From_Dataset_Box::offset size = ", i_offset 00331 00332 ! Open the dataset in the HDF5 file. 00333 CALL h5dopen_f(hid_group_id, s_name, hid_dataset_id, i_err) 00334 CALL CatchHDF5Error("Read_Slab_From_Dataset_Box::h5dopen_f", i_err) 00335 00336 ! Retrieve the total dataspace for the dataset. 00337 CALL h5dget_space_f(hid_dataset_id, hid_dataspace_id, i_err) 00338 CALL CatchHDF5Error("Read_Slab_From_Dataset_Box::h5dget_space_f", i_err) 00339 00340 CALL h5sselect_hyperslab_f(hid_dataspace_id, H5S_SELECT_SET_F, & 00341 ia_data_offset, ia_slab_size, i_err) 00342 CALL CatchHDF5Error("Read_Slab_From_Dataset_Box::h5sselect_hyperslab_f", i_err) 00343 00344 ! Create dataspace for the upcoming hyperslab. 00345 CALL h5screate_simple_f(I_DATASET_RANK, ia_slab_size, hid_memspace_id, i_err) 00346 CALL CatchHDF5Error("Read_Slab_From_Dataset_Box::hscreate_simple_f", i_err) 00347 00348 ! Read data into newly-allocated array. 00349 CALL h5dread_f(hid_dataset_id, hid_box_id, box_data, ia_dims, i_err, & 00350 hid_memspace_id, hid_dataspace_id) 00351 CALL CatchHDF5Error("Read_Slab_From_Dataset_Box::h5dread_f", i_err) 00352 00353 ! Close HDF5 handles. 00354 CALL h5sclose_f(hid_memspace_id, i_err) 00355 CALL CatchHDF5Error("Read_Slab_From_Dataset_Box::h5sclose_f(hid_memspace_id)", i_err) 00356 00357 CALL h5sclose_f(hid_dataspace_id, i_err) 00358 CALL CatchHDF5Error("Read_Slab_From_Dataset_Box::h5sclose_f(hid_dataspace_id)", i_err) 00359 00360 CALL h5dclose_f(hid_dataset_id, i_err) 00361 CALL CatchHDF5Error("Read_Slab_From_Dataset_Box::h5dclose_f(hid_dataset_id)", i_err) 00362 00363 END SUBROUTINE Read_Slab_From_Dataset_Box 00364 00365 00371 SUBROUTINE Read_Slab_From_Dataset_Double(s_name, hid_group_id, dbla_data, i_offset) 00372 00373 IMPLICIT NONE 00374 00375 ! Input parameter declarations 00376 CHARACTER(LEN=*), INTENT(IN) :: s_name ! Name of new dataset. 00377 INTEGER(HID_T), INTENT(IN) :: hid_group_id ! Group identifier. 00378 REAL(KIND=qPrec), POINTER, DIMENSION(:), INTENT(OUT) :: dbla_data ! Data input 00379 INTEGER(HSIZE_T), INTENT(IN) :: i_offset 00380 00381 INTEGER(HID_T) :: hid_dataset_id ! Dataset handle. 00382 INTEGER(HID_T) :: hid_dataspace_id ! Dataspace handle. 00383 INTEGER(HID_T) :: hid_memspace_id ! Dataspace handle for hyperslab reference. 00384 INTEGER(HSIZE_T) :: i_dataset_size ! Size of dataset to be retrieved. 00385 INTEGER(HSIZE_T), DIMENSION(1) :: ia_dims ! Array to hold the size of the input. 00386 INTEGER :: i_err 00387 INTEGER(HSIZE_T), DIMENSION(1) :: ia_data_offset 00388 INTEGER(HSIZE_T), DIMENSION(1) :: ia_slab_size 00389 00390 00391 ! Set the size of the hyperslab. 00392 ia_slab_size(1) = SIZE(dbla_data) 00393 00394 ! Store the size of the offset in a dimension array. 00395 ia_data_offset(1) = i_offset 00396 00397 ! Open the dataset in the HDF5 file. 00398 CALL h5dopen_f(hid_group_id, s_name, hid_dataset_id, i_err) 00399 CALL CatchHDF5Error("Read_Slab_From_Dataset_Double::h5dopen_f", i_err) 00400 00401 00402 ! Retrieve the total dataspace for the dataset. 00403 CALL h5dget_space_f(hid_dataset_id, hid_dataspace_id, i_err) 00404 CALL CatchHDF5Error("Read_Slab_From_Dataset_Double::h5dget_space_f", i_err) 00405 00406 00407 CALL h5sselect_hyperslab_f(hid_dataspace_id, H5S_SELECT_SET_F, & 00408 ia_data_offset, ia_slab_size, i_err) 00409 CALL CatchHDF5Error("Read_Slab_From_Dataset_Double::h5sselect_hyperslab_f", i_err) 00410 00411 00412 ! Create dataspace for the upcoming hyperslab. 00413 CALL h5screate_simple_f(I_DATASET_RANK, ia_slab_size, hid_memspace_id, i_err) 00414 CALL CatchHDF5Error("Read_Slab_From_Dataset_Double::hscreate_simple_f", i_err) 00415 00416 ! Read data into newly-allocated array. 00417 CALL h5dread_f(hid_dataset_id, H5T_NATIVE_DOUBLE, dbla_data, ia_dims, i_err, & 00418 hid_memspace_id, hid_dataspace_id) 00419 CALL CatchHDF5Error("Read_Slab_From_Dataset_Double::h5dread_f", i_err) 00420 00421 00422 ! Close HDF5 handles. 00423 CALL h5sclose_f(hid_memspace_id, i_err) 00424 CALL CatchHDF5Error("Read_Slab_From_Dataset_Double::h5sclose_f(hid_memspace_id)", i_err) 00425 CALL h5sclose_f(hid_dataspace_id, i_err) 00426 CALL CatchHDF5Error("Read_Slab_From_Dataset_Double::h5sclose_f(hid_dataspace_id)", i_err) 00427 CALL h5dclose_f(hid_dataset_id, i_err) 00428 CALL CatchHDF5Error("Read_Slab_From_Dataset_Double::h5dclose_f(hid_dataset_id)", i_err) 00429 00430 END SUBROUTINE Read_Slab_From_Dataset_Double 00431 00437 SUBROUTINE Read_Slab_From_Dataset_Int(s_name, hid_group_id, int_data, i_offset) 00438 00439 IMPLICIT NONE 00440 00441 ! Input parameter declarations 00442 CHARACTER(LEN=*), INTENT(IN) :: s_name ! Name of new dataset. 00443 INTEGER(HID_T), INTENT(IN) :: hid_group_id ! Group identifier. 00444 INTEGER, DIMENSION(:), INTENT(OUT) :: int_data ! Data input 00445 INTEGER(HSIZE_T), INTENT(IN) :: i_offset 00446 00447 INTEGER(HID_T) :: hid_dataset_id ! Dataset handle. 00448 INTEGER(HID_T) :: hid_dataspace_id ! Dataspace handle. 00449 INTEGER(HID_T) :: hid_memspace_id ! Dataspace handle for hyperslab reference. 00450 INTEGER(HSIZE_T) :: i_dataset_size ! Size of dataset to be retrieved. 00451 INTEGER(HSIZE_T), DIMENSION(1) :: ia_dims ! Array to hold the size of the input. 00452 INTEGER :: i_err 00453 INTEGER(HSIZE_T), DIMENSION(1) :: ia_data_offset 00454 INTEGER(HSIZE_T), DIMENSION(1) :: ia_slab_size 00455 00456 00457 ! Set the size of the hyperslab. 00458 ia_slab_size(1) = SIZE(int_data) 00459 00460 ! Store the size of the offset in a dimension array. 00461 ia_data_offset(1) = i_offset 00462 00463 ! Open the dataset in the HDF5 file. 00464 CALL h5dopen_f(hid_group_id, s_name, hid_dataset_id, i_err) 00465 CALL CatchHDF5Error("Read_Slab_From_Dataset_Int::h5dopen_f", i_err) 00466 00467 00468 ! Retrieve the total dataspace for the dataset. 00469 CALL h5dget_space_f(hid_dataset_id, hid_dataspace_id, i_err) 00470 CALL CatchHDF5Error("Read_Slab_From_Dataset_Int::h5dget_space_f", i_err) 00471 00472 00473 CALL h5sselect_hyperslab_f(hid_dataspace_id, H5S_SELECT_SET_F, & 00474 ia_data_offset, ia_slab_size, i_err) 00475 CALL CatchHDF5Error("Read_Slab_From_Dataset_Int::h5sselect_hyperslab_f", i_err) 00476 00477 00478 ! Create dataspace for the upcoming hyperslab. 00479 CALL h5screate_simple_f(I_DATASET_RANK, ia_slab_size, hid_memspace_id, i_err) 00480 CALL CatchHDF5Error("Read_Slab_From_Dataset_Int::hscreate_simple_f", i_err) 00481 00482 00483 00484 ! Read data into newly-allocated array. 00485 CALL h5dread_f(hid_dataset_id, H5T_NATIVE_INTEGER, int_data, ia_dims, i_err, & 00486 hid_memspace_id, hid_dataspace_id) 00487 CALL CatchHDF5Error("Read_Slab_From_Dataset_Int::h5dread_f", i_err) 00488 00489 00490 00491 ! Close HDF5 handles. 00492 CALL h5sclose_f(hid_memspace_id, i_err) 00493 CALL CatchHDF5Error("Read_Slab_From_Dataset_Int::h5sclose_f(hid_memspace_id)", i_err) 00494 CALL h5sclose_f(hid_dataspace_id, i_err) 00495 CALL CatchHDF5Error("Read_Slab_From_Dataset_Int::h5sclose_f(hid_dataspace_id)", i_err) 00496 CALL h5dclose_f(hid_dataset_id, i_err) 00497 CALL CatchHDF5Error("Read_Slab_From_Dataset_Int::h5dclose_f(hid_dataset_id)", i_err) 00498 00499 00500 END SUBROUTINE Read_Slab_From_Dataset_Int 00501 00505 INTEGER FUNCTION Read_HDF5_Attribute_Int(s_name, hid_group_id) 00506 00507 IMPLICIT NONE 00508 00509 ! Input parameter declarations 00510 CHARACTER(LEN=*), INTENT(IN) :: s_name ! Name of new dataset. 00511 INTEGER(HID_T), INTENT(IN) :: hid_group_id ! Group identifier. 00512 00513 ! Variable declarations 00514 INTEGER(HID_T) :: hid_attribute_id ! Dataset handle. 00515 INTEGER(HSIZE_T), DIMENSION(1) :: ia_dims ! Array to hold the size of the input. 00516 INTEGER :: i_output ! Used to store output from HDF5 API. 00517 INTEGER :: i_err 00518 00519 00520 ! Open the attribute in the HDF5 file. 00521 CALL h5aopen_name_f(hid_group_id, s_name, hid_attribute_id, i_err) 00522 00523 ! Die if unable to open dataset. 00524 IF (i_err < 0) THEN 00525 PRINT *, "Read_HDF5_Attribute_Int error ", i_err, ": unable to open ", & 00526 "attribute ", s_name, "." 00527 STOP 00528 END IF 00529 00530 00531 ! Initializing dimension array with size of chombo array. 00532 ia_dims(1) = 1 00533 00534 ! Read data into newly-allocated array. 00535 CALL h5aread_f(hid_attribute_id, H5T_NATIVE_INTEGER, & 00536 i_output, ia_dims, i_err) 00537 00538 ! Die if unable to read data in from dataset. 00539 IF (i_err < 0) THEN 00540 PRINT *, "Read_HDF5_Attribute_Int error ", i_err, ": unable to read value ", & 00541 "from attribute ", s_name, "." 00542 STOP 00543 END IF 00544 00545 ! Close dataset. 00546 CALL h5aclose_f(hid_attribute_id, i_err) 00547 00548 ! Die if unable to close dataset. 00549 IF (i_err < 0) THEN 00550 PRINT *, "Read_HDF5_Attribute_Int error ", i_err, ": unable to close ", & 00551 "attribute ", s_name, "." 00552 STOP 00553 END IF 00554 00555 Read_HDF5_Attribute_Int = i_output 00556 00557 END FUNCTION Read_HDF5_Attribute_Int 00558 00559 00560 00564 REAL(KIND=xPrec) FUNCTION Read_HDF5_Attribute_Double(s_name, hid_group_id) 00565 00566 IMPLICIT NONE 00567 00568 ! Input parameter declarations 00569 CHARACTER(LEN=*), INTENT(IN) :: s_name ! Name of new dataset. 00570 INTEGER(HID_T), INTENT(IN) :: hid_group_id ! Group identifier. 00571 00572 ! Variable declarations 00573 INTEGER(HID_T) :: hid_attribute_id ! Dataset handle. 00574 INTEGER(HSIZE_T), DIMENSION(1) :: ia_dims ! Array to hold the size of the input. 00575 REAL(KIND=xPrec) :: dbl_output ! Used to store output from HDF5 API. 00576 INTEGER :: i_err 00577 00578 00579 ! Open the attribute in the HDF5 file. 00580 CALL h5aopen_name_f(hid_group_id, s_name, hid_attribute_id, i_err) 00581 00582 ! Die if unable to open dataset. 00583 IF (i_err < 0) THEN 00584 PRINT *, "Read_HDF5_Attribute_Double error ", i_err, ": unable to open ", & 00585 "attribute ", s_name, "." 00586 STOP 00587 END IF 00588 00589 00590 ! Initializing dimension array with size of chombo array. 00591 ia_dims(1) = 1 00592 00593 ! Read data into newly-allocated array. 00594 CALL h5aread_f(hid_attribute_id, H5T_NATIVE_DOUBLE, & 00595 dbl_output, ia_dims, i_err) 00596 00597 ! Die if unable to read data in from dataset. 00598 IF (i_err < 0) THEN 00599 PRINT *, "Read_HDF5_Attribute_Double error ", i_err, ": unable to read value ", & 00600 "from attribute ", s_name, "." 00601 STOP 00602 END IF 00603 00604 ! Close dataset. 00605 CALL h5aclose_f(hid_attribute_id, i_err) 00606 00607 ! Die if unable to close dataset. 00608 IF (i_err < 0) THEN 00609 PRINT *, "Read_HDF5_Attribute_Double error ", i_err, ": unable to close ", & 00610 "attribute ", s_name, "." 00611 STOP 00612 END IF 00613 00614 Read_HDF5_Attribute_Double = dbl_output 00615 00616 END FUNCTION Read_HDF5_Attribute_Double 00617 00621 INTEGER FUNCTION IO_GetDatasetElementCount(group_id, s_name) 00622 00623 INTEGER(HID_T) :: group_id 00624 CHARACTER(LEN=*) :: s_name 00625 00626 INTEGER(HSIZE_T), DIMENSION(1) :: dims, maxdims 00627 INTEGER :: i_err 00628 INTEGER(HID_T) :: dataset_id 00629 INTEGER(HID_T) :: dataspace_id 00630 00631 ! Open the dataset in the Chombo file. 00632 CALL h5dopen_f(group_id, s_name, dataset_id, i_err) 00633 CALL CatchHDF5Error("GetDatasetElementCount::h5dopen_f", i_err) 00634 00635 ! Get the associated dataspace. 00636 CALL h5dget_space_f(dataset_id, dataspace_id, i_err) 00637 CALL CatchHDF5Error("GetDatasetElementCount::h5dget_space_f", i_err) 00638 00639 ! Extract the dimensions of said dataspace. 00640 CALL h5sget_simple_extent_dims_f(dataspace_id, dims, maxdims, i_err) 00641 CALL CatchHDF5Error("GetDatasetElementCount::h5sget_simple_extent_dims_f", i_err, I_DATASET_RANK) 00642 00643 ! Close the dataspace. 00644 CALL h5sclose_f(dataspace_id, i_err) 00645 CALL CatchHDF5Error("GetDatasetElementCount::h5sclose_f", i_err) 00646 00647 ! close the dataset. 00648 CALL h5dclose_f(dataset_id, i_err) 00649 CALL CatchHDF5Error("GetDatasetElementCount::h5dclose_f", i_err) 00650 00651 IO_GetDatasetElementCount = dims(1) 00652 00653 END FUNCTION IO_GetDatasetElementCount 00656 00657 !============================= Begin HDF5 Write Functions =============================== 00660 00665 SUBROUTINE CreateHDF5Group(s_name, location_id, hid_group_id) 00666 00667 CHARACTER(LEN=*) :: s_name 00668 INTEGER(HID_T) :: location_id 00669 INTEGER(HID_T) :: hid_group_id 00670 00671 INTEGER :: iErr 00672 00673 ! Test if the given handle is valid. This will fail on handles from unopened files 00674 ! and handles that are already opened. 00675 IF (.NOT. HandleIsValid(hid_group_id)) THEN 00676 PRINT *, "CreateHDF5Group() error: invalid target handle for group '", s_name, "'." 00677 STOP 00678 ELSE IF (.NOT. HandleIsValid(location_id)) THEN 00679 PRINT *, "CreateHDF5Group() error: invalid location handle for group '", s_name, "'." 00680 STOP 00681 END IF 00682 00683 ! Create a new group within the file. 00684 CALL h5gcreate_f(location_id, s_name, hid_group_id, iErr) 00685 CALL CatchHDF5Error("CreateHDF5Group::h5gcreate_f", iErr) 00686 00687 END SUBROUTINE CreateHDF5Group 00688 00693 SUBROUTINE Add_HDF5_Attribute_Int(s_name, hid_group_id, i_value) 00694 00695 IMPLICIT NONE 00696 00697 ! Input parameter declarations 00698 CHARACTER(LEN=*), INTENT(IN) :: s_name ! Name of new attribute. 00699 INTEGER(HID_T), INTENT(IN) :: hid_group_id ! Group identifier. 00700 INTEGER :: i_value ! Data input 00701 00702 ! Variable declarations 00703 INTEGER(HID_T) :: hid_dataspace_id 00704 INTEGER(HID_T) :: hid_attribute_id 00705 00706 INTEGER :: i_err 00707 00708 00709 ! Create a new dataspace. 00710 CALL h5screate_simple_f(HID_ATTRIBUTE_RANK, IA_SCALAR_ATTRIB_DIMS, & 00711 hid_dataspace_id, i_err) 00712 00713 CALL CatchHDF5Error("Add_HDF5_Attribute_Int::h5screate_simple_f", i_err) 00714 00715 ! Create new attribute within dataspace. 00716 CALL h5acreate_f(hid_group_id, s_name, H5T_NATIVE_INTEGER, hid_dataspace_id, & 00717 hid_attribute_id, i_err) 00718 00719 CALL CatchHDF5Error("Add_HDF5_Attribute_Int::h5acreate_f", i_err) 00720 00721 ! Write the datavalue to the dataset. 00722 CALL h5awrite_f(hid_attribute_id, H5T_NATIVE_INTEGER, i_value, & 00723 IA_SCALAR_ATTRIB_DIMS, i_err) 00724 00725 CALL CatchHDF5Error("Add_HDF5_Attribute_Int::h5awrite_f", i_err) 00726 00727 ! Close the attribute object. 00728 CALL h5aclose_f(hid_attribute_id, i_err) 00729 CALL CatchHDF5Error("Add_HDF5_Attribute_Int::h5aclose_f", i_err) 00730 00731 ! Close the dataspace object. 00732 CALL h5sclose_f(hid_dataspace_id, i_err) 00733 CALL CatchHDF5Error("Add_HDF5_Attribute_Int::h5sclose_f", i_err) 00734 00735 END SUBROUTINE Add_HDF5_Attribute_Int 00736 00737 00738 00743 SUBROUTINE Add_HDF5_Attribute_Double(s_name, hid_group_id, dbl_value) 00744 00745 IMPLICIT NONE 00746 00747 ! Input parameter declarations 00748 CHARACTER(LEN=*), INTENT(IN) :: s_name ! Name of new attribute. 00749 INTEGER(HID_T), INTENT(IN) :: hid_group_id ! Group identifier. 00750 REAL(8) :: dbl_value ! Data input 00751 00752 ! Variable declarations 00753 INTEGER(HID_T) :: hid_dataspace_id 00754 INTEGER(HID_T) :: hid_attribute_id 00755 00756 INTEGER :: i_err 00757 00758 00759 00760 ! Create a new dataspace. 00761 CALL h5screate_simple_f(HID_ATTRIBUTE_RANK, IA_SCALAR_ATTRIB_DIMS, & 00762 hid_dataspace_id, i_err) 00763 CALL CatchHDF5Error("Add_HDF5_Attribute_Double::h5screate_simple_f", i_err) 00764 00765 ! Create new attribute within dataspace. 00766 CALL h5acreate_f(hid_group_id, s_name, H5T_NATIVE_DOUBLE, hid_dataspace_id, & 00767 hid_attribute_id, i_err) 00768 CALL CatchHDF5Error("Add_HDF5_Attribute_Double::h5acreate_f", i_err) 00769 00770 ! Write the datavalue to the dataset. 00771 CALL h5awrite_f(hid_attribute_id, H5T_NATIVE_DOUBLE, dbl_value, & 00772 IA_SCALAR_ATTRIB_DIMS, i_err) 00773 CALL CatchHDF5Error("Add_HDF5_Attribute_Double::h5awrite_f", i_err) 00774 00775 ! Close the attribute object. 00776 CALL h5aclose_f(hid_attribute_id, i_err) 00777 CALL CatchHDF5Error("Add_HDF5_Attribute_Double::h5aclose_f", i_err) 00778 00779 ! Close the dataspace object. 00780 CALL h5sclose_f(hid_dataspace_id, i_err) 00781 CALL CatchHDF5Error("Add_HDF5_Attribute_Double::h5sclose_f", i_err) 00782 00783 END SUBROUTINE Add_HDF5_Attribute_Double 00784 00785 00792 SUBROUTINE Add_HDF5_Attribute_IntVector(s_name, hid_type_id, hid_group_id, ia_value, i_size) 00793 00794 IMPLICIT NONE 00795 00796 ! Input parameter declarations 00797 CHARACTER(LEN=*), INTENT(IN) :: s_name ! Name of new attribute. 00798 INTEGER(HID_T), INTENT(IN) :: hid_type_id ! Describes the type of v_value. 00799 INTEGER(HID_T), INTENT(IN) :: hid_group_id ! Group identifier. 00800 INTEGER, DIMENSION(3) :: ia_value ! Data input 00801 INTEGER(size_t), INTENT(IN), OPTIONAL :: i_size ! Size of input type. 00802 00803 ! Variable declarations 00804 INTEGER(HID_T) :: hid_dataspace_id 00805 INTEGER(HID_T) :: hid_attribute_id 00806 00807 INTEGER :: i_err 00808 00809 ! If a size value was passed in, then set the type's size value to that. 00810 IF (PRESENT(i_size)) THEN 00811 CALL h5tset_size_f(hid_type_id, i_size ,i_err) 00812 CALL CatchHDF5Error("Add_HDF5_Attribute_IntVector::h5tset_size_f", i_err) 00813 END IF 00814 00815 ! Create a new dataspace. 00816 CALL h5screate_simple_f(HID_ATTRIBUTE_RANK, IA_SCALAR_ATTRIB_DIMS, & 00817 hid_dataspace_id, i_err) 00818 CALL CatchHDF5Error("Add_HDF5_Attribute_IntVector::h5screate_simple_f", i_err) 00819 00820 ! Create new attribute within dataspace. 00821 CALL h5acreate_f(hid_group_id, s_name, hid_type_id, hid_dataspace_id, & 00822 hid_attribute_id, i_err) 00823 CALL CatchHDF5Error("Add_HDF5_Attribute_IntVector::h5acreate_f", i_err) 00824 00825 ! Write the datavalue to the dataset. 00826 CALL h5awrite_f(hid_attribute_id, hid_type_id, ia_value, & 00827 IA_SCALAR_ATTRIB_DIMS, i_err) 00828 CALL CatchHDF5Error("Add_HDF5_Attribute_IntVector::h5write_f", i_err) 00829 00830 ! Close the attribute object. 00831 CALL h5aclose_f(hid_attribute_id, i_err) 00832 CALL CatchHDF5Error("Add_HDF5_Attribute_IntVector::h5aclose_f", i_err) 00833 00834 ! Close the dataspace object. 00835 CALL h5sclose_f(hid_dataspace_id, i_err) 00836 CALL CatchHDF5Error("Add_HDF5_Attribute_IntVector::h5sclose_f", i_err) 00837 00838 END SUBROUTINE Add_HDF5_Attribute_IntVector 00839 00840 00841 00847 SUBROUTINE Add_HDF5_Attribute_FloatVector(s_name, hid_group_id, dbla_value, i_size) 00848 00849 IMPLICIT NONE 00850 00851 ! Input parameter declarations 00852 CHARACTER(LEN=*), INTENT(IN) :: s_name ! Name of new attribute. 00853 INTEGER(HID_T), INTENT(IN) :: hid_group_id ! Group identifier. 00854 REAL(KIND=qPrec), DIMENSION(3) :: dbla_value ! Data input 00855 INTEGER(size_t), INTENT(IN), OPTIONAL :: i_size ! Size of input type. 00856 00857 ! Variable declarations 00858 INTEGER(HID_T) :: hid_dataspace_id 00859 INTEGER(HID_T) :: hid_attribute_id 00860 00861 INTEGER :: i_err 00862 00863 00864 ! If a size value was passed in, then set the type's size value to that. 00865 IF (PRESENT(i_size)) THEN 00866 CALL h5tset_size_f(hid_floatvect_id, i_size ,i_err) 00867 CALL CatchHDF5Error("Add_HDF5_Attribute_FloatVector(::h5tset_size_f", i_err) 00868 END IF 00869 00870 ! Create a new dataspace. 00871 CALL h5screate_simple_f(HID_ATTRIBUTE_RANK, IA_SCALAR_ATTRIB_DIMS, & 00872 hid_dataspace_id, i_err) 00873 CALL CatchHDF5Error("Add_HDF5_Attribute_FloatVector::h5screate_simple_f", i_err) 00874 00875 ! Create new attribute within dataspace. 00876 CALL h5acreate_f(hid_group_id, s_name, hid_floatvect_id, hid_dataspace_id, & 00877 hid_attribute_id, i_err) 00878 CALL CatchHDF5Error("Add_HDF5_Attribute_FloatVector::h5acreate_f", i_err) 00879 00880 ! Write the datavalue to the dataset. 00881 CALL h5awrite_f(hid_attribute_id, hid_floatvect_id, dbla_value, & 00882 IA_SCALAR_ATTRIB_DIMS, i_err) 00883 CALL CatchHDF5Error("Add_HDF5_Attribute_FloatVector::h5awrite_f", i_err) 00884 00885 ! Close the attribute object. 00886 CALL h5aclose_f(hid_attribute_id, i_err) 00887 CALL CatchHDF5Error("Add_HDF5_Attribute_FloatVector::h5aclose_f", i_err) 00888 00889 ! Close the dataspace object. 00890 CALL h5sclose_f(hid_dataspace_id, i_err) 00891 CALL CatchHDF5Error("Add_HDF5_Attribute_FloatVector::h5sclose_f", i_err) 00892 00893 END SUBROUTINE Add_HDF5_Attribute_FloatVector 00894 00899 SUBROUTINE Add_HDF5_Attribute_String(s_name, hid_group_id, s_value) 00900 00901 IMPLICIT NONE 00902 00903 ! Input parameter declarations 00904 CHARACTER(LEN=*), INTENT(IN) :: s_name ! Name of new attribute. 00905 INTEGER(HID_T), INTENT(IN) :: hid_group_id ! Group identifier. 00906 CHARACTER(LEN=*) :: s_value ! Data input 00907 00908 ! Variable declarations 00909 INTEGER(HID_T) :: hid_type_id 00910 INTEGER(HID_T) :: hid_dataspace_id 00911 INTEGER(HID_T) :: hid_attribute_id 00912 00913 INTEGER :: i_err 00914 00915 ! Copy the H5T_NATIVE_CHARACTER data type to hid_h5string_type_id. 00916 CALL h5tcopy_f(H5T_NATIVE_CHARACTER, hid_type_id, i_err) 00917 CALL CatchHDF5Error("Add_HDF5_Attribute_String::h5tcopy_f", i_err) 00918 00919 ! Set the type's size value. 00920 CALL h5tset_size_f(hid_type_id, INT(LEN(TRIM(s_value)), size_t), i_err) 00921 CALL CatchHDF5Error("Add_HDF5_Attribute_String::h5tset_size_f", i_err) 00922 00923 ! Create a new dataspace. 00924 CALL h5screate_simple_f(HID_ATTRIBUTE_RANK, IA_SCALAR_ATTRIB_DIMS, & 00925 hid_dataspace_id, i_err) 00926 CALL CatchHDF5Error("Add_HDF5_Attribute_Stringh5screate_simple_f", i_err) 00927 00928 ! Create new attribute within dataspace. 00929 CALL h5acreate_f(hid_group_id, s_name, hid_type_id, hid_dataspace_id, & 00930 hid_attribute_id, i_err) 00931 CALL CatchHDF5Error("Add_HDF5_Attribute_String::h5acreate_f", i_err) 00932 00933 ! Write the datavalue to the dataset. 00934 CALL h5awrite_f(hid_attribute_id, hid_type_id, s_value, & 00935 IA_SCALAR_ATTRIB_DIMS, i_err) 00936 CALL CatchHDF5Error("Add_HDF5_Attribute_String::h5awrite_f", i_err) 00937 00938 ! Close the type object. 00939 CALL h5tclose_f(hid_type_id, i_err) 00940 CALL CatchHDF5Error("Add_HDF5_Attribute_String::h5tclose_f", i_err) 00941 00942 ! Close the attribute object. 00943 CALL h5aclose_f(hid_attribute_id, i_err) 00944 CALL CatchHDF5Error("Add_HDF5_Attribute_String::h5aclose_f", i_err) 00945 00946 ! Close the dataspace object. 00947 CALL h5sclose_f(hid_dataspace_id, i_err) 00948 CALL CatchHDF5Error("Add_HDF5_Attribute_String::h5sclose_f", i_err) 00949 00950 END SUBROUTINE Add_HDF5_Attribute_String 00951 00952 00959 SUBROUTINE Add_HDF5_Attribute_Box(s_name, hid_type_id, hid_group_id, ia_value, i_size) 00960 00961 IMPLICIT NONE 00962 00963 ! Input parameter declarations 00964 CHARACTER(LEN=*), INTENT(IN) :: s_name ! Name of new attribute. 00965 INTEGER(HID_T), INTENT(IN) :: hid_type_id ! Used to set the string size. 00966 INTEGER(HID_T), INTENT(IN) :: hid_group_id ! Group identifier. 00967 INTEGER, DIMENSION(I_CHOMBO_BOX_SIZE) :: ia_value ! Data input 00968 INTEGER(size_t), INTENT(IN), OPTIONAL :: i_size ! Size of input type. 00969 00970 ! Variable declarations 00971 INTEGER(HID_T) :: hid_dataspace_id 00972 INTEGER(HID_T) :: hid_attribute_id 00973 00974 INTEGER :: i_err 00975 00976 00977 ! If a size value was passed in, then set the type's size value to that. 00978 IF (PRESENT(i_size)) THEN 00979 CALL h5tset_size_f(hid_type_id, i_size ,i_err) 00980 CALL CatchHDF5Error("Add_HDF5_Attribute_Box::h5tset_size_f", i_err) 00981 END IF 00982 00983 ! Create a new dataspace. 00984 CALL h5screate_simple_f(HID_ATTRIBUTE_RANK, IA_SCALAR_ATTRIB_DIMS, & 00985 hid_dataspace_id, i_err) 00986 CALL CatchHDF5Error("Add_HDF5_Attribute_Box::h5screate_simple_f", i_err) 00987 00988 ! Create new attribute within dataspace. 00989 CALL h5acreate_f(hid_group_id, s_name, hid_type_id, hid_dataspace_id, & 00990 hid_attribute_id, i_err) 00991 CALL CatchHDF5Error("Add_HDF5_Attribute_Box::h5acreate_f", i_err) 00992 00993 ! Write the datavalue to the dataset. 00994 CALL h5awrite_f(hid_attribute_id, hid_type_id, ia_value, & 00995 IA_SCALAR_ATTRIB_DIMS, i_err) 00996 CALL CatchHDF5Error("Add_HDF5_Attribute_Box::h5awrite_f", i_err) 00997 00998 ! Close the attribute object. 00999 CALL h5aclose_f(hid_attribute_id, i_err) 01000 CALL CatchHDF5Error("Add_HDF5_Attribute_Box::h5aclose_f", i_err) 01001 01002 ! Close the dataspace object. 01003 CALL h5sclose_f(hid_dataspace_id, i_err) 01004 CALL CatchHDF5Error("Add_HDF5_Attribute_Box::h5sclose_f", i_err) 01005 01006 END SUBROUTINE Add_HDF5_Attribute_Box 01007 01013 SUBROUTINE Initialize_HDF5_Dataset_Double(s_name, hid_group_id, raw_data_size, lExtensible) 01014 01015 IMPLICIT NONE 01016 01017 ! Input parameter declarations 01018 CHARACTER(LEN=*), INTENT(IN) :: s_name ! Name of new dataset. 01019 INTEGER(HID_T), INTENT(IN) :: hid_group_id ! Group identifier. 01020 INTEGER :: raw_data_size 01021 LOGICAL, OPTIONAL :: lExtensible 01022 ! initialize the dataspace. 01023 01024 INTEGER(HSIZE_T) :: i_dataset_size ! Size of forthcoming data array, used to 01025 INTEGER(HID_T) :: hid_dataspace_id ! Dataspace handle. 01026 INTEGER(HID_T) :: hid_dataset_id ! Dataset handle. 01027 INTEGER(HID_T) :: hid_property_list_id ! Property list handle. 01028 INTEGER(HSIZE_T), DIMENSION(1) :: ia_dataset_dims ! Array to hold the size of the input. 01029 INTEGER(HSIZE_T), DIMENSION(1) :: ia_chunk_dims ! Array to hold the size of the dataset chunks. 01030 INTEGER(HSIZE_T), PARAMETER :: one_hst=1 01031 INTEGER :: i_err 01032 LOGICAL :: l_extend 01033 01034 01035 i_dataset_size = raw_data_size 01036 01037 IF (PRESENT(lExtensible)) THEN 01038 l_extend = lExtensible 01039 ELSE 01040 l_extend = .TRUE. 01041 END IF 01042 01043 ia_dataset_dims(1) = i_dataset_size 01044 01045 ! The dataset chunk size is whatever is smaller--the dataset size, or the maximum dataset 01046 ! chunk size in doubles. Under no circumstances should the chunk size be less than 1. 01047 ia_chunk_dims(1) = MAX(MIN(i_dataset_size, GetChunkSize(H5T_NATIVE_DOUBLE)), INT(1, HSIZE_T)) 01048 01049 ! Create the new dataspace. The H5S_UNLIMITED_F maxdims value is required for extensibility. 01050 IF (l_extend) THEN 01051 CALL h5screate_simple_f(I_DATASET_RANK, ia_dataset_dims, hid_dataspace_id, i_err, & 01052 (/ INT(H5S_UNLIMITED_F, HSIZE_T) /)) 01053 ELSE 01054 CALL h5screate_simple_f(I_DATASET_RANK, ia_dataset_dims, hid_dataspace_id, i_err) 01055 END IF 01056 01057 CALL CatchHDF5Error("Initialize_HDF5_Dataset_Double::h5screate_simple_f", i_err) 01058 01059 ! Create the property list for this dataset. 01060 CALL h5pcreate_f(H5P_DATASET_CREATE_F, hid_property_list_id, i_err) 01061 CALL CatchHDF5Error("Initialize_HDF5_Dataset_Double::h5pcreate_f", i_err) 01062 01063 ! Chunk the dataset for compression. 01064 CALL h5pset_chunk_f(hid_property_list_id, I_DATASET_RANK, ia_chunk_dims, i_err) 01065 CALL CatchHDF5Error("Initialize_HDF5_Dataset_Double::h5pset_chunk_f", i_err) 01066 01067 ! Forces the dataset to allocate the space at creation time. 01068 CALL h5pset_alloc_time_f(hid_property_list_id, H5D_ALLOC_TIME_EARLY_F, i_err) 01069 CALL CatchHDF5Error("Initialize_HDF5_Dataset_Double::h5pset_alloc_time_f", i_err) 01070 01071 ! Set the dataset to use gzip compression (default for deflate option). 01072 ! CALL h5pset_deflate_f(hid_property_list_id,I_DEFLATE_LEVEL,i_err) 01073 ! CALL CatchHDF5Error("Initialize_HDF5_Dataset_Double::h5pset_alloc_time_f", i_err) 01074 01075 ! Create the new dataset. 01076 CALL h5dcreate_f(hid_group_id, s_name, H5T_NATIVE_DOUBLE, & 01077 hid_dataspace_id, hid_dataset_id, i_err, hid_property_list_id) 01078 CALL CatchHDF5Error("Initialize_HDF5_Dataset_Double::h5dcreate_f", i_err) 01079 01080 CALL h5dclose_f(hid_dataset_id, i_err) ! Close the dataset. 01081 CALL CatchHDF5Error("Initialize_HDF5_Dataset_Double::h5dclose_f", i_err) 01082 01083 CALL h5pclose_f(hid_property_list_id, i_err) ! Close the property list. 01084 CALL CatchHDF5Error("Initialize_HDF5_Dataset_Double::h5pclose_f", i_err) 01085 01086 CALL h5sclose_f(hid_dataspace_id, i_err) ! Close the dataspace. 01087 CALL CatchHDF5Error("Initialize_HDF5_Dataset_Double::h5sclose_f", i_err) 01088 01089 END SUBROUTINE Initialize_HDF5_Dataset_Double 01090 01091 01097 SUBROUTINE Write_Slab_To_Dataset_Double(s_dataset_name, hid_group_id, dbla_data, i_offset) 01098 01099 CHARACTER(LEN=*), INTENT(IN) :: s_dataset_name 01100 INTEGER(HID_T) :: hid_group_id 01101 REAL(KIND=qPrec), DIMENSION(:) :: dbla_data ! Data input 01102 INTEGER(HSIZE_T) :: i_offset 01103 01104 INTEGER(HID_T) :: hid_dataset_id 01105 INTEGER(HID_T) :: hid_dataspace_id 01106 INTEGER(HID_T) :: hid_memspace_id 01107 INTEGER(HSIZE_T), DIMENSION(1) :: ia_data_offset 01108 INTEGER(HSIZE_T), DIMENSION(1) :: ia_slab_size 01109 INTEGER(HSIZE_T), DIMENSION(1) :: dims 01110 INTEGER(HSIZE_T), DIMENSION(1) :: maxdims 01111 INTEGER :: i_err 01112 01113 01114 i_err = 0 01115 01116 ! Set the size of the hyperslab. 01117 ia_slab_size(1) = SIZE(dbla_data) 01118 01119 ! Store the size of the offset in a dimension array. 01120 ia_data_offset(1) = i_offset 01121 01122 ! Open the dataset in the HDF5 file. 01123 CALL h5dopen_f(hid_group_id, s_dataset_name, hid_dataset_id, i_err) 01124 CALL CatchHDF5Error("Write_Slab_To_Dataset_Double::h5dopen_f", i_err) 01125 01126 CALL h5dget_space_f(hid_dataset_id, hid_dataspace_id, i_err) 01127 CALL CatchHDF5Error("Write_Slab_To_Dataset_Double::h5dget_space_f", i_err) 01128 01129 ! Extract the dimensions of said dataspace. 01130 CALL h5sget_simple_extent_dims_f(hid_dataspace_id, dims, maxdims, i_err) 01131 CALL CatchHDF5Error("Write_Slab_To_Dataset_Double::h5sget_simple_extent_dims_f", i_err, I_DATASET_RANK) 01132 01133 ! If the dataset is extensible and is about to overflow, then extend the dataset 01134 ! to accomodate the new size. 01135 IF ((maxdims(1) == H5S_UNLIMITED_F) .AND. (ia_slab_size(1) + ia_data_offset(1) > dims(1))) THEN 01136 01137 ! Extend the dataspace. 01138 CALL h5sset_extent_simple_f(hid_dataspace_id, I_DATASET_RANK, & 01139 ia_slab_size + ia_data_offset, maxdims, i_err) 01140 01141 CALL CatchHDF5Error("Write_Slab_To_Dataset_Box::h5sset_extent_f", i_err) 01142 01143 CALL h5dset_extent_f(hid_dataset_id, ia_slab_size + ia_data_offset, i_err) 01144 CALL CatchHDF5Error("Write_Slab_To_Dataset_Box::h5dset_extent_f", i_err) 01145 01146 END IF 01147 01148 ! Extract the dimensions of said dataspace. 01149 CALL h5sget_simple_extent_dims_f(hid_dataspace_id, dims, maxdims, i_err) 01150 CALL CatchHDF5Error("Write_Slab_To_Dataset_Double::h5sget_simple_extent_dims_f", i_err, I_DATASET_RANK) 01151 01152 ! Select hyperslab in the dataset (specifically, the dimensions of the input array). 01153 CALL h5sselect_hyperslab_f(hid_dataspace_id, H5S_SELECT_SET_F, & 01154 ia_data_offset, ia_slab_size, i_err) 01155 CALL CatchHDF5Error("Write_Slab_To_Dataset_Double::h5sselect_hyperslab_f", i_err) 01156 01157 ! Create dataspace for the upcoming hyperslab. 01158 CALL h5screate_simple_f(I_DATASET_RANK, ia_slab_size, hid_memspace_id, i_err) 01159 CALL CatchHDF5Error("Write_Slab_To_Dataset_Double::h5screate_simple_f", i_err) 01160 01161 ! Write the HDF5 data to the dataset. 01162 CALL h5dwrite_f(hid_dataset_id, H5T_NATIVE_DOUBLE, dbla_data, & 01163 ia_slab_size, i_err, hid_memspace_id, hid_dataspace_id) 01164 01165 CALL CatchHDF5Error("Write_Slab_To_Dataset_Double::h5dwrite_f", i_err) 01166 01167 CALL h5sclose_f(hid_memspace_id,i_err) 01168 CALL CatchHDF5Error("Write_Slab_To_Dataset_Double::h5sclose_f",i_err) 01169 01170 ! Close main file memory dataspace. 01171 CALL h5sclose_f(hid_dataspace_id, i_err) 01172 CALL CatchHDF5Error("Write_Slab_To_Dataset_Double::h5sclose_f", i_err) 01173 01174 CALL h5dclose_f(hid_dataset_id, i_err) 01175 CALL CatchHDF5Error("Write_Slab_To_Dataset_Double::h5dclose_f", i_err) 01176 01177 END SUBROUTINE Write_Slab_To_Dataset_Double 01178 01184 SUBROUTINE Initialize_HDF5_Dataset_Box(s_name, hid_group_id, raw_data_size, lExtensible) 01185 01186 IMPLICIT NONE 01187 01188 ! Input parameter declarations 01189 CHARACTER(LEN=*), INTENT(IN) :: s_name ! Name of new dataset. 01190 INTEGER :: raw_data_size 01191 INTEGER(HID_T), INTENT(IN) :: hid_group_id ! Group identifier. 01192 ! initialize the dataspace. 01193 LOGICAL, OPTIONAL :: lExtensible 01194 01195 INTEGER(HSIZE_T) :: i_dataset_size ! Size of forthcoming data array, used to 01196 INTEGER(HID_T) :: hid_dataspace_id ! Dataspace handle. 01197 INTEGER(HID_T) :: hid_dataset_id ! Dataset handle. 01198 INTEGER(HID_T) :: hid_property_list_id ! Property list handle. 01199 INTEGER(HSIZE_T), DIMENSION(1) :: ia_dataset_dims ! Array to hold the size of the input. 01200 INTEGER(HSIZE_T), DIMENSION(1) :: ia_chunk_dims 01201 LOGICAL :: l_extend 01202 INTEGER :: i_err 01203 01204 01205 i_dataset_size = raw_data_size 01206 01207 IF (PRESENT(lExtensible)) THEN 01208 l_extend = lExtensible 01209 ELSE 01210 l_extend = .TRUE. 01211 END IF 01212 01213 ia_dataset_dims(1) = i_dataset_size 01214 01215 ! The dataset chunk size is whatever is smaller--the dataset size, or the maximum dataset 01216 ! chunk size in doubles. Under no circumstances should the chunk size be less than 1. 01217 ia_chunk_dims(1) = MAX(MIN(i_dataset_size, GetChunkSize(hid_box_id)), INT(1, HSIZE_T)) 01218 01219 ! Create the new dataspace. The H5S_UNLIMITED_F maxdims value is required for extensibility. 01220 IF (l_extend) THEN 01221 CALL h5screate_simple_f(I_DATASET_RANK, ia_dataset_dims, hid_dataspace_id, i_err, & 01222 (/ INT(H5S_UNLIMITED_F, HSIZE_T) /)) 01223 ELSE 01224 CALL h5screate_simple_f(I_DATASET_RANK, ia_dataset_dims, hid_dataspace_id, i_err) 01225 END IF 01226 01227 CALL CatchHDF5Error("Initialize_HDF5_Dataset_Box::h5screate_simple_f", i_err) 01228 01229 ! Create the property list for this dataset. 01230 CALL h5pcreate_f(H5P_DATASET_CREATE_F, hid_property_list_id, i_err) 01231 CALL CatchHDF5Error("Initialize_HDF5_Dataset_Box::h5pcreate_f", i_err) 01232 01233 ! Chunk the dataset for compression. 01234 CALL h5pset_chunk_f(hid_property_list_id, I_DATASET_RANK, ia_chunk_dims, i_err) 01235 CALL CatchHDF5Error("Initialize_HDF5_Dataset_Box::h5pset_chunk_f", i_err) 01236 01237 ! Forces the dataset to allocate the space at creation time. 01238 CALL h5pset_alloc_time_f(hid_property_list_id, H5D_ALLOC_TIME_EARLY_F, i_err) 01239 CALL CatchHDF5Error("Initialize_HDF5_Dataset_Box::h5pset_alloc_time_f", i_err) 01240 01241 ! Set the dataset to use gzip compression (default for deflate option). 01242 ! CALL h5pset_deflate_f(hid_property_list_id,I_DEFLATE_LEVEL,i_err) 01243 ! CALL CatchHDF5Error("Initialize_HDF5_Dataset_Box::h5pset_alloc_time_f", i_err) 01244 01245 ! Create the new dataset. 01246 CALL h5dcreate_f(hid_group_id, s_name, hid_box_id, & 01247 hid_dataspace_id, hid_dataset_id, i_err, hid_property_list_id) 01248 CALL CatchHDF5Error("Initialize_HDF5_Dataset_Box::h5dcreate_f", i_err) 01249 01250 CALL h5dclose_f(hid_dataset_id, i_err) ! Close the dataset. 01251 CALL CatchHDF5Error("Initialize_HDF5_Dataset_Box::h5dclose_f", i_err) 01252 01253 CALL h5pclose_f(hid_property_list_id, i_err) ! Close the property list. 01254 CALL CatchHDF5Error("Initialize_HDF5_Dataset_Box::h5pclose_f", i_err) 01255 01256 CALL h5sclose_f(hid_dataspace_id, i_err) ! Close the dataspace. 01257 CALL CatchHDF5Error("Initialize_HDF5_Dataset_Box::h5sclose_f", i_err) 01258 01259 END SUBROUTINE Initialize_HDF5_Dataset_Box 01260 01261 01267 SUBROUTINE Write_Slab_To_Dataset_Box(s_dataset_name, hid_group_id, box_data, i_offset) 01268 01269 CHARACTER(LEN=*), INTENT(IN) :: s_dataset_name 01270 INTEGER(HID_T) :: hid_group_id 01271 INTEGER, DIMENSION(I_CHOMBO_BOX_SIZE) :: box_data ! Data input 01272 INTEGER(HSIZE_T) :: i_offset 01273 01274 INTEGER(HID_T) :: hid_dataset_id 01275 INTEGER(HID_T) :: hid_dataspace_id 01276 INTEGER(HID_T) :: hid_memspace_id 01277 INTEGER(HSIZE_T), DIMENSION(1) :: ia_data_offset 01278 INTEGER(HSIZE_T), DIMENSION(1) :: ia_slab_size 01279 INTEGER(HSIZE_T), DIMENSION(1) :: dims 01280 INTEGER(HSIZE_T), DIMENSION(1) :: maxdims 01281 INTEGER :: i_err 01282 01283 01284 i_err = 0 01285 01286 ! Set the size of the hyperslab. 01287 ia_slab_size(1) = SIZE(box_data) / I_CHOMBO_BOX_SIZE 01288 01289 ! Store the size of the offset in a dimension array. 01290 ia_data_offset(1) = i_offset 01291 01292 ! Open the dataset in the HDF5 file. 01293 CALL h5dopen_f(hid_group_id, s_dataset_name, hid_dataset_id, i_err) 01294 CALL CatchHDF5Error("Write_Slab_To_Dataset_Box::h5dopen_f", i_err) 01295 01296 CALL h5dget_space_f(hid_dataset_id, hid_dataspace_id, i_err) 01297 CALL CatchHDF5Error("Write_Slab_To_Dataset_Box::h5dget_space_f", i_err) 01298 01299 ! Select hyperslab in the dataset (specifically, the dimensions of the input array). 01300 CALL h5sselect_hyperslab_f(hid_dataspace_id, H5S_SELECT_SET_F, & 01301 ia_data_offset, ia_slab_size, i_err) 01302 CALL CatchHDF5Error("Write_Slab_To_Dataset_Box::h5sselect_hyperslab_f", i_err) 01303 01304 ! Create dataspace for the upcoming hyperslab. 01305 CALL h5screate_simple_f(I_DATASET_RANK, ia_slab_size, hid_memspace_id, i_err) 01306 CALL CatchHDF5Error("Write_Slab_To_Dataset_Box::h5screate_simple_f", i_err) 01307 01308 ! Extract the dimensions of said dataspace. 01309 CALL h5sget_simple_extent_dims_f(hid_dataspace_id, dims, maxdims, i_err) 01310 CALL CatchHDF5Error("Write_Slab_To_Dataset_Box::h5sget_simple_extent_dims_f", i_err, I_DATASET_RANK) 01311 01312 ! If the dataset is extensible and is about to overflow, then extend the dataset 01313 ! to accomodate the new size. 01314 IF ((maxdims(1) == H5S_UNLIMITED_F) .AND. (ia_slab_size(1) + ia_data_offset(1) > dims(1))) THEN 01315 01316 ! Extend the dataspace. 01317 CALL h5sset_extent_simple_f(hid_dataspace_id, I_DATASET_RANK, & 01318 ia_slab_size + ia_data_offset, maxdims, i_err) 01319 01320 CALL CatchHDF5Error("Write_Slab_To_Dataset_Box::h5sset_extent_f", i_err) 01321 01322 CALL h5dset_extent_f(hid_dataset_id, ia_slab_size + ia_data_offset, i_err) 01323 CALL CatchHDF5Error("Write_Slab_To_Dataset_Box::h5dset_extent_f", i_err) 01324 01325 END IF 01326 01327 ! Extract the dimensions of said dataspace. 01328 CALL h5sget_simple_extent_dims_f(hid_dataspace_id, dims, maxdims, i_err) 01329 CALL CatchHDF5Error("Write_Slab_To_Dataset_Box::h5sget_simple_extent_dims_f", i_err, I_DATASET_RANK) 01330 01331 ! Write the HDF5 data to the dataset. 01332 CALL h5dwrite_f(hid_dataset_id, hid_box_id, box_data, & 01333 ia_slab_size, i_err, hid_memspace_id, hid_dataspace_id) 01334 CALL CatchHDF5Error("Write_Slab_To_Dataset_Box::h5dwrite_f", i_err) 01335 01336 CALL h5sclose_f(hid_memspace_id,i_err) 01337 CALL CatchHDF5Error("Write_Slab_To_Dataset_Box::h5sclose_f", i_err) 01338 01339 ! Close main file memory dataspace. 01340 CALL h5sclose_f(hid_dataspace_id, i_err) 01341 CALL CatchHDF5Error("Write_Slab_To_Dataset_Box::h5sclose_f", i_err) 01342 01343 CALL h5dclose_f(hid_dataset_id, i_err) 01344 CALL CatchHDF5Error("Write_Slab_To_Dataset_Box::h5dclose_f", i_err) 01345 01346 END SUBROUTINE Write_Slab_To_Dataset_Box 01347 01353 SUBROUTINE Initialize_HDF5_Dataset_Int(s_name, hid_group_id, raw_data_size, lExtensible) 01354 01355 IMPLICIT NONE 01356 01357 ! Input parameter declarations 01358 CHARACTER(LEN=*), INTENT(IN) :: s_name ! Name of new dataset. 01359 INTEGER(HID_T), INTENT(IN) :: hid_group_id ! Group identifier. 01360 INTEGER :: raw_data_size 01361 ! initialize the dataspace. 01362 LOGICAL, OPTIONAL :: lExtensible 01363 01364 INTEGER(HSIZE_T) :: i_dataset_size ! Size of forthcoming data array, used to 01365 INTEGER(HID_T) :: hid_dataspace_id ! Dataspace handle. 01366 INTEGER(HID_T) :: hid_dataset_id ! Dataset handle. 01367 INTEGER(HID_T) :: hid_property_list_id ! Property list handle. 01368 INTEGER(HSIZE_T), DIMENSION(1) :: ia_dataset_dims ! Array to hold the size of the input. 01369 INTEGER(HSIZE_T), DIMENSION(1) :: ia_chunk_dims 01370 LOGICAL :: l_extend 01371 INTEGER :: i_err 01372 01373 01374 i_dataset_size = raw_data_size 01375 01376 IF (PRESENT(lExtensible)) THEN 01377 l_extend = lExtensible 01378 ELSE 01379 l_extend = .TRUE. 01380 END IF 01381 01382 ia_dataset_dims(1) = i_dataset_size 01383 01384 ! The dataset chunk size is whatever is smaller--the dataset size, or the maximum dataset 01385 ! chunk size in doubles. Under no circumstances should the chunk size be less than 1. 01386 ia_chunk_dims(1) = MAX(MIN(i_dataset_size, GetChunkSize(H5T_NATIVE_INTEGER)), INT(1, HSIZE_T)) 01387 01388 ! Create the new dataspace. The H5S_UNLIMITED_F maxdims value is required for extensibility. 01389 IF (l_extend) THEN 01390 CALL h5screate_simple_f(I_DATASET_RANK, ia_dataset_dims, hid_dataspace_id, i_err, & 01391 (/ INT(H5S_UNLIMITED_F, HSIZE_T) /)) 01392 ELSE 01393 CALL h5screate_simple_f(I_DATASET_RANK, ia_dataset_dims, hid_dataspace_id, i_err) 01394 END IF 01395 01396 CALL CatchHDF5Error("Initialize_HDF5_Dataset_Int::h5screate_simple_f", i_err) 01397 01398 ! Create the property list for this dataset. 01399 CALL h5pcreate_f(H5P_DATASET_CREATE_F, hid_property_list_id, i_err) 01400 CALL CatchHDF5Error("Initialize_HDF5_Dataset_Int::h5pcreate_f", i_err) 01401 01402 ! Chunk the dataset for compression. 01403 CALL h5pset_chunk_f(hid_property_list_id, I_DATASET_RANK, ia_chunk_dims, i_err) 01404 CALL CatchHDF5Error("Initialize_HDF5_Dataset_Int::h5pset_chunk_f", i_err) 01405 01406 ! Forces the dataset to allocate the space at creation time. 01407 CALL h5pset_alloc_time_f(hid_property_list_id, H5D_ALLOC_TIME_EARLY_F, i_err) 01408 CALL CatchHDF5Error("Initialize_HDF5_Dataset_Int::h5pset_alloc_time_f", i_err) 01409 01410 ! Set the dataset to use gzip compression (default for deflate option). 01411 ! CALL h5pset_deflate_f(hid_property_list_id,I_DEFLATE_LEVEL,i_err) 01412 ! CALL CatchHDF5Error("Initialize_HDF5_Dataset_Int::h5pset_alloc_time_f", i_err) 01413 01414 ! Create the new dataset. 01415 CALL h5dcreate_f(hid_group_id, s_name, H5T_NATIVE_INTEGER, & 01416 hid_dataspace_id, hid_dataset_id, i_err, hid_property_list_id) 01417 CALL CatchHDF5Error("Initialize_HDF5_Dataset_Int::h5dcreate_f", i_err) 01418 01419 CALL h5dclose_f(hid_dataset_id, i_err) ! Close the dataset. 01420 CALL CatchHDF5Error("Initialize_HDF5_Dataset_Int::h5dclose_f", i_err) 01421 01422 CALL h5pclose_f(hid_property_list_id, i_err) ! Close the property list. 01423 CALL CatchHDF5Error("Initialize_HDF5_Dataset_Int::h5pclose_f", i_err) 01424 01425 CALL h5sclose_f(hid_dataspace_id, i_err) ! Close the dataspace. 01426 CALL CatchHDF5Error("Initialize_HDF5_Dataset_Int::h5sclose_f", i_err) 01427 01428 END SUBROUTINE Initialize_HDF5_Dataset_Int 01429 01435 SUBROUTINE Write_Slab_To_Dataset_Int(s_dataset_name, hid_group_id, int_data, i_offset) 01436 01437 CHARACTER(LEN=*), INTENT(IN) :: s_dataset_name 01438 INTEGER(HID_T) :: hid_group_id 01439 INTEGER, DIMENSION(:) :: int_data ! Data input 01440 INTEGER(HSIZE_T) :: i_offset 01441 01442 INTEGER(HID_T) :: hid_dataset_id 01443 INTEGER(HID_T) :: hid_dataspace_id 01444 INTEGER(HID_T) :: hid_memspace_id 01445 INTEGER(HSIZE_T), DIMENSION(1) :: ia_data_offset 01446 INTEGER(HSIZE_T), DIMENSION(1) :: ia_slab_size 01447 INTEGER(HSIZE_T), DIMENSION(1) :: dims 01448 INTEGER(HSIZE_T), DIMENSION(1) :: maxdims 01449 INTEGER :: i_err 01450 01451 01452 i_err = 0 01453 01454 ! Set the size of the hyperslab. 01455 ia_slab_size(1) = SIZE(int_data) 01456 01457 ! Store the size of the offset in a dimension array. 01458 ia_data_offset(1) = i_offset 01459 01460 ! Open the dataset in the HDF5 file. 01461 CALL h5dopen_f(hid_group_id, s_dataset_name, hid_dataset_id, i_err) 01462 CALL CatchHDF5Error("Write_Slab_To_Dataset_Int::h5dopen_f", i_err) 01463 01464 CALL h5dget_space_f(hid_dataset_id, hid_dataspace_id, i_err) 01465 CALL CatchHDF5Error("Write_Slab_To_Dataset_Int::h5dget_space_f", i_err) 01466 01467 ! Extract the dimensions of said dataspace. 01468 CALL h5sget_simple_extent_dims_f(hid_dataspace_id, dims, maxdims, i_err) 01469 CALL CatchHDF5Error("Write_Slab_To_Dataset_Int::h5sget_simple_extent_dims_f", i_err, I_DATASET_RANK) 01470 01471 ! If the dataset is extensible and is about to overflow, then extend the dataset 01472 ! to accomodate the new size. 01473 IF ((maxdims(1) == H5S_UNLIMITED_F) .AND. (ia_slab_size(1) + ia_data_offset(1) > dims(1))) THEN 01474 01475 ! Extend the dataspace. 01476 CALL h5sset_extent_simple_f(hid_dataspace_id, I_DATASET_RANK, & 01477 ia_slab_size + ia_data_offset, maxdims, i_err) 01478 01479 CALL CatchHDF5Error("Write_Slab_To_Dataset_Box::h5sset_extent_f", i_err) 01480 01481 CALL h5dset_extent_f(hid_dataset_id, ia_slab_size + ia_data_offset, i_err) 01482 CALL CatchHDF5Error("Write_Slab_To_Dataset_Box::h5dset_extent_f", i_err) 01483 01484 END IF 01485 01486 ! Select hyperslab in the dataset (specifically, the dimensions of the input array). 01487 CALL h5sselect_hyperslab_f(hid_dataspace_id, H5S_SELECT_SET_F, & 01488 ia_data_offset, ia_slab_size, i_err) 01489 CALL CatchHDF5Error("Write_Slab_To_Dataset_Int::h5sselect_hyperslab_f", i_err) 01490 01491 ! Create dataspace for the upcoming hyperslab. 01492 CALL h5screate_simple_f(I_DATASET_RANK, ia_slab_size, hid_memspace_id, i_err) 01493 CALL CatchHDF5Error("Write_Slab_To_Dataset_Int::h5screate_simple_f", i_err) 01494 01495 ! Write the HDF5 data to the dataset. 01496 CALL h5dwrite_f(hid_dataset_id, H5T_NATIVE_INTEGER, int_data, & 01497 ia_slab_size, i_err, hid_memspace_id, hid_dataspace_id) 01498 CALL CatchHDF5Error("Write_Slab_To_Dataset_Int::h5dwrite_f", i_err) 01499 01500 CALL h5sclose_f(hid_memspace_id,i_err) 01501 CALL CatchHDF5Error("Write_Slab_To_Dataset_Int::h5sclose_f", i_err) 01502 01503 ! Close main file memory dataspace. 01504 CALL h5sclose_f(hid_dataspace_id, i_err) 01505 CALL CatchHDF5Error("Write_Slab_To_Dataset_Int::h5sclose_f", i_err) 01506 01507 CALL h5dclose_f(hid_dataset_id, i_err) 01508 CALL CatchHDF5Error("Write_Slab_To_Dataset_Int::h5dclose_f", i_err) 01509 01510 END SUBROUTINE Write_Slab_To_Dataset_Int 01511 01514 01517 INTEGER(HSIZE_T) FUNCTION GetChunkSize(hid_type_id) 01518 01519 INTEGER(HID_T) :: hid_type_id 01520 01521 INTEGER(size_t) :: type_size 01522 INTEGER :: i_err 01523 01524 01525 ! Get the type's size in bytes. 01526 CALL h5tget_size_f(hid_type_id, type_size, i_err) 01527 CALL CatchHDF5Error("GetChunkSize::h5tget_size_f", i_err) 01528 01529 ! Obtain the size of the chunk in type elements. 01530 GetChunkSize = FLOOR(I_MAX_CHUNK_SIZE * 1.0 / type_size) 01531 01532 END FUNCTION GetChunkSize 01533 01534 END MODULE HDF5Declarations 01535