TABLE OF CONTENTS


etsf_io_file_check_crystallographic_data

[ Top ] [ etsf_io_file ] [ Methods ]

NAME

etsf_io_file_check_crystallographic_data

FUNCTION

This is a high level routine to inquire a file about a specifications. It returns .true. in lstat if the file is a valid 'crystallographic_data' file. It tests the existence of variables and their definition (type, shape. and dimension names).

INPUTS

OUTPUT

SOURCE

subroutine etsf_io_file_check_crystallographic_data(ncid, lstat, error_data)
  integer, intent(in)                  :: ncid
  logical, intent(out)                 :: lstat
  type(etsf_io_low_error), intent(out) :: error_data

  character(len = *), parameter        :: me = "etsf_io_file_check_crystallographic_data"
  type(etsf_io_low_var_infos)          :: var_infos
  logical                              :: valid
  character(len = etsf_charlen)        :: string_value
  type(etsf_dims)                      :: dims
  type(etsf_split)                     :: split

  ! Read the dimensions
  call etsf_io_dims_get(ncid, dims, lstat, error_data)
  if (.not. lstat) then
     call etsf_io_low_error_update(error_data, me)
     return
  end if

  ! Allocate the split and read it (this will verify variable exist.
  call etsf_io_split_allocate(split, dims)
  call etsf_io_split_get(ncid, split, lstat, error_data)
  if (.not. lstat) then
     call etsf_io_low_error_update(error_data, me)
     return
  end if

  ! Variable primitive_vectors
  write(var_infos%name, "(A)") "primitive_vectors"
  var_infos%nctype  = etsf_io_low_double
  var_infos%ncshape = 2
  allocate(var_infos%ncdimnames(2))
  write(var_infos%ncdimnames(2), "(A)") "number_of_vectors"
  write(var_infos%ncdimnames(1), "(A)") "number_of_cartesian_directions"
  call test_var(ncid, var_infos, lstat, error_data)
  deallocate(var_infos%ncdimnames)
  if (.not. lstat) then
    call etsf_io_split_free(split)
    call etsf_io_low_error_update(error_data, me)
    return
  end if
  
  ! Variable reduced_symmetry_matrices
  write(var_infos%name, "(A)") "reduced_symmetry_matrices"
  var_infos%nctype  = etsf_io_low_integer
  var_infos%ncshape = 3
  allocate(var_infos%ncdimnames(3))
  write(var_infos%ncdimnames(3), "(A)") "number_of_symmetry_operations"
  write(var_infos%ncdimnames(2), "(A)") "number_of_reduced_dimensions"
  write(var_infos%ncdimnames(1), "(A)") "number_of_reduced_dimensions"
  call test_var(ncid, var_infos, lstat, error_data)
  deallocate(var_infos%ncdimnames)
  if (.not. lstat) then
    call etsf_io_split_free(split)
    call etsf_io_low_error_update(error_data, me)
    return
  end if
  
  ! Variable reduced_symmetry_translations
  write(var_infos%name, "(A)") "reduced_symmetry_translations"
  var_infos%nctype  = etsf_io_low_double
  var_infos%ncshape = 2
  allocate(var_infos%ncdimnames(2))
  write(var_infos%ncdimnames(2), "(A)") "number_of_symmetry_operations"
  write(var_infos%ncdimnames(1), "(A)") "number_of_reduced_dimensions"
  call test_var(ncid, var_infos, lstat, error_data)
  deallocate(var_infos%ncdimnames)
  if (.not. lstat) then
    call etsf_io_split_free(split)
    call etsf_io_low_error_update(error_data, me)
    return
  end if
  
  ! Variable space_group
  write(var_infos%name, "(A)") "space_group"
  var_infos%nctype  = etsf_io_low_integer
  var_infos%ncshape = 0
  call test_var(ncid, var_infos, lstat, error_data)
  if (.not. lstat) then
    call etsf_io_split_free(split)
    call etsf_io_low_error_update(error_data, me)
    return
  end if
  
  ! Variable atom_species
  write(var_infos%name, "(A)") "atom_species"
  var_infos%nctype  = etsf_io_low_integer
  var_infos%ncshape = 1
  allocate(var_infos%ncdimnames(1))
  write(var_infos%ncdimnames(1), "(A)") "number_of_atoms"
  call test_var(ncid, var_infos, lstat, error_data)
  deallocate(var_infos%ncdimnames)
  if (.not. lstat) then
    call etsf_io_split_free(split)
    call etsf_io_low_error_update(error_data, me)
    return
  end if
  
  ! Variable reduced_atom_positions
  write(var_infos%name, "(A)") "reduced_atom_positions"
  var_infos%nctype  = etsf_io_low_double
  var_infos%ncshape = 2
  allocate(var_infos%ncdimnames(2))
  write(var_infos%ncdimnames(2), "(A)") "number_of_atoms"
  write(var_infos%ncdimnames(1), "(A)") "number_of_reduced_dimensions"
  call test_var(ncid, var_infos, lstat, error_data)
  deallocate(var_infos%ncdimnames)
  if (.not. lstat) then
    call etsf_io_split_free(split)
    call etsf_io_low_error_update(error_data, me)
    return
  end if
  
  ! Check from a list.
  lstat = .false.
  ! Variable atomic_numbers
  write(var_infos%name, "(A)") "atomic_numbers"
  var_infos%nctype  = etsf_io_low_double
  var_infos%ncshape = 1
  allocate(var_infos%ncdimnames(1))
  write(var_infos%ncdimnames(1), "(A)") "number_of_atom_species"
  call test_var(ncid, var_infos, valid, error_data)
  deallocate(var_infos%ncdimnames)
  if (.not. valid .and. error_data%access_mode_id == ERROR_MODE_SPEC) return
  lstat = lstat .or. valid
  ! Variable atom_species_names
  write(var_infos%name, "(A)") "atom_species_names"
  var_infos%nctype  = etsf_io_low_character
  var_infos%ncshape = 2
  allocate(var_infos%ncdimnames(2))
  write(var_infos%ncdimnames(2), "(A)") "number_of_atom_species"
  write(var_infos%ncdimnames(1), "(A)") "character_string_length"
  call test_var(ncid, var_infos, valid, error_data)
  deallocate(var_infos%ncdimnames)
  if (.not. valid .and. error_data%access_mode_id == ERROR_MODE_SPEC) return
  lstat = lstat .or. valid
  ! Variable chemical_symbols
  write(var_infos%name, "(A)") "chemical_symbols"
  var_infos%nctype  = etsf_io_low_character
  var_infos%ncshape = 2
  allocate(var_infos%ncdimnames(2))
  write(var_infos%ncdimnames(2), "(A)") "number_of_atom_species"
  write(var_infos%ncdimnames(1), "(A)") "symbol_length"
  call test_var(ncid, var_infos, valid, error_data)
  deallocate(var_infos%ncdimnames)
  if (.not. valid .and. error_data%access_mode_id == ERROR_MODE_SPEC) return
  lstat = lstat .or. valid
  if (.not. lstat) then
    call etsf_io_split_free(split)
    call etsf_io_low_error_set(error_data, ERROR_MODE_SPEC, &
                             & ERROR_TYPE_ARG, me, &
   & tgtname = "atomic_numbers, atom_species_names, chemical_symbols", &
                             & errmess = "missing one among the list.")
    return
  end if
  


  ! Deallocate the split data.
  call etsf_io_split_free(split)

  lstat = .true.
end subroutine etsf_io_file_check_crystallographic_data