subroutine computeAreaGrid(grid, petNo, area, regridScheme, rc)
type(ESMF_Grid) :: grid
integer :: petNo
real (ESMF_KIND_R8), pointer :: area(:)
integer :: regridScheme
integer :: rc
type(ESMF_Field) :: areaField
integer :: minIndex(2), maxIndex(2), gridDims(2)
real (ESMF_KIND_R8), pointer :: area2D(:,:)
integer :: i, start, ntiles
integer :: localrc
! Create a field on the grid to hold the areas
areaField=ESMF_FieldCreate(grid, typekind=ESMF_TYPEKIND_R8, &
staggerloc=ESMF_STAGGERLOC_CENTER, name="area", rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (localrc /=ESMF_SUCCESS) then
rc=localrc
return
endif
! compute areas
call ESMF_FieldRegridGetArea(areaField, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (localrc /=ESMF_SUCCESS) then
rc=localrc
return
endif
! Get number of tiles
call ESMF_GridGet(grid, tileCount = ntiles, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! Get size of Grid
call ESMF_GridGet(grid, tile=1, staggerloc=ESMF_STAGGERLOC_CENTER, &
minIndex=minIndex, maxIndex=maxIndex, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (localrc /=ESMF_SUCCESS) then
rc=localrc
return
endif
! Grid size
gridDims(1)=maxIndex(1)-minIndex(1)+1
gridDims(2)=maxIndex(2)-minIndex(2)+1
! Allocate memory for area
allocate(area2D(gridDims(1),gridDims(2)))
! Only do this part on PET 0
if (petNo .eq. 0) then
! Allocate memory for area
allocate(area(gridDims(1)*gridDims(2)*ntiles))
endif
! Get area onto PET 0
start=1
do i=1,ntiles
call ESMF_FieldGather(areaField, farray=area2D, rootPet=0, tile=i, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! copy to 1D array
if (PetNo == 0) then
! flatten area
area(start:start+gridDims(1)*gridDims(2)-1)=RESHAPE(area2D,(/gridDims(1)*gridDims(2)/))
start= start+gridDims(1)*gridDims(2)
endif
enddo
! deallocate memory for 2D area
deallocate(area2D)
end subroutine computeAreaGrid