MOM6
mom_tracer_initialization_from_z Module Reference

Detailed Description

Initializes hydrography from z-coordinate climatology files.

Functions/Subroutines

subroutine, public mom_initialize_tracer_from_z (h, tr, G, GV, US, PF, src_file, src_var_nam, src_var_unit_conversion, src_var_record, homogenize, useALEremapping, remappingScheme, src_var_gridspec)
 Initializes a tracer from a z-space data file. More...
 

Variables

character(len=40) mdl = "MOM_tracer_initialization_from_Z"
 This module's name. More...
 

Function/Subroutine Documentation

◆ mom_initialize_tracer_from_z()

subroutine, public mom_tracer_initialization_from_z::mom_initialize_tracer_from_z ( real, dimension(szi_(g),szj_(g),szk_(g)), intent(in)  h,
real, dimension(:,:,:), pointer  tr,
type(ocean_grid_type), intent(inout)  G,
type(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
type(param_file_type), intent(in)  PF,
character(len=*), intent(in)  src_file,
character(len=*), intent(in)  src_var_nam,
real, intent(in), optional  src_var_unit_conversion,
integer, intent(in), optional  src_var_record,
logical, intent(in), optional  homogenize,
logical, intent(in), optional  useALEremapping,
character(len=*), intent(in), optional  remappingScheme,
character(len=*), intent(in), optional  src_var_gridspec 
)

Initializes a tracer from a z-space data file.

Parameters
[in,out]gOcean grid structure.
[in]gvOcean vertical grid structure.
[in]usA dimensional unit scaling type
[in]hLayer thickness [H ~> m or kg m-2].
trPointer to array to be initialized
[in]pfparameter file
[in]src_filesource filename
[in]src_var_namvariable name in file
[in]src_var_unit_conversionoptional multiplicative unit conversion
[in]src_var_recordrecord to read for multiple time-level files
[in]homogenizeoptionally homogenize to mean value
[in]usealeremappingto remap or not (optional)
[in]remappingschemeremapping scheme to use.
[in]src_var_gridspecSource variable name in a gridspec file. This is not implemented yet.

Definition at line 49 of file MOM_tracer_initialization_from_Z.F90.

49  type(ocean_grid_type), intent(inout) :: G !< Ocean grid structure.
50  type(verticalGrid_type), intent(in) :: GV !< Ocean vertical grid structure.
51  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
52  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
53  intent(in) :: h !< Layer thickness [H ~> m or kg m-2].
54  real, dimension(:,:,:), pointer :: tr !< Pointer to array to be initialized
55  type(param_file_type), intent(in) :: PF !< parameter file
56  character(len=*), intent(in) :: src_file !< source filename
57  character(len=*), intent(in) :: src_var_nam !< variable name in file
58  real, optional, intent(in) :: src_var_unit_conversion !< optional multiplicative unit conversion
59  integer, optional, intent(in) :: src_var_record !< record to read for multiple time-level files
60  logical, optional, intent(in) :: homogenize !< optionally homogenize to mean value
61  logical, optional, intent(in) :: useALEremapping !< to remap or not (optional)
62  character(len=*), optional, intent(in) :: remappingScheme !< remapping scheme to use.
63  character(len=*), optional, intent(in) :: src_var_gridspec !< Source variable name in a gridspec file.
64  !! This is not implemented yet.
65  ! Local variables
66  real :: land_fill = 0.0
67  character(len=200) :: inputdir ! The directory where NetCDF input files are.
68  character(len=200) :: mesg
69  real :: convert
70  integer :: recnum
71  character(len=10) :: remapScheme
72  logical :: homog,useALE
73 
74  ! This include declares and sets the variable "version".
75 # include "version_variable.h"
76  character(len=40) :: mdl = "MOM_initialize_tracers_from_Z" ! This module's name.
77 
78  integer :: is, ie, js, je, nz ! compute domain indices
79  integer :: isd, ied, jsd, jed ! data domain indices
80  integer :: i, j, k, kd
81  real, allocatable, dimension(:,:,:), target :: tr_z, mask_z
82  real, allocatable, dimension(:), target :: z_edges_in, z_in
83 
84  ! Local variables for ALE remapping
85  real, dimension(:,:,:), allocatable :: hSrc ! Source thicknesses [H ~> m or kg m-2].
86  real, dimension(:), allocatable :: h1 ! A 1-d column of source thicknesses [Z ~> m].
87  real :: zTopOfCell, zBottomOfCell, z_bathy ! Heights [Z ~> m].
88  type(remapping_CS) :: remapCS ! Remapping parameters and work arrays
89 
90  real :: missing_value
91  integer :: nPoints
92  integer :: id_clock_routine, id_clock_ALE
93  logical :: reentrant_x, tripolar_n
94 
95  id_clock_routine = cpu_clock_id('(Initialize tracer from Z)', grain=clock_routine)
96  id_clock_ale = cpu_clock_id('(Initialize tracer from Z) ALE', grain=clock_loop)
97 
98  call cpu_clock_begin(id_clock_routine)
99 
100  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
101  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
102 
103  call calltree_enter(trim(mdl)//"(), MOM_state_initialization.F90")
104 
105  call get_param(pf, mdl, "Z_INIT_HOMOGENIZE", homog, &
106  "If True, then horizontally homogenize the interpolated "//&
107  "initial conditions.", default=.false.)
108  call get_param(pf, mdl, "Z_INIT_ALE_REMAPPING", useale, &
109  "If True, then remap straight to model coordinate from file.",&
110  default=.true.)
111  call get_param(pf, mdl, "Z_INIT_REMAPPING_SCHEME", remapscheme, &
112  "The remapping scheme to use if using Z_INIT_ALE_REMAPPING is True.", &
113  default="PLM")
114 
115  ! These are model grid properties, but being applied to the data grid for now.
116  ! need to revisit this (mjh)
117  reentrant_x = .false. ; call get_param(pf, mdl, "REENTRANT_X", reentrant_x,default=.true.)
118  tripolar_n = .false. ; call get_param(pf, mdl, "TRIPOLAR_N", tripolar_n, default=.false.)
119 
120  if (PRESENT(homogenize)) homog=homogenize
121  if (PRESENT(usealeremapping)) useale=usealeremapping
122  if (PRESENT(remappingscheme)) remapscheme=remappingscheme
123  recnum=1
124  if (PRESENT(src_var_record)) recnum = src_var_record
125  convert=1.0
126  if (PRESENT(src_var_unit_conversion)) convert = src_var_unit_conversion
127 
128  call horiz_interp_and_extrap_tracer(src_file, src_var_nam, convert, recnum, &
129  g, tr_z, mask_z, z_in, z_edges_in, missing_value, reentrant_x, tripolar_n, &
130  homog, m_to_z=us%m_to_Z)
131 
132  kd = size(z_edges_in,1)-1
133  call pass_var(tr_z,g%Domain)
134  call pass_var(mask_z,g%Domain)
135 
136 ! Done with horizontal interpolation.
137 ! Now remap to model coordinates
138  if (useale) then
139  call cpu_clock_begin(id_clock_ale)
140  ! First we reserve a work space for reconstructions of the source data
141  allocate( h1(kd) )
142  allocate( hsrc(isd:ied,jsd:jed,kd) )
143  call initialize_remapping( remapcs, remapscheme, boundary_extrapolation=.false. ) ! Data for reconstructions
144  ! Next we initialize the regridding package so that it knows about the target grid
145 
146  do j = js, je ; do i = is, ie
147  if (g%mask2dT(i,j)>0.) then
148  ! Build the source grid
149  ztopofcell = 0. ; zbottomofcell = 0. ; npoints = 0
150  z_bathy = g%bathyT(i,j)
151  do k = 1, kd
152  if (mask_z(i,j,k) > 0.) then
153  zbottomofcell = -min( z_edges_in(k+1), z_bathy )
154  elseif (k>1) then
155  zbottomofcell = -z_bathy
156  endif
157  h1(k) = ztopofcell - zbottomofcell
158  if (h1(k)>0.) npoints = npoints + 1
159  ztopofcell = zbottomofcell ! Bottom becomes top for next value of k
160  enddo
161  h1(kd) = h1(kd) + ( ztopofcell + z_bathy ) ! In case data is deeper than model
162  else
163  tr(i,j,:) = 0.
164  endif ! mask2dT
165  hsrc(i,j,:) = gv%Z_to_H * h1(:)
166  enddo ; enddo
167 
168  call ale_remap_scalar(remapcs, g, gv, kd, hsrc, tr_z, h, tr, all_cells=.false. )
169 
170  deallocate( hsrc )
171  deallocate( h1 )
172 
173  do k=1,nz
174  call mystats(tr(:,:,k), missing_value, is, ie, js, je, k, 'Tracer from ALE()')
175  enddo
176  call cpu_clock_end(id_clock_ale)
177  endif ! useALEremapping
178 
179 ! Fill land values
180  do k=1,nz ; do j=js,je ; do i=is,ie
181  if (tr(i,j,k) == missing_value) then
182  tr(i,j,k)=land_fill
183  endif
184  enddo ; enddo ; enddo
185 
186  call calltree_leave(trim(mdl)//'()')
187  call cpu_clock_end(id_clock_routine)
188 

References mom_ale::ale_remap_scalar(), mom_error_handler::calltree_enter(), mom_error_handler::calltree_leave(), mdl, and mom_horizontal_regridding::mystats().

Here is the call graph for this function:

Variable Documentation

◆ mdl

character(len=40) mom_tracer_initialization_from_z::mdl = "MOM_tracer_initialization_from_Z"
private

This module's name.

Definition at line 41 of file MOM_tracer_initialization_from_Z.F90.

41 character(len=40) :: mdl = "MOM_tracer_initialization_from_Z" !< This module's name.

Referenced by mom_initialize_tracer_from_z().