subroutine f_esmf_regridstorefile(srcField, dstField, fileName, &
srcMaskValues, len1, &
dstMaskValues, len2, &
routehandle, &
regridmethod, &
polemethod, &
regridPoleNPnts, &
linetype, &
normtype, &
unmappedaction, &
ignoreDegenerate, &
createRoutehandle, &
filemode, &
srcFile, &
dstFile, &
srcFileType, &
dstFileType, &
largeFileFlag, &
srcFracField, &
dstFracField, &
rc)
use ESMF_UtilTypesMod
use ESMF_BaseMod
use ESMF_LogErrMod
use ESMF_RHandleMod
use ESMF_FieldRegridMod
use ESMF_FieldMod
use ESMF_FieldCreateMod
use ESMF_FieldGetMod
use ESMF_IOScripMod
use ESMF_GeomMod
use ESMF_GridMod
use ESMF_MeshMod
use ESMF_LocStreamMod
use ESMF_XGridMod
use ESMF_StaggerLocMod
use ESMF_VMMod
use ESMF_UtilRWGMod
implicit none
type(ESMF_Field) :: srcField
type(ESMF_Field) :: dstField
character(*), intent(in) :: fileName
integer :: len1, len2
integer,optional :: srcMaskValues(len1), &
dstMaskValues(len2)
type(ESMF_RouteHandle) :: routehandle
type(ESMF_RegridMethod_Flag) :: regridmethod
type(ESMF_PoleMethod_Flag) :: polemethod
integer :: regridPoleNPnts
type(ESMF_LineType_Flag) :: linetype
type(ESMF_NormType_Flag) :: normtype
type(ESMF_UnmappedAction_Flag) :: unmappedaction
logical :: ignoreDegenerate
logical, optional :: createRoutehandle
type(ESMF_FileMode_Flag), optional :: filemode
character(len=*), optional :: srcFile
character(len=*), optional :: dstFile
type(ESMF_FileFormat_Flag), optional :: srcFileType
type(ESMF_FileFormat_Flag), optional :: dstFileType
logical, optional :: largeFileFlag
real(ESMF_KIND_R8), pointer :: srcArea(:), dstArea(:)
type(ESMF_GeomType_Flag) :: srcgt, dstgt
type(ESMF_TypeKind_Flag) :: srctk, dsttk
type(ESMF_Grid) :: srcgrid, dstgrid
type(ESMF_Mesh) :: srcmesh, dstmesh
integer :: srcslc, dstslc
logical :: ecip
type(ESMF_Field) :: srcFracField
type(ESMF_Field) :: dstFracField
type(ESMF_VM) :: vm
integer :: localPet, petCount
integer :: rc
integer :: localrc
type(ESMF_RouteHandle) :: l_routehandle
real(ESMF_KIND_R8), pointer :: localFactorList(:)
integer(ESMF_KIND_I4), pointer :: localFactorIndexList(:,:)
type(ESMF_FileMode_Flag) :: filemode_local
! initialize return code; assume routine not implemented
rc = ESMF_RC_NOT_IMPL
localrc = ESMF_RC_NOT_IMPL
call ESMF_VMGetCurrent(vm, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! set up local pet info
call ESMF_VMGet(vm, localPet=localPet, petCount=petCount, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (present (createRoutehandle)) then
if (createRoutehandle .eqv. .false.) then
call ESMF_FieldRegridStore(srcField, dstField, &
srcMaskValues=srcMaskValues, &
dstMaskValues=dstMaskValues, &
regridmethod=regridmethod, &
polemethod=polemethod, &
regridPoleNPnts=regridPoleNPnts, &
lineType=linetype, &
normType=normtype, &
unmappedaction=unmappedaction, &
ignoreDegenerate=ignoreDegenerate, &
factorList=localFactorList, &
factorIndexList=localFactorIndexList, &
srcFracField=srcFracField, &
dstFracField=dstFracField, &
rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
else
call ESMF_FieldRegridStore(srcField, dstField, &
srcMaskValues=srcMaskValues, &
dstMaskValues=dstMaskValues, &
regridmethod=regridmethod, &
polemethod=polemethod, &
regridPoleNPnts=regridPoleNPnts, &
lineType=linetype, &
normType=normtype, &
unmappedaction=unmappedaction, &
ignoreDegenerate=ignoreDegenerate, &
routehandle=l_routehandle, &
factorList=localFactorList, &
factorIndexList=localFactorIndexList, &
srcFracField=srcFracField, &
dstFracField=dstFracField, &
rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
else
call ESMF_FieldRegridStore(srcField, dstField, &
srcMaskValues=srcMaskValues, &
dstMaskValues=dstMaskValues, &
regridmethod=regridmethod, &
polemethod=polemethod, &
regridPoleNPnts=regridPoleNPnts, &
lineType=linetype, &
normType=normtype, &
unmappedaction=unmappedaction, &
ignoreDegenerate=ignoreDegenerate, &
routehandle=l_routehandle, &
factorList=localFactorList, &
factorIndexList=localFactorIndexList, &
srcFracField=srcFracField, &
dstFracField=dstFracField, &
rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
! write the weights to file
filemode_local = ESMF_FILEMODE_BASIC
if (present(filemode)) then
filemode_local = filemode
endif
if (filemode_local == ESMF_FILEMODE_BASIC) then
call ESMF_SparseMatrixWrite(localFactorList, localFactorIndexList, &
fileName, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
elseif (filemode_local == ESMF_FILEMODE_WITHAUX) then
! query field for geom type
call ESMF_FieldGet(srcField, geomType=srcgt, typekind=srctk, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! determine which stagger locations are available
srcslc = 0
if (srcgt == ESMF_GEOMTYPE_GRID) then
call ESMF_FieldGet(srcField, grid=srcgrid, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridGetCoord(srcgrid, staggerloc=ESMF_STAGGERLOC_CENTER, &
isPresent=ecip, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (ecip .eqv. .true.) srcslc = 1
call ESMF_GridGetCoord(srcgrid, staggerloc=ESMF_STAGGERLOC_CORNER, &
isPresent=ecip, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (ecip .eqv. .true.) srcslc = 2
else if (srcgt == ESMF_GEOMTYPE_MESH) then
call ESMF_FieldGet(srcField, mesh=srcmesh, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
ecip = .false.
call ESMF_MeshGet(srcmesh, elementCoordsIsPresent=ecip, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
srcslc = 1
if (ecip .eqv. .true.) srcslc = 2
else if (srcgt == ESMF_GEOMTYPE_XGRID) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_NOT_IMPL, &
msg="- xgrid cannot retrieve areas", &
ESMF_CONTEXT, rcToReturn=rc)
endif
! query field for geom type
call ESMF_FieldGet(dstField, geomType=dstgt, typekind=dsttk, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! determine which stagger locations are available
dstslc = 0
if (dstgt == ESMF_GEOMTYPE_GRID) then
call ESMF_FieldGet(dstField, grid=dstgrid, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridGetCoord(dstgrid, staggerloc=ESMF_STAGGERLOC_CENTER, &
isPresent=ecip, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (ecip .eqv. .true.) dstslc = 1
call ESMF_GridGetCoord(dstgrid, staggerloc=ESMF_STAGGERLOC_CORNER, &
isPresent=ecip, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (ecip .eqv. .true.) dstslc = 2
else if (dstgt == ESMF_GEOMTYPE_MESH) then
call ESMF_FieldGet(dstField, mesh=dstmesh, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
ecip = .false.
call ESMF_MeshGet(dstmesh, elementCoordsIsPresent=ecip, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
dstslc = 1
if (ecip .eqv. .true.) dstslc = 2
else if (dstgt == ESMF_GEOMTYPE_XGRID) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_NOT_IMPL, &
msg="- xgrid cannot retrieve areas", &
ESMF_CONTEXT, rcToReturn=rc)
endif
! compute the areas
if (srcslc > 1) then
if (srcgt == ESMF_GEOMTYPE_GRID) then
call computeAreaGrid(srcgrid, localPet, srcarea, 0, localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
else if (srcgt == ESMF_GEOMTYPE_MESH) then
call computeAreaMesh(srcmesh, vm, localPet, petCount, srcarea, localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
endif
if (dstslc > 1) then
if (dstgt == ESMF_GEOMTYPE_GRID) then
call computeAreaGrid(dstgrid, petCount, dstarea, 0, localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
else if (dstgt == ESMF_GEOMTYPE_MESH) then
call computeAreaMesh(dstmesh, vm, localPet, petCount, dstarea, localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
endif
! write the weight file
if (srcslc > 1 .and. dstslc > 1) then
call ESMF_OutputScripWeightFile(fileName, &
localFactorList, localFactorIndexList, &
srcFile=srcFile, dstFile=dstFile, &
srcFileType=srcFileType, &
dstFileType=dstFileType, &
srcArea=srcArea, &
dstArea=dstArea, &
largeFileFlag=largeFileFlag, &
rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
else if (srcslc > 1 .and. dstslc == 1) then
call ESMF_OutputScripWeightFile(fileName, &
localFactorList, localFactorIndexList, &
srcFile=srcFile, dstFile=dstFile, &
srcFileType=srcFileType, &
dstFileType=dstFileType, &
srcArea=srcArea, &
largeFileFlag=largeFileFlag, &
rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
else if (srcslc == 1 .and. dstslc > 1) then
call ESMF_OutputScripWeightFile(fileName, &
localFactorList, localFactorIndexList, &
srcFile=srcFile, dstFile=dstFile, &
srcFileType=srcFileType, &
dstFileType=dstFileType, &
dstArea=dstArea, &
largeFileFlag=largeFileFlag, &
rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
else if (srcslc == 1 .and. dstslc == 1) then
call ESMF_OutputScripWeightFile(fileName, &
localFactorList, localFactorIndexList, &
srcFile=srcFile, dstFile=dstFile, &
srcFileType=srcFileType, &
dstFileType=dstFileType, &
largeFileFlag=largeFileFlag, &
rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
else
call ESMF_LogSetError(rcToCheck=ESMF_RC_VAL_OUTOFRANGE, &
msg="- unrecognized area field options", &
ESMF_CONTEXT, rcToReturn=rc)
endif
else
call ESMF_LogSetError(rcToCheck=ESMF_RC_VAL_OUTOFRANGE, &
msg="- filemode not recognized", &
ESMF_CONTEXT, rcToReturn=rc)
endif
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! because ESMF_RouteHandle.this is private, it cannot be accessed directly
! we use the public interface to do the ptr copy;
! the RouteHandle object returned to the C interface must consist only of
! the 'this' pointer. It must not contain the isInit member.
if (present (createRoutehandle)) then
if (createRoutehandle .eqv. .true.) then
call ESMF_RoutehandleCopyThis(l_routehandle, routehandle, localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
else
call ESMF_RoutehandleCopyThis(l_routehandle, routehandle, localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
deallocate(localFactorList)
deallocate(localFactorIndexList)
rc = ESMF_SUCCESS
end subroutine f_esmf_regridstorefile