22 use gengrid_kinds, only: cl, cs, dbl_kind, real_kind, int_kind
31 use charstrings, only: logmsg, res, dirsrc, dirout, atmres, fv3dir, editsfile
32 use charstrings, only: maskfile, maskname, topofile, toponame, editsfile, staggerlocs, cdate, history
39 real(dbl_kind) :: dxt, dyt
41 real(kind=dbl_kind),
parameter :: pi = 3.14159265358979323846_dbl_kind
42 real(kind=dbl_kind),
parameter :: deg2rad = pi/180.0_dbl_kind
44 real(real_kind),
allocatable,
dimension(:,:) :: ww3dpth
45 integer(int_kind),
allocatable,
dimension(:,:) :: ww3mask
47 character(len=CL) :: fsrc, fdst, fwgt
48 character(len= 2) :: cstagger
50 integer :: rc,ncid,id,xtype
51 integer :: i,j,k,i2,j2
53 integer :: localpet, npet
54 logical :: fexist = .false.
56 type(esmf_regridmethod_flag
) :: method
60 character(len= 6) :: i4fmt =
'(i4.4)'
61 character(len=CS) :: form1
62 character(len=CS) :: form2
63 character(len= 6) :: cnx
69 call esmf_vmgetglobal(vm, rc=rc)
70 call esmf_initialize(vm=vm, logkindflag=esmf_logkind_multi, rc=rc)
71 call esmf_vmget(vm, localpet=localpet, pecount=npet, rc=rc)
72 if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
73 line=__line__, file=__file__)) call esmf_finalize(endflag=esmf_end_abort)
75 if (localpet == 0) mastertask=.true.
77 print *,npet,
' More than one task specified; Aborting '
78 call esmf_finalize(endflag=esmf_end_abort)
88 print
'(a,2i6)',
' output grid requested ',ni,nj
89 print
'(a,2i6)',
' supergrid size used ', nx,ny
90 print
'(a)',
' output grid tag '//trim(res)
91 print
'(a)',
' supergrid source directory '//trim(dirsrc)
92 print
'(a)',
' output grid directory '//trim(dirout)
93 print
'(a)',
' atm resolution '//trim(atmres)
94 print
'(a,i6)',
' fv3 tile grid size ',npx
95 print
'(a)',
' atm mosaic directory '//trim(fv3dir)
96 print
'(a)',
' MOM6 topography file '//trim(topofile)
97 print
'(a)',
' MOM6 edits file '//trim(editsfile)
98 print *,
'editmask flag ',editmask
99 print *,
'debug flag ',debug
100 print *,
'do_postwgts flag ',do_postwgts
106 call esmf_logwrite(
"Starting gen_fixgrid", esmf_logmsg_info)
111 ivertcu = ivertct + 1; jvertcu = jvertct + 0
112 ivertcv = ivertct + 0; jvertcv = jvertct + 1
113 ivertbu = ivertct + 1; jvertbu = jvertct + 1
116 print
'(a8,4i6)',
'iVertCt ',(ivertct(i),i=1,4)
117 print
'(a8,4i6)',
'jVertCt ',(jvertct(i),i=1,4)
119 print
'(a8,4i6)',
'iVertCu ',(ivertcu(i),i=1,4)
120 print
'(a8,4i6)',
'jVertCu ',(jvertcu(i),i=1,4)
122 print
'(a8,4i6)',
'iVertCv ',(ivertcv(i),i=1,4)
123 print
'(a8,4i6)',
'jVertCv ',(jvertcv(i),i=1,4)
125 print
'(a8,4i6)',
'iVertBu ',(ivertbu(i),i=1,4)
126 print
'(a8,4i6)',
'jVertBu ',(jvertbu(i),i=1,4)
130 latct_vert = -9999.0 ; lonct_vert = -9999.0
131 latcu_vert = -9999.0 ; loncu_vert = -9999.0
132 latcv_vert = -9999.0 ; loncv_vert = -9999.0
133 latbu_vert = -9999.0 ; lonbu_vert = -9999.0
139 fsrc = trim(dirsrc)//
'/'//trim(maskfile)
141 rc = nf90_open(fsrc, nf90_nowrite, ncid)
143 print
'(a)',
'reading ocean mask from '//trim(fsrc)
144 if(rc .ne. 0)print
'(a)',
'nf90_open = '//trim(nf90_strerror(rc))
147 wet4 = 0.0; wet8 = 0.0
148 rc = nf90_inq_varid(ncid, trim(maskname), id)
149 rc = nf90_inquire_variable(ncid, id, xtype=xtype)
150 if(xtype .eq. 5)rc = nf90_get_var(ncid, id, wet4)
151 if(xtype .eq. 6)rc = nf90_get_var(ncid, id, wet8)
152 rc = nf90_close(ncid)
154 if(xtype.eq. 6)wet4 =
real(wet8,4)
163 fsrc = trim(dirsrc)//
'/'//trim(topofile)
165 rc = nf90_open(fsrc, nf90_nowrite, ncid)
167 print
'(a)',
'reading ocean topography from '//trim(fsrc)
168 if(rc .ne. 0)print
'(a)',
'nf90_open = '//trim(nf90_strerror(rc))
172 rc = nf90_inq_varid(ncid, trim(toponame), id)
173 rc = nf90_inquire_variable(ncid, id, xtype=xtype)
174 if(xtype .eq. 5)rc = nf90_get_var(ncid, id, dp4)
175 if(xtype .eq. 6)rc = nf90_get_var(ncid, id, dp8)
176 rc = nf90_close(ncid)
178 if(xtype.eq. 6)dp4 =
real(dp8,4)
188 if(trim(editsfile) ==
'none')
then
189 print
'(a)',
'Need a valid editsfile to make mask edits '
193 fsrc = trim(dirsrc)//
'/'//trim(editsfile)
194 fdst = trim(dirout)//
'/'//
'ufs.'//trim(editsfile)
204 fsrc = trim(dirsrc)//
'/'//trim(editsfile)
205 if(editmask)fsrc = trim(dirout)//
'/'//
'ufs.'//trim(editsfile)
212 fsrc = trim(dirsrc)//
'/'//
'ocean_hgrid.nc'
214 rc = nf90_open(fsrc, nf90_nowrite, ncid)
216 print
'(a)',
'reading supergrid from '//trim(fsrc)
217 if(rc .ne. 0)print
'(a)',
'nf90_open = '//trim(nf90_strerror(rc))
220 rc = nf90_inq_varid(ncid,
'x', id)
221 rc = nf90_get_var(ncid, id, x)
223 rc = nf90_inq_varid(ncid,
'y', id)
224 rc = nf90_get_var(ncid, id, y)
226 rc = nf90_inq_varid(ncid,
'dx', id)
227 rc = nf90_get_var(ncid, id, dx)
229 rc = nf90_inq_varid(ncid,
'dy', id)
230 rc = nf90_get_var(ncid, id, dy)
232 rc = nf90_close(ncid)
235 sg_maxlat = maxval(y)
251 ulon(i,j) = x(i2,j2)*deg2rad
252 ulat(i,j) = y(i2,j2)*deg2rad
254 angle(i,j) = -angq(i2,j2)
256 htn(i,j) = (dx(i2-1,j2) + dx(i2,j2))*100._dbl_kind
257 hte(i,j) = (dy(i2,j2-1) + dy(i2,j2))*100._dbl_kind
259 lonbu(i,j) = x(i2,j2)
260 latbu(i,j) = y(i2,j2)
262 lonct(i,j) = x(i2-1,j2-1)
263 loncu(i,j) = x(i2, j2-1)
264 loncv(i,j) = x(i2-1,j2 )
266 latct(i,j) = y(i2-1,j2-1)
267 latcu(i,j) = y(i2, j2-1)
268 latcv(i,j) = y(i2-1,j2 )
270 dxt = dx(i2-1,j2-1) + dx(i2,j2-1)
271 dyt = dy(i2-1,j2-1) + dy(i2-1,j2)
272 areact(i,j) = dxt*dyt
283 print *,
'ANGLET ',minval(anglet),maxval(anglet)
284 print *,
'ANGLE ',minval(angle),maxval(angle)
296 write(logmsg,
'(a,2e12.5)')
'min vals of hte at folds ', &
297 minval(hte(ni/2,:)),minval(hte(ni,:))
298 print
'(a)',trim(logmsg)
302 if(hte(ii,j) .le. 1.0)hte(ii,j) = 0.5*(hte(ii-1,j) + hte(ii+1,j))
304 if(hte(ii,j) .le. 1.0)hte(ii,j) = 0.5*(hte(ii-1,j) + hte( 1,j))
307 write(logmsg,
'(a,2e12.5)')
'min vals of hte at folds ', &
308 minval(hte(ni/2,:)),minval(hte(ni,:))
309 print
'(a)',trim(logmsg)
316 where(lonct .lt. 0.0)lonct = lonct + 360._dbl_kind
317 where(loncu .lt. 0.0)loncu = loncu + 360._dbl_kind
318 where(loncv .lt. 0.0)loncv = loncv + 360._dbl_kind
319 where(lonbu .lt. 0.0)lonbu = lonbu + 360._dbl_kind
330 if(latbu(i,j) .eq. sg_maxlat)ipole(1) = i
333 if(latbu(i,j) .eq. sg_maxlat)ipole(2) = i
336 write(logmsg,
'(a,2i6,2f12.2)')
'poles found at i = ',ipole,latbu(ipole(1),nj), &
338 print
'(a)',trim(logmsg)
344 i2 = ipole(2)+(ipole(1)-i)+1
345 xlonct(i) = lonct(i2,nj)
346 xlatct(i) = latct(i2,nj)
350 i2 = ipole(2)+(ipole(1)-i)
352 xloncu(i) = loncu(i2,nj)
353 xlatcu(i) = latcu(i2,nj)
360 dlatbu(i) = latbu(i,1) + 2.0*(latcu(i,1) - latbu(i,1))
361 dlatcv(i) = latct(i,1) + 2.0*(latct(i,1) - latcv(i,1))
369 call
fill_vertices(2,nj , ivertct,jvertct, latbu,lonbu, latct_vert,lonct_vert)
370 call
fill_bottom(ivertct,jvertct, latbu,lonbu, latct_vert,lonct_vert,dlatbu)
372 call
fill_vertices(2,nj , ivertcu,jvertcu, latcv,loncv, latcu_vert,loncu_vert)
373 call
fill_bottom(ivertcu,jvertcu, latcv,loncv, latcu_vert,loncu_vert,dlatcv)
376 call
fill_vertices(1,nj-1, ivertcv,jvertcv, latcu,loncu, latcv_vert,loncv_vert)
377 call
fill_top(ivertcv,jvertcv, latcu,loncu, latcv_vert,loncv_vert, xlatcu, xloncu)
379 call
fill_vertices(1,nj-1, ivertbu,jvertbu, latct,lonct, latbu_vert,lonbu_vert)
380 call
fill_top(ivertbu,jvertbu, latct,lonct, latbu_vert,lonbu_vert, xlatct, xlonct)
384 if(minval(latct_vert) .lt. -1.e3)stop
385 if(minval(lonct_vert) .lt. -1.e3)stop
386 if(minval(latcu_vert) .lt. -1.e3)stop
387 if(minval(loncu_vert) .lt. -1.e3)stop
388 if(minval(latcv_vert) .lt. -1.e3)stop
389 if(minval(loncv_vert) .lt. -1.e3)stop
390 if(minval(latbu_vert) .lt. -1.e3)stop
391 if(minval(lonbu_vert) .lt. -1.e3)stop
392 deallocate(xlonct, xlatct, xloncu, xlatcu, dlatbu, dlatcv)
399 call date_and_time(date=cdate)
400 history =
'created on '//trim(cdate)//
' from '//trim(fsrc)
403 fdst = trim(dirout)//
'/'//
'tripole.mx'//trim(res)//
'.nc'
407 fdst = trim(dirout)//
'/'//
'grid_cice_NEMS_mx'//trim(res)//
'.nc'
409 deallocate(ulon, ulat, htn, hte)
413 cstagger = trim(staggerlocs(k))
414 fdst = trim(dirout)//
'/'//trim(cstagger)//
'.mx'//trim(res)//
'_SCRIP.nc'
416 logmsg =
'creating SCRIP file '//trim(fdst)
417 print
'(a)',trim(logmsg)
421 deallocate(latcv_vert, loncv_vert)
422 deallocate(latcu_vert, loncu_vert)
423 deallocate(latbu_vert, lonbu_vert)
427 cstagger = trim(staggerlocs(1))
428 fdst= trim(dirout)//
'/'//trim(cstagger)//
'.mx'//trim(res)//
'_SCRIP_land.nc'
430 logmsg =
'creating SCRIP file '//trim(fdst)
431 print
'(a)',trim(logmsg)
434 deallocate(latct_vert, lonct_vert)
442 write(form1,
'(a)')
'('//trim(cnx)//
'f14.8)'
443 write(form2,
'(a)')
'('//trim(cnx)//
'i2)'
445 allocate(ww3mask(1:ni,1:nj)); ww3mask = wet4
446 allocate(ww3dpth(1:ni,1:nj)); ww3dpth = dp4
448 where(latct .ge. maximum_lat)ww3mask = 3
452 open(unit=21,file=trim(dirout)//
'/'//
'ww3.mx'//trim(res)//
'_x.inp',form=
'formatted')
453 open(unit=22,file=trim(dirout)//
'/'//
'ww3.mx'//trim(res)//
'_y.inp',form=
'formatted')
454 open(unit=23,file=trim(dirout)//
'/'//
'ww3.mx'//trim(res)//
'_bottom.inp',form=
'formatted')
455 open(unit=24,file=trim(dirout)//
'/'//
'ww3.mx'//trim(res)//
'_mapsta.inp',form=
'formatted')
457 open(unit=25,file=trim(dirout)//
'/'//
'ww3.mx'//trim(res)//
'_obstr.inp',form=
'formatted')
460 write( 21,trim(form1))lonct(:,j)
461 write( 22,trim(form1))latct(:,j)
464 write( 23,trim(form1))ww3dpth(:,j)
465 write( 24,trim(form2))ww3mask(:,j)
467 write( 25,trim(form2))ww3mask(:,j)*0
468 write( 25,trim(form2))ww3mask(:,j)*0
471 close(21);
close(22);
close(23);
close(24);
close(25)
472 deallocate(ww3mask);
deallocate(ww3dpth)
473 deallocate(wet4, wet8)
481 method=esmf_regridmethod_conserve
482 fsrc = trim(dirout)//
'/'//
'Ct.mx'//trim(res)//
'_SCRIP_land.nc'
483 fdst = trim(fv3dir)//
'/'//trim(atmres)//
'/'//trim(atmres)//
'_mosaic.nc'
484 fwgt = trim(dirout)//
'/'//
'Ct.mx'//trim(res)//
'.to.'//trim(atmres)//
'.nc'
486 logmsg =
'creating weight file '//trim(fwgt)
487 print
'(a)',trim(logmsg)
490 call esmf_regridweightgen(srcfile=trim(fsrc),dstfile=trim(fdst), &
491 weightfile=trim(fwgt), regridmethod=method, &
492 unmappedaction=esmf_unmappedaction_ignore, &
493 ignoredegenerate=.true., &
494 tilefilepath=trim(fv3dir)//
'/'//trim(atmres)//
'/', rc=rc)
495 if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
496 line=__line__, file=__file__)) call esmf_finalize(endflag=esmf_end_abort)
499 logmsg =
'creating mapped ocean mask for '//trim(atmres)
500 print
'(a)',trim(logmsg)
510 if(trim(res) .ne.
'025')
then
511 fsrc = trim(dirout)//
'/'//
'Ct.mx025_SCRIP.nc'
512 inquire(file=trim(fsrc), exist=fexist)
514 method=esmf_regridmethod_nearest_stod
515 fdst = trim(dirout)//
'/'//
'Ct.mx'//trim(res)//
'_SCRIP.nc'
516 fwgt = trim(dirout)//
'/'//
'tripole.mx025.Ct.to.mx'//trim(res)//
'.Ct.neareststod.nc'
518 logmsg =
'creating weight file '//trim(fwgt)
519 print
'(a)',trim(logmsg)
522 call esmf_regridweightgen(srcfile=trim(fsrc),dstfile=trim(fdst), &
523 weightfile=trim(fwgt), regridmethod=method, &
524 ignoredegenerate=.true., &
525 unmappedaction=esmf_unmappedaction_ignore, rc=rc)
526 if (esmf_logfounderror(rctocheck=rc, msg=esmf_logerr_passthru, &
527 line=__line__, file=__file__)) call esmf_finalize(endflag=esmf_end_abort)
530 logmsg =
'ERROR: '//trim(fsrc)//
' is required to generate tripole:triple weights'
531 print
'(a)',trim(logmsg)
547 deallocate(x,y, angq, dx, dy, xsgp1, ysgp1)
548 deallocate(areact, anglet, angle)
549 deallocate(latct, lonct)
550 deallocate(latcv, loncv)
551 deallocate(latcu, loncu)
552 deallocate(latbu, lonbu)
subroutine, public write_scripgrid(fname, cstagger, imask)
Write a SCRIP grid file.
subroutine make_frac_land(src, wgt)
Use ESMF weights to map the ocean land mask to the FV3 tiles and write the mapped mask to 6 tile file...
subroutine allocate_all
Allocate grid variables.
subroutine fill_vertices(jbeg, jend, iVert, jVert, lat, lon, latvert, lonvert)
Fill the vertices for any stagger location between bounding j-rows.
subroutine, public checkseam
Print values across the tripole seam.
subroutine find_ang
Find the rotation angle on center (Ct) grid points.
subroutine make_postwgts
Create the ESMF weights files to remap velocity points from their native stagger location to the cent...
subroutine find_angq
Find the rotation angle on corner grid (Bu) points using the full MOM6 supergrid. ...
subroutine fill_bottom(iVert, jVert, lat, lon, latvert, lonvert, dlat)
Fill the vertices for a stagger location along the bottom j-row.
subroutine, public write_tripolegrid(fname)
Write the tripole grid file.
subroutine, public checkpoint
Print values at specified point.
subroutine, public apply_topoedits(fsrc)
Read the topoedits file and adjust the bathymetry.
subroutine, public add_topoedits(fsrc, fdst)
Read the existing topoedits file, append required topoedits and write a new topoedits file...
program gen_fixgrid
Generate fixed grid files required for coupled model using the MOM6 super grid file and ocean mask fi...
subroutine fill_top(iVert, jVert, lat, lon, latvert, lonvert, xlat, xlon)
Fill the vertices for a stagger location along the top j-row.
subroutine, public checkxlatlon
Print values near the poles and along the domain edges.
subroutine, public write_cicegrid(fname)
Write the CICE6 grid file.