cpld_gridgen  1.7.0
 All Data Structures Files Functions Variables
gen_fixgrid.F90
Go to the documentation of this file.
1 
5 
16 program gen_fixgrid
17 
18  use esmf
19 
20  use grdvars
21  use inputnml
22  use gengrid_kinds, only: cl, cs, dbl_kind, real_kind, int_kind
23  use angles, only: find_angq, find_ang
25  use mapped_mask, only: make_frac_land
26  use postwgts, only: make_postwgts
28  use cicegrid, only: write_cicegrid
29  use scripgrid, only: write_scripgrid
31  use charstrings, only: logmsg, res, dirsrc, dirout, atmres, fv3dir, editsfile
32  use charstrings, only: maskfile, maskname, topofile, toponame, editsfile, staggerlocs, cdate, history
34  use netcdf
35 
36  implicit none
37 
38  ! local variables
39  real(dbl_kind) :: dxt, dyt
40 
41  real(kind=dbl_kind), parameter :: pi = 3.14159265358979323846_dbl_kind
42  real(kind=dbl_kind), parameter :: deg2rad = pi/180.0_dbl_kind
43 
44  real(real_kind), allocatable, dimension(:,:) :: ww3dpth
45  integer(int_kind), allocatable, dimension(:,:) :: ww3mask
46 
47  character(len=CL) :: fsrc, fdst, fwgt
48  character(len= 2) :: cstagger
49 
50  integer :: rc,ncid,id,xtype
51  integer :: i,j,k,i2,j2
52  integer :: ii,jj
53  integer :: localpet, npet
54  logical :: fexist = .false.
55 
56  type(esmf_regridmethod_flag) :: method
57  type(esmf_vm) :: vm
58 
59  !WW3 file format for mod_def generation
60  character(len= 6) :: i4fmt = '(i4.4)'
61  character(len=CS) :: form1
62  character(len=CS) :: form2
63  character(len= 6) :: cnx
64 
65 !-------------------------------------------------------------------------
66 ! Initialize esmf environment.
67 !-------------------------------------------------------------------------
68 
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)
74  mastertask = .false.
75  if (localpet == 0) mastertask=.true.
76  if (npet /= 1) then
77  print *,npet,' More than one task specified; Aborting '
78  call esmf_finalize(endflag=esmf_end_abort)
79  end if
80 
81 !---------------------------------------------------------------------
82 !
83 !---------------------------------------------------------------------
84 
85  call read_inputnml('grid.nml')
86 
87  if(mastertask) then
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
101  print *
102  end if
103 
104  call allocate_all
105 
106  call esmf_logwrite("Starting gen_fixgrid", esmf_logmsg_info)
107 !---------------------------------------------------------------------
108 ! set up the arrays to retrieve the vertices
109 !---------------------------------------------------------------------
110 
111  ivertcu = ivertct + 1; jvertcu = jvertct + 0
112  ivertcv = ivertct + 0; jvertcv = jvertct + 1
113  ivertbu = ivertct + 1; jvertbu = jvertct + 1
114 
115  if(mastertask) then
116  print '(a8,4i6)','iVertCt ',(ivertct(i),i=1,4)
117  print '(a8,4i6)','jVertCt ',(jvertct(i),i=1,4)
118  print *
119  print '(a8,4i6)','iVertCu ',(ivertcu(i),i=1,4)
120  print '(a8,4i6)','jVertCu ',(jvertcu(i),i=1,4)
121  print *
122  print '(a8,4i6)','iVertCv ',(ivertcv(i),i=1,4)
123  print '(a8,4i6)','jVertCv ',(jvertcv(i),i=1,4)
124  print *
125  print '(a8,4i6)','iVertBu ',(ivertbu(i),i=1,4)
126  print '(a8,4i6)','jVertBu ',(jvertbu(i),i=1,4)
127  print *
128  end if
129 
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
134 
135 !---------------------------------------------------------------------
136 ! read the MOM6 land mask
137 !---------------------------------------------------------------------
138 
139  fsrc = trim(dirsrc)//'/'//trim(maskfile)
140 
141  rc = nf90_open(fsrc, nf90_nowrite, ncid)
142  if(mastertask) then
143  print '(a)', 'reading ocean mask from '//trim(fsrc)
144  if(rc .ne. 0)print '(a)', 'nf90_open = '//trim(nf90_strerror(rc))
145  end if
146 
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)
153 
154  if(xtype.eq. 6)wet4 = real(wet8,4)
155 
156  !print *,minval(wet8),maxval(wet8)
157  !print *,minval(wet4),maxval(wet4)
158 
159 !---------------------------------------------------------------------
160 ! read the MOM6 depth file
161 !---------------------------------------------------------------------
162 
163  fsrc = trim(dirsrc)//'/'//trim(topofile)
164 
165  rc = nf90_open(fsrc, nf90_nowrite, ncid)
166  if(mastertask) then
167  print '(a)', 'reading ocean topography from '//trim(fsrc)
168  if(rc .ne. 0)print '(a)', 'nf90_open = '//trim(nf90_strerror(rc))
169  end if
170 
171  dp4 = 0.0; dp8 = 0.0
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)
177 
178  if(xtype.eq. 6)dp4 = real(dp8,4)
179 
180  !print *,minval(dp8),maxval(dp8)
181  !print *,minval(dp4),maxval(dp4)
182 
183  if(editmask)then
184 !---------------------------------------------------------------------
185 ! apply topoedits run time mask changes if required for this config
186 !---------------------------------------------------------------------
187 
188  if(trim(editsfile) == 'none')then
189  print '(a)', 'Need a valid editsfile to make mask edits '
190  stop
191  end if
192 
193  fsrc = trim(dirsrc)//'/'//trim(editsfile)
194  fdst = trim(dirout)//'/'//'ufs.'//trim(editsfile)
195  call add_topoedits(fsrc,fdst)
196  endif
197 
198 !---------------------------------------------------------------------
199 ! MOM6 reads the depth file, applies the topo edits and then adjusts
200 ! depth using masking_depth and min/max depth. This call mimics
201 ! MOM6 routines apply_topography_edits_from_file and limit_topography
202 !---------------------------------------------------------------------
203 
204  fsrc = trim(dirsrc)//'/'//trim(editsfile)
205  if(editmask)fsrc = trim(dirout)//'/'//'ufs.'//trim(editsfile)
206  call apply_topoedits(fsrc)
207 
208 !---------------------------------------------------------------------
209 ! read MOM6 supergrid file
210 !---------------------------------------------------------------------
211 
212  fsrc = trim(dirsrc)//'/'//'ocean_hgrid.nc'
213 
214  rc = nf90_open(fsrc, nf90_nowrite, ncid)
215  if(mastertask) then
216  print '(a)', 'reading supergrid from '//trim(fsrc)
217  if(rc .ne. 0)print '(a)', 'nf90_open = '//trim(nf90_strerror(rc))
218  end if
219 
220  rc = nf90_inq_varid(ncid, 'x', id) !lon
221  rc = nf90_get_var(ncid, id, x)
222 
223  rc = nf90_inq_varid(ncid, 'y', id) !lat
224  rc = nf90_get_var(ncid, id, y)
225 
226  rc = nf90_inq_varid(ncid, 'dx', id)
227  rc = nf90_get_var(ncid, id, dx)
228 
229  rc = nf90_inq_varid(ncid, 'dy', id)
230  rc = nf90_get_var(ncid, id, dy)
231 
232  rc = nf90_close(ncid)
233  !print *,'super grid size ',size(y,1),size(y,2)
234  !print *,'max lat in super grid ',maxval(y)
235  sg_maxlat = maxval(y)
236 
237 !---------------------------------------------------------------------
238 ! find the angle on corners---this requires the supergrid
239 !---------------------------------------------------------------------
240 
241  call find_angq
242 
243 !---------------------------------------------------------------------
244 ! fill grid variables
245 !---------------------------------------------------------------------
246 
247  do j = 1,nj
248  do i = 1,ni
249  i2 = 2*i ; j2 = 2*j
250  !deg->rad
251  ulon(i,j) = x(i2,j2)*deg2rad
252  ulat(i,j) = y(i2,j2)*deg2rad
253  !in rad already
254  angle(i,j) = -angq(i2,j2)
255  !m->cm
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
258  !deg
259  lonbu(i,j) = x(i2,j2)
260  latbu(i,j) = y(i2,j2)
261  !deg
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 )
265  !deg
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 )
269  !m2
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
273  enddo
274  enddo
275 
276 !---------------------------------------------------------------------
277 ! find the angle on centers---this does not requires the supergrid
278 !---------------------------------------------------------------------
279 
280  call find_ang
281 
282  if(mastertask) then
283  print *,'ANGLET ',minval(anglet),maxval(anglet)
284  print *,'ANGLE ',minval(angle),maxval(angle)
285  end if
286 
287 !---------------------------------------------------------------------
288 ! For the 1/4deg grid, hte at j=720 and j = 1440 is identically=0.0 for
289 ! j > 840 (64.0N). These are land points, but since CICE uses hte to
290 ! generate remaining variables, setting them to zero will cause problems
291 ! For 1deg grid, hte at ni/2 and ni are very small O~10-12, so test for
292 ! hte < 1.0
293 !---------------------------------------------------------------------
294 
295  if(mastertask) then
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)
299  end if
300  do j = 1,nj
301  ii = ni/2
302  if(hte(ii,j) .le. 1.0)hte(ii,j) = 0.5*(hte(ii-1,j) + hte(ii+1,j))
303  ii = ni
304  if(hte(ii,j) .le. 1.0)hte(ii,j) = 0.5*(hte(ii-1,j) + hte( 1,j))
305  enddo
306  if(mastertask) then
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)
310  end if
311 
312 !---------------------------------------------------------------------
313 !
314 !---------------------------------------------------------------------
315 
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
320 
321 !---------------------------------------------------------------------
322 ! some basic error checking
323 ! find the i-th index of the poles at j= nj
324 ! the corner points must lie on the pole
325 !---------------------------------------------------------------------
326 
327  ipole = -1
328  j = nj
329  do i = 1,ni/2
330  if(latbu(i,j) .eq. sg_maxlat)ipole(1) = i
331  enddo
332  do i = ni/2+1,ni
333  if(latbu(i,j) .eq. sg_maxlat)ipole(2) = i
334  enddo
335  if(mastertask) then
336  write(logmsg,'(a,2i6,2f12.2)')'poles found at i = ',ipole,latbu(ipole(1),nj), &
337  latbu(ipole(2),nj)
338  print '(a)',trim(logmsg)
339  end if
340 
341  if(mastertask .and. debug)call checkseam
342 
343  do i = 1,ni
344  i2 = ipole(2)+(ipole(1)-i)+1
345  xlonct(i) = lonct(i2,nj)
346  xlatct(i) = latct(i2,nj)
347  enddo
348 
349  do i = 1,ni
350  i2 = ipole(2)+(ipole(1)-i)
351  if(i2 .lt. 1)i2 = ni
352  xloncu(i) = loncu(i2,nj)
353  xlatcu(i) = latcu(i2,nj)
354  enddo
355 
356  if(mastertask .and. debug)call checkxlatlon
357 
358  !approx lat at grid bottom
359  do i = 1,ni
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))
362  enddo
363 
364 !---------------------------------------------------------------------
365 ! fill grid vertices variables
366 !---------------------------------------------------------------------
367 
368  !Ct and Cu grids align in j
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)
371 
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)
374 
375  !Cv and Bu grids align in j
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)
378 
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)
381 
382  if(mastertask .and. debug)call checkpoint
383 
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)
393 
394 !---------------------------------------------------------------------
395 ! write out grid file files
396 !---------------------------------------------------------------------
397 
398  ! create a history attribute
399  call date_and_time(date=cdate)
400  history = 'created on '//trim(cdate)//' from '//trim(fsrc)
401 
402  ! write fix grid
403  fdst = trim(dirout)//'/'//'tripole.mx'//trim(res)//'.nc'
404  call write_tripolegrid(trim(fdst))
405 
406  ! write cice grid
407  fdst = trim(dirout)//'/'//'grid_cice_NEMS_mx'//trim(res)//'.nc'
408  call write_cicegrid(trim(fdst))
409  deallocate(ulon, ulat, htn, hte)
410  ! write scrip grids; only the Ct is required, the remaining
411  ! staggers are used only in the postweights generation
412  do k = 1,nv
413  cstagger = trim(staggerlocs(k))
414  fdst = trim(dirout)//'/'//trim(cstagger)//'.mx'//trim(res)//'_SCRIP.nc'
415  if(mastertask) then
416  logmsg = 'creating SCRIP file '//trim(fdst)
417  print '(a)',trim(logmsg)
418  end if
419  call write_scripgrid(trim(fdst),trim(cstagger))
420  end do
421  deallocate(latcv_vert, loncv_vert)
422  deallocate(latcu_vert, loncu_vert)
423  deallocate(latbu_vert, lonbu_vert)
424 
425  ! write SCRIP file with land mask, used for mapped ocean mask
426  ! and mesh creation
427  cstagger = trim(staggerlocs(1))
428  fdst= trim(dirout)//'/'//trim(cstagger)//'.mx'//trim(res)//'_SCRIP_land.nc'
429  if(mastertask) then
430  logmsg = 'creating SCRIP file '//trim(fdst)
431  print '(a)',trim(logmsg)
432  end if
433  call write_scripgrid(trim(fdst),trim(cstagger),imask=int(wet4))
434  deallocate(latct_vert, lonct_vert)
435 
436 !---------------------------------------------------------------------
437 ! write lat,lon,depth and mask arrays required by ww3 in creating
438 ! mod_def file
439 !---------------------------------------------------------------------
440 
441  write(cnx,i4fmt)nx
442  write(form1,'(a)')'('//trim(cnx)//'f14.8)'
443  write(form2,'(a)')'('//trim(cnx)//'i2)'
444 
445  allocate(ww3mask(1:ni,1:nj)); ww3mask = wet4
446  allocate(ww3dpth(1:ni,1:nj)); ww3dpth = dp4
447 
448  where(latct .ge. maximum_lat)ww3mask = 3
449  !close last row
450  ww3mask(:,nj) = 3
451 
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')
456  ! cice0 .ne. cicen requires obstruction map, should be initialized as zeros (w3grid,ln3032)
457  open(unit=25,file=trim(dirout)//'/'//'ww3.mx'//trim(res)//'_obstr.inp',form='formatted')
458 
459  do j = 1,nj
460  write( 21,trim(form1))lonct(:,j)
461  write( 22,trim(form1))latct(:,j)
462  end do
463  do j = 1,nj
464  write( 23,trim(form1))ww3dpth(:,j)
465  write( 24,trim(form2))ww3mask(:,j)
466  !'obsx' and 'obsy' arrays ???
467  write( 25,trim(form2))ww3mask(:,j)*0
468  write( 25,trim(form2))ww3mask(:,j)*0
469  end do
470 
471  close(21); close(22); close(23); close(24); close(25)
472  deallocate(ww3mask); deallocate(ww3dpth)
473  deallocate(wet4, wet8)
474 
475 !---------------------------------------------------------------------
476 ! use ESMF regridding to produce mapped ocean mask; first generate
477 ! conservative regrid weights from ocean to tiles; then generate the
478 ! tiled files containing the mapped ocean mask
479 !---------------------------------------------------------------------
480 
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'
485  if(mastertask) then
486  logmsg = 'creating weight file '//trim(fwgt)
487  print '(a)',trim(logmsg)
488  end if
489 
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)
497 
498  if(mastertask) then
499  logmsg = 'creating mapped ocean mask for '//trim(atmres)
500  print '(a)',trim(logmsg)
501  end if
502  call make_frac_land(trim(fsrc), trim(fwgt))
503 
504 !---------------------------------------------------------------------
505 ! use ESMF to find the tripole:tripole weights for creation
506 ! of CICE ICs; the source grid is always mx025; don't create this
507 ! file if destination is also mx025
508 !---------------------------------------------------------------------
509 
510  if(trim(res) .ne. '025') then
511  fsrc = trim(dirout)//'/'//'Ct.mx025_SCRIP.nc'
512  inquire(file=trim(fsrc), exist=fexist)
513  if (fexist ) then
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'
517  if(mastertask) then
518  logmsg = 'creating weight file '//trim(fwgt)
519  print '(a)',trim(logmsg)
520  end if
521 
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)
528  else
529  if(mastertask) then
530  logmsg = 'ERROR: '//trim(fsrc)//' is required to generate tripole:triple weights'
531  print '(a)',trim(logmsg)
532  end if
533  stop
534  end if
535  end if
536 
537 !---------------------------------------------------------------------
538 !
539 !---------------------------------------------------------------------
540 
541  if(do_postwgts)call make_postwgts
542 
543 !---------------------------------------------------------------------
544 ! clean up
545 !---------------------------------------------------------------------
546 
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)
553 
554 end program gen_fixgrid
subroutine, public write_scripgrid(fname, cstagger, imask)
Write a SCRIP grid file.
Definition: scripgrid.F90:34
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...
Definition: mapped_mask.F90:26
subroutine allocate_all
Allocate grid variables.
Definition: grdvars.F90:173
subroutine fill_vertices(jbeg, jend, iVert, jVert, lat, lon, latvert, lonvert)
Fill the vertices for any stagger location between bounding j-rows.
Definition: vertices.F90:33
subroutine, public checkseam
Print values across the tripole seam.
Definition: debugprint.F90:27
subroutine find_ang
Find the rotation angle on center (Ct) grid points.
Definition: angles.F90:133
subroutine make_postwgts
Create the ESMF weights files to remap velocity points from their native stagger location to the cent...
Definition: postwgts.F90:24
subroutine find_angq
Find the rotation angle on corner grid (Bu) points using the full MOM6 supergrid. ...
Definition: angles.F90:25
subroutine fill_bottom(iVert, jVert, lat, lon, latvert, lonvert, dlat)
Fill the vertices for a stagger location along the bottom j-row.
Definition: vertices.F90:68
subroutine, public write_tripolegrid(fname)
Write the tripole grid file.
Definition: tripolegrid.F90:33
subroutine, public checkpoint
Print values at specified point.
Definition: debugprint.F90:170
subroutine, public apply_topoedits(fsrc)
Read the topoedits file and adjust the bathymetry.
Definition: topoedits.F90:165
subroutine, public add_topoedits(fsrc, fdst)
Read the existing topoedits file, append required topoedits and write a new topoedits file...
Definition: topoedits.F90:32
program gen_fixgrid
Generate fixed grid files required for coupled model using the MOM6 super grid file and ocean mask fi...
Definition: gen_fixgrid.F90:16
subroutine fill_top(iVert, jVert, lat, lon, latvert, lonvert, xlat, xlon)
Fill the vertices for a stagger location along the top j-row.
Definition: vertices.F90:114
subroutine, public checkxlatlon
Print values near the poles and along the domain edges.
Definition: debugprint.F90:108
subroutine, public write_cicegrid(fname)
Write the CICE6 grid file.
Definition: cicegrid.F90:28
subroutine read_inputnml(fname)
Read input namelist file.
Definition: inputnml.F90:24