TABLE OF CONTENTS


etsf_io_main_get

[ Top ] [ etsf_main ] [ Methods ]

NAME

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

  !Arguments ------------------------------------
  integer, intent(in) :: ncid
  type(etsf_main), 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_main_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_main_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(6))
  
  if (etsf_io_low_var_associated(folder%density)) then
    call etsf_io_low_read_var(ncid, "density", &
                            & folder%density, &
                            & lstat, ncvarid = varid(1), &
                            & error_data = error_data)
    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%exchange_potential)) then
    call etsf_io_low_read_var(ncid, "exchange_potential", &
                            & folder%exchange_potential, &
                            & lstat, ncvarid = varid(2), &
                            & error_data = error_data)
    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%correlation_potential)) then
    call etsf_io_low_read_var(ncid, "correlation_potential", &
                            & folder%correlation_potential, &
                            & lstat, ncvarid = varid(3), &
                            & error_data = error_data)
    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%exchange_correlation_potential)) then
    call etsf_io_low_read_var(ncid, "exchange_correlation_potential", &
                            & folder%exchange_correlation_potential, &
                            & lstat, ncvarid = varid(4), &
                            & error_data = error_data)
    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%coefficients_of_wavefunctions)) then
    allocate(start(6), count(6))
    start(:) = 1
    count(:) = 0
    if (folder%wfs_coeff__spin_access /= etsf_no_sub_access) then
      start(6) = folder%wfs_coeff__spin_access
      count(6) = 1
    end if
    if (folder%wfs_coeff__kpoint_access /= etsf_no_sub_access) then
      start(5) = folder%wfs_coeff__kpoint_access
      count(5) = 1
    end if
    count(4) = folder%wfs_coeff__number_of_states
    if (folder%wfs_coeff__state_access /= etsf_no_sub_access) then
      start(4) = folder%wfs_coeff__state_access
      count(4) = 1
    end if
    count(2) = folder%wfs_coeff__number_of_coefficients
    call etsf_io_low_read_var(ncid, "coefficients_of_wavefunctions", &
                            & folder%coefficients_of_wavefunctions, &
                            & lstat, ncvarid = varid(5), &
                            & 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%real_space_wavefunctions)) then
    allocate(start(8), count(8))
    start(:) = 1
    count(:) = 0
    if (folder%wfs_rsp__spin_access /= etsf_no_sub_access) then
      start(8) = folder%wfs_rsp__spin_access
      count(8) = 1
    end if
    if (folder%wfs_rsp__kpoint_access /= etsf_no_sub_access) then
      start(7) = folder%wfs_rsp__kpoint_access
      count(7) = 1
    end if
    count(6) = folder%wfs_rsp__number_of_states
    if (folder%wfs_rsp__state_access /= etsf_no_sub_access) then
      start(6) = folder%wfs_rsp__state_access
      count(6) = 1
    end if
    call etsf_io_low_read_var(ncid, "real_space_wavefunctions", &
                            & folder%real_space_wavefunctions, &
                            & lstat, ncvarid = varid(6), &
                            & 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
  
  ! Handle all attributes for the group.
  if (etsf_io_low_var_associated(folder%density)) then
    ! Handle the units attribute.
    call etsf_io_low_read_att(ncid, varid(1), &
                            & "units", &
                            & etsf_charlen, folder%density__units, &
                            & lstat, error_data = error_data)
    if (.not. lstat) return
    if (trim(folder%density__units) /= "atomic units") then
      call etsf_io_low_read_att(ncid, varid(1), &
                              & "scale_to_atomic_units", &
                              & folder%density__scale_to_atomic_units, &
                              & lstat, error_data = error_data)
      if (.not. lstat) return
    else
      folder%density__scale_to_atomic_units = 1.0d0
    end if
    if (my_use_atomic_units .and. &
      & folder%density__scale_to_atomic_units /= 1.0d0) then
      call etsf_io_low_var_multiply(folder%density, &
                                  & folder%density__scale_to_atomic_units)
    end if
  end if
  
  if (etsf_io_low_var_associated(folder%exchange_potential)) then
    ! Handle the units attribute.
    call etsf_io_low_read_att(ncid, varid(2), &
                            & "units", &
                            & etsf_charlen, folder%pot_x_only__units, &
                            & lstat, error_data = error_data)
    if (.not. lstat) return
    if (trim(folder%pot_x_only__units) /= "atomic units") then
      call etsf_io_low_read_att(ncid, varid(2), &
                              & "scale_to_atomic_units", &
                              & folder%pot_x_only__scale_to_atomic_units, &
                              & lstat, error_data = error_data)
      if (.not. lstat) return
    else
      folder%pot_x_only__scale_to_atomic_units = 1.0d0
    end if
    if (my_use_atomic_units .and. &
      & folder%pot_x_only__scale_to_atomic_units /= 1.0d0) then
      call etsf_io_low_var_multiply(folder%exchange_potential, &
                                  & folder%pot_x_only__scale_to_atomic_units)
    end if
  end if
  
  if (etsf_io_low_var_associated(folder%correlation_potential)) then
    ! Handle the units attribute.
    call etsf_io_low_read_att(ncid, varid(3), &
                            & "units", &
                            & etsf_charlen, folder%pot_c_only__units, &
                            & lstat, error_data = error_data)
    if (.not. lstat) return
    if (trim(folder%pot_c_only__units) /= "atomic units") then
      call etsf_io_low_read_att(ncid, varid(3), &
                              & "scale_to_atomic_units", &
                              & folder%pot_c_only__scale_to_atomic_units, &
                              & lstat, error_data = error_data)
      if (.not. lstat) return
    else
      folder%pot_c_only__scale_to_atomic_units = 1.0d0
    end if
    if (my_use_atomic_units .and. &
      & folder%pot_c_only__scale_to_atomic_units /= 1.0d0) then
      call etsf_io_low_var_multiply(folder%correlation_potential, &
                                  & folder%pot_c_only__scale_to_atomic_units)
    end if
  end if
  
  if (etsf_io_low_var_associated(folder%exchange_correlation_potential)) then
    ! Handle the units attribute.
    call etsf_io_low_read_att(ncid, varid(4), &
                            & "units", &
                            & etsf_charlen, folder%pot_xc__units, &
                            & lstat, error_data = error_data)
    if (.not. lstat) return
    if (trim(folder%pot_xc__units) /= "atomic units") then
      call etsf_io_low_read_att(ncid, varid(4), &
                              & "scale_to_atomic_units", &
                              & folder%pot_xc__scale_to_atomic_units, &
                              & lstat, error_data = error_data)
      if (.not. lstat) return
    else
      folder%pot_xc__scale_to_atomic_units = 1.0d0
    end if
    if (my_use_atomic_units .and. &
      & folder%pot_xc__scale_to_atomic_units /= 1.0d0) then
      call etsf_io_low_var_multiply(folder%exchange_correlation_potential, &
                                  & folder%pot_xc__scale_to_atomic_units)
    end if
  end if
  
  
  deallocate(varid)

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

end subroutine etsf_io_main_get