diff --git a/diag_manager/diag_data.F90 b/diag_manager/diag_data.F90 index fcf519ff90..6d1ea05694 100644 --- a/diag_manager/diag_data.F90 +++ b/diag_manager/diag_data.F90 @@ -335,6 +335,7 @@ MODULE diag_data_mod character(len=:), allocatable :: att_name !< Name of the attribute contains procedure :: add => fms_add_attribute + procedure :: write_metadata end type fmsDiagAttribute_type ! Include variable "version" to be written to log file. #include @@ -562,8 +563,9 @@ function get_base_second() & res = base_second end function get_base_second + !> @brief Adds an attribute to the attribute type subroutine fms_add_attribute(this, att_name, att_value) - class(fmsDiagAttribute_type), intent(inout) :: this !< Diag attribute type + class(fmsDiagAttribute_type), intent(inout) :: this !< Diag attribute type character(len=*), intent(in) :: att_name !< Name of the attribute class(*), intent(in) :: att_value(:) !< The attribute value to add @@ -589,6 +591,38 @@ subroutine fms_add_attribute(this, att_name, att_value) this%att_value = att_value end select end subroutine fms_add_attribute + + !> @brief Writes out the attributes from an fmsDiagAttribute_type + subroutine write_metadata(this, fileobj, var_name, cell_methods) + class(fmsDiagAttribute_type), intent(inout) :: this !< Diag attribute type + class(FmsNetcdfFile_t), INTENT(INOUT) :: fileobj !< Fms2_io fileobj to write to + character(len=*), intent(in) :: var_name !< The name of the variable to write to + character(len=*), optional, intent(inout) :: cell_methods !< The cell methods attribute + + select type (att_value =>this%att_value) + type is (character(len=*)) + !< If the attribute is cell methods append to the current cell_methods attribute value + !! This will be writen once all of the cell_methods attributes are gathered ... + if (present(cell_methods)) then + if (trim(this%att_name) .eq. "cell_methods") then + cell_methods = trim(cell_methods)//" "//trim(att_value(1)) + return + endif + endif + + call register_variable_attribute(fileobj, var_name, this%att_name, trim(att_value(1)), & + str_len=len_trim(att_value(1))) + type is (real(kind=r8_kind)) + call register_variable_attribute(fileobj, var_name, this%att_name, real(att_value, kind=r8_kind)) + type is (real(kind=r4_kind)) + call register_variable_attribute(fileobj, var_name, this%att_name, real(att_value, kind=r4_kind)) + type is (integer(kind=i4_kind)) + call register_variable_attribute(fileobj, var_name, this%att_name, int(att_value, kind=i4_kind)) + type is (integer(kind=i8_kind)) + call register_variable_attribute(fileobj, var_name, this%att_name, int(att_value, kind=i8_kind)) + end select + + end subroutine write_metadata END MODULE diag_data_mod !> @} ! close documentation grouping diff --git a/diag_manager/diag_manager.F90 b/diag_manager/diag_manager.F90 index fe944297d8..38a32f6663 100644 --- a/diag_manager/diag_manager.F90 +++ b/diag_manager/diag_manager.F90 @@ -4432,6 +4432,11 @@ SUBROUTINE diag_field_add_cell_measures(diag_field_id, area, volume) & 'either area or volume arguments must be present', FATAL ) END IF + if (use_modern_diag) then + call fms_diag_object%fms_diag_field_add_cell_measures(diag_field_id, area, volume) + return + ENDIF + DO j=1, input_fields(diag_field_id)%num_output_fields ind = input_fields(diag_field_id)%output_fields(j) CALL init_field_cell_measures(output_fields(ind), area=area, volume=volume) diff --git a/diag_manager/fms_diag_axis_object.F90 b/diag_manager/fms_diag_axis_object.F90 index f441287e66..343388ae0d 100644 --- a/diag_manager/fms_diag_axis_object.F90 +++ b/diag_manager/fms_diag_axis_object.F90 @@ -103,6 +103,7 @@ module fms_diag_axis_object_mod procedure :: add_structured_axis_ids procedure :: get_structured_axis procedure :: is_unstructured_grid + procedure :: get_edges_id END TYPE fmsDiagAxis_type !> @brief Type to hold the subaxis @@ -152,6 +153,8 @@ module fms_diag_axis_object_mod !! or "UG_DOMAIN") INTEGER , private :: length !< Global axis length INTEGER , private :: direction !< Direction of the axis 0, 1, -1 + INTEGER, ALLOCATABLE, private :: edges_id !< Axis ID for the edges axis + !! This axis will be written to the file CHARACTER(len=:), ALLOCATABLE, private :: edges_name !< Name for the previously defined "edges axis" !! This will be written as an attribute CHARACTER(len=:), ALLOCATABLE, private :: aux !< Auxiliary name, can only be geolon_t @@ -169,7 +172,7 @@ module fms_diag_axis_object_mod PROCEDURE :: add_axis_attribute PROCEDURE :: register => register_diag_axis_obj PROCEDURE :: axis_length => get_axis_length - PROCEDURE :: set_edges_name + PROCEDURE :: set_edges PROCEDURE :: set_axis_id PROCEDURE :: get_compute_domain PROCEDURE :: get_indices @@ -188,7 +191,7 @@ module fms_diag_axis_object_mod !> @brief Initialize the axis subroutine register_diag_axis_obj(this, axis_name, axis_data, units, cart_name, long_name, direction,& & set_name, Domain, Domain2, DomainU, aux, req, tile_count, domain_position, axis_length ) - class(fmsDiagFullAxis_type),INTENT(out) :: this !< Diag_axis obj + class(fmsDiagFullAxis_type),INTENT(inout):: this !< Diag_axis obj CHARACTER(len=*), INTENT(in) :: axis_name !< Name of the axis class(*), INTENT(in) :: axis_data(:) !< Array of coordinate values CHARACTER(len=*), INTENT(in) :: units !< Units for the axis @@ -299,6 +302,9 @@ subroutine write_axis_metadata(this, fileobj, parent_axis) type(fmsDiagFullAxis_type), pointer :: diag_axis !< Local pointer to the diag_axis integer :: type_of_domain !< The type of domain the current axis is in + logical :: is_subaxis !< .true. if the axis is a subaxis + + is_subaxis = .false. select type(this) type is (fmsDiagFullAxis_type) @@ -307,6 +313,7 @@ subroutine write_axis_metadata(this, fileobj, parent_axis) diag_axis => this type_of_domain = this%type_of_domain type is (fmsDiagSubAxis_type) + is_subaxis = .true. axis_name => this%subaxis_name axis_length = this%ending_index - this%starting_index + 1 !< Get all the other information from the parent axis (i.e the cart_name, units, etc) @@ -373,7 +380,7 @@ subroutine write_axis_metadata(this, fileobj, parent_axis) call register_variable_attribute(fileobj, axis_name, "positive", "down", str_len=4) end select - if (allocated(diag_axis%edges_name)) then + if (allocated(diag_axis%edges_name) .and. .not. is_subaxis) then call register_variable_attribute(fileobj, axis_name, "edges", diag_axis%edges_name, & str_len=len_trim(diag_axis%edges_name)) endif @@ -521,6 +528,19 @@ pure function get_structured_axis(this) & end select end function get_structured_axis + + !< @brief Get the edges_id of an axis_object + !! @return The edges_id of an axis object + pure integer function get_edges_id(this) + class(fmsDiagAxis_type), INTENT(in) :: this !< diag_axis obj + + get_edges_id = diag_null + select type (this) + type is (fmsDiagFullAxis_type) + if (allocated(this%edges_id)) get_edges_id = this%edges_id + end select + end function + !> @brief Get the starting and ending indices of the global io domain of the axis subroutine get_global_io_domain(this, global_io_index) class(fmsDiagFullAxis_type), intent(in) :: this !< diag_axis obj @@ -593,13 +613,19 @@ subroutine set_axis_id(this, axis_id) end subroutine set_axis_id - !> @brief Set the name of the edges - subroutine set_edges_name(this, edges_name) - class(fmsDiagFullAxis_type), intent(inout) :: this !< diag_axis obj - CHARACTER(len=*), intent(in) :: edges_name !< Name of the edges + !> @brief Set the name and ids of the edges + subroutine set_edges(this, edges_name, edges_id) + class(fmsDiagFullAxis_type), intent(inout) :: this !< diag_axis obj + CHARACTER(len=*), intent(in) :: edges_name !< Name of the edges + integer, intent(in) :: edges_id !< Axis id of the edges + !< Saving the name and the id of the edges axis because it will make it easier to use + !! downstream (i.e you need the edges name to write the attribute to the current axis, + !! and you need the edges id to add to the diag file object so that you can write the edges + !! to the file) this%edges_name = edges_name - end subroutine + this%edges_id = edges_id + end subroutine set_edges !> @brief Determine if the subRegion is in the current PE. !! If it is, determine the starting and ending indices of the current PE that belong to the subRegion diff --git a/diag_manager/fms_diag_field_object.F90 b/diag_manager/fms_diag_field_object.F90 index 095ca941c3..5ceadc1afa 100644 --- a/diag_manager/fms_diag_field_object.F90 +++ b/diag_manager/fms_diag_field_object.F90 @@ -147,6 +147,8 @@ module fms_diag_field_object_mod procedure :: write_field_metadata procedure :: write_coordinate_attribute procedure :: get_math_needs_to_be_done + procedure :: add_area_volume + procedure :: append_time_cell_methods end type fmsDiagField_type !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! variables !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! type(fmsDiagField_type) :: null_ob @@ -671,11 +673,21 @@ end function get_vartype !> @brief Gets varname !! @return copy of the variable name -pure function get_varname (this) & +pure function get_varname (this, to_write) & result(rslt) - class (fmsDiagField_type), intent(in) :: this !< diag object - character(len=:), allocatable :: rslt - rslt = this%varname + class (fmsDiagField_type), intent(in) :: this !< diag object + logical, optional, intent(in) :: to_write !< .true. if getting the varname that will be writen to the file + character(len=:), allocatable :: rslt + rslt = this%varname + + !< If writing the varname can be the outname which is defined in the yaml + if (present(to_write)) then + if (to_write) then + !TODO this is wrong + rslt = this%diag_field(1)%get_var_outname() + endif + endif + end function get_varname !> @brief Gets longname @@ -1082,13 +1094,15 @@ subroutine write_field_metadata(this, fileobj, file_id, yaml_id, diag_axis, unli class(fmsDiagAxisContainer_type), intent(in) :: diag_axis(:) !< Diag_axis object character(len=*), intent(in) :: unlim_dimname !< The name of the unlimited dimension logical, intent(in) :: is_regional !< Flag indicating if the field is regional - character(len=*), intent(inout) :: cell_measures + character(len=*), intent(in) :: cell_measures !< The cell measures attribute to write type(diagYamlFilesVar_type), pointer :: field_yaml !< pointer to the yaml entry character(len=:), allocatable :: var_name !< Variable name character(len=:), allocatable :: long_name !< Longname to write character(len=:), allocatable :: units !< Units of the field to write character(len=120), allocatable :: dimnames(:) !< Dimension names of the field + character(len=120) :: cell_methods!< Cell methods attribute to write + integer :: i !< For do loops field_yaml => diag_yaml%get_diag_field_from_id(yaml_id) var_name = field_yaml%get_var_outname() @@ -1130,16 +1144,34 @@ subroutine write_field_metadata(this, fileobj, file_id, yaml_id, diag_axis, unli str_len=len_trim(this%get_interp_method())) endif - select case (field_yaml%get_var_reduction()) - case (time_average, time_max, time_min) - call register_variable_attribute(fileobj, var_name, "time_avg_info", & - trim(avg_name)//'_T1,'//trim(avg_name)//'_T2,'//trim(avg_name)//'_DT', & - str_len=len(trim(avg_name)//'_T1,'//trim(avg_name)//'_T2,'//trim(avg_name)//'_DT')) - end select + if (.not. this%static) then + select case (field_yaml%get_var_reduction()) + case (time_average, time_max, time_min, time_diurnal, time_power, time_rms, time_sum) + call register_variable_attribute(fileobj, var_name, "time_avg_info", & + trim(avg_name)//'_T1,'//trim(avg_name)//'_T2,'//trim(avg_name)//'_DT', & + str_len=len(trim(avg_name)//'_T1,'//trim(avg_name)//'_T2,'//trim(avg_name)//'_DT')) + end select + endif - call append_time_cell_measure(cell_measures, field_yaml) - if (trim(cell_measures) .ne. "") & + cell_methods = "" + !< Check if any of the attributes defined via a "diag_field_add_attribute" call + !! are the cell_methods, if so add to the "cell_methods" variable: + do i = 1, this%num_attributes + call this%attributes(i)%write_metadata(fileobj, var_name, & + cell_methods=cell_methods) + enddo + + !< Append the time cell methods based on the variable's reduction + call this%append_time_cell_methods(cell_methods, field_yaml) + if (trim(cell_methods) .ne. "") & call register_variable_attribute(fileobj, var_name, "cell_methods", & + trim(adjustl(cell_methods)), str_len=len_trim(adjustl(cell_methods))) + + !< Write out the cell_measures attribute (i.e Area, Volume) + !! The diag field ids for the Area and Volume are sent in the register call + !! This was defined in file object and passed in here + if (trim(cell_measures) .ne. "") & + call register_variable_attribute(fileobj, var_name, "cell_measures", & trim(adjustl(cell_measures)), str_len=len_trim(adjustl(cell_measures))) !< Write out the standard_name (this was defined in the register call) @@ -1432,30 +1464,65 @@ PURE FUNCTION diag_field_id_from_name(this, module_name, field_name) & endif end function diag_field_id_from_name -!> @brief Append the time cell measured based on the variable's reduction -subroutine append_time_cell_measure(cell_measures, field_yaml) - character(len=*), intent(inout) :: cell_measures !< The cell measures to append to - type(diagYamlFilesVar_type), intent(in) :: field_yaml !< The field's yaml +!> @brief Adds the area and volume id to a field object +subroutine add_area_volume(this, area, volume) + CLASS(fmsDiagField_type), intent(inout) :: this !< The field object + INTEGER, optional, INTENT(in) :: area !< diag ids of area + INTEGER, optional, INTENT(in) :: volume !< diag ids of volume + + if (present(area)) then + if (area > 0) then + this%area = area + else + call mpp_error(FATAL, "diag_field_add_cell_measures: the area id is not valid. "& + &"Verify that the area_id passed in to the field:"//this%varname//& + &" is valid and that the field is registered and in the diag_table.yaml") + endif + endif + + if (present(volume)) then + if (volume > 0) then + this%volume = volume + else + call mpp_error(FATAL, "diag_field_add_cell_measures: the volume id is not valid. "& + &"Verify that the volume_id passed in to the field:"//this%varname//& + &" is valid and that the field is registered and in the diag_table.yaml") + endif + endif + +end subroutine add_area_volume + +!> @brief Append the time cell meathods based on the variable's reduction +subroutine append_time_cell_methods(this, cell_methods, field_yaml) + class (fmsDiagField_type), target, intent(inout) :: this !< diag field + character(len=*), intent(inout) :: cell_methods !< The cell methods var to append to + type(diagYamlFilesVar_type), intent(in) :: field_yaml !< The field's yaml + + if (this%static) then + cell_methods = trim(cell_methods)//" time: point " + return + endif select case (field_yaml%get_var_reduction()) case (time_none) - cell_measures = trim(cell_measures)//" time: point " + cell_methods = trim(cell_methods)//" time: point " case (time_diurnal) - cell_measures = trim(cell_measures)//" time: mean" + cell_methods = trim(cell_methods)//" time: mean" case (time_power) - cell_measures = trim(cell_measures)//" time: mean_pow"//int2str(field_yaml%get_pow_value()) + cell_methods = trim(cell_methods)//" time: mean_pow"//int2str(field_yaml%get_pow_value()) case (time_rms) - cell_measures = trim(cell_measures)//" time: root_mean_square" + cell_methods = trim(cell_methods)//" time: root_mean_square" case (time_max) - cell_measures = trim(cell_measures)//" time: max" + cell_methods = trim(cell_methods)//" time: max" case (time_min) - cell_measures = trim(cell_measures)//" time: min" + cell_methods = trim(cell_methods)//" time: min" case (time_average) - cell_measures = trim(cell_measures)//" time: mean" + cell_methods = trim(cell_methods)//" time: mean" case (time_sum) - cell_measures = trim(cell_measures)//" time: sum" + cell_methods = trim(cell_methods)//" time: sum" end select -end subroutine append_time_cell_measure +end subroutine append_time_cell_methods + !> Dumps any data from a given fmsDiagField_type object subroutine dump_field_obj (this, unit_num) class(fmsDiagField_type), intent(in) :: this diff --git a/diag_manager/fms_diag_file_object.F90 b/diag_manager/fms_diag_file_object.F90 index 12e453fb0f..ef286cb70c 100644 --- a/diag_manager/fms_diag_file_object.F90 +++ b/diag_manager/fms_diag_file_object.F90 @@ -1213,6 +1213,7 @@ subroutine write_axis_metadata(this, diag_axis) integer :: i,k !< For do loops integer :: parent_axis_id !< Id of the parent_axis integer :: structured_ids(2) !< Ids of the uncompress axis + integer :: edges_id !< Id of the axis edge class(fmsDiagAxisContainer_type), pointer :: axis_ptr !< pointer to the axis object currently writing @@ -1234,6 +1235,12 @@ subroutine write_axis_metadata(this, diag_axis) call diag_axis(structured_ids(k))%axis%write_axis_metadata(fileobj) enddo endif + + edges_id = axis_ptr%axis%get_edges_id() + if (edges_id .ne. diag_null) then + if (any(diag_file%axis_ids(1:diag_file%number_of_axis) .eq. edges_id)) return + call diag_axis(edges_id)%axis%write_axis_metadata(fileobj) + endif enddo end subroutine write_axis_metadata @@ -1265,11 +1272,11 @@ subroutine write_field_metadata(this, diag_field, diag_axis) !the file that the fields are in needs to be added cell_measures = "" if (field_ptr%has_area()) then - cell_measures = "area:"//diag_field(field_ptr%get_area())%get_varname() + cell_measures = "area: "//diag_field(field_ptr%get_area())%get_varname(to_write=.true.) endif if (field_ptr%has_volume()) then - cell_measures = trim(cell_measures)//" volume:"//diag_field(field_ptr%get_volume())%get_varname() + cell_measures = trim(cell_measures)//" volume: "//diag_field(field_ptr%get_volume())%get_varname(to_write=.true.) endif call field_ptr%write_field_metadata(fileobj, diag_file%id, diag_file%yaml_ids(i), diag_axis, & diff --git a/diag_manager/fms_diag_object.F90 b/diag_manager/fms_diag_object.F90 index b5ac522847..f895529230 100644 --- a/diag_manager/fms_diag_object.F90 +++ b/diag_manager/fms_diag_object.F90 @@ -78,6 +78,7 @@ module fms_diag_object_mod procedure :: fms_diag_accept_data procedure :: fms_diag_send_complete procedure :: fms_diag_do_io + procedure :: fms_diag_field_add_cell_measures #ifdef use_yaml procedure :: get_diag_buffer #endif @@ -425,7 +426,7 @@ FUNCTION fms_diag_axis_init(this, axis_name, axis_data, units, cart_name, axis_l select type (edges_axis => this%diag_axis(edges)%axis) type is (fmsDiagFullAxis_type) edges_name = edges_axis%get_axis_name() - call axis%set_edges_name(edges_name) + call axis%set_edges(edges_name, edges) end select endif call axis%register(axis_name, axis_data, units, cart_name, long_name=long_name, & @@ -617,6 +618,20 @@ subroutine fms_diag_do_io(this, is_end_of_run) #endif end subroutine fms_diag_do_io +!> @brief Adds the diag ids of the Area and or Volume of the diag_field_object +subroutine fms_diag_field_add_cell_measures(this, diag_field_id, area, volume) + class(fmsDiagObject_type), intent (inout) :: this !< The diag object + integer, intent(in) :: diag_field_id !< diag_field to add the are and volume to + INTEGER, optional, INTENT(in) :: area !< diag ids of area + INTEGER, optional, INTENT(in) :: volume !< diag ids of volume + +#ifndef use_yaml + CALL MPP_ERROR(FATAL,"You can not use the modern diag manager without compiling with -Duse_yaml") +#else + call this%FMS_diag_fields(diag_field_id)%add_area_volume(area, volume) +#endif +end subroutine fms_diag_field_add_cell_measures + !> @brief Add a attribute to the diag_obj using the diag_field_id subroutine fms_diag_field_add_attribute(this, diag_field_id, att_name, att_value) class(fmsDiagObject_type), intent (inout) :: this !< The diag object diff --git a/test_fms/diag_manager/test_modern_diag.F90 b/test_fms/diag_manager/test_modern_diag.F90 index 67000e1ac2..8205b8eee1 100644 --- a/test_fms/diag_manager/test_modern_diag.F90 +++ b/test_fms/diag_manager/test_modern_diag.F90 @@ -26,7 +26,8 @@ program test_modern_diag mpp_get_UG_compute_domain use diag_manager_mod, only: diag_manager_init, diag_manager_end, diag_axis_init, register_diag_field, & diag_axis_add_attribute, diag_field_add_attribute, diag_send_complete, & - diag_manager_set_time_end, send_data, register_static_field + diag_manager_set_time_end, send_data, register_static_field, & + diag_field_add_cell_measures use platform_mod, only: r8_kind, r4_kind use fms_mod, only: fms_init, fms_end use mpp_mod, only: FATAL, mpp_error, mpp_npes, mpp_pe, mpp_root_pe, mpp_broadcast, input_nml_file @@ -178,12 +179,15 @@ program test_modern_diag if (id_var8 .ne. 8) call mpp_error(FATAL, "var8 does not have the expected id") endif +call diag_field_add_cell_measures(id_var6, area=id_var8, volume=id_var8) + call diag_field_add_attribute (id_var1, "some string", "this is a string") call diag_field_add_attribute (id_var1, "integer", 10) call diag_field_add_attribute (id_var1, "1d_integer", (/10, 10/)) call diag_field_add_attribute (id_var1, "real", 10.) call diag_field_add_attribute (id_var2, '1d_real', (/10./)) call diag_field_add_attribute (id_var2, 'formula', 'p(n,k,j,i) = ap(k) + b(k)*ps(n,j,i)') +call diag_field_add_attribute (id_var2, 'cell_methods', 'area: mullions') !! test dump routines !! prints fields from objects for debugging to log if name is provided, othwerise goes to stdout