1! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
2!   Copyright by The HDF Group.                                               *
3!   Copyright by the Board of Trustees of the University of Illinois.         *
4!   All rights reserved.                                                      *
5!                                                                             *
6!   This file is part of HDF5.  The full HDF5 copyright notice, including     *
7!   terms governing use, modification, and redistribution, is contained in    *
8!   the COPYING file, which can be found at the root of the source code       *
9!   distribution tree, or in https://support.hdfgroup.org/ftp/HDF5/releases.  *
10!   If you do not have access to either file, you may request a copy from     *
11!   help@hdfgroup.org.                                                        *
12! * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
13!
14!
15! The following example shows how to write and read to/from an existing dataset.
16! It opens the file created in the previous example, obtains the dataset
17! identifier, writes the data to the dataset in the file,
18! then reads the dataset  to memory.
19!
20! This example is used in the HDF5 Tutorial.
21
22PROGRAM H5_RDWT
23
24  USE HDF5 ! This module contains all necessary modules
25
26  IMPLICIT NONE
27
28  CHARACTER(LEN=8), PARAMETER :: filename = "dsetf.h5" ! File name
29  CHARACTER(LEN=4), PARAMETER :: dsetname = "dset"     ! Dataset name
30
31  INTEGER(HID_T) :: file_id       ! File identifier
32  INTEGER(HID_T) :: dset_id       ! Dataset identifier
33
34  INTEGER     ::   error ! Error flag
35  INTEGER     ::  i, j
36
37  INTEGER, DIMENSION(4,6) :: dset_data, data_out ! Data buffers
38  INTEGER(HSIZE_T), DIMENSION(2) :: data_dims
39
40  !
41  ! Initialize the dset_data array.
42  !
43  DO i = 1, 4
44     DO j = 1, 6
45        dset_data(i,j) = (i-1)*6 + j
46     END DO
47  END DO
48
49  !
50  ! Initialize FORTRAN interface.
51  !
52  CALL h5open_f(error)
53
54  !
55  ! Open an existing file.
56  !
57  CALL h5fopen_f (filename, H5F_ACC_RDWR_F, file_id, error)
58
59  !
60  ! Open an existing dataset.
61  !
62  CALL h5dopen_f(file_id, dsetname, dset_id, error)
63
64  !
65  ! Write the dataset.
66  !
67  data_dims(1) = 4
68  data_dims(2) = 6
69  CALL h5dwrite_f(dset_id, H5T_NATIVE_INTEGER, dset_data, data_dims, error)
70
71  !
72  ! Read the dataset.
73  !
74  CALL h5dread_f(dset_id, H5T_NATIVE_INTEGER, data_out, data_dims, error)
75
76  !
77  ! Close the dataset.
78  !
79  CALL h5dclose_f(dset_id, error)
80
81  !
82  ! Close the file.
83  !
84  CALL h5fclose_f(file_id, error)
85
86  !
87  ! Close FORTRAN interface.
88  !
89  CALL h5close_f(error)
90
91END PROGRAM H5_RDWT
92
93
94
95