MOM6
MOM_horizontal_regridding.F90
Go to the documentation of this file.
1 !> Horizontal interpolation
3 
4 ! This file is part of MOM6. See LICENSE.md for the license.
5 
6 use mom_debugging, only : hchksum
7 use mom_coms, only : max_across_pes, min_across_pes
8 use mom_cpu_clock, only : cpu_clock_id, cpu_clock_begin, cpu_clock_end
9 use mom_cpu_clock, only : clock_routine, clock_loop
10 use mom_domains, only : pass_var, pass_vector, sum_across_pes, broadcast
11 use mom_domains, only : root_pe, to_all, scalar_pair, cgrid_ne, agrid
12 use mom_error_handler, only : mom_mesg, mom_error, fatal, warning, is_root_pe
15 use mom_file_parser, only : log_version
16 use mom_get_input, only : directories
18 use mom_io, only : close_file, fieldtype, file_exists
19 use mom_io, only : open_file, read_data, read_axis_data, single_file, multiple
20 use mom_io, only : slasher, vardesc, write_field
22 use mom_time_manager, only : time_type, get_external_field_size
23 use mom_time_manager, only : init_external_field, time_interp_external
24 use mom_time_manager, only : get_external_field_axes, get_external_field_missing
26 use mpp_io_mod, only : axistype
27 use mpp_domains_mod, only : mpp_global_field, mpp_get_compute_domain
28 use mpp_mod, only : mpp_broadcast,mpp_root_pe,mpp_sync,mpp_sync_self
29 use mpp_mod, only : mpp_max
30 use horiz_interp_mod, only : horiz_interp_new, horiz_interp,horiz_interp_type
31 use horiz_interp_mod, only : horiz_interp_init, horiz_interp_del
32 
33 use mpp_io_mod, only : mpp_get_axis_data
34 use mpp_io_mod, only : mpp_single
35 use netcdf
36 
37 implicit none ; private
38 
39 #include <MOM_memory.h>
40 
42 
43 ! character(len=40) :: mdl = "MOM_horizontal_regridding" ! This module's name.
44 
45 !> Fill grid edges
46 interface fill_boundaries
47  module procedure fill_boundaries_real
48  module procedure fill_boundaries_int
49 end interface
50 
51 !> Extrapolate and interpolate data
55 end interface
56 
57 contains
58 
59 !> Write to the terminal some basic statistics about the k-th level of an array
60 subroutine mystats(array, missing, is, ie, js, je, k, mesg)
61  real, dimension(:,:), intent(in) :: array !< input array (ND)
62  real, intent(in) :: missing !< missing value (ND)
63  !!@{
64  !> Horizontal loop bounds to calculate statistics for
65  integer :: is,ie,js,je
66  !!@}
67  integer :: k !< Level to calculate statistics for
68  character(len=*) :: mesg !< Label to use in message
69  ! Local variables
70  real :: mina, maxa
71  integer :: i,j
72  logical :: found
73  character(len=120) :: lmesg
74  mina = 9.e24 ; maxa = -9.e24 ; found = .false.
75 
76  do j=js,je ; do i=is,ie
77  if (array(i,j) /= array(i,j)) stop 'Nan!'
78  if (abs(array(i,j)-missing) > 1.e-6*abs(missing)) then
79  if (found) then
80  mina = min(mina, array(i,j))
81  maxa = max(maxa, array(i,j))
82  else
83  found = .true.
84  mina = array(i,j)
85  maxa = array(i,j)
86  endif
87  endif
88  enddo ; enddo
89  call min_across_pes(mina)
90  call max_across_pes(maxa)
91  if (is_root_pe()) then
92  write(lmesg(1:120),'(2(a,es12.4),a,i3,x,a)') &
93  'init_from_Z: min=',mina,' max=',maxa,' Level=',k,trim(mesg)
94  call mom_mesg(lmesg,2)
95  endif
96 
97 end subroutine mystats
98 
99 !> Use ICE-9 algorithm to populate points (fill=1) with
100 !! valid data (good=1). If no information is available,
101 !! Then use a previous guess (prev). Optionally (smooth)
102 !! blend the filled points to achieve a more desirable result.
103 subroutine fill_miss_2d(aout, good, fill, prev, G, smooth, num_pass, relc, crit, keep_bug, debug)
104  use mom_coms, only : sum_across_pes
105 
106  type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure.
107  real, dimension(SZI_(G),SZJ_(G)), &
108  intent(inout) :: aout !< The array with missing values to fill
109  real, dimension(SZI_(G),SZJ_(G)), &
110  intent(in) :: good !< Valid data mask for incoming array
111  !! (1==good data; 0==missing data).
112  real, dimension(SZI_(G),SZJ_(G)), &
113  intent(in) :: fill !< Same shape array of points which need
114  !! filling (1==fill;0==dont fill)
115  real, dimension(SZI_(G),SZJ_(G)), &
116  optional, intent(in) :: prev !< First guess where isolated holes exist.
117  logical, optional, intent(in) :: smooth !< If present and true, apply a number of
118  !! Laplacian iterations to the interpolated data
119  integer, optional, intent(in) :: num_pass !< The maximum number of iterations
120  real, optional, intent(in) :: relc !< A relaxation coefficient for Laplacian (ND)
121  real, optional, intent(in) :: crit !< A minimal value for deltas between iterations.
122  logical, optional, intent(in) :: keep_bug !< Use an algorithm with a bug that dates
123  !! to the "sienna" code release.
124  logical, optional, intent(in) :: debug !< If true, write verbose debugging messages.
125 
126  real, dimension(SZI_(G),SZJ_(G)) :: b,r
127  real, dimension(SZI_(G),SZJ_(G)) :: fill_pts, good_, good_new
128 
129  character(len=256) :: mesg ! The text of an error message
130  integer :: i,j,k
131  real :: east,west,north,south,sor
132  real :: ge,gw,gn,gs,ngood
133  logical :: do_smooth,siena_bug
134  real :: nfill, nfill_prev
135  integer, parameter :: num_pass_default = 10000
136  real, parameter :: relc_default = 0.25, crit_default = 1.e-3
137 
138  integer :: npass
139  integer :: is, ie, js, je
140  real :: relax_coeff, acrit, ares
141  logical :: debug_it
142 
143  debug_it=.false.
144  if (PRESENT(debug)) debug_it=debug
145 
146  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
147 
148  npass = num_pass_default
149  if (PRESENT(num_pass)) npass = num_pass
150 
151  relax_coeff = relc_default
152  if (PRESENT(relc)) relax_coeff = relc
153 
154  acrit = crit_default
155  if (PRESENT(crit)) acrit = crit
156 
157  siena_bug=.false.
158  if (PRESENT(keep_bug)) siena_bug = keep_bug
159 
160  do_smooth=.false.
161  if (PRESENT(smooth)) do_smooth=smooth
162 
163  fill_pts(:,:) = fill(:,:)
164 
165  nfill = sum(fill(is:ie,js:je))
166  call sum_across_pes(nfill)
167 
168  nfill_prev = nfill
169  good_(:,:) = good(:,:)
170  r(:,:) = 0.0
171 
172  do while (nfill > 0.0)
173 
174  call pass_var(good_,g%Domain)
175  call pass_var(aout,g%Domain)
176 
177  b(:,:)=aout(:,:)
178  good_new(:,:)=good_(:,:)
179 
180  do j=js,je ; do i=is,ie
181 
182  if (good_(i,j) == 1.0 .or. fill(i,j) == 0.) cycle
183 
184  ge=good_(i+1,j) ; gw=good_(i-1,j)
185  gn=good_(i,j+1) ; gs=good_(i,j-1)
186  east=0.0 ; west=0.0 ; north=0.0 ; south=0.0
187  if (ge == 1.0) east = aout(i+1,j)*ge
188  if (gw == 1.0) west = aout(i-1,j)*gw
189  if (gn == 1.0) north = aout(i,j+1)*gn
190  if (gs == 1.0) south = aout(i,j-1)*gs
191 
192  ngood = ge+gw+gn+gs
193  if (ngood > 0.) then
194  b(i,j)=(east+west+north+south)/ngood
195  !### Replace this with
196  ! b(i,j) = ((east+west) + (north+south))/ngood
197  fill_pts(i,j) = 0.0
198  good_new(i,j) = 1.0
199  endif
200  enddo ; enddo
201 
202  aout(is:ie,js:je) = b(is:ie,js:je)
203  good_(is:ie,js:je) = good_new(is:ie,js:je)
204  nfill_prev = nfill
205  nfill = sum(fill_pts(is:ie,js:je))
206  call sum_across_pes(nfill)
207 
208  if (nfill == nfill_prev .and. PRESENT(prev)) then
209  do j=js,je ; do i=is,ie ; if (fill_pts(i,j) == 1.0) then
210  aout(i,j) = prev(i,j)
211  fill_pts(i,j) = 0.0
212  endif ; enddo ; enddo
213  elseif (nfill == nfill_prev) then
214  call mom_error(warning, &
215  'Unable to fill missing points using either data at the same vertical level from a connected basin'//&
216  'or using a point from a previous vertical level. Make sure that the original data has some valid'//&
217  'data in all basins.', .true.)
218  write(mesg,*) 'nfill=',nfill
219  call mom_error(warning, mesg, .true.)
220  endif
221 
222  nfill = sum(fill_pts(is:ie,js:je))
223  call sum_across_pes(nfill)
224 
225  enddo
226 
227  if (do_smooth) then ; do k=1,npass
228  call pass_var(aout,g%Domain)
229  do j=js,je ; do i=is,ie
230  if (fill(i,j) == 1) then
231  east = max(good(i+1,j),fill(i+1,j)) ; west = max(good(i-1,j),fill(i-1,j))
232  north = max(good(i,j+1),fill(i,j+1)) ; south = max(good(i,j-1),fill(i,j-1))
233  r(i,j) = relax_coeff*(south*aout(i,j-1)+north*aout(i,j+1) + &
234  west*aout(i-1,j)+east*aout(i+1,j) - &
235  (south+north+west+east)*aout(i,j))
236  !### Appropriate parentheses should be added here, but they will change answers.
237  ! r(i,j) = relax_coeff*( ((south*aout(i,j-1) + north*aout(i,j+1)) + &
238  ! (west*aout(i-1,j)+east*aout(i+1,j))) - &
239  ! ((south+north)+(west+east))*aout(i,j) )
240  else
241  r(i,j) = 0.
242  endif
243  enddo ; enddo
244  ares = 0.0
245  do j=js,je ; do i=is,ie
246  aout(i,j) = r(i,j) + aout(i,j)
247  ares = max(ares, abs(r(i,j)))
248  enddo ; enddo
249  call max_across_pes(ares)
250  if (ares <= acrit) exit
251  enddo ; endif
252 
253  do j=js,je ; do i=is,ie
254  if (good_(i,j) == 0.0 .and. fill_pts(i,j) == 1.0) then
255  write(mesg,*) 'In fill_miss, fill, good,i,j= ',fill_pts(i,j),good_(i,j),i,j
256  call mom_error(warning, mesg, .true.)
257  call mom_error(fatal,"MOM_initialize: "// &
258  "fill is true and good is false after fill_miss, how did this happen? ")
259  endif
260  enddo ; enddo
261 
262 end subroutine fill_miss_2d
263 
264 !> Extrapolate and interpolate from a file record
265 subroutine horiz_interp_and_extrap_tracer_record(filename, varnam, conversion, recnum, G, tr_z, &
266  mask_z, z_in, z_edges_in, missing_value, reentrant_x, &
267  tripolar_n, homogenize, m_to_Z)
268 
269  character(len=*), intent(in) :: filename !< Path to file containing tracer to be
270  !! interpolated.
271  character(len=*), intent(in) :: varnam !< Name of tracer in filee.
272  real, intent(in) :: conversion !< Conversion factor for tracer.
273  integer, intent(in) :: recnum !< Record number of tracer to be read.
274  type(ocean_grid_type), intent(inout) :: G !< Grid object
275  real, allocatable, dimension(:,:,:) :: tr_z !< pointer to allocatable tracer array on local
276  !! model grid and input-file vertical levels.
277  real, allocatable, dimension(:,:,:) :: mask_z !< pointer to allocatable tracer mask array on
278  !! local model grid and input-file vertical levels.
279  real, allocatable, dimension(:) :: z_in !< Cell grid values for input data.
280  real, allocatable, dimension(:) :: z_edges_in !< Cell grid edge values for input data.
281  real, intent(out) :: missing_value !< The missing value in the returned array.
282  logical, intent(in) :: reentrant_x !< If true, this grid is reentrant in the x-direction
283  logical, intent(in) :: tripolar_n !< If true, this is a northern tripolar grid
284  logical, optional, intent(in) :: homogenize !< If present and true, horizontally homogenize data
285  !! to produce perfectly "flat" initial conditions
286  real, optional, intent(in) :: m_to_Z !< A conversion factor from meters to the units
287  !! of depth. If missing, G%bathyT must be in m.
288 
289  ! Local variables
290  real, dimension(:,:), allocatable :: tr_in, tr_inp ! A 2-d array for holding input data on
291  ! native horizontal grid and extended grid
292  ! with poles.
293  real, dimension(:,:), allocatable :: mask_in ! A 2-d mask for extended input grid.
294 
295  real :: PI_180
296  integer :: rcode, ncid, varid, ndims, id, jd, kd, jdp
297  integer :: i,j,k
298  integer, dimension(4) :: start, count, dims, dim_id
299  real, dimension(:,:), allocatable :: x_in, y_in
300  real, dimension(:), allocatable :: lon_in, lat_in
301  real, dimension(:), allocatable :: lat_inp, last_row
302  real :: max_lat, min_lat, pole, max_depth, npole
303  real :: roundoff ! The magnitude of roundoff, usually ~2e-16.
304  real :: add_offset, scale_factor
305  logical :: add_np
306  character(len=8) :: laynum
307  type(horiz_interp_type) :: Interp
308  integer :: is, ie, js, je ! compute domain indices
309  integer :: isc,iec,jsc,jec ! global compute domain indices
310  integer :: isg, ieg, jsg, jeg ! global extent
311  integer :: isd, ied, jsd, jed ! data domain indices
312  integer :: id_clock_read
313  character(len=12) :: dim_name(4)
314  logical :: debug=.false.
315  real :: npoints,varAvg
316  real, dimension(SZI_(G),SZJ_(G)) :: lon_out, lat_out, tr_out, mask_out
317  real, dimension(SZI_(G),SZJ_(G)) :: good, fill
318  real, dimension(SZI_(G),SZJ_(G)) :: tr_outf,tr_prev
319  real, dimension(SZI_(G),SZJ_(G)) :: good2,fill2
320  real, dimension(SZI_(G),SZJ_(G)) :: nlevs
321 
322  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
323  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
324  isg = g%isg ; ieg = g%ieg ; jsg = g%jsg ; jeg = g%jeg
325 
326  id_clock_read = cpu_clock_id('(Initialize tracer from Z) read', grain=clock_loop)
327 
328 
329  if (allocated(tr_z)) deallocate(tr_z)
330  if (allocated(mask_z)) deallocate(mask_z)
331  if (allocated(z_edges_in)) deallocate(z_edges_in)
332 
333  pi_180=atan(1.0)/45.
334 
335  ! Open NetCDF file and if present, extract data and spatial coordinate information
336  ! The convention adopted here requires that the data be written in (i,j,k) ordering.
337 
338  call cpu_clock_begin(id_clock_read)
339 
340  rcode = nf90_open(filename, nf90_nowrite, ncid)
341  if (rcode /= 0) call mom_error(fatal,"error opening file "//trim(filename)//&
342  " in hinterp_extrap")
343  rcode = nf90_inq_varid(ncid, varnam, varid)
344  if (rcode /= 0) call mom_error(fatal,"error finding variable "//trim(varnam)//&
345  " in file "//trim(filename)//" in hinterp_extrap")
346 
347  rcode = nf90_inquire_variable(ncid, varid, ndims=ndims, dimids=dims)
348  if (rcode /= 0) call mom_error(fatal,'error inquiring dimensions hinterp_extrap')
349  if (ndims < 3) call mom_error(fatal,"Variable "//trim(varnam)//" in file "// &
350  trim(filename)//" has too few dimensions.")
351 
352  rcode = nf90_inquire_dimension(ncid, dims(1), dim_name(1), len=id)
353  if (rcode /= 0) call mom_error(fatal,"error reading dimension 1 data for "// &
354  trim(varnam)//" in file "// trim(filename)//" in hinterp_extrap")
355  rcode = nf90_inq_varid(ncid, dim_name(1), dim_id(1))
356  if (rcode /= 0) call mom_error(fatal,"error finding variable "//trim(dim_name(1))//&
357  " in file "//trim(filename)//" in hinterp_extrap")
358  rcode = nf90_inquire_dimension(ncid, dims(2), dim_name(2), len=jd)
359  if (rcode /= 0) call mom_error(fatal,"error reading dimension 2 data for "// &
360  trim(varnam)//" in file "// trim(filename)//" in hinterp_extrap")
361  rcode = nf90_inq_varid(ncid, dim_name(2), dim_id(2))
362  if (rcode /= 0) call mom_error(fatal,"error finding variable "//trim(dim_name(2))//&
363  " in file "//trim(filename)//" in hinterp_extrap")
364  rcode = nf90_inquire_dimension(ncid, dims(3), dim_name(3), len=kd)
365  if (rcode /= 0) call mom_error(fatal,"error reading dimension 3 data for "// &
366  trim(varnam)//" in file "// trim(filename)//" in hinterp_extrap")
367  rcode = nf90_inq_varid(ncid, dim_name(3), dim_id(3))
368  if (rcode /= 0) call mom_error(fatal,"error finding variable "//trim(dim_name(3))//&
369  " in file "//trim(filename)//" in hinterp_extrap")
370 
371  missing_value=0.0
372  rcode = nf90_get_att(ncid, varid, "_FillValue", missing_value)
373  if (rcode /= 0) call mom_error(fatal,"error finding missing value for "//&
374  trim(varnam)//" in file "// trim(filename)//" in hinterp_extrap")
375 
376  rcode = nf90_get_att(ncid, varid, "add_offset", add_offset)
377  if (rcode /= 0) add_offset = 0.0
378 
379  rcode = nf90_get_att(ncid, varid, "scale_factor", scale_factor)
380  if (rcode /= 0) scale_factor = 1.0
381 
382  if (allocated(lon_in)) deallocate(lon_in)
383  if (allocated(lat_in)) deallocate(lat_in)
384  if (allocated(z_in)) deallocate(z_in)
385  if (allocated(z_edges_in)) deallocate(z_edges_in)
386  if (allocated(tr_z)) deallocate(tr_z)
387  if (allocated(mask_z)) deallocate(mask_z)
388 
389  allocate(lon_in(id),lat_in(jd),z_in(kd),z_edges_in(kd+1))
390  allocate(tr_z(isd:ied,jsd:jed,kd), mask_z(isd:ied,jsd:jed,kd))
391 
392  start = 1; count = 1; count(1) = id
393  rcode = nf90_get_var(ncid, dim_id(1), lon_in, start, count)
394  if (rcode /= 0) call mom_error(fatal,"error reading dimension 1 values for var_name "// &
395  trim(varnam)//",dim_name "//trim(dim_name(1))//" in file "// trim(filename)//" in hinterp_extrap")
396  start = 1; count = 1; count(1) = jd
397  rcode = nf90_get_var(ncid, dim_id(2), lat_in, start, count)
398  if (rcode /= 0) call mom_error(fatal,"error reading dimension 2 values for var_name "// &
399  trim(varnam)//",dim_name "//trim(dim_name(2))//" in file "// trim(filename)//" in hinterp_extrap")
400  start = 1; count = 1; count(1) = kd
401  rcode = nf90_get_var(ncid, dim_id(3), z_in, start, count)
402  if (rcode /= 0) call mom_error(fatal,"error reading dimension 3 values for var_name "// &
403  trim(varnam//",dim_name "//trim(dim_name(3)))//" in file "// trim(filename)//" in hinterp_extrap")
404 
405  call cpu_clock_end(id_clock_read)
406 
407  if (present(m_to_z)) then ; do k=1,kd ; z_in(k) = m_to_z * z_in(k) ; enddo ; endif
408 
409  ! extrapolate the input data to the north pole using the northerm-most latitude
410 
411  max_lat = maxval(lat_in)
412  add_np=.false.
413  if (max_lat < 90.0) then
414  add_np=.true.
415  jdp=jd+1
416  allocate(lat_inp(jdp))
417  lat_inp(1:jd)=lat_in(:)
418  lat_inp(jd+1)=90.0
419  deallocate(lat_in)
420  allocate(lat_in(1:jdp))
421  lat_in(:)=lat_inp(:)
422  else
423  jdp=jd
424  endif
425 
426  ! construct level cell boundaries as the mid-point between adjacent centers
427 
428  z_edges_in(1) = 0.0
429  do k=2,kd
430  z_edges_in(k)=0.5*(z_in(k-1)+z_in(k))
431  enddo
432  z_edges_in(kd+1)=2.0*z_in(kd) - z_in(kd-1)
433 
434  call horiz_interp_init()
435 
436  lon_in = lon_in*pi_180
437  lat_in = lat_in*pi_180
438  allocate(x_in(id,jdp),y_in(id,jdp))
439  call meshgrid(lon_in,lat_in, x_in, y_in)
440 
441  lon_out(:,:) = g%geoLonT(:,:)*pi_180
442  lat_out(:,:) = g%geoLatT(:,:)*pi_180
443 
444 
445  allocate(tr_in(id,jd)) ; tr_in(:,:)=0.0
446  allocate(tr_inp(id,jdp)) ; tr_inp(:,:)=0.0
447  allocate(mask_in(id,jdp)) ; mask_in(:,:)=0.0
448  allocate(last_row(id)) ; last_row(:)=0.0
449 
450  max_depth = maxval(g%bathyT)
451  call mpp_max(max_depth)
452 
453  if (z_edges_in(kd+1)<max_depth) z_edges_in(kd+1)=max_depth
454 
455  roundoff = 3.0*epsilon(missing_value)
456 
457  ! loop through each data level and interpolate to model grid.
458  ! after interpolating, fill in points which will be needed
459  ! to define the layers
460  do k=1,kd
461  write(laynum,'(I8)') k ; laynum = adjustl(laynum)
462 
463  if (is_root_pe()) then
464  start = 1 ; start(3) = k ; count(:) = 1 ; count(1) = id ; count(2) = jd
465  rcode = nf90_get_var(ncid,varid, tr_in, start, count)
466  if (rcode /= 0) call mom_error(fatal,"hinterp_and_extract_from_Fie: "//&
467  "error reading level "//trim(laynum)//" of variable "//&
468  trim(varnam)//" in file "// trim(filename))
469 
470  if (add_np) then
471  last_row(:)=tr_in(:,jd); pole=0.0;npole=0.0
472  do i=1,id
473  if (abs(tr_in(i,jd)-missing_value) > abs(roundoff*missing_value)) then
474  pole = pole+last_row(i)
475  npole = npole+1.0
476  endif
477  enddo
478  if (npole > 0) then
479  pole=pole/npole
480  else
481  pole=missing_value
482  endif
483  tr_inp(:,1:jd) = tr_in(:,:)
484  tr_inp(:,jdp) = pole
485  else
486  tr_inp(:,:) = tr_in(:,:)
487  endif
488  endif
489 
490  call mpp_sync()
491  call mpp_broadcast(tr_inp, id*jdp, root_pe())
492  call mpp_sync_self()
493 
494  mask_in=0.0
495 
496  do j=1,jdp
497  do i=1,id
498  if (abs(tr_inp(i,j)-missing_value) > abs(roundoff*missing_value)) then
499  mask_in(i,j) = 1.0
500  tr_inp(i,j) = (tr_inp(i,j)*scale_factor+add_offset) * conversion
501  else
502  tr_inp(i,j) = missing_value
503  endif
504  enddo
505  enddo
506 
507 ! call fms routine horiz_interp to interpolate input level data to model horizontal grid
508 
509  if (k == 1) then
510  call horiz_interp_new(interp,x_in,y_in,lon_out(is:ie,js:je),lat_out(is:ie,js:je), &
511  interp_method='bilinear',src_modulo=.true.)
512  endif
513 
514  if (debug) then
515  call mystats(tr_inp,missing_value, is,ie,js,je,k,'Tracer from file')
516  endif
517 
518  tr_out(:,:) = 0.0
519 
520  call horiz_interp(interp,tr_inp,tr_out(is:ie,js:je), missing_value=missing_value, new_missing_handle=.true.)
521 
522  mask_out=1.0
523  do j=js,je
524  do i=is,ie
525  if (abs(tr_out(i,j)-missing_value) < abs(roundoff*missing_value)) mask_out(i,j)=0.
526  enddo
527  enddo
528 
529  fill = 0.0; good = 0.0
530 
531  npoints = 0 ; varavg = 0.
532  do j=js,je
533  do i=is,ie
534  if (mask_out(i,j) < 1.0) then
535  tr_out(i,j)=missing_value
536  else
537  good(i,j)=1.0
538  npoints = npoints + 1
539  varavg = varavg + tr_out(i,j)
540  endif
541  if (g%mask2dT(i,j) == 1.0 .and. z_edges_in(k) <= g%bathyT(i,j) .and. mask_out(i,j) < 1.0) &
542  fill(i,j)=1.0
543  enddo
544  enddo
545  call pass_var(fill,g%Domain)
546  call pass_var(good,g%Domain)
547 
548  if (debug) then
549  call mystats(tr_out,missing_value, is,ie,js,je,k,'variable from horiz_interp()')
550  endif
551 
552  ! Horizontally homogenize data to produce perfectly "flat" initial conditions
553  if (PRESENT(homogenize)) then
554  if (homogenize) then
555  call sum_across_pes(npoints)
556  call sum_across_pes(varavg)
557  if (npoints>0) then
558  varavg = varavg/real(npoints)
559  endif
560  tr_out(:,:) = varavg
561  endif
562  endif
563 
564  ! tr_out contains input z-space data on the model grid with missing values
565  ! now fill in missing values using "ICE-nine" algorithm.
566  tr_outf(:,:) = tr_out(:,:)
567  if (k==1) tr_prev(:,:) = tr_outf(:,:)
568  good2(:,:) = good(:,:)
569  fill2(:,:) = fill(:,:)
570 
571  call fill_miss_2d(tr_outf, good2, fill2, tr_prev, g, smooth=.true.)
572  call mystats(tr_outf, missing_value, is, ie, js, je, k, 'field from fill_miss_2d()')
573 
574  tr_z(:,:,k) = tr_outf(:,:) * g%mask2dT(:,:)
575  mask_z(:,:,k) = good2(:,:) + fill2(:,:)
576 
577  tr_prev(:,:) = tr_z(:,:,k)
578 
579  if (debug) then
580  call hchksum(tr_prev,'field after fill ',g%HI)
581  endif
582 
583  enddo ! kd
584 
586 
587 !> Extrapolate and interpolate using a FMS time interpolation handle
588 subroutine horiz_interp_and_extrap_tracer_fms_id(fms_id, Time, conversion, G, tr_z, mask_z, &
589  z_in, z_edges_in, missing_value, reentrant_x, &
590  tripolar_n, homogenize, spongeOngrid, m_to_Z)
591 
592  integer, intent(in) :: fms_id !< A unique id used by the FMS time interpolator
593  type(time_type), intent(in) :: Time !< A FMS time type
594  real, intent(in) :: conversion !< Conversion factor for tracer.
595  type(ocean_grid_type), intent(inout) :: G !< Grid object
596  real, allocatable, dimension(:,:,:) :: tr_z !< pointer to allocatable tracer array on local
597  !! model grid and native vertical levels.
598  real, allocatable, dimension(:,:,:) :: mask_z !< pointer to allocatable tracer mask array on
599  !! local model grid and native vertical levels.
600  real, allocatable, dimension(:) :: z_in !< Cell grid values for input data.
601  real, allocatable, dimension(:) :: z_edges_in !< Cell grid edge values for input data. (Intent out)
602  real, intent(out) :: missing_value !< The missing value in the returned array.
603  logical, intent(in) :: reentrant_x !< If true, this grid is reentrant in the x-direction
604  logical, intent(in) :: tripolar_n !< If true, this is a northern tripolar grid
605  logical, optional, intent(in) :: homogenize !< If present and true, horizontally homogenize data
606  !! to produce perfectly "flat" initial conditions
607  logical, optional, intent(in) :: spongeOngrid !< If present and true, the sponge data are on the model grid
608  real, optional, intent(in) :: m_to_Z !< A conversion factor from meters to the units
609  !! of depth. If missing, G%bathyT must be in m.
610 
611  ! Local variables
612  real, dimension(:,:), allocatable :: tr_in,tr_inp !< A 2-d array for holding input data on
613  !! native horizontal grid and extended grid
614  !! with poles.
615  real, dimension(:,:,:), allocatable :: data_in !< A buffer for storing the full 3-d time-interpolated array
616  !! on the original grid
617  real, dimension(:,:), allocatable :: mask_in !< A 2-d mask for extended input grid.
618 
619  real :: PI_180
620  integer :: rcode, ncid, varid, ndims, id, jd, kd, jdp
621  integer :: i,j,k
622  integer, dimension(4) :: start, count, dims, dim_id
623  real, dimension(:,:), allocatable :: x_in, y_in
624  real, dimension(:), allocatable :: lon_in, lat_in
625  real, dimension(:), allocatable :: lat_inp, last_row
626  real :: max_lat, min_lat, pole, max_depth, npole
627  real :: roundoff ! The magnitude of roundoff, usually ~2e-16.
628  logical :: add_np
629  character(len=8) :: laynum
630  type(horiz_interp_type) :: Interp
631  type(axistype), dimension(4) :: axes_data
632  integer :: is, ie, js, je ! compute domain indices
633  integer :: isc,iec,jsc,jec ! global compute domain indices
634  integer :: isg, ieg, jsg, jeg ! global extent
635  integer :: isd, ied, jsd, jed ! data domain indices
636  integer :: id_clock_read
637  integer, dimension(4) :: fld_sz
638  character(len=12) :: dim_name(4)
639  logical :: debug=.false.
640  logical :: spongeDataOngrid
641  real :: npoints,varAvg
642  real, dimension(SZI_(G),SZJ_(G)) :: lon_out, lat_out, tr_out, mask_out
643  real, dimension(SZI_(G),SZJ_(G)) :: good, fill
644  real, dimension(SZI_(G),SZJ_(G)) :: tr_outf,tr_prev
645  real, dimension(SZI_(G),SZJ_(G)) :: good2,fill2
646  real, dimension(SZI_(G),SZJ_(G)) :: nlevs
647 
648  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
649  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
650  isg = g%isg ; ieg = g%ieg ; jsg = g%jsg ; jeg = g%jeg
651 
652  id_clock_read = cpu_clock_id('(Initialize tracer from Z) read', grain=clock_loop)
653 
654  pi_180=atan(1.0)/45.
655 
656  ! Open NetCDF file and if present, extract data and spatial coordinate information
657  ! The convention adopted here requires that the data be written in (i,j,k) ordering.
658 
659  call cpu_clock_begin(id_clock_read)
660 
661  fld_sz = get_external_field_size(fms_id)
662  if (allocated(lon_in)) deallocate(lon_in)
663  if (allocated(lat_in)) deallocate(lat_in)
664  if (allocated(z_in)) deallocate(z_in)
665  if (allocated(z_edges_in)) deallocate(z_edges_in)
666  if (allocated(tr_z)) deallocate(tr_z)
667  if (allocated(mask_z)) deallocate(mask_z)
668 
669  axes_data = get_external_field_axes(fms_id)
670 
671  id = fld_sz(1) ; jd = fld_sz(2) ; kd = fld_sz(3)
672 
673  spongedataongrid=.false.
674  if (PRESENT(spongeongrid)) spongedataongrid=spongeongrid
675  if (.not. spongedataongrid) then
676  allocate(lon_in(id),lat_in(jd))
677  call mpp_get_axis_data(axes_data(1), lon_in)
678  call mpp_get_axis_data(axes_data(2), lat_in)
679  endif
680 
681  allocate(z_in(kd),z_edges_in(kd+1))
682 
683  allocate(tr_z(isd:ied,jsd:jed,kd), mask_z(isd:ied,jsd:jed,kd))
684 
685  call mpp_get_axis_data(axes_data(3), z_in)
686 
687  if (present(m_to_z)) then ; do k=1,kd ; z_in(k) = m_to_z * z_in(k) ; enddo ; endif
688 
689  call cpu_clock_end(id_clock_read)
690 
691  missing_value = get_external_field_missing(fms_id)
692 
693  if (.not. spongedataongrid) then
694  ! extrapolate the input data to the north pole using the northerm-most latitude
695  max_lat = maxval(lat_in)
696  add_np=.false.
697  if (max_lat < 90.0) then
698  add_np = .true.
699  jdp = jd+1
700  allocate(lat_inp(jdp))
701  lat_inp(1:jd) = lat_in(:)
702  lat_inp(jd+1) = 90.0
703  deallocate(lat_in)
704  allocate(lat_in(1:jdp))
705  lat_in(:) = lat_inp(:)
706  else
707  jdp=jd
708  endif
709  call horiz_interp_init()
710  lon_in = lon_in*pi_180
711  lat_in = lat_in*pi_180
712  allocate(x_in(id,jdp), y_in(id,jdp))
713  call meshgrid(lon_in, lat_in, x_in, y_in)
714  lon_out(:,:) = g%geoLonT(:,:)*pi_180
715  lat_out(:,:) = g%geoLatT(:,:)*pi_180
716  allocate(data_in(id,jd,kd)) ; data_in(:,:,:)=0.0
717  allocate(tr_in(id,jd)) ; tr_in(:,:)=0.0
718  allocate(tr_inp(id,jdp)) ; tr_inp(:,:)=0.0
719  allocate(mask_in(id,jdp)) ; mask_in(:,:)=0.0
720  allocate(last_row(id)) ; last_row(:)=0.0
721  else
722  allocate(data_in(isd:ied,jsd:jed,kd))
723  endif
724  ! construct level cell boundaries as the mid-point between adjacent centers
725  z_edges_in(1) = 0.0
726  do k=2,kd
727  z_edges_in(k) = 0.5*(z_in(k-1)+z_in(k))
728  enddo
729  z_edges_in(kd+1) = 2.0*z_in(kd) - z_in(kd-1)
730 
731 
732  max_depth = maxval(g%bathyT)
733  call mpp_max(max_depth)
734 
735  if (z_edges_in(kd+1)<max_depth) z_edges_in(kd+1)=max_depth
736 
737  ! roundoff = 3.0*EPSILON(missing_value)
738  roundoff = 1.e-4
739 
740  if (.not.spongedataongrid) then
741  if (is_root_pe()) &
742  call time_interp_external(fms_id, time, data_in, verbose=.true.)
743  ! loop through each data level and interpolate to model grid.
744  ! after interpolating, fill in points which will be needed
745  ! to define the layers
746  do k=1,kd
747  write(laynum,'(I8)') k ; laynum = adjustl(laynum)
748  if (is_root_pe()) then
749  tr_in(1:id,1:jd) = data_in(1:id,1:jd,k)
750  if (add_np) then
751  last_row(:)=tr_in(:,jd); pole=0.0;npole=0.0
752  do i=1,id
753  if (abs(tr_in(i,jd)-missing_value) > abs(roundoff*missing_value)) then
754  pole = pole+last_row(i)
755  npole = npole+1.0
756  endif
757  enddo
758  if (npole > 0) then
759  pole=pole/npole
760  else
761  pole=missing_value
762  endif
763  tr_inp(:,1:jd) = tr_in(:,:)
764  tr_inp(:,jdp) = pole
765  else
766  tr_inp(:,:) = tr_in(:,:)
767  endif
768  endif
769 
770  call mpp_sync()
771  call mpp_broadcast(tr_inp, id*jdp, root_pe())
772  call mpp_sync_self()
773 
774  mask_in=0.0
775 
776  do j=1,jdp ; do i=1,id
777  if (abs(tr_inp(i,j)-missing_value) > abs(roundoff*missing_value)) then
778  mask_in(i,j)=1.0
779  tr_inp(i,j) = tr_inp(i,j) * conversion
780  else
781  tr_inp(i,j) = missing_value
782  endif
783  enddo ; enddo
784 
785  ! call fms routine horiz_interp to interpolate input level data to model horizontal grid
786  if (k == 1) then
787  call horiz_interp_new(interp, x_in, y_in, lon_out(is:ie,js:je), lat_out(is:ie,js:je), &
788  interp_method='bilinear', src_modulo=.true.)
789  endif
790 
791  if (debug) then
792  call mystats(tr_in, missing_value, 1, id, 1, jd, k, 'Tracer from file')
793  endif
794 
795  tr_out(:,:) = 0.0
796 
797  call horiz_interp(interp, tr_inp, tr_out(is:ie,js:je), missing_value=missing_value, &
798  new_missing_handle=.true.)
799 
800  mask_out(:,:) = 1.0
801  do j=js,je ; do i=is,ie
802  if (abs(tr_out(i,j)-missing_value) < abs(roundoff*missing_value)) mask_out(i,j) = 0.
803  enddo ; enddo
804 
805  fill(:,:) = 0.0 ; good(:,:) = 0.0
806 
807  npoints = 0 ; varavg = 0.
808  do j=js,je ; do i=is,ie
809  if (mask_out(i,j) < 1.0) then
810  tr_out(i,j) = missing_value
811  else
812  good(i,j) = 1.0
813  npoints = npoints + 1
814  varavg = varavg + tr_out(i,j)
815  endif
816  if ((g%mask2dT(i,j) == 1.0) .and. (z_edges_in(k) <= g%bathyT(i,j)) .and. &
817  (mask_out(i,j) < 1.0)) &
818  fill(i,j)=1.0
819  enddo ; enddo
820  call pass_var(fill, g%Domain)
821  call pass_var(good, g%Domain)
822 
823  if (debug) then
824  call mystats(tr_out, missing_value, is, ie, js, je, k, 'variable from horiz_interp()')
825  endif
826 
827  ! Horizontally homogenize data to produce perfectly "flat" initial conditions
828  if (PRESENT(homogenize)) then ; if (homogenize) then
829  call sum_across_pes(npoints)
830  call sum_across_pes(varavg)
831  if (npoints>0) then
832  varavg = varavg/real(npoints)
833  endif
834  tr_out(:,:) = varavg
835  endif ; endif
836 
837  ! tr_out contains input z-space data on the model grid with missing values
838  ! now fill in missing values using "ICE-nine" algorithm.
839  tr_outf(:,:) = tr_out(:,:)
840  if (k==1) tr_prev(:,:) = tr_outf(:,:)
841  good2(:,:) = good(:,:)
842  fill2(:,:) = fill(:,:)
843 
844  call fill_miss_2d(tr_outf, good2, fill2, tr_prev, g, smooth=.true.)
845 
846 ! if (debug) then
847 ! call hchksum(tr_outf, 'field from fill_miss_2d ', G%HI)
848 ! endif
849 
850 ! call myStats(tr_outf, missing_value, is, ie, js, je, k, 'field from fill_miss_2d()')
851 
852  tr_z(:,:,k) = tr_outf(:,:)*g%mask2dT(:,:)
853  mask_z(:,:,k) = good2(:,:) + fill2(:,:)
854  tr_prev(:,:) = tr_z(:,:,k)
855 
856  if (debug) then
857  call hchksum(tr_prev,'field after fill ',g%HI)
858  endif
859 
860  enddo ! kd
861  else
862  call time_interp_external(fms_id, time, data_in, verbose=.true.)
863  do k=1,kd
864  do j=js,je
865  do i=is,ie
866  tr_z(i,j,k)=data_in(i,j,k)
867  if (abs(tr_z(i,j,k)-missing_value) < abs(roundoff*missing_value)) mask_z(i,j,k) = 0.
868  enddo
869  enddo
870  enddo
871  endif
872 
874 
875 !> Create a 2d-mesh of grid coordinates from 1-d arrays.
876 subroutine meshgrid(x, y, x_T, y_T)
877  real, dimension(:), intent(in) :: x !< input 1-dimensional vector
878  real, dimension(:), intent(in) :: y !< input 1-dimensional vector
879  real, dimension(size(x,1),size(y,1)), intent(inout) :: x_T !< output 2-dimensional array
880  real, dimension(size(x,1),size(y,1)), intent(inout) :: y_T !< output 2-dimensional array
881 
882  integer :: ni,nj,i,j
883 
884  ni=size(x,1) ; nj=size(y,1)
885 
886  do j=1,nj ; do i=1,ni
887  x_t(i,j) = x(i)
888  enddo ; enddo
889 
890  do j=1,nj ; do i=1,ni
891  y_t(i,j) = y(j)
892  enddo ; enddo
893 
894 end subroutine meshgrid
895 
896 ! None of the subsequent code appears to be used at all.
897 
898 !> Fill grid edges for integer data
899 function fill_boundaries_int(m,cyclic_x,tripolar_n) result(mp)
900  integer, dimension(:,:), intent(in) :: m !< input array (ND)
901  logical, intent(in) :: cyclic_x !< True if domain is zonally re-entrant
902  logical, intent(in) :: tripolar_n !< True if domain has an Arctic fold
903  integer, dimension(0:size(m,1)+1,0:size(m,2)+1) :: mp
904 
905  real, dimension(size(m,1),size(m,2)) :: m_real
906  real, dimension(0:size(m,1)+1,0:size(m,2)+1) :: mp_real
907 
908  m_real = real(m)
909 
910  mp_real = fill_boundaries_real(m_real,cyclic_x,tripolar_n)
911 
912  mp = int(mp_real)
913 
914 end function fill_boundaries_int
915 
916 !> Fill grid edges for real data
917 function fill_boundaries_real(m,cyclic_x,tripolar_n) result(mp)
918  real, dimension(:,:), intent(in) :: m !< input array (ND)
919  logical, intent(in) :: cyclic_x !< True if domain is zonally re-entrant
920  logical, intent(in) :: tripolar_n !< True if domain has an Arctic fold
921  real, dimension(0:size(m,1)+1,0:size(m,2)+1) :: mp
922 
923  integer :: ni,nj,i,j
924 
925  ni=size(m,1); nj=size(m,2)
926 
927  mp(1:ni,1:nj)=m(:,:)
928 
929  if (cyclic_x) then
930  mp(0,1:nj)=m(ni,1:nj)
931  mp(ni+1,1:nj)=m(1,1:nj)
932  else
933  mp(0,1:nj)=m(1,1:nj)
934  mp(ni+1,1:nj)=m(ni,1:nj)
935  endif
936 
937  mp(1:ni,0)=m(1:ni,1)
938  if (tripolar_n) then
939  do i=1,ni
940  mp(i,nj+1)=m(ni-i+1,nj)
941  enddo
942  else
943  mp(1:ni,nj+1)=m(1:ni,nj)
944  endif
945 
946 end function fill_boundaries_real
947 
948 !> Solve del2 (zi) = 0 using successive iterations
949 !! with a 5 point stencil. Only points fill==1 are
950 !! modified. Except where bad==1, information propagates
951 !! isotropically in index space. The resulting solution
952 !! in each region is an approximation to del2(zi)=0 subject to
953 !! boundary conditions along the valid points curve bounding this region.
954 subroutine smooth_heights(zi,fill,bad,sor,niter,cyclic_x, tripolar_n)
955  real, dimension(:,:), intent(inout) :: zi !< input and output array (ND)
956  integer, dimension(size(zi,1),size(zi,2)), intent(in) :: fill !< same shape as zi, 1=fill
957  integer, dimension(size(zi,1),size(zi,2)), intent(in) :: bad !< same shape as zi, 1=bad data
958  real, intent(in) :: sor !< relaxation coefficient (ND)
959  integer, intent(in) :: niter !< maximum number of iterations
960  logical, intent(in) :: cyclic_x !< true if domain is zonally reentrant
961  logical, intent(in) :: tripolar_n !< true if domain has an Arctic fold
962 
963  ! Local variables
964  real, dimension(size(zi,1),size(zi,2)) :: res, m
965  integer, dimension(size(zi,1),size(zi,2),4) :: B
966  real, dimension(0:size(zi,1)+1,0:size(zi,2)+1) :: mp
967  integer, dimension(0:size(zi,1)+1,0:size(zi,2)+1) :: nm
968  integer :: i,j,k,n
969  integer :: ni,nj
970  real :: Isum, bsum
971 
972  ni=size(zi,1) ; nj=size(zi,2)
973 
974 
975  mp(:,:) = fill_boundaries(zi,cyclic_x,tripolar_n)
976 
977  b(:,:,:) = 0.0
978  nm(:,:) = fill_boundaries(bad,cyclic_x,tripolar_n)
979 
980  do j=1,nj
981  do i=1,ni
982  if (fill(i,j) == 1) then
983  b(i,j,1)=1-nm(i+1,j);b(i,j,2)=1-nm(i-1,j)
984  b(i,j,3)=1-nm(i,j+1);b(i,j,4)=1-nm(i,j-1)
985  endif
986  enddo
987  enddo
988 
989  do n=1,niter
990  do j=1,nj
991  do i=1,ni
992  if (fill(i,j) == 1) then
993  bsum = real(b(i,j,1)+b(i,j,2)+b(i,j,3)+b(i,j,4))
994  isum = 1.0/bsum
995  res(i,j)=isum*(b(i,j,1)*mp(i+1,j)+b(i,j,2)*mp(i-1,j)+&
996  b(i,j,3)*mp(i,j+1)+b(i,j,4)*mp(i,j-1)) - mp(i,j)
997  endif
998  enddo
999  enddo
1000  res(:,:)=res(:,:)*sor
1001 
1002  do j=1,nj
1003  do i=1,ni
1004  mp(i,j)=mp(i,j)+res(i,j)
1005  enddo
1006  enddo
1007 
1008  zi(:,:)=mp(1:ni,1:nj)
1009  mp = fill_boundaries(zi,cyclic_x,tripolar_n)
1010  enddo
1011 
1012 end subroutine smooth_heights
1013 
1014 end module mom_horizontal_regridding
mom_time_manager
Wraps the FMS time manager functions.
Definition: MOM_time_manager.F90:2
mom_horizontal_regridding::fill_boundaries_int
integer function, dimension(0:size(m, 1)+1, 0:size(m, 2)+1) fill_boundaries_int(m, cyclic_x, tripolar_n)
Fill grid edges for integer data.
Definition: MOM_horizontal_regridding.F90:900
mom_horizontal_regridding::horiz_interp_and_extrap_tracer
Extrapolate and interpolate data.
Definition: MOM_horizontal_regridding.F90:52
mom_file_parser::log_version
An overloaded interface to log version information about modules.
Definition: MOM_file_parser.F90:109
mom_horizontal_regridding::horiz_interp_and_extrap_tracer_fms_id
subroutine horiz_interp_and_extrap_tracer_fms_id(fms_id, Time, conversion, G, tr_z, mask_z, z_in, z_edges_in, missing_value, reentrant_x, tripolar_n, homogenize, spongeOngrid, m_to_Z)
Extrapolate and interpolate using a FMS time interpolation handle.
Definition: MOM_horizontal_regridding.F90:591
mom_variables::thermo_var_ptrs
Pointers to an assortment of thermodynamic fields that may be available, including potential temperat...
Definition: MOM_variables.F90:78
mom_get_input::directories
Container for paths and parameter file names.
Definition: MOM_get_input.F90:20
mom_string_functions
Handy functions for manipulating strings.
Definition: MOM_string_functions.F90:2
mom_horizontal_regridding::fill_boundaries
Fill grid edges.
Definition: MOM_horizontal_regridding.F90:46
mom_error_handler::mom_mesg
subroutine, public mom_mesg(message, verb, all_print)
This provides a convenient interface for writing an informative comment.
Definition: MOM_error_handler.F90:53
mom_domains::pass_var
Do a halo update on an array.
Definition: MOM_domains.F90:49
mom_file_parser::param_file_type
A structure that can be parsed to read and document run-time parameters.
Definition: MOM_file_parser.F90:54
mom_file_parser::get_param
An overloaded interface to read and log the values of various types of parameters.
Definition: MOM_file_parser.F90:102
mom_horizontal_regridding
Horizontal interpolation.
Definition: MOM_horizontal_regridding.F90:2
mom_io
This module contains I/O framework code.
Definition: MOM_io.F90:2
mom_horizontal_regridding::fill_miss_2d
subroutine fill_miss_2d(aout, good, fill, prev, G, smooth, num_pass, relc, crit, keep_bug, debug)
Use ICE-9 algorithm to populate points (fill=1) with valid data (good=1). If no information is availa...
Definition: MOM_horizontal_regridding.F90:104
mom_string_functions::uppercase
character(len=len(input_string)) function, public uppercase(input_string)
Return a string in which all uppercase letters have been replaced by their lowercase counterparts.
Definition: MOM_string_functions.F90:42
mom_get_input
Reads the only Fortran name list needed to boot-strap the model.
Definition: MOM_get_input.F90:6
mom_io::read_axis_data
subroutine, public read_axis_data(filename, axis_name, var)
Read the data associated with a named axis in a file.
Definition: MOM_io.F90:438
mom_coms
Interfaces to non-domain-oriented communication subroutines, including the MOM6 reproducing sums faci...
Definition: MOM_coms.F90:3
mom_domains::pass_vector
Do a halo update on a pair of arrays representing the two components of a vector.
Definition: MOM_domains.F90:54
mom_domains::to_all
integer, parameter, public to_all
A flag for passing in all directions.
Definition: MOM_domains.F90:132
mom_domains
Describes the decomposed MOM domain and has routines for communications across PEs.
Definition: MOM_domains.F90:2
mom_variables
Provides transparent structures with groups of MOM6 variables and supporting routines.
Definition: MOM_variables.F90:2
mom_cpu_clock
Wraps the MPP cpu clock functions.
Definition: MOM_cpu_clock.F90:2
mom_file_parser
The MOM6 facility to parse input files for runtime parameters.
Definition: MOM_file_parser.F90:2
mom_error_handler::calltree_leave
subroutine, public calltree_leave(mesg)
Writes a message about leaving a subroutine if call tree reporting is active.
Definition: MOM_error_handler.F90:151
mom_horizontal_regridding::meshgrid
subroutine meshgrid(x, y, x_T, y_T)
Create a 2d-mesh of grid coordinates from 1-d arrays.
Definition: MOM_horizontal_regridding.F90:877
mom_horizontal_regridding::horiz_interp_and_extrap_tracer_record
subroutine horiz_interp_and_extrap_tracer_record(filename, varnam, conversion, recnum, G, tr_z, mask_z, z_in, z_edges_in, missing_value, reentrant_x, tripolar_n, homogenize, m_to_Z)
Extrapolate and interpolate from a file record.
Definition: MOM_horizontal_regridding.F90:268
mom_grid
Provides the ocean grid type.
Definition: MOM_grid.F90:2
mom_error_handler::is_root_pe
logical function, public is_root_pe()
This returns .true. if the current PE is the root PE.
Definition: MOM_error_handler.F90:44
mom_horizontal_regridding::mystats
subroutine, public mystats(array, missing, is, ie, js, je, k, mesg)
Write to the terminal some basic statistics about the k-th level of an array.
Definition: MOM_horizontal_regridding.F90:61
mom_debugging
Provides checksumming functions for debugging.
Definition: MOM_debugging.F90:7
mom_io::vardesc
Type for describing a variable, typically a tracer.
Definition: MOM_io.F90:53
mom_io::file_exists
Indicate whether a file exists, perhaps with domain decomposition.
Definition: MOM_io.F90:68
mom_error_handler::mom_error
subroutine, public mom_error(level, message, all_print)
This provides a convenient interface for writing an mpp_error message with run-time filter based on a...
Definition: MOM_error_handler.F90:72
mom_error_handler::calltree_waypoint
subroutine, public calltree_waypoint(mesg, n)
Writes a message about reaching a milestone if call tree reporting is active.
Definition: MOM_error_handler.F90:161
mom_file_parser::log_param
An overloaded interface to log the values of various types of parameters.
Definition: MOM_file_parser.F90:96
mom_horizontal_regridding::fill_boundaries_real
real function, dimension(0:size(m, 1)+1, 0:size(m, 2)+1) fill_boundaries_real(m, cyclic_x, tripolar_n)
Fill grid edges for real data.
Definition: MOM_horizontal_regridding.F90:918
mom_error_handler
Routines for error handling and I/O management.
Definition: MOM_error_handler.F90:2
mom_error_handler::calltree_enter
subroutine, public calltree_enter(mesg, n)
Writes a message about entering a subroutine if call tree reporting is active.
Definition: MOM_error_handler.F90:130
mom_grid::ispointincell
logical function, public ispointincell(G, i, j, x, y)
Returns true if the coordinates (x,y) are within the h-cell (i,j)
Definition: MOM_grid.F90:467
mom_grid::ocean_grid_type
Ocean grid type. See mom_grid for details.
Definition: MOM_grid.F90:26
mom_horizontal_regridding::smooth_heights
subroutine smooth_heights(zi, fill, bad, sor, niter, cyclic_x, tripolar_n)
Solve del2 (zi) = 0 using successive iterations with a 5 point stencil. Only points fill==1 are modif...
Definition: MOM_horizontal_regridding.F90:955
mom_file_parser::read_param
An overloaded interface to read various types of parameters.
Definition: MOM_file_parser.F90:90