FMS  2025.04
Flexible Modeling System
fms_diag_reduction_methods.F90
1 !***********************************************************************
2 !* Apache License 2.0
3 !*
4 !* This file is part of the GFDL Flexible Modeling System (FMS).
5 !*
6 !* Licensed under the Apache License, Version 2.0 (the "License");
7 !* you may not use this file except in compliance with the License.
8 !* You may obtain a copy of the License at
9 !*
10 !* http://www.apache.org/licenses/LICENSE-2.0
11 !*
12 !* FMS is distributed in the hope that it will be useful, but WITHOUT
13 !* WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied;
14 !* without even the implied warranty of MERCHANTABILITY or FITNESS FOR A
15 !* PARTICULAR PURPOSE. See the License for the specific language
16 !* governing permissions and limitations under the License.
17 !***********************************************************************
18 
19 !> @defgroup fms_diag_reduction_methods_mod fms_diag_reduction_methods_mod
20 !> @ingroup diag_manager
21 !! @brief fms_diag_reduction_methods_mod contains routines that are meant to be used for
22 !! error checking and setting up to do the reduction methods
23 
24 !> @file
25 !> @brief File for @ref fms_diag_reduction_methods_mod
26 
27 !> @addtogroup fms_diag_reduction_methods_mod
28 !> @{
29 module fms_diag_reduction_methods_mod
30  use platform_mod, only: r8_kind, r4_kind
31  use fms_diag_bbox_mod, only: fmsdiagibounds_type
32  use fms_string_utils_mod, only: string
33  use diag_data_mod, only: time_diurnal, time_rms
34  use mpp_mod
35  implicit none
36  private
37 
40 
41  !> @brief Does the time_none reduction method. See include/fms_diag_reduction_methods.inc
42  !TODO This needs to be extended to integers
43  interface do_time_none
44  module procedure do_time_none_r4, do_time_none_r8
45  end interface do_time_none
46 
47  !> @brief Does the time_min reduction method. See include/fms_diag_reduction_methods.inc
48  !TODO This needs to be extended to integers
49  interface do_time_min
50  module procedure do_time_min_r4, do_time_min_r8
51  end interface do_time_min
52 
53  !> @brief Does the time_max reduction method. See include/fms_diag_reduction_methods.inc
54  !TODO This needs to be extended to integers
55  interface do_time_max
56  module procedure do_time_max_r4, do_time_max_r8
57  end interface do_time_max
58 
59  !> @brief Sum update updates the buffer for any reductions that involve summation
60  !! (ie. time_sum, avg, rms, pow)
61  !!TODO This needs to be extended to integers
63  module procedure do_time_sum_update_r4, do_time_sum_update_r8
64  end interface
65 
66  !> @brief Finishes a reduction that involves an average
67  !! (ie. time_avg, rms, pow)
68  !! This takes the average at the end of the time step
69  interface time_update_done
70  module procedure sum_update_done_r4, sum_update_done_r8
71  end interface
72 
73  !> @brief Updates the buffer for any reductions that involve summation
74  !! (ie. time_sum, avg, rms, pow)
75  !! In this case the mask is present
76  interface sum_mask
77  module procedure sum_mask_r4, sum_mask_r8
78  end interface
79 
80  !> @brief Updates the buffer for any reductions that involve summation
81  !! (ie. time_sum, avg, rms, pow)
82  !! In this case the mask is present and it varies over time
83  interface sum_mask_variant
84  module procedure sum_mask_variant_r4, sum_mask_variant_r8
85  end interface sum_mask_variant
86 
87  !> @brief Updates the buffer for any reductions that involve summation
88  !! (ie. time_sum, avg, rms, pow)
89  !! In this case the mask is not present
90  interface sum_no_mask
91  module procedure sum_no_mask_r4, sum_no_mask_r8
92  end interface sum_no_mask
93 
94  contains
95 
96  !> @brief Checks improper combinations of is, ie, js, and je.
97  !! @return The error message, empty string if no errors were found
98  !> @note accept_data works in either one or another of two modes.
99  !! 1. Input field is a window (e.g. FMS physics)
100  !! 2. Input field includes halo data
101  !! It cannot handle a window of data that has halos.
102  !! (A field with no windows or halos can be thought of as a special case of either mode.)
103  !! The logic for indexing is quite different for these two modes, but is not clearly separated.
104  !! If both the beggining and ending indices are present, then field is assumed to have halos.
105  !! If only beggining indices are present, then field is assumed to be a window.
106  !> @par
107  !! There are a number of ways a user could mess up this logic, depending on the combination
108  !! of presence/absence of is,ie,js,je. The checks below should catch improper combinations.
109  pure function check_indices_order(is_in, ie_in, js_in, je_in) &
110  result(error_msg)
111  integer, intent(in), optional :: is_in, ie_in, js_in, je_in !< Indices passed to fms_diag_accept_data()
112  character(len=128) :: error_msg !< An error message used only for testing purpose!!!
113 
114  error_msg = ""
115  IF ( PRESENT(ie_in) ) THEN
116  IF ( .NOT.PRESENT(is_in) ) THEN
117  error_msg = 'ie_in present without is_in'
118  return
119  END IF
120  IF ( PRESENT(js_in) .AND. .NOT.PRESENT(je_in) ) THEN
121  error_msg = 'is_in and ie_in present, but js_in present without je_in'
122  return
123  END IF
124  END IF
125 
126  IF ( PRESENT(je_in) ) THEN
127  IF ( .NOT.PRESENT(js_in) ) THEN
128  error_msg = 'je_in present without js_in'
129  return
130  END IF
131  IF ( PRESENT(is_in) .AND. .NOT.PRESENT(ie_in) ) THEN
132  error_msg = 'js_in and je_in present, but is_in present without ie_in'
133  return
134  END IF
135  END IF
136  end function check_indices_order
137 
138  !> @brief Sets the logical mask based on mask or rmask
139  !> @return logical mask
140  function init_mask(rmask, mask, field) &
141  result(oor_mask)
142  LOGICAL, DIMENSION(:,:,:,:), allocatable, INTENT(in) :: mask !< The location of the mask
143  CLASS(*), DIMENSION(:,:,:,:), allocatable, INTENT(in) :: rmask !< The masking values
144  CLASS(*), DIMENSION(:,:,:,:), intent(in) :: field !< Field_data
145 
146  logical, allocatable, dimension(:,:,:,:) :: oor_mask !< mask
147 
148  ALLOCATE(oor_mask(SIZE(field, 1), SIZE(field, 2), SIZE(field, 3), SIZE(field, 4)))
149  oor_mask = .true.
150 
151  if (allocated(mask)) then
152  oor_mask = mask
153  elseif (allocated(rmask)) then
154  select type (rmask)
155  type is (real(kind=r8_kind))
156  WHERE (rmask < 0.5_r8_kind) oor_mask = .false.
157  type is (real(kind=r4_kind))
158  WHERE (rmask < 0.5_r4_kind) oor_mask = .false.
159  end select
160  endif
161 
162  end function init_mask
163 
164  !> @brief Sets the weight based on the weight passed into send_data (1.0_r8_kind if the weight is not passed in)
165  !! The weight will be saved as an r8 and converted to r4 as needed
166  !! @return weight to use when averaging
167  pure function set_weight(weight) &
168  result(out_weight)
169  CLASS(*), INTENT(in), OPTIONAL :: weight !< The weight use when averaging
170 
171  real(kind=r8_kind) :: out_weight
172 
173  out_weight = 1.0_r8_kind
174  if (present(weight)) then
175  select type(weight)
176  type is (real(kind=r8_kind))
177  out_weight = real(weight, kind = r8_kind)
178  type is (real(kind=r4_kind))
179  out_weight = real(weight, kind = r8_kind)
180  end select
181  endif
182  end function set_weight
183 
184 #include "fms_diag_reduction_methods_r4.fh"
185 #include "fms_diag_reduction_methods_r8.fh"
186 
187 end module fms_diag_reduction_methods_mod
188 !> @}
189 ! close documentation grouping
integer, parameter time_diurnal
The reduction method is diurnal.
Definition: diag_data.F90:122
integer, parameter time_rms
The reudction method is root mean square of values.
Definition: diag_data.F90:121
pure real(kind=r8_kind) function, public set_weight(weight)
Sets the weight based on the weight passed into send_data (1.0_r8_kind if the weight is not passed in...
logical function, dimension(:,:,:,:), allocatable, public init_mask(rmask, mask, field)
Sets the logical mask based on mask or rmask.
pure character(len=128) function, public check_indices_order(is_in, ie_in, js_in, je_in)
Checks improper combinations of is, ie, js, and je.
Does the time_max reduction method. See include/fms_diag_reduction_methods.inc.
Does the time_min reduction method. See include/fms_diag_reduction_methods.inc.
Sum update updates the buffer for any reductions that involve summation (ie. time_sum,...
Updates the buffer for any reductions that involve summation (ie. time_sum, avg, rms,...
Updates the buffer for any reductions that involve summation (ie. time_sum, avg, rms,...
Updates the buffer for any reductions that involve summation (ie. time_sum, avg, rms,...
Finishes a reduction that involves an average (ie. time_avg, rms, pow) This takes the average at the ...
character(:) function, allocatable, public string(v, fmt)
Converts a number or a Boolean value to a string.
Does the time_none reduction method. See include/fms_diag_reduction_methods.inc.
Data structure holding a 3D bounding box. It is commonlyused to represent the interval bounds or limi...