cpld_gridgen  1.11.0
All Data Structures Files Functions Variables Pages
postwgts.F90
Go to the documentation of this file.
1 
7 module postwgts
8 
9  use esmf
10 
11  use gengrid_kinds, only : cl,cm,cs
12  use grdvars, only : nv
13  use charstrings, only : dirout, res, staggerlocs, logmsg
14  use netcdf
15 
16  implicit none
17 
18 contains
23 
24  subroutine make_postwgts
25 
26  ! local variables
27  character(len=CL) :: fsrc, fdst, fwgt
28  character(len= 2) :: cstagger
29 
30  character(len=CM), dimension(2) :: methodname = (/'conserve', 'bilinear'/)
31 
32  type(ESMF_RegridMethod_Flag) :: method
33  ! the number of possible destination grids depends on the source grid resolution
34  integer :: k,rc,nd,ndest
35  character(len=CS), allocatable, dimension(:) :: destgrds
36 
37  !---------------------------------------------------------------------
38  ! set the destination grids
39  !---------------------------------------------------------------------
40 
41  if(trim(res) .eq. '400')return
42 
43  if(trim(res) .eq. '100')then
44  ndest = 1
45  allocate(destgrds(ndest))
46  destgrds = (/'1p0 '/)
47  end if
48  if(trim(res) .eq. '050')then
49  ndest = 2
50  allocate(destgrds(ndest))
51  destgrds = (/'1p0 ', '0p5 '/)
52  end if
53  if(trim(res) .eq. '025')then
54  ndest = 3
55  allocate(destgrds(ndest))
56  destgrds = (/'1p0 ', '0p5 ', '0p25'/)
57  end if
58 
59  !---------------------------------------------------------------------
60  ! use ESMF to create the weights for unstaggering the points onto
61  ! the Ct staggers for post; the destination is always Ct
62  !---------------------------------------------------------------------
63 
64  method=esmf_regridmethod_bilinear
65  fdst = trim(dirout)//'/'//'Ct.mx'//trim(res)//'_SCRIP.nc'
66  do k = 2,nv
67  cstagger = trim(staggerlocs(k))
68  fsrc = trim(dirout)//'/'//trim(cstagger)//'.mx'//trim(res)//'_SCRIP.nc'
69  fwgt = trim(dirout)//'/'//'tripole.mx'//trim(res)//'.'//trim(cstagger)//'.to.Ct.bilinear.nc'
70  logmsg = 'creating weight file '//trim(fwgt)
71  print '(a)',trim(logmsg)
72 
73  call esmf_regridweightgen(srcfile=trim(fsrc),dstfile=trim(fdst), &
74  weightfile=trim(fwgt), regridmethod=method, &
75  ignoredegenerate=.true., &
76  unmappedaction=esmf_unmappedaction_ignore, rc=rc)
77  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
78  line=__line__, file=__file__)) call esmf_finalize(endflag=esmf_end_abort)
79  end do
80 
81  !---------------------------------------------------------------------
82  ! use ESMF to create the weights from the Ct tripole to the rectilinear
83  ! grids with conservative and bilinear methods for post; the source
84  ! file is always Ct
85  !---------------------------------------------------------------------
86 
87  do nd = 1,ndest
88  fsrc = trim(dirout)//'/'//'Ct.mx'//trim(res)//'_SCRIP.nc'
89  fdst = trim(dirout)//'/rect.'//trim(destgrds(nd))//'_SCRIP.nc'
90 
91  do k = 1,size(methodname)
92  if(trim(methodname(k)) .eq. 'bilinear')method=esmf_regridmethod_bilinear
93  if(trim(methodname(k)) .eq. 'conserve')method=esmf_regridmethod_conserve
94 
95  fwgt = trim(dirout)//'/'//'tripole.mx'//trim(res)//'.Ct.to.rect.'//trim(destgrds(nd)) &
96  //'.'//trim(methodname(k))//'.nc'
97  logmsg = 'creating weight file '//trim(fwgt)
98  print '(a)',trim(logmsg)
99 
100  call esmf_regridweightgen(srcfile=trim(fsrc),dstfile=trim(fdst), &
101  weightfile=trim(fwgt), regridmethod=method, &
102  ignoredegenerate=.true., &
103  unmappedaction=esmf_unmappedaction_ignore, rc=rc)
104  if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
105  line=__line__, file=__file__)) call esmf_finalize(endflag=esmf_end_abort)
106  end do
107  end do
108 
109  deallocate(destgrds)
110 
111  end subroutine make_postwgts
112 end module postwgts