TABLE OF CONTENTS
etsf_io_data_contents
[ Top ] [ etsf_io_data_group ] [ Methods ]
NAME
etsf_io_data_contents
FUNCTION
High-level routine that get informations from a given @filename. Returned values are the list of dimensions, allocated split definitions (if any), flags for main variables (see FLAGS_VARIABLES) and flags for groups (see FLAGS_GROUPS).
This routine can also be used to get the comprehensive list of read variables with their definitions (name, shape, dimension names...). Use etsf_io_vars_free() to deallocate this list.
INPUTS
- filename = the path to the file to be accessed.
OUTPUT
- dims <type(etsf_dims)> = the dimensions will be read and stored, using etsf_io_dims_get().
- split <type(etsf_split)> = if any, read the split array from the given file and put their values in this argument. If lstat = .true., it may be allocated in output. So, after use, it must be deallocated, using etsf_io_split_free().
- etsf_groups = an integer which is the sum of all present group ids in the read file (see FLAGS_GROUPS).
- etsf_variables <type(etsf_groups_flags)> = an integer for each group detailling which ETSF variables are indeed present in the read file (see FLAGS_VARIABLES).
- lstat = return .true. if all the actions succeed, if not the status of the file is undefined.
- error_data <type(etsf_io_low_error)> = contains the details of the error is @lstat is false.
SIDE EFFECTS
- vars_infos <type(etsf_vars)> = (optional) when reading the file for variable informations, it creates a list of type(etsf_var) that describes all the variables in the file. This list contains non ETSF informations such as variable names, types, shapes, an array of dimension values and an other array of dimension names. It also contains ETSF informations, like group id or if the variable is a split definition. When given, internal pointers are associated in the subroutine. To free them, use etsf_io_vars_free().
NOTES
This file has been automatically generated by the autogen_subroutines.py script. Any change you would bring to it will systematically be overwritten.
SOURCE
subroutine etsf_io_data_contents(filename, dims, split, etsf_groups, etsf_variables, & & lstat, error_data, vars_infos) !Arguments ------------------------------------ character(len=*), intent(in) :: filename type(etsf_dims), intent(out) :: dims type(etsf_split), intent(out) :: split integer, intent(out) :: etsf_groups type(etsf_groups_flags), intent(out) :: etsf_variables logical, intent(out) :: lstat type(etsf_io_low_error), intent(out) :: error_data type(etsf_vars), optional, intent(inout) :: vars_infos !Local variables------------------------------- character(len=*),parameter :: my_name = 'etsf_io_data_contents' integer :: ncid, i type(etsf_vars) :: my_vars_infos logical :: with_dim_name logical :: with_att_name integer :: group_id integer :: var_id logical :: split_id ! ************************************************************************* !DEBUG !write (*,*) 'etsf_io_data_contents : enter' !ENDDEBUG lstat = .false. if (present(vars_infos)) then vars_infos%n_vars = 0 vars_infos%parent => null() with_dim_name = .true. with_att_name = .true. else with_dim_name = .false. with_att_name = .false. end if ! Open file for reading call etsf_io_low_open_read(ncid, trim(filename), & & lstat, error_data = error_data) if (.not. lstat) then call etsf_io_low_error_update(error_data, my_name) return end if ! Get the dimensions. call etsf_io_dims_get(ncid, dims, lstat, error_data) if (.not. lstat) then call etsf_io_low_error_update(error_data, my_name) return end if ! We allocate the split arrays. call etsf_io_split_allocate(split, dims) ! We read the split informations. call etsf_io_split_get(ncid, split, lstat, error_data) if (.not. lstat) then call etsf_io_split_free(split) call etsf_io_low_error_update(error_data, my_name) return end if ! Get all variables definitions. ! It will allocate my_vars_infos%parent array. call etsf_io_low_read_all_var_infos(ncid, my_vars_infos%parent, & & lstat, error_data, with_dim_name = with_dim_name, & & with_att_name = with_att_name) if (.not. lstat) then call etsf_io_split_free(split) call etsf_io_low_error_update(error_data, my_name) return end if etsf_variables%geometry = etsf_geometry_none etsf_variables%electrons = etsf_electrons_none etsf_variables%kpoints = etsf_kpoints_none etsf_variables%basisdata = etsf_basisdata_none etsf_variables%gwdata = etsf_gwdata_none etsf_variables%dielectric = etsf_dielectric_none etsf_variables%main = etsf_main_none if (associated(my_vars_infos%parent)) then my_vars_infos%n_vars = size(my_vars_infos%parent) if (present(vars_infos)) then ! Allocate vars_infos arrays for future use. vars_infos%n_vars = my_vars_infos%n_vars allocate(vars_infos%group(vars_infos%n_vars)) allocate(vars_infos%varid(vars_infos%n_vars)) allocate(vars_infos%split(vars_infos%n_vars)) end if ! get the main_id and the group_id for all variables. do i = 1, my_vars_infos%n_vars, 1 call etsf_io_data_get(group_id, var_id, & & split_id, my_vars_infos%parent(i)%name) select case (group_id) case (etsf_grp_geometry) etsf_variables%geometry = ior(etsf_variables%geometry, var_id) case (etsf_grp_electrons) etsf_variables%electrons = ior(etsf_variables%electrons, var_id) case (etsf_grp_kpoints) etsf_variables%kpoints = ior(etsf_variables%kpoints, var_id) case (etsf_grp_basisdata) etsf_variables%basisdata = ior(etsf_variables%basisdata, var_id) case (etsf_grp_gwdata) etsf_variables%gwdata = ior(etsf_variables%gwdata, var_id) case (etsf_grp_dielectric) etsf_variables%dielectric = ior(etsf_variables%dielectric, var_id) case (etsf_grp_main) etsf_variables%main = ior(etsf_variables%main, var_id) end select etsf_groups = ior(etsf_groups, group_id) if (present(vars_infos)) then ! Update vars_infos arrays. vars_infos%group(i) = group_id vars_infos%varid(i) = var_id vars_infos%split(i) = split_id end if end do end if if (present(vars_infos)) then ! Associate vars_infos%parent to the one computed in my_vars_infos. vars_infos%parent => my_vars_infos%parent else if (associated(my_vars_infos%parent)) then call etsf_io_low_free_all_var_infos(my_vars_infos%parent) end if ! Close file call etsf_io_low_close(ncid, lstat, error_data = error_data) if (.not. lstat) call etsf_io_low_error_update(error_data, my_name) !DEBUG !write (*,*) 'etsf_io_data_contents : exit' !ENDDEBUG end subroutine etsf_io_data_contents