MAPL_CFIOReadBundleRead Subroutine

public subroutine MAPL_CFIOReadBundleRead(MCFIO, tindex, hw, RC)

Arguments

Type IntentOptional Attributes Name
type(MAPL_CFIO), intent(inout) :: MCFIO
integer, intent(inout) :: tindex
integer, intent(in), optional :: hw
integer, intent(out), optional :: RC

Calls

proc~~mapl_cfioreadbundleread~~CallsGraph proc~mapl_cfioreadbundleread MAPL_CFIOReadBundleRead esmf_fieldbundleget esmf_fieldbundleget proc~mapl_cfioreadbundleread->esmf_fieldbundleget esmf_fieldget esmf_fieldget proc~mapl_cfioreadbundleread->esmf_fieldget interface~mapl_am_i_root MAPL_Am_I_Root proc~mapl_cfioreadbundleread->interface~mapl_am_i_root interface~mapl_arrayiscatter MAPL_ArrayIScatter proc~mapl_cfioreadbundleread->interface~mapl_arrayiscatter interface~mapl_assert MAPL_Assert proc~mapl_cfioreadbundleread->interface~mapl_assert nf90_get_var nf90_get_var proc~mapl_cfioreadbundleread->nf90_get_var none~at~128 CFIOCollectionVector%at proc~mapl_cfioreadbundleread->none~at~128 none~find~33 CFIOCollection%find proc~mapl_cfioreadbundleread->none~find~33 none~regrid~12 AbstractRegridder%regrid proc~mapl_cfioreadbundleread->none~regrid~12 none~set_undef_value AbstractRegridder%set_undef_value proc~mapl_cfioreadbundleread->none~set_undef_value proc~mapl_createrequest MAPL_CreateRequest proc~mapl_cfioreadbundleread->proc~mapl_createrequest proc~mapl_gridget MAPL_GridGet proc~mapl_cfioreadbundleread->proc~mapl_gridget proc~mapl_return MAPL_Return proc~mapl_cfioreadbundleread->proc~mapl_return proc~mapl_roundrobinpelist MAPL_RoundRobinPEList proc~mapl_cfioreadbundleread->proc~mapl_roundrobinpelist proc~mapl_verify MAPL_Verify proc~mapl_cfioreadbundleread->proc~mapl_verify

Source Code

  subroutine MAPL_CFIOReadBundleRead(MCFIO,tindex,hw,rc)

    type(MAPL_CFIO  ),               intent(INOUT) :: MCFIO
    integer,                         intent(INOUT) :: tindex
    integer,               optional, intent(IN   ) :: hw
    integer,               optional, intent(  OUT) :: RC

    integer          :: status

    integer   :: nn,l,k,klev,lm,nv,lb
    integer   :: img,jmg,lt,hw_
    integer   :: counts(3)
    logical   :: myGlobal, transAlreadyDone
    real, pointer :: ptr2(:,:)
    real, pointer :: ptr3(:,:,:)
    type(ESMF_Field) :: field
    type(ESMF_CFIO), pointer :: cfiop
    type(CFIOCollection), pointer :: collection
    type(Ptr2Arr), allocatable :: globPtrArr(:)
    type(Ptr2Arr)        :: PtrTypeIn(2)
    type(Ptr2Arr)        :: PtrTypeOut(2)
    integer, allocatable :: varids(:)
    logical, allocatable :: transDone(:)
    integer :: status1,status2
    integer :: alloc_ra

    if (present(hw)) then
       hw_=hw
    else
       hw_=0
    end if

    collection => collections%at(mcfio%collection_id)
    cfiop => collection%find(trim(mcfio%fname), _RC)

    call MAPL_GridGet( MCFIO%GRID, globalCellCountPerDim=COUNTS, RC=STATUS)
    _VERIFY(STATUS)

    call MAPL_RoundRobinPEList(MCFIO%Krank,MCFIO%PartSize,root=MCFIO%ROOT,rc=status)
    _VERIFY(STATUS)

    if (mapL_am_i_root()) write(*,*)"CFIO Read Parallel Processing: ",trim(mcfio%fname)

    call MAPL_RoundRobinPEList(MCFIO%Krank,MCFIO%PartSize,root=MCFIO%ROOT,useFirstRank=.false.,rc=status)
    _VERIFY(STATUS)

    call MAPL_CFIOSetVectorPairs(MCFIO,rc=status)
    _VERIFY(Status)

    img = COUNTS(1)
    jmg = COUNTS(2)

    lt = size(mcfio%reqs)
    allocate(mcfio%buffer(lt),stat=status)
    _VERIFY(status)
    allocate(globPtrArr(lt),stat=status)
    _VERIFY(status)
    allocate(varids(lt),stat=status)
    _VERIFY(status)
    allocate(transDone(lt),source=.false.,stat=status)
    _VERIFY(status)

    nn = 0
    VARS1: do L=1,size(MCFIO%VarDims)

       call ESMF_FieldBundleGet(MCFIO%BUNDLE, trim(MCFIO%VARNAME(L)), field=FIELD,                  RC=STATUS)
       _VERIFY(STATUS)

       CREATE_REQ: if (MCFIO%VarDims(L)==2) then

          nn=nn+1
          varids(nn)=mcfio%varid(L)
          call ESMF_FieldGet(Field,localDE=0, farrayPtr=PTR2, RC=STATUS)
          _VERIFY(STATUS)
          call MAPL_CreateRequest(MCFIO%GRID, MCFIO%Krank(nn), MCFIO%reqs(nn), &
               tag=nn, RequestType=MAPL_IsScatter, DstArray = ptr2 , prepost=.true., hw=hw_, RC=STATUS)

       else if (MCFIO%VarDims(L)==3) then

          call ESMF_FieldGet(Field,localDE=0, farrayPtr=PTR3, RC=STATUS)
          _VERIFY(STATUS)

          do k = 1,MCFIO%lm
             nn=nn+1
             varids(nn)=mcfio%varid(L)
             lb=lbound(ptr3,3)
             call MAPL_CreateRequest(MCFIO%GRID, MCFIO%Krank(nn), MCFIO%reqs(nn), &
                  tag=nn, RequestType=MAPL_IsScatter, DstArray = ptr3(:,:,lb+k-1), prepost=.true.,hw=hw_,RC=STATUS)
          end do

       end if CREATE_REQ

    end do VARS1

    nn = 0
    VARS2: do L=1,size(MCFIO%VarDims)

       if (MCFIO%VarDims(L)==2) then
          LM = 1
       else  if (MCFIO%VarDims(L)==3) then
          LM = MCFIO%lm
       else
          LM = 0
       end if

       do k = 1,lm
          nn=nn+1
          myGlobal = (mcfio%mype == mcfio%krank(nn))
          if (mcfio%kreverse) then
             klev= mcfio%lm-k+1
          else
             klev=k
          end if

          if (myGlobal) then
             nv = mCFIO%pairList(nn)
             alloc_ra = 0
             VECTORTEST: if (nv ==0) then
                alloc_ra = 1
                allocate(mcfio%reqs(nn)%read_array(mcfio%im,mcfio%jm),stat=status)
                _VERIFY(status)

                globPtrArr(nn)%ptr => mcfio%reqs(nn)%read_array
                ptrTypeIn(1)%ptr => globPtrArr(nn)%ptr
                allocate(mcfio%buffer(nn)%ptr(img,jmg),stat=status)
                _VERIFY(STATUS)
                ptrTypeOut(1)%ptr => mcfio%buffer(nn)%ptr
                if (lm == 1) then
                   call readlevel(ptrtypein(1)%ptr,cfiop%fid,varids(nn),cfiop%formatVersion,mcfio%im,mcfio%jm,tindex,rc=status)
                else
                   call readlevel(ptrtypein(1)%ptr,cfiop%fid,varids(nn),cfiop%formatVersion,mcfio%im,mcfio%jm,tindex,lev=klev,rc=status)
                end if
                _VERIFY(STATUS)
                if (mcfio%xshift) call shift180Lon2D_(ptrtypein(1)%ptr)
                call TransAndSave(mcfio,ptrTypein(1:1),ptrtypeout(1:1),mcfio%reqs(nn),.true.,1,hw,rc=status)
                _VERIFY(status)
             else if (nv > 0) then
                ! I am U part of vector
                if (transDone(nn)) then
                   TransAlreadyDone = .true.
                   transDone(nn) = .true.
                   transDone(nv) = .true.
                else
                   TransAlreadyDone = .false.
                   alloc_ra = 2
                   allocate(MCFIO%reqs(nn)%read_array(mcfio%im,mcfio%jm),stat=status)
                   _VERIFY(status)
                   allocate(MCFIO%reqs(nv)%read_array(mcfio%im,mcfio%jm),stat=status)
                   _VERIFY(status)
                   allocate(mcfio%buffer(nn)%ptr(img,jmg),stat=status)
                   _VERIFY(status)
                   allocate(mcfio%buffer(nv)%ptr(img,jmg),stat=status)
                   _VERIFY(status)
                end if
                globPtrArr(nn)%ptr => mcfio%reqs(nn)%read_array
                globPtrArr(nv)%ptr => mcfio%reqs(nv)%read_array
                ptrTypeIn(1)%ptr => globPtrArr(nn)%ptr
                ptrTypeIn(2)%ptr => globPtrArr(nv)%ptr
                ptrTypeOut(1)%ptr => mcfio%buffer(nn)%ptr
                ptrTypeOut(2)%ptr => mcfio%buffer(nv)%ptr
                if (.not.TransAlreadyDone) then
                   if (lm == 1) then
                      call readlevel(ptrtypein(1)%ptr,cfiop%fid,varids(nn),cfiop%formatVersion,mcfio%im,mcfio%jm,tindex,rc=status1)
                      call readlevel(ptrtypein(2)%ptr,cfiop%fid,varids(nv),cfiop%formatVersion,mcfio%im,mcfio%jm,tindex,rc=status2)
                   else
                      call readlevel(ptrtypein(1)%ptr,cfiop%fid,varids(nn),cfiop%formatVersion,mcfio%im,mcfio%jm,tindex,lev=klev,rc=status1)
                      call readlevel(ptrtypein(2)%ptr,cfiop%fid,varids(nv),cfiop%formatVersion,mcfio%im,mcfio%jm,tindex,lev=klev,rc=status2)
                   end if
                   _VERIFY(status1)
                   _VERIFY(status2)
                   if (mcfio%xshift) call shift180Lon2D_(ptrtypein(1)%ptr)
                   if (mcfio%xshift) call shift180Lon2D_(ptrtypein(2)%ptr)
                end if
                call TransAndSave(mcfio,ptrTypein(1:2),ptrtypeout(1:2),mcfio%reqs(nn),.not.TransAlreadyDone,1,hw_,rc=status)
             else
                ! I am V part of vector
                nv = abs(nv)
                if (.not.transDone(nn)) then
                   TransAlreadyDone = .true.
                   transDone(nn)=.true.
                   transDone(nv)=.true.
                else
                   TransAlreadyDone = .false.
                   alloc_ra = 2
                   allocate(MCFIO%reqs(nn)%read_array(mcfio%im,mcfio%jm),stat=status)
                   _VERIFY(status)
                   allocate(MCFIO%reqs(nv)%read_array(mcfio%im,mcfio%jm),stat=status)
                   _VERIFY(status)
                   allocate(mcfio%buffer(nn)%ptr(img,jmg),stat=status)
                   _VERIFY(status)
                   allocate(mcfio%buffer(nv)%ptr(img,jmg),stat=status)
                   _VERIFY(status)
                end if
                globPtrArr(nn)%ptr => mcfio%reqs(nn)%read_array
                globPtrArr(nv)%ptr => mcfio%reqs(nv)%read_array
                ptrTypeIn(1)%ptr => globPtrArr(nv)%ptr
                ptrTypeIn(2)%ptr => globPtrArr(nn)%ptr
                ptrTypeOut(1)%ptr => mcfio%buffer(nv)%ptr
                ptrTypeOut(2)%ptr => mcfio%buffer(nn)%ptr
                if (.not.transAlreadyDone) then
                   if (lm == 1) then
                      call readlevel(ptrtypein(1)%ptr,cfiop%fid,varids(nv),cfiop%formatVersion,mcfio%im,mcfio%jm,tindex,rc=status1)
                      call readlevel(ptrtypein(2)%ptr,cfiop%fid,varids(nn),cfiop%formatVersion,mcfio%im,mcfio%jm,tindex,rc=status2)
                   else
                      call readlevel(ptrtypein(1)%ptr,cfiop%fid,varids(nv),cfiop%formatVersion,mcfio%im,mcfio%jm,tindex,lev=klev,rc=status1)
                      call readlevel(ptrtypein(2)%ptr,cfiop%fid,varids(nn),cfiop%formatVersion,mcfio%im,mcfio%jm,tindex,lev=klev,rc=status2)
                   end if
                   _VERIFY(status1)
                   _VERIFY(status2)
                   if (mcfio%xshift) call shift180Lon2D_(ptrtypein(1)%ptr)
                   if (mcfio%xshift) call shift180Lon2D_(ptrtypein(2)%ptr)
                end if
                call TransAndSave(mcfio,ptrTypein(1:2),ptrtypeout(1:2),mcfio%reqs(nn),.not.TransAlreadyDone,2,hw_,rc=status)
             end if VECTORTEST
             if (alloc_ra > 0) then
                deallocate(mcfio%reqs(nn)%read_array,stat=status)
                _VERIFY(status)
                nullify(mcfio%reqs(nn)%read_array)
                if (alloc_ra > 1) then
                   deallocate(mcfio%reqs(nv)%read_array)
                   nullify(mcfio%reqs(nv)%read_array)
                end if
             end if

          end if
       end do

    end do VARS2

    deallocate(varids,globptrarr,transDone)

    _RETURN(ESMF_SUCCESS)

    contains

  subroutine TransAndSave(mcfio,ptrin,ptrout,req,doTrans,idxOut,hw,rc)
     type(MAPL_CFIO), intent(inout) :: mcfio
     type(Ptr2Arr), intent(inout) :: PtrIn(:)
     type(Ptr2Arr), intent(inout)  :: PtrOut(:)
     type(MAPL_CommRequest), intent(inout)  :: req
     logical, intent(in)  :: doTrans
     integer, intent(in)  :: idxOut
     integer, intent(in)     :: hw
     integer, optional, intent(out) :: rc

     __Iam__('TransAndSave')
     real, pointer  :: gin(:,:)
     real , pointer :: gout(:,:)
     real, dimension(:,:,:), pointer :: uin, uout, vin, vout
     integer :: im, jm
     type(c_ptr) :: cptr

     if (size(ptrin)==1) then

        _ASSERT(idxOut ==1, 'input is scalar, output is not')
        Gin => PtrIn(1)%ptr
        Gout => PtrOut(1)%ptr
        if ( all(shape(gin) == shape(gout)) ) then
           gout = gin
        else
           _ASSERT(associated(mCFIO%regridder), 'mCFIO%regridder is not associated')
           if (mCFIO%regridConservative) then
              call mCFIO%regridder%regrid(gin, gout, rc=status)
              _VERIFY(STATUS)
           else
              call mCFIO%regridder%set_undef_value(MAPL_undef)
              call mCFIO%regridder%regrid(gin,gout,rc=status)
              _VERIFY(STATUS)
           end if
        end if
        if (mcfio%gsiMode) call shift180Lon2D_(gout)
     else

        _ASSERT(size(PtrIn) == 2, 'input is neither a scalar nor a tangent (2d) vector')
        _ASSERT(size(PtrOut) == 2, 'input is a vector, but output is not')
        gout => PtrOut(idxOut)%ptr
        if (doTrans) then

           im = size(PtrIn(1)%ptr,1)
           jm = size(PtrIn(1)%ptr,2)

           cptr = C_loc(PtrIn(1)%ptr(1,1))
           call C_F_pointer (cptr, uin,[im,jm,1])

           cptr = C_loc(PtrIn(2)%ptr(1,1))
           call C_F_pointer (cptr, vin,[im,jm,1])

           im = size(PtrOut(1)%ptr,1)
           jm = size(PtrOut(1)%ptr,2)

           cptr = C_loc(PtrOut(1)%ptr(1,1))
           call C_F_pointer (cptr, uout,[im,jm,1])

           cptr = C_loc(PtrOut(2)%ptr(1,1))
           call C_F_pointer (cptr, vout,[im,jm,1])

           call mCFIO%regridder%set_undef_value(MAPL_undef)
           call mCFIO%regridder%regrid(uin, vin, uout, vout, rotate=mCFIO%doRotate, rc=status)
           _VERIFY(status)

        end if

     end if
     call MAPL_ArrayIScatter(gout,req,hw=hw,rc=status)
     _VERIFY(STATUS)

  end subroutine TransAndSave

    subroutine shift180Lon2D_ ( c  )
    real, intent(inout) :: c(:,:)
    real, allocatable :: cj(:)
    integer :: m(4), n(4), imh, j, im ,jm
    im = size(c,1)
    jm = size(c,2)
    allocate(cj(im))
    imh = nint(im/2.)
    m = [ 1,      imh, 1+imh,    im  ]
    n = [ 1,   im-imh, 1+im-imh, im  ]
    do j = 1, jm
       cj(n(1):n(2)) = c(m(3):m(4),j)
       cj(n(3):n(4)) = c(m(1):m(2),j)
       c(:,j) = cj
    end do
    deallocate(cj)
    return
    end subroutine shift180Lon2D_

    subroutine readlevel(var_array,fid,varid,fformat,im,jm,tindex,lev,rc)
       real, pointer, intent(inout) :: var_array(:,:)
       integer,       intent(in) :: fid
       integer,       intent(in) :: varid
       real,          intent(in) :: fformat
       integer,       intent(in) :: im
       integer,       intent(in) :: jm
       integer,       intent(in) :: tindex
       integer, optional, intent(in) :: lev
       integer, optional, intent(out) :: rc

       integer :: status

       if (present(lev)) then
          if (cfiop%formatVersion > 2.0) then
             status = nf90_get_var(fid,varid,var_array,start=(/1,1,1,lev,tindex/),count=(/im,jm/6,6,1,1/))
          else
             status = nf90_get_var(fid,varid,var_array,start=(/1,1,lev,tindex/),count=(/im,jm,1,1/))
          end if
          _VERIFY(STATUS)
       else
          if (cfiop%formatVersion > 2.0) then
             status = nf90_get_var(fid,varid,var_array,start=(/1,1,1,tindex/),count=(/im,jm/6,6,1/))
          else
             status = nf90_get_var(fid,varid,var_array,start=(/1,1,tindex/),count=(/im,jm,1/))
          end if
          _VERIFY(STATUS)
       end if
       _RETURN(ESMF_SUCCESS)
    end subroutine readlevel

  end subroutine MAPL_CFIOReadBundleRead