TABLE OF CONTENTS


etsf_io_gwdata_get

[ Top ] [ etsf_gwdata ] [ Methods ]

NAME

etsf_io_gwdata_get

FUNCTION

Read an opened ETSF file to get data related to the given group. Only associated pointers of argument @folder will be accessed. If any accessed variable is missing, this routine returns an error (usually an access_mode_id of argument error_data set to ERROR_MODE_INQ). Any other errors implies a return with @lstat = .false..

INPUTS

OUTPUT

SIDE EFFECTS

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_gwdata_get(ncid, folder, lstat, error_data, use_atomic_units)

  !Arguments ------------------------------------
  integer, intent(in) :: ncid
  type(etsf_gwdata), intent(inout) :: folder
  logical, intent(out) :: lstat
  type(etsf_io_low_error), intent(out) :: error_data
  logical, optional, intent(in) :: use_atomic_units

  !Local variables-------------------------------
  character(len = *), parameter :: my_name = 'etsf_io_gwdata_get'
  logical :: my_use_atomic_units
  integer,allocatable :: varid(:)
  integer,allocatable :: start(:)
  integer,allocatable :: count(:)
  integer :: len
  character(etsf_charlen) :: flag


  ! *************************************************************************

!DEBUG
!write (*,*) 'etsf_io_gwdata_get : enter'
!ENDDEBUG

  ! Get values for optional arguments, set default.
  if (present(use_atomic_units)) then
    my_use_atomic_units = use_atomic_units
  else
    my_use_atomic_units = .true.
  end if
  
  allocate(varid(4))
  
  if (etsf_io_low_var_associated(folder%gw_corrections)) then
    allocate(start(4), count(4))
    start(:) = 1
    count(:) = 0
    if (folder%gw_corrections__spin_access /= etsf_no_sub_access) then
      start(4) = folder%gw_corrections__spin_access
      count(4) = 1
    end if
    if (folder%gw_corrections__kpoint_access /= etsf_no_sub_access) then
      start(3) = folder%gw_corrections__kpoint_access
      count(3) = 1
    end if
    count(2) = folder%gw_corrections__number_of_states
    if (folder%gw_corrections__state_access /= etsf_no_sub_access) then
      start(2) = folder%gw_corrections__state_access
      count(2) = 1
    end if
    call etsf_io_low_read_var(ncid, "gw_corrections", &
                            & folder%gw_corrections, &
                            & lstat, ncvarid = varid(1), &
                            & error_data = error_data, start = start, count = count)
    deallocate(start, count)
    if (.not. lstat) then
      call etsf_io_low_error_update(error_data, my_name)
      return
    end if
  end if
  
  if (etsf_io_low_var_associated(folder%kb_formfactor_sign)) then
    allocate(start(3), count(3))
    start(:) = 1
    count(:) = 0
    count(2) = folder%kb_coeff_sig__number_of_angular_momenta
    count(1) = folder%kb_coeff_sig__number_of_projectors
    call etsf_io_low_read_var(ncid, "kb_formfactor_sign", &
                            & folder%kb_formfactor_sign, &
                            & lstat, ncvarid = varid(2), &
                            & error_data = error_data, start = start, count = count)
    deallocate(start, count)
    if (.not. lstat) then
      call etsf_io_low_error_update(error_data, my_name)
      return
    end if
  end if
  
  if (etsf_io_low_var_associated(folder%kb_formfactors)) then
    allocate(start(5), count(5))
    start(:) = 1
    count(:) = 0
    count(4) = folder%kb_coeff__number_of_angular_momenta
    count(3) = folder%kb_coeff__number_of_projectors
    if (folder%kb_coeff__kpoint_access /= etsf_no_sub_access) then
      start(2) = folder%kb_coeff__kpoint_access
      count(2) = 1
    end if
    count(1) = folder%kb_coeff__number_of_coefficients
    call etsf_io_low_read_var(ncid, "kb_formfactors", &
                            & folder%kb_formfactors, &
                            & lstat, ncvarid = varid(3), &
                            & error_data = error_data, start = start, count = count)
    deallocate(start, count)
    if (.not. lstat) then
      call etsf_io_low_error_update(error_data, my_name)
      return
    end if
  end if
  
  if (etsf_io_low_var_associated(folder%kb_formfactor_derivative)) then
    allocate(start(5), count(5))
    start(:) = 1
    count(:) = 0
    count(4) = folder%kb_coeff_der__number_of_angular_momenta
    count(3) = folder%kb_coeff_der__number_of_projectors
    if (folder%kb_coeff_der__kpoint_access /= etsf_no_sub_access) then
      start(2) = folder%kb_coeff_der__kpoint_access
      count(2) = 1
    end if
    count(1) = folder%kb_coeff_der__number_of_coefficients
    call etsf_io_low_read_var(ncid, "kb_formfactor_derivative", &
                            & folder%kb_formfactor_derivative, &
                            & lstat, ncvarid = varid(4), &
                            & error_data = error_data, start = start, count = count)
    deallocate(start, count)
    if (.not. lstat) then
      call etsf_io_low_error_update(error_data, my_name)
      return
    end if
  end if
  
  deallocate(varid)

!DEBUG
!write (*,*) 'etsf_io_gwdata_get : exit'
!ENDDEBUG

end subroutine etsf_io_gwdata_get