subroutine ESMF_OutputScripGridFile(filename, grid, rc)
character(len=*), intent(in) :: filename
type(ESMF_GRID), intent(in) :: grid
integer :: rc
integer :: ncid, dimid, varid, nodedimid
integer :: localrc, ncStatus
integer :: varsize
integer :: ndims, dims(ESMF_MAXDIM)
integer :: PetNo, PetCnt
type(ESMF_STAGGERLOC), allocatable :: Staggers(:)
integer :: staggercnt, decount, xdim, ydim
integer :: londim, londim1, latdim, latdim1, gridid
integer :: rankid, fourid, varid1, varid2, varid3, varid4
integer :: maskid, areaid
integer :: elmtsize, count, i, j, n, nextj
integer, pointer :: minindex(:), maxindex(:)
logical :: xperiod, yperiod
type(ESMF_VM) :: vm
type(ESMF_Grid) :: grid1
type(ESMF_DistGrid) :: distgrid
type(ESMF_Array) :: array
type(ESMF_Array) :: lonarray, latarray
type(ESMF_CoordSys_Flag) :: coordsys
real(ESMF_KIND_R8), pointer::lonArray1d(:), latArray1d(:)
real(ESMF_KIND_R8), pointer::lonArray2d(:,:), latArray2d(:,:)
real(ESMF_KIND_R8), pointer::scripArray(:)
real(ESMF_KIND_R8), pointer::scripArray2(:,:)
integer(ESMF_KIND_I4), pointer :: fptrMask(:,:), scripArrayMask(:)
logical :: hasmask, hasarea
character (len=256) :: errmsg, units
#ifdef ESMF_NETCDF
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
! find out the grid size and other global information
call ESMF_GridGet(grid, dimCount=ndims, localDECount=decount, &
coordDimCount=dims, coordSys=coordsys, distgrid=distgrid, &
staggerlocCount=staggercnt, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! Get list and number of active staggers
allocate(Staggers(staggercnt))
call c_ESMC_gridgetactivestaggers(grid%this, &
staggercnt, Staggers, localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
deallocate(Staggers)
if (ndims /= 2) then
call ESMF_LogSetError(rcToCheck=ESMF_FAILURE, &
msg="- Grid dimension has to be 2 to write out to a SCRIP file", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
if (decount > 1) then
call ESMF_LogSetError(rcToCheck=ESMF_FAILURE, &
msg="- Only one localDE per PET is allowed to write a grid into a SCRIP file", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
allocate(minindex(ndims), maxindex(ndims))
call ESMF_GridGet(grid, tile=1, staggerloc=ESMF_STAGGERLOC_CENTER, minIndex=minindex, &
maxIndex=maxindex,rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! find the size of the grid
londim = maxindex(1)-minindex(1)+1
latdim = maxindex(2)-minindex(2)+1
elmtsize=londim * latdim
if (coordsys == ESMF_COORDSYS_SPH_DEG) then
units="degrees"
elseif (coordsys == ESMF_COORDSYS_SPH_RAD) then
units="radians"
else
!ESMF_COORDSYS_CART not supported -- errors
call ESMF_LogSetError(rcToCheck=ESMF_FAILURE, &
msg="- ESMF_ScripOuputGridFile does not support ESMF_COORDSYS_CART coordinates", &
ESMF_CONTEXT, rcToReturn=rc)
return
endif
if (PetNo==0) then
! Create the GRID file and define dimensions and variables
ncStatus=nf90_create(filename, NF90_CLOBBER, ncid)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
trim(filename),&
rc)) return
ncStatus=nf90_def_dim(ncid,"grid_size",elmtsize, nodedimid)
errmsg = "defining dimension grid_size in "//trim(filename)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus=nf90_def_dim(ncid,"grid_rank", 2, rankid)
errmsg = "defining dimension grid_rank in "//trim(filename)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
if (staggercnt > 1) then
ncStatus=nf90_def_dim(ncid,"grid_corners", 4, fourid)
errmsg = "defining dimension grid_corner in "//trim(filename)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
endif
ncStatus=nf90_def_var(ncid,"grid_dims",NF90_INT, (/rankid/), gridid)
errmsg = "defining variable grid_dims in "//trim(filename)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus=nf90_def_var(ncid,"grid_center_lon",NF90_DOUBLE, (/nodedimid/), varid1)
errmsg = "defining variable grid_center_lon in "//trim(filename)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
! add units attribute based on the COORDSYS
ncStatus=nf90_put_att(ncid,varid1,"units",trim(units))
errmsg = "Adding attribute units for grid_center_lon in "//trim(filename)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus=nf90_def_var(ncid,"grid_center_lat",NF90_DOUBLE, (/nodedimid/), varid2)
errmsg = "defining variable grid_center_lat in "//trim(filename)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus=nf90_put_att(ncid,varid2,"units",trim(units))
errmsg = "Adding attribute units for grid_center_lat in "//trim(filename)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
if (staggercnt > 1) then
ncStatus=nf90_def_var(ncid,"grid_corner_lon",NF90_DOUBLE, (/fourid,nodedimid/), varid3)
errmsg = "defining variable grid_corner_lon in "//trim(filename)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
! add units attribute based on the COORDSYS
ncStatus=nf90_put_att(ncid,varid3,"units",trim(units))
errmsg = "Adding attribute units for grid_corner_lon in "//trim(filename)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus=nf90_def_var(ncid,"grid_corner_lat",NF90_DOUBLE, (/fourid, nodedimid/), varid4)
errmsg = "defining variable grid_corner_lat in "//trim(filename)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
ncStatus=nf90_put_att(ncid,varid4,"units",trim(units))
errmsg = "Adding attribute units for grid_corner_lat in "//trim(filename)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
endif
hasmask = .false.
hasarea = .false.
! If the grid has mask and area, output them too
call ESMF_GridGetItem(grid, ESMF_GRIDITEM_MASK, staggerloc=ESMF_STAGGERLOC_CENTER, &
array=array, rc=localrc)
if (localrc == ESMF_SUCCESS) then
ncStatus=nf90_def_var(ncid, "grid_imask", NF90_INT, (/nodedimid/), maskid)
errmsg = "defining variable grid_imask in "//trim(filename)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
hasmask = .true.
endif
call ESMF_GridGetItem(grid, ESMF_GRIDITEM_AREA, staggerloc=ESMF_STAGGERLOC_CENTER, &
array=array, rc=localrc)
if (localrc == ESMF_SUCCESS) then
ncStatus=nf90_def_var(ncid, "grid_area", NF90_DOUBLE, (/nodedimid/), areaid)
errmsg = "defining variable grid_area in "//trim(filename)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
hasarea = .true.
endif
ncStatus=nf90_enddef(ncid)
errmsg = "nf90_enddef in "//trim(filename)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
! write the grid_dims value
ncStatus = nf90_put_var(ncid,gridid,(/londim, latdim/))
errmsg = "writing grid_dims in "//trim(filename)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
endif
call ESMF_VMBarrier(vm)
! Create a grid at Pet 0 with the same size as the source grid
distgrid = ESMF_DistGridCreate(minIndex=(/1,1/), maxIndex=(/londim, latdim/), &
regDecomp=(/1,1/),rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
grid1 = ESMF_GridCreate(grid, distgrid, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! Get the coordinate data from grid1 at PET 0 to output to the file
if (PetNo == 0) then
call ESMF_GridGetCoord(grid1, 1, staggerloc=ESMF_STAGGERLOC_CENTER, &
array=lonarray, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridGetCoord(grid1, 2, staggerloc=ESMF_STAGGERLOC_CENTER, &
array=latarray, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! Extra the data from array and create 2D lat and lon array to write
! to the file
if (dims(1)==1) then
call ESMF_ArrayGet(lonarray,farrayPtr=lonArray1d, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_ArrayGet(latarray,farrayPtr=latArray1d, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
allocate(scripArray(elmtsize))
do i=1,elmtsize,londim
scripArray(i:i+londim)=lonArray1d
enddo
ncStatus=nf90_put_var(ncid, varid1, scripArray)
errmsg = "Writing variable grid_corner_lon in "//trim(filename)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
j=1
do i=1,latdim
scripArray(j:j+londim)=latArray1d(i)
j=j+londim
enddo
ncStatus=nf90_put_var(ncid, varid2, scripArray)
errmsg = "Writing variable grid_corner_lat in "//trim(filename)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
deallocate(scripArray)
else
call ESMF_ArrayGet(lonarray,farrayPtr=lonArray2d, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_ArrayGet(latarray,farrayPtr=latArray2d, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
allocate(scripArray(elmtsize))
scripArray=reshape(lonArray2d,(/elmtsize/))
ncStatus=nf90_put_var(ncid, varid1, scripArray)
errmsg = "Writing variable grid_corner_lon in "//trim(filename)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
scripArray=reshape(latArray2d,(/elmtsize/))
ncStatus=nf90_put_var(ncid, varid2, scripArray)
errmsg = "Writing variable grid_corner_lat in "//trim(filename)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
deallocate(scripArray)
endif
!call ESMF_ArrayDestroy(lonArray)
!call ESMF_ArrayDestroy(latArray)
! Get corner stagger if exist
if (staggercnt > 1) then
call ESMF_GridGetCoord(grid1, 1, staggerloc=ESMF_STAGGERLOC_CORNER, &
array=lonarray, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_GridGetCoord(grid1, 2, staggerloc=ESMF_STAGGERLOC_CORNER, &
array=latarray, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
! Find the corner stagger array dimension, it may have one periodic dimension
call ESMF_GridGet(grid, tile=1, staggerloc=ESMF_STAGGERLOC_CORNER,&
minIndex=minIndex, maxIndex=maxIndex, rc=localrc)
londim1 = maxIndex(1)-minIndex(1)+1
latdim1 = maxIndex(2)-minIndex(2)+1
xperiod = .false.
yperiod = .false.
if (londim1 == londim) then
!This is periodic dimension, need to wrap around
xperiod = .true.
endif
if (latdim1 == latdim) then
!This is periodic dimension, need to wrap around
yperiod = .true.
endif
deallocate(minindex, maxindex)
!print *, 'Corner stagger dimension ', londim1, latdim1
! Extract the data from array and create 2D lat and lon array to write
! to the file
if (dims(1)==1) then
call ESMF_ArrayGet(lonarray,farrayPtr=lonArray1d, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_ArrayGet(latarray,farrayPtr=latArray1d, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
allocate(scripArray2(4,elmtsize))
count = 1
do j=1,latdim
do i=1,londim-1
scripArray2(1,count)=lonArray1d(i)
scripArray2(2,count)=lonArray1d(i+1)
scripArray2(3,count)=lonArray1d(i+1)
scripArray2(4,count)=lonArray1d(i)
count = count+1
enddo
if (xperiod) then
scripArray2(1,count)=lonArray1d(i)
scripArray2(2,count)=lonArray1d(1)
scripArray2(3,count)=lonArray1d(1)
scripArray2(4,count)=lonArray1d(i)
count = count+1
else
scripArray2(1,count)=lonArray1d(i)
scripArray2(2,count)=lonArray1d(i+1)
scripArray2(3,count)=lonArray1d(i+1)
scripArray2(4,count)=lonArray1d(i)
count = count+1
endif
enddo
ncStatus=nf90_put_var(ncid, varid3, scripArray2)
errmsg = "Writing variable grid_center_lon in "//trim(filename)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
count = 1
do j=1,latdim-1
do i=1,londim
scripArray2(1,count)=latArray1d(j)
scripArray2(2,count)=latArray1d(j)
scripArray2(3,count)=latArray1d(j+1)
scripArray2(4,count)=latArray1d(j+1)
count = count+1
enddo
enddo
if (yperiod) then
do i=1,londim
scripArray2(1,count)=latArray1d(j)
scripArray2(2,count)=latArray1d(j)
scripArray2(3,count)=latArray1d(1)
scripArray2(4,count)=latArray1d(1)
count = count+1
enddo
else
do i=1,londim
scripArray2(1,count)=latArray1d(j)
scripArray2(2,count)=latArray1d(j)
scripArray2(3,count)=latArray1d(j+1)
scripArray2(4,count)=latArray1d(j+1)
count = count+1
enddo
endif
ncStatus=nf90_put_var(ncid, varid4, scripArray2)
errmsg = "Writing variable grid_center_lat in "//trim(filename)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
deallocate(scripArray2)
else !dims==2
call ESMF_ArrayGet(lonarray,farrayPtr=lonArray2d, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
call ESMF_ArrayGet(latarray,farrayPtr=latArray2d, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
allocate(scripArray2(4,elmtsize))
count = 1
do j=1,latdim
if (j==latdim .and. yperiod) then
nextj=1
else
nextj=j+1
endif
do i=1,londim-1
scripArray2(1,count)=lonArray2d(i,j)
scripArray2(2,count)=lonArray2d(i+1,j)
scripArray2(3,count)=lonArray2d(i+1,nextj)
scripArray2(4,count)=lonArray2d(i,nextj)
count = count+1
enddo
if (xperiod) then
scripArray2(1,count)=lonArray2d(i,j)
scripArray2(2,count)=lonArray2d(1,j)
scripArray2(3,count)=lonArray2d(1,nextj)
scripArray2(4,count)=lonArray2d(i,nextj)
count = count+1
else
scripArray2(1,count)=lonArray2d(i,j)
scripArray2(2,count)=lonArray2d(i+1,j)
scripArray2(3,count)=lonArray2d(i+1,nextj)
scripArray2(4,count)=lonArray2d(i,nextj)
count = count+1
endif
enddo
ncStatus=nf90_put_var(ncid, varid3, scripArray2)
errmsg = "Writing variable grid_center_lon in "//trim(filename)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
count = 1
do j=1,latdim
if (j==latdim .and. yperiod) then
nextj=1
else
nextj=j+1
endif
do i=1,londim-1
scripArray2(1,count)=latArray2d(i,j)
scripArray2(2,count)=latArray2d(i+1,j)
scripArray2(3,count)=latArray2d(i+1,nextj)
scripArray2(4,count)=latArray2d(i,nextj)
count = count+1
enddo
if (xperiod) then
scripArray2(1,count)=latArray2d(i,j)
scripArray2(2,count)=latArray2d(1,j)
scripArray2(3,count)=latArray2d(1,nextj)
scripArray2(4,count)=latArray2d(i,nextj)
count = count+1
else
scripArray2(1,count)=latArray2d(i,j)
scripArray2(2,count)=latArray2d(i+1,j)
scripArray2(3,count)=latArray2d(i+1,nextj)
scripArray2(4,count)=latArray2d(i,nextj)
count = count+1
endif
enddo
ncStatus=nf90_put_var(ncid, varid4, scripArray2)
errmsg = "Writing variable grid_center_lat in "//trim(filename)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
deallocate(scripArray2)
endif !staggercnt > 1
if (hasmask) then
call ESMF_GridGetItem(grid1, ESMF_GRIDITEM_MASK, farrayPtr=fptrMask, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
allocate(scripArrayMask(elmtsize))
scripArrayMask=reshape(fptrMask,(/elmtsize/))
ncStatus=nf90_put_var(ncid, maskid, scripArrayMask)
errmsg = "Writing variable grid_imask in "//trim(filename)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
deallocate(scripArrayMask)
endif
if (hasarea) then
call ESMF_GridGetItem(grid1, ESMF_GRIDITEM_AREA, farrayPtr=scripArray2, rc=localrc)
if (ESMF_LogFoundError(localrc, ESMF_ERR_PASSTHRU, &
ESMF_CONTEXT, rcToReturn=rc)) return
allocate(scripArray(elmtsize))
scripArray=reshape(fptrMask,(/elmtsize/))
ncStatus=nf90_put_var(ncid, areaid, scripArray)
errmsg = "Writing variable grid_area in "//trim(filename)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
errmsg,&
rc)) return
deallocate(scripArray)
endif
endif !PetNo==0
!call ESMF_ArrayDestroy(lonArray)
!call ESMF_ArrayDestroy(latArray)
ncStatus = nf90_close(ncid)
if (CDFCheckError (ncStatus, &
ESMF_METHOD, &
ESMF_SRCLINE,&
trim(filename),&
rc)) return
endif !PET==0
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_OutputScripGridFile