ESMFL_RegridStore Subroutine

public subroutine ESMFL_RegridStore(srcFLD, SRCgrid2D, dstFLD, DSTgrid2D, vm, rh, rc)

Given a srcFLD and its associated 3dGrid and a dstFLD and its associated 3DGrid, the subroutine ESMFL_RegridStore creates their corresponding 2DGrids and a 2D routehandle.

History

  • 17Oct2005 Cruz Initial code.

Arguments

Type IntentOptional Attributes Name
type(ESMF_Field), intent(inout) :: srcFLD
type(ESMF_Grid), intent(out) :: SRCgrid2D
type(ESMF_Field), intent(inout) :: dstFLD
type(ESMF_Grid), intent(out) :: DSTgrid2D
type(ESMF_VM), intent(in) :: vm
type(ESMF_RouteHandle), intent(inout) :: rh
integer, intent(out), optional :: rc

Source Code

   subroutine ESMFL_RegridStore (srcFLD, SRCgrid2D, dstFLD, DSTgrid2D, &
                                 vm, rh, rc)
!
   implicit NONE

! !ARGUMENTS:

   type(ESMF_Field), intent(inout)       :: srcFLD
   type(ESMF_Field), intent(inout)       :: dstFLD
   type(ESMF_Grid), intent(out)          :: SRCgrid2D
   type(ESMF_Grid), intent(out)          :: DSTgrid2D
   type(ESMF_RouteHandle), intent(inout) :: rh
   type(ESMF_VM), intent(in)             :: vm  ! should be intent IN
   integer, optional, intent(OUT)        :: rc
!
!-------------------------------------------------------------------------

#if 0
! local vars

   type(ESMF_DELayout)     :: layout
   type(ESMF_Field)        :: dstFld2D
   type(ESMF_Field)        :: srcFld2D
   type(ESMF_ArraySpec)    :: arrayspec
   type(ESMF_Array)        :: dstARR
   type(ESMF_Array)        :: srcARR
   type(ESMF_FieldDataMap) :: dmap
   type(ESMF_Grid)         :: grid3D
   real(kind=REAL32), pointer        :: Sptr2d(:,:)
   real(kind=REAL32), pointer        :: Dptr2d(:,:)
   real(ESMF_KIND_R8)      :: deltaX, deltaY
   real(ESMF_KIND_R8)      :: min(2), max(2)
   integer :: status, rank
   integer :: sCPD(3), dCPD(3) ! localCellCountPerDim
   integer :: gccpd(3)         ! globalCellCountPerDim

   type(ESMF_AxisIndex), dimension(:,:), pointer :: AI
   integer, allocatable, dimension(:)            :: ims, jms
   integer                                   :: nDEs, de
   integer                                   :: NX, NY
   integer                                   :: IX, JY
   integer                                   :: gridRank
   integer                                   :: decpd(7)
#endif

! start

   _UNUSED_DUMMY(srcFLD)
   _UNUSED_DUMMY(SRCgrid2D)
   _UNUSED_DUMMY(dstFLD)
   _UNUSED_DUMMY(DSTgrid2d)
   _UNUSED_DUMMY(rh)
   _UNUSED_DUMMY(vm)

#if 0
   ! can get the rank from either srcFLD or dstFLD
   call ESMF_FieldGetArray(srcFLD, srcARR, rc=status)
   _VERIFY(status)
   call ESMF_ArrayGet(srcARR, RANK=rank, rc=status)
   _VERIFY(status)
   ! datamap - with rank=2!
   call ESMF_FieldDataMapSetDefault(dmap, 2, rc=status)
   _VERIFY(STATUS)

   ! Create 2D grids and fields

   ! SOURCE FIELD

   call ESMF_FieldGet(srcFLD, grid=grid3D, rc=status)
   _VERIFY(status)

   call ESMF_GridGet(GRID3d, dimCount=gridRank, rc=STATUS)
   _VERIFY(STATUS)

   call ESMF_GridGetDELayout(grid3D, layout, rc=status)
   _VERIFY(STATUS)

   call ESMF_DELayoutGet(layout, deCount=nDEs, deCountPerDim = decpd, rc=status)

   NX = decpd(1)
   NY = decpd(2)

   allocate (ims(NX), jms(NY), stat=status)
   _VERIFY(STATUS)

   call ESMF_GridGet(grid3D, &
                     horzRelLoc=ESMF_CELL_CENTER, vertRelLoc=ESMF_CELL_CELL, &
                     globalCellCountPerDim=gccpd, &
                     minGlobalCoordPerDim=min, &
                     maxGlobalCoordPerDim=max, &
                     rc=status)
   _VERIFY(status)

#if 0
   ! kludge min,max DEGREE/RADIANS conversion
   if (min(1) < 0.0) then
      min(1) = min(1) * 180./MAPL_PI
      min(2) = min(2) * 180./MAPL_PI
      max(1) = max(1) * 180./MAPL_PI
      max(2) = max(2) * 180./MAPL_PI
   end if
   ! overide min,max - just create a simple 2-D grid
   min(1) = 0.0
   min(2) = 0.0
   max(1) = 360.0
   max(2) = 180.0
#endif

   ! this is probably incorrect unless parent grid is XYuni
   !SRCGrid2D = ESMF_GridCreateHorzXYUni( &
   !                           counts=gccpd(1:2), &
   !                           minGlobalCoordPerDim=min, &
   !                           maxGlobalCoordPerDim=max, &
   !                           horzStagger=ESMF_GRID_HORZ_STAGGER_A, &
   !                           periodic=(/ESMF_TRUE, ESMF_FALSE/),     &
   !                           name="SRC 2D grid", rc=status)
   !_VERIFY(status)
   ! instead use the following...
   deltaX = 2.0*MAPL_PI/gccpd(1)
   deltaY = MAPL_PI/(gccpd(2)-1)
   SRCGrid2D = ESMF_GridCreateHorzLatLonUni(         &
         counts = gccpd(1:2),        &
         minGlobalCoordPerDim=min(1:2),     &
         deltaPerDim=(/deltaX, deltaY /),        &
         horzStagger=ESMF_Grid_Horz_Stagger_A,   &
         periodic=(/ESMF_TRUE, ESMF_FALSE/),     &
         name='SRC 2D grid', rc=status)
   _VERIFY(STATUS)

   allocate (AI(nDEs,gridRank), stat=status)
   _VERIFY(STATUS)

   if (gridRank == 3) then
      call ESMF_GridGetAllAxisIndex(grid3d, &
                                    horzRelLoc=ESMF_CELL_CENTER, &
                                    vertRelLoc=ESMF_CELL_CELL, &
                                    globalAI=AI, rc=status)
   else
      call ESMF_GridGetAllAxisIndex(grid3d, &
                                    horzRelLoc=ESMF_CELL_CENTER, &
                                    globalAI=AI, rc=status)
   end if
   _VERIFY(STATUS)

   JY = 1
   DO IX = 1, NX

!ALT: this is a workaround to compute deId from Layout coords
      de = (jy-1)*NX + ix

      ims(IX) = AI(de,1)%max - AI(de,1)%min + 1
   END DO

   IX = 1
   DO JY = 1, NY

!ALT: same workaround as above
      de = (jy-1)*NX + ix

      jms(JY) = AI(de,2)%max - AI(de,2)%min + 1
   END DO

   deallocate(AI)

   if (verbose .and. MAPL_AM_I_Root()) then
      print *, 'ims=', ims
      print *, 'jms=', jms
   end if

   call ESMF_GridDistribute(SRCGrid2D, delayout=layout, &
                            countsPerDEDim1=ims,        &
                            countsPerDEDim2=jms,        &
                            rc=status)
   _VERIFY(status)

   deallocate(jms, ims)

   call ESMF_GridGetDELocalInfo(SRCGrid2D, &
         horzRelLoc=ESMF_CELL_CENTER, &
!         vertRelLoc=ESMF_CELL_CELL, &
         localCellCountPerDim=sCPD,RC=STATUS)
   _VERIFY(STATUS)
   allocate(Sptr2d(1:sCPD(1),1:sCPD(2)), STAT=STATUS)
   srcArr = ESMF_ArrayCreate(Sptr2d, ESMF_DATA_REF, RC=STATUS)
   _VERIFY(STATUS)
   srcFLD2D = ESMF_FieldCreate(SRCGrid2D, srcARR,    &
                 horzRelloc = ESMF_CELL_CENTER,         &
                 datamap=dmap,                       &
                 haloWidth=0, &
                 name   = "PS", rc=status )  ! can be any name, all we want
   _VERIFY(STATUS)                           ! is the route handle

   ! DESTINATION FIELD

   call ESMF_FieldGet(dstFLD, grid=grid3D, rc=status)
   _VERIFY(status)

   call ESMF_GridGetDELayout(grid3D, layout, rc=status)
   _VERIFY(STATUS)

   call ESMF_DELayoutGet(layout, deCount=nDEs, deCountPerDim = decpd, rc=status)

   NX = decpd(1)
   NY = decpd(2)

   allocate (ims(NX), jms(NY), stat=status)
   _VERIFY(STATUS)

   call ESMF_GridGet(grid3D, &
                     horzRelLoc=ESMF_CELL_CENTER, vertRelLoc=ESMF_CELL_CELL, &
                     globalCellCountPerDim=gccpd, &
                     minGlobalCoordPerDim=min, &
                     maxGlobalCoordPerDim=max, &
                     rc=status)
   _VERIFY(status)

   ! this is probably incorrect unless parent grid is XYuni
   !DSTGrid2D = ESMF_GridCreateHorzXYUni( &
   !                           counts=gccpd(1:2), &
   !                           minGlobalCoordPerDim=min, &
   !                           maxGlobalCoordPerDim=max, &
   !                           horzStagger=ESMF_GRID_HORZ_STAGGER_A, &
   !                           periodic=(/ESMF_TRUE, ESMF_FALSE/),     &
   !                           name="DST 2D grid", rc=status)
   !_VERIFY(status)
   ! instead use the following ...
   deltaX = 2.0*MAPL_PI/gccpd(1)
   deltaY = MAPL_PI/(gccpd(2)-1)
   DSTGrid2D = ESMF_GridCreateHorzLatLonUni(         &
         counts = gccpd(1:2),        &
         minGlobalCoordPerDim=min(1:2),     &
         deltaPerDim=(/deltaX, deltaY /),        &
         horzStagger=ESMF_Grid_Horz_Stagger_A,   &
         periodic=(/ESMF_TRUE, ESMF_FALSE/),     &
         name='DST 2D grid', rc=status)
   _VERIFY(STATUS)

   allocate (AI(nDEs,gridRank), stat=status)
   _VERIFY(STATUS)

   if (gridRank == 3) then
      call ESMF_GridGetAllAxisIndex(grid3d, &
                                    horzRelLoc=ESMF_CELL_CENTER, &
                                    vertRelLoc=ESMF_CELL_CELL, &
                                    globalAI=AI, rc=status)
   else
      call ESMF_GridGetAllAxisIndex(grid3d, &
                                    horzRelLoc=ESMF_CELL_CENTER, &
                                    globalAI=AI, rc=status)
   end if
   _VERIFY(STATUS)

   JY = 1
   DO IX = 1, NX

!ALT: this is a workaround to compute deId from Layout coords
      de = (jy-1)*NX + ix

      ims(IX) = AI(de,1)%max - AI(de,1)%min + 1
   END DO

   IX = 1
   DO JY = 1, NY

!ALT: same workaround as above
      de = (jy-1)*NX + ix

      jms(JY) = AI(de,2)%max - AI(de,2)%min + 1
   END DO

   deallocate(AI)

   if (verbose .and. MAPL_AM_I_Root()) then
      print *, 'ims=', ims
      print *, 'jms=', jms
   end if

   call ESMF_GridDistribute(DSTGrid2D, delayout=layout, &
                            countsPerDEDim1=ims,        &
                            countsPerDEDim2=jms,        &
                            rc=status)
   _VERIFY(status)

   deallocate(jms, ims)

   call ESMF_GridGetDELocalInfo(DSTGrid2D, &
         horzRelLoc=ESMF_CELL_CENTER, &
!         vertRelLoc=ESMF_CELL_CELL, &
         localCellCountPerDim=dCPD,RC=STATUS)
   _VERIFY(STATUS)
   allocate(Dptr2d(1:dCPD(1),1:dCPD(2)), STAT=STATUS)
   _VERIFY(STATUS)
   dstArr = ESMF_ArrayCreate(Dptr2d, ESMF_DATA_REF, RC=STATUS)
   _VERIFY(STATUS)
   dstFLD2D = ESMF_FieldCreate(DSTGrid2D, dstARR,    &
                 horzRelloc = ESMF_CELL_CENTER,         &
                 datamap=dmap,                       &
                 haloWidth=0, &
                 name   = "PS", rc=status )  ! can be any name, all we want
   _VERIFY(STATUS)                           ! is the correct route handle

   call ESMF_FIELDRegridStore(srcFLD2D, dstFLD2D, vm, rh,  &
                             regridmethod=ESMF_REGRID_METHOD_BILINEAR, &
                             rc=status)
   _VERIFY(status)

   deallocate(Sptr2d, stat=status)
   _VERIFY(STATUS)
   deallocate(Dptr2d, stat=status)
   _VERIFY(STATUS)
#endif

   if (present(rc)) then
      rc = 0
   end if

   end subroutine ESMFL_RegridStore