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