subroutine IORead2D(vm, field, &
minIndexPDe, maxIndexPDe, minIndexPTile, maxIndexPTile, keywordEnforcer, &
fileName, iofmt, localDe, ncid, rc)
type(ESMF_VM), intent(in) :: vm
type(ESMF_Field), intent(in) :: field
integer, dimension(:), intent(in) :: minIndexPDe
integer, dimension(:), intent(in) :: maxIndexPDe
integer, dimension(:), intent(in) :: minIndexPTile
integer, dimension(:), intent(in) :: maxIndexPTile
type(ESMF_KeywordEnforcer), optional:: keywordEnforcer ! must use keywords below
character(len=*), intent(in), optional :: fileName
type(ESMF_IOFmt_flag), intent(in), optional :: iofmt
integer, intent(in), optional :: localDe
integer, intent(in), optional :: ncid
integer, intent(out), optional :: rc
! -- local variables
integer :: localrc
integer :: ilen, jlen, lbuf, lde, localpe, rank
integer :: lncid, varId, ncStatus, ndims, xtype
integer, dimension(2) :: elb, eub
integer(ESMF_KIND_I4), dimension(:), allocatable :: bcstbuf_i4
integer(ESMF_KIND_I4), dimension(:,:), allocatable :: buf_i4
integer(ESMF_KIND_I4), dimension(:,:), pointer :: fp_i4 => null()
real(ESMF_KIND_R4), dimension(:), allocatable :: bcstbuf_r4
real(ESMF_KIND_R4), dimension(:,:), allocatable :: buf_r4
real(ESMF_KIND_R4), dimension(:,:), pointer :: fp_r4 => null()
real(ESMF_KIND_R8), dimension(:), allocatable :: bcstbuf_r8
real(ESMF_KIND_R8), dimension(:,:), allocatable :: buf_r8
real(ESMF_KIND_R8), dimension(:,:), pointer :: fp_r8 => null()
character(len=ESMF_MAXSTR) :: fieldName, dataSetName
type(ESMF_TypeKind_Flag) :: typekind
! -- begin
if (present(rc)) rc = ESMF_SUCCESS
lde = 0
if (present(localDe)) lde = localDe
call ESMF_FieldGet(field, name=fieldName, rank=rank, &
typekind=typekind, rc=localrc)
if (ESMF_LogFoundError(rcToCheck=localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return ! bail out
ilen = maxIndexPTile(1)-minIndexPTile(1)+1
jlen = maxIndexPTile(2)-minIndexPTile(2)+1
lbuf = ilen * jlen
if (typekind == ESMF_TYPEKIND_I4) then
allocate(buf_i4(minIndexPTile(1):maxIndexPTile(1), &
minIndexPTile(2):maxIndexPTile(2)), stat=localrc)
if (ESMF_LogFoundAllocError(statusToCheck=localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return ! bail out
buf_i4 = 0_ESMF_KIND_I4
else if (typekind == ESMF_TYPEKIND_R4) then
allocate(buf_r4(minIndexPTile(1):maxIndexPTile(1), &
minIndexPTile(2):maxIndexPTile(2)), stat=localrc)
if (ESMF_LogFoundAllocError(statusToCheck=localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return ! bail out
buf_r4 = 0._ESMF_KIND_R4
else if (typekind == ESMF_TYPEKIND_R8) then
allocate(buf_r8(minIndexPTile(1):maxIndexPTile(1), &
minIndexPTile(2):maxIndexPTile(2)), stat=localrc)
if (ESMF_LogFoundAllocError(statusToCheck=localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return ! bail out
buf_r8 = 0._ESMF_KIND_R8
end if
call ESMF_VMGet(vm, localPet=localpe, rc=localrc)
if (ESMF_LogFoundError(rcToCheck=localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return ! bail out
if (localpe == 0) then
if (iofmt == ESMF_IOFMT_NETCDF) then
#ifdef ESMF_NETCDF
lncid = 0
dataSetName = "NetCDF data set"
if (present(ncid)) then
lncid = ncid
else if (present(fileName)) then
dataSetName = trim(dataSetName) // " " // trim(fileName)
ncStatus = nf90_open(trim(fileName), NF90_NOWRITE, lncid)
if (ESMF_LogFoundNetCDFError(ncerrToCheck=ncStatus, &
msg="Field "//trim(fieldName)//" not defined in "//trim(dataSetName), &
ESMF_CONTEXT, rcToReturn=rc)) return ! bail out
end if
! -- add data
if (lncid /= 0) then
ncStatus = nf90_inq_varid(lncid, trim(fieldName), varId)
if (ESMF_LogFoundNetCDFError(ncerrToCheck=ncStatus, &
msg="Field "//trim(fieldName)//" not defined in "//trim(dataSetName), &
ESMF_CONTEXT, rcToReturn=rc)) return ! bail out
ncStatus = nf90_inquire_variable(lncid, varId, xtype=xtype, ndims=ndims)
if (ESMF_LogFoundNetCDFError(ncerrToCheck=ncStatus, &
msg="Error inquiring variable "//trim(fieldName)//" in "//trim(dataSetName), &
ESMF_CONTEXT, rcToReturn=rc)) return ! bail out
if (ndims /= rank) then
call ESMF_LogSetError(ESMF_RC_NOT_IMPL, &
msg="Variable "//trim(fieldName)//" has different rank than Field", &
ESMF_CONTEXT, rcToReturn=rc)
return ! bail out
end if
if (typekind == ESMF_TYPEKIND_I4) then
if (xtype /= NF90_INT) then
call ESMF_LogSetError(ESMF_RC_NOT_IMPL, &
msg="Variable "//trim(fieldName)//" has different typekind than Field", &
ESMF_CONTEXT, rcToReturn=rc)
return ! bail out
end if
ncStatus = nf90_get_var(lncid, varId, buf_i4)
if (ESMF_LogFoundNetCDFError(ncerrToCheck=ncStatus, &
msg="Error reading field "//trim(fieldName)//" from "//trim(dataSetName), &
ESMF_CONTEXT, rcToReturn=rc)) return ! bail out
else if (typekind == ESMF_TYPEKIND_R4) then
if (xtype /= NF90_FLOAT) then
call ESMF_LogSetError(ESMF_RC_NOT_IMPL, &
msg="Variable "//trim(fieldName)//" has different typekind than Field", &
ESMF_CONTEXT, rcToReturn=rc)
return ! bail out
end if
ncStatus = nf90_get_var(lncid, varId, buf_r4)
if (ESMF_LogFoundNetCDFError(ncerrToCheck=ncStatus, &
msg="Error reading field "//trim(fieldName)//" from "//trim(dataSetName), &
ESMF_CONTEXT, rcToReturn=rc)) return ! bail out
else if (typekind == ESMF_TYPEKIND_R8) then
if (xtype /= NF90_DOUBLE) then
call ESMF_LogSetError(ESMF_RC_NOT_IMPL, &
msg="Variable "//trim(fieldName)//" has different typekind than Field", &
ESMF_CONTEXT, rcToReturn=rc)
return ! bail out
end if
ncStatus = nf90_get_var(lncid, varId, buf_r8)
if (ESMF_LogFoundNetCDFError(ncerrToCheck=ncStatus, &
msg="Error reading field "//trim(fieldName)//" from "//trim(dataSetName), &
ESMF_CONTEXT, rcToReturn=rc)) return ! bail out
else
call ESMF_LogSetError(ESMF_RC_NOT_IMPL, &
msg="Field: "//trim(fieldName)//" - typekind not supported", &
ESMF_CONTEXT, rcToReturn=rc)
return ! bail out
end if
end if
#else
call ESMF_LogSetError(rcToCheck=ESMF_RC_LIB_NOT_PRESENT, &
msg="- ESMF_NETCDF not defined when lib was compiled", &
ESMF_CONTEXT, rcToReturn=rc)
#endif
end if
end if
if (typekind == ESMF_TYPEKIND_I4) then
allocate(bcstbuf_i4(lbuf), stat=localrc)
if (ESMF_LogFoundAllocError(statusToCheck=localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return ! bail out
bcstbuf_i4 = reshape(buf_i4, (/lbuf/))
call ESMF_VMBroadcast(vm, bcstbuf_i4, lbuf, 0, rc=localrc)
if (ESMF_LogFoundError(rcToCheck=localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return ! bail out
call ESMF_FieldGet(field, localDe=lde, farrayPtr=fp_i4, &
exclusiveLBound=elb, exclusiveUBound=eub, rc=localrc)
if (ESMF_LogFoundError(rcToCheck=localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return ! bail out
buf_i4 = reshape(bcstbuf_i4, (/ilen,jlen/))
fp_i4(elb(1):eub(1),elb(2):eub(2)) = buf_i4(minIndexPDe(1):maxIndexPDe(1), &
minIndexPDe(2):maxIndexPDe(2))
deallocate(buf_i4, bcstbuf_i4, stat=localrc)
if (ESMF_LogFoundDeallocError(statusToCheck=localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return ! bail out
else if (typekind == ESMF_TYPEKIND_R4) then
allocate(bcstbuf_r4(lbuf), stat=localrc)
if (ESMF_LogFoundAllocError(statusToCheck=localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return ! bail out
bcstbuf_r4 = reshape(buf_r4, (/lbuf/))
call ESMF_VMBroadcast(vm, bcstbuf_r4, lbuf, 0, rc=localrc)
if (ESMF_LogFoundError(rcToCheck=localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return ! bail out
call ESMF_FieldGet(field, localDe=lde, farrayPtr=fp_r4, &
exclusiveLBound=elb, exclusiveUBound=eub, rc=localrc)
if (ESMF_LogFoundError(rcToCheck=localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return ! bail out
buf_r4 = reshape(bcstbuf_r4, (/ilen,jlen/))
fp_r4(elb(1):eub(1),elb(2):eub(2)) = buf_r4(minIndexPDe(1):maxIndexPDe(1), &
minIndexPDe(2):maxIndexPDe(2))
deallocate(buf_r4, bcstbuf_r4, stat=localrc)
if (ESMF_LogFoundDeallocError(statusToCheck=localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return ! bail out
else if (typekind == ESMF_TYPEKIND_R8) then
allocate(bcstbuf_r8(lbuf), stat=localrc)
if (ESMF_LogFoundAllocError(statusToCheck=localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return ! bail out
bcstbuf_r8 = reshape(buf_r8, (/lbuf/))
call ESMF_VMBroadcast(vm, bcstbuf_r8, lbuf, 0, rc=localrc)
if (ESMF_LogFoundError(rcToCheck=localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return ! bail out
call ESMF_FieldGet(field, localDe=lde, farrayPtr=fp_r8, &
exclusiveLBound=elb, exclusiveUBound=eub, rc=localrc)
if (ESMF_LogFoundError(rcToCheck=localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return ! bail out
buf_r8 = reshape(bcstbuf_r8, (/ilen,jlen/))
fp_r8(elb(1):eub(1),elb(2):eub(2)) = buf_r8(minIndexPDe(1):maxIndexPDe(1), &
minIndexPDe(2):maxIndexPDe(2))
deallocate(buf_r8, bcstbuf_r8, stat=localrc)
if (ESMF_LogFoundDeallocError(statusToCheck=localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return ! bail out
end if
end subroutine IORead2D