diag_manager/fms_diag_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 module fms_diag_object_mod
  20 use mpp_mod, only: fatal, note, warning, mpp_error, mpp_pe, mpp_root_pe, stdout
  21 use diag_data_mod,  only: diag_null, diag_not_found, diag_not_registered, diag_registered_id, &
  22                          &DIAG_FIELD_NOT_FOUND, diag_not_registered, max_axes, TWO_D_DOMAIN, &
  23                          &get_base_time, NULL_AXIS_ID, get_var_type, diag_not_registered, &
  24                          &time_none, time_max, time_min, time_sum, time_average, time_diurnal, &
  25                          &time_power, time_rms, r8
  26
  27   USE time_manager_mod, ONLY: set_time, set_date, get_time, time_type, OPERATOR(>=), OPERATOR(>),&
  28        & OPERATOR(<), OPERATOR(==), OPERATOR(/=), OPERATOR(/), OPERATOR(+), ASSIGNMENT(=), get_date, &
  29        & get_ticks_per_second
  30 #ifdef use_yaml
  31 use fms_diag_file_object_mod, only: fmsDiagFileContainer_type, fmsDiagFile_type, fms_diag_files_object_init
  32 use fms_diag_field_object_mod, only: fmsDiagField_type, fms_diag_fields_object_init, get_default_missing_value, &
  33                                      check_for_slices
  34 use fms_diag_yaml_mod, only: diag_yaml_object_init, diag_yaml_object_end, find_diag_field, &
  35                            & get_diag_files_id, diag_yaml, get_diag_field_ids, DiagYamlFilesVar_type
  36 use fms_diag_axis_object_mod, only: fms_diag_axis_object_init, fmsDiagAxis_type, fmsDiagSubAxis_type, &
  37                                    &diagDomain_t, get_domain_and_domain_type, diagDomain2d_t, &
  38                                    &fmsDiagAxisContainer_type, fms_diag_axis_object_end, fmsDiagFullAxis_type, &
  39                                    &parse_compress_att, get_axis_id_from_name
  40 use fms_diag_output_buffer_mod
  41 use fms_mod, only: fms_error_handler
  42 use fms_diag_reduction_methods_mod, only: check_indices_order, init_mask, set_weight
  43 use constants_mod, only: SECONDS_PER_DAY
  44 #endif
  45 USE fms_diag_bbox_mod, ONLY: fmsDiagIbounds_type, determine_if_block_is_in_region
  46 #if defined(_OPENMP)
  47 use omp_lib
  48 #endif
  49 use mpp_domains_mod, only: domain1d, domain2d, domainUG, null_domain2d
  50 use fms_string_utils_mod, only: string
  51 use platform_mod
  52 implicit none
  53 private
  54
  55 type fmsDiagObject_type
  56 !TODO add container arrays
  57 #ifdef use_yaml
  58 private
  59 !TODO: Remove FMS prefix from variables in this type
  60   class(fmsDiagFileContainer_type), allocatable :: FMS_diag_files (:) !< array of diag files
  61   class(fmsDiagField_type), allocatable :: FMS_diag_fields(:) !< Array of diag fields
  62   type(fmsDiagOutputBuffer_type), allocatable :: FMS_diag_output_buffers(:) !< array of output buffer objects
  63                                                                        !! one for each variable in the diag_table.yaml
  64   integer, private :: registered_buffers = 0 !< number of registered buffers, per dimension
  65   class(fmsDiagAxisContainer_type), allocatable :: diag_axis(:) !< Array of diag_axis
  66   type(time_type)  :: current_model_time !< The current model time
  67   integer, private :: registered_variables !< Number of registered variables
  68   integer, private :: registered_axis !< Number of registered axis
  69   logical, private :: initialized=.false. !< True if the fmsDiagObject is initialized
  70   logical, private :: files_initialized=.false. !< True if the fmsDiagObject is initialized
  71   logical, private :: fields_initialized=.false. !< True if the fmsDiagObject is initialized
  72   logical, private :: buffers_initialized=.false. !< True if the fmsDiagObject is initialized
  73   logical, private :: axes_initialized=.false. !< True if the fmsDiagObject is initialized
  74 #endif
  75   contains
  76     procedure :: init => fms_diag_object_init
  77     procedure :: diag_end => fms_diag_object_end
  78     procedure :: fms_register_diag_field_scalar
  79     procedure :: fms_register_diag_field_array
  80     procedure :: fms_register_static_field
  81     procedure :: fms_diag_axis_init
  82     procedure :: register => fms_register_diag_field_obj !! Merely initialize fields.
  83     procedure :: fms_diag_field_add_attribute
  84     procedure :: fms_diag_axis_add_attribute
  85     procedure :: fms_get_domain2d
  86     procedure :: fms_get_axis_length
  87     procedure :: fms_get_diag_field_id_from_name
  88     procedure :: fms_get_field_name_from_id
  89     procedure :: fms_get_axis_name_from_id
  90     procedure :: fms_diag_accept_data
  91     procedure :: fms_diag_send_complete
  92     procedure :: do_buffer_math
  93     procedure :: fms_diag_do_io
  94     procedure :: fms_diag_do_reduction
  95     procedure :: fms_diag_field_add_cell_measures
  96     procedure :: allocate_diag_field_output_buffers
  97     procedure :: fms_diag_compare_window
  98     procedure :: update_current_model_time
  99 #ifdef use_yaml
 100     procedure :: get_diag_buffer
 101 #endif
 102 end type fmsDiagObject_type
 103
 104 type (fmsDiagObject_type), target :: fms_diag_object
 105
 106 public :: fms_register_diag_field_obj
 107 public :: fms_register_diag_field_scalar
 108 public :: fms_register_diag_field_array
 109 public :: fms_register_static_field
 110 public :: fms_diag_field_add_attribute
 111 public :: fms_get_diag_field_id_from_name
 112 public :: fms_diag_object
 113 public :: fmsDiagObject_type
 114 integer, private :: registered_variables !< Number of registered variables
 115 public :: dump_diag_obj
 116
 117 contains
 118
 119 !> @brief Initiliazes the  fms_diag_object.
 120 !! Reads the diag_table.yaml and fills in the yaml object
 121 !! Allocates the diag manager object arrays for files, fields, and buffers
 122 !! Initializes variables
 123 subroutine fms_diag_object_init (this,diag_subset_output)
 124  class(fmsDiagObject_type) :: this !< Diag mediator/controller object
 125  integer :: diag_subset_output !< Subset of the diag output?
 126 #ifdef use_yaml
 127  if (this%initialized) return
 128
 129 ! allocate(diag_objs(get_num_unique_fields()))
 130   CALL diag_yaml_object_init(diag_subset_output)
 131   this%axes_initialized = fms_diag_axis_object_init(this%diag_axis)
 132   this%files_initialized = fms_diag_files_object_init(this%FMS_diag_files)
 133   this%fields_initialized = fms_diag_fields_object_init(this%FMS_diag_fields)
 134   this%buffers_initialized =fms_diag_output_buffer_init(this%FMS_diag_output_buffers,SIZE(diag_yaml%get_diag_fields()))
 135   this%registered_variables = 0
 136   this%registered_axis = 0
 137   this%current_model_time = get_base_time()
 138   this%initialized = .true.
 139 #else
 140   call mpp_error("fms_diag_object_init",&
 141     "You must compile with -Duse_yaml to use the option use_modern_diag", FATAL)
 142 #endif
 143 end subroutine fms_diag_object_init
 144
 145 !> \description Loops through all files and does one final write.
 146 !! Closes all files
 147 !! Deallocates all buffers, fields, and files
 148 !! Uninitializes the fms_diag_object
 149 subroutine fms_diag_object_end (this, time)
 150   class(fmsDiagObject_type) :: this
 151   TYPE(time_type), INTENT(in) :: time
 152
 153   integer                   :: i
 154 #ifdef use_yaml
 155   !TODO: loop through files and force write
 156   if (.not. this%initialized) return
 157
 158   call this%do_buffer_math()
 159   call this%fms_diag_do_io(end_time=time)
 160   !TODO: Deallocate diag object arrays and clean up all memory
 161   do i=1, size(this%FMS_diag_output_buffers)
 162     call this%FMS_diag_output_buffers(i)%flush_buffer()
 163   enddo
 164   deallocate(this%FMS_diag_output_buffers)
 165   this%axes_initialized = fms_diag_axis_object_end(this%diag_axis)
 166   this%initialized = .false.
 167   call diag_yaml_object_end
 168 #else
 169   call mpp_error(FATAL, "You can not call fms_diag_object%end without yaml")
 170 #endif
 171 end subroutine fms_diag_object_end
 172
 173 !> @brief Registers a field.
 174 !! @description This to avoid having duplicate code in each of the _scalar, _array and _static register calls
 175 !! @return field index to be used in subsequent calls to send_data or DIAG_FIELD_NOT_FOUND if the field is not
 176 !! in the diag_table.yaml
 177 integer function fms_register_diag_field_obj &
 178        (this, modname, varname, axes, init_time, &
 179        longname, units, missing_value, varRange, mask_variant, standname, &
 180        do_not_log, err_msg, interp_method, tile_count, area, volume, realm, static, &
 181        multiple_send_data)
 182
 183  class(fmsDiagObject_type),TARGET,INTENT(inout):: this       !< Diaj_obj to fill
 184  CHARACTER(len=*),               INTENT(in)    :: modname               !< The module name
 185  CHARACTER(len=*),               INTENT(in)    :: varname               !< The variable name
 186  TYPE(time_type),  OPTIONAL,     INTENT(in)    :: init_time             !< Initial time
 187  INTEGER, TARGET,  OPTIONAL,     INTENT(in)    :: axes(:)               !< The axes indicies
 188  CHARACTER(len=*), OPTIONAL,     INTENT(in)    :: longname              !< THe variables long name
 189  CHARACTER(len=*), OPTIONAL,     INTENT(in)    :: units                 !< The units of the variables
 190  CHARACTER(len=*), OPTIONAL,     INTENT(in)    :: standname             !< The variables stanard name
 191  class(*),         OPTIONAL,     INTENT(in)    :: missing_value         !< Missing value to add as a attribute
 192  class(*),         OPTIONAL,     INTENT(in)    :: varRANGE(2)           !< Range to add as a attribute
 193  LOGICAL,          OPTIONAL,     INTENT(in)    :: mask_variant          !< .True. if mask changes over time
 194  LOGICAL,          OPTIONAL,     INTENT(in)    :: do_not_log            !< if TRUE, field info is not logged
 195  CHARACTER(len=*), OPTIONAL,     INTENT(out)   :: err_msg               !< Error message to be passed back up
 196  CHARACTER(len=*), OPTIONAL,     INTENT(in)    :: interp_method         !< The interp method to be used when
 197                                                                         !! regridding the field in post-processing.
 198                                                                         !! Valid options are "conserve_order1",
 199                                                                         !! "conserve_order2", and "none".
 200  INTEGER,          OPTIONAL,     INTENT(in)    :: tile_count            !< the number of tiles
 201  INTEGER,          OPTIONAL,     INTENT(in)    :: area                  !< diag_field_id of the cell area field
 202  INTEGER,          OPTIONAL,     INTENT(in)    :: volume                !< diag_field_id of the cell volume field
 203  CHARACTER(len=*), OPTIONAL,     INTENT(in)    :: realm                 !< String to set as the value to the
 204                                                                         !! modeling_realm attribute
 205  LOGICAL,          OPTIONAL,     INTENT(in)    :: static                !< True if the variable is static
 206  LOGICAL,          OPTIONAL,     INTENT(in)    :: multiple_send_data    !< .True. if send data is called, multiple
 207                                                                         !! times for the same time
 208
 209 #ifdef use_yaml
 210
 211  class (fmsDiagFile_type), pointer :: fileptr !< Pointer to the diag_file
 212  class (fmsDiagField_type), pointer :: fieldptr !< Pointer to the diag_field
 213  class (fmsDiagOutputBuffer_type), pointer :: bufferptr !< Pointer to the output buffer
 214  class (diagYamlFilesVar_type), pointer :: yamlfptr !< Pointer to yaml object to get the reduction method
 215  integer, allocatable :: file_ids(:) !< The file IDs for this variable
 216  integer :: i !< For do loops
 217  integer, allocatable :: diag_field_indices(:) !< indices where the field was found in the yaml
 218 #endif
 219 #ifndef use_yaml
 220 fms_register_diag_field_obj = DIAG_FIELD_NOT_FOUND
 221 CALL MPP_ERROR(FATAL,"You can not use the modern diag manager without compiling with -Duse_yaml")
 222 #else
 223  diag_field_indices = find_diag_field(varname, modname)
 224  if (diag_field_indices(1) .eq. diag_null) then
 225     !< The field was not found in the table, so return DIAG_FIELD_NOT_FOUND
 226     fms_register_diag_field_obj = DIAG_FIELD_NOT_FOUND
 227     deallocate(diag_field_indices)
 228     return
 229   endif
 230
 231   this%registered_variables = this%registered_variables + 1
 232   fms_register_diag_field_obj = this%registered_variables
 233
 234   call this%FMS_diag_fields(this%registered_variables)%&
 235     &setID(this%registered_variables)
 236
 237 !> Use pointers for convenience
 238   fieldptr => this%FMS_diag_fields(this%registered_variables)
 239 !> Get the file IDs from the field indicies from the yaml
 240   file_ids = get_diag_files_id(diag_field_indices)
 241   call fieldptr%set_file_ids(file_ids)
 242
 243 !> Allocate and initialize member buffer_allocated of this field
 244   fieldptr%buffer_allocated = .false.
 245   fieldptr%buffer_ids = get_diag_field_ids(diag_field_indices)
 246
 247 !> Register the data for the field
 248   call fieldptr%register(modname, varname, diag_field_indices, this%diag_axis, &
 249        axes=axes, longname=longname, units=units, missing_value=missing_value, varRange= varRange, &
 250        mask_variant= mask_variant, standname=standname, do_not_log=do_not_log, err_msg=err_msg, &
 251        interp_method=interp_method, tile_count=tile_count, area=area, volume=volume, realm=realm, &
 252        static=static, multiple_send_data=multiple_send_data)
 253
 254 !> Add the axis information, initial time, and field IDs to the files
 255   if (present(axes) .and. present(init_time)) then
 256     do i = 1, size(file_ids)
 257      fileptr => this%FMS_diag_files(file_ids(i))%FMS_diag_file
 258      call fileptr%add_field_and_yaml_id(fieldptr%get_id(), diag_field_indices(i))
 259      call fileptr%add_buffer_id(fieldptr%buffer_ids(i))
 260      call fileptr%set_file_domain(fieldptr%get_domain(), fieldptr%get_type_of_domain())
 261      call fileptr%init_diurnal_axis(this%diag_axis, this%registered_axis, diag_field_indices(i))
 262      call fileptr%add_axes(axes, this%diag_axis, this%registered_axis, diag_field_indices(i), &
 263        fieldptr%buffer_ids(i), this%FMS_diag_output_buffers)
 264      call fileptr%add_start_time(init_time, this%current_model_time)
 265      call fileptr%set_file_time_ops (fieldptr%diag_field(i), fieldptr%is_static())
 266     enddo
 267   elseif (present(axes)) then !only axes present
 268     do i = 1, size(file_ids)
 269      fileptr => this%FMS_diag_files(file_ids(i))%FMS_diag_file
 270      call fileptr%add_field_and_yaml_id(fieldptr%get_id(), diag_field_indices(i))
 271      call fileptr%add_buffer_id(fieldptr%buffer_ids(i))
 272      call fileptr%init_diurnal_axis(this%diag_axis, this%registered_axis, diag_field_indices(i))
 273      call fileptr%set_file_domain(fieldptr%get_domain(), fieldptr%get_type_of_domain())
 274      call fileptr%add_axes(axes, this%diag_axis, this%registered_axis, diag_field_indices(i), &
 275        fieldptr%buffer_ids(i), this%FMS_diag_output_buffers)
 276      call fileptr%set_file_time_ops (fieldptr%diag_field(i), fieldptr%is_static())
 277     enddo
 278   elseif (present(init_time)) then !only inti time present
 279     do i = 1, size(file_ids)
 280      fileptr => this%FMS_diag_files(file_ids(i))%FMS_diag_file
 281      call fileptr%add_field_and_yaml_id(fieldptr%get_id(), diag_field_indices(i))
 282      call fileptr%add_buffer_id(fieldptr%buffer_ids(i))
 283      call fileptr%add_start_time(init_time, this%current_model_time)
 284      call fileptr%set_file_time_ops (fieldptr%diag_field(i), fieldptr%is_static())
 285     enddo
 286   else !no axis or init time present
 287     do i = 1, size(file_ids)
 288      fileptr => this%FMS_diag_files(file_ids(i))%FMS_diag_file
 289      call fileptr%add_field_and_yaml_id(fieldptr%get_id(), diag_field_indices(i))
 290      call fileptr%add_buffer_id(fieldptr%buffer_ids(i))
 291      call fileptr%set_file_time_ops (fieldptr%diag_field(i), fieldptr%is_static())
 292     enddo
 293   endif
 294
 295   !> Initialize buffer_ids of this field with the diag_field_indices(diag_field_indices)
 296 !! of the sorted variable list
 297   do i = 1, size(fieldptr%buffer_ids)
 298     bufferptr => this%FMS_diag_output_buffers(fieldptr%buffer_ids(i))
 299     call bufferptr%set_field_id(this%registered_variables)
 300     call bufferptr%set_yaml_id(fieldptr%buffer_ids(i))
 301     ! check if diurnal reduction for this buffer and if so set the diurnal sample size
 302     yamlfptr => diag_yaml%diag_fields(fieldptr%buffer_ids(i))
 303     if( yamlfptr%get_var_reduction() .eq. time_diurnal) then
 304       call bufferptr%set_diurnal_sample_size(yamlfptr%get_n_diurnal())
 305     endif
 306     call bufferptr%init_buffer_time(init_time)
 307     call bufferptr%set_next_output(this%FMS_diag_files(file_ids(i))%get_next_output(), fieldptr%is_static())
 308   enddo
 309
 310   nullify (fileptr)
 311   nullify (fieldptr)
 312   deallocate(diag_field_indices)
 313 #endif
 314 end function fms_register_diag_field_obj
 315
 316 !> @brief Registers a scalar field
 317 !! @return field index to be used in subsequent calls to send_data or DIAG_FIELD_NOT_FOUND if the field is not
 318 !! in the diag_table.yaml
 319 INTEGER FUNCTION fms_register_diag_field_scalar(this,module_name, field_name, init_time, &
 320        & long_name, units, missing_value, var_range, standard_name, do_not_log, err_msg,&
 321        & area, volume, realm, multiple_send_data)
 322     class(fmsDiagObject_type),TARGET,INTENT(inout):: this       !< Diaj_obj to fill
 323     CHARACTER(len=*),           INTENT(in) :: module_name   !< Module where the field comes from
 324     CHARACTER(len=*),           INTENT(in) :: field_name    !< Name of the field
 325     TYPE(time_type),  OPTIONAL, INTENT(in) :: init_time     !< Time to start writing data from
 326     CHARACTER(len=*), OPTIONAL, INTENT(in) :: long_name     !< Long_name to add as a variable attribute
 327     CHARACTER(len=*), OPTIONAL, INTENT(in) :: units         !< Units to add as a variable_attribute
 328     CHARACTER(len=*), OPTIONAL, INTENT(in) :: standard_name !< Standard_name to name the variable in the file
 329     CLASS(*),         OPTIONAL, INTENT(in) :: missing_value !< Missing value to add as a variable attribute
 330     CLASS(*),         OPTIONAL, INTENT(in) :: var_range(:)  !< Range to add a variable attribute
 331     LOGICAL,          OPTIONAL, INTENT(in) :: do_not_log    !< If TRUE, field information is not logged
 332     CHARACTER(len=*), OPTIONAL, INTENT(out):: err_msg       !< Error_msg from call
 333     INTEGER,          OPTIONAL, INTENT(in) :: area          !< Id of the area field
 334     INTEGER,          OPTIONAL, INTENT(in) :: volume        !< Id of the volume field
 335     CHARACTER(len=*), OPTIONAL, INTENT(in) :: realm         !< String to set as the modeling_realm attribute
 336     LOGICAL,          OPTIONAL, INTENT(in) :: multiple_send_data !< .True. if send data is called, multiple times
 337                                                                  !! for the same time
 338
 339 #ifndef use_yaml
 340 fms_register_diag_field_scalar=DIAG_FIELD_NOT_FOUND
 341 CALL MPP_ERROR(FATAL,"You can not use the modern diag manager without compiling with -Duse_yaml")
 342 #else
 343     fms_register_diag_field_scalar = this%register(&
 344       & module_name, field_name, init_time=init_time, &
 345       & longname=long_name, units=units, missing_value=missing_value, varrange=var_range, &
 346       & standname=standard_name, do_not_log=do_not_log, err_msg=err_msg, &
 347       & area=area, volume=volume, realm=realm, multiple_send_data=multiple_send_data)
 348 #endif
 349 end function fms_register_diag_field_scalar
 350
 351 !> @brief Registers an array field
 352 !! @return field index to be used in subsequent calls to send_data or DIAG_FIELD_NOT_FOUND if the field is not
 353 !! in the diag_table.yaml
 354 INTEGER FUNCTION fms_register_diag_field_array(this, module_name, field_name, axes, init_time, &
 355        & long_name, units, missing_value, var_range, mask_variant, standard_name, verbose,&
 356        & do_not_log, err_msg, interp_method, tile_count, area, volume, realm, &
 357        & multiple_send_data)
 358     class(fmsDiagObject_type),TARGET,INTENT(inout):: this       !< Diaj_obj to fill
 359     CHARACTER(len=*),           INTENT(in) :: module_name   !< Module where the field comes from
 360     CHARACTER(len=*),           INTENT(in) :: field_name    !< Name of the field
 361     INTEGER,                    INTENT(in) :: axes(:)       !< Ids corresponding to the variable axis
 362     TYPE(time_type),  OPTIONAL, INTENT(in) :: init_time     !< Time to start writing data from
 363     CHARACTER(len=*), OPTIONAL, INTENT(in) :: long_name     !< Long_name to add as a variable attribute
 364     CHARACTER(len=*), OPTIONAL, INTENT(in) :: units         !< Units to add as a variable_attribute
 365     CLASS(*),         OPTIONAL, INTENT(in) :: missing_value !< Missing value to add as a variable attribute
 366     CLASS(*),         OPTIONAL, INTENT(in) :: var_range(:)  !< Range to add a variable attribute
 367     LOGICAL,          OPTIONAL, INTENT(in) :: mask_variant  !< .True. if mask changes over time
 368     CHARACTER(len=*), OPTIONAL, INTENT(in) :: standard_name !< Standard_name to name the variable in the file
 369     LOGICAL,          OPTIONAL, INTENT(in) :: verbose       !< Print more information
 370     LOGICAL,          OPTIONAL, INTENT(in) :: do_not_log    !< If TRUE, field information is not logged
 371     CHARACTER(len=*), OPTIONAL, INTENT(out):: err_msg       !< Error_msg from call
 372     CHARACTER(len=*), OPTIONAL, INTENT(in) :: interp_method !< The interp method to be used when
 373                                                             !! regridding the field in post-processing.
 374                                                             !! Valid options are "conserve_order1",
 375                                                             !! "conserve_order2", and "none".
 376     INTEGER,          OPTIONAL, INTENT(in) :: tile_count    !< The current tile number
 377     INTEGER,          OPTIONAL, INTENT(in) :: area          !< Id of the area field
 378     INTEGER,          OPTIONAL, INTENT(in) :: volume        !< Id of the volume field
 379     CHARACTER(len=*), OPTIONAL, INTENT(in) :: realm         !< String to set as the modeling_realm attribute
 380     LOGICAL,          OPTIONAL, INTENT(in) :: multiple_send_data !< .True. if send data is called, multiple times
 381                                                                  !! for the same time
 382
 383
 384 #ifndef use_yaml
 385 fms_register_diag_field_array=DIAG_FIELD_NOT_FOUND
 386 CALL MPP_ERROR(FATAL,"You can not use the modern diag manager without compiling with -Duse_yaml")
 387 #else
 388     fms_register_diag_field_array = this%register( &
 389       & module_name, field_name, init_time=init_time, &
 390       & axes=axes, longname=long_name, units=units, missing_value=missing_value, varrange=var_range, &
 391       & mask_variant=mask_variant, standname=standard_name, do_not_log=do_not_log, err_msg=err_msg, &
 392       & interp_method=interp_method, tile_count=tile_count, area=area, volume=volume, realm=realm, &
 393       & multiple_send_data=multiple_send_data)
 394 #endif
 395 end function fms_register_diag_field_array
 396
 397 !> @brief Return field index for subsequent call to send_data.
 398 !! @return field index to be used in subsequent calls to send_data or DIAG_FIELD_NOT_FOUND if the field is not
 399 !! in the diag_table.yaml
 400 INTEGER FUNCTION fms_register_static_field(this, module_name, field_name, axes, long_name, units,&
 401        & missing_value, range, mask_variant, standard_name, DYNAMIC, do_not_log, interp_method,&
 402        & tile_count, area, volume, realm)
 403     class(fmsDiagObject_type),TARGET,INTENT(inout):: this       !< Diaj_obj to fill
 404     CHARACTER(len=*),                         INTENT(in) :: module_name   !< Name of the module, the field is on
 405     CHARACTER(len=*),                         INTENT(in) :: field_name    !< Name of the field
 406     INTEGER,          DIMENSION(:),           INTENT(in) :: axes          !< Axes_id of the field
 407     CHARACTER(len=*),               OPTIONAL, INTENT(in) :: long_name     !< Longname to be added as a attribute
 408     CHARACTER(len=*),               OPTIONAL, INTENT(in) :: units         !< Units to be added as a attribute
 409     CHARACTER(len=*),               OPTIONAL, INTENT(in) :: standard_name !< Standard name to be added as a attribute
 410     CLASS(*),                       OPTIONAL, INTENT(in) :: missing_value !< Missing value to be added as a attribute
 411     CLASS(*),                       OPTIONAL, INTENT(in) :: range(:)      !< Range to be added as a attribute
 412     LOGICAL,                        OPTIONAL, INTENT(in) :: mask_variant  !< .True. if mask changes over time
 413     LOGICAL,                        OPTIONAL, INTENT(in) :: DYNAMIC       !< Flag indicating if the field is dynamic
 414     LOGICAL,                        OPTIONAL, INTENT(in) :: do_not_log    !< if TRUE, field information is not logged
 415     CHARACTER(len=*),               OPTIONAL, INTENT(in) :: interp_method !< The interp method to be used when
 416                                                                           !! regridding the field in post-processing
 417                                                                           !! Valid options are "conserve_order1",
 418                                                                           !! "conserve_order2", and "none".
 419     INTEGER,                        OPTIONAL, INTENT(in) :: tile_count    !! Number of tiles
 420     INTEGER,                        OPTIONAL, INTENT(in) :: area          !< Field ID for the area field associated
 421                                                                           !! with this field
 422     INTEGER,                        OPTIONAL, INTENT(in) :: volume        !< Field ID for the volume field associated
 423                                                                           !! with this field
 424     CHARACTER(len=*),               OPTIONAL, INTENT(in) :: realm         !< String to set as the value to the
 425                                                                           !! modeling_realm attribute
 426
 427 #ifndef use_yaml
 428 fms_register_static_field=DIAG_FIELD_NOT_FOUND
 429 CALL MPP_ERROR(FATAL,"You can not use the modern diag manager without compiling with -Duse_yaml")
 430 #else
 431   !TODO The register_static_field interface does not have the capabiliy to register a variable as a "scalar"
 432   !     since the axes argument is required, this forced model code to pass in a null_axis_id as an argument
 433   if (size(axes) .eq. 1 .and. axes(1) .eq. null_axis_id) then
 434     ! If they are passing in the null_axis_ids, ignore the `axes` argument
 435     fms_register_static_field = this%register( &
 436       & module_name, field_name, &
 437       & longname=long_name, units=units, missing_value=missing_value, varrange=range, &
 438       & mask_variant=mask_variant, do_not_log=do_not_log, interp_method=interp_method, tile_count=tile_count, &
 439       & standname=standard_name, area=area, volume=volume, realm=realm, &
 440       & static=.true.)
 441   else
 442     fms_register_static_field = this%register( &
 443       & module_name, field_name, axes=axes, &
 444       & longname=long_name, units=units, missing_value=missing_value, varrange=range, &
 445       & mask_variant=mask_variant, do_not_log=do_not_log, interp_method=interp_method, tile_count=tile_count, &
 446       & standname=standard_name, area=area, volume=volume, realm=realm, &
 447       & static=.true.)
 448   endif
 449 #endif
 450 end function fms_register_static_field
 451
 452 !> @brief Wrapper for the register_diag_axis subroutine. This is needed to keep the diag_axis_init
 453 !! interface the same
 454 !> @return Axis id
 455 FUNCTION fms_diag_axis_init(this, axis_name, axis_data, units, cart_name, axis_length, long_name, direction,&
 456   & set_name, edges, Domain, Domain2, DomainU, aux, req, tile_count, domain_position ) &
 457   & result(id)
 458
 459   class(fmsDiagObject_type),TARGET,INTENT(inout):: this       !< Diaj_obj to fill
 460   CHARACTER(len=*),   INTENT(in)           :: axis_name       !< Name of the axis
 461   CLASS(*),           INTENT(in)           :: axis_data(:)    !< Array of coordinate values
 462   CHARACTER(len=*),   INTENT(in)           :: units           !< Units for the axis
 463   CHARACTER(len=1),   INTENT(in)           :: cart_name       !< Cartesian axis ("X", "Y", "Z", "T", "U", "N")
 464   integer,            intent(in)           :: axis_length     !< The length of the axis size(axis_data(:))
 465   CHARACTER(len=*),   INTENT(in), OPTIONAL :: long_name       !< Long name for the axis.
 466   CHARACTER(len=*),   INTENT(in), OPTIONAL :: set_name        !< Name of the parent axis, if it is a subaxis
 467   INTEGER,            INTENT(in), OPTIONAL :: direction       !< Indicates the direction of the axis
 468   INTEGER,            INTENT(in), OPTIONAL :: edges           !< Axis ID for the previously defined "edges axis"
 469   TYPE(domain1d),     INTENT(in), OPTIONAL :: Domain          !< 1D domain
 470   TYPE(domain2d),     INTENT(in), OPTIONAL :: Domain2         !< 2D domain
 471   TYPE(domainUG),     INTENT(in), OPTIONAL :: DomainU         !< Unstructured domain
 472   CHARACTER(len=*),   INTENT(in), OPTIONAL :: aux             !< Auxiliary name, can only be <TT>geolon_t</TT>
 473                                                                !! or <TT>geolat_t</TT>
 474   CHARACTER(len=*),   INTENT(in), OPTIONAL :: req             !< Required field names.
 475   INTEGER,            INTENT(in), OPTIONAL :: tile_count      !< Number of tiles
 476   INTEGER,            INTENT(in), OPTIONAL :: domain_position !< Domain position, "NORTH" or "EAST"
 477   integer :: id
 478
 479 #ifndef use_yaml
 480 id = diag_null
 481 CALL MPP_ERROR(FATAL,"You can not use the modern diag manager without compiling with -Duse_yaml")
 482 #else
 483   CHARACTER(len=:),   ALLOCATABLE :: edges_name !< Name of the edges
 484
 485   this%registered_axis = this%registered_axis + 1
 486
 487   if (this%registered_axis > max_axes) call mpp_error(FATAL, &
 488     &"diag_axis_init: max_axes exceeded, increase via diag_manager_nml")
 489
 490   allocate(fmsDiagFullAxis_type :: this%diag_axis(this%registered_axis)%axis)
 491
 492   select type (axis => this%diag_axis(this%registered_axis)%axis )
 493   type is (fmsDiagFullAxis_type)
 494     if(present(edges)) then
 495       if (edges < 0 .or. edges > this%registered_axis) &
 496         call mpp_error(FATAL, "diag_axit_init: The edge axis has not been defined. "&
 497                                "Call diag_axis_init for the edge axis first")
 498       select type (edges_axis => this%diag_axis(edges)%axis)
 499       type is (fmsDiagFullAxis_type)
 500         edges_name = edges_axis%get_axis_name()
 501         call axis%set_edges(edges_name, edges)
 502       end select
 503     endif
 504     call axis%register(axis_name, axis_data, units, cart_name, long_name=long_name, &
 505       & direction=direction, set_name=set_name, Domain=Domain, Domain2=Domain2, DomainU=DomainU, aux=aux, &
 506       & req=req, tile_count=tile_count, domain_position=domain_position, axis_length=axis_length)
 507
 508     id = this%registered_axis
 509     call axis%set_axis_id(id)
 510   end select
 511 #endif
 512 end function fms_diag_axis_init
 513
 514 !> Accepts data from the send_data functions.  If this is in an openmp region with more than
 515 !! one thread, the data is buffered in the field object and processed later.  If only a single thread
 516 !! is being used, then the processing can be done and stored in the buffer object.  The hope is that
 517 !! the increase in memory footprint related to buffering can be handled by the shared memory of the
 518 !! multithreaded case.
 519 !! \note If some of the diag manager is offloaded in the future, then it should be treated similarly
 520 !! to the multi-threaded option for processing later
 521 logical function fms_diag_accept_data (this, diag_field_id, field_data, mask, rmask, &
 522                                        time, is_in, js_in, ks_in, &
 523                                        ie_in, je_in, ke_in, weight, err_msg)
 524   class(fmsDiagObject_type),TARGET,      INTENT(inout)          :: this          !< Diaj_obj to fill
 525   INTEGER,                               INTENT(in)             :: diag_field_id !< The ID of the diag field
 526   CLASS(*), DIMENSION(:,:,:,:),          INTENT(in)             :: field_data    !< The data for the diag_field
 527   LOGICAL,  allocatable,                 INTENT(in)             :: mask(:,:,:,:) !< Logical mask indicating the grid
 528                                                                                  !! points to mask (null if no mask)
 529   CLASS(*), allocatable,                 INTENT(in)             :: rmask(:,:,:,:)!< real mask indicating the grid
 530                                                                                  !! points to mask (null if no mask)
 531   CLASS(*),                              INTENT(in),   OPTIONAL :: weight        !< The weight used for averaging
 532   TYPE (time_type),                      INTENT(in),   OPTIONAL :: time          !< The current time
 533   INTEGER,                               INTENT(in),   OPTIONAL :: is_in, js_in, ks_in !< Starting indices
 534   INTEGER,                               INTENT(in),   OPTIONAL :: ie_in, je_in, ke_in !< Ending indices
 535   CHARACTER(len=*),                      INTENT(out),  OPTIONAL :: err_msg       !< An error message returned
 536
 537   integer                                  :: is, js, ks      !< Starting indicies of the field_data
 538   integer                                  :: ie, je, ke      !< Ending indicies of the field_data
 539   integer                                  :: omp_num_threads !< Number of openmp threads
 540   integer                                  :: omp_level       !< The openmp active level
 541   logical                                  :: buffer_the_data !< True if the user selects to buffer the data and run
 542                                                               !! the calculationslater.  \note This is experimental
 543   character(len=128)                       :: error_string    !< Store error text
 544   logical                                  :: data_buffer_is_allocated !< .true. if the data buffer is allocated
 545   character(len=256)                       :: field_info      !< String holding info about the field to append to the
 546                                                               !! error message
 547   logical, allocatable, dimension(:,:,:,:) :: oor_mask        !< Out of range mask
 548   real(kind=r8_kind)                       :: field_weight    !< Weight to use when averaging (it will be converted
 549                                                               !! based on the type of field_data when doing the math)
 550   type(fmsDiagIbounds_type)                :: bounds          !< Bounds (starting ending indices) for the field
 551   logical                                  :: has_halos       !< .True. if field_data contains halos
 552   logical                                  :: using_blocking  !< .True. if field_data is passed in blocks
 553 #ifndef use_yaml
 554 CALL MPP_ERROR(FATAL,"You can not use the modern diag manager without compiling with -Duse_yaml")
 555 #else
 556
 557   !TODO this%FMS_diag_fields(diag_field_id) should be a pointer!
 558   field_info = " Check send data call for field:"//trim(this%FMS_diag_fields(diag_field_id)%get_varname())//&
 559     " and module:"//trim(this%FMS_diag_fields(diag_field_id)%get_modname())
 560
 561   !< Check if time should be present for this field
 562   if (.not.this%FMS_diag_fields(diag_field_id)%is_static() .and. .not.present(time)) &
 563     call mpp_error(FATAL, "Time must be present if the field is not static. "//trim(field_info))
 564
 565   !< Set the field_weight. If "weight" is not present it will be set to 1.0_r8_kind
 566   field_weight = set_weight(weight)
 567
 568   !< Check that the indices are present in the correct combination
 569   error_string = check_indices_order(is_in, ie_in, js_in, je_in)
 570   if (trim(error_string) .ne. "") call mpp_error(FATAL, trim(error_string)//". "//trim(field_info))
 571
 572   using_blocking = .false.
 573   if ((present(is_in) .and. .not. present(ie_in)) .or. (present(js_in) .and. .not. present(je_in))) &
 574     using_blocking = .true.
 575
 576   has_halos = .false.
 577   if ((present(is_in) .and. present(ie_in)) .or. (present(js_in) .and. present(je_in))) &
 578     has_halos = .true.
 579
 580   !< If the field has `mask_variant=.true.`, check that mask OR rmask are present
 581   if (this%FMS_diag_fields(diag_field_id)%is_mask_variant()) then
 582     if (.not. allocated(mask) .and. .not. allocated(rmask)) call mpp_error(FATAL, &
 583       "The field was registered with mask_variant, but mask or rmask are not present in the send_data call. "//&
 584       trim(field_info))
 585   endif
 586
 587   !< Check that mask and rmask are not both present
 588   if (allocated(mask) .and. allocated(rmask)) call mpp_error(FATAL, &
 589     "mask and rmask are both present in the send_data call. "//&
 590     trim(field_info))
 591
 592   !< Create the oor_mask based on the "mask" and "rmask" arguments
 593   oor_mask = init_mask(rmask, mask, field_data)
 594
 595   !> Does the user want to push off calculations until send_diag_complete?
 596   buffer_the_data = .false.
 597
 598   !> initialize the number of threads and level to be 0
 599   omp_num_threads = 0
 600   omp_level = 0
 601 #if defined(_OPENMP)
 602   omp_num_threads = omp_get_num_threads()
 603   omp_level = omp_get_level()
 604   buffer_the_data = (omp_num_threads > 1 .AND. omp_level > 0)
 605 #endif
 606
 607   !> Calculate the i,j,k start and end
 608   ! If is, js, or ks not present default them to 1
 609   is = 1
 610   js = 1
 611   ks = 1
 612   IF ( PRESENT(is_in) ) is = is_in
 613   IF ( PRESENT(js_in) ) js = js_in
 614   IF ( PRESENT(ks_in) ) ks = ks_in
 615   ie = is+SIZE(field_data, 1)-1
 616   je = js+SIZE(field_data, 2)-1
 617   ke = ks+SIZE(field_data, 3)-1
 618   IF ( PRESENT(ie_in) ) ie = ie_in
 619   IF ( PRESENT(je_in) ) je = je_in
 620   IF ( PRESENT(ke_in) ) ke = ke_in
 621
 622   if (.not. buffer_the_data .and. using_blocking) then
 623     ! If running with only 1 thread and using blocking, check if the data was sent in blocks
 624     ! if it is, then buffer the data
 625     buffer_the_data = check_for_slices(this%FMS_diag_fields(diag_field_id), this%diag_axis, &
 626       shape(field_data))
 627   endif
 628
 629   !< If send data is called multiple times, buffer the data
 630   !! This is so that the other reduction methods work and just averaging
 631   if (this%FMS_diag_fields(diag_field_id)%get_multiple_send_data()) &
 632     buffer_the_data = .true.
 633
 634   !If this is true, buffer data
 635   main_if: if (buffer_the_data) then
 636 !> Only 1 thread allocates the output buffer and sets set_math_needs_to_be_done
 637 !$omp critical
 638
 639     if (present(time)) call this%update_current_model_time(time)
 640
 641     !< These set_* calls need to be done inside an omp_critical to avoid any race conditions
 642     !! and allocation issues
 643     if(has_halos) call this%FMS_diag_fields(diag_field_id)%set_halo_present()
 644
 645     !< Set the variable type based off passed in field data
 646     if(.not. this%FMS_diag_fields(diag_field_id)%has_vartype()) &
 647       call this%FMS_diag_fields(diag_field_id)%set_type(field_data(1,1,1,1))
 648
 649     if (allocated(mask) .or. allocated(rmask)) then
 650       call this%FMS_diag_fields(diag_field_id)%set_var_is_masked(.True.)
 651     else
 652       call this%FMS_diag_fields(diag_field_id)%set_var_is_masked(.False.)
 653     endif
 654
 655     if (.not. this%FMS_diag_fields(diag_field_id)%is_data_buffer_allocated()) then
 656       data_buffer_is_allocated = &
 657         this%FMS_diag_fields(diag_field_id)%allocate_data_buffer(field_data, this%diag_axis)
 658       if(.not. this%FMS_diag_fields(diag_field_id)%has_mask_allocated()) &
 659         call this%FMS_diag_fields(diag_field_id)%allocate_mask(oor_mask, this%diag_axis)
 660     endif
 661     call this%FMS_diag_fields(diag_field_id)%set_send_data_time(time)
 662     call this%FMS_diag_fields(diag_field_id)%set_data_buffer_is_allocated(.TRUE.)
 663     call this%FMS_diag_fields(diag_field_id)%set_math_needs_to_be_done(.TRUE.)
 664 !$omp end critical
 665     call this%FMS_diag_fields(diag_field_id)%set_data_buffer(field_data, oor_mask, field_weight, &
 666                                                              is, js, ks, ie, je, ke)
 667     fms_diag_accept_data = .TRUE.
 668     return
 669   else
 670     if (present(time)) call this%update_current_model_time(time)
 671
 672     !< At this point if we are no longer in an openmp region or running with 1 thread
 673     !! so it is safe to have these set_* calls
 674     if(has_halos) call this%FMS_diag_fields(diag_field_id)%set_halo_present()
 675
 676     !< Set the variable type based off passed in field data
 677     if(.not. this%FMS_diag_fields(diag_field_id)%has_vartype()) &
 678       call this%FMS_diag_fields(diag_field_id)%set_type(field_data(1,1,1,1))
 679
 680     if (allocated(mask) .or. allocated(rmask)) then
 681       call this%FMS_diag_fields(diag_field_id)%set_var_is_masked(.True.)
 682     else
 683       call this%FMS_diag_fields(diag_field_id)%set_var_is_masked(.False.)
 684     endif
 685
 686     error_string = bounds%set_bounds(field_data, is, ie, js, je, ks, ke, has_halos)
 687     if (trim(error_string) .ne. "") call mpp_error(FATAL, trim(error_string)//". "//trim(field_info))
 688
 689     call this%allocate_diag_field_output_buffers(field_data, diag_field_id)
 690     error_string = this%fms_diag_do_reduction(field_data, diag_field_id, oor_mask, field_weight, &
 691       bounds, using_blocking, Time=Time)
 692     if (trim(error_string) .ne. "") call mpp_error(FATAL, trim(error_string)//". "//trim(field_info))
 693     call this%FMS_diag_fields(diag_field_id)%set_math_needs_to_be_done(.FALSE.)
 694     if(.not. this%FMS_diag_fields(diag_field_id)%has_mask_allocated()) &
 695       call this%FMS_diag_fields(diag_field_id)%allocate_mask(oor_mask)
 696     call this%FMS_diag_fields(diag_field_id)%set_mask(oor_mask, field_info)
 697     return
 698   end if main_if
 699   !> Return false if nothing is done
 700   fms_diag_accept_data = .FALSE.
 701   return
 702 #endif
 703 end function fms_diag_accept_data
 704
 705 !< @brief Do the math for all the buffers
 706 subroutine do_buffer_math(this)
 707   class(fmsDiagObject_type), target, intent (inout) :: this      !< The diag object
 708
 709 #ifdef use_yaml
 710   integer :: i !< For do loops
 711   integer :: ifile !< For file loops
 712   integer :: ifield !< For field loops
 713
 714   class(fmsDiagFileContainer_type), pointer :: diag_file !< Pointer to this%FMS_diag_files(i) (for convenience
 715   class(fmsDiagField_type), pointer :: diag_field !< Pointer to this%FMS_diag_files(i)%diag_field(j)
 716   logical :: math !< True if the math functions need to be called using the data buffer,
 717   !! False if the math functions were done in accept_data
 718   integer, dimension(:), allocatable :: file_field_ids !< Array of field IDs for a file
 719   class(*), pointer :: input_data_buffer(:,:,:,:)
 720   character(len=128) :: error_string
 721   type(fmsDiagIbounds_type) :: bounds
 722   integer, dimension(:), allocatable :: file_ids !< Array of file IDs for a field
 723   logical, parameter :: DEBUG_SC = .false. !< turn on output for debugging
 724
 725     !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
 726   !! In the future, this may be parallelized for offloading
 727   ! loop through each field
 728   field_loop: do ifield = 1, size(this%FMS_diag_fields)
 729     diag_field => this%FMS_diag_fields(ifield)
 730     if(.not. diag_field%is_registered()) cycle
 731     if(DEBUG_SC) call mpp_error(NOTE, "fms_diag_send_complete:: var: "//diag_field%get_varname())
 732     ! get files the field is in
 733     allocate (file_ids(size(diag_field%get_file_ids() )))
 734     file_ids = diag_field%get_file_ids()
 735     math = diag_field%get_math_needs_to_be_done()
 736     ! if doing math loop through each file for given field
 737     doing_math: if (size(file_ids) .ge. 1 .and. math) then
 738       ! Check if buffer alloc'd
 739       has_input_buff: if (diag_field%has_input_data_buffer()) then
 740         call diag_field%prepare_data_buffer()
 741         input_data_buffer => diag_field%get_data_buffer()
 742         ! reset bounds, allocate output buffer, and update it with reduction
 743         call bounds%reset_bounds_from_array_4D(input_data_buffer)
 744         call this%allocate_diag_field_output_buffers(input_data_buffer, ifield)
 745         error_string = this%fms_diag_do_reduction(input_data_buffer, ifield, &
 746                               diag_field%get_mask(), diag_field%get_weight(), &
 747                               bounds, .False., Time=diag_field%get_send_data_time())
 748         call diag_field%init_data_buffer()
 749         if (trim(error_string) .ne. "") call mpp_error(FATAL, "Field:"//trim(diag_field%get_varname()//&
 750                                                        " -"//trim(error_string)))
 751       else
 752         call mpp_error(FATAL, "diag_send_complete:: no input buffer allocated for field"//diag_field%get_longname())
 753       endif has_input_buff
 754     endif doing_math
 755     call diag_field%set_math_needs_to_be_done(.False.)
 756     !> Clean up, clean up, everybody do your share
 757     if (allocated(file_ids)) deallocate(file_ids)
 758     if (associated(diag_field)) nullify(diag_field)
 759   enddo field_loop
 760 #endif
 761 end subroutine do_buffer_math
 762
 763 !> @brief Loops through all the files, open the file, writes out axis and
 764 !! variable metadata and data when necessary.
 765 subroutine fms_diag_send_complete(this, time_step)
 766   class(fmsDiagObject_type), target, intent (inout) :: this      !< The diag object
 767   TYPE (time_type),                  INTENT(in)     :: time_step !< The time_step
 768
 769 #ifndef use_yaml
 770 CALL MPP_ERROR(FATAL,"You can not use the modern diag manager without compiling with -Duse_yaml")
 771 #else
 772   call this%do_buffer_math()
 773   call this%fms_diag_do_io()
 774 #endif
 775
 776 end subroutine fms_diag_send_complete
 777
 778 !> @brief Loops through all the files, open the file, writes out axis and
 779 !! variable metadata and data when necessary.
 780 !! TODO: passing in the saved mask from the field obj to diag_reduction_done_wrapper
 781 !! for performance
 782 subroutine fms_diag_do_io(this, end_time)
 783   class(fmsDiagObject_type), target, intent(inout)  :: this          !< The diag object
 784   type(time_type), optional, target, intent(in)     :: end_time      !< the model end_time
 785 #ifdef use_yaml
 786   integer :: i !< For do loops
 787   class(fmsDiagFileContainer_type), pointer :: diag_file !< Pointer to this%FMS_diag_files(i) (for convenience)
 788   class(fmsDiagOutputBuffer_type), pointer  :: diag_buff !< pointer to output buffers iterated in buff_loop
 789   class(fmsDiagField_type), pointer         :: diag_field !< pointer to output buffers iterated in buff_loop
 790   class(DiagYamlFilesVar_type), pointer     :: field_yaml !< Pointer to a field from yaml fields
 791   TYPE (time_type),                 pointer :: model_time!< The current model time
 792   integer, allocatable                      :: buff_ids(:) !< ids for output buffers to loop through
 793   integer                                   :: ibuff !< buffer index
 794   logical :: file_is_opened_this_time_step !< True if the file was opened in this time_step
 795                                            !! If true the metadata will need to be written
 796   logical :: force_write !< force the last write if at end of run
 797   logical :: finish_writing !< true if finished writing for all the fields
 798   logical :: has_mask !< whether we have a mask
 799   logical, parameter :: DEBUG_REDUCT = .false. !< enables debugging output
 800   class(*), allocatable :: missing_val !< netcdf missing value for a given field
 801   real(r8_kind) :: mval !< r8 copy of missing value
 802   character(len=128) :: error_string !< outputted error string from reducti
 803
 804   force_write = .false.
 805   if (present (end_time)) then
 806     force_write = .true.
 807     model_time => end_time
 808   else
 809     model_time => this%current_model_time
 810   endif
 811
 812   do i = 1, size(this%FMS_diag_files)
 813     diag_file => this%FMS_diag_files(i)
 814
 815     !< Go away if the file is a subregional file and the current PE does not have any data for it
 816     if (.not. diag_file%writing_on_this_pe()) cycle
 817
 818     call diag_file%open_diag_file(model_time, file_is_opened_this_time_step)
 819     if (file_is_opened_this_time_step) then
 820       call diag_file%write_global_metadata()
 821       call diag_file%write_axis_metadata(this%diag_axis)
 822       call diag_file%write_time_metadata()
 823       call diag_file%write_field_metadata(this%FMS_diag_fields, this%diag_axis)
 824       call diag_file%write_axis_data(this%diag_axis)
 825       call diag_file%increase_unlim_dimension_level()
 826     endif
 827
 828     finish_writing = diag_file%is_time_to_write(model_time, this%FMS_diag_output_buffers)
 829
 830     ! finish reduction method if its time to write
 831     buff_ids = diag_file%FMS_diag_file%get_buffer_ids()
 832     ! loop through the buffers and finish reduction if needed
 833     buff_loop: do ibuff=1, SIZE(buff_ids)
 834       diag_buff => this%FMS_diag_output_buffers(buff_ids(ibuff))
 835       field_yaml => diag_yaml%diag_fields(diag_buff%get_yaml_id())
 836       diag_field => this%FMS_diag_fields(diag_buff%get_field_id())
 837
 838       ! Go away if there is no data to write
 839       if (.not. diag_buff%is_there_data_to_write()) cycle
 840
 841       if ( diag_buff%is_time_to_finish_reduction(end_time)) then
 842         ! sets missing value
 843         mval = diag_field%find_missing_value(missing_val)
 844         ! time_average and greater values all involve averaging so need to be "finished" before written
 845         if( field_yaml%has_var_reduction()) then
 846           if( field_yaml%get_var_reduction() .ge. time_average) then
 847             if(DEBUG_REDUCT)call mpp_error(NOTE, "fms_diag_do_io:: finishing reduction for "//diag_field%get_longname())
 848             error_string = diag_buff%diag_reduction_done_wrapper( &
 849                                     field_yaml%get_var_reduction(), &
 850                                    mval, diag_field%get_var_is_masked(), diag_field%get_mask_variant())
 851           endif
 852         endif
 853         call diag_file%write_field_data(diag_field, diag_buff)
 854         call diag_buff%set_next_output(diag_file%get_next_next_output())
 855       endif
 856       nullify(diag_buff)
 857       nullify(field_yaml)
 858     enddo buff_loop
 859     deallocate(buff_ids)
 860
 861     if (finish_writing) then
 862       call diag_file%write_time_data()
 863       call diag_file%update_next_write(model_time)
 864       call diag_file%update_current_new_file_freq_index(model_time)
 865       call diag_file%increase_unlim_dimension_level()
 866       if (diag_file%is_time_to_close_file(model_time)) call diag_file%close_diag_file(this%FMS_diag_output_buffers, &
 867         diag_fields = this%FMS_diag_fields)
 868     else if (force_write) then
 869       call diag_file%write_time_data()
 870       call diag_file%close_diag_file(this%FMS_diag_output_buffers, diag_fields = this%FMS_diag_fields)
 871     endif
 872   enddo
 873 #endif
 874 end subroutine fms_diag_do_io
 875
 876 !> @brief Computes average, min, max, rms error, etc.
 877 !! based on the specified reduction method for the field.
 878 !> @return Empty string if successful, error message if it fails
 879 function fms_diag_do_reduction(this, field_data, diag_field_id, oor_mask, weight, &
 880   bounds, using_blocking, time) &
 881   result(error_msg)
 882   class(fmsDiagObject_type), intent(inout), target:: this                !< Diag Object
 883   class(*),                  intent(in)           :: field_data(:,:,:,:) !< Field data
 884   integer,                   intent(in)           :: diag_field_id       !< ID of the input field
 885   logical,                   intent(in), target   :: oor_mask(:,:,:,:)   !< mask
 886   real(kind=r8_kind),        intent(in)           :: weight              !< Must be a updated weight
 887   type(fmsDiagIbounds_type), intent(in)           :: bounds              !< Bounds for the field
 888   logical,                   intent(in)           :: using_blocking      !< .True. if field data is passed
 889                                                                          !! in blocks
 890   type(time_type),           intent(in), optional :: time                !< Current time
 891
 892   character(len=150)         :: error_msg          !< Error message to check
 893   !TODO Mostly everything
 894 #ifdef use_yaml
 895   type(fmsDiagField_type),          pointer :: field_ptr      !< Pointer to the field's object
 896   type(fmsDiagOutputBuffer_type),   pointer :: buffer_ptr     !< Pointer to the field's buffer
 897   class(fmsDiagFileContainer_type), pointer :: file_ptr       !< Pointer to the field's file
 898   type(diagYamlFilesVar_type),      pointer :: field_yaml_ptr !< Pointer to the field's yaml
 899
 900   integer                   :: reduction_method   !< Integer representing a reduction method
 901   integer                   :: ids                !< For looping through buffer ids
 902   integer                   :: buffer_id          !< Id of the buffer
 903   integer                   :: file_id            !< File id
 904   integer, pointer          :: axis_ids(:)        !< Axis ids for the buffer
 905   logical                   :: is_subregional     !< .True. if the buffer is subregional
 906   logical                   :: reduced_k_range    !< .True. is the field is only outputing a section
 907                                                   !! of the z dimension
 908   type(fmsDiagIbounds_type) :: bounds_in          !< Starting and ending indices of the input field_data
 909   type(fmsDiagIbounds_type) :: bounds_out         !< Starting and ending indices of the output buffer
 910   integer                   :: i                  !< For looping through axid ids
 911   integer                   :: sindex             !< Starting index of a subregion
 912   integer                   :: eindex             !< Ending index of a subregion
 913   integer                   :: compute_idx(2)     !< Starting and Ending of the compute domain
 914   character(len=1)          :: cart_axis          !< Cartesian axis of the axis
 915   logical                   :: block_in_subregion !< .True. if the current block is part of the subregion
 916   integer                   :: starting           !< Starting index of the subregion relative to the compute domain
 917   integer                   :: ending             !< Ending index of the subregion relative to the compute domain
 918   real(kind=r8_kind)        :: missing_value      !< Missing_value for data points that are masked
 919                                                   !! This will obtained as r8 and converted to the right type as
 920                                                   !! needed. This is to avoid yet another select type ...
 921
 922   !TODO mostly everything
 923   field_ptr => this%FMS_diag_fields(diag_field_id)
 924   if (field_ptr%has_missing_value()) then
 925     select type (missing_val => field_ptr%get_missing_value(r8))
 926     type is (real(kind=r8_kind))
 927       missing_value = missing_val
 928     class default
 929       call mpp_error(FATAl, "The missing value for the field:"//trim(field_ptr%get_varname())//&
 930         &" was not allocated to the correct type. This shouldn't have happened")
 931     end select
 932   else
 933     select type (missing_val => get_default_missing_value(r8))
 934     type is (real(kind=r8_kind))
 935       missing_value = missing_val
 936     class default
 937       call mpp_error(FATAl, "The missing value for the field:"//trim(field_ptr%get_varname())//&
 938         &" was not allocated to the correct type. This shouldn't have happened")
 939     end select
 940   endif
 941
 942   buffer_loop: do ids = 1, size(field_ptr%buffer_ids)
 943     error_msg = ""
 944     buffer_id = this%FMS_diag_fields(diag_field_id)%buffer_ids(ids)
 945     file_id = this%FMS_diag_fields(diag_field_id)%file_ids(ids)
 946
 947     !< Gather all the objects needed for the buffer
 948     field_yaml_ptr => field_ptr%diag_field(ids)
 949     buffer_ptr     => this%FMS_diag_output_buffers(buffer_id)
 950     file_ptr       => this%FMS_diag_files(file_id)
 951
 952     !< Go away if the file is a subregional file and the current PE does not have any data for it
 953     if (.not. file_ptr%writing_on_this_pe()) cycle
 954
 955     !< Go away if finished doing math for this buffer
 956     if (buffer_ptr%is_done_with_math()) cycle
 957
 958     bounds_out = bounds
 959     if (.not. using_blocking) then
 960       !< Set output bounds to start at 1:size(buffer_ptr%buffer)
 961       call bounds_out%reset_bounds_from_array_4D(buffer_ptr%buffer(:,:,:,:,1))
 962     endif
 963
 964     bounds_in = bounds
 965     if (.not. bounds%has_halos) then
 966       !< If field_data does not contain halos, set bounds_in to start at 1:size(field_data)
 967       call bounds_in%reset_bounds_from_array_4D(field_data)
 968     endif
 969
 970     is_subregional = file_ptr%is_regional()
 971     reduced_k_range = field_yaml_ptr%has_var_zbounds()
 972
 973     !< Reset the bounds based on the reduced k range and subregional
 974     is_subregional_reduced_k_range: if (is_subregional .or. reduced_k_range) then
 975       call buffer_ptr%get_axis_ids(axis_ids)
 976       block_in_subregion = .true.
 977       axis_loops: do i = 1, size(axis_ids)
 978         !< Move on if the block does not have any data for the subregion
 979         if (.not. block_in_subregion) cycle
 980
 981         select type (diag_axis => this%diag_axis(axis_ids(i))%axis)
 982         type is (fmsDiagSubAxis_type)
 983           sindex = diag_axis%get_starting_index()
 984           eindex = diag_axis%get_ending_index()
 985           compute_idx = diag_axis%get_compute_indices()
 986           starting=sindex-compute_idx(1)+1
 987           ending=eindex-compute_idx(1)+1
 988           if (using_blocking) then
 989             block_in_subregion = determine_if_block_is_in_region(starting, ending, bounds, i)
 990             if (.not. block_in_subregion) cycle
 991
 992             !< Set bounds_in so that you can the correct section of the data for the block (starting at 1)
 993             call bounds_in%rebase_input(bounds, starting, ending, i)
 994
 995             !< Set bounds_out to be the correct section relative to the block starting and ending indices
 996             call bounds_out%rebase_output(starting, ending, i)
 997           else
 998             !< Set bounds_in so that only the subregion section of the data will be used (starting at 1)
 999             call bounds_in%update_index(starting, ending, i, .false.)
1000
1001             !< Set bounds_out to 1:size(subregion) for the PE
1002             call bounds_out%update_index(1, ending-starting+1, i, .true.)
1003           endif
1004         end select
1005       enddo axis_loops
1006       nullify(axis_ids)
1007       !< Move on to the next buffer if the block does not have any data for the subregion
1008       if (.not. block_in_subregion) cycle
1009     endif is_subregional_reduced_k_range
1010
1011     !< Determine the reduction method for the buffer
1012     reduction_method = field_yaml_ptr%get_var_reduction()
1013     if (present(time)) call buffer_ptr%update_buffer_time(time)
1014     call buffer_ptr%set_send_data_called()
1015     select case(reduction_method)
1016     case (time_none)
1017       error_msg = buffer_ptr%do_time_none_wrapper(field_data, oor_mask, field_ptr%get_var_is_masked(), &
1018         bounds_in, bounds_out, missing_value)
1019       if (trim(error_msg) .ne. "") then
1020         return
1021       endif
1022     case (time_min)
1023       error_msg = buffer_ptr%do_time_min_wrapper(field_data, oor_mask, field_ptr%get_var_is_masked(), &
1024         bounds_in, bounds_out, missing_value)
1025       if (trim(error_msg) .ne. "") then
1026         return
1027       endif
1028     case (time_max)
1029       error_msg = buffer_ptr%do_time_max_wrapper(field_data, oor_mask, field_ptr%get_var_is_masked(), &
1030         bounds_in, bounds_out, missing_value)
1031       if (trim(error_msg) .ne. "") then
1032         return
1033       endif
1034     case (time_sum)
1035       error_msg = buffer_ptr%do_time_sum_wrapper(field_data, oor_mask, field_ptr%get_var_is_masked(), &
1036         field_ptr%get_mask_variant(), bounds_in, bounds_out, missing_value, field_ptr%has_missing_value())
1037       if (trim(error_msg) .ne. "") then
1038         return
1039       endif
1040     case (time_average)
1041       error_msg = buffer_ptr%do_time_sum_wrapper(field_data, oor_mask, field_ptr%get_var_is_masked(), &
1042         field_ptr%get_mask_variant(), bounds_in, bounds_out, missing_value, field_ptr%has_missing_value())
1043       if (trim(error_msg) .ne. "") then
1044         return
1045       endif
1046     case (time_power)
1047       error_msg = buffer_ptr%do_time_sum_wrapper(field_data, oor_mask, field_ptr%get_var_is_masked(), &
1048         field_ptr%get_mask_variant(), bounds_in, bounds_out, missing_value, field_ptr%has_missing_value(), &
1049         pow_value=field_yaml_ptr%get_pow_value())
1050       if (trim(error_msg) .ne. "") then
1051         return
1052       endif
1053     case (time_rms)
1054       error_msg = buffer_ptr%do_time_sum_wrapper(field_data, oor_mask, field_ptr%get_var_is_masked(), &
1055         field_ptr%get_mask_variant(), bounds_in, bounds_out, missing_value, field_ptr%has_missing_value(), &
1056         pow_value = 2)
1057       if (trim(error_msg) .ne. "") then
1058         return
1059       endif
1060     case (time_diurnal)
1061       if(.not. present(time)) call mpp_error(FATAL, &
1062                             "fms_diag_do_reduction:: time must be present when using diurnal reductions")
1063       ! sets the diurnal index for reduction within the buffer object
1064       call buffer_ptr%set_diurnal_section_index(time)
1065       error_msg = buffer_ptr%do_time_sum_wrapper(field_data, oor_mask, field_ptr%get_var_is_masked(), &
1066         field_ptr%get_mask_variant(), bounds_in, bounds_out, missing_value, field_ptr%has_missing_value())
1067       if (trim(error_msg) .ne. "") then
1068         return
1069       endif
1070     case default
1071       error_msg = "The reduction method is not supported. "//&
1072         "Only none, min, max, sum, average, power, rms, and diurnal are supported."
1073     end select
1074
1075     if (field_ptr%is_static() .or. file_ptr%FMS_diag_file%is_done_writing_data()) then
1076       call buffer_ptr%set_done_with_math()
1077     endif
1078   enddo buffer_loop
1079 #else
1080   error_msg = ""
1081   CALL MPP_ERROR(FATAL,"You can not use the modern diag manager without compiling with -Duse_yaml")
1082 #endif
1083 end function fms_diag_do_reduction
1084
1085 !> @brief Adds the diag ids of the Area and or Volume of the diag_field_object
1086 subroutine fms_diag_field_add_cell_measures(this, diag_field_id, area, volume)
1087   class(fmsDiagObject_type), intent (inout) :: this          !< The diag object
1088   integer,                   intent(in)     :: diag_field_id !< diag_field to add the are and volume to
1089   INTEGER, optional,         INTENT(in)     :: area          !< diag ids of area
1090   INTEGER, optional,         INTENT(in)     :: volume        !< diag ids of volume
1091
1092 #ifndef use_yaml
1093   CALL MPP_ERROR(FATAL,"You can not use the modern diag manager without compiling with -Duse_yaml")
1094 #else
1095   call this%FMS_diag_fields(diag_field_id)%add_area_volume(area, volume)
1096 #endif
1097 end subroutine fms_diag_field_add_cell_measures
1098
1099 !> @brief Add a attribute to the diag_obj using the diag_field_id
1100 subroutine fms_diag_field_add_attribute(this, diag_field_id, att_name, att_value)
1101   class(fmsDiagObject_type), intent (inout) :: this !< The diag object
1102   integer,          intent(in) :: diag_field_id      !< Id of the axis to add the attribute to
1103   character(len=*), intent(in) :: att_name     !< Name of the attribute
1104   class(*),         intent(in) :: att_value(:) !< The attribute value to add
1105 #ifndef use_yaml
1106 CALL MPP_ERROR(FATAL,"You can not use the modern diag manager without compiling with -Duse_yaml")
1107 #else
1108 !TODO: Value for diag not found
1109   if ( diag_field_id .LE. 0 ) THEN
1110     RETURN
1111   else
1112     if (this%FMS_diag_fields(diag_field_id)%is_registered() ) &
1113       call this%FMS_diag_fields(diag_field_id)%add_attribute(att_name, att_value)
1114   endif
1115 #endif
1116 end subroutine fms_diag_field_add_attribute
1117
1118 !> @brief Add an attribute to an axis
1119 subroutine fms_diag_axis_add_attribute(this, axis_id, att_name, att_value)
1120   class(fmsDiagObject_type), intent (inout) :: this !< The diag object
1121   integer,          intent(in) :: axis_id      !< Id of the axis to add the attribute to
1122   character(len=*), intent(in) :: att_name     !< Name of the attribute
1123   class(*),         intent(in) :: att_value(:) !< The attribute value to add
1124
1125   character(len=20) :: axis_names(2) !< Names of the uncompress axis
1126   character(len=20) :: set_name      !< Name of the axis set
1127   integer           :: uncmx_ids(2)  !< Ids of the uncompress axis
1128   integer           :: j             !< For do loops
1129 #ifndef use_yaml
1130 CALL MPP_ERROR(FATAL,"You can not use the modern diag manager without compiling with -Duse_yaml")
1131 #else
1132   if (axis_id < 0 .and. axis_id > this%registered_axis) &
1133     call mpp_error(FATAL, "diag_axis_add_attribute: The axis_id is not valid")
1134
1135   select type (axis => this%diag_axis(axis_id)%axis)
1136   type is (fmsDiagFullAxis_type)
1137     call axis%add_axis_attribute(att_name, att_value)
1138
1139     !! Axis that are in the "unstructured" domain require a "compress" attribute for the
1140     !! combiner and PP. This attribute is passed in via a diag_axis_add_attribute call in the model code
1141     !! The compress attribute indicates the names of the axis that were compressed
1142     !! For example grid_index:compress = "grid_yt grid_xt"
1143     !! The metadata and the data for these axis also needs to be written to the file
1144     if (trim(att_name) .eq. "compress") then
1145       !< If the attribute is the "compress" attribute, get the axis names,
1146       !! and the ids of the axis and add it to the axis object so it can be written to netcdf files
1147       !! that use this axis
1148       axis_names = parse_compress_att(att_value)
1149       set_name = ""
1150       if (axis%has_set_name()) set_name = axis%get_set_name()
1151       do j = 1, size(axis_names)
1152         uncmx_ids(j) = get_axis_id_from_name(axis_names(j), this%diag_axis, this%registered_axis, set_name)
1153         if (uncmx_ids(j) .eq. diag_null) call mpp_error(FATAL, &
1154           &"Error parsing the compress attribute for axis: "//trim(axis%get_axis_name())//&
1155           &". Be sure that the axes in the compress attribute are registered")
1156       enddo
1157       call axis%add_structured_axis_ids(uncmx_ids)
1158     endif
1159   end select
1160 #endif
1161 end subroutine fms_diag_axis_add_attribute
1162
1163 !> \brief Gets the field_name from the diag_field
1164 !> \returns a copy of the field_name
1165 function fms_get_field_name_from_id (this, field_id) &
1166   result(field_name)
1167
1168   class(fmsDiagObject_type), intent (in) :: this     !< The diag object, the caller
1169   integer,                   intent (in) :: field_id !< Field id to get the name for
1170   character(len=:), allocatable :: field_name
1171 #ifndef use_yaml
1172   CALL MPP_ERROR(FATAL,"You can not use the modern diag manager without compiling with -Duse_yaml")
1173 #else
1174   field_name = this%FMS_diag_fields(field_id)%get_varname()
1175 #endif
1176 end function fms_get_field_name_from_id
1177
1178 !> \brief Gets the diag field ID from the module name and field name.
1179 !> \returns a copy of the ID of the diag field or DIAG_FIELD_NOT_FOUND if the field is not registered
1180 FUNCTION fms_get_diag_field_id_from_name(this, module_name, field_name) &
1181   result(diag_field_id)
1182   class(fmsDiagObject_type), intent (in) :: this !< The diag object, the caller
1183   CHARACTER(len=*), INTENT(in) :: module_name !< Module name that registered the variable
1184   CHARACTER(len=*), INTENT(in) :: field_name !< Variable name
1185   integer :: diag_field_id
1186
1187 #ifdef use_yaml
1188   integer              :: i                     !< For looping
1189   integer, allocatable :: diag_field_indices(:) !< indices where the field was found in the yaml
1190
1191   diag_field_id = DIAG_FIELD_NOT_FOUND
1192
1193   !> Loop through fields to find it.
1194   do i=1, this%registered_variables
1195     !< Check if the field was registered, if it was return the diag_field_id
1196     diag_field_id = this%FMS_diag_fields(i)%id_from_name(module_name, field_name)
1197     if(diag_field_id .ne. DIAG_FIELD_NOT_FOUND) return
1198   enddo
1199
1200   !< Check if the field is in the diag_table.yaml. If it is, return DIAG_FIELD_NOT_REGISTERED
1201   !! Otherwsie it will return DIAG_FIELD_NOT_FOUND
1202   diag_field_indices = find_diag_field(field_name, module_name)
1203   if (diag_field_indices(1) .ne. diag_null) then
1204     diag_field_id = DIAG_NOT_REGISTERED
1205   endif
1206   deallocate(diag_field_indices)
1207 #else
1208   diag_field_id = DIAG_FIELD_NOT_FOUND
1209   CALL MPP_ERROR(FATAL,"You can not use the modern diag manager without compiling with -Duse_yaml")
1210 #endif
1211 END FUNCTION fms_get_diag_field_id_from_name
1212
1213 #ifdef use_yaml
1214 !> returns the buffer object for the given id
1215 !! actual data comes from %get_buffer_data() on the returned object
1216 function get_diag_buffer(this, bufferid) &
1217 result(rslt)
1218   class(fmsDiagObject_type), intent(in) :: this
1219   integer, intent(in)                   :: bufferid
1220   type(fmsDiagOutputBuffer_type),allocatable:: rslt
1221   if( (bufferid .gt. UBOUND(this%FMS_diag_output_buffers, 1)) .or. &
1222       (bufferid .lt. LBOUND(this%FMS_diag_output_buffers, 1))) &
1223     call mpp_error(FATAL, 'get_diag_bufer: invalid bufferid given')
1224   rslt = this%FMS_diag_output_buffers(bufferid)
1225 end function
1226 #endif
1227
1228 !> @brief Return the 2D domain for the axis IDs given.
1229 !! @return 2D domain for the axis IDs given
1230 type(domain2d) FUNCTION fms_get_domain2d(this, ids)
1231   class(fmsDiagObject_type), intent (in) :: this !< The diag object
1232   INTEGER, DIMENSION(:), INTENT(in) :: ids !< Axis IDs.
1233
1234 #ifndef use_yaml
1235 CALL MPP_ERROR(FATAL,"You can not use the modern diag manager without compiling with -Duse_yaml")
1236 fms_get_domain2d = null_domain2d
1237 #else
1238   INTEGER                      :: type_of_domain !< The type of domain
1239   CLASS(diagDomain_t), POINTER :: domain         !< Diag Domain pointer
1240
1241   call get_domain_and_domain_type(fms_diag_object%diag_axis, ids, type_of_domain, domain, "get_domain2d")
1242   if (type_of_domain .ne. TWO_D_DOMAIN) &
1243     call mpp_error(FATAL, 'diag_axis_mod::get_domain2d- The axis do not correspond to a 2d Domain')
1244   select type(domain)
1245   type is (diagDomain2d_t)
1246     fms_get_domain2d = domain%domain2
1247   end select
1248 #endif
1249 END FUNCTION fms_get_domain2d
1250
1251  !> @brief Gets the length of the axis based on the axis_id
1252  !> @return Axis_length
1253  integer function fms_get_axis_length(this, axis_id)
1254   class(fmsDiagObject_type), intent (in) :: this !< The diag object
1255   INTEGER, INTENT(in) :: axis_id !< Axis ID of the axis to the length of
1256
1257 #ifndef use_yaml
1258 CALL MPP_ERROR(FATAL,"You can not use the modern diag manager without compiling with -Duse_yaml")
1259 fms_get_axis_length = 0
1260 #else
1261 fms_get_axis_length = 0
1262
1263   if (axis_id < 0 .and. axis_id > this%registered_axis) &
1264     call mpp_error(FATAL, "fms_get_axis_length: The axis_id is not valid")
1265
1266   select type (axis => this%diag_axis(axis_id)%axis)
1267   type is (fmsDiagFullAxis_type)
1268     fms_get_axis_length = axis%axis_length()
1269   type is (fmsDiagSubAxis_type)
1270     fms_get_axis_length = axis%axis_length()
1271   end select
1272 #endif
1273 end function fms_get_axis_length
1274
1275 !> @brief Gets the name of the axis based on the axis_id
1276  !> @return The axis_name
1277 function fms_get_axis_name_from_id (this, axis_id) &
1278 result(axis_name)
1279   class(fmsDiagObject_type), intent (in) :: this !< The diag object
1280   INTEGER, INTENT(in) :: axis_id !< Axis ID of the axis to the length of
1281
1282   character (len=:), allocatable :: axis_name
1283
1284 #ifndef use_yaml
1285 CALL MPP_ERROR(FATAL,"You can not use the modern diag manager without compiling with -Duse_yaml")
1286 axis_name=" "
1287 #else
1288     if (axis_id < 0 .and. axis_id > this%registered_axis) &
1289       call mpp_error(FATAL, "fms_get_axis_length: The axis_id is not valid")
1290
1291     !! if its a scalar (null axis id) just returns n/a since no axis is defined
1292     if (axis_id .eq. NULL_AXIS_ID) then
1293       allocate(character(len=3) :: axis_name)
1294       axis_name = "n/a"
1295       return
1296     endif
1297
1298
1299     select type (axis => this%diag_axis(axis_id)%axis)
1300     type is (fmsDiagFullAxis_type)
1301       axis_name = axis%get_axis_name()
1302     end select
1303 #endif
1304 end function fms_get_axis_name_from_id
1305
1306 !> Dumps as much data as it can from the fmsDiagObject_type.
1307 !! Will dump any fields and files as well (see d)
1308 subroutine dump_diag_obj( filename )
1309   character(len=*), intent(in), optional :: filename !< optional filename to print to,
1310                                             !! otherwise prints to stdout
1311 #ifdef use_yaml
1312   !type(fmsDiagObject_type) :: diag_obj
1313   type(fmsDiagFile_type), pointer :: fileptr !<  pointer for traversing file list
1314   type(fmsDiagField_type), pointer :: fieldptr !<  pointer for traversing field list
1315   integer :: i !< do loops
1316   integer :: unit_num !< unit num of opened log file or stdout
1317
1318   if( present(filename) ) then
1319     open(newunit=unit_num, file=trim(filename), action='WRITE')
1320   else
1321     unit_num = stdout()
1322   endif
1323   if( mpp_pe() .eq. mpp_root_pe()) then
1324     write(unit_num, *) '********** dumping diag object ***********'
1325     write(unit_num, *) 'registered_variables:', fms_diag_object%registered_variables
1326     write(unit_num, *) 'registered_axis:', fms_diag_object%registered_axis
1327     write(unit_num, *) 'initialized:', fms_diag_object%initialized
1328     write(unit_num, *) 'files_initialized:', fms_diag_object%files_initialized
1329     write(unit_num, *) 'fields_initialized:', fms_diag_object%fields_initialized
1330     write(unit_num, *) 'buffers_initialized:', fms_diag_object%buffers_initialized
1331     write(unit_num, *) 'axes_initialized:', fms_diag_object%axes_initialized
1332     write(unit_num, *) 'Files:'
1333     if( fms_diag_object%files_initialized ) then
1334       do i=1, SIZE(fms_diag_object%FMS_diag_files)
1335         write(unit_num, *) 'File num:', i
1336         fileptr => fms_diag_object%FMS_diag_files(i)%FMS_diag_file
1337         call fileptr%dump_file_obj(unit_num)
1338       enddo
1339     else
1340       write(unit_num, *) 'files not initialized'
1341     endif
1342     if( fms_diag_object%fields_initialized) then
1343       do i=1, SIZE(fms_diag_object%FMS_diag_fields)
1344         write(unit_num, *) 'Field num:', i
1345         fieldptr => fms_diag_object%FMS_diag_fields(i)
1346         call fieldptr%dump_field_obj(unit_num)
1347       enddo
1348     else
1349       write(unit_num, *) 'fields not initialized'
1350     endif
1351     if( present(filename) ) close(unit_num)
1352   endif
1353 #else
1354   call mpp_error( FATAL, "You can not use the modern diag manager without compiling with -Duse_yaml")
1355 #endif
1356 end subroutine
1357
1358 !> @brief Allocates the output buffers of the fields corresponding to the registered variable
1359 !! Input arguments are the field and its ID passed to routine fms_diag_accept_data()
1360 subroutine allocate_diag_field_output_buffers(this, field_data, field_id)
1361   class(fmsDiagObject_type), target, intent(inout) :: this !< diag object
1362   class(*), dimension(:,:,:,:), intent(in) :: field_data !< field data
1363   integer, intent(in) :: field_id !< Id of the field data
1364 #ifdef use_yaml
1365   integer :: ndims !< Number of dimensions in the input field data
1366   integer :: buffer_id !< Buffer index of FMS_diag_buffers
1367   integer :: num_diurnal_samples !< Number of diurnal samples from diag_yaml
1368   integer :: axes_length(4) !< Length of each axis
1369   integer :: i, j !< For looping
1370   class(fmsDiagOutputBuffer_type), pointer :: ptr_diag_buffer_obj !< Pointer to the buffer class
1371   class(DiagYamlFilesVar_type), pointer :: ptr_diag_field_yaml !< Pointer to a field from yaml fields
1372   integer, pointer :: axis_ids(:) !< Pointer to indices of axes of the field variable
1373   integer :: var_type !< Stores type of the field data (r4, r8, i4, i8, and string) represented as an integer.
1374   character(len=:), allocatable :: var_name !< Field name to initialize output buffers
1375   logical :: is_scalar !< Flag indicating that the variable is a scalar
1376   integer :: yaml_id !< Yaml id for the buffer
1377   integer :: file_id !< File id for the buffer
1378
1379   if (this%FMS_diag_fields(field_id)%buffer_allocated) return
1380
1381   ! Determine the type of the field data
1382   var_type = get_var_type(field_data(1, 1, 1, 1))
1383
1384   ! Get variable/field name
1385   var_name = this%FMS_diag_fields(field_id)%get_varname()
1386
1387   ! Determine dimensions of the field
1388   is_scalar = this%FMS_diag_fields(field_id)%is_scalar()
1389
1390   ! Loop over a number of fields/buffers where this variable occurs
1391   do i = 1, size(this%FMS_diag_fields(field_id)%buffer_ids)
1392     buffer_id = this%FMS_diag_fields(field_id)%buffer_ids(i)
1393     file_id = this%FMS_diag_fields(field_id)%file_ids(i)
1394
1395     !< Go away if the file is a subregional file and the current PE does not have any data for it
1396     if (.not. this%FMS_diag_files(file_id)%writing_on_this_pe()) cycle
1397
1398     ndims = 0
1399     if (.not. is_scalar) then
1400       call this%FMS_diag_output_buffers(buffer_id)%get_axis_ids(axis_ids)
1401       ndims = size(axis_ids)
1402     endif
1403
1404     yaml_id = this%FMS_diag_output_buffers(buffer_id)%get_yaml_id()
1405
1406     ptr_diag_field_yaml => diag_yaml%diag_fields(yaml_id)
1407     num_diurnal_samples = ptr_diag_field_yaml%get_n_diurnal() !< Get number of diurnal samples
1408
1409     axes_length = 1
1410     do j = 1, ndims
1411       axes_length(j) = this%fms_get_axis_length(axis_ids(j))
1412     enddo
1413
1414     if (num_diurnal_samples .ne. 0) then
1415       ndims = ndims + 1 !< Add one more dimension for the diurnal axis
1416     endif
1417
1418     ptr_diag_buffer_obj => this%FMS_diag_output_buffers(buffer_id)
1419     call ptr_diag_buffer_obj%allocate_buffer(field_data(1, 1, 1, 1), ndims, axes_length(1:4), &
1420       this%FMS_diag_fields(field_id)%get_mask_variant(), var_name, num_diurnal_samples)
1421     call ptr_diag_buffer_obj%initialize_buffer(ptr_diag_field_yaml%get_var_reduction(), var_name)
1422
1423   enddo
1424   nullify(axis_ids)
1425
1426   this%FMS_diag_fields(field_id)%buffer_allocated = .true.
1427 #else
1428   call mpp_error( FATAL, "allocate_diag_field_output_buffers: "//&
1429     "you can not use the modern diag manager without compiling with -Duse_yaml")
1430 #endif
1431 end subroutine allocate_diag_field_output_buffers
1432
1433 !> @brief Determines if the window defined by the input bounds is a physics window.
1434 !> @return TRUE if the window size is less then the actual field size else FALSE.
1435 function fms_diag_compare_window(this, field, field_id, &
1436   is_in, ie_in, js_in, je_in, ks_in, ke_in) result(is_phys_win)
1437   class(fmsDiagObject_type), intent(in) :: this !< Diag Object
1438   class(*), intent(in) :: field(:,:,:,:) !< Field data
1439   integer, intent(in) :: field_id !< ID of the input field
1440   integer, intent(in) :: is_in, js_in !< Starting field indices for the first 2 dimensions;
1441                                       !< pass reconditioned indices fis and fjs
1442                                       !< which are computed elsewhere.
1443   integer, intent(in) :: ie_in, je_in !< Ending field indices for the first 2 dimensions;
1444                                       !< pass reconditioned indices fie and fje
1445                                       !< which are computed elsewhere.
1446   integer, intent(in) :: ks_in, ke_in !< Starting and ending indices of the field in 3rd dimension
1447   logical :: is_phys_win !< Return flag
1448 #ifdef use_yaml
1449   integer, pointer :: axis_ids(:)
1450   integer :: total_elements
1451   integer :: i !< For do loop
1452   integer :: field_size
1453   integer, allocatable :: field_shape(:) !< Shape of the field data
1454   integer :: window_size
1455
1456   !> Determine shape of the field defined by the input bounds
1457   field_shape = shape(field(is_in:ie_in, js_in:je_in, ks_in:ke_in, :))
1458
1459   window_size = field_shape(1) * field_shape(2) * field_shape(3)
1460
1461   total_elements = 1
1462   axis_ids => this%FMS_diag_fields(field_id)%get_axis_id()
1463   do i=1, size(axis_ids)
1464     total_elements = total_elements * this%fms_get_axis_length(axis_ids(i))
1465   enddo
1466
1467   if (total_elements > window_size) then
1468     is_phys_win = .true.
1469   else
1470     is_phys_win = .false.
1471   end if
1472 #else
1473   is_phys_win = .false.
1474   call mpp_error( FATAL, "fms_diag_compare_window: "//&
1475     "you can not use the modern diag manager without compiling with -Duse_yaml")
1476 #endif
1477 end function fms_diag_compare_window
1478
1479 !> @brief Update the current model time in the diag object
1480 subroutine update_current_model_time(this, time)
1481   class(fmsDiagObject_type), intent(inout) :: this !< Diag Object
1482   type(time_type),           intent(in)    :: time !< Current diag manager time
1483 #ifdef use_yaml
1484   if(time > this%current_model_time) this%current_model_time = time
1485 #endif
1486 end subroutine update_current_model_time
1487
1488 end module fms_diag_object_mod