subroutine MyGridCreate_ ( vm, grid2d, grid3d, rc)
type (ESMF_VM), intent(IN ) :: VM
type(ESMF_GRID), intent(INOUT) :: grid2d, grid3d
integer, optional, intent(OUT) :: rc
! Local vars
integer :: status
character(len=ESMF_MAXSTR), parameter :: IAm='MyGridCreate'
integer :: LM
integer :: L
integer :: NX, NY
integer, allocatable :: IMXY(:), JMXY(:)
character(len=ESMF_MAXSTR) :: gridname
real(ESMF_KIND_R8) :: minCoord(3)
real(ESMF_KIND_R8) :: deltaX, deltaY, deltaZ
real :: LON0, LAT0
real :: pi, d2r
! grid create
lm = KM_WORLD ! no. vertical layers
nx = 2
ny = 2
pi = 4.0 * atan ( 1.0 )
d2r = pi / 180.
LON0 = -180 * d2r
LAT0 = -90.0 * d2r
! Get the IMXY vector
! -------------------
allocate( imxy(0:nx-1) )
call MAPL_GET_LOCAL_DIMS ( IM_WORLD, imxy, nx )
! Get the JMXY vector
! -------------------
allocate( jmxy(0:ny-1) )
call MAPL_GET_LOCAL_DIMS ( JM_WORLD, jmxy, ny )
deltaX = 2.0*pi/IM_WORLD
deltaY = pi/(JM_WORLD-1)
deltaZ = 1.0
! Define South-West Corner of First Grid-Box
! ------------------------------------------
minCoord(1) = LON0 - deltaX/2
minCoord(2) = LAT0 - deltaY/2
minCoord(3) = deltaZ/2.
layout = ESMF_DELayoutCreate(vm, deCountList=(/NX, NY/), rc=status)
_VERIFY(STATUS)
if ( MAPL_Am_I_Root() ) then
print *
print *, ' Testing CFIO Read Array'
print *, ' -----------------------'
print *
print *, 'Decomposition Information:'
print *, ' X AXIS: nx, imxy = ', nx, imxy
print *, ' Y AXIS: ny, jmxy = ', ny, jmxy
print *, ' Z AXIS: nz = ', lm, ' (not decomposed)'
end if
! 2D grid
! -------
grid2d = ESMF_GridCreateHorzLatLonUni( &
counts = (/IM_WORLD, JM_WORLD/), &
minGlobalCoordPerDim=minCoord(1:2), &
deltaPerDim=(/deltaX, deltaY /), &
horzStagger=ESMF_Grid_Horz_Stagger_A, &
periodic=(/ESMF_TRUE, ESMF_FALSE/), &
name='Grid2d', rc=status)
_VERIFY(STATUS)
call ESMF_GridDistribute(grid2d, &
deLayout=layout, &
countsPerDEDim1=imxy, &
countsPerDEDim2=jmxy, &
rc=status)
_VERIFY(STATUS)
! 3D Grid
! -------
grid3d = ESMF_GridCreateHorzLatLonUni( &
counts = (/IM_WORLD, JM_WORLD/), &
minGlobalCoordPerDim=minCoord(1:2), &
deltaPerDim=(/deltaX, deltaY /), &
horzStagger=ESMF_Grid_Horz_Stagger_A, &
periodic=(/ESMF_TRUE, ESMF_FALSE/), &
name='Grid3d', rc=status)
_VERIFY(STATUS)
call ESMF_GridAddVertHeight(grid3d, &
delta=(/(deltaZ, L=1,LM) /), &
rc=status)
_VERIFY(STATUS)
call ESMF_GridDistribute(grid3d, &
deLayout=layout, &
countsPerDEDim1=imxy, &
countsPerDEDim2=jmxy, &
rc=status)
_VERIFY(STATUS)
call MAPL_GridGetLatLons ( grid3d, lons, lats )
if ( MAPL_Am_I_Root() ) then
print *
print *, 'Grid points:'
print *, ' - Longitudes: '
write(*,100) lons(:)
print *, ' - Latitudes:'
write(*,100) lats(:)
print *
100 format(( ' ',8F8.2))
endif
deallocate(imxy)
deallocate(jmxy)
_RETURN(STATUS)
end subroutine MyGridCreate_