subroutine test_regridNearestLocStreamToLocStream(rc)
integer, intent(out) :: rc
logical :: correct
integer :: localrc
type(ESMF_Field) :: srcField
type(ESMF_Field) :: dstField
type(ESMF_RouteHandle) :: routeHandle
type(ESMF_ArraySpec) :: arrayspec
type(ESMF_VM) :: vm
real(ESMF_KIND_R8), pointer :: lon(:),lat(:)
real(ESMF_KIND_R8), pointer :: temperature(:)
real(ESMF_KIND_R8), pointer :: farrayPtr1D(:)
integer :: numLocationsOnThisPet,i
type(ESMF_LocStream) :: srcLocStream,dstLocStream
integer :: localPet, petCount
! init success flag
correct=.true.
rc=ESMF_SUCCESS
! get global VM
call ESMF_VMGetGlobal(vm, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_VMGet(vm, localPet=localPet, petCount=petCount, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! If we don't have 1 or 4 PETS then exit successfully
if ((petCount .ne. 1) .and. (petCount .ne. 4)) then
print*,'ERROR: test must be run using exactly 1 or 4 PETS - detected ',petCount
rc=ESMF_FAILURE
return
endif
! Setup Src LocStream
numLocationsOnThisPet=6
!-------------------------------------------------------------------
! Allocate and set example Field data
!-------------------------------------------------------------------
allocate(temperature(numLocationsOnThisPet))
do i=1,numLocationsOnThisPet
temperature(i)=80.0+i
enddo
!-------------------------------------------------------------------
! Create the LocStream: Allocate space for the LocStream object,
! define the number and distribution of the locations.
!-------------------------------------------------------------------
srcLocStream=ESMF_LocStreamCreate(name="Equatorial Measurements", &
localCount=numLocationsOnThisPet, &
coordSys=ESMF_COORDSYS_SPH_DEG, &
rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble creating locStream'
rc=ESMF_FAILURE
return
endif
!-------------------------------------------------------------------
! Add key data (internally allocating memory).
!-------------------------------------------------------------------
call ESMF_LocStreamAddKey(srcLocStream, &
keyName="ESMF:Lat", &
KeyTypeKind=ESMF_TYPEKIND_R8, &
keyUnits="Degrees", &
keyLongName="Latitude", rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble adding LocStream key for Lat'
rc=ESMF_FAILURE
return
endif
call ESMF_LocStreamAddKey(srcLocStream, &
keyName="ESMF:Lon", &
KeyTypeKind=ESMF_TYPEKIND_R8, &
keyUnits="Degrees", &
keyLongName="Longitude", rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble adding LocStream key for Lon'
rc=ESMF_FAILURE
return
endif
!-------------------------------------------------------------------
! Get key data.
!-------------------------------------------------------------------
call ESMF_LocStreamGetKey(srcLocStream, &
localDE=0, &
keyName="ESMF:Lat", &
farray=lat, &
rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble getting LocStream key for Latitude'
rc=ESMF_FAILURE
return
endif
call ESMF_LocStreamGetKey(srcLocStream, &
localDE=0, &
keyName="ESMF:Lon", &
farray=lon, &
rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble getting LocStream key for Longitude'
rc=ESMF_FAILURE
return
endif
!-------------------------------------------------------------------
! Set key data.
!-------------------------------------------------------------------
do i=1,numLocationsOnThisPet
lon(i)=(i-1)*360.0/numLocationsOnThisPet
lat(i)=0.0
enddo
!-------------------------------------------------------------------
! Create a Field on the Location Stream. In this case the
! Field is created from a user array, but any of the other
! Field create methods (e.g. from ArraySpec) would also apply.
!-------------------------------------------------------------------
srcField=ESMF_FieldCreate(srcLocStream, &
temperature, &
name="temperature", &
rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble creating field on locStream'
rc=ESMF_FAILURE
return
endif
! setup Dst locStream
!-------------------------------------------------------------------
! Create the LocStream: Allocate space for the LocStream object,
! define the number and distribution of the locations.
!-------------------------------------------------------------------
dstLocStream=ESMF_LocStreamCreate(name="Equatorial Measurements", &
localCount=numLocationsOnThisPet, &
coordSys=ESMF_COORDSYS_SPH_DEG, &
rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble creating locStream'
rc=ESMF_FAILURE
return
endif
!-------------------------------------------------------------------
! Add key data (internally allocating memory).
!-------------------------------------------------------------------
call ESMF_LocStreamAddKey(dstLocStream, &
keyName="ESMF:Lat", &
KeyTypeKind=ESMF_TYPEKIND_R8, &
keyUnits="Degrees", &
keyLongName="Latitude", rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble adding LocStream key for Lat'
rc=ESMF_FAILURE
return
endif
call ESMF_LocStreamAddKey(dstLocStream, &
keyName="ESMF:Lon", &
KeyTypeKind=ESMF_TYPEKIND_R8, &
keyUnits="Degrees", &
keyLongName="Longitude", rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble adding LocStream key for Lon'
rc=ESMF_FAILURE
return
endif
!-------------------------------------------------------------------
! Get key data.
!-------------------------------------------------------------------
call ESMF_LocStreamGetKey(dstLocStream, &
localDE=0, &
keyName="ESMF:Lat", &
farray=lat, &
rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble getting LocStream key for Latitude'
rc=ESMF_FAILURE
return
endif
call ESMF_LocStreamGetKey(dstLocStream, &
localDE=0, &
keyName="ESMF:Lon", &
farray=lon, &
rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble getting LocStream key for Longitude'
rc=ESMF_FAILURE
return
endif
!-------------------------------------------------------------------
! Set key data.
!-------------------------------------------------------------------
!setting the coordinates with offset just under limit before nearest neighbor
!regridding would shift to next point
do i=1,numLocationsOnThisPet
lon(i)=(i-1)*360.0/numLocationsOnThisPet+29.999
lat(i)=0.0
enddo
!-------------------------------------------------------------------
! Create a Field on the Location Stream. In this case the
! Field is created from a user array, but any of the other
! Field create methods (e.g. from ArraySpec) would also apply.
!-------------------------------------------------------------------
! Create dest field
call ESMF_ArraySpecSet(arrayspec, 1, ESMF_TYPEKIND_R8, rc=rc)
dstField = ESMF_FieldCreate(dstLocStream, arrayspec, &
name="dest", rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! clear destination Field
! Should only be 1 localDE
call ESMF_FieldGet(dstField, 0, farrayPtr1D, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
farrayPtr1D=0.0
rc=ESMF_SUCCESS
!!! Regrid forward from one locstream to another
! Regrid store
call ESMF_FieldRegridStore( &
srcField, &
dstField=dstField, &
routeHandle=routeHandle, &
regridmethod=ESMF_REGRIDMETHOD_NEAREST_STOD, &
rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! Do regrid
call ESMF_FieldRegrid(srcField, dstField, routeHandle, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
call ESMF_FieldRegridRelease(routeHandle, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
!check error
call ESMF_FieldGet(dstField, 0, farrayPtr1D, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
! loop through nodes and make sure interpolated values are reasonable
! interpolated point should always map to same index
do i=1,numLocationsOnThisPet
if ( abs( farrayPtr1D(i)-temperature(i) ) > 0.0001) then
correct=.false.
write(*,*) localPet,i,"::",farrayPtr1D(i),temperature(i)
endif
enddo
deallocate(temperature)
! Destroy the Fields
call ESMF_FieldDestroy(srcField, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
call ESMF_FieldDestroy(dstField, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
rc=ESMF_FAILURE
return
endif
!destroy locStream objects
call ESMF_LocStreamDestroy(srcLocStream, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble destroying location stream'
rc=ESMF_FAILURE
return
endif
call ESMF_LocStreamDestroy(dstLocStream, rc=localrc)
if (localrc /=ESMF_SUCCESS) then
print*,'ERROR: trouble destroying location stream'
rc=ESMF_FAILURE
return
endif
! return answer based on correct flag
if (correct) then
rc=ESMF_SUCCESS
else
rc=ESMF_FAILURE
endif
end subroutine test_regridNearestLocStreamToLocStream