FMS  2025.03
Flexible Modeling System
fms_diag_file_object.F90
1 !***********************************************************************
2 !* GNU Lesser General Public License
3 !*
4 !* This file is part of the GFDL Flexible Modeling System (FMS).
5 !*
6 !* FMS is free software: you can redistribute it and/or modify it under
7 !* the terms of the GNU Lesser General Public License as published by
8 !* the Free Software Foundation, either version 3 of the License, or (at
9 !* your option) any later version.
10 !*
11 !* FMS is distributed in the hope that it will be useful, but WITHOUT
12 !* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
13 !* FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
14 !* for more details.
15 !*
16 !* You should have received a copy of the GNU Lesser General Public
17 !* License along with FMS. If not, see <http://www.gnu.org/licenses/>.
18 !***********************************************************************
19 !> @defgroup fms_diag_output_yaml_mod fms_diag_output_yaml_mod
20 !> @ingroup diag_manager
21 !! @brief fms_diag_file_object_mod handles the file objects data, functions, and subroutines.
22 !! @author Tom Robinson
23 !! @description The fmsDiagFile_type contains the information for each history file to be written. It has
24 !! a pointer to the information from the diag yaml, additional metadata that comes from the model, and a
25 !! list of the variables and their variable IDs that are in the file.
26 module fms_diag_file_object_mod
27 #ifdef use_yaml
28 use fms2_io_mod, only: fmsnetcdffile_t, fmsnetcdfunstructureddomainfile_t, fmsnetcdfdomainfile_t, &
29  get_instance_filename, open_file, close_file, get_mosaic_tile_file, unlimited, &
30  register_axis, register_field, register_variable_attribute, write_data, &
31  dimension_exists, register_global_attribute, flush_file
32 use diag_data_mod, only: diag_null, no_domain, max_axes, sub_regional, get_base_time, diag_not_registered, &
33  two_d_domain, ug_domain, prepend_date, diag_days, very_large_file_freq, &
37  middle_time, begin_time, end_time, max_str_len, index_gridtype, latlon_gridtype, &
38  null_gridtype, flush_nc_files, diag_init_time
39 use time_manager_mod, only: time_type, operator(>), operator(/=), operator(==), get_date, get_calendar_type, &
40  valid_calendar_types, operator(>=), date_to_string, &
41  OPERATOR(/), OPERATOR(+), operator(<)
42 use fms_diag_time_utils_mod, only: diag_time_inc, get_time_string, get_date_dif
43 use fms_diag_yaml_mod, only: diag_yaml, diagyamlobject_type, diagyamlfiles_type, subregion_type, diagyamlfilesvar_type
44 use fms_diag_axis_object_mod, only: diagdomain_t, get_domain_and_domain_type, fmsdiagaxis_type, &
49 use fms_diag_field_object_mod, only: fmsdiagfield_type
51 use mpp_mod, only: mpp_get_current_pelist, mpp_npes, mpp_root_pe, mpp_pe, mpp_error, fatal, stdout, &
52  uppercase, lowercase, note, mpp_max
53 use platform_mod, only: fms_file_len
54 use mpp_domains_mod, only: mpp_get_ntile_count, mpp_get_ug_domain_ntiles, mpp_get_io_domain_layout, &
55  mpp_get_io_domain_ug_layout
56 implicit none
57 private
58 
60 public :: fmsdiagfile_type, fms_diag_files_object_init, fms_diag_files_object_initialized
61 
62 logical :: fms_diag_files_object_initialized = .false.
63 
64 integer, parameter :: var_string_len = 25
65 
67  private
68  integer :: id !< The number associated with this file in the larger array of files
69  TYPE(time_type) :: model_time !< The last time data was sent for any of the buffers in this file object
70  TYPE(time_type) :: start_time !< The start time for the file
71  TYPE(time_type) :: last_output !< Time of the last time output was writen
72  TYPE(time_type) :: next_output !< Time of the next write
73  TYPE(time_type) :: next_next_output !< Time of the next next write
74  TYPE(time_type) :: no_more_data !< Time to stop receiving data for this file
75  logical :: done_writing_data!< .True. if finished writing data
76 
77  !< This will be used when using the new_file_freq keys in the diag_table.yaml
78  TYPE(time_type) :: next_close !< Time to close the file
79  logical :: is_file_open !< .True. if the file is opened
80 
81  class(fmsnetcdffile_t), allocatable :: fms2io_fileobj !< fms2_io file object for this history file
82  type(diagyamlfiles_type), pointer :: diag_yaml_file => null() !< Pointer to the diag_yaml_file data
83  integer :: type_of_domain !< The type of domain to use to open the file
84  !! NO_DOMAIN, TWO_D_DOMAIN, UG_DOMAIN, SUB_REGIONAL
85  class(diagdomain_t), pointer :: domain !< The domain to use,
86  !! null if NO_DOMAIN or SUB_REGIONAL
87  character(len=:) , dimension(:), allocatable :: file_metadata_from_model !< File metadata that comes from
88  !! the model.
89  integer, dimension(:), allocatable :: field_ids !< Variable IDs corresponding to file_varlist
90  integer, dimension(:), allocatable :: yaml_ids !< IDs corresponding to the yaml field section
91  logical, dimension(:), private, allocatable :: field_registered !< Array corresponding to `field_ids`, .true.
92  !! if the variable has been registered and
93  !! `field_id` has been set for the variable
94  integer, allocatable :: num_registered_fields !< The number of fields registered
95  !! to the file
96  integer, dimension(:), allocatable :: axis_ids !< Array of axis ids in the file
97  integer :: number_of_axis !< Number of axis in the file
98  integer, dimension(:), allocatable :: buffer_ids !< array of buffer ids associated with the file
99  integer :: number_of_buffers !< Number of buffers that have been added to the file
100  logical, allocatable :: time_ops !< .True. if file contains variables that are time_min, time_max, time_average,
101  !! or time_sum
102  integer :: unlim_dimension_level !< The unlimited dimension level currently being written
103  integer :: num_time_levels !< The number of time levels that were actually written to the file
104  logical :: data_has_been_written !< .True. if data has been written for the current unlimited dimension level
105  logical :: is_static !< .True. if the frequency is -1
106  integer :: nz_subaxis !< The number of Z axis currently added to the file
107 
108  contains
109  procedure, public :: add_field_and_yaml_id
110  procedure, public :: add_buffer_id
111  procedure, public :: is_field_registered
112  procedure, public :: init_diurnal_axis
113  procedure, public :: has_file_metadata_from_model
114  procedure, public :: has_fileobj
115  procedure, public :: has_diag_yaml_file
116  procedure, public :: set_domain_from_axis
117  procedure, public :: set_file_domain
118  procedure, public :: add_axes
119  procedure, public :: add_new_axis
120  procedure, public :: update_write_on_this_pe
121  procedure, public :: get_write_on_this_pe
122  procedure, public :: does_axis_exist
123  procedure, public :: define_new_subaxis
124  procedure, public :: add_start_time
125  procedure, public :: set_file_time_ops
126  procedure, public :: has_field_ids
127  procedure, public :: get_time_ops
128  procedure, public :: get_id
129 ! TODO procedure, public :: get_fileobj ! TODO
130 ! TODO procedure, public :: get_diag_yaml_file ! TODO
131  procedure, public :: get_file_metadata_from_model
132  procedure, public :: get_field_ids
133 ! The following fuctions come will use the yaml inquiry functions
134  procedure, public :: get_file_fname
135  procedure, public :: get_file_frequnit
136  procedure, public :: get_file_freq
137  procedure, public :: get_file_timeunit
138  procedure, public :: get_file_unlimdim
139  procedure, public :: get_file_sub_region
140  procedure, public :: get_file_sub_region_grid_type
141  procedure, public :: get_file_new_file_freq
142  procedure, public :: get_filename_time
143  procedure, public :: get_file_new_file_freq_units
144  procedure, public :: get_file_start_time
145  procedure, public :: get_file_duration
146  procedure, public :: get_file_duration_units
147  procedure, public :: get_file_varlist
148  procedure, public :: get_file_global_meta
149  procedure, public :: is_done_writing_data
150  procedure, public :: has_file_fname
151  procedure, public :: has_file_frequnit
152  procedure, public :: has_file_freq
153  procedure, public :: has_file_timeunit
154  procedure, public :: has_file_unlimdim
155  procedure, public :: has_file_sub_region
156  procedure, public :: has_file_new_file_freq
157  procedure, public :: has_file_new_file_freq_units
158  procedure, public :: has_file_start_time
159  procedure, public :: has_file_duration
160  procedure, public :: has_file_duration_units
161  procedure, public :: has_file_varlist
162  procedure, public :: has_file_global_meta
163  procedure, public :: dump_file_obj
164  procedure, public :: get_buffer_ids
165  procedure, public :: get_number_of_buffers
166  procedure, public :: has_send_data_been_called
167  procedure, public :: check_buffer_times
168 end type fmsdiagfile_type
169 
171  integer, dimension(:), allocatable :: sub_axis_ids !< Array of axis ids in the file
172  logical :: write_on_this_pe !< Flag indicating if the subregion is on the current PE
173  logical :: is_subaxis_defined !< Flag indicating if the subaxes have already been defined
174 end type subregionalfile_type
175 
176 !> \brief A container for fmsDiagFile_type. This is used to create the array of files
178  class(fmsdiagfile_type),allocatable :: fms_diag_file !< The individual file object
179 
180  contains
181  procedure :: is_regional
182  procedure :: is_file_static
183  procedure :: open_diag_file
184  procedure :: write_global_metadata
185  procedure :: write_time_metadata
186  procedure :: write_field_data
187  procedure :: write_axis_metadata
188  procedure :: write_field_metadata
189  procedure :: write_axis_data
190  procedure :: writing_on_this_pe
191  procedure :: check_file_times
192  procedure :: is_time_to_close_file
193  procedure :: write_time_data
194  procedure :: update_next_write
195  procedure :: prepare_for_force_write
196  procedure :: init_unlim_dim
197  procedure :: update_current_new_file_freq_index
198  procedure :: get_unlim_dimension_level
199  procedure :: get_num_time_levels
200  procedure :: get_num_tiles
201  procedure :: get_ndistributedfiles
202  procedure :: flush_diag_file
203  procedure :: get_next_output
204  procedure :: get_next_next_output
205  procedure :: close_diag_file
206  procedure :: set_model_time
207  procedure :: get_model_time
208  procedure :: time_to_start_doing_math
210 
211 !type(fmsDiagFile_type), dimension (:), allocatable, target :: FMS_diag_file !< The array of diag files
212 !class(fmsDiagFileContainer_type),dimension (:), allocatable, target :: FMS_diag_file
213 
214 contains
215 
216 !< @brief Allocates the number of files and sets an ID based for each file
217 !! @return true if there are files allocated in the YAML object
218 logical function fms_diag_files_object_init (files_array)
219  class(fmsdiagfilecontainer_type), allocatable, target, intent(inout) :: files_array (:) !< array of diag files
220  class(fmsdiagfile_type), pointer :: obj => null() !< Pointer for each member of the array
221  integer :: nfiles !< Number of files in the diag yaml
222  integer :: i !< Looping iterator
223 
224  if (diag_yaml%has_diag_files()) then
225  nfiles = diag_yaml%size_diag_files()
226  allocate (files_array(nfiles))
227  set_ids_loop: do i= 1,nfiles
228  !> If the file has a sub_regional, define it as one and allocate the sub_axis_ids array.
229  !! This will be set in a add_axes
230  if (diag_yaml%diag_files(i)%has_file_sub_region()) then
231  allocate(subregionalfile_type :: files_array(i)%FMS_diag_file)
232  obj => files_array(i)%FMS_diag_file
233  select type (obj)
234  type is (subregionalfile_type)
235  allocate(obj%sub_axis_ids(max_axes))
236  obj%sub_axis_ids = diag_null
237  obj%write_on_this_pe = .true.
238  obj%is_subaxis_defined = .false.
239  obj%number_of_axis = 0
240  end select
241  else
242  allocate(fmsdiagfile_type::files_array(i)%FMS_diag_file)
243  obj => files_array(i)%FMS_diag_file
244  endif
245  !!
246  obj%diag_yaml_file => diag_yaml%diag_files(i)
247  obj%id = i
248  allocate(obj%field_ids(diag_yaml%diag_files(i)%size_file_varlist()))
249  allocate(obj%buffer_ids(diag_yaml%diag_files(i)%size_file_varlist()))
250  allocate(obj%yaml_ids(diag_yaml%diag_files(i)%size_file_varlist()))
251  allocate(obj%field_registered(diag_yaml%diag_files(i)%size_file_varlist()))
252  !! Initialize the integer arrays
253  obj%field_ids = diag_not_registered
254  obj%yaml_ids = diag_not_registered
255  obj%buffer_ids = diag_not_registered
256  obj%field_registered = .false.
257  obj%num_registered_fields = 0
258  obj%number_of_buffers = 0
259 
260  !> These will be set in a set_file_domain
261  obj%type_of_domain = no_domain
262  obj%domain => null()
263 
264  !> This will be set in a add_axes
265  allocate(obj%axis_ids(max_axes))
266  obj%number_of_axis = 0
267 
268  !> Set the start_time of the file to the base_time and set up the *_output variables
269  obj%done_writing_data = .false.
270 
271  !! Set this to the time passed in to diag_manager_init
272  !! This will be the base_time if nothing was passed in
273  !! This time is appended to the filename if the prepend_date namelist is .True.
274  if (obj%has_file_start_time()) then
275  obj%start_time = obj%get_file_start_time()
276  else
277  obj%start_time = diag_init_time
278  endif
279  obj%last_output = obj%start_time
280  obj%model_time = diag_init_time
281  obj%next_output = diag_time_inc(obj%start_time, obj%get_file_freq(), obj%get_file_frequnit())
282  obj%next_next_output = diag_time_inc(obj%next_output, obj%get_file_freq(), obj%get_file_frequnit())
283 
284  if (obj%has_file_new_file_freq()) then
285  obj%next_close = diag_time_inc(obj%start_time, obj%get_file_new_file_freq(), &
286  obj%get_file_new_file_freq_units())
287  else
288  if (obj%has_file_duration()) then
289  obj%next_close = diag_time_inc(obj%start_time, obj%get_file_duration(), &
290  obj%get_file_duration_units())
291  else
292  obj%next_close = diag_time_inc(obj%start_time, very_large_file_freq, diag_days)
293  endif
294  endif
295  obj%is_file_open = .false.
296 
297  if(obj%has_file_duration()) then
298  obj%no_more_data = diag_time_inc(obj%start_time, obj%get_file_duration(), &
299  obj%get_file_duration_units())
300  else
301  obj%no_more_data = diag_time_inc(obj%start_time, very_large_file_freq, diag_days)
302  endif
303 
304  obj%unlim_dimension_level = 0
305  obj%num_time_levels = 0
306  obj%is_static = obj%get_file_freq() .eq. -1
307  obj%nz_subaxis = 0
308 
309  nullify(obj)
310  enddo set_ids_loop
311  fms_diag_files_object_init = .true.
312  else
313  fms_diag_files_object_init = .false.
314 ! mpp_error("fms_diag_files_object_init: The diag_table.yaml file has not been correctly parsed.",&
315 ! FATAL)
316  endif
317 end function fms_diag_files_object_init
318 
319 !< @brief Determine if the field corresponding to the field_id was registered to the file
320 !! @return .True. if the field was registed to the file
321 pure logical function is_field_registered(this, field_id)
322  class(fmsdiagfile_type), intent(in) :: this !< The file object
323  integer, intent(in) :: field_id !< Id of the field to check
324 
325  is_field_registered = this%field_registered(field_id)
326 end function is_field_registered
327 
328 !> \brief Adds a field and yaml ID to the file
329 subroutine add_field_and_yaml_id (this, new_field_id, yaml_id)
330  class(fmsdiagfile_type), intent(inout) :: this !< The file object
331  integer, intent(in) :: new_field_id !< The field ID to be added to field_ids
332  integer, intent(in) :: yaml_id !< The yaml_id
333 
334  this%num_registered_fields = this%num_registered_fields + 1
335  if (this%num_registered_fields .le. size(this%field_ids)) then
336  this%field_ids( this%num_registered_fields ) = new_field_id
337  this%yaml_ids( this%num_registered_fields ) = yaml_id
338  this%field_registered( this%num_registered_fields ) = .true.
339  else
340  call mpp_error(fatal, "The file: "//this%get_file_fname()//" has already been assigned its maximum "//&
341  "number of fields.")
342  endif
343 end subroutine add_field_and_yaml_id
344 
345 !> \brief Adds a buffer_id to the file object
346 subroutine add_buffer_id (this, buffer_id)
347  class(fmsdiagfile_type), intent(inout) :: this !< The file object
348  integer, intent(in) :: buffer_id !< Buffer id to add to the file
349 
350  this%number_of_buffers = this%number_of_buffers + 1
351  this%buffer_ids(this%number_of_buffers) = buffer_id
352 
353 end subroutine add_buffer_id
354 
355 !> \brief Initializes a diurnal axis for a fileobj
356 !! \note This is going to be called for every variable in the file, if the variable is not a diurnal variable
357 !! it will do nothing. It only defined a diurnal axis once.
358 subroutine init_diurnal_axis(this, diag_axis, naxis, yaml_id)
359  class(fmsdiagfile_type), intent(inout) :: this !< The file object
360  class(fmsdiagaxiscontainer_type), intent(inout) :: diag_axis(:) !< Array of diag_axis object
361  integer, intent(inout) :: naxis !< Number of diag_axis that heve been defined
362  integer, intent(in) :: yaml_id !< The ID to the variable's yaml
363 
364  integer :: i !< For do loops
365  type(diagyamlfilesvar_type), pointer :: field_yaml !< pointer to the yaml entry
366 
367  field_yaml => diag_yaml%get_diag_field_from_id(yaml_id)
368 
369  !< Go away if the file does not need a diurnal axis
370  if (.not. field_yaml%has_n_diurnal()) return
371 
372  !< Check if the diurnal axis is already defined for this number of diurnal samples
373  do i = 1, this%number_of_axis
374  select type(axis=>diag_axis(this%axis_ids(i))%axis)
375  type is (fmsdiagdiurnalaxis_type)
376  if(field_yaml%get_n_diurnal() .eq. axis%get_diurnal_axis_samples()) return
377  end select
378  end do
379 
380  !< If it is not already defined, define it
381  call define_diurnal_axis(diag_axis, naxis, field_yaml%get_n_diurnal(), .true.)
382  call define_diurnal_axis(diag_axis, naxis, field_yaml%get_n_diurnal(), .false.)
383 
384  !< Add it to the list of axis for the file
385  this%number_of_axis = this%number_of_axis + 1
386  this%axis_ids(this%number_of_axis) = naxis !< This is the diurnal axis edges
387 
388  this%number_of_axis = this%number_of_axis + 1
389  this%axis_ids(this%number_of_axis) = naxis - 1 !< This the diurnal axis
390 
391 end subroutine init_diurnal_axis
392 
393 !> \brief Set the time_ops variable in the diag_file object
394 subroutine set_file_time_ops(this, VarYaml, is_static)
395  class(fmsdiagfile_type), intent(inout) :: this !< The file object
396  type (diagyamlfilesvar_type), intent(in) :: varyaml !< The variable's yaml file
397  logical, intent(in) :: is_static !< Flag indicating if variable is static
398  integer, allocatable :: var_reduct !< temp to hold enumerated reduction type
399 
400  !< Go away if the file is static
401  if (this%is_static) return
402  if (is_static) return
403 
404  ! Set time_ops the first time this subroutine it is called
405  if (.not. allocated(this%time_ops)) then
406  var_reduct = varyaml%get_var_reduction()
407 
408  select case (var_reduct)
410  this%time_ops = .true.
411  case (time_none)
412  this%time_ops = .false.
413  end select
414 
415  return
416  endif
417 
418  if (this%time_ops) then
419  if (varyaml%get_var_reduction() .eq. time_none) &
420  call mpp_error(fatal, "The file: "//this%get_file_fname()//&
421  " has variables that are time averaged and instantaneous")
422  else
423  if (varyaml%get_var_reduction() .ne. time_none) &
424  call mpp_error(fatal, "The file: "//this%get_file_fname()//&
425  " has variables that are time averaged and instantaneous")
426  endif
427 end subroutine set_file_time_ops
428 
429 !> \brief Logical function to determine if the variable file_metadata_from_model has been allocated or associated
430 !! \return .True. if file_metadata_from_model exists .False. if file_metadata_from_model has not been set
431 pure logical function has_file_metadata_from_model (this)
432  class(fmsdiagfile_type), intent(in) :: this !< The file object
433  has_file_metadata_from_model = allocated(this%file_metadata_from_model)
434 end function has_file_metadata_from_model
435 
436 !> \brief Logical function to determine if the variable fileobj has been allocated or associated
437 !! \return .True. if fileobj exists .False. if fileobj has not been set
438 pure logical function has_fileobj (this)
439  class(fmsdiagfile_type), intent(in) :: this !< The file object
440  has_fileobj = allocated(this%fms2io_fileobj)
441 end function has_fileobj
442 
443 !> \brief Logical function to determine if the variable diag_yaml_file has been allocated or associated
444 !! \return .True. if diag_yaml_file exists .False. if diag_yaml has not been set
445 pure logical function has_diag_yaml_file (this)
446  class(fmsdiagfile_type), intent(in) :: this !< The file object
447  has_diag_yaml_file = associated(this%diag_yaml_file)
448 end function has_diag_yaml_file
449 
450 !> \brief Get the time to use to determine the filename, if using a wildcard file name (i.e ocn%4yr%2mo%2dy%2hr)
451 !! \return The time to use when determining the filename
452 function get_filename_time(this) &
453  result(res)
454  class(fmsdiagfile_type), intent(in) :: this !< The file object
455  type(time_type) :: res
456  type(time_type) :: file_end_time
457 
458  file_end_time = this%next_close
459  if (this%next_close > this%no_more_data) file_end_time = this%no_more_data
460  select case (this%diag_yaml_file%get_filename_time())
461  case (begin_time)
462  res = this%last_output
463  case (middle_time)
464  res = (this%last_output + file_end_time )/2
465  case (end_time)
466  res = file_end_time
467  end select
468 end function get_filename_time
469 
470 !> \brief Logical function to determine if the variable field_ids has been allocated or associated
471 !! \return .True. if field_ids exists .False. if field_ids has not been set
472 pure logical function has_field_ids (this)
473  class(fmsdiagfile_type), intent(in) :: this !< The file object
474  has_field_ids = allocated(this%field_ids)
475 end function has_field_ids
476 
477 !> \brief Returns a copy of the value of id
478 !! \return A copy of id
479 pure function get_id (this) result (res)
480  class(fmsdiagfile_type), intent(in) :: this !< The file object
481  integer :: res
482  res = this%id
483 end function get_id
484 
485 !> \brief Returns a copy of the value of time_ops
486 !! \return A copy of time_ops
487 pure function get_time_ops (this) result (res)
488  class(fmsdiagfile_type), intent(in) :: this !< The file object
489  logical :: res
490 
491  if (.not. allocated(this%time_ops)) then
492  res = .false.
493  else
494  res = this%time_ops
495  endif
496 end function get_time_ops
497 
498 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
499 !! TODO
500 !> \brief Returns a copy of the value of fileobj
501 !! \return A copy of fileobj
502 !pure function get_fileobj (obj) result (res)
503 ! class(fmsDiagFile_type), intent(in) :: obj !< The file object
504 ! class(FmsNetcdfFile_t) :: res
505 ! res = obj%fileobj
506 !end function get_fileobj
507 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
508 !! TODO
509 !!> \brief Returns a copy of the value of diag_yaml_file
510 !!! \return A copy of diag_yaml_file
511 !pure function get_diag_yaml_file (obj) result (res)
512 ! class(fmsDiagFile_type), intent(in) :: obj !< The file object
513 ! type(diagYamlFiles_type) :: res
514 ! res = obj%diag_yaml_file
515 !end function get_diag_yaml_file
516 
517 !> \brief Returns a copy of the value of file_metadata_from_model
518 !! \return A copy of file_metadata_from_model
519 pure function get_file_metadata_from_model (this) result (res)
520  class(fmsdiagfile_type), intent(in) :: this !< The file object
521  character(len=:), dimension(:), allocatable :: res
522  res = this%file_metadata_from_model
523 end function get_file_metadata_from_model
524 
525 !> \brief Returns a copy of the value of field_ids
526 !! \return A copy of field_ids
527 pure function get_field_ids (this) result (res)
528  class(fmsdiagfile_type), intent(in) :: this !< The file object
529  integer, dimension(:), allocatable :: res
530  allocate(res(size(this%field_ids)))
531  res = this%field_ids
532 end function get_field_ids
533 
534 !!!!!!!!! Functions from diag_yaml_file
535 !> \brief Returns a copy of file_fname from the yaml object
536 !! \return Copy of file_fname
537 pure function get_file_fname (this) result(res)
538  class(fmsdiagfile_type), intent(in) :: this !< The file object
539  character (len=:), allocatable :: res
540  res = this%diag_yaml_file%get_file_fname()
541 end function get_file_fname
542 
543 !> \brief Returns a copy of file_frequnit from the yaml object
544 !! \return Copy of file_frequnit
545 pure function get_file_frequnit (this) result(res)
546  class(fmsdiagfile_type), intent(in) :: this !< The file object
547  integer :: res
548  res = this%diag_yaml_file%get_file_frequnit()
549 end function get_file_frequnit
550 
551 !> \brief Returns a copy of file_freq from the yaml object
552 !! \return Copy of file_freq
553 pure function get_file_freq (this) result(res)
554  class(fmsdiagfile_type), intent(in) :: this !< The file object
555  integer :: res
556  res = this%diag_yaml_file%get_file_freq()
557 end function get_file_freq
558 
559 !> \brief Returns a copy of file_timeunit from the yaml object
560 !! \return Copy of file_timeunit
561 pure function get_file_timeunit (this) result(res)
562  class(fmsdiagfile_type), intent(in) :: this !< The file object
563  integer :: res
564  res = this%diag_yaml_file%get_file_timeunit()
565 end function get_file_timeunit
566 
567 !> \brief Returns a copy of file_unlimdim from the yaml object
568 !! \return Copy of file_unlimdim
569 pure function get_file_unlimdim (this) result(res)
570  class(fmsdiagfile_type), intent(in) :: this !< The file object
571  character (len=:), allocatable :: res
572  res = this%diag_yaml_file%get_file_unlimdim()
573 end function get_file_unlimdim
574 
575 !> \brief Returns a copy of file_sub_region from the yaml object
576 !! \return Pointer to file_sub_region
577 function get_file_sub_region (obj) result(res)
578  class(fmsdiagfile_type), target, intent(in) :: obj !< The file object
579  type(subregion_type), pointer :: res
580  res => obj%diag_yaml_file%get_file_sub_region()
581 end function get_file_sub_region
582 
583 !< @brief Query for the subregion grid type (latlon or index)
584 !! @return Pointer to subregion grid type
585 function get_file_sub_region_grid_type(this) &
586  result(res)
587  class(fmsdiagfile_type), intent(in) :: this !< Diag file object
588  integer :: res
589 
590  type(subregion_type), pointer :: subregion !< Subregion type
591 
592  if(this%diag_yaml_file%has_file_sub_region()) then
593  subregion => this%diag_yaml_file%get_file_sub_region()
594  res = subregion%grid_type
595  else
596  res = null_gridtype
597  endif
598 end function get_file_sub_region_grid_type
599 
600 !> \brief Returns a copy of file_new_file_freq from the yaml object
601 !! \return Copy of file_new_file_freq
602 pure function get_file_new_file_freq (this) result(res)
603  class(fmsdiagfile_type), intent(in) :: this !< The file object
604  integer :: res
605  res = this%diag_yaml_file%get_file_new_file_freq()
606 end function get_file_new_file_freq
607 
608 !> \brief Returns a copy of file_new_file_freq_units from the yaml object
609 !! \return Copy of file_new_file_freq_units
610 pure function get_file_new_file_freq_units (this) result(res)
611  class(fmsdiagfile_type), intent(in) :: this !< The file object
612  integer :: res
613  res = this%diag_yaml_file%get_file_new_file_freq_units()
614 end function get_file_new_file_freq_units
615 
616 !> \brief Returns a copy of file_start_time from the yaml object
617 !! \return Copy of file_start_time
618 pure function get_file_start_time (this) result(res)
619  class(fmsdiagfile_type), intent(in) :: this !< The file object
620  type(time_type):: res
621  res = this%diag_yaml_file%get_file_start_time()
622 end function get_file_start_time
623 
624 !> \brief Returns a copy of file_duration from the yaml object
625 !! \return Copy of file_duration
626 pure function get_file_duration (this) result(res)
627  class(fmsdiagfile_type), intent(in) :: this !< The file object
628  integer :: res
629  res = this%diag_yaml_file%get_file_duration()
630 end function get_file_duration
631 
632 !> \brief Returns a copy of file_duration_units from the yaml object
633 !! \return Copy of file_duration_units
634 pure function get_file_duration_units (this) result(res)
635  class(fmsdiagfile_type), intent(in) :: this !< The file object
636  integer :: res
637  res = this%diag_yaml_file%get_file_duration_units()
638 end function get_file_duration_units
639 
640 !> \brief Returns a copy of file_varlist from the yaml object
641 !! \return Copy of file_varlist
642 pure function get_file_varlist (this) result(res)
643  class(fmsdiagfile_type), intent(in) :: this !< The file object
644  character (len=:), allocatable, dimension(:) :: res
645  res = this%diag_yaml_file%get_file_varlist()
646 end function get_file_varlist
647 
648 !> \brief Returns a copy of file_global_meta from the yaml object
649 !! \return Copy of file_global_meta
650 pure function get_file_global_meta (this) result(res)
651  class(fmsdiagfile_type), intent(in) :: this !< The file object
652  character (len=MAX_STR_LEN), allocatable, dimension(:,:) :: res
653  res = this%diag_yaml_file%get_file_global_meta()
654 end function get_file_global_meta
655 
656 !> \brief Determines if done writing data
657 !! \return .True. if done writing data
658 pure function is_done_writing_data (this) result(res)
659  class(fmsdiagfile_type), intent(in) :: this !< The file object
660  logical :: res
661  res = this%done_writing_data
662  if (this%is_file_open) res = .false.
663 end function is_done_writing_data
664 
665 !> \brief Checks if file_fname is allocated in the yaml object
666 !! \return true if file_fname is allocated
667 pure function has_file_fname (this) result(res)
668  class(fmsdiagfile_type), intent(in) :: this !< The file object
669  logical :: res
670  res = this%diag_yaml_file%has_file_fname()
671 end function has_file_fname
672 
673 !> \brief Checks if file_frequnit is allocated in the yaml object
674 !! \return true if file_frequnit is allocated
675 pure function has_file_frequnit (this) result(res)
676  class(fmsdiagfile_type), intent(in) :: this !< The file object
677  logical :: res
678  res = this%diag_yaml_file%has_file_frequnit()
679 end function has_file_frequnit
680 
681 !> \brief Checks if file_freq is allocated in the yaml object
682 !! \return true if file_freq is allocated
683 pure function has_file_freq (this) result(res)
684  class(fmsdiagfile_type), intent(in) :: this !< The file object
685  logical :: res
686  res = this%diag_yaml_file%has_file_freq()
687 end function has_file_freq
688 
689 !> \brief Checks if file_timeunit is allocated in the yaml object
690 !! \return true if file_timeunit is allocated
691 pure function has_file_timeunit (this) result(res)
692  class(fmsdiagfile_type), intent(in) :: this !< The file object
693  logical :: res
694  res = this%diag_yaml_file%has_file_timeunit()
695 end function has_file_timeunit
696 
697 !> \brief Checks if file_unlimdim is allocated in the yaml object
698 !! \return true if file_unlimdim is allocated
699 pure function has_file_unlimdim (this) result(res)
700  class(fmsdiagfile_type), intent(in) :: this !< The file object
701  logical :: res
702  res = this%diag_yaml_file%has_file_unlimdim()
703 end function has_file_unlimdim
704 
705 !> \brief Checks if file_sub_region is allocated in the yaml object
706 !! \return true if file_sub_region is allocated
707 pure function has_file_sub_region (this) result(res)
708  class(fmsdiagfile_type), intent(in) :: this !< The file object
709  logical :: res
710  res = this%diag_yaml_file%has_file_sub_region()
711 end function has_file_sub_region
712 
713 !> \brief Checks if file_new_file_freq is allocated in the yaml object
714 !! \return true if file_new_file_freq is allocated
715 pure function has_file_new_file_freq (this) result(res)
716  class(fmsdiagfile_type), intent(in) :: this !< The file object
717  logical :: res
718  res = this%diag_yaml_file%has_file_new_file_freq()
719 end function has_file_new_file_freq
720 
721 !> \brief Checks if file_new_file_freq_units is allocated in the yaml object
722 !! \return true if file_new_file_freq_units is allocated
723 pure function has_file_new_file_freq_units (this) result(res)
724  class(fmsdiagfile_type), intent(in) :: this !< The file object
725  logical :: res
726  res = this%diag_yaml_file%has_file_new_file_freq_units()
727 end function has_file_new_file_freq_units
728 
729 !> \brief Checks if file_start_time is allocated in the yaml object
730 !! \return true if file_start_time is allocated
731 pure function has_file_start_time (this) result(res)
732  class(fmsdiagfile_type), intent(in) :: this !< The file object
733  logical :: res
734  res = this%diag_yaml_file%has_file_start_time()
735 end function has_file_start_time
736 
737 !> \brief Checks if file_duration is allocated in the yaml object
738 !! \return true if file_duration is allocated
739 pure function has_file_duration (this) result(res)
740  class(fmsdiagfile_type), intent(in) :: this !< The file object
741  logical :: res
742  res = this%diag_yaml_file%has_file_duration()
743 end function has_file_duration
744 
745 !> \brief Checks if file_duration_units is allocated in the yaml object
746 !! \return true if file_duration_units is allocated
747 pure function has_file_duration_units (this) result(res)
748  class(fmsdiagfile_type), intent(in) :: this !< The file object
749  logical :: res
750  res = this%diag_yaml_file%has_file_duration_units()
751 end function has_file_duration_units
752 
753 !> \brief Checks if file_varlist is allocated in the yaml object
754 !! \return true if file_varlist is allocated
755 pure function has_file_varlist (this) result(res)
756  class(fmsdiagfile_type), intent(in) :: this !< The file object
757  logical :: res
758  res = this%diag_yaml_file%has_file_varlist()
759 end function has_file_varlist
760 
761 !> \brief Checks if file_global_meta is allocated in the yaml object
762 !! \return true if file_global_meta is allocated
763 pure function has_file_global_meta (this) result(res)
764  class(fmsdiagfile_type), intent(in) :: this !< The file object
765  logical :: res
766  res = this%diag_yaml_file%has_file_global_meta()
767 end function has_file_global_meta
768 
769 !> @brief Sets the domain and type of domain from the axis IDs
770 subroutine set_domain_from_axis(this, diag_axis, axes)
771  class(fmsdiagfile_type), intent(inout) :: this !< The file object
772  class(fmsdiagaxiscontainer_type), intent(in) :: diag_axis(:) !< Array of diag_axis
773  integer, intent(in) :: axes (:)
774 
775  call get_domain_and_domain_type(diag_axis, axes, this%type_of_domain, this%domain, this%get_file_fname())
776 end subroutine set_domain_from_axis
777 
778 !> @brief Set the domain and the type_of_domain for a file
779 !> @details This subroutine is going to be called once by every variable in the file
780 !! in register_diag_field. It will update the domain and the type_of_domain if needed and verify that
781 !! all the variables are in the same domain
782 subroutine set_file_domain(this, domain, type_of_domain)
783  class(fmsdiagfile_type), intent(inout) :: this !< The file object
784  integer, INTENT(in) :: type_of_domain !< fileobj_type to use
785  CLASS(diagdomain_t), INTENT(in), pointer :: domain !< Domain
786 
787  if (type_of_domain .ne. this%type_of_domain) then
788  !! If the current type_of_domain in the file obj is not the same as the variable calling this subroutine
789 
790  if (type_of_domain .eq. no_domain .or. this%type_of_domain .eq. no_domain) then
791  !! If they are not the same then one of them can be NO_DOMAIN
792  !! (i.e a file can have variables that are not domain decomposed and variables that are)
793 
794  if (type_of_domain .ne. no_domain) then
795  !! Update the file's type_of_domain and domain if needed
796  this%type_of_domain = type_of_domain
797  this%domain => domain
798  endif
799 
800  else
801  !! If they are not the same and of them is not NO_DOMAIN, then crash because the variables don't have the
802  !! same domain (i.e a file has a variable is that in a 2D domain and one that is in a UG domain)
803 
804  call mpp_error(fatal, "The file: "//this%get_file_fname()//" has variables that are not in the same domain")
805  endif
806  endif
807 
808 end subroutine set_file_domain
809 
810 !> @brief Loops through a variable's axis_ids and adds them to the FMSDiagFile object if they don't exist
811 subroutine add_axes(this, axis_ids, diag_axis, naxis, yaml_id, buffer_id, output_buffers)
812  class(fmsdiagfile_type), intent(inout) :: this !< The file object
813  integer, INTENT(in) :: axis_ids(:) !< Array of axes_ids
814  class(fmsdiagaxiscontainer_type), intent(inout) :: diag_axis(:) !< Diag_axis object
815  integer, intent(inout) :: naxis !< Number of axis that have been
816  !! registered
817  integer, intent(in) :: yaml_id !< Yaml id of the field section for
818  !! this var
819  integer, intent(in) :: buffer_id !< ID of the buffer
820  type(fmsdiagoutputbuffer_type), intent(inout) :: output_buffers(:) !< Array of output buffers
821 
822  type(diagyamlfilesvar_type), pointer :: field_yaml !< pointer to the yaml entry
823 
824  integer :: i, j !< For do loops
825  logical :: is_cube_sphere !< Flag indicating if the file's domain is a cubesphere
826  logical :: axis_found !< Flag indicating that the axis was already to the file obj
827  integer, allocatable :: var_axis_ids(:) !< Array of the variable's axis ids
828  integer :: x_y_axis_id(2) !< Ids of the x and y axis
829  integer :: x_or_y !< integer indicating if the axis is x or y
830  logical :: is_x_or_y !< flag indicating if the axis is x or y
831  integer :: subregion_gridtype !< The type of the subregion (latlon or index)
832  logical :: write_on_this_pe !< Flag indicating if the current pe is in the subregion
833 
834  is_cube_sphere = .false.
835  subregion_gridtype = this%get_file_sub_region_grid_type()
836 
837  field_yaml => diag_yaml%get_diag_field_from_id(yaml_id)
838 
839  !< Created a copy here, because if the variable has a z subaxis var_axis_ids will be modified in
840  !! `create_new_z_subaxis` to contain the id of the new z subaxis instead of the parent axis,
841  !! which will be added to the list of axis in the file object (axis_ids is intent(in),
842  !! which is why the copy was needed)
843  var_axis_ids = axis_ids
844 
845  if (field_yaml%has_var_zbounds()) then
846  call create_new_z_subaxis(field_yaml%get_var_zbounds(), var_axis_ids, diag_axis, naxis, &
847  this%axis_ids, this%number_of_axis, this%nz_subaxis)
848  endif
849 
850  select type(this)
851  type is (subregionalfile_type)
852  if (associated(this%domain)) then
853  if (this%domain%get_ntiles() .eq. 6) is_cube_sphere = .true.
854  endif
855  if (.not. this%get_write_on_this_pe()) return
856  subaxis_defined: if (this%is_subaxis_defined) then
857  do i = 1, size(var_axis_ids)
858  select type (parent_axis => diag_axis(var_axis_ids(i))%axis)
859  type is (fmsdiagfullaxis_type)
860  axis_found = .false.
861  is_x_or_y = parent_axis%is_x_or_y_axis()
862  do j = 1, this%number_of_axis
863  if (is_x_or_y) then
864  if(is_parent_axis(this%axis_ids(j), var_axis_ids(i), diag_axis)) then
865  axis_found = .true.
866  var_axis_ids(i) = this%axis_ids(j) !Set the var_axis_id to the sub axis_id
867  cycle
868  endif
869  elseif (var_axis_ids(i) .eq. this%axis_ids(j)) then
870  axis_found = .true.
871  endif
872  enddo
873 
874  if (.not. axis_found) then
875  if (is_x_or_y) then
876  if (subregion_gridtype .eq. latlon_gridtype .and. is_cube_sphere) &
877  call mpp_error(fatal, "If using the cube sphere and defining the subregion with latlon "//&
878  "the variable need to have the same x and y axis. Please check the variables in the file "//&
879  trim(this%get_file_fname())//" or use indices to define the subregion.")
880 
881  select case (subregion_gridtype)
882  case (index_gridtype)
883  call define_new_subaxis_index(parent_axis, this%get_file_sub_region(), diag_axis, naxis, &
884  i, write_on_this_pe)
885  case (latlon_gridtype)
886  call define_new_subaxis_latlon(diag_axis, var_axis_ids(i:i), naxis, this%get_file_sub_region(), &
887  .false., write_on_this_pe)
888  end select
889  call this%update_write_on_this_pe(write_on_this_pe)
890  if (.not. this%get_write_on_this_pe()) cycle
891  call this%add_new_axis(naxis)
892  var_axis_ids(i) = naxis
893  else
894  call this%add_new_axis(var_axis_ids(i))
895  endif
896  endif
897  type is (fmsdiagsubaxis_type)
898  axis_found = this%does_axis_exist(var_axis_ids(i))
899  if (.not. axis_found) call this%add_new_axis(var_axis_ids(i))
900  end select
901  enddo
902  else
903  x_y_axis_id = diag_null
904  do i = 1, size(var_axis_ids)
905  select type (parent_axis => diag_axis(var_axis_ids(i))%axis)
906  type is (fmsdiagfullaxis_type)
907  if (.not. parent_axis%is_x_or_y_axis(x_or_y)) then
908  axis_found = this%does_axis_exist(var_axis_ids(i))
909  if (.not. axis_found) call this%add_new_axis(var_axis_ids(i))
910  else
911  x_y_axis_id(x_or_y) = var_axis_ids(i)
912  endif
913  type is (fmsdiagsubaxis_type)
914  axis_found = this%does_axis_exist(var_axis_ids(i))
915  if (.not. axis_found) call this%add_new_axis(var_axis_ids(i))
916  end select
917  enddo
918 
919  call this%define_new_subaxis(var_axis_ids, x_y_axis_id, is_cube_sphere, diag_axis, naxis)
920  this%is_subaxis_defined = .true.
921  endif subaxis_defined
922  type is (fmsdiagfile_type)
923  do i = 1, size(var_axis_ids)
924  axis_found = this%does_axis_exist(var_axis_ids(i))
925  if (.not. axis_found) call this%add_new_axis(var_axis_ids(i))
926  enddo
927  end select
928  !> Add the axis to the buffer object
929  call output_buffers(buffer_id)%add_axis_ids(var_axis_ids)
930 end subroutine add_axes
931 
932 !> @brief Adds a new axis the list of axis in the diag file object
933 subroutine add_new_axis(this, var_axis_id)
934  class(fmsdiagfile_type), intent(inout) :: this !< The file object
935  integer, intent(in) :: var_axis_id !< Axis id of the variable
936 
937  this%number_of_axis = this%number_of_axis + 1
938  this%axis_ids(this%number_of_axis) = var_axis_id
939 end subroutine add_new_axis
940 
941 !> @brief This updates write on this pe
942 subroutine update_write_on_this_pe(this, write_on_this_pe)
943  class(fmsdiagfile_type), intent(inout) :: this !< The file object
944  logical, intent(in) :: write_on_this_pe !< .True. if the current PE is in
945  !! subregion
946 
947  select type (this)
948  type is (subregionalfile_type)
949  if (this%write_on_this_pe) this%write_on_this_pe = write_on_this_pe
950  end select
951 end subroutine update_write_on_this_pe
952 
953 !> @brief Query for the write_on_this_pe member of the diag file object
954 !! @return the write_on_this_pe member of the diag file object
955 function get_write_on_this_pe(this) &
956  result(rslt)
957  class(fmsdiagfile_type), intent(inout) :: this !< The file object
958  logical :: rslt
959  rslt = .true.
960  select type (this)
961  type is (subregionalfile_type)
962  rslt= this%write_on_this_pe
963  end select
964 end function get_write_on_this_pe
965 
966 !< @brief Determine if an axis is already in the list of axis for a diag file
967 !! @return .True. if the axis is already in the list of axis for a diag file
968 function does_axis_exist(this, var_axis_id) &
969  result(rslt)
970  class(fmsdiagfile_type), intent(inout) :: this !< The file object
971  integer, intent(in) :: var_axis_id !< Variable axis id to check
972 
973  logical :: rslt
974  integer :: j !< For do loops
975 
976  rslt = .false.
977  do j = 1, this%number_of_axis
978  !> Check if the axis already exists, move on
979  if (var_axis_id .eq. this%axis_ids(j)) then
980  rslt = .true.
981  return
982  endif
983  enddo
984 end function
985 
986 !> @brief Define a new sub axis
987 subroutine define_new_subaxis(this, var_axis_ids, x_y_axis_id, is_cube_sphere, diag_axis, naxis)
988  class(fmsdiagfile_type), intent(inout) :: this !< The file object
989  integer, INTENT(inout) :: var_axis_ids(:) !< Original variable axis ids
990  integer, INTENT(in) :: x_y_axis_id(:) !< The ids of the x and y axis
991  logical, intent(in) :: is_cube_sphere !< .True. if the axis is in the cubesphere
992  integer, intent(inout) :: naxis !< Number of axis current registered
993  class(fmsdiagaxiscontainer_type), intent(inout) :: diag_axis(:) !< Diag_axis object
994 
995  logical :: write_on_this_pe !< .True. if the current PE is in the subregion
996  integer :: i, j !< For do loop
997 
998  select case (this%get_file_sub_region_grid_type())
999  case(latlon_gridtype)
1000  call define_new_subaxis_latlon(diag_axis, x_y_axis_id, naxis, this%get_file_sub_region(), is_cube_sphere, &
1001  write_on_this_pe)
1002  call this%update_write_on_this_pe(write_on_this_pe)
1003  if (.not. this%get_write_on_this_pe()) return
1004  call this%add_new_axis(naxis)
1005  call this%add_new_axis(naxis-1)
1006  do j = 1, size(var_axis_ids)
1007  if (x_y_axis_id(1) .eq. var_axis_ids(j)) var_axis_ids(j) = naxis - 1
1008  if (x_y_axis_id(2) .eq. var_axis_ids(j)) var_axis_ids(j) = naxis
1009  enddo
1010  case (index_gridtype)
1011  do i = 1, size(x_y_axis_id)
1012  select type (parent_axis => diag_axis(x_y_axis_id(i))%axis)
1013  type is (fmsdiagfullaxis_type)
1014  call define_new_subaxis_index(parent_axis, this%get_file_sub_region(), diag_axis, naxis, i, &
1015  write_on_this_pe)
1016  call this%update_write_on_this_pe(write_on_this_pe)
1017  if (.not. this%get_write_on_this_pe()) return
1018  call this%add_new_axis(naxis)
1019  do j = 1, size(var_axis_ids)
1020  if (x_y_axis_id(i) .eq. var_axis_ids(j)) var_axis_ids(j) = naxis
1021  enddo
1022  end select
1023  enddo
1024  end select
1025 end subroutine define_new_subaxis
1026 
1027 !> @brief adds the start time to the fileobj
1028 !! @note This should be called from the register field calls. It can be called multiple times (one for each variable)
1029 !! So it needs to make sure that the start_time is the same for each variable. The initial value is the base_time
1030 subroutine add_start_time(this, start_time)
1031  class(fmsdiagfile_type), intent(inout) :: this !< The file object
1032  TYPE(time_type), intent(in) :: start_time !< Start time passed into register_diag_field
1033 
1034  !< If the start_time sent in is equal to the diag_init_time return because
1035  !! this%start_time was already set to the diag_init_time
1036  if (start_time .eq. diag_init_time) return
1037 
1038  !< If the start_time sent is is greater than or equal to the start time already
1039  !! in the diag file obj return because either this%start_time was already updated
1040  !! or the file has start_time defined in the yaml
1041  if (this%start_time >= start_time) return
1042 
1043  if (this%start_time .ne. diag_init_time) then
1044  !> If the this%start_time is not equal to the diag_init_time from the diag_table
1045  !! this%start_time was already updated so make sure it is the same for the current variable
1046  !! or error out
1047  if (this%start_time .ne. start_time)&
1048  call mpp_error(fatal, "The variables associated with the file:"//this%get_file_fname()//" have &
1049  &different start_time")
1050  else
1051  !> If the this%start_time is equal to the diag_init_time,
1052  !! simply update it with the start_time and set up the *_output variables
1053  this%model_time = start_time
1054  this%start_time = start_time
1055  this%last_output = start_time
1056  this%next_output = diag_time_inc(start_time, this%get_file_freq(), this%get_file_frequnit())
1057  this%next_next_output = diag_time_inc(this%next_output, this%get_file_freq(), this%get_file_frequnit())
1058  if (this%has_file_new_file_freq()) then
1059  this%next_close = diag_time_inc(this%start_time, this%get_file_new_file_freq(), &
1060  this%get_file_new_file_freq_units())
1061  else
1062  if (this%is_static) then
1063  ! If the file is static, set the close time to be equal to the start_time, so that it can be closed
1064  ! after the first write!
1065  this%next_close = this%start_time
1066  this%next_next_output = diag_time_inc(this%start_time, very_large_file_freq, diag_days)
1067  else
1068  this%next_close = diag_time_inc(this%start_time, very_large_file_freq, diag_days)
1069  endif
1070  endif
1071 
1072  if(this%has_file_duration()) then
1073  this%no_more_data = diag_time_inc(this%start_time, this%get_file_duration(), &
1074  this%get_file_duration_units())
1075  else
1076  this%no_more_data = diag_time_inc(this%start_time, very_large_file_freq, diag_days)
1077  endif
1078 
1079  endif
1080 
1081 end subroutine
1082 
1083 !> writes out internal values for fmsDiagFile_type object
1084 subroutine dump_file_obj(this, unit_num)
1085  class(fmsdiagfile_type), intent(in) :: this !< the file object
1086  integer, intent(in) :: unit_num !< passed in from dump_diag_obj
1087  !! will either be for new log file or stdout
1088  write( unit_num, *) 'file id:', this%id
1089  write( unit_num, *) 'start time:', date_to_string(this%start_time)
1090  write( unit_num, *) 'last_output', date_to_string(this%last_output)
1091  write( unit_num, *) 'next_output', date_to_string(this%next_output)
1092  write( unit_num, *)'next_next_output', date_to_string(this%next_next_output)
1093  write( unit_num, *)'next_close', date_to_string(this%next_close)
1094 
1095  if( allocated(this%fms2io_fileobj)) write( unit_num, *)'fileobj path', this%fms2io_fileobj%path
1096 
1097  write( unit_num, *)'type_of_domain', this%type_of_domain
1098  if( allocated(this%file_metadata_from_model)) write( unit_num, *) 'file_metadata_from_model', &
1099  this%file_metadata_from_model
1100  if( allocated(this%field_ids)) write( unit_num, *)'field_ids', this%field_ids
1101  if( allocated(this%field_registered)) write( unit_num, *)'field_registered', this%field_registered
1102  if( allocated(this%num_registered_fields)) write( unit_num, *)'num_registered_fields', this%num_registered_fields
1103  if( allocated(this%axis_ids)) write( unit_num, *)'axis_ids', this%axis_ids(1:this%number_of_axis)
1104 
1105 end subroutine
1106 
1107 !> @brief Determine if a file is regional
1108 !! @return Flag indicating if the file is regional or not
1109 logical pure function is_regional(this)
1110  class(fmsdiagfilecontainer_type), intent(in) :: this !< The file object
1111 
1112  select type (wut=>this%FMS_diag_file)
1113  type is (subregionalfile_type)
1114  is_regional = .true.
1115  type is (fmsdiagfile_type)
1116  is_regional = .false.
1117  end select
1118 
1119 end function is_regional
1120 
1121 !> @brief Determine if a file is static
1122 !! @return Flag indicating if the file is static or not
1123 logical pure function is_file_static(this)
1124 class(fmsdiagfilecontainer_type), intent(in) :: this !< The file object
1125 
1126 is_file_static = .false.
1127 
1128 select type (fileptr=>this%FMS_diag_file)
1129 type is (fmsdiagfile_type)
1130  is_file_static = fileptr%is_static
1131 end select
1132 
1133 end function is_file_static
1134 
1135 !< @brief Opens the diag_file if it is time to do so
1136 subroutine open_diag_file(this, time_step, file_is_opened)
1137  class(fmsdiagfilecontainer_type), intent(inout), target :: this !< The file object
1138  TYPE(time_type), intent(in) :: time_step !< Current model step time
1139  logical, intent(out) :: file_is_opened !< .true. if the file was opened in this
1140  !! time
1141 
1142  class(fmsdiagfile_type), pointer :: diag_file !< Diag_file object to open
1143  class(diagdomain_t), pointer :: domain !< The domain used in the file
1144  character(len=:), allocatable :: diag_file_name !< The file name as defined in the yaml
1145  character(len=FMS_FILE_LEN) :: base_name !< The file name as defined in the yaml
1146  !! without the wildcard definition
1147  character(len=FMS_FILE_LEN) :: file_name !< The file name as it will be written to disk
1148  character(len=FMS_FILE_LEN) :: temp_name !< Temp variable to store the file_name
1149  character(len=128) :: start_date !< The start_time as a string that will be added to
1150  !! the begining of the filename (start_date.filename)
1151  character(len=128) :: suffix !< The current time as a string that will be added to
1152  !! the end of filename
1153  integer :: pos !< Index of the filename with the first "%" in the file name
1154  INTEGER :: year !< The year of the start_date
1155  INTEGER :: month !< The month of the start_date
1156  INTEGER :: day !< The day of the start_date
1157  INTEGER :: hour !< The hour of the start_date
1158  INTEGER :: minute !< The minute of the start_date
1159  INTEGER :: second !< The second of the start_date
1160  character(len=4) :: mype_string !< The pe as a string
1161  logical :: is_regional !< Flag indicating if the file is regional
1162  integer, allocatable :: pes(:) !< Array of the pes in the current pelist
1163 
1164  diag_file => this%FMS_diag_file
1165  domain => diag_file%domain
1166 
1167  file_is_opened = .false.
1168  !< Go away if it the file is already open
1169  if (diag_file%is_file_open) return
1170 
1171  is_regional = .false.
1172  !< Figure out what fms2io_fileobj to use!
1173  if (.not. allocated(diag_file%fms2io_fileobj)) then
1174  select type (diag_file)
1175  type is (subregionalfile_type)
1176  !< In this case each PE is going to write its own file
1177  allocate(fmsnetcdffile_t :: diag_file%fms2io_fileobj)
1178  is_regional = .true.
1179  type is (fmsdiagfile_type)
1180  !< Use the type_of_domain to get the correct fms2io_fileobj
1181  select case (diag_file%type_of_domain)
1182  case (no_domain)
1183  allocate(fmsnetcdffile_t :: diag_file%fms2io_fileobj)
1184  case (two_d_domain)
1185  allocate(fmsnetcdfdomainfile_t :: diag_file%fms2io_fileobj)
1186  case (ug_domain)
1187  allocate(fmsnetcdfunstructureddomainfile_t :: diag_file%fms2io_fileobj)
1188  end select
1189  end select
1190  endif
1191 
1192  !< Figure out what to name of the file
1193  diag_file_name = diag_file%get_file_fname()
1194 
1195  !< If using the new_file_freq figure out what the name is based on the current time
1196  if (diag_file%has_file_new_file_freq()) then
1197  !< If using a wildcard file name (i.e ocn%4yr%2mo%2dy%2hr), get the basename (i.e ocn)
1198  pos = index(diag_file_name, '%')
1199  if (pos > 0) base_name = diag_file_name(1:pos-1)
1200  suffix = get_time_string(diag_file_name, diag_file%get_filename_time())
1201  base_name = trim(base_name)//trim(suffix)
1202  else
1203  base_name = trim(diag_file_name)
1204  endif
1205 
1206  !< Add the ens number to the file name (if it exists)
1207  file_name = trim(base_name)
1208  call get_instance_filename(base_name, file_name)
1209 
1210  !< Prepend the file start_time to the file name if prepend_date == .TRUE. in
1211  !! the namelist
1212  IF ( prepend_date ) THEN
1213  call get_date(diag_file%start_time, year, month, day, hour, minute, second)
1214  write (start_date, '(1I20.4, 2I2.2)') year, month, day
1215 
1216  file_name = trim(adjustl(start_date))//'.'//trim(file_name)
1217  END IF
1218 
1219  file_name = trim(file_name)//".nc"
1220 
1221  !< If this is a regional file add the PE and the tile_number to the filename
1222  if (is_regional) then
1223  !< Get the pe number that will be appended to the end of the file
1224  write(mype_string,'(I0.4)') mpp_pe()
1225 
1226  !< Add the tile number if appropriate
1227  select type (domain)
1228  type is (diagdomain2d_t)
1229  temp_name = file_name
1230  call get_mosaic_tile_file(temp_name, file_name, .true., domain%domain2)
1231  end select
1232 
1233  file_name = trim(file_name)//"."//trim(mype_string)
1234  endif
1235 
1236  !< Open the file!
1237  select type (fms2io_fileobj => diag_file%fms2io_fileobj)
1238  type is (fmsnetcdffile_t)
1239  if (is_regional) then
1240  if (.not. open_file(fms2io_fileobj, file_name, "overwrite", pelist=(/mpp_pe()/))) &
1241  &call mpp_error(fatal, "Error opening the file:"//file_name)
1242  call register_global_attribute(fms2io_fileobj, "is_subregional", "True", str_len=4)
1243  else
1244  allocate(pes(mpp_npes()))
1245  call mpp_get_current_pelist(pes)
1246 
1247  if (.not. open_file(fms2io_fileobj, file_name, "overwrite", pelist=pes)) &
1248  &call mpp_error(fatal, "Error opening the file:"//file_name)
1249  endif
1250  type is (fmsnetcdfdomainfile_t)
1251  select type (domain)
1252  type is (diagdomain2d_t)
1253  if (.not. open_file(fms2io_fileobj, file_name, "overwrite", domain%Domain2)) &
1254  &call mpp_error(fatal, "Error opening the file:"//file_name)
1255  end select
1256  type is (fmsnetcdfunstructureddomainfile_t)
1257  select type (domain)
1258  type is (diagdomainug_t)
1259  if (.not. open_file(fms2io_fileobj, file_name, "overwrite", domain%DomainUG)) &
1260  &call mpp_error(fatal, "Error opening the file:"//file_name)
1261  end select
1262  end select
1263 
1264  file_is_opened = .true.
1265  diag_file%is_file_open = file_is_opened
1266  domain => null()
1267  diag_file => null()
1268 end subroutine open_diag_file
1269 
1270 !< @brief Write global attributes in the diag_file
1271 subroutine write_global_metadata(this)
1272  class(fmsdiagfilecontainer_type), intent(inout), target :: this !< The file object
1273 
1274  class(fmsnetcdffile_t), pointer :: fms2io_fileobj !< The fileobj to write to
1275  integer :: i !< For do loops
1276  character (len=MAX_STR_LEN), allocatable :: yaml_file_attributes(:,:) !< Global attributes defined in the yaml
1277 
1278  type(diagyamlfiles_type), pointer :: diag_file_yaml !< The diag_file yaml
1279 
1280  diag_file_yaml => this%FMS_diag_file%diag_yaml_file
1281  fms2io_fileobj => this%FMS_diag_file%fms2io_fileobj
1282 
1283  if (diag_file_yaml%has_file_global_meta()) then
1284  yaml_file_attributes = diag_file_yaml%get_file_global_meta()
1285  do i = 1, size(yaml_file_attributes,1)
1286  call register_global_attribute(fms2io_fileobj, trim(yaml_file_attributes(i,1)), &
1287  trim(yaml_file_attributes(i,2)), str_len=len_trim(yaml_file_attributes(i,2)))
1288  enddo
1289  deallocate(yaml_file_attributes)
1290  endif
1291 end subroutine write_global_metadata
1292 
1293 !< @brief Writes a variable's metadata in the netcdf file
1294 subroutine write_var_metadata(fms2io_fileobj, variable_name, dimensions, long_name, units)
1295  class(fmsnetcdffile_t), intent(inout) :: fms2io_fileobj !< The file object to write into
1296  character(len=*) , intent(in) :: variable_name !< The name of the time variables
1297  character(len=*) , intent(in) :: dimensions(:) !< The dimensions of the variable
1298  character(len=*) , intent(in) :: long_name !< The long_name of the variable
1299  character(len=*) , intent(in) :: units !< The units of the variable
1300 
1301  call register_field(fms2io_fileobj, variable_name, pack_size_str, dimensions)
1302  call register_variable_attribute(fms2io_fileobj, variable_name, "long_name", &
1303  trim(long_name), str_len=len_trim(long_name))
1304  if (trim(units) .ne. no_units) &
1305  call register_variable_attribute(fms2io_fileobj, variable_name, "units", &
1306  trim(units), str_len=len_trim(units))
1307 end subroutine write_var_metadata
1308 
1309 !> \brief Write the time metadata to the diag file
1310 subroutine write_time_metadata(this)
1311  class(fmsdiagfilecontainer_type), intent(inout), target :: this !< The file object
1312 
1313  class(fmsdiagfile_type), pointer :: diag_file !< Diag_file object to open
1314  class(fmsnetcdffile_t), pointer :: fms2io_fileobj !< The fileobj to write to
1315  character(len=50) :: time_units_str !< Time units written as a string
1316  character(len=50) :: calendar !< The calendar name
1317 
1318  character(len=:), allocatable :: time_var_name !< The name of the time variable as it is defined in the yaml
1319  character(len=50) :: dimensions(2) !< Array of dimensions names for the variable
1320 
1321  diag_file => this%FMS_diag_file
1322  fms2io_fileobj => diag_file%fms2io_fileobj
1323 
1324  time_var_name = diag_file%get_file_unlimdim()
1325  call register_axis(fms2io_fileobj, time_var_name, unlimited)
1326 
1327  WRITE(time_units_str, 11) &
1328  trim(time_unit_list(diag_file%get_file_timeunit())), get_base_year(),&
1329  & get_base_month(), get_base_day(), get_base_hour(), get_base_minute(), get_base_second()
1330 11 FORMAT(a, ' since ', i4.4, '-', i2.2, '-', i2.2, ' ', i2.2, ':', i2.2, ':', i2.2)
1331 
1332  dimensions(1) = "nv"
1333  dimensions(2) = trim(time_var_name)
1334 
1335  call write_var_metadata(fms2io_fileobj, time_var_name, dimensions(2:2), &
1336  time_var_name, time_units_str)
1337 
1338  !< Add additional variables to the time variable
1339  call register_variable_attribute(fms2io_fileobj, time_var_name, "axis", "T", str_len=1 )
1340 
1341  !TODO no need to have both attributes, probably?
1342  calendar = valid_calendar_types(get_calendar_type())
1343  call register_variable_attribute(fms2io_fileobj, time_var_name, "calendar_type", &
1344  uppercase(trim(calendar)), str_len=len_trim(calendar))
1345  call register_variable_attribute(fms2io_fileobj, time_var_name, "calendar", &
1346  lowercase(trim(calendar)), str_len=len_trim(calendar))
1347 
1348  if (diag_file%get_time_ops()) then
1349  call register_variable_attribute(fms2io_fileobj, time_var_name, "bounds", &
1350  trim(time_var_name)//"_bnds", str_len=len_trim(time_var_name//"_bnds"))
1351 
1352  !< It is possible that the "nv" "axis" was registered via "diag_axis_init" call
1353  !! so only adding it if it doesn't exist already
1354  if ( .not. dimension_exists(fms2io_fileobj, "nv")) then
1355  call register_axis(fms2io_fileobj, "nv", 2) !< Time bounds need a vertex number
1356  call write_var_metadata(fms2io_fileobj, "nv", dimensions(1:1), &
1357  "vertex number", no_units)
1358  endif
1359  call write_var_metadata(fms2io_fileobj, time_var_name//"_bnds", dimensions, &
1360  trim(time_var_name)//" axis boundaries", time_units_str)
1361  endif
1362 
1363 end subroutine write_time_metadata
1364 
1365 !> \brief Write out the field data to the file
1366 subroutine write_field_data(this, field_obj, buffer_obj, unlim_dim_was_increased)
1367  class(fmsdiagfilecontainer_type), intent(in), target :: this !< The diag file object to write to
1368  type(fmsdiagfield_type), intent(in), target :: field_obj !< The field object to write from
1369  type(fmsdiagoutputbuffer_type), intent(inout), target :: buffer_obj !< The buffer object with the data
1370  logical, intent(inout) :: unlim_dim_was_increased
1371 
1372  class(fmsdiagfile_type), pointer :: diag_file !< Diag_file object to open
1373  class(fmsnetcdffile_t), pointer :: fms2io_fileobj !< Fileobj to write to
1374  logical :: has_diurnal !< indicates if theres a diurnal axis to adjust for
1375 
1376  diag_file => this%FMS_diag_file
1377  fms2io_fileobj => diag_file%fms2io_fileobj
1378 
1379  !< Increase the unlim dimension index for the output buffer and update the output buffer for the file
1380  !! if haven't already
1381  call buffer_obj%increase_unlim_dim()
1382  if (buffer_obj%get_unlim_dim() > diag_file%unlim_dimension_level) then
1383  diag_file%unlim_dimension_level = buffer_obj%get_unlim_dim()
1384  unlim_dim_was_increased = .true.
1385  endif
1386 
1387  !TODO This may be offloaded in the future
1388  if (diag_file%is_static) then
1389  !< Here the file is static so there is no need for the unlimited dimension
1390  !! as a variables are static
1391  call buffer_obj%write_buffer(fms2io_fileobj)
1392  diag_file%data_has_been_written = .true.
1393  else
1394  if (field_obj%is_static()) then
1395  !< If the variable is static, only write it the first time
1396  if (buffer_obj%get_unlim_dim() .eq. 1) then
1397  call buffer_obj%write_buffer(fms2io_fileobj)
1398  diag_file%data_has_been_written = .true.
1399  endif
1400  else
1401  if (unlim_dim_was_increased) diag_file%data_has_been_written = .true.
1402  has_diurnal = buffer_obj%get_diurnal_sample_size() .gt. 1
1403  call buffer_obj%write_buffer(fms2io_fileobj, &
1404  unlim_dim_level=buffer_obj%get_unlim_dim(), is_diurnal=has_diurnal)
1405  endif
1406  endif
1407 
1408 end subroutine write_field_data
1409 
1410 !> \brief Determine if it is time to close the file
1411 !! \return .True. if it is time to close the file
1412 logical function is_time_to_close_file (this, time_step, force_close)
1413  class(fmsdiagfilecontainer_type), intent(in), target :: this !< The file object
1414  TYPE(time_type), intent(in) :: time_step !< Current model step time
1415  logical, intent(in) :: force_close !< if .true. return true
1416 
1417  if (force_close .or. this%FMS_diag_file%done_writing_data) then
1418  is_time_to_close_file = .true.
1419  elseif (time_step >= this%FMS_diag_file%next_close) then
1420  is_time_to_close_file = .true.
1421  else
1422  if (this%FMS_diag_file%is_static) then
1423  is_time_to_close_file = .true.
1424  else
1425  is_time_to_close_file = .false.
1426  endif
1427  endif
1428 end function
1429 
1430 !> \brief Determine if it is time to start doing mathz
1431 !! \return .True. if it is time to start doing mathz
1432 logical function time_to_start_doing_math (this)
1433  class(fmsdiagfilecontainer_type), intent(in), target :: this !< The file object
1434  time_to_start_doing_math = .false.
1435  if (this%FMS_diag_file%model_time >= this%FMS_diag_file%start_time) then
1436  time_to_start_doing_math = .true.
1437  endif
1438 end function
1439 
1440 !> \brief Determine if it is time to "write" to the file
1441 subroutine check_file_times(this, time_step, output_buffers, diag_fields, do_not_write)
1442  class(fmsdiagfilecontainer_type), intent(inout), target :: this !< The file object
1443  TYPE(time_type), intent(in) :: time_step !< Current model step time
1444  type(fmsdiagoutputbuffer_type), intent(in) :: output_buffers(:) !< Array of output buffer.
1445  !! This is needed for error messages!
1446  type(fmsdiagfield_type), intent(in) :: diag_fields(:) !< Array of diag_fields objects
1447  logical, intent(out) :: do_not_write !< .True. only if this is not a new
1448  !! time step and you are writing
1449  !! at every time step
1450 
1451  do_not_write = .false.
1452  if (time_step > this%FMS_diag_file%next_output) then
1453  if (this%FMS_diag_file%is_static) return
1454  if (time_step > this%FMS_diag_file%next_next_output) then
1455  if (this%FMS_diag_file%get_file_freq() .eq. 0) then
1456  !! If the diag file is being written at every time step
1457  if (time_step .ne. this%FMS_diag_file%next_output) then
1458  !! Only write and update the next_output if it is a new time
1459  call this%FMS_diag_file%check_buffer_times(output_buffers, diag_fields)
1460  this%FMS_diag_file%next_output = time_step
1461  this%FMS_diag_file%next_next_output = time_step
1462  endif
1463  elseif (this%FMS_diag_file%num_registered_fields .eq. 0) then
1464  !! If no variables have been registered, write a dummy time dimension for the first level
1465  !! At least one time level is needed for the combiner to work ...
1466  if (this%FMS_diag_file%unlim_dimension_level .eq. 0) then
1467  call mpp_error(note, this%FMS_diag_file%get_file_fname()//&
1468  ": diag_manager_mod: This file does not have any variables registered. Fill values will be written")
1469  this%FMS_diag_file%data_has_been_written = .true.
1470  this%FMS_diag_file%unlim_dimension_level = 1
1471  endif
1472  else
1473  !! Only fail if send data has actually been called for at least one variable
1474  if (this%FMS_diag_file%has_send_data_been_called(output_buffers, .false.)) &
1475  call mpp_error(fatal, this%FMS_diag_file%get_file_fname()//&
1476  ": diag_manager_mod: You skipped a time_step. Be sure that diag_send_complete is called at every "//&
1477  "time_step needed by the file.")
1478  endif
1479  endif
1480  else
1481  if(this%FMS_diag_file%get_file_freq() .eq. 0) then
1482  do_not_write = .true.
1483  endif
1484  endif
1485 end subroutine check_file_times
1486 
1487 !> \brief Determine if the current PE has data to write
1488 logical function writing_on_this_pe(this)
1489  class(fmsdiagfilecontainer_type), intent(in), target :: this !< The file object
1490 
1491  select type(diag_file => this%FMS_diag_file)
1492  type is (subregionalfile_type)
1493  writing_on_this_pe = diag_file%write_on_this_pe
1494  class default
1495  writing_on_this_pe = .true.
1496  end select
1497 
1498 end function
1499 
1500 !> \brief Write out the time data to the file
1501 subroutine write_time_data(this)
1502  class(fmsdiagfilecontainer_type), intent(in), target :: this !< The file object
1503 
1504  real :: dif !< The time as a real number
1505  class(fmsdiagfile_type), pointer :: diag_file !< Diag_file object to open
1506  class(fmsnetcdffile_t), pointer :: fms2io_fileobj !< The fileobj to write to
1507  TYPE(time_type) :: middle_time !< The middle time of the averaging period
1508 
1509  real :: t1 !< The beginning time of the averaging period
1510  real :: t2 !< The ending time of the averaging period
1511 
1512  diag_file => this%FMS_diag_file
1513  fms2io_fileobj => diag_file%fms2io_fileobj
1514 
1515  !< If data has not been written for the current unlimited dimension leave the subroutine
1516  if (.not. diag_file%data_has_been_written) return
1517 
1518  if (diag_file%get_time_ops()) then
1519  middle_time = (diag_file%last_output+diag_file%next_output)/2
1520  dif = get_date_dif(middle_time, get_base_time(), diag_file%get_file_timeunit())
1521  else
1522  dif = get_date_dif(diag_file%next_output, get_base_time(), diag_file%get_file_timeunit())
1523  endif
1524 
1525  call write_data(fms2io_fileobj, diag_file%get_file_unlimdim(), dif, &
1526  unlim_dim_level=diag_file%unlim_dimension_level)
1527 
1528  if (diag_file%get_time_ops()) then
1529  t1 = get_date_dif(diag_file%last_output, get_base_time(), diag_file%get_file_timeunit())
1530  t2 = get_date_dif(diag_file%next_output, get_base_time(), diag_file%get_file_timeunit())
1531 
1532  call write_data(fms2io_fileobj, trim(diag_file%get_file_unlimdim())//"_bnds", &
1533  (/t1, t2/), unlim_dim_level=diag_file%unlim_dimension_level)
1534 
1535  if (diag_file%unlim_dimension_level .eq. 1) then
1536  call write_data(fms2io_fileobj, "nv", (/1, 2/))
1537  endif
1538  endif
1539 
1540  diag_file%data_has_been_written = .false.
1541 end subroutine write_time_data
1542 
1543 !> \brief Updates the current_new_file_freq_index if using a new_file_freq
1544 subroutine update_current_new_file_freq_index(this, time_step)
1545  class(fmsdiagfilecontainer_type), intent(inout), target :: this !< The file object
1546  TYPE(time_type), intent(in) :: time_step !< Current model step time
1547 
1548  class(fmsdiagfile_type), pointer :: diag_file !< Diag_file object to open
1549 
1550  diag_file => this%FMS_diag_file
1551 
1552  if (time_step >= diag_file%no_more_data) then
1553  call diag_file%diag_yaml_file%increase_new_file_freq_index()
1554 
1555  if (diag_file%has_file_duration()) then
1556  diag_file%no_more_data = diag_time_inc(diag_file%no_more_data, diag_file%get_file_duration(), &
1557  diag_file%get_file_duration_units())
1558  else
1559  !< At this point you are done writing data
1560  diag_file%done_writing_data = .true.
1561  diag_file%no_more_data = diag_time_inc(diag_file%no_more_data, very_large_file_freq, diag_days)
1562  diag_file%next_output = diag_file%no_more_data
1563  diag_file%next_next_output = diag_file%no_more_data
1564  diag_file%last_output = diag_file%no_more_data
1565  endif
1566  endif
1567 
1568  if (diag_file%is_static) diag_file%done_writing_data = .true.
1569 end subroutine update_current_new_file_freq_index
1570 
1571 !> \brief Set up the next_output and next_next_output variable in a file obj
1572 subroutine update_next_write(this, time_step)
1573  class(fmsdiagfilecontainer_type), intent(in), target :: this !< The file object
1574  TYPE(time_type), intent(in) :: time_step !< Current model step time
1575 
1576  class(fmsdiagfile_type), pointer :: diag_file !< Diag_file object to open
1577 
1578  diag_file => this%FMS_diag_file
1579  if (diag_file%is_static) then
1580  diag_file%last_output = diag_file%next_output
1581  diag_file%next_output = diag_time_inc(diag_file%next_output, very_large_file_freq, diag_days)
1582  diag_file%next_next_output = diag_time_inc(diag_file%next_output, very_large_file_freq, diag_days)
1583  else
1584  diag_file%last_output = diag_file%next_output
1585  diag_file%next_output = diag_time_inc(diag_file%next_output, diag_file%get_file_freq(), &
1586  diag_file%get_file_frequnit())
1587  diag_file%next_next_output = diag_time_inc(diag_file%next_output, diag_file%get_file_freq(), &
1588  diag_file%get_file_frequnit())
1589  endif
1590 
1591 end subroutine update_next_write
1592 
1593 !> \brief Prepare the diag file for the force_write
1594 subroutine prepare_for_force_write(this)
1595  class(fmsdiagfilecontainer_type), intent(inout), target :: this !< The file object
1596 
1597  if (this%FMS_diag_file%unlim_dimension_level .eq. 0) then
1598  this%FMS_diag_file%unlim_dimension_level = 1
1599  this%FMS_diag_file%data_has_been_written = .true.
1600  endif
1601 end subroutine prepare_for_force_write
1602 
1603 !> \brief Initialize the unlim dimension in the file and in its buffer objects to 0
1604 subroutine init_unlim_dim(this, output_buffers)
1605  class(fmsdiagfilecontainer_type), intent(inout), target :: this !< The file object
1606  type(fmsdiagoutputbuffer_type), intent(in), target :: output_buffers(:) !< Array of output buffer.
1607 
1608  class(fmsdiagfile_type), pointer :: diag_file !< Diag_file object
1609  type(fmsdiagoutputbuffer_type), pointer :: output_buffer_obj !< Buffer object
1610  integer :: i !< For looping through buffers
1611 
1612  diag_file => this%FMS_diag_file
1613  diag_file%unlim_dimension_level = 0
1614  do i = 1, diag_file%number_of_buffers
1615  output_buffer_obj => output_buffers(diag_file%buffer_ids(i))
1616  call output_buffer_obj%init_buffer_unlim_dim()
1617  enddo
1618 end subroutine init_unlim_dim
1619 
1620 !> \brief Get the number of time levels that were actually written to the file
1621 !! \return Number of time levels that were actually written to the file
1622 function get_num_time_levels(this) &
1623 result(res)
1624  class(fmsdiagfilecontainer_type), intent(in), target :: this !< The file object
1625  integer :: res
1626 
1627  if (this%is_regional()) then
1628  !! If this is a subregional file, num_time_levels will be set to 0 for
1629  !! all PEs that are not in the subregion. If the root pe is not in the subregion
1630  !! then the num_time_levels will not be correct, so this is just getting the max
1631  !! from all PEs
1632  res = this%FMS_diag_file%num_time_levels
1633  call mpp_max(res)
1634  else
1635  res = this%FMS_diag_file%num_time_levels
1636  endif
1637 end function
1638 
1639 !> \brief Get the number of tiles in the file's domain
1640 !! \return Number of tiles in the file's domain
1641 function get_num_tiles(this) &
1642 result(res)
1643  class(fmsdiagfilecontainer_type), intent(in), target :: this !< The file object
1644  integer :: res
1645 
1646  select case(this%FMS_diag_file%type_of_domain)
1647  case (two_d_domain, ug_domain)
1648  select type(domain => this%FMS_diag_file%domain)
1649  type is (diagdomain2d_t)
1650  res = mpp_get_ntile_count(domain%Domain2)
1651  type is (diagdomainug_t)
1652  res = mpp_get_ug_domain_ntiles(domain%DomainUG)
1653  end select
1654  case default
1655  res = 1
1656  end select
1657 end function get_num_tiles
1658 
1659 !> \brief Get number of distributed files that were written
1660 !! \return The number of distributed files that were written
1661 function get_ndistributedfiles(this) &
1662 result(res)
1663  class(fmsdiagfilecontainer_type), intent(in), target :: this !< The file object
1664  integer :: res
1665  integer :: io_layout(2) !< The io_layout
1666 
1667  select case(this%FMS_diag_file%type_of_domain)
1668  case (two_d_domain, ug_domain)
1669  select type(domain => this%FMS_diag_file%domain)
1670  type is (diagdomain2d_t)
1671  io_layout = mpp_get_io_domain_layout(domain%Domain2)
1672  res = io_layout(1) * io_layout(2)
1673  type is (diagdomainug_t)
1674  res = mpp_get_io_domain_ug_layout(domain%DomainUG)
1675  end select
1676  case default
1677  res = 1
1678  end select
1679 end function get_ndistributedfiles
1680 
1681 !> \brief Get the unlimited dimension level that is in the file
1682 !! \return The unlimited dimension
1683 pure function get_unlim_dimension_level(this) &
1684 result(res)
1685  class(fmsdiagfilecontainer_type), intent(in), target :: this !< The file object
1686  integer :: res
1687 
1688  res = this%FMS_diag_file%unlim_dimension_level
1689 end function
1690 
1691 !> \brief Flushes the netcdf file to disk if flush_nc_files is set to .True. in the namelist
1692 subroutine flush_diag_file(this)
1693  class(fmsdiagfilecontainer_type), intent(inout), target :: this !< The file object
1694 
1695  if (flush_nc_files) then
1696  call flush_file(this%FMS_diag_file%fms2io_fileobj)
1697  endif
1698 end subroutine flush_diag_file
1699 
1700 !> \brief Get the next_output for the file object
1701 !! \return The next_output
1702 pure function get_next_output(this) &
1703 result(res)
1704  class(fmsdiagfilecontainer_type), intent(in), target :: this !< The file object
1705  type(time_type) :: res
1706 
1707  res = this%FMS_diag_file%next_output
1708 end function get_next_output
1709 
1710 !> \brief Get the next_output for the file object
1711 !! \return The next_output
1712 pure function get_next_next_output(this) &
1713 result(res)
1714  class(fmsdiagfilecontainer_type), intent(in), target :: this !< The file object
1715  type(time_type) :: res
1716 
1717  res = this%FMS_diag_file%next_next_output
1718  if (this%FMS_diag_file%is_static) then
1719  res = this%FMS_diag_file%no_more_data
1720  endif
1721 end function get_next_next_output
1722 
1723 !< @brief Writes the axis metadata for the file
1724 subroutine write_axis_metadata(this, diag_axis)
1725  class(fmsdiagfilecontainer_type), intent(inout), target :: this !< The file object
1726  class(fmsdiagaxiscontainer_type), intent(in), target :: diag_axis(:) !< Diag_axis object
1727 
1728  class(fmsdiagfile_type), pointer :: diag_file !< Diag_file object to open
1729  class(fmsnetcdffile_t), pointer :: fms2io_fileobj !< The fileobj to write to
1730  integer :: i,k !< For do loops
1731  integer :: parent_axis_id !< Id of the parent_axis
1732  integer :: structured_ids(2) !< Ids of the uncompress axis
1733  integer :: edges_id !< Id of the axis edge
1734 
1735  class(fmsdiagaxiscontainer_type), pointer :: axis_ptr !< pointer to the axis object currently writing
1736  logical :: edges_in_file !< .true. if the edges are already in the file
1737 
1738  diag_file => this%FMS_diag_file
1739  fms2io_fileobj => diag_file%fms2io_fileobj
1740 
1741  do i = 1, diag_file%number_of_axis
1742  edges_in_file = .false.
1743  axis_ptr => diag_axis(diag_file%axis_ids(i))
1744  parent_axis_id = axis_ptr%axis%get_parent_axis_id()
1745 
1746  edges_id = axis_ptr%axis%get_edges_id()
1747  if (edges_id .ne. diag_null) then
1748  !< write the edges if is not in the list of axis in the file, otherwrise ignore
1749  if (any(diag_file%axis_ids(1:diag_file%number_of_axis) .eq. edges_id)) then
1750  edges_in_file = .true.
1751  else
1752  call diag_axis(edges_id)%axis%write_axis_metadata(fms2io_fileobj, .true.)
1753  call diag_file%add_new_axis(edges_id)
1754  endif
1755  endif
1756 
1757  if (parent_axis_id .eq. diag_null) then
1758  call axis_ptr%axis%write_axis_metadata(fms2io_fileobj, edges_in_file)
1759  else
1760  call axis_ptr%axis%write_axis_metadata(fms2io_fileobj, edges_in_file, diag_axis(parent_axis_id)%axis)
1761  endif
1762 
1763  if (axis_ptr%axis%is_unstructured_grid()) then
1764  structured_ids = axis_ptr%axis%get_structured_axis()
1765  do k = 1, size(structured_ids)
1766  call diag_axis(structured_ids(k))%axis%write_axis_metadata(fms2io_fileobj, .false.)
1767  enddo
1768  endif
1769 
1770  enddo
1771 
1772 end subroutine write_axis_metadata
1773 
1774 !< @brief Writes the field metadata for the file
1775 subroutine write_field_metadata(this, diag_field, diag_axis)
1776  class(fmsdiagfilecontainer_type), intent(inout), target :: this !< The file object
1777  class(fmsdiagfield_type) , intent(inout), target :: diag_field(:) !<
1778  class(fmsdiagaxiscontainer_type), intent(in) :: diag_axis(:) !< Diag_axis object
1779 
1780  class(fmsnetcdffile_t), pointer :: fms2io_fileobj !< The fileobj to write to
1781  class(fmsdiagfile_type), pointer :: diag_file !< Diag_file object to open
1782  class(fmsdiagfield_type), pointer :: field_ptr !< diag_field(diag_file%field_ids(i)), for convenience
1783 
1784  integer :: i !< For do loops
1785  logical :: is_regional !< Flag indicating if the field is in a regional file
1786  character(len=255) :: cell_measures !< cell_measures attributes for the field
1787  logical :: need_associated_files !< .True. if the 'associated_files' global attribute is needed
1788  character(len=FMS_FILE_LEN) :: associated_files !< Associated files attribute to add
1789 
1790  is_regional = this%is_regional()
1791 
1792  diag_file => this%FMS_diag_file
1793  fms2io_fileobj => diag_file%fms2io_fileobj
1794 
1795  associated_files = ""
1796  need_associated_files = .false.
1797  do i = 1, size(diag_file%field_ids)
1798  if (.not. diag_file%field_registered(i)) cycle !TODO do something else here
1799  field_ptr => diag_field(diag_file%field_ids(i))
1800 
1801  cell_measures = ""
1802  if (field_ptr%has_area()) then
1803  cell_measures = "area: "//diag_field(field_ptr%get_area())%get_varname(to_write=.true.)
1804 
1805  !! Determine if the area field is already in the file. If it is not create the "associated_files" attribute
1806  !! which contains the file name of the file the area field is in. This is needed for PP/fregrid.
1807  if (.not. diag_field(field_ptr%get_area())%is_variable_in_file(diag_file%id)) then
1808  need_associated_files = .true.
1809  call diag_field(field_ptr%get_area())%generate_associated_files_att(associated_files, diag_file%start_time)
1810  endif
1811  endif
1812 
1813  if (field_ptr%has_volume()) then
1814  cell_measures = trim(cell_measures)//" volume: "//diag_field(field_ptr%get_volume())%get_varname(to_write=.true.)
1815 
1816  !! Determine if the volume field is already in the file. If it is not create the "associated_files" attribute
1817  !! which contains the file name of the file the volume field is in. This is needed for PP/fregrid.
1818  if (.not. diag_field(field_ptr%get_volume())%is_variable_in_file(diag_file%id)) then
1819  need_associated_files = .true.
1820  call diag_field(field_ptr%get_volume())%generate_associated_files_att(associated_files, diag_file%start_time)
1821  endif
1822  endif
1823 
1824  call field_ptr%write_field_metadata(fms2io_fileobj, diag_file%id, diag_file%yaml_ids(i), diag_axis, &
1825  this%FMS_diag_file%get_file_unlimdim(), is_regional, cell_measures)
1826  enddo
1827 
1828  if (need_associated_files) &
1829  call register_global_attribute(fms2io_fileobj, "associated_files", trim(adjustl(associated_files)), &
1830  str_len=len_trim(adjustl(associated_files)))
1831 
1832 end subroutine write_field_metadata
1833 
1834 !< @brief Writes the axis data for the file
1835 subroutine write_axis_data(this, diag_axis)
1836  class(fmsdiagfilecontainer_type), intent(inout), target :: this !< The file object
1837  class(fmsdiagaxiscontainer_type), intent(in) :: diag_axis(:) !< Diag_axis object
1838 
1839  class(fmsdiagfile_type), pointer :: diag_file !< Diag_file object to open
1840  class(fmsnetcdffile_t), pointer :: fms2io_fileobj !< The fileobj to write to
1841  integer :: i, k !< For do loops
1842  integer :: j !< diag_file%axis_ids(i) (for less typing)
1843  integer :: parent_axis_id !< Id of the parent_axis
1844  integer :: structured_ids(2) !< Ids of the uncompress axis
1845 
1846  diag_file => this%FMS_diag_file
1847  fms2io_fileobj => diag_file%fms2io_fileobj
1848 
1849  do i = 1, diag_file%number_of_axis
1850  j = diag_file%axis_ids(i)
1851  parent_axis_id = diag_axis(j)%axis%get_parent_axis_id()
1852  if (parent_axis_id .eq. diag_null) then
1853  call diag_axis(j)%axis%write_axis_data(fms2io_fileobj)
1854  else
1855  call diag_axis(j)%axis%write_axis_data(fms2io_fileobj, diag_axis(parent_axis_id)%axis)
1856  endif
1857 
1858  if (diag_axis(j)%axis%is_unstructured_grid()) then
1859  structured_ids = diag_axis(j)%axis%get_structured_axis()
1860  do k = 1, size(structured_ids)
1861  call diag_axis(structured_ids(k))%axis%write_axis_data(fms2io_fileobj)
1862  enddo
1863  endif
1864  enddo
1865 
1866 end subroutine write_axis_data
1867 
1868 !< @brief Closes the diag_file
1869 subroutine close_diag_file(this, output_buffers, model_end_time, diag_fields)
1870  class(fmsdiagfilecontainer_type), intent(inout), target :: this !< The file object
1871  type(fmsdiagoutputbuffer_type), intent(in) :: output_buffers(:) !< Array of output buffers
1872  !! This is needed for error checking
1873  type(time_type), intent(in) :: model_end_time !< Time that simulation ends
1874  type(fmsdiagfield_type), intent(in), optional :: diag_fields(:) !< Array of diag fields
1875  !! This is needed for error checking
1876 
1877  if (.not. this%FMS_diag_file%is_file_open) return
1878 
1879  !< The select types are needed here because otherwise the code will go to the
1880  !! wrong close_file routine and things will not close propertly
1881  select type( fms2io_fileobj => this%FMS_diag_file%fms2io_fileobj)
1882  type is (fmsnetcdfdomainfile_t)
1883  call close_file(fms2io_fileobj)
1884  type is (fmsnetcdffile_t)
1885  call close_file(fms2io_fileobj)
1886  type is (fmsnetcdfunstructureddomainfile_t)
1887  call close_file(fms2io_fileobj)
1888  end select
1889 
1890  !< Keep track of the number of time levels that were written to the file
1891  !! If the file is using the new_file_freq key, it will be closing the file multiple
1892  !! time, so this ensures that we keep a running count
1893  !! This is going be written to the output yaml, after all the files are closed
1894  this%FMS_diag_file%num_time_levels = this%FMS_diag_file%num_time_levels + &
1895  this%FMS_diag_file%unlim_dimension_level
1896 
1897  !< Reset the unlimited dimension level back to 0, in case the fms2io_fileobj is re-used
1898  this%FMS_diag_file%unlim_dimension_level = 0
1899  this%FMS_diag_file%is_file_open = .false.
1900 
1901  if (this%FMS_diag_file%has_file_new_file_freq()) then
1902  this%FMS_diag_file%next_close = diag_time_inc(this%FMS_diag_file%next_close, &
1903  this%FMS_diag_file%get_file_new_file_freq(), &
1904  this%FMS_diag_file%get_file_new_file_freq_units())
1905  else
1906  this%FMS_diag_file%next_close = model_end_time
1907  endif
1908 
1909  if (this%FMS_diag_file%model_time >= model_end_time) &
1910  this%FMS_diag_file%done_writing_data = .true.
1911  if (this%FMS_diag_file%has_send_data_been_called(output_buffers, .true., diag_fields)) return
1912 end subroutine close_diag_file
1913 
1914 !> \brief Set the model time for the diag file object
1915 subroutine set_model_time(this, model_time)
1916  class(fmsdiagfilecontainer_type), intent(inout) :: this !< The file object
1917  type(time_type), intent(in) :: model_time !< Model time to add
1918 
1919  if (model_time > this%FMS_diag_file%model_time) this%FMS_diag_file%model_time = model_time
1920 end subroutine
1921 
1922 !> \brief Get the model time from the file object
1923 !! \result A pointer to the model time
1924 function get_model_time(this) &
1925  result(rslt)
1926  class(fmsdiagfilecontainer_type), intent(inout), target :: this !< The file object
1927  type(time_type), pointer :: rslt
1928 
1929  rslt => this%FMS_diag_file%model_time
1930 end function get_model_time
1931 
1932 !> \brief Gets the buffer_id list from the file object
1933 pure function get_buffer_ids (this)
1934  class(fmsdiagfile_type), intent(in) :: this !< The file object
1935  integer, allocatable :: get_buffer_ids(:) !< returned buffer ids for this file
1936 
1937  allocate(get_buffer_ids(this%number_of_buffers))
1938  get_buffer_ids = this%buffer_ids(1:this%number_of_buffers)
1939 end function get_buffer_ids
1940 
1941 !> Gets the stored number of buffers from the file object
1942 pure function get_number_of_buffers(this)
1943  class(fmsdiagfile_type), intent(in) :: this !< file object
1944  integer :: get_number_of_buffers !< returned number of buffers
1945  get_number_of_buffers = this%number_of_buffers
1946 end function get_number_of_buffers
1947 
1948 !> Check to ensure that send_data was called at the time step for every output buffer in the file
1949 !! This is only needed when you are output data at every time step
1950 subroutine check_buffer_times(this, output_buffers, diag_fields)
1951  class(fmsdiagfile_type), intent(in) :: this !< file object
1952  type(fmsdiagoutputbuffer_type), intent(in), target :: output_buffers(:) !< Array of output buffers
1953  type(fmsdiagfield_type), intent(in) :: diag_fields(:) !< Array of diag_fields
1954 
1955  integer :: i !< For do loop
1956  type(time_type) :: current_buffer_time !< The buffer time for the current buffer in the do loop
1957  character(len=:), allocatable :: field_name !< The field name (for error messages)
1958  logical :: buffer_time_set !< .True. if current_buffer_time has been set
1959  type(fmsdiagoutputbuffer_type), pointer :: output_buffer_obj !< Pointer to the output buffer
1960 
1961  buffer_time_set = .false.
1962  do i = 1, this%number_of_buffers
1963  output_buffer_obj => output_buffers(this%buffer_ids(i))
1964  if (diag_fields(output_buffer_obj%get_field_id())%is_static()) cycle
1965  if (.not. buffer_time_set) then
1966  current_buffer_time = output_buffer_obj%get_buffer_time()
1967  field_name = output_buffer_obj%get_buffer_name()
1968  buffer_time_set = .true.
1969  else
1970  if (current_buffer_time .ne. output_buffer_obj%get_buffer_time()) &
1971  call mpp_error(fatal, "Send data has not been called at the same time steps for the fields:"//&
1972  field_name//" and "//output_buffer_obj%get_buffer_name()//&
1973  " in file:"//this%get_file_fname())
1974  endif
1975  enddo
1976 end subroutine
1977 
1978 !> @brief Determine if send_data has been called for any fields in the file. Prints out warnings, if indicated
1979 !! @return .True. if send_data has been called for any fields in the file
1980 function has_send_data_been_called(this, output_buffers, print_warnings, diag_fields) &
1981 result(rslt)
1982  class(fmsdiagfile_type), intent(in) :: this !< file object
1983  type(fmsdiagoutputbuffer_type), intent(in), target :: output_buffers(:) !< Array of output buffers
1984  logical, intent(in) :: print_warnings !< .True. if printing warnings
1985  type(fmsdiagfield_type), intent(in), optional :: diag_fields(:) !< Array of diag fields
1986 
1987  logical :: rslt
1988  integer :: i !< For do loops
1989  integer :: field_id !< Field id
1990 
1991  rslt = .false.
1992 
1993  if (print_warnings) then
1994  do i = 1, this%number_of_buffers
1995  if (.not. output_buffers(this%buffer_ids(i))%is_there_data_to_write()) then
1996  field_id = output_buffers(this%buffer_ids(i))%get_field_id()
1997  call mpp_error(note, "Send data was never called for field:"//&
1998  trim(diag_fields(field_id)%get_varname())//" mod: "//trim(diag_fields(field_id)%get_modname())//&
1999  " in file: "//trim(this%get_file_fname())//". Writing FILL VALUES!")
2000  endif
2001  enddo
2002  else
2003  do i = 1, this%number_of_buffers
2004  if (output_buffers(this%buffer_ids(i))%is_there_data_to_write()) then
2005  rslt = .true.
2006  return
2007  endif
2008  enddo
2009  endif
2010 end function has_send_data_been_called
2011 #endif
2012 end module fms_diag_file_object_mod
character(len=8) no_units
String indicating that the variable has no units.
Definition: diag_data.F90:125
integer, parameter sub_regional
This is a file with a sub_region use the FmsNetcdfFile_t fileobj.
Definition: diag_data.F90:104
integer, parameter max_str_len
Max length for a string.
Definition: diag_data.F90:129
integer function get_base_minute()
gets the module variable base_minute
Definition: diag_data.F90:551
integer function get_base_year()
gets the module variable base_year
Definition: diag_data.F90:519
integer function get_base_hour()
gets the module variable base_hour
Definition: diag_data.F90:543
logical flush_nc_files
Control if diag_manager will force a flush of the netCDF file on each write. Note: changing this to ....
Definition: diag_data.F90:365
integer, parameter no_domain
Use the FmsNetcdfFile_t fileobj.
Definition: diag_data.F90:101
character(len=7) avg_name
Name of the average fields.
Definition: diag_data.F90:124
integer, parameter end_time
Use the end of the time average bounds.
Definition: diag_data.F90:128
character(len=6) pack_size_str
Pack size as a string to be used in fms2_io register call set to "double" or "float".
Definition: diag_data.F90:409
type(time_type) function get_base_time()
gets the module variable base_time
Definition: diag_data.F90:511
integer max_axes
Maximum number of independent axes.
Definition: diag_data.F90:361
integer function get_base_day()
gets the module variable base_day
Definition: diag_data.F90:535
type(time_type) diag_init_time
Time diag_manager_init called. If init_time not included in diag_manager_init call,...
Definition: diag_data.F90:417
integer, parameter time_min
The reduction method is min value.
Definition: diag_data.F90:117
integer, parameter ug_domain
Use the FmsNetcdfUnstructuredDomainFile_t fileobj.
Definition: diag_data.F90:103
integer, parameter time_diurnal
The reduction method is diurnal.
Definition: diag_data.F90:122
integer, parameter time_power
The reduction method is average with exponents.
Definition: diag_data.F90:123
logical prepend_date
Should the history file have the start date prepended to the file name. .TRUE. is only supported if t...
Definition: diag_data.F90:386
integer, parameter begin_time
Use the begining of the time average bounds.
Definition: diag_data.F90:126
integer, parameter time_average
The reduction method is average of values.
Definition: diag_data.F90:120
integer function get_base_month()
gets the module variable base_month
Definition: diag_data.F90:527
integer, parameter time_sum
The reduction method is sum of values.
Definition: diag_data.F90:119
integer, parameter time_rms
The reudction method is root mean square of values.
Definition: diag_data.F90:121
integer, parameter middle_time
Use the middle of the time average bounds.
Definition: diag_data.F90:127
integer, parameter time_none
There is no reduction method.
Definition: diag_data.F90:116
integer function get_base_second()
gets the module variable base_second
Definition: diag_data.F90:559
integer, parameter time_max
The reduction method is max value.
Definition: diag_data.F90:118
integer, parameter two_d_domain
Use the FmsNetcdfDomainFile_t fileobj.
Definition: diag_data.F90:102
Close a netcdf or domain file opened with open_file or open_virtual_file.
Definition: fms2_io.F90:166
Opens a given netcdf or domain file.
Definition: fms2_io.F90:122
Add a dimension to a given file.
Definition: fms2_io.F90:190
Defines a new field within the given file Example usage:
Definition: fms2_io.F90:212
Write data to a defined field within a file Example usage:
Definition: fms2_io.F90:262
character(len=128) function, public get_time_string(filename, current_time)
This function determines a string based on current time. This string is used as suffix in output file...
real function, public get_date_dif(t2, t1, units)
Return the difference between two times in units.
subroutine, public define_new_subaxis_latlon(diag_axis, axis_ids, naxis, subRegion, is_cube_sphere, write_on_this_pe)
Fill in the subaxis object for a subRegion defined by lat lon.
subroutine, public define_new_subaxis_index(parent_axis, subRegion, diag_axis, naxis, is_x_or_y, write_on_this_pe)
Fill in the subaxis object for a subRegion defined by index.
subroutine, public define_diurnal_axis(diag_axis, naxis, n_diurnal_samples, is_edges)
Defined a new diurnal axis.
subroutine, public create_new_z_subaxis(zbounds, var_axis_ids, diag_axis, naxis, file_axis_id, nfile_axis, nz_subaxis)
Creates a new z subaxis to use.
logical function, public is_parent_axis(axis_id, parent_axis_id, diag_axis)
Determine if the diag_axis(parent_axis_id) is the parent of diag_axis(axis_id)
subroutine, public get_domain_and_domain_type(diag_axis, axis_id, domain_type, domain, var_name)
Loop through a variable's axis_id to determine and return the domain type and domain to use.
Object that holds the information of the diag_yaml.
integer function mpp_get_ntile_count(domain)
Returns number of tiles in mosaic.
integer function, dimension(2) mpp_get_io_domain_layout(domain)
Set user stack size.
integer function stdout()
This function returns the current standard fortran unit numbers for output.
Definition: mpp_util.inc:43
integer function mpp_npes()
Returns processor count for current pelist.
Definition: mpp_util.inc:421
integer function mpp_pe()
Returns processor ID.
Definition: mpp_util.inc:407
Error handler.
Definition: mpp.F90:382
Reduction operations. Find the max of scalar a from the PEs in pelist result is also automatically br...
Definition: mpp.F90:538
character(len=15) function, public date_to_string(time, err_msg)
Get the a character string that represents the time. The format will be yyyymmdd.hhmmss.
character(len=24) function, public valid_calendar_types(ncal, err_msg)
Returns a character string that describes the calendar type corresponding to the input integer.
subroutine, public get_date(time, year, month, day, hour, minute, second, tick, err_msg)
Gets the date for different calendar types. Given a time_interval, returns the corresponding date und...
integer function, public get_calendar_type()
Returns default calendar type for mapping from time to date.
Type to represent amounts of time. Implemented as seconds and days to allow for larger intervals.
Contains buffer types and routines for the diag manager.
Type to hold the domain info for an axis This type was created to avoid having to send in "Domain",...
Type to hold the unstructured domain.
Type to hold the diagnostic axis description.
Type to hold the diag_axis (either subaxis or a full axis)
Type to hold the diagnostic axis description.
Object that holds all variable information.
A container for fmsDiagFile_type. This is used to create the array of files.
type to hold the diag_file information
type to hold the info a diag_field
type to hold the sub region information about a file