MOM6
MOM_diag_remap.F90
Go to the documentation of this file.
1 !> provides runtime remapping of diagnostics to z star, sigma and
2 !! rho vertical coordinates.
3 !!
4 !! The diag_remap_ctrl type represents a remapping of diagnostics to a particular
5 !! vertical coordinate. The module is used by the diag mediator module in the
6 !! following way:
7 !! 1. diag_remap_init() is called to initialize a diag_remap_ctrl instance.
8 !! 2. diag_remap_configure_axes() is called to read the configuration file and set up the
9 !! vertical coordinate / axes definitions.
10 !! 3. diag_remap_get_axes_info() returns information needed for the diag mediator to
11 !! define new axes for the remapped diagnostics.
12 !! 4. diag_remap_update() is called periodically (whenever h, T or S change) to either
13 !! create or update the target remapping grids.
14 !! 5. diag_remap_do_remap() is called from within a diag post() to do the remapping before
15 !! the diagnostic is written out.
16 
17 
18 ! NOTE: In the following functions, the fields are passed using 1-based
19 ! indexing, which requires special handling within the grid index loops.
20 !
21 ! * diag_remap_do_remap
22 ! * vertically_reintegrate_diag_field
23 ! * vertically_interpolate_diag_field
24 ! * horizontally_average_diag_field
25 !
26 ! Symmetric grids add an additional row of western and southern points to u-
27 ! and v-grids. Non-symmetric grids are 1-based and symmetric grids are
28 ! zero-based, allowing the same expressions to be used when accessing the
29 ! fields. But if u- or v-points become 1-indexed, as in these functions, then
30 ! the stencils must be re-assessed.
31 !
32 ! For interpolation between h and u grids, we use the following relations:
33 !
34 ! h->u: f_u[ig] = 0.5 * (f_h[ ig ] + f_h[ig+1])
35 ! f_u[i1] = 0.5 * (f_h[i1-1] + f_h[ i1 ])
36 !
37 ! u->h: f_h[ig] = 0.5 * (f_u[ig-1] + f_u[ ig ])
38 ! f_h[i1] = 0.5 * (f_u[ i1 ] + f_u[i1+1])
39 !
40 ! where ig is the grid index and i1 is the 1-based index. That is, a 1-based
41 ! u-point is ahead of its matching h-point in non-symmetric mode, but behind
42 ! its matching h-point in non-symmetric mode.
43 !
44 ! We can combine these expressions by applying to ig a -1 shift on u-grids and
45 ! a +1 shift on h-grids in symmetric mode.
46 !
47 ! We do not adjust the h-point indices, since they are assumed to be 1-based.
48 ! This is only correct when global indexing is disabled. If global indexing is
49 ! enabled, then all indices will need to be defined relative to the data
50 ! domain.
51 !
52 ! Finally, note that the mask input fields are pointers to arrays which are
53 ! zero-indexed, and do not need any corrections over grid index loops.
54 
55 
57 
58 ! This file is part of MOM6. See LICENSE.md for the license.
59 
60 use mom_coms, only : reproducing_sum
61 use mom_error_handler, only : mom_error, fatal, assert, warning
64 use mom_io, only : slasher, mom_read_data
65 use mom_io, only : file_exists, field_size
67 use mom_grid, only : ocean_grid_type
70 use mom_eos, only : eos_type
77 use regrid_consts, only : coordinatemode
80 use coord_rho, only : build_rho_column
81 
82 use diag_axis_mod, only : get_diag_axis_name
83 use diag_manager_mod, only : diag_axis_init
84 
86 implicit none ; private
87 
88 public diag_remap_ctrl
97 
98 !> Represents remapping of diagnostics to a particular vertical coordinate.
99 !!
100 !! There is one of these types for each vertical coordinate. The vertical axes
101 !! of a diagnostic will reference an instance of this type indicating how (or
102 !! if) the diagnostic should be vertically remapped when being posted.
104  logical :: configured = .false. !< Whether vertical coordinate has been configured
105  logical :: initialized = .false. !< Whether remappping initialized
106  logical :: used = .false. !< Whether this coordinate actually gets used.
107  integer :: vertical_coord = 0 !< The vertical coordinate that we remap to
108  character(len=10) :: vertical_coord_name ='' !< The coordinate name as understood by ALE
109  character(len=16) :: diag_coord_name = '' !< A name for the purpose of run-time parameters
110  character(len=8) :: diag_module_suffix = '' !< The suffix for the module to appear in diag_table
111  type(remapping_cs) :: remap_cs !< Remapping control structure use for this axes
112  type(regridding_cs) :: regrid_cs !< Regridding control structure that defines the coordiantes for this axes
113  integer :: nz = 0 !< Number of vertical levels used for remapping
114  real, dimension(:,:,:), allocatable :: h !< Remap grid thicknesses
115  real, dimension(:), allocatable :: dz !< Nominal layer thicknesses
116  integer :: interface_axes_id = 0 !< Vertical axes id for remapping at interfaces
117  integer :: layer_axes_id = 0 !< Vertical axes id for remapping on layers
118 end type diag_remap_ctrl
119 
120 contains
121 
122 !> Initialize a diagnostic remapping type with the given vertical coordinate.
123 subroutine diag_remap_init(remap_cs, coord_tuple)
124  type(diag_remap_ctrl), intent(inout) :: remap_cs !< Diag remapping control structure
125  character(len=*), intent(in) :: coord_tuple !< A string in form of
126  !! MODULE_SUFFIX PARAMETER_SUFFIX COORDINATE_NAME
127 
128  remap_cs%diag_module_suffix = trim(extractword(coord_tuple, 1))
129  remap_cs%diag_coord_name = trim(extractword(coord_tuple, 2))
130  remap_cs%vertical_coord_name = trim(extractword(coord_tuple, 3))
131  remap_cs%vertical_coord = coordinatemode(remap_cs%vertical_coord_name)
132  remap_cs%configured = .false.
133  remap_cs%initialized = .false.
134  remap_cs%used = .false.
135  remap_cs%nz = 0
136 
137 end subroutine diag_remap_init
138 
139 !> De-init a diagnostic remapping type.
140 !! Free allocated memory.
141 subroutine diag_remap_end(remap_cs)
142  type(diag_remap_ctrl), intent(inout) :: remap_cs !< Diag remapping control structure
143 
144  if (allocated(remap_cs%h)) deallocate(remap_cs%h)
145  if (allocated(remap_cs%dz)) deallocate(remap_cs%dz)
146  remap_cs%configured = .false.
147  remap_cs%initialized = .false.
148  remap_cs%used = .false.
149  remap_cs%nz = 0
150 
151 end subroutine diag_remap_end
152 
153 !> Inform that all diagnostics have been registered.
154 !! If _set_active() has not been called on the remapping control structure
155 !! will be disabled. This saves time in the case that a vertical coordinate was
156 !! configured but no diagnostics which use the coordinate appeared in the
157 !! diag_table.
158 subroutine diag_remap_diag_registration_closed(remap_cs)
159  type(diag_remap_ctrl), intent(inout) :: remap_cs !< Diag remapping control structure
160 
161  if (.not. remap_cs%used) then
162  call diag_remap_end(remap_cs)
163  endif
164 
166 
167 !> Indicate that this remapping type is actually used by the diag manager.
168 !! If this is never called then the type will be disabled to save time.
169 !! See further explanation with diag_remap_registration_closed.
170 subroutine diag_remap_set_active(remap_cs)
171  type(diag_remap_ctrl), intent(inout) :: remap_cs !< Diag remapping control structure
172 
173  remap_cs%used = .true.
174 
175 end subroutine diag_remap_set_active
176 
177 !> Configure the vertical axes for a diagnostic remapping control structure.
178 !! Reads a configuration parameters to determine coordinate generation.
179 subroutine diag_remap_configure_axes(remap_cs, GV, US, param_file)
180  type(diag_remap_ctrl), intent(inout) :: remap_cs !< Diag remap control structure
181  type(verticalgrid_type), intent(in) :: gv !< ocean vertical grid structure
182  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
183  type(param_file_type), intent(in) :: param_file !< Parameter file structure
184  ! Local variables
185  integer :: nzi(4), nzl(4), k
186  character(len=200) :: inputdir, string, filename, int_varname, layer_varname
187  character(len=40) :: mod = "MOM_diag_remap" ! This module's name.
188  character(len=8) :: units, expected_units
189  character(len=34) :: longname, string2
190 
191  character(len=256) :: err_msg
192  logical :: ierr
193 
194  real, allocatable, dimension(:) :: interfaces, layers
195 
196  call initialize_regridding(remap_cs%regrid_cs, gv, us, gv%max_depth, param_file, mod, &
197  trim(remap_cs%vertical_coord_name), "DIAG_COORD", trim(remap_cs%diag_coord_name))
198  call set_regrid_params(remap_cs%regrid_cs, min_thickness=0., integrate_downward_for_e=.false.)
199 
200  remap_cs%nz = get_regrid_size(remap_cs%regrid_cs)
201 
202  if (remap_cs%vertical_coord == coordinatemode('SIGMA')) then
203  units = 'nondim'
204  longname = 'Fraction'
205  elseif (remap_cs%vertical_coord == coordinatemode('RHO')) then
206  units = 'kg m-3'
207  longname = 'Target Potential Density'
208  else
209  units = 'meters'
210  longname = 'Depth'
211  endif
212 
213  ! Make axes objects
214  allocate(interfaces(remap_cs%nz+1))
215  allocate(layers(remap_cs%nz))
216 
217  interfaces(:) = getcoordinateinterfaces(remap_cs%regrid_cs)
218  layers(:) = 0.5 * ( interfaces(1:remap_cs%nz) + interfaces(2:remap_cs%nz+1) )
219 
220  remap_cs%interface_axes_id = diag_axis_init(lowercase(trim(remap_cs%diag_coord_name))//'_i', &
221  interfaces, trim(units), 'z', &
222  trim(longname)//' at interface', direction=-1)
223  remap_cs%layer_axes_id = diag_axis_init(lowercase(trim(remap_cs%diag_coord_name))//'_l', &
224  layers, trim(units), 'z', &
225  trim(longname)//' at cell center', direction=-1, &
226  edges=remap_cs%interface_axes_id)
227 
228  ! Axes have now been configured.
229  remap_cs%configured = .true.
230 
231  deallocate(interfaces)
232  deallocate(layers)
233 
234 end subroutine diag_remap_configure_axes
235 
236 !> Get layer and interface axes ids for this coordinate
237 !! Needed when defining axes groups.
238 subroutine diag_remap_get_axes_info(remap_cs, nz, id_layer, id_interface)
239  type(diag_remap_ctrl), intent(in) :: remap_cs !< Diagnostic coordinate control structure
240  integer, intent(out) :: nz !< Number of vertical levels for the coordinate
241  integer, intent(out) :: id_layer !< 1D-axes id for layer points
242  integer, intent(out) :: id_interface !< 1D-axes id for interface points
243 
244  nz = remap_cs%nz
245  id_layer = remap_cs%layer_axes_id
246  id_interface = remap_cs%interface_axes_id
247 
248 end subroutine diag_remap_get_axes_info
249 
250 
251 !> Whether or not the axes for this vertical coordinated has been configured.
252 !! Configuration is complete when diag_remap_configure_axes() has been
253 !! successfully called.
254 function diag_remap_axes_configured(remap_cs)
255  type(diag_remap_ctrl), intent(in) :: remap_cs !< Diagnostic coordinate control structure
256  logical :: diag_remap_axes_configured
257 
258  diag_remap_axes_configured = remap_cs%configured
259 
260 end function
261 
262 !> Build/update target vertical grids for diagnostic remapping.
263 !! \note The target grids need to be updated whenever sea surface
264 !! height or layer thicknesses changes. In the case of density-based
265 !! coordinates then technically we should also regenerate the
266 !! target grid whenever T/S change.
267 subroutine diag_remap_update(remap_cs, G, GV, US, h, T, S, eqn_of_state)
268  type(diag_remap_ctrl), intent(inout) :: remap_cs !< Diagnostic coordinate control structure
269  type(ocean_grid_type), pointer :: g !< The ocean's grid type
270  type(verticalgrid_type), intent(in) :: gv !< ocean vertical grid structure
271  type(unit_scale_type), intent(in) :: us !< A dimensional unit scaling type
272  real, dimension(:, :, :), intent(in) :: h !< New thickness
273  real, dimension(:, :, :), intent(in) :: t !< New T
274  real, dimension(:, :, :), intent(in) :: s !< New S
275  type(eos_type), pointer :: eqn_of_state !< A pointer to the equation of state
276 
277  ! Local variables
278  real, dimension(remap_cs%nz + 1) :: zinterfaces
279  real :: h_neglect, h_neglect_edge
280  integer :: i, j, k, nz
281 
282  ! Note that coordinateMode('LAYER') is never 'configured' so will
283  ! always return here.
284  if (.not. remap_cs%configured) then
285  return
286  endif
287 
288  !### Try replacing both of these with GV%H_subroundoff
289  if (gv%Boussinesq) then
290  h_neglect = gv%m_to_H*1.0e-30 ; h_neglect_edge = gv%m_to_H*1.0e-10
291  else
292  h_neglect = gv%kg_m2_to_H*1.0e-30 ; h_neglect_edge = gv%kg_m2_to_H*1.0e-10
293  endif
294  nz = remap_cs%nz
295 
296  if (.not. remap_cs%initialized) then
297  ! Initialize remapping and regridding on the first call
298  call initialize_remapping(remap_cs%remap_cs, 'PPM_IH4', boundary_extrapolation=.false.)
299  allocate(remap_cs%h(g%isd:g%ied,g%jsd:g%jed, nz))
300  remap_cs%initialized = .true.
301  endif
302 
303  ! Calculate remapping thicknesses for different target grids based on
304  ! nominal/target interface locations. This happens for every call on the
305  ! assumption that h, T, S has changed.
306  do j=g%jsc-1, g%jec+1 ; do i=g%isc-1, g%iec+1
307  if (g%mask2dT(i,j)==0.) then
308  remap_cs%h(i,j,:) = 0.
309  cycle
310  endif
311 
312  if (remap_cs%vertical_coord == coordinatemode('ZSTAR')) then
313  call build_zstar_column(get_zlike_cs(remap_cs%regrid_cs), &
314  gv%Z_to_H*g%bathyT(i,j), sum(h(i,j,:)), &
315  zinterfaces, zscale=gv%Z_to_H)
316  elseif (remap_cs%vertical_coord == coordinatemode('SIGMA')) then
317  call build_sigma_column(get_sigma_cs(remap_cs%regrid_cs), &
318  gv%Z_to_H*g%bathyT(i,j), sum(h(i,j,:)), zinterfaces)
319  elseif (remap_cs%vertical_coord == coordinatemode('RHO')) then
320  call build_rho_column(get_rho_cs(remap_cs%regrid_cs), g%ke, &
321  us%Z_to_m*g%bathyT(i,j), h(i,j,:), t(i,j,:), s(i,j,:), &
322  eqn_of_state, zinterfaces, h_neglect, h_neglect_edge)
323  elseif (remap_cs%vertical_coord == coordinatemode('SLIGHT')) then
324 ! call build_slight_column(remap_cs%regrid_cs,remap_cs%remap_cs, nz, &
325 ! US%Z_to_m*G%bathyT(i,j), sum(h(i,j,:)), zInterfaces)
326  call mom_error(fatal,"diag_remap_update: SLIGHT coordinate not coded for diagnostics yet!")
327  elseif (remap_cs%vertical_coord == coordinatemode('HYCOM1')) then
328 ! call build_hycom1_column(remap_cs%regrid_cs, nz, &
329 ! US%Z_to_m*G%bathyT(i,j), sum(h(i,j,:)), zInterfaces)
330  call mom_error(fatal,"diag_remap_update: HYCOM1 coordinate not coded for diagnostics yet!")
331  endif
332  remap_cs%h(i,j,:) = zinterfaces(1:nz) - zinterfaces(2:nz+1)
333  enddo ; enddo
334 
335 end subroutine diag_remap_update
336 
337 !> Remap diagnostic field to alternative vertical grid.
338 subroutine diag_remap_do_remap(remap_cs, G, GV, h, staggered_in_x, staggered_in_y, &
339  mask, missing_value, field, remapped_field)
340  type(diag_remap_ctrl), intent(in) :: remap_cs !< Diagnostic coodinate control structure
341  type(ocean_grid_type), intent(in) :: g !< Ocean grid structure
342  type(verticalgrid_type), intent(in) :: gv !< ocean vertical grid structure
343  real, dimension(:,:,:), intent(in) :: h !< The current thicknesses
344  logical, intent(in) :: staggered_in_x !< True is the x-axis location is at u or q points
345  logical, intent(in) :: staggered_in_y !< True is the y-axis location is at v or q points
346  real, dimension(:,:,:), pointer :: mask !< A mask for the field
347  real, intent(in) :: missing_value !< A missing_value to assign land/vanished points
348  real, dimension(:,:,:), intent(in) :: field(:,:,:) !< The diagnostic field to be remapped
349  real, dimension(:,:,:), intent(inout) :: remapped_field !< Field remapped to new coordinate
350  ! Local variables
351  real, dimension(remap_cs%nz) :: h_dest
352  real, dimension(size(h,3)) :: h_src
353  real :: h_neglect, h_neglect_edge
354  integer :: nz_src, nz_dest
355  integer :: i, j, k !< Grid index
356  integer :: i1, j1 !< 1-based index
357  integer :: i_lo, i_hi, j_lo, j_hi !< (uv->h) interpolation indices
358  integer :: shift !< Symmetric offset for 1-based indexing
359 
360  call assert(remap_cs%initialized, 'diag_remap_do_remap: remap_cs not initialized.')
361  call assert(size(field, 3) == size(h, 3), &
362  'diag_remap_do_remap: Remap field and thickness z-axes do not match.')
363 
364  !### Try replacing both of these with GV%H_subroundoff
365  if (gv%Boussinesq) then
366  h_neglect = gv%m_to_H*1.0e-30 ; h_neglect_edge = gv%m_to_H*1.0e-10
367  else
368  h_neglect = gv%kg_m2_to_H*1.0e-30 ; h_neglect_edge = gv%kg_m2_to_H*1.0e-10
369  endif
370 
371  nz_src = size(field,3)
372  nz_dest = remap_cs%nz
373  remapped_field(:,:,:) = 0.
374 
375  ! Symmetric grid offset under 1-based indexing; see header for details.
376  shift = 0; if (g%symmetric) shift = 1
377 
378  if (staggered_in_x .and. .not. staggered_in_y) then
379  ! U-points
380  do j=g%jsc, g%jec
381  do i=g%iscB, g%iecB
382  i1 = i - g%isdB + 1
383  i_lo = i1 - shift; i_hi = i_lo + 1
384  if (associated(mask)) then
385  if (mask(i,j,1) == 0.) cycle
386  endif
387  h_src(:) = 0.5 * (h(i_lo,j,:) + h(i_hi,j,:))
388  h_dest(:) = 0.5 * (remap_cs%h(i_lo,j,:) + remap_cs%h(i_hi,j,:))
389  call remapping_core_h(remap_cs%remap_cs, &
390  nz_src, h_src(:), field(i1,j,:), &
391  nz_dest, h_dest(:), remapped_field(i1,j,:), &
392  h_neglect, h_neglect_edge)
393  enddo
394  enddo
395  elseif (staggered_in_y .and. .not. staggered_in_x) then
396  ! V-points
397  do j=g%jscB, g%jecB
398  j1 = j - g%jsdB + 1
399  j_lo = j1 - shift; j_hi = j_lo + 1
400  do i=g%isc, g%iec
401  if (associated(mask)) then
402  if (mask(i,j,1) == 0.) cycle
403  endif
404  h_src(:) = 0.5 * (h(i,j_lo,:) + h(i,j_hi,:))
405  h_dest(:) = 0.5 * (remap_cs%h(i,j_lo,:) + remap_cs%h(i,j_hi,:))
406  call remapping_core_h(remap_cs%remap_cs, &
407  nz_src, h_src(:), field(i,j1,:), &
408  nz_dest, h_dest(:), remapped_field(i,j1,:), &
409  h_neglect, h_neglect_edge)
410  enddo
411  enddo
412  elseif ((.not. staggered_in_x) .and. (.not. staggered_in_y)) then
413  ! H-points
414  do j=g%jsc, g%jec
415  do i=g%isc, g%iec
416  if (associated(mask)) then
417  if (mask(i,j,1) == 0.) cycle
418  endif
419  h_src(:) = h(i,j,:)
420  h_dest(:) = remap_cs%h(i,j,:)
421  call remapping_core_h(remap_cs%remap_cs, &
422  nz_src, h_src(:), field(i,j,:), &
423  nz_dest, h_dest(:), remapped_field(i,j,:), &
424  h_neglect, h_neglect_edge)
425  enddo
426  enddo
427  else
428  call assert(.false., 'diag_remap_do_remap: Unsupported axis combination')
429  endif
430 
431 end subroutine diag_remap_do_remap
432 
433 !> Calculate masks for target grid
434 subroutine diag_remap_calc_hmask(remap_cs, G, mask)
435  type(diag_remap_ctrl), intent(in) :: remap_cs !< Diagnostic coodinate control structure
436  type(ocean_grid_type), intent(in) :: g !< Ocean grid structure
437  real, dimension(:,:,:), intent(out) :: mask !< h-point mask for target grid
438  ! Local variables
439  real, dimension(remap_cs%nz) :: h_dest
440  integer :: i, j, k
441  logical :: mask_vanished_layers
442  real :: h_tot, h_err
443 
444  call assert(remap_cs%initialized, 'diag_remap_calc_hmask: remap_cs not initialized.')
445 
446  ! Only z*-like diagnostic coordinates should have a 3d mask
447  mask_vanished_layers = (remap_cs%vertical_coord == coordinatemode('ZSTAR'))
448  mask(:,:,:) = 0.
449 
450  do j=g%jsc-1, g%jec+1 ; do i=g%isc-1, g%iec+1
451  if (g%mask2dT(i,j)>0.) then
452  if (mask_vanished_layers) then
453  h_dest(:) = remap_cs%h(i,j,:)
454  h_tot = 0.
455  h_err = 0.
456  do k=1, remap_cs%nz
457  h_tot = h_tot + h_dest(k)
458  ! This is an overestimate of how thick a vanished layer might be, that
459  ! appears due to round-off.
460  h_err = h_err + epsilon(h_tot) * h_tot
461  ! Mask out vanished layers
462  if (h_dest(k)<=8.*h_err) then
463  mask(i,j,k) = 0.
464  else
465  mask(i,j,k) = 1.
466  endif
467  enddo
468  else ! all layers might contain data
469  mask(i,j,:) = 1.
470  endif
471  endif
472  enddo ; enddo
473 
474 end subroutine diag_remap_calc_hmask
475 
476 !> Vertically re-grid an already vertically-integrated diagnostic field to alternative vertical grid.
477 subroutine vertically_reintegrate_diag_field(remap_cs, G, h, staggered_in_x, staggered_in_y, &
478  mask, missing_value, field, reintegrated_field)
479  type(diag_remap_ctrl), intent(in) :: remap_cs !< Diagnostic coodinate control structure
480  type(ocean_grid_type), intent(in) :: g !< Ocean grid structure
481  real, dimension(:,:,:), intent(in) :: h !< The current thicknesses
482  logical, intent(in) :: staggered_in_x !< True is the x-axis location is at u or q points
483  logical, intent(in) :: staggered_in_y !< True is the y-axis location is at v or q points
484  real, dimension(:,:,:), pointer :: mask !< A mask for the field
485  real, intent(in) :: missing_value !< A missing_value to assign land/vanished points
486  real, dimension(:,:,:), intent(in) :: field !< The diagnostic field to be remapped
487  real, dimension(:,:,:), intent(inout) :: reintegrated_field !< Field argument remapped to alternative coordinate
488  ! Local variables
489  real, dimension(remap_cs%nz) :: h_dest
490  real, dimension(size(h,3)) :: h_src
491  integer :: nz_src, nz_dest
492  integer :: i, j, k !< Grid index
493  integer :: i1, j1 !< 1-based index
494  integer :: i_lo, i_hi, j_lo, j_hi !< (uv->h) interpolation indices
495  integer :: shift !< Symmetric offset for 1-based indexing
496 
497  call assert(remap_cs%initialized, 'vertically_reintegrate_diag_field: remap_cs not initialized.')
498  call assert(size(field, 3) == size(h, 3), &
499  'vertically_reintegrate_diag_field: Remap field and thickness z-axes do not match.')
500 
501  nz_src = size(field,3)
502  nz_dest = remap_cs%nz
503  reintegrated_field(:,:,:) = 0.
504 
505  ! Symmetric grid offset under 1-based indexing; see header for details.
506  shift = 0; if (g%symmetric) shift = 1
507 
508  if (staggered_in_x .and. .not. staggered_in_y) then
509  ! U-points
510  do j=g%jsc, g%jec
511  do i=g%iscB, g%iecB
512  i1 = i - g%isdB + 1
513  i_lo = i1 - shift; i_hi = i_lo + 1
514  if (associated(mask)) then
515  if (mask(i,j,1) == 0.) cycle
516  endif
517  h_src(:) = 0.5 * (h(i_lo,j,:) + h(i_hi,j,:))
518  h_dest(:) = 0.5 * (remap_cs%h(i_lo,j,:) + remap_cs%h(i_hi,j,:))
519  call reintegrate_column(nz_src, h_src, field(i1,j,:), &
520  nz_dest, h_dest, 0., reintegrated_field(i1,j,:))
521  enddo
522  enddo
523  elseif (staggered_in_y .and. .not. staggered_in_x) then
524  ! V-points
525  do j=g%jscB, g%jecB
526  j1 = j - g%jsdB + 1
527  j_lo = j1 - shift; j_hi = j_lo + 1
528  do i=g%isc, g%iec
529  if (associated(mask)) then
530  if (mask(i,j,1) == 0.) cycle
531  endif
532  h_src(:) = 0.5 * (h(i,j_lo,:) + h(i,j_hi,:))
533  h_dest(:) = 0.5 * (remap_cs%h(i,j_lo,:) + remap_cs%h(i,j_hi,:))
534  call reintegrate_column(nz_src, h_src, field(i,j1,:), &
535  nz_dest, h_dest, 0., reintegrated_field(i,j1,:))
536  enddo
537  enddo
538  elseif ((.not. staggered_in_x) .and. (.not. staggered_in_y)) then
539  ! H-points
540  do j=g%jsc, g%jec
541  do i=g%isc, g%iec
542  if (associated(mask)) then
543  if (mask(i,j,1) == 0.) cycle
544  endif
545  h_src(:) = h(i,j,:)
546  h_dest(:) = remap_cs%h(i,j,:)
547  call reintegrate_column(nz_src, h_src, field(i,j,:), &
548  nz_dest, h_dest, 0., reintegrated_field(i,j,:))
549  enddo
550  enddo
551  else
552  call assert(.false., 'vertically_reintegrate_diag_field: Q point remapping is not coded yet.')
553  endif
554 
556 
557 !> Vertically interpolate diagnostic field to alternative vertical grid.
558 subroutine vertically_interpolate_diag_field(remap_cs, G, h, staggered_in_x, staggered_in_y, &
559  mask, missing_value, field, interpolated_field)
560  type(diag_remap_ctrl), intent(in) :: remap_cs !< Diagnostic coodinate control structure
561  type(ocean_grid_type), intent(in) :: g !< Ocean grid structure
562  real, dimension(:,:,:), intent(in) :: h !< The current thicknesses
563  logical, intent(in) :: staggered_in_x !< True is the x-axis location is at u or q points
564  logical, intent(in) :: staggered_in_y !< True is the y-axis location is at v or q points
565  real, dimension(:,:,:), pointer :: mask !< A mask for the field
566  real, intent(in) :: missing_value !< A missing_value to assign land/vanished points
567  real, dimension(:,:,:), intent(in) :: field !< The diagnostic field to be remapped
568  real, dimension(:,:,:), intent(inout) :: interpolated_field !< Field argument remapped to alternative coordinate
569  ! Local variables
570  real, dimension(remap_cs%nz) :: h_dest
571  real, dimension(size(h,3)) :: h_src
572  integer :: nz_src, nz_dest
573  integer :: i, j, k !< Grid index
574  integer :: i1, j1 !< 1-based index
575  integer :: i_lo, i_hi, j_lo, j_hi !< (uv->h) interpolation indices
576  integer :: shift !< Symmetric offset for 1-based indexing
577 
578  call assert(remap_cs%initialized, 'vertically_interpolate_diag_field: remap_cs not initialized.')
579  call assert(size(field, 3) == size(h, 3)+1, &
580  'vertically_interpolate_diag_field: Remap field and thickness z-axes do not match.')
581 
582  interpolated_field(:,:,:) = 0.
583 
584  nz_src = size(h,3)
585  nz_dest = remap_cs%nz
586 
587  ! Symmetric grid offset under 1-based indexing; see header for details.
588  shift = 0; if (g%symmetric) shift = 1
589 
590  if (staggered_in_x .and. .not. staggered_in_y) then
591  ! U-points
592  do j=g%jsc, g%jec
593  do i=g%iscB, g%iecB
594  i1 = i - g%isdB + 1
595  i_lo = i1 - shift; i_hi = i_lo + 1
596  if (associated(mask)) then
597  if (mask(i,j,1) == 0.) cycle
598  endif
599  h_src(:) = 0.5 * (h(i_lo,j,:) + h(i_hi,j,:))
600  h_dest(:) = 0.5 * (remap_cs%h(i_lo,j,:) + remap_cs%h(i_hi,j,:))
601  call interpolate_column(nz_src, h_src, field(i1,j,:), &
602  nz_dest, h_dest, 0., interpolated_field(i1,j,:))
603  enddo
604  enddo
605  elseif (staggered_in_y .and. .not. staggered_in_x) then
606  ! V-points
607  do j=g%jscB, g%jecB
608  j1 = j - g%jsdB + 1
609  j_lo = j1 - shift; j_hi = j_lo + 1
610  do i=g%isc, g%iec
611  if (associated(mask)) then
612  if (mask(i,j,1) == 0.) cycle
613  endif
614  h_src(:) = 0.5 * (h(i,j_lo,:) + h(i,j_hi,:))
615  h_dest(:) = 0.5 * (remap_cs%h(i,j_lo,:) + remap_cs%h(i,j_hi,:))
616  call interpolate_column(nz_src, h_src, field(i,j1,:), &
617  nz_dest, h_dest, 0., interpolated_field(i,j1,:))
618  enddo
619  enddo
620  elseif ((.not. staggered_in_x) .and. (.not. staggered_in_y)) then
621  ! H-points
622  do j=g%jsc, g%jec
623  do i=g%isc, g%iec
624  if (associated(mask)) then
625  if (mask(i,j,1) == 0.) cycle
626  endif
627  h_src(:) = h(i,j,:)
628  h_dest(:) = remap_cs%h(i,j,:)
629  call interpolate_column(nz_src, h_src, field(i,j,:), &
630  nz_dest, h_dest, 0., interpolated_field(i,j,:))
631  enddo
632  enddo
633  else
634  call assert(.false., 'vertically_interpolate_diag_field: Q point remapping is not coded yet.')
635  endif
636 
638 
639 !> Horizontally average field
640 subroutine horizontally_average_diag_field(G, GV, h, staggered_in_x, staggered_in_y, &
641  is_layer, is_extensive, &
642  missing_value, field, averaged_field, &
643  averaged_mask)
644  type(ocean_grid_type), intent(in) :: g !< Ocean grid structure
645  type(verticalgrid_type), intent(in) :: gv !< The ocean vertical grid structure
646  real, dimension(:,:,:), intent(in) :: h !< The current thicknesses
647  logical, intent(in) :: staggered_in_x !< True if the x-axis location is at u or q points
648  logical, intent(in) :: staggered_in_y !< True if the y-axis location is at v or q points
649  logical, intent(in) :: is_layer !< True if the z-axis location is at h points
650  logical, intent(in) :: is_extensive !< True if the z-direction is spatially integrated (over layers)
651  real, intent(in) :: missing_value !< A missing_value to assign land/vanished points
652  real, dimension(:,:,:), intent(in) :: field !< The diagnostic field to be remapped
653  real, dimension(:), intent(inout) :: averaged_field !< Field argument horizontally averaged
654  logical, dimension(:), intent(inout) :: averaged_mask !< Mask for horizontally averaged field
655 
656  ! Local variables
657  real, dimension(G%isc:G%iec, G%jsc:G%jec, size(field,3)) :: volume, stuff
658  real, dimension(size(field, 3)) :: vol_sum, stuff_sum ! nz+1 is needed for interface averages
659  real :: height
660  integer :: i, j, k, nz
661  integer :: i1, j1 !< 1-based index
662 
663  nz = size(field, 3)
664 
665  ! TODO: These averages could potentially be modified to use the function in
666  ! the MOM_spatial_means module.
667  ! NOTE: Reproducible sums must be computed in the original MKS units
668 
669  if (staggered_in_x .and. .not. staggered_in_y) then
670  if (is_layer) then
671  ! U-points
672  do k=1,nz
673  vol_sum(k) = 0.
674  stuff_sum(k) = 0.
675  if (is_extensive) then
676  do j=g%jsc, g%jec ; do i=g%isc, g%iec
677  i1 = i - g%isdB + 1
678  volume(i,j,k) = (g%US%L_to_m**2 * g%areaCu(i,j)) * g%mask2dCu(i,j)
679  stuff(i,j,k) = volume(i,j,k) * field(i1,j,k)
680  enddo ; enddo
681  else ! Intensive
682  do j=g%jsc, g%jec ; do i=g%isc, g%iec
683  i1 = i - g%isdB + 1
684  height = 0.5 * (h(i,j,k) + h(i+1,j,k))
685  volume(i,j,k) = (g%US%L_to_m**2 * g%areaCu(i,j)) &
686  * (gv%H_to_m * height) * g%mask2dCu(i,j)
687  stuff(i,j,k) = volume(i,j,k) * field(i1,j,k)
688  enddo ; enddo
689  endif
690  enddo
691  else ! Interface
692  do k=1,nz
693  do j=g%jsc, g%jec ; do i=g%isc, g%iec
694  i1 = i - g%isdB + 1
695  volume(i,j,k) = (g%US%L_to_m**2 * g%areaCu(i,j)) * g%mask2dCu(i,j)
696  stuff(i,j,k) = volume(i,j,k) * field(i1,j,k)
697  enddo ; enddo
698  enddo
699  endif
700  elseif (staggered_in_y .and. .not. staggered_in_x) then
701  if (is_layer) then
702  ! V-points
703  do k=1,nz
704  if (is_extensive) then
705  do j=g%jsc, g%jec ; do i=g%isc, g%iec
706  j1 = j - g%jsdB + 1
707  volume(i,j,k) = (g%US%L_to_m**2 * g%areaCv(i,j)) * g%mask2dCv(i,j)
708  stuff(i,j,k) = volume(i,j,k) * field(i,j1,k)
709  enddo ; enddo
710  else ! Intensive
711  do j=g%jsc, g%jec ; do i=g%isc, g%iec
712  j1 = j - g%jsdB + 1
713  height = 0.5 * (h(i,j,k) + h(i,j+1,k))
714  volume(i,j,k) = (g%US%L_to_m**2 * g%areaCv(i,j)) &
715  * (gv%H_to_m * height) * g%mask2dCv(i,j)
716  stuff(i,j,k) = volume(i,j,k) * field(i,j1,k)
717  enddo ; enddo
718  endif
719  enddo
720  else ! Interface
721  do k=1,nz
722  do j=g%jsc, g%jec ; do i=g%isc, g%iec
723  j1 = j - g%jsdB + 1
724  volume(i,j,k) = (g%US%L_to_m**2 * g%areaCv(i,j)) * g%mask2dCv(i,j)
725  stuff(i,j,k) = volume(i,j,k) * field(i,j1,k)
726  enddo ; enddo
727  enddo
728  endif
729  elseif ((.not. staggered_in_x) .and. (.not. staggered_in_y)) then
730  if (is_layer) then
731  ! H-points
732  do k=1,nz
733  if (is_extensive) then
734  do j=g%jsc, g%jec ; do i=g%isc, g%iec
735  if (h(i,j,k) > 0.) then
736  volume(i,j,k) = (g%US%L_to_m**2 * g%areaT(i,j)) * g%mask2dT(i,j)
737  stuff(i,j,k) = volume(i,j,k) * field(i,j,k)
738  else
739  volume(i,j,k) = 0.
740  stuff(i,j,k) = 0.
741  endif
742  enddo ; enddo
743  else ! Intensive
744  do j=g%jsc, g%jec ; do i=g%isc, g%iec
745  volume(i,j,k) = (g%US%L_to_m**2 * g%areaT(i,j)) &
746  * (gv%H_to_m * h(i,j,k)) * g%mask2dT(i,j)
747  stuff(i,j,k) = volume(i,j,k) * field(i,j,k)
748  enddo ; enddo
749  endif
750  enddo
751  else ! Interface
752  do k=1,nz
753  do j=g%jsc, g%jec ; do i=g%isc, g%iec
754  volume(i,j,k) = (g%US%L_to_m**2 * g%areaT(i,j)) * g%mask2dT(i,j)
755  stuff(i,j,k) = volume(i,j,k) * field(i,j,k)
756  enddo ; enddo
757  enddo
758  endif
759  else
760  call assert(.false., 'horizontally_average_diag_field: Q point averaging is not coded yet.')
761  endif
762 
763  do k = 1,nz
764  vol_sum(k) = reproducing_sum(volume(:,:,k))
765  stuff_sum(k) = reproducing_sum(stuff(:,:,k))
766  enddo
767 
768  averaged_mask(:) = .true.
769  do k=1,nz
770  if (vol_sum(k) > 0.) then
771  averaged_field(k) = stuff_sum(k) / vol_sum(k)
772  else
773  averaged_field(k) = 0.
774  averaged_mask(k) = .false.
775  endif
776  enddo
777 
778 end subroutine horizontally_average_diag_field
779 
780 end module mom_diag_remap
mom_diag_remap::diag_remap_diag_registration_closed
subroutine, public diag_remap_diag_registration_closed(remap_cs)
Inform that all diagnostics have been registered. If _set_active() has not been called on the remappi...
Definition: MOM_diag_remap.F90:159
mom_regridding::get_sigma_cs
type(sigma_cs) function, public get_sigma_cs(CS)
This returns a copy of the sigma_CS stored in the regridding control structure.
Definition: MOM_regridding.F90:2326
coord_rho
Regrid columns for the continuous isopycnal (rho) coordinate.
Definition: coord_rho.F90:2
mom_diag_remap::vertically_reintegrate_diag_field
subroutine, public vertically_reintegrate_diag_field(remap_cs, G, h, staggered_in_x, staggered_in_y, mask, missing_value, field, reintegrated_field)
Vertically re-grid an already vertically-integrated diagnostic field to alternative vertical grid.
Definition: MOM_diag_remap.F90:479
mom_regridding::set_regrid_params
subroutine, public set_regrid_params(CS, boundary_extrapolation, min_thickness, old_grid_weight, interp_scheme, depth_of_time_filter_shallow, depth_of_time_filter_deep, compress_fraction, dz_min_surface, nz_fixed_surface, Rho_ML_avg_depth, nlay_ML_to_interior, fix_haloclines, halocline_filt_len, halocline_strat_tol, integrate_downward_for_e, adaptTimeRatio, adaptZoom, adaptZoomCoeff, adaptBuoyCoeff, adaptAlpha, adaptDoMin)
Can be used to set any of the parameters for MOM_regridding.
Definition: MOM_regridding.F90:2218
mom_diag_vkernels
Provides kernels for single-column interpolation, re-integration (re-mapping of integrated quantities...
Definition: MOM_diag_vkernels.F90:3
mom_regridding::initialize_regridding
subroutine, public initialize_regridding(CS, GV, US, max_depth, param_file, mdl, coord_mode, param_prefix, param_suffix)
Initialization and configures a regridding control structure based on customizable run-time parameter...
Definition: MOM_regridding.F90:177
mom_regridding::regridding_cs
Regridding control structure.
Definition: MOM_regridding.F90:45
mom_diag_vkernels::interpolate_column
subroutine, public interpolate_column(nsrc, h_src, u_src, ndest, h_dest, missing_value, u_dest)
Linearly interpolate interface data, u_src, from grid h_src to a grid h_dest.
Definition: MOM_diag_vkernels.F90:17
mom_remapping::remapping_core_h
subroutine, public remapping_core_h(CS, n0, h0, u0, n1, h1, u1, h_neglect, h_neglect_edge)
Remaps column of values u0 on grid h0 to grid h1 assuming the top edge is aligned.
Definition: MOM_remapping.F90:188
mom_verticalgrid
Provides a transparent vertical ocean grid type and supporting routines.
Definition: MOM_verticalGrid.F90:2
mom_diag_remap::diag_remap_end
subroutine, public diag_remap_end(remap_cs)
De-init a diagnostic remapping type. Free allocated memory.
Definition: MOM_diag_remap.F90:142
mom_diag_remap::diag_remap_init
subroutine, public diag_remap_init(remap_cs, coord_tuple)
Initialize a diagnostic remapping type with the given vertical coordinate.
Definition: MOM_diag_remap.F90:124
mom_string_functions
Handy functions for manipulating strings.
Definition: MOM_string_functions.F90:2
coord_zlike::build_zstar_column
subroutine, public build_zstar_column(CS, depth, total_thickness, zInterface, z_rigid_top, eta_orig, zScale)
Builds a z* coordinate with a minimum thickness.
Definition: coord_zlike.F90:65
mom_eos
Provides subroutines for quantities specific to the equation of state.
Definition: MOM_EOS.F90:2
mom_remapping::remapping_cs
Container for remapping parameters.
Definition: MOM_remapping.F90:22
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_diag_remap::diag_remap_do_remap
subroutine, public diag_remap_do_remap(remap_cs, G, GV, h, staggered_in_x, staggered_in_y, mask, missing_value, field, remapped_field)
Remap diagnostic field to alternative vertical grid.
Definition: MOM_diag_remap.F90:340
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_io
This module contains I/O framework code.
Definition: MOM_io.F90:2
mom_debugging::check_column_integrals
logical function, public check_column_integrals(nk_1, field_1, nk_2, field_2, missing_value)
Returns false if the column integrals of two given quantities are within roundoff of each other.
Definition: MOM_debugging.F90:823
mom_diag_remap::diag_remap_get_axes_info
subroutine, public diag_remap_get_axes_info(remap_cs, nz, id_layer, id_interface)
Get layer and interface axes ids for this coordinate Needed when defining axes groups.
Definition: MOM_diag_remap.F90:239
mom_unit_scaling::unit_scale_type
Describes various unit conversion factors.
Definition: MOM_unit_scaling.F90:14
mom_coms
Interfaces to non-domain-oriented communication subroutines, including the MOM6 reproducing sums faci...
Definition: MOM_coms.F90:3
regrid_consts
Contains constants for interpreting input parameters that control regridding.
Definition: regrid_consts.F90:2
mom_remapping
Provides column-wise vertical remapping functions.
Definition: MOM_remapping.F90:2
mom_regridding::getcoordinateinterfaces
real function, dimension(cs%nk+1), public getcoordinateinterfaces(CS, undo_scaling)
Query the target coordinate interface positions.
Definition: MOM_regridding.F90:2113
mom_eos::eos_type
A control structure for the equation of state.
Definition: MOM_EOS.F90:91
mom_regridding::get_regrid_size
integer function, public get_regrid_size(CS)
Returns the number of levels/layers in the regridding control structure.
Definition: MOM_regridding.F90:2310
mom_verticalgrid::verticalgrid_type
Describes the vertical ocean grid, including unit conversion factors.
Definition: MOM_verticalGrid.F90:24
coord_zlike
Regrid columns for a z-like coordinate (z-star, z-level)
Definition: coord_zlike.F90:2
coord_sigma::build_sigma_column
subroutine, public build_sigma_column(CS, depth, totalThickness, zInterface)
Build a sigma coordinate column.
Definition: coord_sigma.F90:64
mom_io::mom_read_data
Read a data field from a file.
Definition: MOM_io.F90:74
mom_file_parser
The MOM6 facility to parse input files for runtime parameters.
Definition: MOM_file_parser.F90:2
mom_diag_remap
provides runtime remapping of diagnostics to z star, sigma and rho vertical coordinates.
Definition: MOM_diag_remap.F90:56
mom_string_functions::extractword
character(len=120) function, public extractword(string, n)
Returns the string corresponding to the nth word in the argument or "" if the string is not long enou...
Definition: MOM_string_functions.F90:196
mom_grid
Provides the ocean grid type.
Definition: MOM_grid.F90:2
mom_unit_scaling
Provides a transparent unit rescaling type to facilitate dimensional consistency testing.
Definition: MOM_unit_scaling.F90:2
regrid_consts::coordinatemode
integer function coordinatemode(string)
Parse a string parameter specifying the coordinate mode and return the appropriate enumerated integer...
Definition: regrid_consts.F90:54
mom_error_handler::assert
subroutine, public assert(logical_arg, msg)
Issues a FATAL error if the assertion fails, i.e. the first argument is false.
Definition: MOM_error_handler.F90:182
mom_diag_remap::horizontally_average_diag_field
subroutine, public horizontally_average_diag_field(G, GV, h, staggered_in_x, staggered_in_y, is_layer, is_extensive, missing_value, field, averaged_field, averaged_mask)
Horizontally average field.
Definition: MOM_diag_remap.F90:644
mom_regridding::get_rho_cs
type(rho_cs) function, public get_rho_cs(CS)
This returns a copy of the rho_CS stored in the regridding control structure.
Definition: MOM_regridding.F90:2334
mom_debugging
Provides checksumming functions for debugging.
Definition: MOM_debugging.F90:7
mom_io::file_exists
Indicate whether a file exists, perhaps with domain decomposition.
Definition: MOM_io.F90:68
coord_rho::build_rho_column
subroutine, public build_rho_column(CS, nz, depth, h, T, S, eqn_of_state, z_interface, h_neglect, h_neglect_edge)
Build a rho coordinate column.
Definition: coord_rho.F90:102
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_diag_remap::vertically_interpolate_diag_field
subroutine, public vertically_interpolate_diag_field(remap_cs, G, h, staggered_in_x, staggered_in_y, mask, missing_value, field, interpolated_field)
Vertically interpolate diagnostic field to alternative vertical grid.
Definition: MOM_diag_remap.F90:560
mom_diag_remap::diag_remap_calc_hmask
subroutine, public diag_remap_calc_hmask(remap_cs, G, mask)
Calculate masks for target grid.
Definition: MOM_diag_remap.F90:435
mom_string_functions::lowercase
character(len=len(input_string)) function, public lowercase(input_string)
Return a string in which all uppercase letters have been replaced by their lowercase counterparts.
Definition: MOM_string_functions.F90:24
mom_coms::reproducing_sum
Find an accurate and order-invariant sum of distributed 2d or 3d fields.
Definition: MOM_coms.F90:54
mom_diag_remap::diag_remap_configure_axes
subroutine, public diag_remap_configure_axes(remap_cs, GV, US, param_file)
Configure the vertical axes for a diagnostic remapping control structure. Reads a configuration param...
Definition: MOM_diag_remap.F90:180
mom_regridding
Generates vertical grids as part of the ALE algorithm.
Definition: MOM_regridding.F90:2
mom_file_parser::log_param
An overloaded interface to log the values of various types of parameters.
Definition: MOM_file_parser.F90:96
mom_diag_remap::diag_remap_set_active
subroutine, public diag_remap_set_active(remap_cs)
Indicate that this remapping type is actually used by the diag manager. If this is never called then ...
Definition: MOM_diag_remap.F90:171
mom_diag_remap::diag_remap_update
subroutine, public diag_remap_update(remap_cs, G, GV, US, h, T, S, eqn_of_state)
Build/update target vertical grids for diagnostic remapping.
Definition: MOM_diag_remap.F90:268
mom_diag_remap::diag_remap_ctrl
Represents remapping of diagnostics to a particular vertical coordinate.
Definition: MOM_diag_remap.F90:103
mom_regridding::get_zlike_cs
type(zlike_cs) function, public get_zlike_cs(CS)
This returns a copy of the zlike_CS stored in the regridding control structure.
Definition: MOM_regridding.F90:2318
mom_diag_remap::diag_remap_axes_configured
logical function, public diag_remap_axes_configured(remap_cs)
Whether or not the axes for this vertical coordinated has been configured. Configuration is complete ...
Definition: MOM_diag_remap.F90:255
mom_remapping::initialize_remapping
subroutine, public initialize_remapping(CS, remapping_scheme, boundary_extrapolation, check_reconstruction, check_remapping, force_bounds_in_subcell, answers_2018)
Constructor for remapping control structure.
Definition: MOM_remapping.F90:1547
coord_sigma
Regrid columns for the sigma coordinate.
Definition: coord_sigma.F90:2
mom_error_handler
Routines for error handling and I/O management.
Definition: MOM_error_handler.F90:2
mom_grid::ocean_grid_type
Ocean grid type. See mom_grid for details.
Definition: MOM_grid.F90:26
mom_diag_vkernels::reintegrate_column
subroutine, public reintegrate_column(nsrc, h_src, uh_src, ndest, h_dest, missing_value, uh_dest)
Conservatively calculate integrated data, uh_dest, on grid h_dest, from layer-integrated data,...
Definition: MOM_diag_vkernels.F90:92