subroutine ESMF_OutputScripWeightFile (wgtFile, factorList, factorIndexList, &
srcFile, dstFile, srcFileType, dstFileType,&
title, method, normType, srcArea, dstArea, &
srcFrac, dstFrac, largeFileFlag, netcdf4FileFlag, &
srcmeshname, dstmeshname, &
srcMissingValue, dstMissingValue, &
srcvarname, dstvarname, &
useSrcCorner, useDstCorner, &
srccoordnames, dstcoordnames, tileFilePath, rc)
!
! !ARGUMENTS:
character(len=*), intent(in) :: wgtFile
real(ESMF_KIND_R8) , intent(in) :: factorList(:)
integer(ESMF_KIND_I4) , intent(in), target :: factorIndexList(:,:)
character(len=*), optional, intent(in) :: srcFile
character(len=*), optional, intent(in) :: dstFile
type(ESMF_FileFormat_Flag), optional, intent(in) :: srcFileType
type(ESMF_FileFormat_Flag), optional, intent(in) :: dstFileType
character(len=*), optional, intent(in) :: title
type(ESMF_RegridMethod_Flag), optional, intent(in) :: method
type(ESMF_NormType_Flag), intent(in), optional :: normType
real(ESMF_KIND_R8),optional, intent(in) :: srcArea(:),dstArea(:)
real(ESMF_KIND_R8),optional, intent(in) :: srcFrac(:), dstFrac(:)
logical, optional, intent(in) :: largeFileFlag, netcdf4FileFlag
character(len=*), optional, intent(in) :: srcmeshname, dstmeshname
logical, optional, intent(in) :: srcMissingValue, dstMissingValue
character(len=*), optional, intent(in) :: srcvarname, dstvarname
character(len=*), optional, intent(in) :: srccoordnames(:), dstcoordnames(:)
logical, optional, intent(in) :: useSrcCorner, useDstCorner
character(len=*), optional, intent(in) :: tileFilePath
integer, optional :: rc
type(ESMF_VM):: vm
integer :: PetNo, PetCnt
integer :: total, localCount(1)
integer :: ncid, ncid1
integer :: ncStatus
integer :: status
integer :: i,j, k, start
integer :: srcDim, dstDim
integer:: naDimId, nbDimId, nsDimId, srankDimId, drankDimId, varId, varId1, varId2
integer :: nvaDimId, nvbDimId
integer :: src_grid_rank, dst_grid_rank, src_grid_corner, dst_grid_corner
integer, pointer :: src_grid_dims(:), dst_grid_dims(:)
integer :: src_ndims, dst_ndims
integer :: src_coordids(2), dst_coordids(2)
integer :: src_dimids(2), dst_dimids(2)
real(ESMF_KIND_R8), pointer :: coords(:),area(:),frac(:)
real(ESMF_KIND_R8), pointer :: latBuffer(:),lonBuffer(:)
real(ESMF_KIND_R8), pointer :: latBuffer1D(:),lonBuffer1D(:)
real(ESMF_KIND_R8), pointer :: latBuffer2(:,:),lonBuffer2(:,:)
real(ESMF_KIND_R8), pointer :: cornerlon2D(:,:), cornerlat2D(:,:)
real(ESMF_KIND_R8), pointer :: cornerlon3D(:,:,:), cornerlat3D(:,:,:)
real(ESMF_KIND_R8), pointer :: weightbuf(:), varBuffer(:,:)
real(ESMF_KIND_R8), pointer :: varBuffer1D(:)
real(ESMF_KIND_R8) :: missing_value
integer(ESMF_KIND_I4), pointer:: indexbuf(:), next(:)
integer(ESMF_KIND_I4), pointer:: mask(:)
integer(ESMF_KIND_I4), pointer:: allCounts(:)
character(len=256) :: titlelocal, norm, map_method
character(len=256) :: esmf_regrid_method,conventions
type(ESMF_RegridMethod_Flag) :: methodlocal
character(len=80) :: srcunits, dstunits
integer :: maxcount
type(ESMF_FileFormat_Flag) :: srcFileTypeLocal, dstFileTypeLocal
integer :: srcNodeDim, dstNodeDim, srcCoordDim, dstCoordDim
integer, parameter :: nf90_noerror = 0
character(len=256) :: errmsg
character(len=20) :: varStr
type(ESMF_Logical) :: largeFileFlaglocal
type(ESMF_Logical) :: netcdf4FileFlaglocal
logical :: faceCoordFlag
logical :: srchasbound, dsthasbound
logical :: useSrcCornerlocal, useDstCornerlocal
type(ESMF_NormType_Flag):: localNormType
logical :: src_has_area, dst_has_area
type(ESMF_Mosaic) :: srcmosaic, dstmosaic
character(len=ESMF_MAXPATHLEN) :: tempname
integer :: totalsize, totallen
integer :: meshId
integer :: memstat
integer :: dim
type(ESMF_CoordSys_Flag):: coordsys
#ifdef ESMF_NETCDF
! write out the indices and weights table sequentially to the output file
! first find out the starting index of my portion of table
! Global reduce
src_has_area = .false.
dst_has_area = .false.
coordsys = ESMF_COORDSYS_SPH_DEG
if (present(srcFileType)) then
srcFileTypeLocal = srcFileType
else
srcFileTypeLocal = ESMF_FILEFORMAT_SCRIP
endif
if (present(dstFileType)) then
dstFileTypeLocal = dstFileType
else
dstFileTypeLocal = ESMF_FILEFORMAT_SCRIP
endif
if (present(largeFileFlag)) then
largeFileFlaglocal = largeFileFlag
else
largeFileFlaglocal = .false.
endif
if (present(netcdf4FileFlag)) then
netcdf4FileFlaglocal = netcdf4FileFlag
else
netcdf4FileFlaglocal = .false.
endif
! Handle optional normType argument
if (present(normType)) then
localNormType=normType
else
localNormType=ESMF_NORMTYPE_DSTAREA
endif
if (present(useSrcCorner)) then
useSrcCornerlocal=useSrcCorner
else
useSrcCornerlocal=.false.
endif
if (present(useDstCorner)) then
useDstCornerlocal=useDstCorner
else
useDstCornerlocal=.false.
endif
call ESMF_VMGetCurrent(vm, rc=rc)
if (rc /= ESMF_SUCCESS) return
! set up local pet info
call ESMF_VMGet(vm, localPet=PetNo, petCount=PetCnt, rc=rc)
if (rc /= ESMF_SUCCESS) return
localCount(1)=size(factorList,1)
allocate(allCounts(PetCnt), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_VMAllGather(vm,localCount,allCounts,1,rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! calculate the size of the global weight table
total = 0
do i=1,PetCnt
total=allCounts(i)+total
end do
!print *, PetNo, 'local count ', localCount(1), AllCounts(PetNo+1), total
!Read the variables from the input grid files at PET0
if (PetNo == 0) then
! Check if srcFile and dstFile exists
if (.not. present(srcFile)) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_BAD, &
msg="- The srcFile argument does not exist on PET0 ", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
if (.not. present(dstFile)) then
call ESMF_LogSetError(rcToCheck=ESMF_RC_ARG_BAD, &
msg="- The dstFile argument does not exist on PET0 ", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
! Create output file and create dimensions and variables
call c_nc_create(wgtFile, NF90_CLOBBER, &
largeFileFlaglocal, netcdf4FileFlaglocal, ncid, status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! global variables
if (present(title)) then
titlelocal = trim(title)
else
titlelocal = "ESMF Offline Regridding Weight Generator"
endif
! Regrid method
if (present(method)) then
methodlocal = method
if (methodlocal%regridmethod == ESMF_REGRIDMETHOD_BILINEAR%regridmethod) then
map_method = "Bilinear remapping"
esmf_regrid_method = "Bilinear"
elseif (methodlocal%regridmethod == ESMF_REGRIDMETHOD_PATCH%regridmethod) then
!scrip_test does not recognize patch remapping
map_method = "Bilinear remapping"
esmf_regrid_method = "Higher-order Patch"
elseif (methodlocal%regridmethod == ESMF_REGRIDMETHOD_CONSERVE%regridmethod) then
map_method = "Conservative remapping"
esmf_regrid_method = "First-order Conservative"
elseif (methodlocal%regridmethod == ESMF_REGRIDMETHOD_CONSERVE_2ND%regridmethod) then
map_method = "Conservative remapping"
esmf_regrid_method = "Second-order Conservative"
elseif (methodlocal%regridmethod == ESMF_REGRIDMETHOD_NEAREST_STOD%regridmethod) then
map_method = "Bilinear remapping"
esmf_regrid_method = "Nearest source to destination"
elseif (methodlocal%regridmethod == ESMF_REGRIDMETHOD_NEAREST_DTOS%regridmethod) then
map_method = "Bilinear remapping"
esmf_regrid_method = "Nearest destination to source"
else
!report error
call ESMF_LogSetError(rcToCheck=ESMF_FAILURE, &
msg="- regridmethod not recongized", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
else
methodlocal = ESMF_REGRIDMETHOD_BILINEAR
map_method = "Bilinear remapping"
esmf_regrid_method = "Bilinear"
endif
! Norm type
if (methodlocal%regridmethod == ESMF_REGRIDMETHOD_CONSERVE%regridmethod .or. &
methodlocal%regridmethod == ESMF_REGRIDMETHOD_CONSERVE_2ND%regridmethod) then
if (localNormType .eq. ESMF_NORMTYPE_DSTAREA) then
norm = "destarea"
elseif (localNormType .eq. ESMF_NORMTYPE_FRACAREA) then
norm = "fracarea"
else
norm = "unknown"
endif
else
! For regrid methods other than conservative, the normalization type is irrelevant
norm = "N/A"
end if
conventions = "NCAR-CSM"
ncStatus = nf90_put_att(ncid, NF90_GLOBAL, "title", trim(titlelocal))
errmsg = "Attribute title in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus = nf90_put_att(ncid, NF90_GLOBAL, "normalization", trim(norm))
errmsg = "Attribute normalization in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus = nf90_put_att(ncid, NF90_GLOBAL, "map_method", trim(map_method))
errmsg = "Attribute map_method in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus = nf90_put_att(ncid, NF90_GLOBAL, "ESMF_regrid_method", &
trim(esmf_regrid_method))
errmsg = "Attribute esmf_regrid_method in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus = nf90_put_att(ncid, NF90_GLOBAL, "conventions", trim(conventions))
errmsg = "Attribute conventions in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus = nf90_put_att(ncid, NF90_GLOBAL, "domain_a", trim(srcFile))
errmsg = "Attribute domain_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus = nf90_put_att(ncid, NF90_GLOBAL, "domain_b", trim(dstFile))
errmsg = "Attribute domain_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus = nf90_put_att(ncid, NF90_GLOBAL, "grid_file_src", trim(srcFile))
errmsg = "Attribute grid_file_src in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus = nf90_put_att(ncid, NF90_GLOBAL, "grid_file_dst", trim(dstFile))
errmsg = "Attribute grid_file_dst in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
! Output version information
#ifdef ESMF_VERSION_STRING_GIT
ncStatus = nf90_put_att(ncid, NF90_GLOBAL, "ESMF_version", ESMF_VERSION_STRING_GIT)
#else
ncStatus = nf90_put_att(ncid, NF90_GLOBAL, "ESMF_version", "(No version information available.)")
#endif
errmsg = "Attribute CVS_revision in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
! Get Source Grid dimension and variables
if (srcFileTypeLocal == ESMF_FILEFORMAT_SCRIP) then
allocate(src_grid_dims(2), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_ScripInq(srcFile, grid_rank=src_grid_rank, grid_size=srcDim, &
grid_dims=src_grid_dims, grid_corners=src_grid_corner, &
has_area=src_has_area, rc=status)
! The grid_dims for an unstructured grie (grid_rank = 1) is not used
! by SCRIP, thus some of the SCRIP files did not set this value correctly.
! This causes the scrip_test generating wrong output data, i.e.
! the grid coordinates will not be written in the weight file
if (src_grid_rank == 1) src_grid_dims(1) = srcDim
call ESMF_ScripInqUnits(srcFile,units = srcunits, rc=status)
else if (srcFileTypeLocal == ESMF_FILEFORMAT_GRIDSPEC) then
allocate(src_grid_dims(2), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridspecInq(srcFile, src_ndims, src_grid_dims, &
dimids = src_dimids, coordids = src_coordids, &
coord_names = srccoordnames, hasbound=srchasbound, &
units=srcunits,rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
src_grid_rank = 2
srcDim = src_grid_dims(1)*src_grid_dims(2)
if (src_ndims == 1) then
src_grid_corner = 2
else
src_grid_corner = 4
endif
else if (srcFileTypeLocal == ESMF_FILEFORMAT_TILE) then
allocate(src_grid_dims(2), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
! this returns the size of the center stagger, not the supergrid
call ESMF_GridSpecQueryTileSize(srcfile, src_grid_dims(1), src_grid_dims(2), &
units=srcunits, rc=status)
if (ESMF_LogFoundError(status, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
src_grid_rank = 2
srcDim = src_grid_dims(1)*src_grid_dims(2)
src_grid_corner = 4
else if (srcFileTypeLocal == ESMF_FILEFORMAT_ESMFMESH) then
! If bilinear, we have to switch node and elment, so the nodeCount became srcDim and
! elementCount becomes srcNodeDim. Hard code src_grid_corner to 3. The xv_a and xv_b
! will be empty
if (useSrcCornerlocal .and. &
(methodlocal%regridmethod ==ESMF_REGRIDMETHOD_BILINEAR%regridmethod .or. &
methodlocal%regridmethod ==ESMF_REGRIDMETHOD_PATCH%regridmethod .or. &
methodlocal%regridmethod ==ESMF_REGRIDMETHOD_NEAREST_STOD%regridmethod .or. &
methodlocal%regridmethod ==ESMF_REGRIDMETHOD_NEAREST_DTOS%regridmethod)) then
call ESMF_EsmfInq(srcFile, nodeCount=srcDim, &
coordDim = srcCoordDim, elementCount=srcNodeDim, rc=status)
src_grid_corner =3
else
call ESMF_EsmfInq(srcFile, elementCount=srcDim, maxNodePElement=src_grid_corner, &
coordDim = srcCoordDim, nodeCount=srcNodeDim, haveArea=src_has_area, rc=status)
endif
call ESMF_EsmfInqUnits(srcFile,units = srcunits, rc=status)
allocate(src_grid_dims(1), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
src_grid_dims(1)=srcDim
src_grid_rank = 1
else if (srcFileTypeLocal == ESMF_FILEFORMAT_UGRID) then
if (useSrcCornerlocal .and. &
(methodlocal%regridmethod ==ESMF_REGRIDMETHOD_BILINEAR%regridmethod .or. &
methodlocal%regridmethod ==ESMF_REGRIDMETHOD_PATCH%regridmethod .or. &
methodlocal%regridmethod ==ESMF_REGRIDMETHOD_NEAREST_STOD%regridmethod .or. &
methodlocal%regridmethod ==ESMF_REGRIDMETHOD_NEAREST_DTOS%regridmethod)) then
call ESMF_UGridInq(srcFile, srcmeshname, nodeCount=srcDim, &
nodeCoordDim=dim, units=srcunits, rc=status)
! If it is 1D network topology, there is no corner coordinates
if (dim==1) then
src_grid_corner = 0
else
src_grid_corner =3
endif
else
call ESMF_UGridInq(srcFile, srcmeshname, elementCount=srcDim, &
maxNodePElement=src_grid_corner, units=srcunits, &
nodeCount = srcNodeDim, rc=status)
endif
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
srcCoordDim = 2
allocate(src_grid_dims(1), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
src_grid_dims(1)=srcDim
src_grid_rank = 1
else if (srcFileTypeLocal == ESMF_FILEFORMAT_MOSAIC) then
call ESMF_GridSpecReadMosaic(srcFile, srcmosaic, tileFilePath=tileFilePath, rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
srcCoordDim = 2
src_grid_rank = 1
src_grid_corner = 0
allocate(src_grid_dims(1))
srcDim = srcmosaic%nx * srcmosaic%ny * srcmosaic%ntiles
src_grid_dims(1)=srcDim
srcunits = 'degrees'
endif
if (dstFileTypelocal == ESMF_FILEFORMAT_SCRIP) then
allocate(dst_grid_dims(2), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_ScripInq(dstFile, grid_rank=dst_grid_rank, grid_size=dstDim, &
grid_dims=dst_grid_dims, grid_corners=dst_grid_corner, &
has_area=dst_has_area, rc=status)
! The grid_dims for an unstructured grie (grid_rank = 1) is not used
! by SCRIP, thus some of the SCRIP files did not set this value correctly.
! This causes the scrip_test generating wrong output data, i.e.
! the grid coordinates will not be written in the weight file
if (dst_grid_rank == 1) dst_grid_dims(1) = dstDim
call ESMF_ScripInqUnits(dstFile,units = dstunits, rc=status)
else if (dstFileTypeLocal == ESMF_FILEFORMAT_GRIDSPEC) then
allocate(dst_grid_dims(2), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridspecInq(dstFile, dst_ndims, dst_grid_dims, &
dimids = dst_dimids, coordids = dst_coordids, &
coord_names = dstcoordnames, hasbound=dsthasbound, &
units=dstunits, rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
dst_grid_rank = 2
dstDim = dst_grid_dims(1)*dst_grid_dims(2)
if (dst_ndims == 1) then
dst_grid_corner = 2
else
dst_grid_corner = 4
endif
else if (dstFileTypeLocal == ESMF_FILEFORMAT_TILE) then
allocate(dst_grid_dims(2), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
! this returns the size of the center stagger, not the supergrid
call ESMF_GridSpecQueryTileSize(dstfile, dst_grid_dims(1),dst_grid_dims(2), &
units=dstunits, rc=status)
if (ESMF_LogFoundError(status, &
ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
dst_grid_rank = 2
dstDim = dst_grid_dims(1)*dst_grid_dims(2)
dst_grid_corner = 4
else if (dstFileTypeLocal == ESMF_FILEFORMAT_ESMFMESH) then
! If bilinear, we have to switch node and elment, so the nodeCount became dstDim and
! elementCount becomes dstNodeDim. Hard code dst_grid_corner to 3. The xv_a and xv_b
! will be empty
if (useDstCornerlocal .and. &
(methodlocal%regridmethod ==ESMF_REGRIDMETHOD_BILINEAR%regridmethod .or. &
methodlocal%regridmethod ==ESMF_REGRIDMETHOD_PATCH%regridmethod .or. &
methodlocal%regridmethod ==ESMF_REGRIDMETHOD_NEAREST_STOD%regridmethod .or. &
methodlocal%regridmethod ==ESMF_REGRIDMETHOD_NEAREST_DTOS%regridmethod)) then
call ESMF_EsmfInq(dstFile, nodeCount=dstDim, &
coordDim = dstCoordDim, elementCount=dstNodeDim, rc=status)
dst_grid_corner =3
else
call ESMF_EsmfInq(dstFile, elementCount=dstDim, maxNodePElement=dst_grid_corner, &
coordDim = dstCoordDim, nodeCount=dstNodeDim, haveArea=dst_has_area, rc=status)
endif
call ESMF_EsmfInqUnits(dstFile,units = dstunits, rc=status)
allocate(dst_grid_dims(1), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
dst_grid_dims(1)=dstDim
dst_grid_rank = 1
else if (dstFileTypeLocal == ESMF_FILEFORMAT_UGRID) then
if (useDstCornerlocal .and. &
(methodlocal%regridmethod ==ESMF_REGRIDMETHOD_BILINEAR%regridmethod .or. &
methodlocal%regridmethod ==ESMF_REGRIDMETHOD_PATCH%regridmethod .or. &
methodlocal%regridmethod ==ESMF_REGRIDMETHOD_NEAREST_STOD%regridmethod .or. &
methodlocal%regridmethod ==ESMF_REGRIDMETHOD_NEAREST_DTOS%regridmethod)) then
call ESMF_UGridInq(dstFile, dstmeshname, nodeCount=dstDim, &
nodeCoordDim=dim, units=dstunits, rc=status)
if (dim==1) then
dst_grid_corner = 0
else
dst_grid_corner = 3
endif
else
call ESMF_UGridInq(dstFile, dstmeshname, elementCount=dstDim, &
maxNodePElement=dst_grid_corner, units=dstunits, &
nodeCount = dstNodeDim, rc=status)
endif
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
dstCoordDim = 2
allocate(dst_grid_dims(1), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
dst_grid_dims(1)=dstDim
dst_grid_rank = 1
else if (dstFileTypeLocal == ESMF_FILEFORMAT_MOSAIC) then
call ESMF_GridSpecReadMosaic(dstFile, dstmosaic, tileFilePath=tileFilePath, rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
dstCoordDim = 2
dst_grid_rank = 1
dst_grid_corner = 0
allocate(dst_grid_dims(1), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
dstDim = dstmosaic%nx * dstmosaic%ny * dstmosaic%ntiles
dst_grid_dims(1)=dstDim
dstunits = 'degrees'
endif
! define dimensions
ncStatus = nf90_def_dim(ncid,"n_a",srcDim, naDimId)
errmsg = "Dimension n_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus = nf90_def_dim(ncid,"n_b",dstDim, nbDimId)
errmsg = "Dimension n_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus = nf90_def_dim(ncid,"n_s",total, nsDimId)
errmsg = "Dimension n_s in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
! if src_grid_corner < 0 -- the input file is a ESMF Mesh file with ragged array
! do not define nv_a, xv_a and yv_a
if (src_grid_corner > 0) then
! define max number of vertices
ncStatus = nf90_def_dim(ncid,"nv_a",src_grid_corner, nvaDimId)
errmsg = "Dimension nv_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
endif
if (dst_grid_corner > 0) then
ncStatus = nf90_def_dim(ncid,"nv_b",dst_grid_corner, nvbDimId)
errmsg = "Dimension nv_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
endif
! define max number of vertices
ncStatus = nf90_def_dim(ncid,"num_wgts",1, VarId)
errmsg = "Dimension num_wgts in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
! define grid ranks
ncStatus = nf90_def_dim(ncid,"src_grid_rank",src_grid_rank, srankDimId)
errmsg = "Dimension src_grid_rank in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus = nf90_def_dim(ncid,"dst_grid_rank",dst_grid_rank, drankDimId)
errmsg = "Dimension dst_grid_rank in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
! define variables
! Grid Dims
ncStatus = nf90_def_var(ncid,"src_grid_dims",NF90_INT, (/srankDimId/), VarId)
errmsg = "Variable src_grid_dims in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus = nf90_def_var(ncid,"dst_grid_dims",NF90_INT, (/drankDimId/), VarId)
errmsg = "Variable dst_grid_dims in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
! yc_a: source vertex coordinate (latitude)
ncStatus = nf90_def_var(ncid,"yc_a",NF90_DOUBLE, (/naDimId/), VarId)
errmsg = "Variable yc_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus = nf90_put_att(ncid, VarId, "units", trim(srcunits))
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
! yc_b: destination vertex coordinate (latitude)
ncStatus = nf90_def_var(ncid,"yc_b",NF90_DOUBLE, (/nbDimId/), VarId)
errmsg = "Variable yc_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus = nf90_put_att(ncid, VarId, "units", trim(dstunits))
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
! xc_a: source vertex coordinate (longitude)
ncStatus = nf90_def_var(ncid,"xc_a",NF90_DOUBLE, (/naDimId/), VarId)
errmsg = "Variable xc_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus = nf90_put_att(ncid, VarId, "units", trim(srcunits))
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
! xc_b: dest. vertex coordinate (longitude)
ncStatus = nf90_def_var(ncid,"xc_b",NF90_DOUBLE, (/nbDimId/), VarId)
errmsg = "Variable xc_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus = nf90_put_att(ncid, VarId, "units", trim(dstunits))
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
if (src_grid_corner > 0) then
! yv_a: source corner coordinate (latitude)
ncStatus = nf90_def_var(ncid,"yv_a",NF90_DOUBLE, (/nvaDimId,naDimId/), VarId)
errmsg = "Variable yv_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus = nf90_put_att(ncid, VarId, "units", trim(srcunits))
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
! xv_a: source corner coordinate (longitude)
ncStatus = nf90_def_var(ncid,"xv_a",NF90_DOUBLE, (/nvaDimId,naDimId/), VarId)
errmsg = "Variable xv_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus = nf90_put_att(ncid, VarId, "units", trim(srcunits))
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
endif
if (dst_grid_corner > 0) then
! yv_b: source corner coordinate (latitude)
ncStatus = nf90_def_var(ncid,"yv_b",NF90_DOUBLE, (/nvbDimId,nbDimId/), VarId)
errmsg = "Variable yv_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus = nf90_put_att(ncid, VarId, "units", trim(dstunits))
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
! xv_b: source corner coordinate (longitude)
ncStatus = nf90_def_var(ncid,"xv_b",NF90_DOUBLE, (/nvbDimId,nbDimId/), VarId)
errmsg = "Variable xv_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus = nf90_put_att(ncid, VarId, "units", trim(dstunits))
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
endif
! mask_a
ncStatus = nf90_def_var(ncid,"mask_a",NF90_INT, (/naDimId/), VarId)
errmsg = "Variable mask_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus = nf90_put_att(ncid, VarId, "units", "unitless")
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
! mask_b
ncStatus = nf90_def_var(ncid,"mask_b",NF90_INT, (/nbDimId/), VarId)
errmsg = "Variable mask_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus = nf90_put_att(ncid, VarId, "units", "unitless")
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
! area_a
ncStatus = nf90_def_var(ncid,"area_a",NF90_DOUBLE, (/naDimId/), VarId)
errmsg = "Variable area_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus = nf90_put_att(ncid, VarId, "units", "square radians")
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
! area_b
ncStatus = nf90_def_var(ncid,"area_b",NF90_DOUBLE, (/nbDimId/), VarId)
errmsg = "Variable area_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus = nf90_put_att(ncid, VarId, "units", "square radians")
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
! frac_a
ncStatus = nf90_def_var(ncid,"frac_a",NF90_DOUBLE, (/naDimId/), VarId)
errmsg = "Variable frac_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus = nf90_put_att(ncid, VarId, "units", "unitless")
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
! frac_b
ncStatus = nf90_def_var(ncid,"frac_b",NF90_DOUBLE, (/nbDimId/), VarId)
errmsg = "Variable frac_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus = nf90_put_att(ncid, VarId, "units", "unitless")
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
! col: sparse matrix weight table
ncStatus = nf90_def_var(ncid,"col",NF90_INT, (/nsDimId/), VarId)
errmsg = "Variable col in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
! row: sparse matrix weight table
ncStatus = nf90_def_var(ncid,"row",NF90_INT, (/nsDimId/), VarId)
errmsg = "Variable row in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
! S: sparse matrix weight table
ncStatus = nf90_def_var(ncid,"S",NF90_DOUBLE, (/nsDimId/), VarId)
errmsg = "Variable S in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus=nf90_enddef(ncid)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
trim(wgtFile),&
rc)) return
! write src_grid_dims and dst_grid_dims
ncStatus=nf90_inq_varid(ncid,"src_grid_dims",VarId)
ncStatus=nf90_put_var(ncid,VarId, src_grid_dims(1:src_grid_rank))
errmsg = "Variable src_grid_dims in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus=nf90_inq_varid(ncid,"dst_grid_dims",VarId)
ncStatus=nf90_put_var(ncid,VarId, dst_grid_dims(1:dst_grid_rank))
errmsg = "Variable dst_grid_dims in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
if (srcFileTypeLocal == ESMF_FILEFORMAT_SCRIP) then
! Read the srcGrid variables and write them out
allocate(latBuffer(srcDim), lonBuffer(srcDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_ScripGetVar(srcFile, grid_center_lon=lonBuffer, &
grid_center_lat=latBuffer, rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
ncStatus=nf90_inq_varid(ncid,"xc_a",VarId)
ncStatus=nf90_put_var(ncid,VarId, lonBuffer)
errmsg = "Variable xc_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus=nf90_inq_varid(ncid,"yc_a",VarId)
ncStatus=nf90_put_var(ncid,VarId, latBuffer)
errmsg = "Variable yc_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
deallocate(latBuffer, lonBuffer, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
! Write xv_a, yv_a
allocate(latBuffer2(src_grid_corner,srcDim),lonBuffer2(src_grid_corner,srcDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_ScripGetVar(srcFile, grid_corner_lon=lonBuffer2, &
grid_corner_lat=latBuffer2, rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
ncStatus=nf90_inq_varid(ncid,"xv_a",VarId)
ncStatus=nf90_put_var(ncid,VarId, lonBuffer2)
errmsg = "Variable xv_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus=nf90_inq_varid(ncid,"yv_a",VarId)
ncStatus=nf90_put_var(ncid,VarId, latBuffer2)
errmsg = "Variable yv_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
deallocate(latBuffer2, lonBuffer2, stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
allocate(mask(srcDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_ScripGetVar(srcFile, grid_imask=mask, rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
ncStatus=nf90_inq_varid(ncid,"mask_a",VarId)
ncStatus=nf90_put_var(ncid,VarId, mask)
errmsg = "Variable mask_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
deallocate(mask, stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (src_has_area .and. .not. present(srcArea)) then
allocate(area(srcDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_ScripGetVar(srcFile, grid_area=area, rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
ncStatus=nf90_inq_varid(ncid,"area_a",VarId)
ncStatus=nf90_put_var(ncid,VarId, area)
errmsg = "Variable area_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
deallocate(area, stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
else if (srcFileTypeLocal == ESMF_FILEFORMAT_GRIDSPEC) then
allocate(lonBuffer(srcDim), latBuffer(srcDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (src_ndims == 1) then
allocate(lonBuffer1D(src_grid_dims(1)), latBuffer1D(src_grid_dims(2)), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (srchasbound) then
allocate(cornerlon2D(src_grid_corner,src_grid_dims(1)), &
cornerlat2D(src_grid_corner,src_grid_dims(2)), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridspecGetVar1D(srcFile, src_coordids, lonBuffer1D, latBuffer1D, &
cornerlon = cornerlon2D, cornerlat=cornerlat2D, rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
else
call ESMF_GridspecGetVar1D(srcFile, src_coordids, lonBuffer1D, latBuffer1D, &
rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
k=1
do j=1,src_grid_dims(2)
do i=1,src_grid_dims(1)
lonBuffer(k)=lonBuffer1D(i)
latBuffer(k)=latBuffer1D(j)
k=k+1
enddo
enddo
deallocate(lonBuffer1D, latBuffer1D, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (srchasbound) then
allocate(lonBuffer2(src_grid_corner, srcDim),latBuffer2(src_grid_corner,srcDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
k=1
do j=1,src_grid_dims(2)
do i=1,src_grid_dims(1)
lonBuffer2(:,k)=cornerlon2D(:,i)
latBuffer2(:,k)=cornerlat2D(:,j)
k=k+1
enddo
enddo
endif
else
allocate(lonBuffer2(src_grid_dims(1),src_grid_dims(2)), &
latBuffer2(src_grid_dims(1),src_grid_dims(2)), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (srchasbound) then
allocate(cornerlon3D(src_grid_corner,src_grid_dims(1), src_grid_dims(2)),&
cornerlat3D(src_grid_corner,src_grid_dims(1), src_grid_dims(2)), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridspecGetVar2D(srcFile, src_coordids, lonBuffer2, latBuffer2, &
cornerlon=cornerlon3D, cornerlat=cornerlat3D, rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
else
call ESMF_GridspecGetVar2D(srcFile, src_coordids, lonBuffer2, latBuffer2, &
rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
lonBuffer = reshape(lonBuffer2, (/srcDim/))
latBuffer = reshape(latBuffer2, (/srcDim/))
deallocate(lonBuffer2, latBuffer2, stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (srchasbound) then
allocate(lonBuffer2(src_grid_corner, srcDim),latBuffer2(src_grid_corner,srcDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
lonBuffer2=reshape(cornerlon3D, (/src_grid_corner, srcDim/))
latBuffer2=reshape(cornerlat3D, (/src_grid_corner, srcDim/))
deallocate(cornerlon3D, cornerlat3D, stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
endif
ncStatus=nf90_inq_varid(ncid,"xc_a",VarId)
ncStatus=nf90_put_var(ncid,VarId, lonBuffer)
errmsg = "Variable xc_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus=nf90_inq_varid(ncid,"yc_a",VarId)
ncStatus=nf90_put_var(ncid,VarId, latBuffer)
errmsg = "Variable yc_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
deallocate(latBuffer, lonBuffer, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
! Write xv_a, yv_a
if (srchasbound) then
ncStatus=nf90_inq_varid(ncid,"xv_a",VarId)
ncStatus=nf90_put_var(ncid,VarId, lonBuffer2)
errmsg = "Variable xv_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus=nf90_inq_varid(ncid,"yv_a",VarId)
ncStatus=nf90_put_var(ncid,VarId, latBuffer2)
errmsg = "Variable yv_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
deallocate(lonBuffer2, latBuffer2, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
! Mask
allocate(mask(srcDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
mask(:)=1
if (present(srcMissingValue)) then
if (srcMissingValue) then
allocate(varBuffer(src_grid_dims(1),src_grid_dims(2)), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridspecGetVarByName(srcFile, srcvarname, src_dimids, &
varBuffer, missing_value = missing_value, &
rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
k=1
do i=1,size(varBuffer,2)
do j=1,size(varBuffer,1)
if (varBuffer(j,i) == missing_value) mask(k)=0
k=k+1
enddo
enddo
deallocate(varBuffer, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
endif
ncStatus=nf90_inq_varid(ncid,"mask_a",VarId)
ncStatus=nf90_put_var(ncid,VarId, mask)
errmsg = "Variable mask_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
deallocate(mask, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
else if (srcFileTypeLocal == ESMF_FILEFORMAT_TILE) then
allocate(lonBuffer2(src_grid_dims(1),src_grid_dims(2)), &
latBuffer2(src_grid_dims(1), src_grid_dims(2)), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
allocate(cornerlon2D(src_grid_dims(1)+1,src_grid_dims(2)+1), &
cornerlat2D(src_grid_dims(1)+1, src_grid_dims(2)+1), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridspecReadTile(srcFile, src_grid_dims(1), src_grid_dims(2), &
lonBuffer2, latBuffer2, cornerLon=cornerlon2D, cornerLat=cornerlat2D, &
rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
allocate(lonBuffer(srcDim), latBuffer(srcDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
lonBuffer = reshape(lonBuffer2, (/srcDim/))
latBuffer = reshape(latBuffer2, (/srcDim/))
deallocate(lonBuffer2, latBuffer2, stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
allocate(lonBuffer2(src_grid_corner, srcDim),latBuffer2(src_grid_corner,srcDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
k=1
do j=1,src_grid_dims(2)
do i=1,src_grid_dims(1)
lonBuffer2(1,k)=cornerlon2D(i,j)
lonBuffer2(2,k)=cornerlon2D(i,j+1)
lonBuffer2(3,k)=cornerlon2D(i+1,j+1)
lonBuffer2(4,k)=cornerlon2D(i+1,j)
latBuffer2(1,k)=cornerlat2D(i,j)
latBuffer2(2,k)=cornerlat2D(i,j+1)
latBuffer2(3,k)=cornerlat2D(i+1,j+1)
latBuffer2(4,k)=cornerlat2D(i+1,j)
k=k+1
enddo
enddo
deallocate(cornerlon2D, cornerlat2D, stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
ncStatus=nf90_inq_varid(ncid,"xc_a",VarId)
ncStatus=nf90_put_var(ncid,VarId, lonBuffer)
errmsg = "Variable xc_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus=nf90_inq_varid(ncid,"yc_a",VarId)
ncStatus=nf90_put_var(ncid,VarId, latBuffer)
errmsg = "Variable yc_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
deallocate(latBuffer, lonBuffer, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
! Write xv_a, yv_a
ncStatus=nf90_inq_varid(ncid,"xv_a",VarId)
ncStatus=nf90_put_var(ncid,VarId, lonBuffer2)
errmsg = "Variable xv_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus=nf90_inq_varid(ncid,"yv_a",VarId)
ncStatus=nf90_put_var(ncid,VarId, latBuffer2)
errmsg = "Variable yv_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
deallocate(lonBuffer2, latBuffer2, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
else if (srcFileTypeLocal == ESMF_FILEFORMAT_ESMFMESH) then
! ESMF unstructured grid
ncStatus=nf90_open(srcFile,NF90_NOWRITE,ncid1)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
trim(srcFile),&
rc)) return
if (useSrcCornerlocal .and. (methodlocal%regridmethod ==ESMF_REGRIDMETHOD_BILINEAR%regridmethod &
.or. methodlocal%regridmethod ==ESMF_REGRIDMETHOD_PATCH%regridmethod &
.or. methodlocal%regridmethod ==ESMF_REGRIDMETHOD_NEAREST_STOD%regridmethod &
.or. methodlocal%regridmethod ==ESMF_REGRIDMETHOD_NEAREST_DTOS%regridmethod)) then
! check if centerCoords exit
ncStatus=nf90_inq_varid(ncid1,"nodeCoords",VarId)
varStr = "nodeCoords"
src_has_area = .false.
else
ncStatus=nf90_inq_varid(ncid1,"centerCoords",VarId)
varStr = "centerCoords"
endif
if (ncStatus /= nf90_noerror) then
write(*,*) "Warning: "//trim(varStr)// &
" not present in src grid file, so not outputting xc_a and yc_a to weight file."
write(*,*)
else
allocate(latBuffer2(srcCoordDim, srcDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
allocate(latBuffer(srcDim), lonBuffer(srcDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
ncStatus=nf90_get_var(ncid1,VarId, latBuffer2)
errmsg = "Variable "//trim(varStr)//" in "//trim(srcFile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
do i=1,srcDim
lonBuffer(i)=latBuffer2(1,i)
latBuffer(i)=latBuffer2(2,i)
enddo
ncStatus=nf90_inq_varid(ncid,"xc_a",VarId)
ncStatus=nf90_put_var(ncid,VarId, lonBuffer)
errmsg = "Variable xc_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus=nf90_inq_varid(ncid,"yc_a",VarId)
ncStatus=nf90_put_var(ncid,VarId, latBuffer)
errmsg = "Variable yc_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
deallocate(latBuffer, lonBuffer, latBuffer2, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
! only write out xv_a and yv_a when the regrid method is conserve
if (.not. useSrcCornerlocal .or. &
methodlocal%regridmethod == ESMF_REGRIDMETHOD_CONSERVE%regridmethod .or. &
methodlocal%regridmethod ==ESMF_REGRIDMETHOD_CONSERVE_2ND%regridmethod) then
! output xv_a and yv_a is harder, we have to read in the nodeCoords and
! elementConn and construct the the latitudes and longitudes for
! all the corner vertices
if (src_grid_corner > 0) then
allocate(latBuffer2(src_grid_corner,srcDim),lonBuffer2(src_grid_corner,srcDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_EsmfGetVerts(ncid1, srcFile, srcDim, src_grid_corner, srcNodeDim, &
latBuffer2, lonBuffer2,status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
ncStatus=nf90_inq_varid(ncid,"xv_a",VarId)
ncStatus=nf90_put_var(ncid,VarId, lonBuffer2)
errmsg = "Variable xv_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE, errmsg, &
rc)) return
ncStatus=nf90_inq_varid(ncid,"yv_a",VarId)
ncStatus=nf90_put_var(ncid,VarId, latBuffer2)
errmsg = "Variable yv_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
deallocate(latBuffer2, lonBuffer2, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
endif
allocate(mask(srcDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (.not. useSrcCornerlocal .and. &
(methodlocal%regridmethod ==ESMF_REGRIDMETHOD_CONSERVE%regridmethod .or. &
methodlocal%regridmethod ==ESMF_REGRIDMETHOD_CONSERVE_2ND%regridmethod)) then
ncStatus=nf90_inq_varid(ncid1,"elementMask",VarId)
if (ncStatus /= nf90_noerror) then
write(*,*) "Warning: elementMask"// &
" not present in src grid file, so setting mask_a=1 in weight file."
write(*,*)
mask = 1
else
ncStatus=nf90_get_var(ncid1,VarId, mask)
errmsg = "Variable elementMask in "//trim(srcFile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
end if
else
mask = 1
endif
ncStatus=nf90_inq_varid(ncid,"mask_a",VarId)
ncStatus=nf90_put_var(ncid,VarId, mask)
errmsg = "Variable mask_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
deallocate(mask, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (src_has_area .and. .not. present(srcArea)) then
allocate(area(srcDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
ncStatus=nf90_inq_varid(ncid1,"elementArea",VarId)
errmsg = "Variable elementArea in "//trim(srcFile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
ncStatus=nf90_get_var(ncid1,VarId, area)
errmsg = "Variable elementArea in "//trim(srcFile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
ncStatus=nf90_inq_varid(ncid,"area_a",VarId)
ncStatus=nf90_put_var(ncid,VarId, area)
errmsg = "Variable area_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
deallocate(area, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
ncStatus=nf90_close(ncid1)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,trim(srcFile),&
rc)) return
else if (srcFileTypeLocal == ESMF_FILEFORMAT_UGRID) then
call ESMF_UGridInq(srcfile, srcmeshname, meshId=meshId, faceCoordFlag=faceCoordFlag)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (.not. useSrcCornerlocal .or. &
methodlocal%regridmethod ==ESMF_REGRIDMETHOD_CONSERVE%regridmethod .or. &
methodlocal%regridmethod ==ESMF_REGRIDMETHOD_CONSERVE_2ND%regridmethod) then
! check if faceCoords exit
allocate(latBuffer2(src_grid_corner,srcDim),&
lonBuffer2(src_grid_corner,srcDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (faceCoordFlag) then
allocate(latBuffer(srcDim), lonBuffer(srcDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_UGridGetVar(srcfile, meshId, &
faceXcoords=lonBuffer, faceYcoords=latBuffer, &
faceNodeConnX=lonBuffer2, faceNodeConnY=latBuffer2, rc=status)
else
write(*,*) "Warning: face coordinates not present in src grid file,"// &
" so not outputting xc_a and yc_a to weight file."
write(*,*)
call ESMF_UGridGetVar(srcfile, meshId, &
faceNodeConnX=lonBuffer2, faceNodeConnY=latBuffer2, rc=status)
endif
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (faceCoordFlag) then
ncStatus=nf90_inq_varid(ncid,"xc_a",VarId)
ncStatus=nf90_put_var(ncid,VarId, lonBuffer)
errmsg = "Variable xc_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
ncStatus=nf90_inq_varid(ncid,"yc_a",VarId)
ncStatus=nf90_put_var(ncid,VarId, latBuffer)
errmsg = "Variable yc_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
deallocate(latBuffer, lonBuffer, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
ncStatus=nf90_inq_varid(ncid,"xv_a",VarId)
ncStatus=nf90_put_var(ncid,VarId, lonBuffer2)
errmsg = "Variable xv_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
ncStatus=nf90_inq_varid(ncid,"yv_a",VarId)
ncStatus=nf90_put_var(ncid,VarId, latBuffer2)
errmsg = "Variable yv_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
deallocate(latBuffer2, lonBuffer2, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
else
allocate(latBuffer(srcDim), lonBuffer(srcDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_UGridGetVar(srcfile, meshId, &
nodeXcoords=lonBuffer, nodeYcoords=latBuffer, rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
ncStatus=nf90_inq_varid(ncid,"xc_a",VarId)
ncStatus=nf90_put_var(ncid,VarId, lonBuffer)
errmsg = "Variable xc_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
ncStatus=nf90_inq_varid(ncid,"yc_a",VarId)
ncStatus=nf90_put_var(ncid,VarId, latBuffer)
errmsg = "Variable yc_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
deallocate(latBuffer, lonBuffer, stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
! Write out mask
allocate(mask(srcDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
mask(:)=1
if (present(srcMissingValue)) then
if (srcMissingValue) then
allocate(varBuffer1D(srcDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_UgridGetVarByName(srcFile, srcvarname, &
varBuffer1D, missingvalue = missing_value, &
rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
do j=1,size(varBuffer1D)
if (varBuffer1D(j) == missing_value) mask(j)=0
enddo
deallocate(varBuffer1D, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
endif
ncStatus=nf90_inq_varid(ncid,"mask_a",VarId)
ncStatus=nf90_put_var(ncid,VarId, mask)
errmsg = "Variable mask_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
deallocate(mask, stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
else if (srcFileTypeLocal == ESMF_FILEFORMAT_MOSAIC) then
!read coordinates from the tile files
allocate(latBuffer2(srcmosaic%nx, srcmosaic%ny), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
allocate(lonBuffer2(srcmosaic%nx, srcmosaic%ny), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
totalsize = srcmosaic%nx*srcmosaic%ny
allocate(varBuffer1D(totalsize), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
ncStatus=nf90_inq_varid(ncid,"xc_a",VarId1)
ncStatus=nf90_inq_varid(ncid,"yc_a",VarId2)
errmsg = "Variable xc_a or yc_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
do k=1,srcmosaic%ntiles
totallen = len_trim(srcmosaic%filenames(k))+len_trim(srcmosaic%tileDirectory)
tempname = trim(srcmosaic%tileDirectory)//trim(srcmosaic%filenames(k))
call ESMF_GridspecReadTile(trim(tempname),srcmosaic%nx, srcmosaic%ny, &
centerlon=lonBuffer2, centerlat=latBuffer2, rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
varBuffer1D=reshape(lonBuffer2,(/totalsize/))
start=(k-1)*totalsize+1
ncStatus = nf90_put_var(ncid, VarId1,varBuffer1D, (/start/), (/totalsize/))
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
"writing xc_a in weight file", &
rc)) return
varBuffer1D=reshape(latBuffer2,(/totalsize/))
ncStatus = nf90_put_var(ncid, VarId2, varBuffer1D, (/start/), (/totalsize/))
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
"writing xc_a in weight file", &
rc)) return
enddo
deallocate(lonBuffer2, latBuffer2, varBuffer1D, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
! Read the dstGrid variables and write them out
if (dstFileTypeLocal == ESMF_FILEFORMAT_SCRIP) then
allocate(latBuffer(dstDim), lonBuffer(dstDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_ScripGetVar(dstFile, grid_center_lon=lonBuffer, &
grid_center_lat=latBuffer, rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
ncStatus=nf90_inq_varid(ncid,"xc_b",VarId)
ncStatus=nf90_put_var(ncid,VarId, lonBuffer)
errmsg = "Variable xc_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
ncStatus=nf90_inq_varid(ncid,"yc_b",VarId)
ncStatus=nf90_put_var(ncid,VarId, latBuffer)
errmsg = "Variable yc_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
deallocate(latBuffer, lonBuffer, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
allocate(latBuffer2(dst_grid_corner,dstDim),lonBuffer2(dst_grid_corner,dstDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_ScripGetVar(dstFile, grid_corner_lon=lonBuffer2, &
grid_corner_lat=latBuffer2, rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
ncStatus=nf90_inq_varid(ncid,"xv_b",VarId)
ncStatus=nf90_put_var(ncid,VarId, lonBuffer2)
errmsg = "Variable xv_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
ncStatus=nf90_inq_varid(ncid,"yv_b",VarId)
ncStatus=nf90_put_var(ncid,VarId, latBuffer2)
errmsg = "Variable yv_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
deallocate(latBuffer2, lonBuffer2, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
allocate(mask(dstDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_ScripGetVar(dstFile, grid_imask=mask, rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
ncStatus=nf90_inq_varid(ncid,"mask_b",VarId)
ncStatus=nf90_put_var(ncid,VarId, mask)
errmsg = "Variable mask_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
deallocate(mask, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (dst_has_area .and. .not. present(dstArea)) then
allocate(area(dstDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_ScripGetVar(dstFile, grid_area=area, rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
ncStatus=nf90_inq_varid(ncid,"area_b",VarId)
ncStatus=nf90_put_var(ncid,VarId, area)
errmsg = "Variable area_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
deallocate(area, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
else if (dstFileTypeLocal == ESMF_FILEFORMAT_GRIDSPEC) then
allocate(lonBuffer(dstDim), latBuffer(dstDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
! check if bound variables exist or not
if (dst_ndims == 1) then
allocate(lonBuffer1D(dst_grid_dims(1)), latBuffer1D(dst_grid_dims(2)), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (dsthasbound) then
allocate(cornerlon2D(dst_grid_corner,dst_grid_dims(1)), &
cornerlat2D(dst_grid_corner,dst_grid_dims(2)), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridspecGetVar1D(dstFile, dst_coordids, lonBuffer1D, latBuffer1D,&
cornerlon=cornerlon2D, cornerlat=cornerlat2D, rc=status)
else
call ESMF_GridspecGetVar1D(dstFile, dst_coordids, lonBuffer1D, latBuffer1D,&
rc=status)
endif
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
k=1
do j=1,dst_grid_dims(2)
do i=1,dst_grid_dims(1)
lonBuffer(k)=lonBuffer1D(i)
latBuffer(k)=latBuffer1D(j)
k=k+1
enddo
enddo
deallocate(lonBuffer1D, latBuffer1D, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (dsthasbound) then
allocate(lonBuffer2(dst_grid_corner, dstDim),latBuffer2(dst_grid_corner,dstDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
k=1
do j=1,dst_grid_dims(2)
do i=1,dst_grid_dims(1)
lonBuffer2(:,k)=cornerlon2D(:,i)
latBuffer2(:,k)=cornerlat2D(:,j)
k=k+1
enddo
enddo
deallocate(cornerlon2D, cornerlat2D, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
else
allocate(lonBuffer2(dst_grid_dims(1),dst_grid_dims(2)), &
latBuffer2(dst_grid_dims(1),dst_grid_dims(2)), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (dsthasbound) then
allocate(cornerlon3D(dst_grid_corner,dst_grid_dims(1), dst_grid_dims(2)),&
cornerlat3D(dst_grid_corner,dst_grid_dims(1), dst_grid_dims(2)), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridspecGetVar2D(dstFile, dst_coordids, lonBuffer2, latBuffer2, &
cornerlon=cornerlon3D, cornerlat=cornerlat3D, rc=status)
else
call ESMF_GridspecGetVar2D(dstFile, dst_coordids, lonBuffer2, latBuffer2, &
rc=status)
endif
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
lonBuffer = reshape(lonBuffer2, (/dstDim/))
latBuffer = reshape(latBuffer2, (/dstDim/))
deallocate(lonBuffer2, latBuffer2, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (dsthasbound) then
allocate(lonBuffer2(dst_grid_corner, dstDim),latBuffer2(dst_grid_corner,dstDim), &
stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
lonBuffer2=reshape(cornerlon3D, (/dst_grid_corner, dstDim/))
latBuffer2=reshape(cornerlat3D, (/dst_grid_corner, dstDim/))
deallocate(cornerlon3D, cornerlat3D, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
endif
ncStatus=nf90_inq_varid(ncid,"xc_b",VarId)
ncStatus=nf90_put_var(ncid,VarId, lonBuffer)
errmsg = "Variable xc_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus=nf90_inq_varid(ncid,"yc_b",VarId)
ncStatus=nf90_put_var(ncid,VarId, latBuffer)
errmsg = "Variable yc_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
deallocate(latBuffer, lonBuffer, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
! Write xv_b, yv_b
if (dsthasbound) then
ncStatus=nf90_inq_varid(ncid,"xv_b",VarId)
ncStatus=nf90_put_var(ncid,VarId, lonBuffer2)
errmsg = "Variable xv_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus=nf90_inq_varid(ncid,"yv_b",VarId)
ncStatus=nf90_put_var(ncid,VarId, latBuffer2)
errmsg = "Variable yv_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
deallocate(lonBuffer2, latBuffer2, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
! Mask
allocate(mask(dstDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
mask(:)=1
if (present(dstMissingValue)) then
if (dstMissingValue) then
allocate(varBuffer(dst_grid_dims(1),dst_grid_dims(2)), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridspecGetVarByName(dstFile, dstvarname, dst_dimids, &
varBuffer, missing_value = missing_value, &
rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
k=1
do i=1,size(varBuffer,2)
do j=1,size(varBuffer,1)
if (varBuffer(j,i) == missing_value) mask(k)=0
k=k+1
enddo
enddo
deallocate(varBuffer, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
endif
ncStatus=nf90_inq_varid(ncid,"mask_b",VarId)
ncStatus=nf90_put_var(ncid,VarId, mask)
errmsg = "Variable mask_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
deallocate(mask, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
else if (dstFileTypeLocal == ESMF_FILEFORMAT_TILE) then
allocate(lonBuffer2(dst_grid_dims(1),dst_grid_dims(2)), &
latBuffer2(dst_grid_dims(1), dst_grid_dims(2)), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
allocate(cornerlon2D(dst_grid_dims(1)+1,dst_grid_dims(2)+1), &
cornerlat2D(dst_grid_dims(1)+1, dst_grid_dims(2)+1), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridspecReadTile(dstFile, dst_grid_dims(1), dst_grid_dims(2), &
lonBuffer2, latBuffer2, cornerLon=cornerlon2D, cornerLat=cornerlat2D, &
rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
allocate(lonBuffer(dstDim), latBuffer(dstDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
lonBuffer = reshape(lonBuffer2, (/dstDim/))
latBuffer = reshape(latBuffer2, (/dstDim/))
deallocate(lonBuffer2, latBuffer2, stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
allocate(lonBuffer2(dst_grid_corner, dstDim),latBuffer2(dst_grid_corner,dstDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
k=1
do j=1,dst_grid_dims(2)
do i=1,dst_grid_dims(1)
lonBuffer2(1,k)=cornerlon2D(i,j)
lonBuffer2(2,k)=cornerlon2D(i,j+1)
lonBuffer2(3,k)=cornerlon2D(i+1,j+1)
lonBuffer2(4,k)=cornerlon2D(i+1,j)
latBuffer2(1,k)=cornerlat2D(i,j)
latBuffer2(2,k)=cornerlat2D(i,j+1)
latBuffer2(3,k)=cornerlat2D(i+1,j+1)
latBuffer2(4,k)=cornerlat2D(i+1,j)
k=k+1
enddo
enddo
deallocate(cornerlon2D, cornerlat2D, stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
ncStatus=nf90_inq_varid(ncid,"xc_b",VarId)
ncStatus=nf90_put_var(ncid,VarId, lonBuffer)
errmsg = "Variable xc_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus=nf90_inq_varid(ncid,"yc_b",VarId)
ncStatus=nf90_put_var(ncid,VarId, latBuffer)
errmsg = "Variable yc_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
deallocate(latBuffer, lonBuffer, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
! Write xv_b, yv_b
ncStatus=nf90_inq_varid(ncid,"xv_b",VarId)
ncStatus=nf90_put_var(ncid,VarId, lonBuffer2)
errmsg = "Variable xv_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus=nf90_inq_varid(ncid,"yv_b",VarId)
ncStatus=nf90_put_var(ncid,VarId, latBuffer2)
errmsg = "Variable yv_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
deallocate(lonBuffer2, latBuffer2, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
else if (dstFileTypeLocal == ESMF_FILEFORMAT_ESMFMESH) then
ncStatus=nf90_open(dstFile,NF90_NOWRITE,ncid1)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
! only write out xv_a and yv_a when the regrid method is conserve
if (useDstCornerlocal .and. (methodlocal%regridmethod ==ESMF_REGRIDMETHOD_BILINEAR%regridmethod &
.or. methodlocal%regridmethod ==ESMF_REGRIDMETHOD_PATCH%regridmethod &
.or. methodlocal%regridmethod ==ESMF_REGRIDMETHOD_NEAREST_STOD%regridmethod &
.or. methodlocal%regridmethod ==ESMF_REGRIDMETHOD_NEAREST_DTOS%regridmethod)) then
! check if centerCoords exit
ncStatus=nf90_inq_varid(ncid1,"nodeCoords",VarId)
varStr = "nodeCoords"
dst_has_area = .false.
else
ncStatus=nf90_inq_varid(ncid1,"centerCoords",VarId)
varStr = "centerCoords"
endif
if (ncStatus /= nf90_noerror) then
write(*,*) "Warning: "//trim(varStr)// &
" not present in dst grid file, so not outputting xc_b and yc_b to weight file."
write(*,*)
else
allocate(latBuffer2(dstCoordDim, dstDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
allocate(latBuffer(dstDim), lonBuffer(dstDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
ncStatus=nf90_get_var(ncid1,VarId, latBuffer2)
errmsg = "Variable "//varStr//" in "//trim(dstFile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
do i=1,dstDim
lonBuffer(i)=latBuffer2(1,i)
latBuffer(i)=latBuffer2(2,i)
enddo
ncStatus=nf90_inq_varid(ncid,"xc_b",VarId)
ncStatus=nf90_put_var(ncid,VarId, lonBuffer)
errmsg = "Variable xc_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
ncStatus=nf90_inq_varid(ncid,"yc_b",VarId)
ncStatus=nf90_put_var(ncid,VarId, latBuffer)
errmsg = "Variable yc_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
deallocate(latBuffer, lonBuffer, latBuffer2, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
! output xv_b and yv_b is harder, we have to read in the nodeCoords and
! elementConn and construct the the latitudes and longitudes for
! all the corner vertices
if (.not. useDstCornerlocal .or. &
methodlocal%regridmethod ==ESMF_REGRIDMETHOD_CONSERVE%regridmethod .or. &
methodlocal%regridmethod ==ESMF_REGRIDMETHOD_CONSERVE_2ND%regridmethod) then
if (dst_grid_corner > 0) then
allocate(latBuffer2(dst_grid_corner,dstDim),lonBuffer2(dst_grid_corner,dstDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_EsmfGetVerts(ncid1, dstFile, dstDim, dst_grid_corner, dstNodeDim, &
latBuffer2, lonBuffer2, status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
ncStatus=nf90_inq_varid(ncid,"xv_b",VarId)
ncStatus=nf90_put_var(ncid,VarId, lonBuffer2)
errmsg = "Variable xv_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
ncStatus=nf90_inq_varid(ncid,"yv_b",VarId)
ncStatus=nf90_put_var(ncid,VarId, latBuffer2)
errmsg = "Variable yv_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
deallocate(latBuffer2, lonBuffer2, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
endif
! Write mask_b
allocate(mask(dstDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (.not. useDstCornerlocal .or. &
methodlocal%regridmethod ==ESMF_REGRIDMETHOD_CONSERVE%regridmethod .or. &
methodlocal%regridmethod ==ESMF_REGRIDMETHOD_CONSERVE_2ND%regridmethod) then
ncStatus=nf90_inq_varid(ncid1,"elementMask",VarId)
if (ncStatus /= nf90_noerror) then
write(*,*) "Warning: elementMask"// &
" not present in dst grid file, so setting mask_b=1 in weight file."
write(*,*)
mask = 1
else
ncStatus=nf90_get_var(ncid1,VarId, mask)
errmsg = "Variable elementMask in "//trim(dstFile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
end if
else
! make mask all 1
mask=1
endif
ncStatus=nf90_inq_varid(ncid,"mask_b",VarId)
ncStatus=nf90_put_var(ncid,VarId, mask)
errmsg = "Variable mask_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
deallocate(mask, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (dst_has_area .and. .not. present(dstArea)) then
allocate(area(dstDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
ncStatus=nf90_inq_varid(ncid1,"elementArea",VarId)
errmsg = "Variable elementArea in "//trim(dstFile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
ncStatus=nf90_get_var(ncid1,VarId, area)
errmsg = "Variable elementArea in "//trim(dstFile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
ncStatus=nf90_inq_varid(ncid,"area_b",VarId)
ncStatus=nf90_put_var(ncid,VarId, area)
errmsg = "Variable mask_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
deallocate(area, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
ncStatus=nf90_close(ncid1)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,trim(dstFile),&
rc)) return
else if (dstFileTypeLocal == ESMF_FILEFORMAT_UGRID) then
call ESMF_UGridInq(dstfile, dstmeshname, meshId=meshId, faceCoordFlag=faceCoordFlag)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (.not. useDstCornerlocal .or. &
methodlocal%regridmethod ==ESMF_REGRIDMETHOD_CONSERVE%regridmethod .or. &
methodlocal%regridmethod ==ESMF_REGRIDMETHOD_CONSERVE_2ND%regridmethod) then
! check if faceCoords exit
allocate(latBuffer2(dst_grid_corner,dstDim),&
lonBuffer2(dst_grid_corner,dstDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (faceCoordFlag) then
allocate(latBuffer(dstDim), lonBuffer(dstDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_UGridGetVar(dstfile, meshId, &
faceXcoords=lonBuffer, faceYcoords=latBuffer, &
faceNodeConnX=lonBuffer2, faceNodeConnY=latBuffer2, rc=status)
else
write(*,*) "Warning: face coordinates not present in dst grid file,"// &
" so not outputting xc_a and yc_a to weight file."
write(*,*)
call ESMF_UGridGetVar(dstfile, meshId, &
faceNodeConnX=lonBuffer2, faceNodeConnY=latBuffer2, rc=status)
endif
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (faceCoordFlag) then
ncStatus=nf90_inq_varid(ncid,"xc_b",VarId)
ncStatus=nf90_put_var(ncid,VarId, lonBuffer)
errmsg = "Variable xc_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
ncStatus=nf90_inq_varid(ncid,"yc_b",VarId)
ncStatus=nf90_put_var(ncid,VarId, latBuffer)
errmsg = "Variable yc_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
deallocate(latBuffer, lonBuffer, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
ncStatus=nf90_inq_varid(ncid,"xv_b",VarId)
ncStatus=nf90_put_var(ncid,VarId, lonBuffer2)
errmsg = "Variable xv_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
ncStatus=nf90_inq_varid(ncid,"yv_b",VarId)
ncStatus=nf90_put_var(ncid,VarId, latBuffer2)
errmsg = "Variable yv_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
deallocate(latBuffer2, lonBuffer2, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
else
allocate(latBuffer(dstDim), lonBuffer(dstDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_UGridGetVar(dstfile, meshId, &
nodeXcoords=lonBuffer, nodeYcoords=latBuffer, rc=rc)
ncStatus=nf90_inq_varid(ncid,"xc_b",VarId)
ncStatus=nf90_put_var(ncid,VarId, lonBuffer)
errmsg = "Variable xc_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
ncStatus=nf90_inq_varid(ncid,"yc_b",VarId)
ncStatus=nf90_put_var(ncid,VarId, latBuffer)
errmsg = "Variable yc_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
deallocate(latBuffer, lonBuffer, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
! Write out mask
allocate(mask(dstDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
mask(:)=1
if (present(dstMissingValue)) then
if (dstMissingValue) then
allocate(varBuffer1D(dstDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_UgridGetVarByName(dstFile, dstvarname, &
varBuffer1D, missingvalue = missing_value, &
rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
do j=1,size(varBuffer1D)
if (varBuffer1D(j) == missing_value) mask(j)=0
enddo
deallocate(varBuffer1D, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
endif
ncStatus=nf90_inq_varid(ncid,"mask_b",VarId)
ncStatus=nf90_put_var(ncid,VarId, mask)
errmsg = "Variable mask_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
deallocate(mask, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
else if (dstFileTypeLocal == ESMF_FILEFORMAT_MOSAIC) then
!read coordinates from the tile files
allocate(latBuffer2(dstmosaic%nx, dstmosaic%ny), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
allocate(lonBuffer2(dstmosaic%nx, dstmosaic%ny), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
totalsize = dstmosaic%nx*dstmosaic%ny
allocate(varBuffer1D(totalsize), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
ncStatus=nf90_inq_varid(ncid,"xc_b",VarId1)
ncStatus=nf90_inq_varid(ncid,"yc_b",VarId2)
errmsg = "Variable xc_b or yc_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
do k=1,dstmosaic%ntiles
totallen = len_trim(dstmosaic%filenames(k))+len_trim(dstmosaic%tileDirectory)
tempname = trim(dstmosaic%tileDirectory)//trim(dstmosaic%filenames(k))
call ESMF_GridspecReadTile(trim(tempname),dstmosaic%nx, dstmosaic%ny, &
centerlon=lonBuffer2, centerlat=latBuffer2, rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
varBuffer1D=reshape(lonBuffer2,(/totalsize/))
start=(k-1)*totalsize+1
ncStatus = nf90_put_var(ncid, VarId1, varBuffer1D, (/start/), (/totalsize/))
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
"writing xc_b in weight file", &
rc)) return
varBuffer1D=reshape(latBuffer2,(/totalsize/))
ncStatus = nf90_put_var(ncid, VarId2,varBuffer1D, (/start/), (/totalsize/))
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
"writing yc_b in weight file", &
rc)) return
enddo
deallocate(lonBuffer2, latBuffer2, varBuffer1D, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
! Write area_a
ncStatus=nf90_inq_varid(ncid,"area_a",VarId)
if (present(srcArea)) then
ncStatus=nf90_put_var(ncid,VarId, srcArea)
errmsg = "Variable area_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
else if (.not. src_has_area) then
! Just set these to 0.0, because not provided
allocate(area(srcDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
area=0.0
ncStatus=nf90_put_var(ncid,VarId, area)
errmsg = "Variable area_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
deallocate(area, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
! Write area_b
ncStatus=nf90_inq_varid(ncid,"area_b",VarId)
if (present(dstArea)) then
ncStatus=nf90_put_var(ncid,VarId, dstArea)
errmsg = "Variable area_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
else if (.not. dst_has_area) then
! Just set these to 0.0, because not provided
allocate(area(dstDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
area=0.0
ncStatus=nf90_put_var(ncid,VarId, area)
errmsg = "Variable area_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
deallocate(area, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
! Write frac_a
ncStatus=nf90_inq_varid(ncid,"frac_a",VarId)
if (present(srcFrac)) then
ncStatus=nf90_put_var(ncid,VarId, srcFrac)
else
allocate(frac(srcDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
frac=0.0
ncStatus=nf90_put_var(ncid,VarId, frac)
deallocate(frac, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
errmsg = "Variable frac_a in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
! Write frac_b
ncStatus=nf90_inq_varid(ncid,"frac_b",VarId)
if (present(dstFrac)) then
ncStatus=nf90_put_var(ncid,VarId, dstFrac)
else
allocate(frac(dstDim), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
frac=1.0
ncStatus=nf90_put_var(ncid,VarId, frac)
deallocate(frac, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif
errmsg = "Variable frac_b in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
deallocate(src_grid_dims, dst_grid_dims, stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
endif ! PetNo==0
! Block all other PETs until the NetCDF file has been created
call ESMF_VMBarrier(vm)
! find the max of allCounts(i) and allocate colrow
maxcount=0
do i=1,PetCnt
if (allCounts(i) > maxcount) maxcount = allCounts(i)
enddo
if (PetNo == 0) then
! First write out its own weight and indices, then receive the data from
! other PETs and write them out
start = 1
! allocate indexbuf and weightbuf to receive data from other PETs
allocate(indexbuf(maxcount*2), weightbuf(maxcount), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
do i=1, PetCnt
! write the local factorList and factorIndexList first
localCount(1)=allCounts(i)
if (i==1) then
!do j=1,localCount(1)
! indexbuf(j) = factorIndexList(1,j)
!enddo
next => factorIndexList(1,:)
ncStatus=nf90_inq_varid(ncid,"col",VarId)
ncStatus=nf90_put_var(ncid,VarId, next,(/start/),localCount)
errmsg = "Variable col in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
!do j=1,localCount(1)
! indexbuf(j) = factorIndexList(2,j)
!enddo
next => factorIndexList(2,:)
ncStatus=nf90_inq_varid(ncid,"row",VarId)
ncStatus=nf90_put_var(ncid,VarId, next,(/start/),localCount)
errmsg = "Variable row in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
ncStatus=nf90_inq_varid(ncid,"S",VarId)
ncStatus=nf90_put_var(ncid,VarId, factorList, (/start/),localCount)
errmsg = "Variable S in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
else
if (localCount(1) > 0) then
! receive the factorList and factorIndexList
call ESMF_VMRecv(vm, indexbuf, localCount(1)*2, i-1, rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_VMRecv(vm, weightbuf, localCount(1), i-1, rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
ncStatus=nf90_inq_varid(ncid,"col",VarId)
ncStatus=nf90_put_var(ncid,VarId, indexbuf,(/start/),localCount)
errmsg = "Variable col in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
next => indexbuf(localCount(1)+1:localCount(1)*2)
ncStatus=nf90_inq_varid(ncid,"row",VarId)
ncStatus=nf90_put_var(ncid,VarId, next ,(/start/),localCount)
errmsg = "Variable row in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
ncStatus=nf90_inq_varid(ncid,"S",VarId)
ncStatus=nf90_put_var(ncid,VarId, weightbuf, (/start/),localCount)
errmsg = "Variable S in "//trim(wgtfile)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,errmsg,&
rc)) return
end if
end if
start = start + localCount(1)
end do
else
allocate(indexbuf(localcount(1)*2), stat=memstat)
if (ESMF_LogFoundAllocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (localcount(1) > 0) then
do j=1,localCount(1)
indexbuf(j) = factorIndexList(1,j)
indexbuf(j+localCount(1)) = factorIndexList(2,j)
enddo
! a non-root PET, send the results to PET 0
call ESMF_VMSend(vm, indexbuf, localCount(1)*2, 0, rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_VMSend(vm, factorList, localCount(1), 0, rc=status)
if (ESMF_LogFoundError(status, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
end if
end if
deallocate(allCounts, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (PetNo == 0) then
ncStatus = nf90_close(ncid)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE, trim(wgtfile),&
rc)) return
deallocate(weightbuf, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
end if
call ESMF_VMBarrier(vm)
deallocate(indexbuf, stat=memstat)
if (ESMF_LogFoundDeallocError(memstat, &
ESMF_CONTEXT, rcToReturn=rc)) return
if (present(rc)) rc = ESMF_SUCCESS
return
#else
if (ESMF_LogFoundError(ESMF_RC_LIB_NOT_PRESENT, &
msg="- ESMF_NETCDF not defined when lib was compiled", &
ESMF_CONTEXT, rcToReturn=rc)) return
#endif
end subroutine ESMF_OutputScripWeightFile