subroutine ESMF_RegridWeightGenDG(srcFile, dstFile, regridRouteHandle, &
keywordEnforcer, srcElementDistgrid, dstElementDistgrid, &
srcNodalDistgrid, dstNodalDistgrid, &
weightFile, regridmethod, lineType, normType, &
extrapMethod, extrapNumSrcPnts, extrapDistExponent, extrapNumLevels,&
unmappedaction, ignoreDegenerate, useUserAreaFlag, &
largefileFlag, netcdf4fileFlag, &
weightOnlyFlag, verboseFlag, rc)
! !ARGUMENTS:
character(len=*), intent(in) :: srcFile
character(len=*), intent(in) :: dstFile
type(ESMF_RouteHandle), intent(out) :: regridRouteHandle
type(ESMF_KeywordEnforcer), optional:: keywordEnforcer ! must use keywords below
type(ESMF_DistGrid), intent(in), optional :: srcElementDistgrid
type(ESMF_DistGrid), intent(in), optional :: dstElementDistgrid
character(len=*), intent(in), optional :: weightFile
type(ESMF_DistGrid), intent(in), optional :: srcNodalDistgrid
type(ESMF_DistGrid), intent(in), optional :: dstNodalDistgrid
type(ESMF_RegridMethod_Flag), intent(in), optional :: regridmethod
type(ESMF_LineType_Flag), intent(in), optional :: lineType
type(ESMF_NormType_Flag), intent(in), optional :: normType
type(ESMF_ExtrapMethod_Flag), intent(in), optional :: extrapMethod
integer, intent(in), optional :: extrapNumSrcPnts
real, intent(in), optional :: extrapDistExponent
integer, intent(in), optional :: extrapNumLevels
type(ESMF_UnmappedAction_Flag),intent(in), optional :: unmappedaction
logical, intent(in), optional :: ignoreDegenerate
logical, intent(in), optional :: useUserAreaFlag
logical, intent(in), optional :: largefileFlag
logical, intent(in), optional :: netcdf4fileFlag
logical, intent(in), optional :: weightOnlyFlag
logical, intent(in), optional :: verboseFlag
integer, intent(out), optional :: rc
! !DESCRIPTION:
! This subroutine does online regridding weight generation from files with user specified distribution.
! The main differences between this API and the one in \ref{api:esmf_regridweightgenfile} are listed below:
! \begin{itemize}
! \item The input grids are always represented as {\tt ESMF\_Mesh} whether they are logically rectangular or unstructured.
! \item The input grids will be decomposed using a user-specified distribution instead of a fixed decomposition in the
! other subroutine if {\tt srcElementDistgrid} and {\tt dstElementDistgrid} are specified.
! \item The source and destination grid files have to be in the SCRIP grid file format.
! \item This subroutine has one additional required argument {\tt regridRouteHandle} and four additional optional
! arguments: {\tt srcElementDistgrid}, {\tt dstElementDistgrid}, {\tt srcNodelDistgrid} and {\tt dstNodalDistgrid}.
! These four arguments are of type {\tt ESMF\_DistGrid}, they are used to define the distribution of the source
! and destination grid elements and nodes. The output {\tt regridRouteHandle} allows users to regrid the field
! values later in the application.
! \item The {\tt weightFile} argument is optional. When it is given, a weightfile will be generated as well.
! \end{itemize}
! \smallskip
!
! The arguments are:
! \begin{description}
! \item [srcFile]
! The source grid file name in SCRIP grid file format
! \item [dstFile]
! The destination grid file name in SCRIP grid file format
! \item [regridRouteHandle]
! The regrid RouteHandle returned by {\tt ESMF\_FieldRegridStore()}
! \item [srcElementDistgrid]
! An optional distGrid that specifies the distribution of the source grid's elements. If not
! specified, a system-defined block decomposition is used.
! \item [dstElementDistgrid]
! An optional distGrid that specifies the distribution of the destination grid's elements. If
! not specified, a system-defined block decomposition is used.
! \item [weightFile]
! The interpolation weight file name. If present, an output weight file will be generated.
! \item [srcNodalDistgrid]
! An optional distGrid that specifies the distribution of the source grid's nodes
! \item [dstNodalDistgrid]
! An optional distGrid that specifies the distribution of the destination grid's nodes
! \item [{[regridmethod]}]
! The type of interpolation. Please see Section~\ref{opt:regridmethod}
! for a list of valid options. If not specified, defaults to
! {\tt ESMF\_REGRIDMETHOD\_BILINEAR}.
! \item [{[lineType]}]
! This argument controls the path of the line which connects two points on a sphere surface. This in
! turn controls the path along which distances are calculated and the shape of the edges that make
! up a cell. Both of these quantities can influence how interpolation weights are calculated.
! As would be expected, this argument is only applicable when {\tt srcField} and {\tt dstField} are
! built on grids which lie on the surface of a sphere. Section~\ref{opt:lineType} shows a
! list of valid options for this argument. If not specified, the default depends on the
! regrid method. Section~\ref{opt:lineType} has the defaults by line type. Figure~\ref{line_type_support} shows
! which line types are supported for each regrid method as well as showing the default line type by regrid method.
! \item [{[normType]}]
! This argument controls the type of normalization used when generating conservative weights. This option
! only applies to weights generated with {\tt regridmethod=ESMF\_REGRIDMETHOD\_CONSERVE}. Please see
! Section~\ref{opt:normType} for a
! list of valid options. If not specified {\tt normType} defaults to {\tt ESMF\_NORMTYPE\_DSTAREA}.
! \item [{[extrapMethod]}]
! The type of extrapolation. Please see Section~\ref{opt:extrapmethod}
! for a list of valid options. If not specified, defaults to
! {\tt ESMF\_EXTRAPMETHOD\_NONE}.
! \item [{[extrapNumSrcPnts]}]
! The number of source points to use for the extrapolation methods that use more than one source point
! (e.g. {\tt ESMF\_EXTRAPMETHOD\_NEAREST\_IDAVG}). If not specified, defaults to 8..
! \item [{[extrapDistExponent]}]
! The exponent to raise the distance to when calculating weights for
! the {\tt ESMF\_EXTRAPMETHOD\_NEAREST\_IDAVG} extrapolation method. A higher value reduces the influence
! of more distant points. If not specified, defaults to 2.0.
! \item [{[unmappedaction]}]
! Specifies what should happen if there are destination points that
! can't be mapped to a source cell. Please see Section~\ref{const:unmappedaction} for a
! list of valid options. If not specified, {\tt unmappedaction} defaults to {\tt ESMF\_UNMAPPEDACTION\_ERROR}.
! \item [{[ignoreDegenerate]}]
! Ignore degenerate cells when checking the input Grids or Meshes for errors. If this is set to true, then the
! regridding proceeds, but degenerate cells will be skipped. If set to false, a degenerate cell produces an error.
! If not specified, {\tt ignoreDegenerate} defaults to false.
! \item [{[useUserAreaFlag]}]
! If .TRUE., the element area values defined in the grid files are used.
! Only the SCRIP and ESMF format grid files have user specified areas. This flag
! is only used for conservative regridding. The default is .FALSE.
! \item [{[largefileFlag]}]
! If .TRUE., the output weight file is in NetCDF 64bit offset format.
! The default is .FALSE.
! \item [{[netcdf4fileFlag]}]
! If .TRUE., the output weight file is in NetCDF4 file format.
! The default is .FALSE.
! \item [{[weightOnlyFlag]}]
! If .TRUE., the output weight file only contains factorList and factorIndexList.
! The default is .FALSE.
! \item [{[verboseFlag]}]
! If .TRUE., it will print summary information about the regrid parameters,
! default to .FALSE.
! \item [{[rc]}]
! Return code; equals {\tt ESMF\_SUCCESS} if there are no errors.
! \end{description}
!EOP
type(ESMF_RegridMethod_Flag) :: localRegridMethod
logical :: localUserAreaFlag
logical :: localLargefileFlag
logical :: localNetcdf4fileFlag
logical :: localWeightOnlyFlag
logical :: localVerboseFlag
integer :: localrc
type(ESMF_VM) :: vm
integer :: PetNo, PetCnt
type(ESMF_Mesh) :: srcMesh, dstMesh
type(ESMF_Field) :: srcField, dstField
type(ESMF_Field) :: srcFracField, dstFracField
type(ESMF_ArraySpec) :: arrayspec
integer(ESMF_KIND_I4) :: maskvals(1)
integer(ESMF_KIND_I4), pointer:: factorIndexList(:,:)
real(ESMF_KIND_R8), pointer :: factorList(:)
integer :: ind
logical :: convert3D
logical :: isConserve
logical :: convertToDual
type(ESMF_MeshLoc) :: meshloc
character(len=256) :: methodStr
real(ESMF_KIND_R8), pointer :: srcArea(:)
real(ESMF_KIND_R8), pointer :: dstArea(:)
real(ESMF_KIND_R8), pointer :: dstFrac(:), srcFrac(:)
integer :: regridScheme
logical :: wasCompacted
integer(ESMF_KIND_I4), pointer:: compactedFactorIndexList(:,:)
real(ESMF_KIND_R8), pointer :: compactedFactorList(:)
type(ESMF_UnmappedAction_Flag) :: localUnmappedaction
character(len=256) :: argStr
!real(ESMF_KIND_R8) :: starttime, endtime
type(ESMF_LineType_Flag) :: localLineType
type(ESMF_NormType_Flag):: localNormType
logical :: localIgnoreDegenerate
#ifdef ESMF_NETCDF
!------------------------------------------------------------------------
! get global vm information
!
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=PetNo, petCount=PetCnt, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! Default values
localRegridMethod = ESMF_REGRIDMETHOD_BILINEAR
localVerboseFlag = .false.
localLargeFileFlag = .false.
localNetcdf4FileFlag = .false.
localWeightOnlyFlag = .false.
localUserAreaflag = .false.
localIgnoreDegenerate = .false.
if (present(regridMethod)) then
localRegridMethod = regridMethod
endif
if (present(unmappedaction)) then
localUnmappedaction = unmappedaction
else
localUnmappedaction = ESMF_UNMAPPEDACTION_ERROR
endif
if (present(ignoreDegenerate)) then
localIgnoreDegenerate = ignoreDegenerate
endif
if (present(largefileFlag)) then
localLargeFileFlag = largefileFlag
endif
if (present(netcdf4fileFlag)) then
localNetcdf4FileFlag = netcdf4fileFlag
endif
if (present(weightOnlyFlag)) then
localWeightOnlyFlag = weightOnlyFlag
endif
if (present(useUserAreaFlag)) then
localUserAreaFlag = useUserAreaFlag
endif
if (present(verboseFlag)) then
localVerboseFlag = verboseFlag
endif
! Handle optional normType argument
if (present(normType)) then
localNormType=normType
else
localNormType=ESMF_NORMTYPE_DSTAREA
endif
! Handle optional lineType argument
if (present(lineType)) then
localLineType=lineType
else
if (localRegridMethod == ESMF_REGRIDMETHOD_CONSERVE) then
localLineType=ESMF_LINETYPE_GREAT_CIRCLE
else if (localRegridMethod == ESMF_REGRIDMETHOD_CONSERVE_2ND) then
localLineType=ESMF_LINETYPE_GREAT_CIRCLE
else
localLineType=ESMF_LINETYPE_CART
endif
endif
! user area only needed for conservative regridding
if (localUserAreaFlag .and. .not. ((localRegridMethod == ESMF_REGRIDMETHOD_CONSERVE) .or. &
(localRegridMethod == ESMF_REGRIDMETHOD_CONSERVE_2ND))) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg = " user defined area is only used for the conservative regridding", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
! Print the regrid options
if (localVerboseFlag .and. PetNo == 0) then
print *, "Starting weight generation with these inputs: "
print *, " Source File: ", trim(srcfile)
print *, " Destination File: ", trim(dstfile)
if (present(weightFile)) then
print *, " Weight File: ", trim(weightFile)
if (localWeightOnlyFlag) then
print *, " only output weights in the weight file"
endif
endif
if (localRegridMethod == ESMF_REGRIDMETHOD_BILINEAR) then
print *, " Regrid Method: bilinear"
elseif (localRegridMethod == ESMF_REGRIDMETHOD_CONSERVE) then
print *, " Regrid Method: conserve"
elseif (localRegridMethod == ESMF_REGRIDMETHOD_CONSERVE_2ND) then
print *, " Regrid Method: conserve2nd"
elseif (localRegridMethod == ESMF_REGRIDMETHOD_PATCH) then
print *, " Regrid Method: patch"
elseif (localRegridMethod == ESMF_REGRIDMETHOD_NEAREST_STOD) then
print *, " Regrid Method: nearest source to destination"
elseif (localRegridMethod == ESMF_REGRIDMETHOD_NEAREST_DTOS) then
print *, " Regrid Method: nearest destination to source"
endif
if (localUnmappedaction .eq. ESMF_UNMAPPEDACTION_IGNORE) then
print *, " Ignore unmapped destination points"
endif
if (localLargeFileFlag) then
print *, " Output weight file in 64bit offset NetCDF file format"
endif
if (localNetcdf4FileFlag) then
print *, " Output weight file in NetCDF4 file format"
endif
if (localUserAreaFlag) then
print *, " Use user defined cell area for both the source and destination grids"
endif
if (localLineType .eq. ESMF_LINETYPE_CART) then
print *, " Line Type: cartesian"
elseif (localLineType .eq. ESMF_LINETYPE_GREAT_CIRCLE) then
print *, " Line Type: greatcircle"
endif
if (localNormType .eq. ESMF_NORMTYPE_DSTAREA) then
print *, " Norm Type: dstarea"
elseif (localNormType .eq. ESMF_NORMTYPE_FRACAREA) then
print *, " Norm Type: fracarea"
endif
if (present(extrapMethod)) then
if (extrapMethod%extrapmethod .eq. &
ESMF_EXTRAPMETHOD_NONE%extrapmethod) then
print *, " Extrap. Method: none"
else if (extrapMethod%extrapmethod .eq. &
ESMF_EXTRAPMETHOD_NEAREST_STOD%extrapmethod) then
print *, " Extrap. Method: neareststod"
else if (extrapMethod%extrapmethod .eq. &
ESMF_EXTRAPMETHOD_NEAREST_IDAVG%extrapmethod) then
print *, " Extrap. Method: nearestidavg"
if (present(extrapNumSrcPnts)) then
print '(a,i0)', " Extrap. Number of Source Points: ",extrapNumSrcPnts
else
print '(a,i0)', " Extrap. Number of Source Points: ",8
endif
if (present(extrapDistExponent)) then
print *, " Extrap. Dist. Exponent: ",extrapDistExponent
else
print *, " Extrap. Dist. Exponent: ",2.0
endif
else
print *, " Extrap. Method: unknown"
endif
else
print *, " Extrap. Method: none"
endif
write(*,*)
endif
! Set flags according to the regrid method
if ((localRegridMethod == ESMF_REGRIDMETHOD_CONSERVE) .or. &
(localRegridMethod == ESMF_REGRIDMETHOD_CONSERVE_2ND)) then
isConserve=.true.
convertToDual=.false.
meshloc=ESMF_MESHLOC_ELEMENT
else
isConserve=.false.
convertToDual=.true.
meshloc=ESMF_MESHLOC_NODE
endif
regridScheme = ESMF_REGRID_SCHEME_FULL3D
maskvals(1) = 0
!Read in the srcfile and create the corresponding ESMF_Mesh object
if (present(srcElementDistgrid) .and. present(srcNodalDistgrid)) then
srcMesh = ESMF_MeshCreate(srcfile, ESMF_FILEFORMAT_SCRIP, &
convertToDual=convertToDual, addUserArea=localUserAreaFlag, &
elementDistgrid=srcElementDistgrid, nodalDistgrid=srcNodalDistgrid, &
rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
elseif (present(srcElementDistgrid)) then
srcMesh = ESMF_MeshCreate(srcfile, ESMF_FILEFORMAT_SCRIP, &
convertToDual=convertToDual, addUserArea=localUserAreaFlag, &
elementDistgrid=srcElementDistgrid, &
rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
elseif (present(srcNodalDistgrid)) then
srcMesh = ESMF_MeshCreate(srcfile, ESMF_FILEFORMAT_SCRIP, &
convertToDual=convertToDual, addUserArea=localUserAreaFlag, &
nodalDistgrid=srcNodalDistgrid, &
rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
else
srcMesh = ESMF_MeshCreate(srcfile, ESMF_FILEFORMAT_SCRIP, &
convertToDual=convertToDual, addUserArea=localUserAreaFlag, &
rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
!call ESMF_MeshWrite(srcMesh, "srcMesh", rc)
call ESMF_ArraySpecSet(arrayspec, 1, ESMF_TYPEKIND_R8, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
srcField=ESMF_FieldCreate(srcMesh,arrayspec,meshloc=meshloc,rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
!Read in the dstfile and create the corresponding ESMF object (either
! ESMF_Grid or ESMF_Mesh)
if (present(dstElementDistgrid) .and. present(dstNodalDistgrid)) then
dstMesh = ESMF_MeshCreate(dstfile, ESMF_FILEFORMAT_SCRIP, &
convertToDual=convertToDual, addUserArea=localUserAreaFlag, &
elementDistgrid=dstElementDistgrid, nodalDistgrid=dstNodalDistgrid, &
rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
elseif (present(dstElementDistgrid)) then
dstMesh = ESMF_MeshCreate(dstfile, ESMF_FILEFORMAT_SCRIP, &
convertToDual=convertToDual, addUserArea=localUserAreaFlag, &
elementDistgrid=dstElementDistgrid, &
rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
elseif (present(dstNodalDistgrid)) then
dstMesh = ESMF_MeshCreate(dstfile, ESMF_FILEFORMAT_SCRIP, &
convertToDual=convertToDual, addUserArea=localUserAreaFlag, &
nodalDistgrid=dstNodalDistgrid, &
rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
else
dstMesh = ESMF_MeshCreate(dstfile, ESMF_FILEFORMAT_SCRIP, &
convertToDual=convertToDual, addUserArea=localUserAreaFlag, &
rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
call ESMF_ArraySpecSet(arrayspec, 1, ESMF_TYPEKIND_R8, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
dstField=ESMF_FieldCreate(dstMesh,arrayspec,meshloc=meshloc,rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! Create Frac Fields if conservative
if (isConserve) then
call ESMF_ArraySpecSet(arrayspec, 1, ESMF_TYPEKIND_R8, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
srcFracField=ESMF_FieldCreate(srcMesh,arrayspec,meshloc=meshloc,rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_ArraySpecSet(arrayspec, 1, ESMF_TYPEKIND_R8, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
dstFracField=ESMF_FieldCreate(dstMesh,arrayspec,meshloc=meshloc,rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
!call ESMF_VMBarrier(vm)
!call ESMF_VMWtime(starttime, rc=localrc)
if (localRegridMethod == ESMF_REGRIDMETHOD_BILINEAR) then
methodStr = "Bilinear remapping"
else if (localRegridMethod == ESMF_REGRIDMETHOD_PATCH) then
methodStr = "Bilinear remapping" ! SCRIP doesn't recognize Patch
else if (localRegridMethod == ESMF_REGRIDMETHOD_NEAREST_STOD) then
methodStr = "Bilinear remapping" ! SCRIP doesn't recognize Nearest neighbor
else if (localRegridMethod == ESMF_REGRIDMETHOD_NEAREST_DTOS) then
methodStr = "Bilinear remapping" ! SCRIP doesn't recognize Nearest neighbor
else if (localRegridMethod == ESMF_REGRIDMETHOD_CONSERVE) then
methodStr = "Conservative remapping"
else if (localRegridMethod == ESMF_REGRIDMETHOD_CONSERVE_2ND) then
methodStr = "Conservative remapping" ! SCRIP doesn't recognize 2nd order
else ! nothing recognizable so report error
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_WRONG, &
msg="unrecognized RegridMethod", &
ESMF_CONTEXT, rcToReturn=rc)
endif
if (present(weightFile)) then
call ESMF_FieldRegridStore(srcField=srcField, dstField=dstField, &
srcMaskValues = maskvals, dstMaskValues = maskvals, &
unmappedaction=localUnmappedaction, &
ignoreDegenerate=localIgnoreDegenerate, &
routehandle=regridRouteHandle, &
factorIndexList=factorIndexList, factorList=factorList, &
srcFracField=srcFracField, dstFracField=dstFracField, &
regridmethod = localRegridMethod, &
polemethod = ESMF_POLEMETHOD_NONE, &
lineType=locallineType, &
normType=localNormType, &
extrapMethod=extrapMethod, &
extrapNumSrcPnts=extrapNumSrcPnts, &
extrapDistExponent=extrapDistExponent, &
extrapNumLevels=extrapNumLevels, &
rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
else
call ESMF_FieldRegridStore(srcField=srcField, dstField=dstField, &
srcMaskValues = maskvals, dstMaskValues = maskvals, &
unmappedaction=localUnmappedaction, &
ignoreDegenerate=localIgnoreDegenerate, &
routehandle=regridRouteHandle, &
srcFracField=srcFracField, dstFracField=dstFracField, &
regridmethod = localRegridMethod, &
polemethod = ESMF_POLEMETHOD_NONE, &
lineType=locallineType, &
normType=localNormType, &
extrapMethod=extrapMethod, &
extrapNumSrcPnts=extrapNumSrcPnts, &
extrapDistExponent=extrapDistExponent, &
extrapNumLevels=extrapNumLevels, &
rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
! Only compute area, fraction and output weight file is weightFile is present
if (present(weightFile) .and. .not. localWeightOnlyFlag) then
! Compute areas if conservative
! Area only valid on PET 0 right now, when parallel Array
! write works, then make area io parallel
if (isConserve) then
call computeRedistAreaMesh(srcMesh, vm, petNo, petCnt, srcArea, localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call computeRedistAreaMesh(dstMesh, vm, petNo, petCnt, dstArea, localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
call compactMatrix(factorList, factorIndexList, &
wasCompacted, &
compactedFactorList, compactedFactorIndexList, &
localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! If the list was compacted get rid of the old lists and
! point to the new lists
if (wasCompacted) then
deallocate(factorList)
factorList=>compactedFactorList
deallocate(factorIndexList)
factorIndexList=>compactedFactorIndexList
endif
! Computer fraction if bilinear
! src fraction is always 0
! destination fraction depends on the src mask, dst mask, and the weight
if ((localRegridMethod /= ESMF_REGRIDMETHOD_CONSERVE) .and. &
(localRegridMethod /= ESMF_REGRIDMETHOD_CONSERVE_2ND)) then
call computeFracMesh(dstMesh, vm, factorIndexList, dstFrac, localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
else
call gatherRedistFracFieldMesh(srcMesh, vm, srcFracField, petNo, petCnt, &
srcFrac, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call gatherRedistFracFieldMesh(dstMesh, vm, dstFracField, petNo, petCnt, &
dstFrac, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
!! Write the weight table into a SCRIP format NetCDF file
if (PetNo == 0) then
if (isConserve) then
call ESMF_OutputScripWeightFile(weightFile, factorList, factorIndexList, &
srcFile=srcfile, dstFile=dstfile, method = localRegridMethod, &
srcArea=srcArea, dstArea=dstArea, srcFrac=srcFrac, &
normType=localNormType, &
dstFrac=dstFrac, largeFileFlag=localLargeFileFlag, &
netcdf4FileFlag = localNetcdf4FileFlag, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
else
call ESMF_OutputScripWeightFile(weightFile, factorList, factorIndexList, &
srcFile=srcfile, dstFile=dstfile, &
normType=localNormType, &
method = localRegridMethod, dstFrac=dstFrac, &
largeFileFlag=localLargeFileFlag, &
netcdf4FileFlag = localNetcdf4FileFlag, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
else
call ESMF_OutputScripWeightFile(weightFile, factorList, factorIndexList, &
normType=localNormType, rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
!call ESMF_VMBarrier(vm)
!call ESMF_VMWtime(endtime, rc=localrc)
! Get rid of conservative arrays
if (isConserve) then
if (PetNo == 0) then
deallocate(srcArea)
deallocate(dstArea)
endif
endif
if (PetNo == 0) then
if (isConserve) deallocate(srcFrac)
deallocate(dstFrac)
endif
else if (present(weightFile)) then !localWeightOnlyFlag = .TRUE.
call compactMatrix(factorList, factorIndexList, &
wasCompacted, &
compactedFactorList, compactedFactorIndexList, &
localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! If the list was compacted get rid of the old lists and
! point to the new lists
if (wasCompacted) then
deallocate(factorList)
factorList=>compactedFactorList
deallocate(factorIndexList)
factorIndexList=>compactedFactorIndexList
endif
! write simple weight file
call ESMF_OutputSimpleWeightFile(weightFile, factorList, factorIndexList, &
title = "ESMF Regrid Weight Generator", &
method = localRegridMethod, &
largeFileFlag=localLargeFileFlag, &
netcdf4FileFlag = localNetcdf4FileFlag, &
rc=localrc)
if (ESMF_LogFoundError(localrc, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
! clean up
call ESMF_FieldDestroy(srcField)
call ESMF_FieldDestroy(dstField)
if (isConserve) then
call ESMF_FieldDestroy(srcFracField)
call ESMF_FieldDestroy(dstFracField)
endif
! ESMF_MeshDestory() will destroy the distgrid passed in as input argument
! Work Around: use ESMF_MeshFreeMemory() instead
! call ESMF_MeshDestroy(srcMesh)
! call ESMF_MeshDestroy(dstMesh)
call ESMF_MeshFreeMemory(srcMesh)
call ESMF_MeshFreeMemory(dstMesh)
if (present(weightFile)) then
deallocate(factorList, factorIndexList)
endif
rc = ESMF_SUCCESS
return
#else
call ESMF_LogSetError(rcToCheck=ESMF_RC_LIB_NOT_PRESENT, &
msg="- ESMF_NETCDF not defined when lib was compiled", &
ESMF_CONTEXT, rcToReturn=rc)
return
#endif
end subroutine ESMF_RegridWeightGenDG