MAPL_IOChangeRes Subroutine

public subroutine MAPL_IOChangeRes(cfIn, cfOut, dimNames, dimSizes, rc)

Arguments

Type IntentOptional Attributes Name
type(FileMetadata), intent(inout) :: cfIn
type(FileMetadata), intent(inout) :: cfOut
character(len=*) :: dimNames(:)
integer, intent(in) :: dimSizes(:)
integer, intent(out), optional :: rc

Calls

proc~~mapl_iochangeres~~CallsGraph proc~mapl_iochangeres MAPL_IOChangeRes at at proc~mapl_iochangeres->at begin begin proc~mapl_iochangeres->begin first first proc~mapl_iochangeres->first insert insert proc~mapl_iochangeres->insert interface~mapl_range MAPL_Range proc~mapl_iochangeres->interface~mapl_range next next proc~mapl_iochangeres->next none~first~238 StringVariableMapIterator%first proc~mapl_iochangeres->none~first~238 none~ftn_begin~35 StringVariableMap%ftn_begin proc~mapl_iochangeres->none~ftn_begin~35 none~ftn_end~35 StringVariableMap%ftn_end proc~mapl_iochangeres->none~ftn_end~35 none~get_coordinate_data CoordinateVariable%get_coordinate_data proc~mapl_iochangeres->none~get_coordinate_data none~get_coordinate_variable~2 FileMetadata%get_coordinate_variable proc~mapl_iochangeres->none~get_coordinate_variable~2 none~get_dimensions~2 FileMetadata%get_dimensions proc~mapl_iochangeres->none~get_dimensions~2 none~get_variables~2 FileMetadata%get_variables proc~mapl_iochangeres->none~get_variables~2 none~modify_dimension FileMetadata%modify_dimension proc~mapl_iochangeres->none~modify_dimension none~next~96 StringVariableMapIterator%next proc~mapl_iochangeres->none~next~96 none~replace_coordinate_data CoordinateVariable%replace_coordinate_data proc~mapl_iochangeres->none~replace_coordinate_data proc~mapl_return MAPL_Return proc~mapl_iochangeres->proc~mapl_return proc~mapl_verify MAPL_Verify proc~mapl_iochangeres->proc~mapl_verify none~of~182 StringVariableMapIterator%of none~first~238->none~of~182 to_node to_node none~ftn_end~35->to_node none~get_coordinate_data->proc~mapl_return none~get_coordinate_variable~2->proc~mapl_return none~get_coordinate_variable~2->proc~mapl_verify interface~mapl_assert MAPL_Assert none~get_coordinate_variable~2->interface~mapl_assert none~at~212 StringVariableMap%at none~get_coordinate_variable~2->none~at~212 none~get_variables~2->proc~mapl_return none~modify_dimension->proc~mapl_return set set none~modify_dimension->set none~next~96->to_node none~replace_coordinate_data->proc~mapl_return proc~mapl_return->at proc~mapl_return->insert proc~mapl_throw_exception MAPL_throw_exception proc~mapl_return->proc~mapl_throw_exception proc~mapl_verify->proc~mapl_throw_exception none~at_rc~10 StringVariableMap%at_rc none~at~212->none~at_rc~10 none~find~51 StringVariableMap%find none~at_rc~10->none~find~51

Source Code

  subroutine MAPL_IOChangeRes(cfIn,cfOut,dimNames,dimSizes,rc)
  type(FileMetadata), intent(inout) :: cfIn
  type(Filemetadata), intent(inout) :: cfOut
  character(len=*) :: dimNames(:)
  integer, intent(in) :: dimSizes(:)
  integer, intent(out), optional :: rc

  integer :: status
  type(StringIntegerMap) :: newDims
  integer :: i

  do i=1,size(dimNames)
     call newDims%insert(trim(dimNames(i)),dimSizes(i))
  enddo

  cfOut = cfIn
  call modify_grid_dimensions(rc=status)
  _VERIFY(status)
  call modify_coordinate_vars(rc=status)

  _RETURN(ESMF_SUCCESS)

  contains

      subroutine modify_grid_dimensions(rc)
         integer, optional, intent(out) :: rc
         integer :: status
         type(StringIntegerMap), pointer :: dims
         type(StringIntegerMapIterator) :: iter
         character(len=:), pointer :: name
         integer, pointer :: newExtent => null()

         dims => cfIn%get_dimensions()

         iter = dims%begin()
         do while (iter /= dims%end())
            name => iter%first()
            newExtent => newDims%at(trim(name))
            if (associated(newExtent)) then
               call cfOut%modify_dimension(trim(name),newExtent,rc=status)
               nullify(newExtent)
            end if
            call iter%next()
         enddo

         _RETURN(ESMF_SUCCESS)

      end subroutine modify_grid_dimensions

      subroutine modify_coordinate_vars(rc)
         integer, optional, intent(out) :: rc

         integer :: status
         type(StringVariableMap), pointer :: vars
         type(StringVariableMapIterator) :: iter
         type(CoordinateVariable), pointer :: cvar
         character(len=:), pointer :: name
         real(kind=REAL32) :: r32_x1,r32_x0
         real(kind=REAL64) :: r64_x1,r64_x0
         real(kind=REAL32), allocatable :: var32(:)
         real(kind=REAL64), allocatable :: var64(:)
         integer, pointer :: newExtent => null()
         class(*), pointer :: dim_var_values(:)
         class(*), allocatable :: coordinate_data(:)

          vars => cfIn%get_variables(_RC)

         iter = vars%ftn_begin()
         do while (iter /= vars%ftn_end())
            call iter%next()

            name => iter%first()
            newExtent => newDims%at(trim(name))
            if (associated(newExtent)) then
               cvar => cfOut%get_coordinate_variable(trim(name),rc=status)
               if (status==ESMF_SUCCESS) then
                  dim_var_values => cvar%get_coordinate_data()
                  select type(q => dim_var_values)
                  type is (real(REAL32))
                     r32_x0=1.0d0
                     r32_x1=dble(newExtent)
                     var32 = MAPL_Range(r32_x0,r32_x1,newExtent)
                     allocate(coordinate_data,source=var32)
                     call cvar%replace_coordinate_data(coordinate_data)
                     deallocate(coordinate_data,var32)
                  type is (real(REAL64))
                     r64_x0=1.0d0
                     r64_x1=dble(newExtent)
                     var64 = MAPL_Range(r64_x0,r64_x1,newExtent)
                     allocate(coordinate_data,source=var64)
                     call cvar%replace_coordinate_data(coordinate_data)
                     deallocate(coordinate_data,var64)
                  class default
                     status = ESMF_FAILURE
                  end select

               end if

               nullify(newExtent)
            end if
         enddo

         _RETURN(ESMF_SUCCESS)

      end subroutine modify_coordinate_vars

  end subroutine MAPL_IOChangeRes