MOM6
mom_io Module Reference

Detailed Description

This module contains I/O framework code.

This file contains a number of subroutines that manipulate NetCDF files and handle input and output of fields. These subroutines, along with their purpose, are:

  • create_file: create a new file and set up structures that are needed for subsequent output and write out the coordinates.
  • reopen_file: reopen an existing file for writing and set up structures that are needed for subsequent output.
  • open_input_file: open the indicated file for reading only.
  • close_file: close an open file.
  • synch_file: flush the buffers, completing all pending output.
  • write_field: write a field to an open file.
  • write_time: write a value of the time axis to an open file.
  • read_data: read a variable from an open file.
  • read_time: read a time from an open file.
  • name_output_file: provide a name for an output file based on a name root and the time of the output.
  • find_input_file: find a file that has been previously written by MOM and named by name_output_file and open it for reading.
  • handle_error: write an error code and quit.

Data Types

interface  file_exists
 Indicate whether a file exists, perhaps with domain decomposition. More...
 
interface  mom_read_data
 Read a data field from a file. More...
 
interface  mom_read_vector
 Read a pair of data fields representing the two components of a vector from a file. More...
 
type  vardesc
 Type for describing a variable, typically a tracer. More...
 

Functions/Subroutines

subroutine, public create_file (unit, filename, vars, novars, fields, threading, timeunit, G, dG, GV, checksums)
 Routine creates a new NetCDF file. It also sets up structures that describe this file and variables that will later be written to this file. Type for describing a variable, typically a tracer. More...
 
subroutine, public reopen_file (unit, filename, vars, novars, fields, threading, timeunit, G, dG, GV)
 This routine opens an existing NetCDF file for output. If it does not find the file, a new file is created. It also sets up structures that describe this file and the variables that will later be written to this file. More...
 
subroutine, public read_axis_data (filename, axis_name, var)
 Read the data associated with a named axis in a file. More...
 
integer function, public num_timelevels (filename, varname, min_dims)
 This function determines how many time levels a variable has. More...
 
type(vardesc) function, public var_desc (name, units, longname, hor_grid, z_grid, t_grid, cmor_field_name, cmor_units, cmor_longname, conversion, caller)
 Returns a vardesc type whose elements have been filled with the provided fields. The argument name is required, while the others are optional and have default values that are empty strings or are appropriate for a 3-d tracer field at the tracer cell centers. More...
 
subroutine, public modify_vardesc (vd, name, units, longname, hor_grid, z_grid, t_grid, cmor_field_name, cmor_units, cmor_longname, conversion, caller)
 This routine modifies the named elements of a vardesc type. All arguments are optional, except the vardesc type to be modified. More...
 
character(len=len(longname)) function, public cmor_long_std (longname)
 This function returns the CMOR standard name given a CMOR longname, based on the standard pattern of character conversions. More...
 
subroutine, public query_vardesc (vd, name, units, longname, hor_grid, z_grid, t_grid, cmor_field_name, cmor_units, cmor_longname, conversion, caller)
 This routine queries vardesc. More...
 
subroutine safe_string_copy (str1, str2, fieldnm, caller)
 Copies a string. More...
 
character(len=len(name)) function, public ensembler (name, ens_no_in)
 Returns a name with "%#E" or "%E" replaced with the ensemble member number. More...
 
logical function mom_file_exists (filename, MOM_Domain)
 Returns true if the named file or its domain-decomposed variant exists. More...
 
logical function fms_file_exists (filename, domain, no_domain)
 Returns true if the named file or its domain-decomposed variant exists. More...
 
subroutine mom_read_data_1d (filename, fieldname, data, timelevel, scale)
 This function uses the fms_io function read_data to read 1-D data field named "fieldname" from file "filename". More...
 
subroutine mom_read_data_2d (filename, fieldname, data, MOM_Domain, timelevel, position, scale)
 This function uses the fms_io function read_data to read a distributed 2-D data field named "fieldname" from file "filename". Valid values for "position" include CORNER, CENTER, EAST_FACE and NORTH_FACE. More...
 
subroutine mom_read_data_3d (filename, fieldname, data, MOM_Domain, timelevel, position, scale)
 This function uses the fms_io function read_data to read a distributed 3-D data field named "fieldname" from file "filename". Valid values for "position" include CORNER, CENTER, EAST_FACE and NORTH_FACE. More...
 
subroutine mom_read_data_4d (filename, fieldname, data, MOM_Domain, timelevel, position, scale)
 This function uses the fms_io function read_data to read a distributed 4-D data field named "fieldname" from file "filename". Valid values for "position" include CORNER, CENTER, EAST_FACE and NORTH_FACE. More...
 
subroutine mom_read_vector_2d (filename, u_fieldname, v_fieldname, u_data, v_data, MOM_Domain, timelevel, stagger, scalar_pair, scale)
 This function uses the fms_io function read_data to read a pair of distributed 2-D data fields with names given by "[uv]_fieldname" from file "filename". Valid values for "stagger" include CGRID_NE, BGRID_NE, and AGRID. More...
 
subroutine mom_read_vector_3d (filename, u_fieldname, v_fieldname, u_data, v_data, MOM_Domain, timelevel, stagger, scalar_pair, scale)
 This function uses the fms_io function read_data to read a pair of distributed 3-D data fields with names given by "[uv]_fieldname" from file "filename". Valid values for "stagger" include CGRID_NE, BGRID_NE, and AGRID. More...
 
subroutine, public mom_io_init (param_file)
 Initialize the MOM_io module. More...
 

Function/Subroutine Documentation

◆ cmor_long_std()

character(len=len(longname)) function, public mom_io::cmor_long_std ( character(len=*), intent(in)  longname)

This function returns the CMOR standard name given a CMOR longname, based on the standard pattern of character conversions.

Parameters
[in]longnameThe CMOR longname being converted
Returns
The CMOR standard name generated from longname

Definition at line 683 of file MOM_io.F90.

683  character(len=*), intent(in) :: longname !< The CMOR longname being converted
684  character(len=len(longname)) :: std_name !< The CMOR standard name generated from longname
685 
686  integer :: k
687 
688  std_name = lowercase(longname)
689 
690  do k=1, len_trim(std_name)
691  if (std_name(k:k) == ' ') std_name(k:k) = '_'
692  enddo
693 

References mom_string_functions::lowercase().

Referenced by mom_tracer_registry::register_tracer_diagnostics().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ create_file()

subroutine, public mom_io::create_file ( integer, intent(out)  unit,
character(len=*), intent(in)  filename,
type(vardesc), dimension(:), intent(in)  vars,
integer, intent(in)  novars,
type(fieldtype), dimension(:), intent(inout)  fields,
integer, intent(in), optional  threading,
real, intent(in), optional  timeunit,
type(ocean_grid_type), intent(in), optional  G,
type(dyn_horgrid_type), intent(in), optional  dG,
type(verticalgrid_type), intent(in), optional  GV,
integer(kind=8), dimension(:,:), intent(in), optional  checksums 
)

Routine creates a new NetCDF file. It also sets up structures that describe this file and variables that will later be written to this file. Type for describing a variable, typically a tracer.

Parameters
[out]unitunit id of an open file or -1 on a nonwriting PE with single file output
[in]filenamefull path to the file to create
[in]varsstructures describing fields written to filename
[in]novarsnumber of fields written to filename
[in,out]fieldsarray of fieldtypes for each variable
[in]threadingSINGLE_FILE or MULTIPLE
[in]timeunitlength of the units for time [s]. The default value is 86400.0, for 1 day.
[in]gocean horizontal grid structure; G or dG is required if the new file uses any horizontal grid axes.
[in]dgdynamic horizontal grid structure; G or dG is required if the new file uses any horizontal grid axes.
[in]gvocean vertical grid structure, which is required if the new file uses any vertical grid axes.
[in]checksumschecksums of vars

Definition at line 93 of file MOM_io.F90.

93  integer, intent(out) :: unit !< unit id of an open file or -1 on a
94  !! nonwriting PE with single file output
95  character(len=*), intent(in) :: filename !< full path to the file to create
96  type(vardesc), intent(in) :: vars(:) !< structures describing fields written to filename
97  integer, intent(in) :: novars !< number of fields written to filename
98  type(fieldtype), intent(inout) :: fields(:) !< array of fieldtypes for each variable
99  integer, optional, intent(in) :: threading !< SINGLE_FILE or MULTIPLE
100  real, optional, intent(in) :: timeunit !< length of the units for time [s]. The
101  !! default value is 86400.0, for 1 day.
102  type(ocean_grid_type), optional, intent(in) :: G !< ocean horizontal grid structure; G or dG
103  !! is required if the new file uses any
104  !! horizontal grid axes.
105  type(dyn_horgrid_type), optional, intent(in) :: dG !< dynamic horizontal grid structure; G or dG
106  !! is required if the new file uses any
107  !! horizontal grid axes.
108  type(verticalGrid_type), optional, intent(in) :: GV !< ocean vertical grid structure, which is
109  !! required if the new file uses any
110  !! vertical grid axes.
111  integer(kind=8), optional, intent(in) :: checksums(:,:) !< checksums of vars
112 
113  logical :: use_lath, use_lonh, use_latq, use_lonq, use_time
114  logical :: use_layer, use_int, use_periodic
115  logical :: one_file, domain_set
116  type(axistype) :: axis_lath, axis_latq, axis_lonh, axis_lonq
117  type(axistype) :: axis_layer, axis_int, axis_time, axis_periodic
118  type(axistype) :: axes(4)
119  type(MOM_domain_type), pointer :: Domain => null()
120  type(domain1d) :: x_domain, y_domain
121  integer :: numaxes, pack, thread, k
122  integer :: isg, ieg, jsg, jeg, IsgB, IegB, JsgB, JegB
123  integer :: var_periods, num_periods=0
124  real, dimension(:), allocatable :: period_val
125  real, pointer, dimension(:) :: &
126  gridLatT => null(), & ! The latitude or longitude of T or B points for
127  gridlatb => null(), & ! the purpose of labeling the output axes.
128  gridlont => null(), gridlonb => null()
129  character(len=40) :: time_units, x_axis_units, y_axis_units
130  character(len=8) :: t_grid, t_grid_read
131 
132  use_lath = .false. ; use_lonh = .false.
133  use_latq = .false. ; use_lonq = .false.
134  use_time = .false. ; use_periodic = .false.
135  use_layer = .false. ; use_int = .false.
136 
137  thread = single_file
138  if (PRESENT(threading)) thread = threading
139 
140  domain_set = .false.
141  if (present(g)) then
142  domain_set = .true. ; domain => g%Domain
143  gridlatt => g%gridLatT ; gridlatb => g%gridLatB
144  gridlont => g%gridLonT ; gridlonb => g%gridLonB
145  x_axis_units = g%x_axis_units ; y_axis_units = g%y_axis_units
146  isg = g%isg ; ieg = g%ieg ; jsg = g%jsg ; jeg = g%jeg
147  isgb = g%IsgB ; iegb = g%IegB ; jsgb = g%JsgB ; jegb = g%JegB
148  elseif (present(dg)) then
149  domain_set = .true. ; domain => dg%Domain
150  gridlatt => dg%gridLatT ; gridlatb => dg%gridLatB
151  gridlont => dg%gridLonT ; gridlonb => dg%gridLonB
152  x_axis_units = dg%x_axis_units ; y_axis_units = dg%y_axis_units
153  isg = dg%isg ; ieg = dg%ieg ; jsg = dg%jsg ; jeg = dg%jeg
154  isgb = dg%IsgB ; iegb = dg%IegB ; jsgb = dg%JsgB ; jegb = dg%JegB
155  endif
156 
157  one_file = .true.
158  if (domain_set) one_file = (thread == single_file)
159 
160  if (one_file) then
161  call open_file(unit, filename, mpp_overwr, mpp_netcdf, threading=thread)
162  else
163  call open_file(unit, filename, mpp_overwr, mpp_netcdf, domain=domain%mpp_domain)
164  endif
165 
166 ! Define the coordinates.
167  do k=1,novars
168  select case (vars(k)%hor_grid)
169  case ('h') ; use_lath = .true. ; use_lonh = .true.
170  case ('q') ; use_latq = .true. ; use_lonq = .true.
171  case ('u') ; use_lath = .true. ; use_lonq = .true.
172  case ('v') ; use_latq = .true. ; use_lonh = .true.
173  case ('T') ; use_lath = .true. ; use_lonh = .true.
174  case ('Bu') ; use_latq = .true. ; use_lonq = .true.
175  case ('Cu') ; use_lath = .true. ; use_lonq = .true.
176  case ('Cv') ; use_latq = .true. ; use_lonh = .true.
177  case ('1') ! Do nothing.
178  case default
179  call mom_error(warning, "MOM_io create_file: "//trim(vars(k)%name)//&
180  " has unrecognized hor_grid "//trim(vars(k)%hor_grid))
181  end select
182  select case (vars(k)%z_grid)
183  case ('L') ; use_layer = .true.
184  case ('i') ; use_int = .true.
185  case ('1') ! Do nothing.
186  case default
187  call mom_error(fatal, "MOM_io create_file: "//trim(vars(k)%name)//&
188  " has unrecognized z_grid "//trim(vars(k)%z_grid))
189  end select
190  t_grid = adjustl(vars(k)%t_grid)
191  select case (t_grid(1:1))
192  case ('s', 'a', 'm') ; use_time = .true.
193  case ('p') ; use_periodic = .true.
194  if (len_trim(t_grid(2:8)) <= 0) call mom_error(fatal, &
195  "MOM_io create_file: No periodic axis length was specified in "//&
196  trim(vars(k)%t_grid) // " in the periodic axes of variable "//&
197  trim(vars(k)%name)//" in file "//trim(filename))
198  var_periods = -9999999
199  t_grid_read = adjustl(t_grid(2:8))
200  read(t_grid_read,*) var_periods
201  if (var_periods == -9999999) call mom_error(fatal, &
202  "MOM_io create_file: Failed to read the number of periods from "//&
203  trim(vars(k)%t_grid) // " in the periodic axes of variable "//&
204  trim(vars(k)%name)//" in file "//trim(filename))
205  if (var_periods < 1) call mom_error(fatal, "MOM_io create_file: "//&
206  "variable "//trim(vars(k)%name)//" in file "//trim(filename)//&
207  " uses a periodic time axis, and must have a positive "//&
208  "value for the number of periods in "//vars(k)%t_grid )
209  if ((num_periods > 0) .and. (var_periods /= num_periods)) &
210  call mom_error(fatal, "MOM_io create_file: "//&
211  "Only one value of the number of periods can be used in the "//&
212  "create_file call for file "//trim(filename)//". The second is "//&
213  "variable "//trim(vars(k)%name)//" with t_grid "//vars(k)%t_grid )
214 
215  num_periods = var_periods
216  case ('1') ! Do nothing.
217  case default
218  call mom_error(warning, "MOM_io create_file: "//trim(vars(k)%name)//&
219  " has unrecognized t_grid "//trim(vars(k)%t_grid))
220  end select
221  enddo
222 
223  if ((use_lath .or. use_lonh .or. use_latq .or. use_lonq)) then
224  if (.not.domain_set) call mom_error(fatal, "create_file: "//&
225  "An ocean_grid_type or dyn_horgrid_type is required to create a file with a horizontal coordinate.")
226 
227  call mpp_get_domain_components(domain%mpp_domain, x_domain, y_domain)
228  endif
229  if ((use_layer .or. use_int) .and. .not.present(gv)) call mom_error(fatal, &
230  "create_file: A vertical grid type is required to create a file with a vertical coordinate.")
231 
232 ! Specify all optional arguments to mpp_write_meta: name, units, longname, cartesian, calendar, sense,
233 ! domain, data, min). Otherwise if optional arguments are added to mpp_write_meta the compiler may
234 ! (and in case of GNU does) get confused and crash.
235  if (use_lath) &
236  call mpp_write_meta(unit, axis_lath, name="lath", units=y_axis_units, longname="Latitude", &
237  cartesian='Y', domain = y_domain, data=gridlatt(jsg:jeg))
238 
239  if (use_lonh) &
240  call mpp_write_meta(unit, axis_lonh, name="lonh", units=x_axis_units, longname="Longitude", &
241  cartesian='X', domain = x_domain, data=gridlont(isg:ieg))
242 
243  if (use_latq) &
244  call mpp_write_meta(unit, axis_latq, name="latq", units=y_axis_units, longname="Latitude", &
245  cartesian='Y', domain = y_domain, data=gridlatb(jsgb:jegb))
246 
247  if (use_lonq) &
248  call mpp_write_meta(unit, axis_lonq, name="lonq", units=x_axis_units, longname="Longitude", &
249  cartesian='X', domain = x_domain, data=gridlonb(isgb:iegb))
250 
251  if (use_layer) &
252  call mpp_write_meta(unit, axis_layer, name="Layer", units=trim(gv%zAxisUnits), &
253  longname="Layer "//trim(gv%zAxisLongName), cartesian='Z', &
254  sense=1, data=gv%sLayer(1:gv%ke))
255 
256  if (use_int) &
257  call mpp_write_meta(unit, axis_int, name="Interface", units=trim(gv%zAxisUnits), &
258  longname="Interface "//trim(gv%zAxisLongName), cartesian='Z', &
259  sense=1, data=gv%sInterface(1:gv%ke+1))
260 
261  if (use_time) then ; if (present(timeunit)) then
262  ! Set appropriate units, depending on the value.
263  if (timeunit < 0.0) then
264  time_units = "days" ! The default value.
265  elseif ((timeunit >= 0.99) .and. (timeunit < 1.01)) then
266  time_units = "seconds"
267  elseif ((timeunit >= 3599.0) .and. (timeunit < 3601.0)) then
268  time_units = "hours"
269  elseif ((timeunit >= 86399.0) .and. (timeunit < 86401.0)) then
270  time_units = "days"
271  elseif ((timeunit >= 3.0e7) .and. (timeunit < 3.2e7)) then
272  time_units = "years"
273  else
274  write(time_units,'(es8.2," s")') timeunit
275  endif
276 
277  call mpp_write_meta(unit, axis_time, name="Time", units=time_units, longname="Time", cartesian='T')
278  else
279  call mpp_write_meta(unit, axis_time, name="Time", units="days", longname="Time",cartesian= 'T')
280  endif ; endif
281 
282  if (use_periodic) then
283  if (num_periods <= 1) call mom_error(fatal, "MOM_io create_file: "//&
284  "num_periods for file "//trim(filename)//" must be at least 1.")
285  ! Define a periodic axis with unit labels.
286  allocate(period_val(num_periods))
287  do k=1,num_periods ; period_val(k) = real(k) ; enddo
288  call mpp_write_meta(unit, axis_periodic, name="Period", units="nondimensional", &
289  longname="Periods for cyclical varaiables", cartesian= 't', data=period_val)
290  deallocate(period_val)
291  endif
292 
293  do k=1,novars
294  numaxes = 0
295  select case (vars(k)%hor_grid)
296  case ('h') ; numaxes = 2 ; axes(1) = axis_lonh ; axes(2) = axis_lath
297  case ('q') ; numaxes = 2 ; axes(1) = axis_lonq ; axes(2) = axis_latq
298  case ('u') ; numaxes = 2 ; axes(1) = axis_lonq ; axes(2) = axis_lath
299  case ('v') ; numaxes = 2 ; axes(1) = axis_lonh ; axes(2) = axis_latq
300  case ('T') ; numaxes = 2 ; axes(1) = axis_lonh ; axes(2) = axis_lath
301  case ('Bu') ; numaxes = 2 ; axes(1) = axis_lonq ; axes(2) = axis_latq
302  case ('Cu') ; numaxes = 2 ; axes(1) = axis_lonq ; axes(2) = axis_lath
303  case ('Cv') ; numaxes = 2 ; axes(1) = axis_lonh ; axes(2) = axis_latq
304  case ('1') ! Do nothing.
305  case default
306  call mom_error(warning, "MOM_io create_file: "//trim(vars(k)%name)//&
307  " has unrecognized hor_grid "//trim(vars(k)%hor_grid))
308  end select
309  select case (vars(k)%z_grid)
310  case ('L') ; numaxes = numaxes+1 ; axes(numaxes) = axis_layer
311  case ('i') ; numaxes = numaxes+1 ; axes(numaxes) = axis_int
312  case ('1') ! Do nothing.
313  case default
314  call mom_error(fatal, "MOM_io create_file: "//trim(vars(k)%name)//&
315  " has unrecognized z_grid "//trim(vars(k)%z_grid))
316  end select
317  t_grid = adjustl(vars(k)%t_grid)
318  select case (t_grid(1:1))
319  case ('s', 'a', 'm') ; numaxes = numaxes+1 ; axes(numaxes) = axis_time
320  case ('p') ; numaxes = numaxes+1 ; axes(numaxes) = axis_periodic
321  case ('1') ! Do nothing.
322  case default
323  call mom_error(warning, "MOM_io create_file: "//trim(vars(k)%name)//&
324  " has unrecognized t_grid "//trim(vars(k)%t_grid))
325  end select
326  pack = 1
327 
328  if (present(checksums)) then
329  call mpp_write_meta(unit, fields(k), axes(1:numaxes), vars(k)%name, vars(k)%units, &
330  vars(k)%longname, pack = pack, checksum=checksums(k,:))
331  else
332  call mpp_write_meta(unit, fields(k), axes(1:numaxes), vars(k)%name, vars(k)%units, &
333  vars(k)%longname, pack = pack)
334  endif
335  enddo
336 
337  if (use_lath) call write_field(unit, axis_lath)
338  if (use_latq) call write_field(unit, axis_latq)
339  if (use_lonh) call write_field(unit, axis_lonh)
340  if (use_lonq) call write_field(unit, axis_lonq)
341  if (use_layer) call write_field(unit, axis_layer)
342  if (use_int) call write_field(unit, axis_int)
343  if (use_periodic) call write_field(unit, axis_periodic)
344 

References mom_error_handler::mom_error().

Referenced by mom_ale::ale_writecoordinatefile(), and reopen_file().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ ensembler()

character(len=len(name)) function, public mom_io::ensembler ( character(len=*), intent(in)  name,
integer, intent(in), optional  ens_no_in 
)

Returns a name with "%#E" or "%E" replaced with the ensemble member number.

Parameters
[in]nameThe name to be modified
[in]ens_no_inThe number of the current ensemble member
Returns
The name encoded with the ensemble number

Definition at line 763 of file MOM_io.F90.

763  character(len=*), intent(in) :: name !< The name to be modified
764  integer, optional, intent(in) :: ens_no_in !< The number of the current ensemble member
765  character(len=len(name)) :: en_nm !< The name encoded with the ensemble number
766 
767  ! This function replaces "%#E" or "%E" with the ensemble number anywhere it
768  ! occurs in name, with %E using 4 or 6 digits (depending on the ensemble size)
769  ! and %#E using # digits, where # is a number from 1 to 9.
770 
771  character(len=len(name)) :: tmp
772  character(10) :: ens_num_char
773  character(3) :: code_str
774  integer :: ens_no
775  integer :: n, is, ie
776 
777  en_nm = trim(name)
778  if (index(name,"%") == 0) return
779 
780  if (present(ens_no_in)) then
781  ens_no = ens_no_in
782  else
783  ens_no = get_ensemble_id()
784  endif
785 
786  write(ens_num_char, '(I10)') ens_no ; ens_num_char = adjustl(ens_num_char)
787  do
788  is = index(en_nm,"%E")
789  if (is == 0) exit
790  if (len(en_nm) < len(trim(en_nm)) + len(trim(ens_num_char)) - 2) &
791  call mom_error(fatal, "MOM_io ensembler: name "//trim(name)// &
792  " is not long enough for %E expansion for ens_no "//trim(ens_num_char))
793  tmp = en_nm(1:is-1)//trim(ens_num_char)//trim(en_nm(is+2:))
794  en_nm = tmp
795  enddo
796 
797  if (index(name,"%") == 0) return
798 
799  write(ens_num_char, '(I10.10)') ens_no
800  do n=1,9 ; do
801  write(code_str, '("%",I1,"E")') n
802 
803  is = index(en_nm,code_str)
804  if (is == 0) exit
805  if (ens_no < 10**n) then
806  if (len(en_nm) < len(trim(en_nm)) + n-3) call mom_error(fatal, &
807  "MOM_io ensembler: name "//trim(name)//" is not long enough for %E expansion.")
808  tmp = en_nm(1:is-1)//trim(ens_num_char(11-n:10))//trim(en_nm(is+3:))
809  else
810  call mom_error(fatal, "MOM_io ensembler: Ensemble number is too large "//&
811  "to be encoded with "//code_str//" in "//trim(name))
812  endif
813  en_nm = tmp
814  enddo ; enddo
815 

References mom_error_handler::mom_error().

Here is the call graph for this function:

◆ fms_file_exists()

logical function mom_io::fms_file_exists ( character(len=*), intent(in)  filename,
type(domain2d), intent(in), optional  domain,
logical, intent(in), optional  no_domain 
)
private

Returns true if the named file or its domain-decomposed variant exists.

Parameters
[in]filenameThe name of the file being inquired about
[in]domainThe mpp domain2d that describes the decomposition
[in]no_domainThis file does not use domain decomposition

Definition at line 835 of file MOM_io.F90.

835  character(len=*), intent(in) :: filename !< The name of the file being inquired about
836  type(domain2d), optional, intent(in) :: domain !< The mpp domain2d that describes the decomposition
837  logical, optional, intent(in) :: no_domain !< This file does not use domain decomposition
838 ! This function uses the fms_io function file_exist to determine whether
839 ! a named file (or its decomposed variant) exists.
840 
841  logical :: FMS_file_exists
842 
843  fms_file_exists = file_exist(filename, domain, no_domain)
844 

◆ modify_vardesc()

subroutine, public mom_io::modify_vardesc ( type(vardesc), intent(inout)  vd,
character(len=*), intent(in), optional  name,
character(len=*), intent(in), optional  units,
character(len=*), intent(in), optional  longname,
character(len=*), intent(in), optional  hor_grid,
character(len=*), intent(in), optional  z_grid,
character(len=*), intent(in), optional  t_grid,
character(len=*), intent(in), optional  cmor_field_name,
character(len=*), intent(in), optional  cmor_units,
character(len=*), intent(in), optional  cmor_longname,
real, intent(in), optional  conversion,
character(len=*), intent(in), optional  caller 
)

This routine modifies the named elements of a vardesc type. All arguments are optional, except the vardesc type to be modified.

Parameters
[in,out]vdvardesc type that is modified
[in]namename of variable
[in]unitsunits of variable
[in]longnamelong name of variable
[in]hor_gridhorizonal staggering of variable
[in]z_gridvertical staggering of variable
[in]t_gridtime description: s, p, or 1
[in]cmor_field_nameCMOR name
[in]cmor_unitsCMOR physical dimensions of variable
[in]cmor_longnameCMOR long name
[in]conversionfor unit conversions, such as needed to convert from intensive to extensive
[in]callercalling routine?

Definition at line 640 of file MOM_io.F90.

640  type(vardesc), intent(inout) :: vd !< vardesc type that is modified
641  character(len=*), optional, intent(in) :: name !< name of variable
642  character(len=*), optional, intent(in) :: units !< units of variable
643  character(len=*), optional, intent(in) :: longname !< long name of variable
644  character(len=*), optional, intent(in) :: hor_grid !< horizonal staggering of variable
645  character(len=*), optional, intent(in) :: z_grid !< vertical staggering of variable
646  character(len=*), optional, intent(in) :: t_grid !< time description: s, p, or 1
647  character(len=*), optional, intent(in) :: cmor_field_name !< CMOR name
648  character(len=*), optional, intent(in) :: cmor_units !< CMOR physical dimensions of variable
649  character(len=*), optional, intent(in) :: cmor_longname !< CMOR long name
650  real , optional, intent(in) :: conversion !< for unit conversions, such as needed
651  !! to convert from intensive to extensive
652  character(len=*), optional, intent(in) :: caller !< calling routine?
653 
654  character(len=120) :: cllr
655  cllr = "mod_vardesc"
656  if (present(caller)) cllr = trim(caller)
657 
658  if (present(name)) call safe_string_copy(name, vd%name, "vd%name", cllr)
659 
660  if (present(longname)) call safe_string_copy(longname, vd%longname, &
661  "vd%longname of "//trim(vd%name), cllr)
662  if (present(units)) call safe_string_copy(units, vd%units, &
663  "vd%units of "//trim(vd%name), cllr)
664  if (present(hor_grid)) call safe_string_copy(hor_grid, vd%hor_grid, &
665  "vd%hor_grid of "//trim(vd%name), cllr)
666  if (present(z_grid)) call safe_string_copy(z_grid, vd%z_grid, &
667  "vd%z_grid of "//trim(vd%name), cllr)
668  if (present(t_grid)) call safe_string_copy(t_grid, vd%t_grid, &
669  "vd%t_grid of "//trim(vd%name), cllr)
670 
671  if (present(cmor_field_name)) call safe_string_copy(cmor_field_name, vd%cmor_field_name, &
672  "vd%cmor_field_name of "//trim(vd%name), cllr)
673  if (present(cmor_units)) call safe_string_copy(cmor_units, vd%cmor_units, &
674  "vd%cmor_units of "//trim(vd%name), cllr)
675  if (present(cmor_longname)) call safe_string_copy(cmor_longname, vd%cmor_longname, &
676  "vd%cmor_longname of "//trim(vd%name), cllr)
677 

References safe_string_copy().

Referenced by var_desc().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ mom_file_exists()

logical function mom_io::mom_file_exists ( character(len=*), intent(in)  filename,
type(mom_domain_type), intent(in)  MOM_Domain 
)
private

Returns true if the named file or its domain-decomposed variant exists.

Parameters
[in]filenameThe name of the file being inquired about
[in]mom_domainThe MOM_Domain that describes the decomposition

Definition at line 821 of file MOM_io.F90.

821  character(len=*), intent(in) :: filename !< The name of the file being inquired about
822  type(MOM_domain_type), intent(in) :: MOM_Domain !< The MOM_Domain that describes the decomposition
823 
824 ! This function uses the fms_io function file_exist to determine whether
825 ! a named file (or its decomposed variant) exists.
826 
827  logical :: MOM_file_exists
828 
829  mom_file_exists = file_exist(filename, mom_domain%mpp_domain)
830 

◆ mom_io_init()

subroutine, public mom_io::mom_io_init ( type(param_file_type), intent(in)  param_file)

Initialize the MOM_io module.

Parameters
[in]param_filestructure indicating the open file to parse for model parameter values.

Definition at line 1043 of file MOM_io.F90.

1043  type(param_file_type), intent(in) :: param_file !< structure indicating the open file to
1044  !! parse for model parameter values.
1045 
1046 ! This include declares and sets the variable "version".
1047 #include "version_variable.h"
1048  character(len=40) :: mdl = "MOM_io" ! This module's name.
1049 
1050  call log_version(param_file, mdl, version)
1051 

◆ mom_read_data_1d()

subroutine mom_io::mom_read_data_1d ( character(len=*), intent(in)  filename,
character(len=*), intent(in)  fieldname,
real, dimension(:), intent(inout)  data,
integer, intent(in), optional  timelevel,
real, intent(in), optional  scale 
)
private

This function uses the fms_io function read_data to read 1-D data field named "fieldname" from file "filename".

Parameters
[in]filenameThe name of the file to read
[in]fieldnameThe variable name of the data in the file
[in,out]dataThe 1-dimensional array into which the data
[in]timelevelThe time level in the file to read
[in]scaleA scaling factor that the field is multiplied by before they are returned.

Definition at line 850 of file MOM_io.F90.

850  character(len=*), intent(in) :: filename !< The name of the file to read
851  character(len=*), intent(in) :: fieldname !< The variable name of the data in the file
852  real, dimension(:), intent(inout) :: data !< The 1-dimensional array into which the data
853  integer, optional, intent(in) :: timelevel !< The time level in the file to read
854  real, optional, intent(in) :: scale !< A scaling factor that the field is multiplied
855  !! by before they are returned.
856 
857  call read_data(filename, fieldname, data, timelevel=timelevel, no_domain=.true.)
858 
859  if (present(scale)) then ; if (scale /= 1.0) then
860  data(:) = scale*data(:)
861  endif ; endif
862 

◆ mom_read_data_2d()

subroutine mom_io::mom_read_data_2d ( character(len=*), intent(in)  filename,
character(len=*), intent(in)  fieldname,
real, dimension(:,:), intent(inout)  data,
type(mom_domain_type), intent(in)  MOM_Domain,
integer, intent(in), optional  timelevel,
integer, intent(in), optional  position,
real, intent(in), optional  scale 
)
private

This function uses the fms_io function read_data to read a distributed 2-D data field named "fieldname" from file "filename". Valid values for "position" include CORNER, CENTER, EAST_FACE and NORTH_FACE.

Parameters
[in]filenameThe name of the file to read
[in]fieldnameThe variable name of the data in the file
[in,out]dataThe 2-dimensional array into which the data should be read
[in]mom_domainThe MOM_Domain that describes the decomposition
[in]timelevelThe time level in the file to read
[in]positionA flag indicating where this data is located
[in]scaleA scaling factor that the field is multiplied by before they are returned.

Definition at line 870 of file MOM_io.F90.

870  character(len=*), intent(in) :: filename !< The name of the file to read
871  character(len=*), intent(in) :: fieldname !< The variable name of the data in the file
872  real, dimension(:,:), intent(inout) :: data !< The 2-dimensional array into which the data
873  !! should be read
874  type(MOM_domain_type), intent(in) :: MOM_Domain !< The MOM_Domain that describes the decomposition
875  integer, optional, intent(in) :: timelevel !< The time level in the file to read
876  integer, optional, intent(in) :: position !< A flag indicating where this data is located
877  real, optional, intent(in) :: scale !< A scaling factor that the field is multiplied
878  !! by before they are returned.
879 
880  integer :: is, ie, js, je
881 
882  call read_data(filename, fieldname, data, mom_domain%mpp_domain, &
883  timelevel=timelevel, position=position)
884 
885  if (present(scale)) then ; if (scale /= 1.0) then
886  call get_simple_array_i_ind(mom_domain, size(data,1), is, ie)
887  call get_simple_array_j_ind(mom_domain, size(data,2), js, je)
888  data(is:ie,js:je) = scale*data(is:ie,js:je)
889  endif ; endif
890 

References mom_domains::get_simple_array_i_ind(), and mom_domains::get_simple_array_j_ind().

Here is the call graph for this function:

◆ mom_read_data_3d()

subroutine mom_io::mom_read_data_3d ( character(len=*), intent(in)  filename,
character(len=*), intent(in)  fieldname,
real, dimension(:,:,:), intent(inout)  data,
type(mom_domain_type), intent(in)  MOM_Domain,
integer, intent(in), optional  timelevel,
integer, intent(in), optional  position,
real, intent(in), optional  scale 
)
private

This function uses the fms_io function read_data to read a distributed 3-D data field named "fieldname" from file "filename". Valid values for "position" include CORNER, CENTER, EAST_FACE and NORTH_FACE.

Parameters
[in]filenameThe name of the file to read
[in]fieldnameThe variable name of the data in the file
[in,out]dataThe 3-dimensional array into which the data should be read
[in]mom_domainThe MOM_Domain that describes the decomposition
[in]timelevelThe time level in the file to read
[in]positionA flag indicating where this data is located
[in]scaleA scaling factor that the field is multiplied by before they are returned.

Definition at line 898 of file MOM_io.F90.

898  character(len=*), intent(in) :: filename !< The name of the file to read
899  character(len=*), intent(in) :: fieldname !< The variable name of the data in the file
900  real, dimension(:,:,:), intent(inout) :: data !< The 3-dimensional array into which the data
901  !! should be read
902  type(MOM_domain_type), intent(in) :: MOM_Domain !< The MOM_Domain that describes the decomposition
903  integer, optional, intent(in) :: timelevel !< The time level in the file to read
904  integer, optional, intent(in) :: position !< A flag indicating where this data is located
905  real, optional, intent(in) :: scale !< A scaling factor that the field is multiplied
906  !! by before they are returned.
907 
908  integer :: is, ie, js, je
909 
910  call read_data(filename, fieldname, data, mom_domain%mpp_domain, &
911  timelevel=timelevel, position=position)
912 
913  if (present(scale)) then ; if (scale /= 1.0) then
914  call get_simple_array_i_ind(mom_domain, size(data,1), is, ie)
915  call get_simple_array_j_ind(mom_domain, size(data,2), js, je)
916  data(is:ie,js:je,:) = scale*data(is:ie,js:je,:)
917  endif ; endif
918 

References mom_domains::get_simple_array_i_ind(), and mom_domains::get_simple_array_j_ind().

Here is the call graph for this function:

◆ mom_read_data_4d()

subroutine mom_io::mom_read_data_4d ( character(len=*), intent(in)  filename,
character(len=*), intent(in)  fieldname,
real, dimension(:,:,:,:), intent(inout)  data,
type(mom_domain_type), intent(in)  MOM_Domain,
integer, intent(in), optional  timelevel,
integer, intent(in), optional  position,
real, intent(in), optional  scale 
)
private

This function uses the fms_io function read_data to read a distributed 4-D data field named "fieldname" from file "filename". Valid values for "position" include CORNER, CENTER, EAST_FACE and NORTH_FACE.

Parameters
[in]filenameThe name of the file to read
[in]fieldnameThe variable name of the data in the file
[in,out]dataThe 4-dimensional array into which the data should be read
[in]mom_domainThe MOM_Domain that describes the decomposition
[in]timelevelThe time level in the file to read
[in]positionA flag indicating where this data is located
[in]scaleA scaling factor that the field is multiplied by before they are returned.

Definition at line 926 of file MOM_io.F90.

926  character(len=*), intent(in) :: filename !< The name of the file to read
927  character(len=*), intent(in) :: fieldname !< The variable name of the data in the file
928  real, dimension(:,:,:,:), intent(inout) :: data !< The 4-dimensional array into which the data
929  !! should be read
930  type(MOM_domain_type), intent(in) :: MOM_Domain !< The MOM_Domain that describes the decomposition
931  integer, optional, intent(in) :: timelevel !< The time level in the file to read
932  integer, optional, intent(in) :: position !< A flag indicating where this data is located
933  real, optional, intent(in) :: scale !< A scaling factor that the field is multiplied
934  !! by before they are returned.
935 
936  integer :: is, ie, js, je
937 
938  call read_data(filename, fieldname, data, mom_domain%mpp_domain, &
939  timelevel=timelevel, position=position)
940 
941  if (present(scale)) then ; if (scale /= 1.0) then
942  call get_simple_array_i_ind(mom_domain, size(data,1), is, ie)
943  call get_simple_array_j_ind(mom_domain, size(data,2), js, je)
944  data(is:ie,js:je,:,:) = scale*data(is:ie,js:je,:,:)
945  endif ; endif
946 

References mom_domains::get_simple_array_i_ind(), and mom_domains::get_simple_array_j_ind().

Here is the call graph for this function:

◆ mom_read_vector_2d()

subroutine mom_io::mom_read_vector_2d ( character(len=*), intent(in)  filename,
character(len=*), intent(in)  u_fieldname,
character(len=*), intent(in)  v_fieldname,
real, dimension(:,:), intent(inout)  u_data,
real, dimension(:,:), intent(inout)  v_data,
type(mom_domain_type), intent(in)  MOM_Domain,
integer, intent(in), optional  timelevel,
integer, intent(in), optional  stagger,
logical, intent(in), optional  scalar_pair,
real, intent(in), optional  scale 
)
private

This function uses the fms_io function read_data to read a pair of distributed 2-D data fields with names given by "[uv]_fieldname" from file "filename". Valid values for "stagger" include CGRID_NE, BGRID_NE, and AGRID.

Parameters
[in]filenameThe name of the file to read
[in]u_fieldnameThe variable name of the u data in the file
[in]v_fieldnameThe variable name of the v data in the file
[in,out]u_dataThe 2 dimensional array into which the u-component of the data should be read
[in,out]v_dataThe 2 dimensional array into which the v-component of the data should be read
[in]mom_domainThe MOM_Domain that describes the decomposition
[in]timelevelThe time level in the file to read
[in]staggerA flag indicating where this vector is discretized
[in]scalar_pairIf true, a pair of scalars are to be read.cretized
[in]scaleA scaling factor that the fields are multiplied by before they are returned.

Definition at line 955 of file MOM_io.F90.

955  character(len=*), intent(in) :: filename !< The name of the file to read
956  character(len=*), intent(in) :: u_fieldname !< The variable name of the u data in the file
957  character(len=*), intent(in) :: v_fieldname !< The variable name of the v data in the file
958  real, dimension(:,:), intent(inout) :: u_data !< The 2 dimensional array into which the
959  !! u-component of the data should be read
960  real, dimension(:,:), intent(inout) :: v_data !< The 2 dimensional array into which the
961  !! v-component of the data should be read
962  type(MOM_domain_type), intent(in) :: MOM_Domain !< The MOM_Domain that describes the decomposition
963  integer, optional, intent(in) :: timelevel !< The time level in the file to read
964  integer, optional, intent(in) :: stagger !< A flag indicating where this vector is discretized
965  logical, optional, intent(in) :: scalar_pair !< If true, a pair of scalars are to be read.cretized
966  real, optional, intent(in) :: scale !< A scaling factor that the fields are multiplied
967  !! by before they are returned.
968  integer :: is, ie, js, je
969  integer :: u_pos, v_pos
970 
971  u_pos = east_face ; v_pos = north_face
972  if (present(stagger)) then
973  if (stagger == cgrid_ne) then ; u_pos = east_face ; v_pos = north_face
974  elseif (stagger == bgrid_ne) then ; u_pos = corner ; v_pos = corner
975  elseif (stagger == agrid) then ; u_pos = center ; v_pos = center ; endif
976  endif
977 
978  call read_data(filename, u_fieldname, u_data, mom_domain%mpp_domain, &
979  timelevel=timelevel, position=u_pos)
980  call read_data(filename, v_fieldname, v_data, mom_domain%mpp_domain, &
981  timelevel=timelevel, position=v_pos)
982 
983  if (present(scale)) then ; if (scale /= 1.0) then
984  call get_simple_array_i_ind(mom_domain, size(u_data,1), is, ie)
985  call get_simple_array_j_ind(mom_domain, size(u_data,2), js, je)
986  u_data(is:ie,js:je) = scale*u_data(is:ie,js:je)
987  call get_simple_array_i_ind(mom_domain, size(v_data,1), is, ie)
988  call get_simple_array_j_ind(mom_domain, size(v_data,2), js, je)
989  v_data(is:ie,js:je) = scale*v_data(is:ie,js:je)
990  endif ; endif
991 

References mom_domains::get_simple_array_i_ind(), and mom_domains::get_simple_array_j_ind().

Here is the call graph for this function:

◆ mom_read_vector_3d()

subroutine mom_io::mom_read_vector_3d ( character(len=*), intent(in)  filename,
character(len=*), intent(in)  u_fieldname,
character(len=*), intent(in)  v_fieldname,
real, dimension(:,:,:), intent(inout)  u_data,
real, dimension(:,:,:), intent(inout)  v_data,
type(mom_domain_type), intent(in)  MOM_Domain,
integer, intent(in), optional  timelevel,
integer, intent(in), optional  stagger,
logical, intent(in), optional  scalar_pair,
real, intent(in), optional  scale 
)
private

This function uses the fms_io function read_data to read a pair of distributed 3-D data fields with names given by "[uv]_fieldname" from file "filename". Valid values for "stagger" include CGRID_NE, BGRID_NE, and AGRID.

Parameters
[in]filenameThe name of the file to read
[in]u_fieldnameThe variable name of the u data in the file
[in]v_fieldnameThe variable name of the v data in the file
[in,out]u_dataThe 3 dimensional array into which the u-component of the data should be read
[in,out]v_dataThe 3 dimensional array into which the v-component of the data should be read
[in]mom_domainThe MOM_Domain that describes the decomposition
[in]timelevelThe time level in the file to read
[in]staggerA flag indicating where this vector is discretized
[in]scalar_pairIf true, a pair of scalars are to be read.cretized
[in]scaleA scaling factor that the fields are multiplied by before they are returned.

Definition at line 1000 of file MOM_io.F90.

1000  character(len=*), intent(in) :: filename !< The name of the file to read
1001  character(len=*), intent(in) :: u_fieldname !< The variable name of the u data in the file
1002  character(len=*), intent(in) :: v_fieldname !< The variable name of the v data in the file
1003  real, dimension(:,:,:), intent(inout) :: u_data !< The 3 dimensional array into which the
1004  !! u-component of the data should be read
1005  real, dimension(:,:,:), intent(inout) :: v_data !< The 3 dimensional array into which the
1006  !! v-component of the data should be read
1007  type(MOM_domain_type), intent(in) :: MOM_Domain !< The MOM_Domain that describes the decomposition
1008  integer, optional, intent(in) :: timelevel !< The time level in the file to read
1009  integer, optional, intent(in) :: stagger !< A flag indicating where this vector is discretized
1010  logical, optional, intent(in) :: scalar_pair !< If true, a pair of scalars are to be read.cretized
1011  real, optional, intent(in) :: scale !< A scaling factor that the fields are multiplied
1012  !! by before they are returned.
1013 
1014  integer :: is, ie, js, je
1015  integer :: u_pos, v_pos
1016 
1017  u_pos = east_face ; v_pos = north_face
1018  if (present(stagger)) then
1019  if (stagger == cgrid_ne) then ; u_pos = east_face ; v_pos = north_face
1020  elseif (stagger == bgrid_ne) then ; u_pos = corner ; v_pos = corner
1021  elseif (stagger == agrid) then ; u_pos = center ; v_pos = center ; endif
1022  endif
1023 
1024  call read_data(filename, u_fieldname, u_data, mom_domain%mpp_domain, &
1025  timelevel=timelevel, position=u_pos)
1026  call read_data(filename, v_fieldname, v_data, mom_domain%mpp_domain, &
1027  timelevel=timelevel, position=v_pos)
1028 
1029  if (present(scale)) then ; if (scale /= 1.0) then
1030  call get_simple_array_i_ind(mom_domain, size(u_data,1), is, ie)
1031  call get_simple_array_j_ind(mom_domain, size(u_data,2), js, je)
1032  u_data(is:ie,js:je,:) = scale*u_data(is:ie,js:je,:)
1033  call get_simple_array_i_ind(mom_domain, size(v_data,1), is, ie)
1034  call get_simple_array_j_ind(mom_domain, size(v_data,2), js, je)
1035  v_data(is:ie,js:je,:) = scale*v_data(is:ie,js:je,:)
1036  endif ; endif
1037 

References mom_domains::get_simple_array_i_ind(), and mom_domains::get_simple_array_j_ind().

Here is the call graph for this function:

◆ num_timelevels()

integer function, public mom_io::num_timelevels ( character(len=*), intent(in)  filename,
character(len=*), intent(in)  varname,
integer, intent(in), optional  min_dims 
)

This function determines how many time levels a variable has.

Parameters
[in]filenamename of the file to read
[in]varnamevariable whose number of time levels are to be returned
[in]min_dimsThe minimum number of dimensions a variable must have if it has a time dimension. If the variable has 1 less dimension than this, then 0 is returned.
Returns
number of time levels varname has in filename

Definition at line 478 of file MOM_io.F90.

478  character(len=*), intent(in) :: filename !< name of the file to read
479  character(len=*), intent(in) :: varname !< variable whose number of time levels
480  !! are to be returned
481  integer, optional, intent(in) :: min_dims !< The minimum number of dimensions a variable must have
482  !! if it has a time dimension. If the variable has 1 less
483  !! dimension than this, then 0 is returned.
484  integer :: n_time !< number of time levels varname has in filename
485 
486  logical :: found
487  character(len=200) :: msg
488  character(len=nf90_max_name) :: name
489  integer :: ncid, nvars, status, varid, ndims, n
490  integer, allocatable :: varids(:)
491  integer, dimension(nf90_max_var_dims) :: dimids
492 
493  n_time = -1
494  found = .false.
495 
496  status = nf90_open(filename, nf90_nowrite, ncid)
497  if (status /= nf90_noerr) then
498  call mom_error(warning,"num_timelevels: "//&
499  " Difficulties opening "//trim(filename)//" - "//&
500  trim(nf90_strerror(status)))
501  return
502  endif
503 
504  status = nf90_inquire(ncid, nvariables=nvars)
505  if (status /= nf90_noerr) then
506  call mom_error(warning,"num_timelevels: "//&
507  " Difficulties getting the number of variables in file "//&
508  trim(filename)//" - "//trim(nf90_strerror(status)))
509  return
510  endif
511 
512  if (nvars < 1) then
513  call mom_error(warning,"num_timelevels: "//&
514  " There appear not to be any variables in "//trim(filename))
515  return
516  endif
517 
518 
519  allocate(varids(nvars))
520 
521  status = nf90_inq_varids(ncid, nvars, varids)
522  if (status /= nf90_noerr) then
523  call mom_error(warning,"num_timelevels: "//&
524  " Difficulties getting the variable IDs in file "//&
525  trim(filename)//" - "//trim(nf90_strerror(status)))
526  deallocate(varids) ; return
527  endif
528 
529  do n = 1,nvars
530  status = nf90_inquire_variable(ncid, varids(n), name=name)
531  if (status /= nf90_noerr) then
532  call mom_error(warning,"num_timelevels: "//&
533  " Difficulties getting a variable name in file "//&
534  trim(filename)//" - "//trim(nf90_strerror(status)))
535  endif
536 
537  if (trim(lowercase(name)) == trim(lowercase(varname))) then
538  if (found) then
539  call mom_error(warning,"num_timelevels: "//&
540  " Two variables match the case-insensitive name "//trim(varname)//&
541  " in file "//trim(filename)//" - "//trim(nf90_strerror(status)))
542  else
543  varid = varids(n) ; found = .true.
544  endif
545  endif
546  enddo
547 
548  deallocate(varids)
549 
550  if (.not.found) then
551  call mom_error(warning,"num_timelevels: "//&
552  " variable "//trim(varname)//" was not found in file "//&
553  trim(filename))
554  return
555  endif
556 
557  status = nf90_inquire_variable(ncid, varid, ndims = ndims)
558  if (status /= nf90_noerr) then
559  call mom_error(warning,"num_timelevels: "//&
560  trim(nf90_strerror(status))//" Getting number of dimensions of "//&
561  trim(varname)//" in "//trim(filename))
562  return
563  endif
564 
565  if (present(min_dims)) then
566  if (ndims < min_dims-1) then
567  write(msg, '(I3)') min_dims
568  call mom_error(warning, "num_timelevels: variable "//trim(varname)//&
569  " in file "//trim(filename)//" has fewer than min_dims = "//trim(msg)//&
570  " dimensions.")
571  elseif (ndims == min_dims - 1) then
572  n_time = 0 ; return
573  endif
574  endif
575 
576  status = nf90_inquire_variable(ncid, varid, dimids = dimids(1:ndims))
577  if (status /= nf90_noerr) then
578  call mom_error(warning,"num_timelevels: "//&
579  trim(nf90_strerror(status))//" Getting last dimension ID for "//&
580  trim(varname)//" in "//trim(filename))
581  return
582  endif
583 
584  status = nf90_inquire_dimension(ncid, dimids(ndims), len=n_time)
585  if (status /= nf90_noerr) call mom_error(warning,"num_timelevels: "//&
586  trim(nf90_strerror(status))//" Getting number of time levels of "//&
587  trim(varname)//" in "//trim(filename))
588 
589  return
590 

References mom_string_functions::lowercase(), and mom_error_handler::mom_error().

Referenced by mom_surface_forcing::surface_forcing_init().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ query_vardesc()

subroutine, public mom_io::query_vardesc ( type(vardesc), intent(in)  vd,
character(len=*), intent(out), optional  name,
character(len=*), intent(out), optional  units,
character(len=*), intent(out), optional  longname,
character(len=*), intent(out), optional  hor_grid,
character(len=*), intent(out), optional  z_grid,
character(len=*), intent(out), optional  t_grid,
character(len=*), intent(out), optional  cmor_field_name,
character(len=*), intent(out), optional  cmor_units,
character(len=*), intent(out), optional  cmor_longname,
real, intent(out), optional  conversion,
character(len=*), intent(in), optional  caller 
)

This routine queries vardesc.

Parameters
[in]vdvardesc type that is queried
[out]namename of variable
[out]unitsunits of variable
[out]longnamelong name of variable
[out]hor_gridhoriz staggering of variable
[out]z_gridvert staggering of variable
[out]t_gridtime description: s, p, or 1
[out]cmor_field_nameCMOR name
[out]cmor_unitsCMOR physical dimensions of variable
[out]cmor_longnameCMOR long name
[out]conversionfor unit conversions, such as needed to convert from intensive to extensive
[in]callercalling routine?

Definition at line 699 of file MOM_io.F90.

699  type(vardesc), intent(in) :: vd !< vardesc type that is queried
700  character(len=*), optional, intent(out) :: name !< name of variable
701  character(len=*), optional, intent(out) :: units !< units of variable
702  character(len=*), optional, intent(out) :: longname !< long name of variable
703  character(len=*), optional, intent(out) :: hor_grid !< horiz staggering of variable
704  character(len=*), optional, intent(out) :: z_grid !< vert staggering of variable
705  character(len=*), optional, intent(out) :: t_grid !< time description: s, p, or 1
706  character(len=*), optional, intent(out) :: cmor_field_name !< CMOR name
707  character(len=*), optional, intent(out) :: cmor_units !< CMOR physical dimensions of variable
708  character(len=*), optional, intent(out) :: cmor_longname !< CMOR long name
709  real , optional, intent(out) :: conversion !< for unit conversions, such as needed to
710  !! convert from intensive to extensive
711  character(len=*), optional, intent(in) :: caller !< calling routine?
712 
713 
714  character(len=120) :: cllr
715  cllr = "mod_vardesc"
716  if (present(caller)) cllr = trim(caller)
717 
718  if (present(name)) call safe_string_copy(vd%name, name, &
719  "vd%name of "//trim(vd%name), cllr)
720  if (present(longname)) call safe_string_copy(vd%longname, longname, &
721  "vd%longname of "//trim(vd%name), cllr)
722  if (present(units)) call safe_string_copy(vd%units, units, &
723  "vd%units of "//trim(vd%name), cllr)
724  if (present(hor_grid)) call safe_string_copy(vd%hor_grid, hor_grid, &
725  "vd%hor_grid of "//trim(vd%name), cllr)
726  if (present(z_grid)) call safe_string_copy(vd%z_grid, z_grid, &
727  "vd%z_grid of "//trim(vd%name), cllr)
728  if (present(t_grid)) call safe_string_copy(vd%t_grid, t_grid, &
729  "vd%t_grid of "//trim(vd%name), cllr)
730 
731  if (present(cmor_field_name)) call safe_string_copy(vd%cmor_field_name, cmor_field_name, &
732  "vd%cmor_field_name of "//trim(vd%name), cllr)
733  if (present(cmor_units)) call safe_string_copy(vd%cmor_units, cmor_units, &
734  "vd%cmor_units of "//trim(vd%name), cllr)
735  if (present(cmor_longname)) call safe_string_copy(vd%cmor_longname, cmor_longname, &
736  "vd%cmor_longname of "//trim(vd%name), cllr)
737 

References safe_string_copy().

Referenced by advection_test_tracer::advection_test_stock(), boundary_impulse_tracer::boundary_impulse_stock(), regional_dyes::dye_stock(), ideal_age_example::ideal_age_stock(), advection_test_tracer::initialize_advection_test_tracer(), boundary_impulse_tracer::initialize_boundary_impulse_tracer(), dome_tracer::initialize_dome_tracer(), dyed_obc_tracer::initialize_dyed_obc_tracer(), ideal_age_example::initialize_ideal_age_tracer(), isomip_tracer::initialize_isomip_tracer(), oil_tracer::initialize_oil_tracer(), pseudo_salt_tracer::initialize_pseudo_salt_tracer(), rgc_tracer::initialize_rgc_tracer(), mom_ocmip2_cfc::ocmip2_cfc_stock(), oil_tracer::oil_stock(), pseudo_salt_tracer::pseudo_salt_stock(), boundary_impulse_tracer::register_boundary_impulse_tracer(), regional_dyes::register_dye_tracer(), ideal_age_example::register_ideal_age_tracer(), oil_tracer::register_oil_tracer(), pseudo_salt_tracer::register_pseudo_salt_tracer(), mom_tracer_registry::register_tracer(), user_tracer_example::user_initialize_tracer(), and user_tracer_example::user_tracer_stock().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ read_axis_data()

subroutine, public mom_io::read_axis_data ( character(len=*), intent(in)  filename,
character(len=*), intent(in)  axis_name,
real, dimension(:), intent(out)  var 
)

Read the data associated with a named axis in a file.

Parameters
[in]filenameName of the file to read
[in]axis_nameName of the axis to read
[out]varThe axis location data

Definition at line 438 of file MOM_io.F90.

438  character(len=*), intent(in) :: filename !< Name of the file to read
439  character(len=*), intent(in) :: axis_name !< Name of the axis to read
440  real, dimension(:), intent(out) :: var !< The axis location data
441 
442  integer :: i,len,unit, ndim, nvar, natt, ntime
443  logical :: axis_found
444  type(axistype), allocatable :: axes(:)
445  type(axistype) :: time_axis
446  character(len=32) :: name, units
447 
448  call open_file(unit, trim(filename), action=mpp_rdonly, form=mpp_netcdf, &
449  threading=mpp_multi, fileset=single_file)
450 
451 !Find the number of variables (nvar) in this file
452  call mpp_get_info(unit, ndim, nvar, natt, ntime)
453 ! -------------------------------------------------------------------
454 ! Allocate space for the number of axes in the data file.
455 ! -------------------------------------------------------------------
456  allocate(axes(ndim))
457  call mpp_get_axes(unit, axes, time_axis)
458 
459  axis_found = .false.
460  do i = 1, ndim
461  call mpp_get_atts(axes(i), name=name,len=len,units=units)
462  if (name == axis_name) then
463  axis_found = .true.
464  call get_axis_data(axes(i),var)
465  exit
466  endif
467  enddo
468 
469  if (.not.axis_found) call mom_error(fatal, "MOM_io read_axis_data: "//&
470  "Unable to find axis "//trim(axis_name)//" in file "//trim(filename))
471 
472  deallocate(axes)
473 

References mom_error_handler::mom_error().

Here is the call graph for this function:

◆ reopen_file()

subroutine, public mom_io::reopen_file ( integer, intent(out)  unit,
character(len=*), intent(in)  filename,
type(vardesc), dimension(:), intent(in)  vars,
integer, intent(in)  novars,
type(fieldtype), dimension(:), intent(inout)  fields,
integer, intent(in), optional  threading,
real, intent(in), optional  timeunit,
type(ocean_grid_type), intent(in), optional  G,
type(dyn_horgrid_type), intent(in), optional  dG,
type(verticalgrid_type), intent(in), optional  GV 
)

This routine opens an existing NetCDF file for output. If it does not find the file, a new file is created. It also sets up structures that describe this file and the variables that will later be written to this file.

Parameters
[out]unitunit id of an open file or -1 on a nonwriting PE with single file output
[in]filenamefull path to the file to create
[in]varsstructures describing fields written to filename
[in]novarsnumber of fields written to filename
[in,out]fieldsarray of fieldtypes for each variable
[in]threadingSINGLE_FILE or MULTIPLE
[in]timeunitlength of the units for time [s]. The default value is 86400.0, for 1 day.
[in]gocean horizontal grid structure; G or dG is required if a new file uses any horizontal grid axes.
[in]dgdynamic horizontal grid structure; G or dG is required if a new file uses any horizontal grid axes.
[in]gvocean vertical grid structure, which is required if a new file uses any vertical grid axes.

Definition at line 353 of file MOM_io.F90.

353  integer, intent(out) :: unit !< unit id of an open file or -1 on a
354  !! nonwriting PE with single file output
355  character(len=*), intent(in) :: filename !< full path to the file to create
356  type(vardesc), intent(in) :: vars(:) !< structures describing fields written to filename
357  integer, intent(in) :: novars !< number of fields written to filename
358  type(fieldtype), intent(inout) :: fields(:) !< array of fieldtypes for each variable
359  integer, optional, intent(in) :: threading !< SINGLE_FILE or MULTIPLE
360  real, optional, intent(in) :: timeunit !< length of the units for time [s]. The
361  !! default value is 86400.0, for 1 day.
362  type(ocean_grid_type), optional, intent(in) :: G !< ocean horizontal grid structure; G or dG
363  !! is required if a new file uses any
364  !! horizontal grid axes.
365  type(dyn_horgrid_type), optional, intent(in) :: dG !< dynamic horizontal grid structure; G or dG
366  !! is required if a new file uses any
367  !! horizontal grid axes.
368  type(verticalGrid_type), optional, intent(in) :: GV !< ocean vertical grid structure, which is
369  !! required if a new file uses any
370  !! vertical grid axes.
371 
372  type(MOM_domain_type), pointer :: Domain => null()
373  character(len=200) :: check_name, mesg
374  integer :: length, ndim, nvar, natt, ntime, thread
375  logical :: exists, one_file, domain_set
376 
377  thread = single_file
378  if (PRESENT(threading)) thread = threading
379 
380  check_name = filename
381  length = len(trim(check_name))
382  if (check_name(length-2:length) /= ".nc") check_name = trim(check_name)//".nc"
383  if (thread /= single_file) check_name = trim(check_name)//".0000"
384 
385  inquire(file=check_name,exist=exists)
386 
387  if (.not.exists) then
388  call create_file(unit, filename, vars, novars, fields, threading, timeunit, &
389  g=g, dg=dg, gv=gv)
390  else
391 
392  domain_set = .false.
393  if (present(g)) then
394  domain_set = .true. ; domain => g%Domain
395  elseif (present(dg)) then
396  domain_set = .true. ; domain => dg%Domain
397  endif
398 
399  one_file = .true.
400  if (domain_set) one_file = (thread == single_file)
401 
402  if (one_file) then
403  call open_file(unit, filename, mpp_append, mpp_netcdf, threading=thread)
404  else
405  call open_file(unit, filename, mpp_append, mpp_netcdf, domain=domain%mpp_domain)
406  endif
407  if (unit < 0) return
408 
409  call mpp_get_info(unit, ndim, nvar, natt, ntime)
410 
411  if (nvar == -1) then
412  write (mesg,*) "Reopening file ",trim(filename)," apparently had ",nvar,&
413  " variables. Clobbering and creating file with ",novars," instead."
414  call mom_error(warning,"MOM_io: "//mesg)
415  call create_file(unit, filename, vars, novars, fields, threading, timeunit, g=g, gv=gv)
416  elseif (nvar /= novars) then
417  write (mesg,*) "Reopening file ",trim(filename)," with ",novars,&
418  " variables instead of ",nvar,"."
419  call mom_error(fatal,"MOM_io: "//mesg)
420  endif
421 
422  if (nvar>0) call mpp_get_fields(unit,fields(1:nvar))
423 
424  ! Check the field names...
425 ! do i=1,nvar
426 ! call mpp_get_field_atts(fields(i),name)
427 ! !if (trim(name) /= trim(vars%name) then
428 ! !write (mesg,'("Reopening file ",a," variable ",a," is called ",a,".")',&
429 ! ! filename,vars%name,name)
430 ! !call MOM_error(NOTE,"MOM_io: "//mesg)
431 ! enddo
432  endif
433 

References create_file(), and mom_error_handler::mom_error().

Here is the call graph for this function:

◆ safe_string_copy()

subroutine mom_io::safe_string_copy ( character(len=*), intent(in)  str1,
character(len=*), intent(out)  str2,
character(len=*), intent(in), optional  fieldnm,
character(len=*), intent(in), optional  caller 
)
private

Copies a string.

Parameters
[in]str1The string being copied
[out]str2The string being copied into
[in]fieldnmThe name of the field for error messages
[in]callerThe calling routine for error messages

Definition at line 743 of file MOM_io.F90.

743  character(len=*), intent(in) :: str1 !< The string being copied
744  character(len=*), intent(out) :: str2 !< The string being copied into
745  character(len=*), optional, intent(in) :: fieldnm !< The name of the field for error messages
746  character(len=*), optional, intent(in) :: caller !< The calling routine for error messages
747 
748  if (len(trim(str1)) > len(str2)) then
749  if (present(fieldnm) .and. present(caller)) then
750  call mom_error(fatal, trim(caller)//" attempted to copy the overly long"//&
751  " string "//trim(str1)//" into "//trim(fieldnm))
752  else
753  call mom_error(fatal, "safe_string_copy: The string "//trim(str1)//&
754  " is longer than its intended target.")
755  endif
756  endif
757  str2 = trim(str1)

References mom_error_handler::mom_error().

Referenced by modify_vardesc(), query_vardesc(), and var_desc().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ var_desc()

type(vardesc) function, public mom_io::var_desc ( character(len=*), intent(in)  name,
character(len=*), intent(in), optional  units,
character(len=*), intent(in), optional  longname,
character(len=*), intent(in), optional  hor_grid,
character(len=*), intent(in), optional  z_grid,
character(len=*), intent(in), optional  t_grid,
character(len=*), intent(in), optional  cmor_field_name,
character(len=*), intent(in), optional  cmor_units,
character(len=*), intent(in), optional  cmor_longname,
real, intent(in), optional  conversion,
character(len=*), intent(in), optional  caller 
)

Returns a vardesc type whose elements have been filled with the provided fields. The argument name is required, while the others are optional and have default values that are empty strings or are appropriate for a 3-d tracer field at the tracer cell centers.

Parameters
[in]namevariable name
[in]unitsvariable units
[in]longnamevariable long name
[in]hor_gridvariable horizonal staggering
[in]z_gridvariable vertical staggering
[in]t_gridtime description: s, p, or 1
[in]cmor_field_nameCMOR name
[in]cmor_unitsCMOR physical dimensions of variable
[in]cmor_longnameCMOR long name
[in]conversionfor unit conversions, such as needed to convert from intensive to extensive
[in]callercalling routine?
Returns
vardesc type that is created

Definition at line 600 of file MOM_io.F90.

600  character(len=*), intent(in) :: name !< variable name
601  character(len=*), optional, intent(in) :: units !< variable units
602  character(len=*), optional, intent(in) :: longname !< variable long name
603  character(len=*), optional, intent(in) :: hor_grid !< variable horizonal staggering
604  character(len=*), optional, intent(in) :: z_grid !< variable vertical staggering
605  character(len=*), optional, intent(in) :: t_grid !< time description: s, p, or 1
606  character(len=*), optional, intent(in) :: cmor_field_name !< CMOR name
607  character(len=*), optional, intent(in) :: cmor_units !< CMOR physical dimensions of variable
608  character(len=*), optional, intent(in) :: cmor_longname !< CMOR long name
609  real , optional, intent(in) :: conversion !< for unit conversions, such as needed to
610  !! convert from intensive to extensive
611  character(len=*), optional, intent(in) :: caller !< calling routine?
612  type(vardesc) :: vd !< vardesc type that is created
613 
614  character(len=120) :: cllr
615  cllr = "var_desc"
616  if (present(caller)) cllr = trim(caller)
617 
618  call safe_string_copy(name, vd%name, "vd%name", cllr)
619 
620  vd%longname = "" ; vd%units = ""
621  vd%hor_grid = 'h' ; vd%z_grid = 'L' ; vd%t_grid = 's'
622 
623  vd%cmor_field_name = ""
624  vd%cmor_units = ""
625  vd%cmor_longname = ""
626  vd%conversion = 1.0
627 
628  call modify_vardesc(vd, units=units, longname=longname, hor_grid=hor_grid, &
629  z_grid=z_grid, t_grid=t_grid, &
630  cmor_field_name=cmor_field_name,cmor_units=cmor_units, &
631  cmor_longname=cmor_longname, conversion=conversion, caller=cllr)
632 

References modify_vardesc(), and safe_string_copy().

Referenced by mom_meke::meke_alloc_register_restart(), mom_mixed_layer_restrat::mixedlayer_restrat_register_restarts(), advection_test_tracer::register_advection_test_tracer(), mom_barotropic::register_barotropic_restarts(), boundary_impulse_tracer::register_boundary_impulse_tracer(), mom_controlled_forcing::register_ctrl_forcing_restarts(), dome_tracer::register_dome_tracer(), regional_dyes::register_dye_tracer(), dyed_obc_tracer::register_dyed_obc_tracer(), ideal_age_example::register_ideal_age_tracer(), isomip_tracer::register_isomip_tracer(), mom_ocmip2_cfc::register_ocmip2_cfc(), oil_tracer::register_oil_tracer(), pseudo_salt_tracer::register_pseudo_salt_tracer(), mom_dynamics_split_rk2::register_restarts_dyn_split_rk2(), rgc_tracer::register_rgc_tracer(), user_tracer_example::user_register_tracer_example(), and mom_coord_initialization::write_vertgrid_file().

Here is the call graph for this function:
Here is the caller graph for this function: