FMS  2025.04
Flexible Modeling System
diag_util.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 !> @defgroup diag_util_mod diag_util_mod
19 !> @ingroup diag_manager
20 !! @brief Functions and subroutines necessary for the <TT>diag_manager_mod</TT>.
21 !! @author Seth Underwood
22 
23 MODULE diag_util_mod
24 
25 use platform_mod
26 use,intrinsic :: iso_fortran_env, only: real128
27 use,intrinsic :: iso_c_binding, only: c_double,c_float,c_int64_t, &
28  c_int32_t,c_int16_t,c_intptr_t
29 
30  ! <FUTURE>
31  ! Make an interface <TT>check_bounds_are_exact</TT> for the subroutines <TT>check_bounds_are_exact_static</TT>
32  ! and
33  ! <TT>check_bounds_are_exact_dynamic</TT>.
34  ! <PRE>
35  ! INTERFACE check_bounds_are_exact
36  ! MODULE PROCEDURE check_bounds_are_exact_static
37  ! MODULE PROCEDURE check_bounds_are_exact_dynamic
38  ! END INTERFACE check_bounds_are_exact
39  ! </PRE>
40  ! </FUTURE>
41 
42  USE diag_data_mod, ONLY: output_fields, input_fields, files, do_diag_field_log, diag_log_unit,&
43  & very_large_axis_length, time_zero, very_large_file_freq, end_of_run, every_time,&
44  & diag_seconds, diag_minutes, diag_hours, diag_days, diag_months, diag_years, get_base_time,&
48  & mix_snapshot_average_fields, global_descriptor, cmor_missing_value, use_cmor, pack_size,&
52  USE diag_data_mod, ONLY: fileobju, fileobj, fnum_for_domain, fileobjnd
57  USE diag_output_mod, ONLY: diag_output_init, write_axis_meta_data,&
59  USE diag_output_mod, ONLY: diag_field_write, diag_write_time !<fms2_io
60  USE diag_grid_mod, ONLY: get_local_indexes
61  USE fms_diag_time_utils_mod, ONLY: diag_time_inc, get_time_string, get_date_dif
62  USE fms_mod, ONLY: error_mesg, fatal, warning, note, mpp_pe, mpp_root_pe, lowercase, fms_error_handler,&
63  & string, write_version_number
64  USE mpp_domains_mod,ONLY: domain1d, domain2d, mpp_get_compute_domain, null_domain1d, null_domain2d,&
65  & OPERATOR(.NE.), OPERATOR(.EQ.), mpp_modify_domain, mpp_get_domain_components,&
67  & domainug, null_domainug
68  USE time_manager_mod,ONLY: time_type, OPERATOR(==), OPERATOR(>), no_calendar, increment_date,&
70  & OPERATOR(<), OPERATOR(>=), OPERATOR(<=), OPERATOR(==)
71  USE mpp_mod, ONLY: mpp_npes
72  USE constants_mod, ONLY: seconds_per_day, seconds_per_hour, seconds_per_minute
73  USE fms2_io_mod
74  USE fms_diag_bbox_mod, ONLY: fmsdiagibounds_type
75  USE netcdf, ONLY: nf90_char
76 
77  IMPLICIT NONE
78  PRIVATE
86 
87 
88  !> @brief Prepend a value to a string attribute in the output field or output file.
89  !> @ingroup diag_util_mod
91  MODULE PROCEDURE prepend_attribute_field
92  MODULE PROCEDURE prepend_attribute_file
93  END INTERFACE prepend_attribute
94 
95  !> @brief Allocates the atttype in out_file.
96  !> @ingroup diag_util_mod
97  INTERFACE attribute_init
98  MODULE PROCEDURE attribute_init_field
99  MODULE PROCEDURE attribute_init_file
100  END INTERFACE attribute_init
101 
103  module procedure fms_diag_check_out_of_bounds_r4
104  module procedure fms_diag_check_out_of_bounds_r8
105  END INTERFACE fms_diag_check_out_of_bounds
106 
107 
108 !> @addtogroup diag_util_mod
109 !> @{
110 
111  ! Include variable "version" to be written to log file.
112 #include <file_version.h>
113 
114  LOGICAL :: module_initialized = .false.
115 
116  character(len=1), public :: field_log_separator = '|' !< separator used for csv-style log of registered fields
117  !! set by nml in diag_manager init
118 
119 
120 CONTAINS
121 
122  !> @brief Write the version number of this file to the log file.
123  SUBROUTINE diag_util_init()
124  IF (module_initialized) THEN
125  RETURN
126  END IF
127 
128  ! Write version number out to log file
129  call write_version_number("DIAG_UTIL_MOD", version)
130  END SUBROUTINE diag_util_init
131 
132  !> @brief Get the size, start, and end indices for output fields.
133  SUBROUTINE get_subfield_size(axes, outnum)
134  INTEGER, INTENT(in) :: axes(:) !< axes of the input_field
135  INTEGER, INTENT(in) :: outnum !< position in array output_fields
136 
137  REAL, ALLOCATABLE :: global_lat(:), global_lon(:), global_depth(:)
138  INTEGER :: global_axis_size, global_axis_sizey
139  INTEGER :: i,xbegin,xend,ybegin,yend,xbegin_l,xend_l,ybegin_l,yend_l
140  CHARACTER(len=1) :: cart
141  TYPE(domain2d) :: domain2, domain2_new
142  TYPE(domain1d) :: domain1, domain1x, domain1y
143  REAL :: start(3) !< start coordinates in 3 axes
144  REAL :: end(3) !< end coordinates in 3 axes
145  INTEGER :: gstart_indx(3)!< global start indices of output domain in 3 axes
146  INTEGER :: gend_indx(3) !< global end indices of output domain in 3 axes
147  REAL, ALLOCATABLE :: subaxis_x(:) !< containing local coordinates in x,y,z axes
148  REAL, ALLOCATABLE :: subaxis_y(:) !< containing local coordinates in x,y,z axes
149  REAL, ALLOCATABLE :: subaxis_z(:) !< containing local coordinates in x,y,z axes
150  CHARACTER(len=128) :: msg
151  INTEGER :: ishift, jshift
152  INTEGER :: grv !< Value used to determine if the region defined in the diag_table is for the whole
153  !! axis, or a sub-axis
154  CHARACTER(len=128), DIMENSION(2) :: axis_domain_name
155 
156  !initilization for local output
157  ! initially out of (lat/lon/depth) range
158  start = -1.e10
159  end = -1.e10
160  gstart_indx = -1
161  gend_indx = -1
162 
163  ! get the value to compare to determine if writing full axis data
164  IF ( region_out_use_alt_value ) THEN
165  grv = glo_reg_val_alt
166  ELSE
167  grv = glo_reg_val
168  END IF
169 
170  ! get axis data (lat, lon, depth) and indices
171  start = output_fields(outnum)%output_grid%start
172  end = output_fields(outnum)%output_grid%end
173 
174  CALL get_diag_axis_domain_name(axes(1), axis_domain_name(1))
175  CALL get_diag_axis_domain_name(axes(2), axis_domain_name(2))
176 
177  IF ( index(lowercase(axis_domain_name(1)), 'cubed') == 0 .AND. &
178  & index(lowercase(axis_domain_name(2)), 'cubed') == 0 ) THEN
179  DO i = 1, SIZE(axes(:))
180  global_axis_size = get_axis_global_length(axes(i))
181  output_fields(outnum)%output_grid%subaxes(i) = -1
182  CALL get_diag_axis_cart(axes(i), cart)
183  SELECT CASE(cart)
184  CASE ('X')
185  ! <ERROR STATUS="FATAL">wrong order of axes. X should come first.</ERROR>
186  IF( i.NE.1 ) CALL error_mesg('diag_util_mod::get_subfield_size',&
187  & 'wrong order of axes, X should come first',fatal)
188  ALLOCATE(global_lon(global_axis_size))
189  CALL get_diag_axis_data(axes(i),global_lon)
190  IF( int(start(i)) == grv .AND. int(end(i)) == grv ) THEN
191  gstart_indx(i) = 1
192  gend_indx(i) = global_axis_size
193  output_fields(outnum)%output_grid%subaxes(i) = axes(i)
194  ELSE
195  gstart_indx(i) = get_index(start(i),global_lon)
196  gend_indx(i) = get_index(end(i),global_lon)
197  END IF
198  ALLOCATE(subaxis_x(gstart_indx(i):gend_indx(i)))
199  subaxis_x=global_lon(gstart_indx(i):gend_indx(i))
200  CASE ('Y')
201  ! <ERROR STATUS="FATAL">wrong order of axes, Y should come second.</ERROR>
202  IF( i.NE.2 ) CALL error_mesg('diag_util_mod::get_subfield_size',&
203  & 'wrong order of axes, Y should come second',fatal)
204  ALLOCATE(global_lat(global_axis_size))
205  CALL get_diag_axis_data(axes(i),global_lat)
206  IF( int(start(i)) == grv .AND. int(end(i)) == grv ) THEN
207  gstart_indx(i) = 1
208  gend_indx(i) = global_axis_size
209  output_fields(outnum)%output_grid%subaxes(i) = axes(i)
210  ELSE
211  gstart_indx(i) = get_index(start(i),global_lat)
212  gend_indx(i) = get_index(end(i),global_lat)
213  END IF
214  ALLOCATE(subaxis_y(gstart_indx(i):gend_indx(i)))
215  subaxis_y=global_lat(gstart_indx(i):gend_indx(i))
216  CASE ('Z')
217  ! <ERROR STATUS="FATAL">wrong values in vertical axis of region</ERROR>
218  IF ( start(i)*end(i)<0. ) CALL error_mesg('diag_util_mod::get_subfield_size',&
219  & 'wrong values in vertical axis of region',fatal)
220  IF ( start(i)>=0. .AND. end(i)>0. ) THEN
221  ALLOCATE(global_depth(global_axis_size))
222  CALL get_diag_axis_data(axes(i),global_depth)
223  gstart_indx(i) = get_index(start(i),global_depth)
224  gend_indx(i) = get_index(end(i),global_depth)
225  ALLOCATE(subaxis_z(gstart_indx(i):gend_indx(i)))
226  subaxis_z=global_depth(gstart_indx(i):gend_indx(i))
227  output_fields(outnum)%output_grid%subaxes(i) =&
228  & diag_subaxes_init(axes(i),subaxis_z, gstart_indx(i),gend_indx(i))
229  DEALLOCATE(subaxis_z,global_depth)
230  ELSE ! regional vertical axis is the same as global vertical axis
231  gstart_indx(i) = 1
232  gend_indx(i) = global_axis_size
233  output_fields(outnum)%output_grid%subaxes(i) = axes(i)
234  ! <ERROR STATUS="FATAL">i should equal 3 for z axis</ERROR>
235  IF( i /= 3 ) CALL error_mesg('diag_util_mod::get_subfield_size',&
236  & 'i should equal 3 for z axis', fatal)
237  END IF
238  CASE default
239  ! <ERROR STATUS="FATAL">Wrong axis_cart</ERROR>
240  CALL error_mesg('diag_util_mod::get_subfield_size', 'Wrong axis_cart', fatal)
241  END SELECT
242  END DO
243 
244  DO i = 1, SIZE(axes(:))
245  IF ( gstart_indx(i) == -1 .OR. gend_indx(i) == -1 ) THEN
246  ! <ERROR STATUS="FATAL">
247  ! can not find gstart_indx/gend_indx for <output_fields(outnum)%output_name>,
248  ! check region bounds for axis <i>.
249  ! </ERROR>
250  WRITE(msg,'(A,I2)') ' check region bounds for axis ', i
251  CALL error_mesg('diag_util_mod::get_subfield_size', 'can not find gstart_indx/gend_indx for '&
252  & //trim(output_fields(outnum)%output_name)//','//trim(msg), fatal)
253  END IF
254  END DO
255  ELSE ! cubed sphere
256  ! get the i and j start and end indexes
257  CALL get_local_indexes(lonstart=start(1), lonend=end(1), &
258  & latstart=start(2), latend=end(2), &
259  & istart=gstart_indx(1), iend=gend_indx(1), &
260  & jstart=gstart_indx(2), jend=gend_indx(2))
261  global_axis_size = get_axis_global_length(axes(1))
262  ALLOCATE(global_lon(global_axis_size))
263  global_axis_sizey = get_axis_global_length(axes(2))
264  ALLOCATE(global_lat(global_axis_sizey))
265  CALL get_diag_axis_data(axes(1),global_lon)
266  CALL get_diag_axis_data(axes(2),global_lat)
267 
268  !Potential fix for out-of-bounds error for global_lon and global_lat.
269  IF ((gstart_indx(1) .GT. 0 .AND. gstart_indx(2) .GT. 0) .AND. &
270  (gstart_indx(1) .LE. global_axis_size .AND. gstart_indx(2) .LE. global_axis_sizey) .AND. &
271  (gend_indx(1) .GT. 0 .AND. gend_indx(2) .GT. 0) .AND. &
272  (gend_indx(1) .LE. global_axis_size .AND. gend_indx(2) .LE. global_axis_sizey)) THEN
273  ALLOCATE(subaxis_x(gstart_indx(1):gend_indx(1)))
274  ALLOCATE(subaxis_y(gstart_indx(2):gend_indx(2)))
275  subaxis_x=global_lon(gstart_indx(1):gend_indx(1))
276  subaxis_y=global_lat(gstart_indx(2):gend_indx(2))
277  END IF
278 
279  ! Now deal with the Z component
280  IF ( SIZE(axes(:)) > 2 ) THEN
281  global_axis_size = get_axis_global_length(axes(3))
282  output_fields(outnum)%output_grid%subaxes(3) = -1
283  CALL get_diag_axis_cart(axes(3), cart)
284  ! <ERROR STATUS="FATAL">
285  ! axis(3) should be Z-axis
286  ! </ERROR>
287  IF ( lowercase(cart) /= 'z' ) CALL error_mesg('diag_util_mod::get_subfield_size', &
288  &'axis(3) should be Z-axis', fatal)
289  ! <ERROR STATUS="FATAL">
290  ! wrong values in vertical axis of region
291  ! </ERROR>
292  IF ( start(3)*end(3)<0. ) CALL error_mesg('diag_util_mod::get_subfield_size',&
293  & 'wrong values in vertical axis of region',fatal)
294  IF ( start(3)>=0. .AND. end(3)>0. ) THEN
295  ALLOCATE(global_depth(global_axis_size))
296  CALL get_diag_axis_data(axes(3),global_depth)
297  gstart_indx(3) = get_index(start(3),global_depth)
298  IF( start(3) == 0.0 ) gstart_indx(3) = 1
299  gend_indx(3) = get_index(end(3),global_depth)
300  IF( start(3) >= maxval(global_depth) ) gstart_indx(3)= global_axis_size
301  IF( end(3) >= maxval(global_depth) ) gend_indx(3) = global_axis_size
302 
303  ALLOCATE(subaxis_z(gstart_indx(3):gend_indx(3)))
304  subaxis_z=global_depth(gstart_indx(3):gend_indx(3))
305  output_fields(outnum)%output_grid%subaxes(3) =&
306  & diag_subaxes_init(axes(3),subaxis_z, gstart_indx(3),gend_indx(3))
307  DEALLOCATE(subaxis_z,global_depth)
308  ELSE ! regional vertical axis is the same as global vertical axis
309  gstart_indx(3) = 1
310  gend_indx(3) = global_axis_size
311  output_fields(outnum)%output_grid%subaxes(3) = axes(3)
312  END IF
313  END IF
314  END IF
315 
316  ! get domain and compute_domain(xbegin,xend,ybegin,yend)
317  xbegin = -1
318  xend = -1
319  ybegin = -1
320  yend = -1
321 
322  domain2 = get_domain2d(axes)
323  IF ( domain2 .NE. null_domain2d ) THEN
324  CALL mpp_get_compute_domain(domain2, xbegin, xend, ybegin, yend)
325  CALL mpp_get_domain_components(domain2, domain1x, domain1y)
326  ELSE
327  DO i = 1, min(SIZE(axes(:)),2)
328  domain1 = get_domain1d(axes(i))
329  IF ( domain1 .NE. null_domain1d ) THEN
330  CALL get_diag_axis_cart(axes(i),cart)
331  SELECT CASE(cart)
332  CASE ('X')
333  domain1x = get_domain1d(axes(i))
334  CALL mpp_get_compute_domain(domain1x, xbegin, xend)
335  CASE ('Y')
336  domain1y = get_domain1d(axes(i))
337  CALL mpp_get_compute_domain(domain1y, ybegin, yend)
338  CASE default ! do nothing here
339  END SELECT
340  ELSE
341  ! <ERROR STATUS="FATAL">No domain available</ERROR>
342  CALL error_mesg('diag_util_mod::get_subfield_size', 'NO domain available', fatal)
343  END IF
344  END DO
345  END IF
346 
347  CALL get_axes_shift(axes, ishift, jshift)
348  xend = xend+ishift
349  yend = yend+jshift
350 
351  IF ( xbegin == -1 .OR. xend == -1 .OR. ybegin == -1 .OR. yend == -1 ) THEN
352  ! <ERROR STATUS="FATAL">wrong compute domain indices</ERROR>
353  CALL error_mesg('diag_util_mod::get_subfield_size', 'wrong compute domain indices',fatal)
354  END IF
355 
356  ! get the area containing BOTH compute domain AND local output area
357  IF( gstart_indx(1) > xend .OR. xbegin > gend_indx(1) ) THEN
358  output_fields(outnum)%output_grid%l_start_indx(1) = -1
359  output_fields(outnum)%output_grid%l_end_indx(1) = -1
360  output_fields(outnum)%need_compute = .false. ! not involved
361  ELSEIF ( gstart_indx(2) > yend .OR. ybegin > gend_indx(2) ) THEN
362  output_fields(outnum)%output_grid%l_start_indx(2) = -1
363  output_fields(outnum)%output_grid%l_end_indx(2) = -1
364  output_fields(outnum)%need_compute = .false. ! not involved
365  ELSE
366  output_fields(outnum)%output_grid%l_start_indx(1) = max(xbegin, gstart_indx(1))
367  output_fields(outnum)%output_grid%l_start_indx(2) = max(ybegin, gstart_indx(2))
368  output_fields(outnum)%output_grid%l_end_indx(1) = min(xend, gend_indx(1))
369  output_fields(outnum)%output_grid%l_end_indx(2) = min(yend, gend_indx(2))
370  output_fields(outnum)%need_compute = .true. ! involved in local output
371  END IF
372 
373  IF ( output_fields(outnum)%need_compute ) THEN
374  ! need to modify domain1d and domain2d for subaxes
375  xbegin_l = output_fields(outnum)%output_grid%l_start_indx(1)
376  xend_l = output_fields(outnum)%output_grid%l_end_indx(1)
377  ybegin_l = output_fields(outnum)%output_grid%l_start_indx(2)
378  yend_l = output_fields(outnum)%output_grid%l_end_indx(2)
379  CALL mpp_modify_domain(domain2, domain2_new, xbegin_l,xend_l, ybegin_l,yend_l,&
380  & gstart_indx(1),gend_indx(1), gstart_indx(2),gend_indx(2))
381 
382  output_fields(outnum)%output_grid%subaxes(1) =&
383  & diag_subaxes_init(axes(1),subaxis_x, gstart_indx(1),gend_indx(1),domain2_new)
384  output_fields(outnum)%output_grid%subaxes(2) =&
385  & diag_subaxes_init(axes(2),subaxis_y, gstart_indx(2),gend_indx(2),domain2_new)
386  DO i = 1, SIZE(axes(:))
387  IF ( output_fields(outnum)%output_grid%subaxes(i) == -1 ) THEN
388  ! <ERROR STATUS="FATAL">
389  ! <output_fields(outnum)%output_name> error at i = <i>
390  ! </ERROR>
391  WRITE(msg,'(a,"/",I4)') 'at i = ',i
392  CALL error_mesg('diag_util_mod::get_subfield_size '//trim(output_fields(outnum)%output_name),&
393  'error '//trim(msg), fatal)
394  END IF
395  END DO
396 
397  ! local start index should start from 1
398  output_fields(outnum)%output_grid%l_start_indx(1) = max(xbegin, gstart_indx(1)) - xbegin + 1
399  output_fields(outnum)%output_grid%l_start_indx(2) = max(ybegin, gstart_indx(2)) - ybegin + 1
400  output_fields(outnum)%output_grid%l_end_indx(1) = min(xend, gend_indx(1)) - xbegin + 1
401  output_fields(outnum)%output_grid%l_end_indx(2) = min(yend, gend_indx(2)) - ybegin + 1
402  IF ( SIZE(axes(:))>2 ) THEN
403  output_fields(outnum)%output_grid%l_start_indx(3) = gstart_indx(3)
404  output_fields(outnum)%output_grid%l_end_indx(3) = gend_indx(3)
405  ELSE
406  output_fields(outnum)%output_grid%l_start_indx(3) = 1
407  output_fields(outnum)%output_grid%l_end_indx(3) = 1
408  END IF
409  END IF
410  IF ( ALLOCATED(subaxis_x) ) DEALLOCATE(subaxis_x, global_lon)
411  IF ( ALLOCATED(subaxis_y) ) DEALLOCATE(subaxis_y, global_lat)
412  END SUBROUTINE get_subfield_size
413 
414  !> @brief Get size, start and end indices for output fields.
415  SUBROUTINE get_subfield_vert_size(axes, outnum)
416  INTEGER, DIMENSION(:), INTENT(in) :: axes !< axes of the input_field
417  INTEGER, INTENT(in) :: outnum !< position in array output_fields
418 
419  REAL, DIMENSION(3) :: start !< start and end coordinates in 3 axes
420  REAL, DIMENSION(3) :: end !< start and end coordinates in 3 axes
421  REAL, ALLOCATABLE, DIMENSION(:) :: global_depth
422  REAL, ALLOCATABLE, DIMENSION(:) :: subaxis_z !< containing local coordinates in x,y,z axes
423  INTEGER :: i, global_axis_size
424  INTEGER, DIMENSION(3) :: gstart_indx !< global start and end indices of output domain in 3 axes
425  INTEGER, DIMENSION(3) :: gend_indx !< global start and end indices of output domain in 3 axes
426  CHARACTER(len=1) :: cart
427  CHARACTER(len=128) :: msg
428 !----------
429 !ug support
430  integer :: vert_dim_num
431 !----------
432 
433  !initilization for local output
434  start = -1.e10
435  end = -1.e10 ! initially out of (lat/lon/depth) range
436  gstart_indx = -1
437  gend_indx=-1
438 
439  ! get axis data (lat, lon, depth) and indices
440  start= output_fields(outnum)%output_grid%start
441  end = output_fields(outnum)%output_grid%end
442 
443 !----------
444 !ug support
445  vert_dim_num = 3
446 !----------
447  DO i = 1, SIZE(axes(:))
448  global_axis_size = get_axis_global_length(axes(i))
449  output_fields(outnum)%output_grid%subaxes(i) = -1
450  CALL get_diag_axis_cart(axes(i), cart)
451  SELECT CASE(cart)
452  CASE ('X')
453  ! <ERROR STATUS="FATAL">wrong order of axes, X should come first</ERROR>
454  IF ( i.NE.1 ) CALL error_mesg('diag_util_mod::get_subfield_vert_size',&
455  & 'wrong order of axes, X should come first',fatal)
456  gstart_indx(i) = 1
457  gend_indx(i) = global_axis_size
458  output_fields(outnum)%output_grid%subaxes(i) = axes(i)
459  CASE ('Y')
460  ! <ERROR STATUS="FATAL">wrong order of axes, Y should come second</ERROR>
461  IF( i.NE.2 ) CALL error_mesg('diag_util_mod::get_subfield_vert_size',&
462  & 'wrong order of axes, Y should come second',fatal)
463  gstart_indx(i) = 1
464  gend_indx(i) = global_axis_size
465  output_fields(outnum)%output_grid%subaxes(i) = axes(i)
466 !----------
467 !ug support
468  case ("U")
469  if (i .ne. 1) then
470  call error_mesg("diag_util_mod::get_subfield_vert_size", &
471  "the unstructured axis must be the first dimension.", &
472  fatal)
473  endif
474  gstart_indx(i) = 1
475  gend_indx(i) = global_axis_size
476  output_fields(outnum)%output_grid%subaxes(i) = axes(i)
477  vert_dim_num = 2
478  start(vert_dim_num) = start(3)
479  end(vert_dim_num) = end(3)
480 !----------
481  CASE ('Z')
482 !----------
483 !ug support
484  if (i .ne. vert_dim_num) then
485  call error_mesg("diag_util_mod::get_subfield_vert_size",&
486  "i should equal vert_dim_num for z axis", &
487  fatal)
488  endif
489 !----------
490  ! <ERROR STATUS="FATAL">wrong values in vertical axis of region</ERROR>
491  IF( start(i)*end(i) < 0. ) CALL error_mesg('diag_util_mod::get_subfield_vert_size',&
492  & 'wrong values in vertical axis of region',fatal)
493  IF( start(i) >= 0. .AND. end(i) > 0. ) THEN
494  ALLOCATE(global_depth(global_axis_size))
495  CALL get_diag_axis_data(axes(i),global_depth)
496  gstart_indx(i) = get_index(start(i),global_depth)
497  IF( start(i) == 0.0 ) gstart_indx(i) = 1
498 
499  gend_indx(i) = get_index(end(i),global_depth)
500  IF( start(i) >= maxval(global_depth) ) gstart_indx(i)= global_axis_size
501  IF( end(i) >= maxval(global_depth) ) gend_indx(i) = global_axis_size
502 
503  ALLOCATE(subaxis_z(gstart_indx(i):gend_indx(i)))
504  subaxis_z=global_depth(gstart_indx(i):gend_indx(i))
505  output_fields(outnum)%output_grid%subaxes(i) = &
506  diag_subaxes_init(axes(i),subaxis_z, gstart_indx(i),gend_indx(i))
507  DEALLOCATE(subaxis_z,global_depth)
508  ELSE ! vertical axis is the same as global vertical axis
509  gstart_indx(i) = 1
510  gend_indx(i) = global_axis_size
511  output_fields(outnum)%output_grid%subaxes(i) = axes(i)
512  END IF
513  CASE default
514  ! <ERROR STATUS="FATAL">Wrong axis_cart</ERROR>
515  CALL error_mesg('diag_util_mod::get_subfield_vert_size', 'Wrong axis_cart', fatal)
516  END SELECT
517  END DO
518 
519  DO i = 1,SIZE(axes(:))
520  IF ( gstart_indx(i) == -1 .OR. gend_indx(i) == -1 ) THEN
521  ! <ERROR STATUS="FATAL">
522  ! can not find gstart_indx/gend_indx for <output_fields(outnum)%output_name>
523  ! check region bounds for axis
524  ! </ERROR>
525  WRITE(msg,'(A,I2)') ' check region bounds for axis ', i
526  CALL error_mesg('diag_util_mod::get_subfield_vert_size', 'can not find gstart_indx/gend_indx for '&
527  & //trim(output_fields(outnum)%output_name)//','//trim(msg), fatal)
528  END IF
529  END DO
530 
531  DO i= 1, 2
532  output_fields(outnum)%output_grid%l_start_indx(i) = gstart_indx(i)
533  output_fields(outnum)%output_grid%l_end_indx(i) = gend_indx(i)
534  END DO
535 
536  IF( SIZE(axes(:)) > 2 ) THEN
537  output_fields(outnum)%output_grid%l_start_indx(3) = gstart_indx(3)
538  output_fields(outnum)%output_grid%l_end_indx(3) = gend_indx(3)
539  ELSE
540  output_fields(outnum)%output_grid%l_start_indx(3) = 1
541  output_fields(outnum)%output_grid%l_end_indx(3) = 1
542  END IF
543  END SUBROUTINE get_subfield_vert_size
544 
545  !> @brief Find index <TT>i</TT> of array such that <TT>array(i)</TT> is closest to number.
546  INTEGER FUNCTION get_index(number, array)
547  REAL, INTENT(in) :: number !<
548  REAL, INTENT(in), DIMENSION(:) :: array !<
549 
550  INTEGER :: i, n
551  LOGICAL :: found
552 
553  n = SIZE(array(:))
554  ! check if array is monotonous
555  DO i = 2, n-1
556  IF( (array(i-1)<array(i).AND.array(i)>array(i+1)) .OR. (array(i-1)>array(i).AND.array(i)<array(i+1))) THEN
557  ! <ERROR STATUS="FATAL">array NOT monotonously ordered</ERROR>
558  CALL error_mesg('diag_util_mod::get_index', 'array NOT monotonously ordered',fatal)
559  END IF
560  END DO
561  get_index = -1
562  found = .false.
563  ! search in increasing array
564  DO i = 1, n-1
565  IF ( (array(i)<=number).AND.(array(i+1)>= number) ) THEN
566  IF( number - array(i) <= array(i+1) - number ) THEN
567  get_index = i
568  found=.true.
569  ELSE
570  get_index = i+1
571  found=.true.
572  ENDIF
573  EXIT
574  END IF
575  END DO
576  ! if not found, search in decreasing array
577  IF( .NOT.found ) THEN
578  DO i = 1, n-1
579  IF ( (array(i)>=number).AND.(array(i+1)<= number) ) THEN
580  IF ( array(i)-number <= number-array(i+1) ) THEN
581  get_index = i
582  found = .true.
583  ELSE
584  get_index = i+1
585  found = .true.
586  END IF
587  EXIT
588  END IF
589  END DO
590  END IF
591  ! if still not found, is it less than the first element
592  ! or greater than last element? (Increasing Array)
593  ! But it must be within 2x the axis spacing
594  ! i.e. array(1)-(array(3)-array(1)).LT.number .AND. or 2*array(1)-array(3).LT.number
595  IF ( .NOT. found ) THEN
596  IF ( 2*array(1)-array(3).LT.number .AND. number.LT.array(1) ) THEN
597  get_index = 1
598  found = .true.
599  ELSE IF ( array(n).LT.number .AND. number.LT.2*array(n)-array(n-2) ) THEN
600  get_index = n
601  found = .true.
602  ELSE
603  found = .false.
604  END IF
605  END IF
606 
607  ! if still not found, is it greater than the first element
608  ! or less than the last element? (Decreasing Array)
609  ! But it must be within 2x the axis spacing (see above)
610  IF ( .NOT. found ) THEN
611  IF ( 2*array(1)-array(3).GT.number .AND. number.GT.array(1) ) THEN
612  get_index = 1
613  found = .true.
614  ELSE IF ( array(n).GT.number .AND. number.GT.2*array(n)-array(n-2) ) THEN
615  get_index = n
616  found = .true.
617  ELSE
618  found = .false.
619  END IF
620  END IF
621  END FUNCTION get_index
622 
623  !> @brief Writes brief diagnostic field info to the log file.
624  !! @details If the <TT>do_diag_field_log</TT> namelist parameter is .TRUE.,
625  !! then a line briefly describing diagnostic field is added to
626  !! the log file. Normally users should not call this subroutine
627  !! directly, since it is called by register_static_field and
628  !! register_diag_field if do_not_log is not set to .TRUE.. It is
629  !! used, however, in LM3 to avoid excessive logs due to the
630  !! number of fields registered for each of the tile types. LM3
631  !! code uses a do_not_log parameter in the registration calls,
632  !! and subsequently calls this subroutine to log field information
633  !! under a generic name.
634  SUBROUTINE log_diag_field_info(module_name, field_name, axes, long_name, units,&
635  & missing_value, range, dynamic )
636  CHARACTER(len=*), INTENT(in) :: module_name !< Module name
637  CHARACTER(len=*), INTENT(in) :: field_name !< Field name
638  INTEGER, DIMENSION(:), INTENT(in) :: axes !< Axis IDs
639  CHARACTER(len=*), OPTIONAL, INTENT(in) :: long_name !< Long name for field.
640  CHARACTER(len=*), OPTIONAL, INTENT(in) :: units !< Unit of field.
641  CLASS(*), OPTIONAL, INTENT(in) :: missing_value !< Missing value value.
642  CLASS(*), DIMENSION(:), OPTIONAL, INTENT(IN) :: range !< Valid range of values for field.
643  LOGICAL, OPTIONAL, INTENT(in) :: dynamic !< <TT>.TRUE.</TT> if field is not static.
644 
645  ! ---- local vars
646  CHARACTER(len=256) :: lmodule, lfield, lname, lunits
647  CHARACTER(len=64) :: lmissval, lmin, lmax
648  CHARACTER(len=8) :: numaxis, timeaxis
649  CHARACTER(len=1) :: sep = '|'
650  CHARACTER(len=256) :: axis_name, axes_list
651  INTEGER :: i
652  REAL :: missing_value_use !< Local copy of missing_value
653  REAL, DIMENSION(2) :: range_use !< Local copy of range
654 
655  IF ( .NOT.do_diag_field_log ) RETURN
656  IF ( mpp_pe().NE.mpp_root_pe() ) RETURN
657 
658  ! Fatal error if range is present and its extent is not 2.
659  IF ( PRESENT(range) ) THEN
660  IF ( SIZE(range) .NE. 2 ) THEN
661  CALL error_mesg('diag_util_mod::fms_log_field_info', 'extent of range should be 2', fatal)
662  END IF
663  END IF
664 
665  lmodule = trim(module_name)
666  lfield = trim(field_name)
667 
668  IF ( PRESENT(long_name) ) THEN
669  lname = trim(long_name)
670  ELSE
671  lname = ''
672  END IF
673 
674  IF ( PRESENT(units) ) THEN
675  lunits = trim(units)
676  ELSE
677  lunits = ''
678  END IF
679 
680  WRITE (numaxis,'(i1)') SIZE(axes)
681 
682  IF (PRESENT(missing_value)) THEN
683  IF ( use_cmor ) THEN
684  WRITE (lmissval,*) cmor_missing_value
685  ELSE
686  SELECT TYPE (missing_value)
687  TYPE IS (real(kind=r4_kind))
688  missing_value_use = missing_value
689  TYPE IS (real(kind=r8_kind))
690  missing_value_use = real(missing_value)
691  CLASS DEFAULT
692  CALL error_mesg ('diag_util_mod::log_diag_field_info',&
693  & 'The missing_value is not one of the supported types of real(kind=4) or real(kind=8)', fatal)
694  END SELECT
695  WRITE (lmissval,*) missing_value_use
696  END IF
697  ELSE
698  lmissval = ''
699  ENDIF
700 
701  IF ( PRESENT(range) ) THEN
702  SELECT TYPE (range)
703  TYPE IS (real(kind=r4_kind))
704  range_use = range
705  TYPE IS (real(kind=r8_kind))
706  range_use = real(range)
707  CLASS DEFAULT
708  CALL error_mesg ('diag_util_mod::log_diag_field_info',&
709  & 'The range is not one of the supported types of real(kind=4) or real(kind=8)', fatal)
710  END SELECT
711  WRITE (lmin,*) range_use(1)
712  WRITE (lmax,*) range_use(2)
713  ELSE
714  lmin = ''
715  lmax = ''
716  END IF
717 
718  IF ( PRESENT(dynamic) ) THEN
719  IF (dynamic) THEN
720  timeaxis = 'T'
721  ELSE
722  timeaxis = 'F'
723  END IF
724  ELSE
725  timeaxis = ''
726  END IF
727 
728  axes_list=''
729  DO i = 1, SIZE(axes)
730  CALL get_diag_axis_name(axes(i),axis_name)
731  IF ( trim(axes_list) /= '' ) axes_list = trim(axes_list)//','
732  axes_list = trim(axes_list)//trim(axis_name)
733  END DO
734 
735  WRITE (diag_log_unit,'(777a)') &
736  & trim(lmodule), field_log_separator, trim(lfield), field_log_separator, trim(lname), field_log_separator,&
737  & trim(lunits), field_log_separator, trim(numaxis), field_log_separator, trim(timeaxis), field_log_separator,&
738  & trim(lmissval), field_log_separator, trim(lmin), field_log_separator, trim(lmax), field_log_separator,&
739  & trim(axes_list)
740  END SUBROUTINE log_diag_field_info
741 
742 
743 
744  !> @brief Update the <TT>output_fields</TT> x, y, and z min and max boundaries (array indices)
745  !! with the six specified bounds values.
746  SUBROUTINE update_bounds(out_num, lower_i, upper_i, lower_j, upper_j, lower_k, upper_k)
747  INTEGER, INTENT(in) :: out_num !< output field ID
748  INTEGER, INTENT(in) :: lower_i !< Lower i bound.
749  INTEGER, INTENT(in) :: upper_i !< Upper i bound.
750  INTEGER, INTENT(in) :: lower_j !< Lower j bound.
751  INTEGER, INTENT(in) :: upper_j !< Upper j bound.
752  INTEGER, INTENT(in) :: lower_k !< Lower k bound.
753  INTEGER, INTENT(in) :: upper_k !< Upper k bound.
754  CALL output_fields(out_num)%buff_bounds%update_bounds &
755  & ( lower_i, upper_i, lower_j, upper_j, lower_k, upper_k )
756  END SUBROUTINE update_bounds
757 
758 
759 
760  !> @brief Compares the bounding indices of an array specified in "current_bounds"
761 !! to the corresponding lower and upper bounds specified in "bounds"
762 !! Comparison is done by the two user specified input functions lowerb_comp and upperb_comp.
763 !! If any compariosn function returns true, then, after filling error_str, this routine also returns
764 !! true. The suplied comparison functions should return true for errors : for indices out of bounds,
765 !! or indices are not equal when expected to be equal.
766 LOGICAL FUNCTION compare_buffer_bounds_to_size(current_bounds, bounds, error_str, lowerb_comp, upperb_comp)
767  TYPE (fmsdiagibounds_type), INTENT(in) :: current_bounds !<A bounding box holding the current bounds
768  !! of an array.
769  TYPE (fmsdiagibounds_type), INTENT(in):: bounds !< The bounding box to check against.
770  CHARACTER(*), INTENT(out) :: error_str !< The return status, which is set to non-empty message
771  !! if the check fails.
772 
773  !> @brief Interface lowerb_comp should be used for comparison to lower bounds of buffer.
774  INTERFACE
775  LOGICAL FUNCTION lowerb_comp(a , b)
776  INTEGER, INTENT(IN) :: a !< One of the two args that are to be compared to each other.
777  INTEGER, INTENT(IN) :: b !< One of the two args that are to be compared to each other.
778  END FUNCTION lowerb_comp
779  END INTERFACE
780 
781  !> @brief Interface lowerb_comp should be used for comparison to upper bounds of buffer.
782  INTERFACE
783  LOGICAL FUNCTION upperb_comp(a, b)
784  INTEGER, INTENT(IN) :: a !< One of the two args that are to be compared to each other.
785  INTEGER, INTENT(IN) :: b !< One of the two args that are to be compared to each other.
786  END FUNCTION upperb_comp
787  END INTERFACE
788 
790 
791  IF (lowerb_comp( bounds%get_imin() , current_bounds%get_imin()) .OR. &
792  upperb_comp( bounds%get_imax() , current_bounds%get_imax()).OR.&
793  lowerb_comp( bounds%get_jmin() , current_bounds%get_jmin()) .OR.&
794  upperb_comp( bounds%get_jmax() , current_bounds%get_jmax()) .OR.&
795  lowerb_comp( bounds%get_kmin() , current_bounds%get_kmin()) .OR.&
796  upperb_comp( bounds%get_kmax() , current_bounds%get_kmax())) THEN
798  error_str ='Buffer bounds= : , : , : Actual bounds= : , : , : '
799  WRITE(error_str(15:17),'(i3)') current_bounds%get_imin()
800  WRITE(error_str(19:21),'(i3)') current_bounds%get_imax()
801  WRITE(error_str(23:25),'(i3)') current_bounds%get_jmin()
802  WRITE(error_str(27:29),'(i3)') current_bounds%get_jmax()
803  WRITE(error_str(31:33),'(i3)') current_bounds%get_kmin()
804  WRITE(error_str(35:37),'(i3)') current_bounds%get_kmax()
805  WRITE(error_str(54:56),'(i3)') bounds%get_imin()
806  WRITE(error_str(58:60),'(i3)') bounds%get_imax()
807  WRITE(error_str(62:64),'(i3)') bounds%get_jmin()
808  WRITE(error_str(66:68),'(i3)') bounds%get_jmax()
809  WRITE(error_str(70:72),'(i3)') bounds%get_kmin()
810  WRITE(error_str(74:76),'(i3)') bounds%get_kmax()
811  ELSE
813  error_str = ''
814  END IF
816 
817 !> @brief return true iff a<b.
818 LOGICAL FUNCTION a_lessthan_b(a , b)
819  INTEGER, INTENT(IN) :: a !< The first of the two integer args that are to be compared to each other.
820  INTEGER, INTENT(IN) :: b !< The first of the two integer args that are to be compared to each other.
821  a_lessthan_b = a < b
822 END FUNCTION a_lessthan_b
823 
824 !> @brief return true iff a>b.
825 LOGICAL FUNCTION a_greaterthan_b(a, b)
826  INTEGER, INTENT(IN) :: a !< The first of the two integer args that are to be compared to each other.
827  INTEGER, INTENT(IN) :: b !< The first of the two integer args that are to be compared to each other.
828  a_greaterthan_b = a > b
829 END FUNCTION a_greaterthan_b
830 
831 !> @brief return true iff a /= b
832 LOGICAL FUNCTION a_noteq_b(a, b)
833 INTEGER, INTENT(IN) :: a !< The first of the two integer args that are to be compared to each other.
834 INTEGER, INTENT(IN) :: b !< The first of the two integer args that are to be compared to each other.
835 a_noteq_b = a /= b
836 END FUNCTION a_noteq_b
837 
838  !> @brief Checks if the array indices for <TT>output_fields(out_num)</TT> are outside the
839  !! <TT>output_fields(out_num)%buffer</TT> upper and lower bounds.
840  !! If there is an error then error message will be filled.
841 SUBROUTINE check_out_of_bounds(out_num, diag_field_id, err_msg)
842  INTEGER, INTENT(in) :: out_num !< Output field ID number.
843  INTEGER, INTENT(in) :: diag_field_id !< Input field ID number.
844  CHARACTER(len=*), INTENT(out) :: err_msg !< Return status of <TT>check_out_of_bounds</TT>. An empty
845  !! error string indicates the x, y, and z indices are not outside the
846 
847  CHARACTER(len=128) :: error_string1, error_string2
848  LOGICAL :: out_of_bounds = .true.
849  TYPE (fmsdiagibounds_type) :: array_bounds
850  associate(buff_bounds => output_fields(out_num)%buff_bounds)
851 
852  CALL array_bounds%reset_bounds_from_array_4D(output_fields(out_num)%buffer)
853 
854  out_of_bounds = compare_buffer_bounds_to_size(array_bounds, buff_bounds, &
855  & error_string2, a_lessthan_b, a_greaterthan_b)
856 
857  IF (out_of_bounds .EQV. .true.) THEN
858  WRITE(error_string1,'(a,"/",a)') trim(input_fields(diag_field_id)%module_name), &
859  & trim(output_fields(out_num)%output_name)
860  err_msg = 'module/output_field='//trim(error_string1)//&
861  & ' Bounds of buffer exceeded. '//trim(error_string2)
862  ! imax, imin, etc need to be reset in case the program is not terminated.
863  call buff_bounds%reset(very_large_axis_length, 0)
864  ELSE
865  err_msg = ''
866  END IF
867  end associate
868 END SUBROUTINE check_out_of_bounds
869 
870  !> @brief Checks if the array indices for <TT>output_fields(out_num)</TT> are outside the
871  !! <TT>output_fields(out_num)%buffer</TT> upper and lower bounds.
872  !! If there is an error then error message will be filled.
873 SUBROUTINE fms_diag_check_out_of_bounds_r4(ofb, bounds, output_name, module_name, err_msg)
874  REAL(kind=r4_kind), INTENT (in), DIMENSION(:,:,:,:,:) :: ofb !< The output field buffer to check
875  TYPE (fmsDiagIbounds_type), INTENT(inout) :: bounds !< The bounding box to check against
876  CHARACTER(:), ALLOCATABLE, INTENT(in) :: output_name !< output name for placing in error message
877  CHARACTER(:), ALLOCATABLE, INTENT(in) :: module_name !< module name for placing in error message
878  CHARACTER(len=*), INTENT(inout) :: err_msg !< Return status of <TT>check_out_of_bounds</TT>. An empty
879  !! error string indicates the x, y, and z indices are not outside the
880 
881  CHARACTER(len=128) :: error_string1, error_string2
882  LOGICAL :: out_of_bounds = .true.
883  TYPE (fmsDiagIbounds_type) :: array_bounds
884 
885  CALL array_bounds%reset_bounds_from_array_5D(ofb)
886 
887  out_of_bounds = compare_buffer_bounds_to_size(array_bounds, bounds, &
888  & error_string2, a_lessthan_b, a_greaterthan_b)
889 
890  IF (out_of_bounds .EQV. .true.) THEN
891  WRITE(error_string1,'(a,"/",a)') trim(module_name), trim(output_name)
892  err_msg = 'module/output_field='//trim(error_string1)//&
893  & ' Bounds of buffer exceeded. '//trim(error_string2)
894  ! imax, imin, etc need to be reset in case the program is not terminated.
895  call bounds%reset(very_large_axis_length,0)
896  ELSE
897  err_msg = ''
898  END IF
899 END SUBROUTINE fms_diag_check_out_of_bounds_r4
900 
901  !> @brief Checks if the array indices for output_field buffer (ofb) are outside the
902  !! are outside the bounding box (bounds).
903  !! If there is an error then error message will be filled.
904 
905 SUBROUTINE fms_diag_check_out_of_bounds_r8(ofb, bounds, output_name, module_name, err_msg)
906  REAL(kind=r8_kind), INTENT (in), DIMENSION(:,:,:,:,:) :: ofb !< The output field buffer to check
907  TYPE (fmsDiagIbounds_type), INTENT(inout) :: bounds !< The bounding box to check against
908  CHARACTER(:), ALLOCATABLE, INTENT(in) :: output_name !< output name for placing in error message
909  CHARACTER(:), ALLOCATABLE, INTENT(in) :: module_name !< module name for placing in error message
910  CHARACTER(len=*), INTENT(out) :: err_msg !< Return status of <TT>check_out_of_bounds</TT>. An empty
911  !! error string indicates the x, y, and z indices are not outside the
912 
913  CHARACTER(len=128) :: error_string1, error_string2
914  LOGICAL :: out_of_bounds = .true.
915  TYPE (fmsDiagIbounds_type) :: array_bounds !<A bounding box holdstore the current bounds ofb
916 
917  CALL array_bounds%reset_bounds_from_array_5D(ofb)
918 
919  out_of_bounds = compare_buffer_bounds_to_size(array_bounds, bounds, &
920  & error_string2, a_lessthan_b, a_greaterthan_b)
921 
922  IF (out_of_bounds .EQV. .true.) THEN
923  WRITE(error_string1,'(a,"/",a)') trim(module_name), trim(output_name)
924  err_msg = 'module/output_field='//trim(error_string1)//&
925  & ' Bounds of buffer exceeded. '//trim(error_string2)
926  ! imax, imin, etc need to be reset in case the program is not terminated.
927  call bounds%reset(very_large_axis_length,0)
928  ELSE
929  err_msg = ''
930  END IF
931 END SUBROUTINE fms_diag_check_out_of_bounds_r8
932 
933 
934 !> @brief Checks that array indices specified in the bounding box "current_bounds"
935 !! are identical to those in the bounding box "bounds" match exactly. The check
936 !! occurs only when the time changed.
937 !! If there is an error then error message will be filled.
938 SUBROUTINE fms_diag_check_bounds_are_exact_dynamic(current_bounds, bounds, output_name, module_name, &
939  & Time, field_prev_Time, err_msg)
940  TYPE (fmsdiagibounds_type), INTENT(in) :: current_bounds !<A bounding box holding the current bounds
941  !! of an array.
942  TYPE (fmsdiagibounds_type), INTENT(inout) :: bounds !<The bounding box to check against
943  CHARACTER(:), ALLOCATABLE, INTENT(in) :: output_name !< output name for placing in error message
944  CHARACTER(:), ALLOCATABLE, INTENT(in) :: module_name !< module name for placing in error message
945  TYPE(time_type), INTENT(in) :: time !< Time to use in check. The check is only performed if
946  !! <TT>output_fields(out_num)%Time_of_prev_field_data</TT> is not
947  !! equal to <TT>Time</TT> or <TT>Time_zero</TT>.
948  TYPE(time_type), INTENT(inout) :: field_prev_time !< <TT>output_fields(out_num)%Time_of_prev_field_data</TT>
949  CHARACTER(len=*), INTENT(out) :: err_msg !< Return status of <TT>check_bounds_are_exact_dynamic</TT>.
950  !! An empty error string indicates the x, y, and z indices are
951  !! equal to the buffer array boundaries.
952 
953  CHARACTER(len=128) :: error_string1, error_string2
954  LOGICAL :: do_check
955  LOGICAL :: lims_not_exact
956 
957  err_msg = ''
958 
959  ! Check bounds only when the value of Time changes. When windows are used,
960  ! a change in Time indicates that a new loop through the windows has begun,
961  ! so a check of the previous loop can be done.
962  IF ( time == field_prev_time ) THEN
963  do_check = .false.
964  ELSE
965  IF ( field_prev_time == time_zero ) THEN
966  ! It may or may not be OK to check, I don't know how to tell.
967  ! Check will be done on subsequent calls anyway.
968  do_check = .false.
969  ELSE
970  do_check = .true.
971  END IF
972  field_prev_time = time
973  END IF
974 
975  IF ( do_check ) THEN
976  lims_not_exact = compare_buffer_bounds_to_size(current_bounds, bounds, &
977  & error_string2, a_noteq_b, a_noteq_b)
978  IF( lims_not_exact .eqv. .true.) THEN
979  WRITE(error_string1,'(a,"/",a)') trim(module_name), trim(output_name)
980  err_msg = trim(error_string1)//' Bounds of data do not match those of buffer. '//trim(error_string2)
981  END IF
982  call bounds%reset(very_large_axis_length, 0)
983  END IF
985 
986 
987 !> @brief This is an adaptor to the check_bounds_are_exact_dynamic_modern function to
988 !! maintain an interface servicing the legacy diag_manager.
989 SUBROUTINE check_bounds_are_exact_dynamic(out_num, diag_field_id, Time, err_msg)
990  INTEGER, INTENT(in) :: out_num !< Output field ID number.
991  INTEGER, INTENT(in) :: diag_field_id !< Input field ID number.
992  TYPE(time_type), INTENT(in) :: time !< Time to use in check. The check is only performed if
993  !! <TT>output_fields(out_num)%Time_of_prev_field_data</TT> is not
994  !! equal to <TT>Time</TT> or <TT>Time_zero</TT>.
995  CHARACTER(len=*), INTENT(out) :: err_msg !< Return status of <TT>check_bounds_are_exact_dynamic</TT>.
996  !! An empty error string indicates the x, y, and z indices are
997  !! equal to the buffer array boundaries.
998  CHARACTER(:), ALLOCATABLE :: output_name !< output name for placing in error message
999  CHARACTER(:), ALLOCATABLE :: module_name !< module name for placing in error message
1000  TYPE (fmsdiagibounds_type) :: current_bounds !< a bounding box to store the current bounds of the array.
1001 
1002  output_name = output_fields(out_num)%output_name
1003  module_name = input_fields(diag_field_id)%module_name
1004 
1005  CALL current_bounds%reset_bounds_from_array_4D(output_fields(out_num)%buffer)
1006 
1007  CALL fms_diag_check_bounds_are_exact_dynamic(current_bounds, output_fields(out_num)%buff_bounds, &
1008  & output_name, module_name, &
1009  & time, output_fields(out_num)%Time_of_prev_field_data, err_msg)
1010 
1011 END SUBROUTINE check_bounds_are_exact_dynamic
1012 
1013 
1014  !> @brief Check if the array indices for <TT>output_fields(out_num)</TT> are equal to the
1015  !! <TT>output_fields(out_num)%buffer</TT> upper and lower bounds.
1016  SUBROUTINE check_bounds_are_exact_static(out_num, diag_field_id, err_msg)
1017  INTEGER, INTENT(in) :: out_num !< Output field ID
1018  INTEGER, INTENT(in) :: diag_field_id !< Input field ID.
1019  CHARACTER(len=*), INTENT(out) :: err_msg !< The return status, which is set to non-empty message
1020  !! if the check fails.
1021  CHARACTER(:), ALLOCATABLE :: output_name !< output name for placing in error message
1022  CHARACTER(:), ALLOCATABLE :: module_name !< output name for placing in error message
1023  TYPE (fmsdiagibounds_type) :: current_bounds !< a bounding box to store the current bounds of the array.
1024 
1025  output_name = output_fields(out_num)%output_name
1026  module_name = input_fields(diag_field_id)%module_name
1027 
1028  CALL current_bounds%reset_bounds_from_array_4D(output_fields(out_num)%buffer)
1029 
1030  CALL fms_diag_check_bounds_are_exact_static(current_bounds, output_fields(out_num)%buff_bounds, &
1031  & output_name, module_name, err_msg)
1032  END SUBROUTINE check_bounds_are_exact_static
1033 
1034 
1035  !> @brief Check if the array indices specified in the bounding box "current_bounds" are equal to those
1036  !! specified in the bounding box "bounds" output_fields are equal to the buffer upper and lower bounds.
1037  !! If there is an error then error message will be filled.
1038  SUBROUTINE fms_diag_check_bounds_are_exact_static(current_bounds, bounds, output_name, module_name, err_msg)
1039  TYPE (fmsdiagibounds_type), INTENT(in) :: current_bounds !<A bounding box holding the current bounds
1040  !! of the array.
1041  TYPE (fmsdiagibounds_type), INTENT(inout) :: bounds !<The bounding box to check against
1042  CHARACTER(:), ALLOCATABLE, INTENT(in) :: output_name !< output name for placing in error message
1043  CHARACTER(:), ALLOCATABLE, INTENT(in) :: module_name !< module name for placing in error message
1044  CHARACTER(len=*), INTENT(out) :: err_msg !< The return status, which is set to non-empty message
1045  !! if the check fails.
1046 
1047  CHARACTER(len=128) :: error_string1, error_string2
1048  LOGICAL :: lims_not_exact
1049 
1050  err_msg = ''
1051  lims_not_exact = compare_buffer_bounds_to_size(current_bounds, bounds, &
1052  & error_string2, a_noteq_b, a_noteq_b)
1053  IF( lims_not_exact .eqv. .true.) THEN
1054  WRITE(error_string1,'(a,"/",a)') trim(module_name), trim(output_name)
1055  err_msg = trim(error_string1)//' Bounds of data do not match those of buffer. '//trim(error_string2)
1056  END IF
1057  call bounds%reset(very_large_axis_length, 0)
1059 
1060 
1061  !> @brief Initialize the output file.
1062  SUBROUTINE init_file(name, output_freq, output_units, format, time_units, long_name, tile_count,&
1063  & new_file_freq, new_file_freq_units, start_time, file_duration, file_duration_units, filename_time_bounds)
1064  CHARACTER(len=*), INTENT(in) :: name !< File name.
1065  CHARACTER(len=*), INTENT(in) :: long_name !< Long name for time axis.
1066  INTEGER, INTENT(in) :: output_freq !< How often data is to be written to the file.
1067  INTEGER, INTENT(in) :: output_units !< The output frequency unit. (MIN, HOURS, DAYS, etc.)
1068  INTEGER, INTENT(in) :: format !< Number type/kind the data is to be written out to the file.
1069  INTEGER, INTENT(in) :: time_units !< Time axis units.
1070  INTEGER, INTENT(in) :: tile_count !< Tile number.
1071  INTEGER, INTENT(in), OPTIONAL :: new_file_freq !< How often a new file is to be created.
1072  INTEGER, INTENT(in), OPTIONAL :: new_file_freq_units !< The new file frequency unit. (MIN, HOURS, DAYS, etc.)</IN>
1073  INTEGER, INTENT(in), OPTIONAL :: file_duration !< How long file is to be used.
1074  INTEGER, INTENT(in), OPTIONAL :: file_duration_units !< File duration unit. (MIN, HOURS, DAYS, etc.)
1075  TYPE(time_type), INTENT(in), OPTIONAL :: start_time !< Time when the file is to start
1076  CHARACTER(len=*), INTENT(in), OPTIONAL :: filename_time_bounds !< Set time bound file name to
1077  ! use "begin" "middle" or "end"
1078  INTEGER :: new_file_freq1, new_file_freq_units1
1079  INTEGER :: file_duration1, file_duration_units1
1080  INTEGER :: n
1081  LOGICAL :: same_file_err !< .FALSE. indicates that if the file name had
1082  !! previously been registered, this new file
1083  !! contained differences from the previous.
1084  REAL, DIMENSION(1) :: tdata
1085  CHARACTER(len=128) :: time_units_str
1086 
1087  ! Check if this file has already been defined
1088  same_file_err=.false. ! To indicate that if this file was previously defined
1089  ! no differences in this registration was detected.
1090  DO n=1,num_files
1091  IF ( trim(files(n)%name) == trim(name) ) THEN
1092  ! File is defined, check if all inputs are the same
1093  ! Start with the required parameters
1094  IF ( files(n)%output_freq.NE.output_freq .OR.&
1095  & files(n)%output_units.NE.output_units .OR.&
1096  & files(n)%format.NE.format .OR.&
1097  & files(n)%time_units.NE.time_units .OR.&
1098  & trim(files(n)%long_name).NE.trim(long_name) .OR.&
1099  & files(n)%tile_count.NE.tile_count ) THEN
1100  same_file_err=.true.
1101  END IF
1102 
1103  ! Now check the OPTIONAL parameters
1104  IF ( PRESENT(new_file_freq) ) THEN
1105  IF ( files(n)%new_file_freq.NE.new_file_freq ) THEN
1106  same_file_err=.true.
1107  END IF
1108  END IF
1109 
1110  IF ( PRESENT(new_file_freq_units) ) THEN
1111  IF ( files(n)%new_file_freq_units.NE.new_file_freq_units ) THEN
1112  same_file_err=.true.
1113  END IF
1114  END IF
1115 
1116  IF ( PRESENT(start_time) ) THEN
1117  IF ( files(n)%start_time==start_time ) THEN
1118  same_file_err=.true.
1119  END IF
1120  END IF
1121 
1122  IF ( PRESENT(file_duration) ) THEN
1123  IF ( files(n)%duration.NE.file_duration) THEN
1124  same_file_err=.true.
1125  END IF
1126  END IF
1127 
1128  IF ( PRESENT(file_duration_units) ) THEN
1129  IF ( files(n)%duration_units.NE.file_duration_units ) THEN
1130  same_file_err=.true.
1131  END IF
1132  END IF
1133 
1134  ! If the same file was defined twice, simply return, else FATAL
1135  IF ( same_file_err ) THEN
1136  ! Something in this file is not identical to the previously defined
1137  ! file of the same name. FATAL
1138  CALL error_mesg('diag_util_mod::init_file',&
1139  & 'The file "'//trim(name)//'" is defined multiple times in&
1140  & the diag_table.', fatal)
1141  ELSE
1142  ! Issue a note that the same file is defined multiple times
1143  CALL error_mesg('diag_util_mod::init_file',&
1144  & 'The file "'//trim(name)//'" is defined multiple times in&
1145  & the diag_table.', note)
1146  ! Return to the calling function
1147  RETURN
1148  END IF
1149  END IF
1150  END DO
1151 
1152  ! Get a number for this file
1153  num_files = num_files + 1
1154  IF ( num_files >= max_files ) THEN
1155  ! <ERROR STATUS="FATAL">
1156  ! max_files exceeded, increase max_files via the max_files variable
1157  ! in the namelist diag_manager_nml.
1158  ! </ERROR>
1159  CALL error_mesg('diag_util_mod::init_file',&
1160  & ' max_files exceeded, increase max_files via the max_files variable&
1161  & in the namelist diag_manager_nml.', fatal)
1162  END IF
1163 
1164  IF ( PRESENT(new_file_freq) ) THEN
1165  new_file_freq1 = new_file_freq
1166  ELSE
1167  new_file_freq1 = very_large_file_freq
1168  END IF
1169 
1170  IF ( PRESENT(new_file_freq_units) ) THEN
1171  new_file_freq_units1 = new_file_freq_units
1172  ELSE IF ( get_calendar_type() == no_calendar ) THEN
1173  new_file_freq_units1 = diag_days
1174  ELSE
1175  new_file_freq_units1 = diag_years
1176  END IF
1177 
1178  IF ( PRESENT(file_duration) ) THEN
1179  file_duration1 = file_duration
1180  ELSE
1181  file_duration1 = new_file_freq1
1182  END IF
1183 
1184  IF ( PRESENT(file_duration_units) ) THEN
1185  file_duration_units1 = file_duration_units
1186  ELSE
1187  file_duration_units1 = new_file_freq_units1
1188  END IF
1189 
1190  files(num_files)%tile_count = tile_count
1191  files(num_files)%name = trim(name)
1192  files(num_files)%output_freq = output_freq
1193  files(num_files)%output_units = output_units
1194  files(num_files)%format = FORMAT
1195  files(num_files)%time_units = time_units
1196  files(num_files)%long_name = trim(long_name)
1197  files(num_files)%num_fields = 0
1198  files(num_files)%local = .false.
1199  files(num_files)%last_flush = get_base_time()
1200  files(num_files)%file_unit = -1
1201  files(num_files)%new_file_freq = new_file_freq1
1202  files(num_files)%new_file_freq_units = new_file_freq_units1
1203  files(num_files)%duration = file_duration1
1204  files(num_files)%duration_units = file_duration_units1
1205 !> Initialize the times to 0
1206  files(num_files)%rtime_current = -1.0
1207  files(num_files)%time_index = 0
1208  files(num_files)%filename_time_bounds = filename_time_bounds
1209 
1210  IF ( PRESENT(start_time) ) THEN
1211  files(num_files)%start_time = start_time
1212  ELSE
1213  files(num_files)%start_time = get_base_time()
1214  END IF
1215  files(num_files)%next_open=diag_time_inc(files(num_files)%start_time,new_file_freq1,new_file_freq_units1)
1216  files(num_files)%close_time = diag_time_inc(files(num_files)%start_time,file_duration1, file_duration_units1)
1217  IF ( files(num_files)%close_time>files(num_files)%next_open ) THEN
1218  ! <ERROR STATUS="FATAL">
1219  ! close time GREATER than next_open time, check file duration,
1220  ! file frequency in <files(num_files)%name>
1221  ! </ERROR>
1222  CALL error_mesg('diag_util_mod::init_file', 'close time GREATER than next_open time, check file duration,&
1223  & file frequency in '//files(num_files)%name, fatal)
1224  END IF
1225 
1226  ! add time_axis_id and time_bounds_id here
1227  WRITE(time_units_str, 11) trim(time_unit_list(files(num_files)%time_units)), get_base_year(),&
1228  & get_base_month(), get_base_day(), get_base_hour(), get_base_minute(), get_base_second()
1229 11 FORMAT(a, ' since ', i4.4, '-', i2.2, '-', i2.2, ' ', i2.2, ':', i2.2, ':', i2.2)
1230  files(num_files)%time_axis_id = diag_axis_init(trim(long_name), tdata, time_units_str, 'T',&
1231  & trim(long_name) , set_name=trim(name) )
1232  !---- register axis for storing time boundaries
1233  files(num_files)%time_bounds_id = diag_axis_init( 'nv',(/1.,2./),'none','N','vertex number',&
1234  & set_name='nv')
1235  END SUBROUTINE init_file
1236 
1237  !> @brief Synchronize the file's start and close times with the model start and end times.
1238  !! @details <TT>sync_file_times</TT> checks to see if the file start time is less than the
1239  !! model's init time (passed in as the only argument). If it is less, then the
1240  !! both the file start time and end time are synchronized using the passed in initial time
1241  !! and the duration as calculated by the <TT>diag_time_inc</TT> function. <TT>sync_file_times</TT>
1242  !! will also increase the <TT>next_open</TT> until it is greater than the init_time.
1243  SUBROUTINE sync_file_times(file_id, init_time, err_msg)
1244  INTEGER, INTENT(in) :: file_id !< The file ID
1245  TYPE(time_type), INTENT(in) :: init_time !< Initial time use for the synchronization.
1246  CHARACTER(len=*), OPTIONAL, INTENT(out) :: err_msg !< Return error message
1247 
1248  CHARACTER(len=128) :: msg
1249 
1250  IF ( PRESENT(err_msg) ) err_msg = ''
1251 
1252  IF ( files(file_id)%start_time < init_time ) THEN
1253  ! Sync the start_time of the file with the initial time of the model
1254  files(file_id)%start_time = init_time
1255  ! Sync the file's close time also
1256  files(file_id)%close_time = diag_time_inc(files(file_id)%start_time,&
1257  & files(file_id)%duration, files(file_id)%duration_units)
1258  END IF
1259 
1260  ! Need to increase next_open until it is greate than init_time
1261  DO WHILE ( files(file_id)%next_open <= init_time )
1262  files(file_id)%next_open = diag_time_inc(files(file_id)%next_open,&
1263  & files(file_id)%new_file_freq, files(file_id)%new_file_freq_units, err_msg=msg)
1264  IF ( msg /= '' ) THEN
1265  IF ( fms_error_handler('diag_util_mod::sync_file_times',&
1266  & ' file='//trim(files(file_id)%name)//': '//trim(msg), err_msg) ) RETURN
1267  END IF
1268  END DO
1269  END SUBROUTINE sync_file_times
1270 
1271  !> @brief Return the file number for file name and tile.
1272  !! @return Integer find_file
1273  INTEGER FUNCTION find_file(name, tile_count)
1274  INTEGER, INTENT(in) :: tile_count !< Tile number.
1275  CHARACTER(len=*), INTENT(in) :: name !< File name.
1276 
1277  INTEGER :: i
1278 
1279  find_file = -1
1280  DO i = 1, num_files
1281  IF( trim(files(i)%name) == trim(name) .AND. tile_count == files(i)%tile_count ) THEN
1282  find_file = i
1283  RETURN
1284  END IF
1285  END DO
1286  END FUNCTION find_file
1287 
1288  !> @brief Return the field number for the given module name, field name, and tile number.
1289  !! @return Integer find_input_field
1290  INTEGER FUNCTION find_input_field(module_name, field_name, tile_count)
1291  CHARACTER(len=*), INTENT(in) :: module_name !< Module name.
1292  CHARACTER(len=*), INTENT(in) :: field_name !< field name.
1293  INTEGER, INTENT(in) :: tile_count !< Tile number.
1294 
1295  INTEGER :: i
1296 
1297  find_input_field = diag_field_not_found ! Default return value if not found.
1298  DO i = 1, num_input_fields
1299  IF(tile_count == input_fields(i)%tile_count .AND.&
1300  & trim(input_fields(i)%module_name) == trim(module_name) .AND.&
1301  & lowercase(trim(input_fields(i)%field_name)) == lowercase(trim(field_name))) THEN
1302  find_input_field = i
1303  RETURN
1304  END IF
1305  END DO
1306  END FUNCTION find_input_field
1307 
1308  !> @brief Initialize the input field.
1309  SUBROUTINE init_input_field(module_name, field_name, tile_count)
1310  CHARACTER(len=*), INTENT(in) :: module_name !< Module name.
1311  CHARACTER(len=*), INTENT(in) :: field_name !< Input field name.
1312  INTEGER, INTENT(in) :: tile_count !< Tile number.
1313 
1314  ! Get a number for this input_field if not already set up
1315  IF ( find_input_field(module_name, field_name, tile_count) < 0 ) THEN
1316  num_input_fields = num_input_fields + 1
1317  IF ( num_input_fields > max_input_fields ) THEN
1318  ! <ERROR STATUS="FATAL">max_input_fields exceeded, increase it via diag_manager_nml</ERROR>
1319  CALL error_mesg('diag_util_mod::init_input_field',&
1320  & 'max_input_fields exceeded, increase it via diag_manager_nml', fatal)
1321  END IF
1322  ELSE
1323  ! If this is already initialized do not need to do anything
1324  RETURN
1325  END IF
1326 
1327  input_fields(num_input_fields)%module_name = trim(module_name)
1328  input_fields(num_input_fields)%field_name = trim(field_name)
1329  input_fields(num_input_fields)%num_output_fields = 0
1330  ! Set flag that this field has not been registered
1331  input_fields(num_input_fields)%register = .false.
1332  input_fields(num_input_fields)%local = .false.
1333  input_fields(num_input_fields)%standard_name = 'none'
1334  input_fields(num_input_fields)%tile_count = tile_count
1335  input_fields(num_input_fields)%numthreads = 1
1336  input_fields(num_input_fields)%active_omp_level = 0
1337  input_fields(num_input_fields)%time = time_zero
1338  END SUBROUTINE init_input_field
1339 
1340  !> @brief Initialize the output field.
1341  SUBROUTINE init_output_field(module_name, field_name, output_name, output_file,&
1342  & time_method, pack, tile_count, local_coord)
1343  CHARACTER(len=*), INTENT(in) :: module_name !< Module name.
1344  CHARACTER(len=*), INTENT(in) :: field_name !< Output field name.
1345  CHARACTER(len=*), INTENT(in) :: output_name !< Output name written to file.
1346  CHARACTER(len=*), INTENT(in) :: output_file !< File where field should be written.
1347  CHARACTER(len=*), INTENT(in) :: time_method !< Data reduction method.
1348  !! See <LINK SRC="diag_manager.html">diag_manager_mod</LINK>
1349  !! for valid methods.</IN>
1350  INTEGER, INTENT(in) :: pack !< Packing method.
1351  INTEGER, INTENT(in) :: tile_count !< Tile number.
1352  TYPE(coord_type), INTENT(in), OPTIONAL :: local_coord !< Region to be written.
1353  !! If missing, then all data to be written.</IN>
1354  INTEGER :: out_num, in_num, file_num, file_num_tile1
1355  INTEGER :: num_fields, i, method_selected, l1
1356  INTEGER :: ioerror
1357  REAL :: pow_value
1358  INTEGER :: grv !< Value used to determine if the region defined in the diag_table is for the whole
1359  !! axis, or a sub-axis
1360  CHARACTER(len=128) :: error_msg
1361  CHARACTER(len=50) :: t_method
1362  character(len=256) :: tmp_name
1363 
1364  ! Value to use to determine if a region is to be output on the full axis, or sub-axis
1365  ! get the value to compare to determine if writing full axis data
1366  IF ( region_out_use_alt_value ) THEN
1367  grv = glo_reg_val_alt
1368  ELSE
1369  grv = glo_reg_val
1370  END IF
1371 
1372 
1373  ! Get a number for this output field
1374  num_output_fields = num_output_fields + 1
1375  IF ( num_output_fields > max_output_fields ) THEN
1376  ! <ERROR STATUS="FATAL">max_output_fields = <max_output_fields> exceeded. Increase via diag_manager_nml</ERROR>
1377  WRITE (unit=error_msg,fmt=*) max_output_fields
1378  CALL error_mesg('diag_util_mod::init_output_field', 'max_output_fields = '//trim(error_msg)//' exceeded.&
1379  & Increase via diag_manager_nml', fatal)
1380  END IF
1381  out_num = num_output_fields
1382 
1383  ! First, find the index to the associated input field
1384  in_num = find_input_field(module_name, field_name, tile_count)
1385  IF ( in_num < 0 ) THEN
1386  IF ( tile_count > 1 ) THEN
1387  WRITE (error_msg,'(A,"/",A,"/",A)') trim(module_name),trim(field_name),&
1388  & "tile_count="//trim(string(tile_count))
1389  ELSE
1390  WRITE (error_msg,'(A,"/",A)') trim(module_name),trim(field_name)
1391  END IF
1392  ! <ERROR STATUS="FATAL">module_name/field_name <module_name>/<field_name>[/tile_count=<tile_count>]
1393  ! NOT registered</ERROR>
1394  CALL error_mesg('diag_util_mod::init_output_field',&
1395  & 'module_name/field_name '//trim(error_msg)//' NOT registered', fatal)
1396  END IF
1397 
1398  ! Add this output field into the list for this input field
1399  input_fields(in_num)%num_output_fields =&
1400  & input_fields(in_num)%num_output_fields + 1
1401  IF ( input_fields(in_num)%num_output_fields > max_out_per_in_field ) THEN
1402  ! <ERROR STATUS="FATAL">
1403  ! MAX_OUT_PER_IN_FIELD = <MAX_OUT_PER_IN_FIELD> exceeded for <module_name>/<field_name>,
1404  ! increase MAX_OUT_PER_IN_FIELD
1405  ! in the diag_manager_nml namelist.
1406  ! </ERROR>
1407  WRITE (unit=error_msg,fmt=*) max_out_per_in_field
1408  CALL error_mesg('diag_util_mod::init_output_field',&
1409  & 'MAX_OUT_PER_IN_FIELD exceeded for '//trim(module_name)//"/"//trim(field_name)//&
1410  &', increase MAX_OUT_PER_IN_FIELD in the diag_manager_nml namelist', fatal)
1411  END IF
1412  input_fields(in_num)%output_fields(input_fields(in_num)%num_output_fields) = out_num
1413 
1414  ! Also put pointer to input field in this output field
1415  output_fields(out_num)%input_field = in_num
1416 
1417  ! Next, find the number for the corresponding file
1418  IF ( trim(output_file).EQ.'null' ) THEN
1419  file_num = max_files
1420  ELSE
1421  file_num = find_file(output_file, 1)
1422  IF ( file_num < 0 ) THEN
1423  ! <ERROR STATUS="FATAL">
1424  ! file <file_name> is NOT found in the diag_table.
1425  ! </ERROR>
1426  CALL error_mesg('diag_util_mod::init_output_field', 'file '&
1427  & //trim(output_file)//' is NOT found in the diag_table', fatal)
1428  END IF
1429  IF ( tile_count > 1 ) THEN
1430  file_num_tile1 = file_num
1431  file_num = find_file(output_file, tile_count)
1432  IF(file_num < 0) THEN
1433  CALL init_file(files(file_num_tile1)%name, files(file_num_tile1)%output_freq,&
1434  & files(file_num_tile1)%output_units, files(file_num_tile1)%format,&
1435  & files(file_num_tile1)%time_units, files(file_num_tile1)%long_name,&
1436  & tile_count, files(file_num_tile1)%new_file_freq,&
1437  & files(file_num_tile1)%new_file_freq_units, files(file_num_tile1)%start_time,&
1438  & files(file_num_tile1)%duration, files(file_num_tile1)%duration_units )
1439  file_num = find_file(output_file, tile_count)
1440  IF ( file_num < 0 ) THEN
1441  ! <ERROR STATUS="FATAL">
1442  ! file <output_file> is not initialized for tile_count = <tile_count>
1443  ! </ERROR>
1444  CALL error_mesg('diag_util_mod::init_output_field', 'file '//trim(output_file)//&
1445  & ' is not initialized for tile_count = '//trim(string(tile_count)), fatal)
1446  END IF
1447  END IF
1448  END IF
1449  END IF
1450 
1451  ! Insert this field into list for this file
1452  files(file_num)%num_fields = files(file_num)%num_fields + 1
1453  IF ( files(file_num)%num_fields > max_fields_per_file ) THEN
1454  WRITE (unit=error_msg, fmt=*) max_fields_per_file
1455  ! <ERROR STATUS="FATAL">
1456  ! MAX_FIELDS_PER_FILE = <MAX_FIELDS_PER_FILE> exceeded. Increase MAX_FIELDS_PER_FILE in diag_data.F90.
1457  ! </ERROR>
1458  CALL error_mesg('diag_util_mod::init_output_field',&
1459  & 'MAX_FIELDS_PER_FILE = '//trim(error_msg)// &
1460  & ' exceeded. Increase MAX_FIELDS_PER_FILE in diag_data.F90.', fatal)
1461  END IF
1462  num_fields = files(file_num)%num_fields
1463  files(file_num)%fields(num_fields) = out_num
1464 
1465  ! Set the file for this output field
1466  output_fields(out_num)%output_file = file_num
1467 
1468  ! Enter the other data for this output field
1469  output_fields(out_num)%output_name = trim(output_name)
1470  output_fields(out_num)%pack = pack
1471  output_fields(out_num)%pow_value = 1
1472  output_fields(out_num)%num_axes = 0
1473  output_fields(out_num)%total_elements = 0
1474  output_fields(out_num)%region_elements = 0
1475 
1476  call output_fields(out_num)%buff_bounds%reset(very_large_axis_length, 0)
1477 
1478  ! initialize the size of the diurnal axis to 1
1479  output_fields(out_num)%n_diurnal_samples = 1
1480 
1481  ! Initialize all time method to false
1482  method_selected = 0
1483  output_fields(out_num)%time_average = .false.
1484  output_fields(out_num)%time_rms = .false.
1485  output_fields(out_num)%time_min = .false.
1486  output_fields(out_num)%time_max = .false.
1487  output_fields(out_num)%time_sum = .false.
1488  output_fields(out_num)%time_ops = .false.
1489  output_fields(out_num)%written_once = .false.
1490 
1491  t_method = lowercase(time_method)
1492  ! cannot time average fields output every time
1493  IF ( files(file_num)%output_freq == every_time ) THEN
1494  output_fields(out_num)%time_average = .false.
1495  method_selected = method_selected+1
1496  t_method = 'point'
1497  ELSEIF ( index(t_method,'diurnal') == 1 ) THEN
1498  ! get the integer number from the t_method
1499  READ (unit=t_method(8:len_trim(t_method)), fmt=*, iostat=ioerror) output_fields(out_num)%n_diurnal_samples
1500  IF ( ioerror /= 0 ) THEN
1501  ! <ERROR STATUS="FATAL">
1502  ! could not find integer number of diurnal samples in string "<t_method>"
1503  ! </ERROR>
1504  CALL error_mesg('diag_util_mod::init_output_field',&
1505  & 'could not find integer number of diurnal samples in string "' //trim(t_method)//'"', fatal)
1506  ELSE IF ( output_fields(out_num)%n_diurnal_samples <= 0 ) THEN
1507  ! <ERROR STATUS="FATAL">
1508  ! The integer value of diurnal samples must be greater than zero.
1509  ! </ERROR>
1510  CALL error_mesg('diag_util_mod::init_output_field',&
1511  & 'The integer value of diurnal samples must be greater than zero.', fatal)
1512  END IF
1513  output_fields(out_num)%time_average = .true.
1514  method_selected = method_selected+1
1515  t_method='mean'
1516  ELSEIF ( index(t_method,'pow') == 1 ) THEN
1517  ! Get the integer number from the t_method
1518  READ (unit=t_method(4:len_trim(t_method)), fmt=*, iostat=ioerror) pow_value
1519  IF ( ioerror /= 0 .OR. output_fields(out_num)%pow_value < 1 .OR. floor(pow_value) /= ceiling(pow_value) ) THEN
1520  ! <ERROR STATUS="FATAL">
1521  ! Invalid power number in time operation "<t_method>". Must be a positive integer.
1522  ! </ERROR>
1523  CALL error_mesg('diag_util_mod::init_output_field',&
1524  & 'Invalid power number in time operation "'//trim(t_method)//'". Must be a positive integer', fatal)
1525  END IF
1526  output_fields(out_num)%pow_value = int(pow_value)
1527  output_fields(out_num)%time_average = .true.
1528  method_selected = method_selected+1
1529  t_method = 'mean_pow('//t_method(4:len_trim(t_method))//')'
1530  ELSE
1531  SELECT CASE(trim(t_method))
1532  CASE ( '.true.', 'mean', 'average', 'avg' )
1533  output_fields(out_num)%time_average = .true.
1534  method_selected = method_selected+1
1535  t_method = 'mean'
1536  CASE ( 'rms' )
1537  output_fields(out_num)%time_average = .true.
1538  output_fields(out_num)%time_rms = .true.
1539  output_fields(out_num)%pow_value = 2.0
1540  method_selected = method_selected+1
1541  t_method = 'root_mean_square'
1542  CASE ( '.false.', 'none', 'point' )
1543  output_fields(out_num)%time_average = .false.
1544  method_selected = method_selected+1
1545  t_method = 'point'
1546  CASE ( 'maximum', 'max' )
1547  output_fields(out_num)%time_max = .true.
1548  l1 = len_trim(output_fields(out_num)%output_name)
1549  if (l1 .ge. 3) then
1550  tmp_name = trim(adjustl(output_fields(out_num)%output_name(l1-2:l1)))
1551  IF (lowercase(trim(tmp_name)) /= 'max' ) then
1552  output_fields(out_num)%output_name = trim(output_name)//'_max'
1553  endif
1554  endif
1555  method_selected = method_selected+1
1556  t_method = 'max'
1557  CASE ( 'minimum', 'min' )
1558  output_fields(out_num)%time_min = .true.
1559  l1 = len_trim(output_fields(out_num)%output_name)
1560  if (l1 .ge. 3) then
1561  tmp_name = trim(adjustl(output_fields(out_num)%output_name(l1-2:l1)))
1562  IF (lowercase(trim(tmp_name)) /= 'min' ) then
1563  output_fields(out_num)%output_name = trim(output_name)//'_min'
1564  endif
1565  endif
1566  method_selected = method_selected+1
1567  t_method = 'min'
1568  CASE ( 'sum', 'cumsum' )
1569  output_fields(out_num)%time_sum = .true.
1570  l1 = len_trim(output_fields(out_num)%output_name)
1571  IF ( output_fields(out_num)%output_name(l1-2:l1) /= 'sum' )&
1572  & output_fields(out_num)%output_name = trim(output_name)//'_sum'
1573  method_selected = method_selected+1
1574  t_method = 'sum'
1575  END SELECT
1576  END IF
1577 
1578  ! reconcile logical flags
1579  output_fields(out_num)%time_ops = output_fields(out_num)%time_min.OR.output_fields(out_num)%time_max&
1580  & .OR.output_fields(out_num)%time_average .OR. output_fields(out_num)%time_sum
1581 
1582  output_fields(out_num)%phys_window = .false.
1583  ! need to initialize grid_type = -1(start, end, l_start_indx,l_end_indx etc...)
1584  IF ( PRESENT(local_coord) ) THEN
1585  input_fields(in_num)%local = .true.
1586  input_fields(in_num)%local_coord = local_coord
1587  IF ( int(local_coord%xbegin) == grv .AND. int(local_coord%xend) == grv .AND.&
1588  & int(local_coord%ybegin) == grv .AND. int(local_coord%yend) == grv ) THEN
1589  output_fields(out_num)%local_output = .false.
1590  output_fields(out_num)%need_compute = .false.
1591  output_fields(out_num)%reduced_k_range = .true.
1592  ELSE
1593  output_fields(out_num)%local_output = .true.
1594  output_fields(out_num)%need_compute = .false.
1595  output_fields(out_num)%reduced_k_range = .false.
1596  END IF
1597 
1598  output_fields(out_num)%output_grid%start(1) = local_coord%xbegin
1599  output_fields(out_num)%output_grid%start(2) = local_coord%ybegin
1600  output_fields(out_num)%output_grid%start(3) = local_coord%zbegin
1601  output_fields(out_num)%output_grid%end(1) = local_coord%xend
1602  output_fields(out_num)%output_grid%end(2) = local_coord%yend
1603  output_fields(out_num)%output_grid%end(3) = local_coord%zend
1604  DO i = 1, 3
1605  output_fields(out_num)%output_grid%l_start_indx(i) = -1
1606  output_fields(out_num)%output_grid%l_end_indx(i) = -1
1607  output_fields(out_num)%output_grid%subaxes(i) = -1
1608  END DO
1609  ELSE
1610  output_fields(out_num)%local_output = .false.
1611  output_fields(out_num)%need_compute = .false.
1612  output_fields(out_num)%reduced_k_range = .false.
1613  END IF
1614 
1615  ! <ERROR STATUS="FATAL">
1616  ! improper time method in diag_table for output field <output_name>
1617  ! </ERROR>
1618  IF ( method_selected /= 1 ) CALL error_mesg('diag_util_mod::init_output_field',&
1619  &'improper time method in diag_table for output field:'//trim(output_name),fatal)
1620 
1621  output_fields(out_num)%time_method = trim(t_method)
1622 
1623  ! allocate counters: NOTE that for simplicity we always allocate them, even
1624  ! if they are superceeded by 4D "counter" array. This isn't most memory
1625  ! efficient, approach, but probably tolerable since they are so small anyway
1626  ALLOCATE(output_fields(out_num)%count_0d(output_fields(out_num)%n_diurnal_samples))
1627  ALLOCATE(output_fields(out_num)%num_elements(output_fields(out_num)%n_diurnal_samples))
1628  output_fields(out_num)%count_0d(:) = 0
1629  output_fields(out_num)%num_elements(:) = 0
1630  output_fields(out_num)%num_attributes = 0
1631  END SUBROUTINE init_output_field
1632 
1633  !> @brief Open file for output, and write the meta data.
1634  SUBROUTINE opening_file(file, time, filename_time)
1635  ! WARNING: Assumes that all data structures are fully initialized
1636  INTEGER, INTENT(in) :: file !< File ID.
1637  TYPE(time_type), INTENT(in) :: time !< Time for the file time stamp
1638  TYPE(time_type), INTENT(in), optional :: filename_time !< Time used in setting the filename when
1639  !! writing periodic files
1640 
1641  TYPE(time_type) :: fname_time !< Time used in setting the filename when writing periodic files
1642  REAL, DIMENSION(2) :: open_file_data
1643  INTEGER :: j, field_num, input_field_num, num_axes, k
1644  INTEGER :: field_num1
1645  INTEGER :: position
1646  INTEGER :: dir, edges
1647  INTEGER :: year, month, day, hour, minute, second
1648  INTEGER, DIMENSION(1) :: time_axis_id, time_bounds_id
1649  ! size of this axes array must be at least max num. of
1650  ! axes per field + 2; the last two elements are for time
1651  ! and time bounds dimensions
1652  INTEGER, DIMENSION(6) :: axes
1653  INTEGER, ALLOCATABLE :: axesc(:) ! indices if compressed axes associated with the field
1654  LOGICAL :: time_ops, aux_present, match_aux_name, req_present, match_req_fields
1655  CHARACTER(len=7) :: avg_name = 'average'
1656  CHARACTER(len=128) :: time_units, timeb_units, avg, error_string, aux_name, req_fields, fieldname
1657  CHARACTER(len=FMS_FILE_LEN) :: filename
1658  CHARACTER(len=128) :: suffix, base_name
1659  CHARACTER(len=32) :: time_name, timeb_name,time_longname, timeb_longname, cart_name
1660  CHARACTER(len=FMS_FILE_LEN) :: fname
1661  CHARACTER(len=24) :: start_date
1662  TYPE(domain1d) :: domain
1663  TYPE(domain2d) :: domain2
1664  TYPE(domainug) :: domainU
1665  INTEGER :: is, ie, last, ind
1666  class(fmsnetcdffile_t), pointer :: fileob
1667  integer :: actual_num_axes !< The actual number of axes to write including time
1668 
1669  aux_present = .false.
1670  match_aux_name = .false.
1671  req_present = .false.
1672  match_req_fields = .false.
1673 
1674  ! Here is where time_units string must be set up; time since base date
1675  WRITE (time_units, 11) trim(time_unit_list(files(file)%time_units)), get_base_year(),&
1676  & get_base_month(), get_base_day(), get_base_hour(), get_base_minute(), get_base_second()
1677 11 FORMAT(a, ' since ', i4.4, '-', i2.2, '-', i2.2, ' ', i2.2, ':', i2.2, ':', i2.2)
1678  base_name = files(file)%name
1679  IF ( files(file)%new_file_freq < very_large_file_freq ) THEN
1680  position = index(files(file)%name, '%')
1681  IF ( position > 0 ) THEN
1682  base_name = base_name(1:position-1)
1683  ELSE
1684  ! <ERROR STATUS="FATAL">
1685  ! filename <files(file)%name> does not contain % for time stamp string
1686  ! </ERROR>
1687  CALL error_mesg('diag_util_mod::opening_file',&
1688  & 'file name '//trim(files(file)%name)//' does not contain % for time stamp string', fatal)
1689  END IF
1690  if (present(filename_time)) then
1691  fname_time = filename_time
1692  else
1693  fname_time = time
1694  endif
1695  suffix = get_time_string(files(file)%name, fname_time)
1696  ELSE
1697  suffix = ' '
1698  END IF
1699 
1700  ! Add ensemble ID to filename
1701  fname=base_name
1702  call get_instance_filename(fname, base_name)
1703 
1704  ! Set the filename
1705  filename = trim(base_name)//trim(suffix)
1706 
1707  ! prepend the file start date if prepend_date == .TRUE.
1708  IF ( prepend_date ) THEN
1709  call get_date(diag_init_time, year, month, day, hour, minute, second)
1710  write (start_date, '(1I20.4, 2I2.2)') year, month, day
1711 
1712  filename = trim(adjustl(start_date))//'.'//trim(filename)
1713  END IF
1714 
1715  ! Loop through all fields with this file to output axes
1716  ! JWD: This is a klooge; need something more robust
1717  domain2 = null_domain2d
1718  domainu = null_domainug
1719  DO j = 1, files(file)%num_fields
1720  field_num = files(file)%fields(j)
1721  if (output_fields(field_num)%local_output .AND. .NOT. output_fields(field_num)%need_compute) cycle
1722  num_axes = output_fields(field_num)%num_axes
1723  IF ( num_axes > 1 ) THEN
1724  domain2 = get_domain2d( output_fields(field_num)%axes(1:num_axes) )
1725  domainu = get_domainug( output_fields(field_num)%axes(1) )
1726  IF ( domain2 .NE. null_domain2d ) EXIT
1727  ELSEIF (num_axes == 1) THEN
1728  if (domainu .EQ. null_domainug) then
1729  domainu = get_domainug( output_fields(field_num)%axes(num_axes) )
1730  endif
1731  END IF
1732  END DO
1733 
1734  IF (domainu .NE. null_domainug .AND. domain2 .NE. null_domain2d) THEN
1735  CALL error_mesg('diag_util_mod::opening_file', &
1736  'Domain2 and DomainU are somehow both set.', &
1737  fatal)
1738  ENDIF
1739 
1740  IF ( allocated(files(file)%attributes) ) THEN
1741  CALL diag_output_init(filename, global_descriptor,&
1742  & files(file)%file_unit, domain2, domainu,&
1743  & fileobj(file),fileobju(file), fileobjnd(file), fnum_for_domain(file),&
1744  & attributes=files(file)%attributes(1:files(file)%num_attributes))
1745  ELSE
1746  CALL diag_output_init(filename, global_descriptor,&
1747  & files(file)%file_unit, domain2,domainu, &
1748  & fileobj(file),fileobju(file),fileobjnd(file),fnum_for_domain(file))
1749  END IF
1750  !> update fnum_for_domain with the correct domain
1751  files(file)%bytes_written = 0
1752  ! Does this file contain time_average fields?
1753  time_ops = .false.
1754  DO j = 1, files(file)%num_fields
1755  field_num = files(file)%fields(j)
1756  IF ( output_fields(field_num)%time_ops ) THEN
1757  time_ops = .true.
1758  EXIT
1759  END IF
1760  END DO
1761  ! Loop through all fields with this file to output axes
1762  DO j = 1, files(file)%num_fields
1763  field_num = files(file)%fields(j)
1764  input_field_num = output_fields(field_num)%input_field
1765  IF (.NOT.input_fields(input_field_num)%register) THEN
1766  WRITE (error_string,'(A,"/",A)') trim(input_fields(input_field_num)%module_name),&
1767  & trim(input_fields(input_field_num)%field_name)
1768  IF(mpp_pe() .EQ. mpp_root_pe()) THEN
1769  CALL error_mesg('diag_util_mod::opening_file',&
1770  & 'module/field_name ('//trim(error_string)//') NOT registered', warning)
1771  END IF
1772  cycle
1773  END IF
1774  if (output_fields(field_num)%local_output .AND. .NOT. output_fields(field_num)%need_compute) cycle
1775 
1776  ! Put the time axis in the axis field
1777  num_axes = output_fields(field_num)%num_axes
1778  axes(1:num_axes) = output_fields(field_num)%axes(1:num_axes)
1779  ! make sure that axis_id are not -1
1780  DO k = 1, num_axes
1781  IF ( axes(k) < 0 ) THEN
1782  WRITE(error_string,'(a)') output_fields(field_num)%output_name
1783  ! <ERROR STATUS="FATAL">
1784  ! ouptut_name <output_fields(field_num)%output_name> has axis_id = -1
1785  ! </ERROR>
1786  CALL error_mesg('diag_util_mod::opening_file','output_name '//trim(error_string)//&
1787  & ' has axis_id = -1', fatal)
1788  END IF
1789  END DO
1790  ! check if aux is present in any axes
1791  IF ( .NOT.aux_present ) THEN
1792  DO k = 1, num_axes
1793  aux_name = get_axis_aux(axes(k))
1794  IF ( trim(aux_name) /= 'none' ) THEN
1795  aux_present = .true.
1796  EXIT
1797  END IF
1798  END DO
1799  END IF
1800  ! check if required fields are present in any axes
1801  IF ( .NOT.req_present ) THEN
1802  DO k = 1, num_axes
1803  req_fields = get_axis_reqfld(axes(k))
1804  IF ( trim(req_fields) /= 'none' ) THEN
1805  CALL error_mesg('diag_util_mod::opening_file','required fields found: '//&
1806  &trim(req_fields)//' in file '//trim(files(file)%name),note)
1807  req_present = .true.
1808  EXIT
1809  END IF
1810  END DO
1811  END IF
1812 
1813  axes(num_axes + 1) = files(file)%time_axis_id
1814 !> Allocate the is_time_axis_registered field and set it to false for the first trip
1815  if (.not. allocated(files(file)%is_time_axis_registered)) then
1816  allocate(files(file)%is_time_axis_registered)
1817  files(file)%is_time_axis_registered = .false.
1818  endif
1819  if (time_ops) then
1820  !< If the file contains time_average fields write the "time" and "nv" dimension
1821  actual_num_axes = num_axes + 2
1822  axes(num_axes + 2) = files(file)%time_bounds_id
1823  else
1824  !< If the file doesn't contain time_average fields write the "time" dimension
1825  actual_num_axes = num_axes + 1
1826  endif
1827 
1828  if (fnum_for_domain(file) == "2d") then
1829  CALL write_axis_meta_data(files(file)%file_unit, axes(1:actual_num_axes),fileobj(file), time_ops=time_ops, &
1830  time_axis_registered=files(file)%is_time_axis_registered)
1831  elseif (fnum_for_domain(file) == "nd") then
1832  CALL write_axis_meta_data(files(file)%file_unit, axes(1:actual_num_axes),fileobjnd(file), time_ops=time_ops,&
1833  time_axis_registered=files(file)%is_time_axis_registered)
1834  elseif (fnum_for_domain(file) == "ug") then
1835  CALL write_axis_meta_data(files(file)%file_unit, axes(1:actual_num_axes),fileobju(file), time_ops=time_ops, &
1836  time_axis_registered=files(file)%is_time_axis_registered)
1837  endif
1838 
1839  ! write metadata for axes used in compression-by-gathering, e.g. for unstructured
1840  ! grid
1841  DO k = 1, num_axes
1842  IF (axis_is_compressed(axes(k))) THEN
1843  CALL get_compressed_axes_ids(axes(k), axesc) ! returns allocatable array
1844  if (fnum_for_domain(file) == "ug") then
1845  CALL write_axis_meta_data(files(file)%file_unit, axesc,fileobju(file), &
1846  time_axis_registered=files(file)%is_time_axis_registered)
1847  else
1848  CALL error_mesg('diag_util_mod::opening_file::'//trim(filename), "Compressed "//&
1849  "dimensions are only allowed with axis in the unstructured dimension", fatal)
1850  endif
1851  DEALLOCATE(axesc)
1852  ENDIF
1853  ENDDO
1854  END DO
1855 
1856  ! Looking for the first NON-static field in a file
1857  field_num1 = files(file)%fields(1)
1858  DO j = 1, files(file)%num_fields
1859  field_num = files(file)%fields(j)
1860  IF ( output_fields(field_num)%time_ops ) THEN
1861  field_num1 = field_num
1862  EXIT
1863  END IF
1864  END DO
1865  nfields_loop: DO j = 1, files(file)%num_fields
1866  field_num = files(file)%fields(j)
1867  input_field_num = output_fields(field_num)%input_field
1868  IF (.NOT.input_fields(input_field_num)%register) cycle
1869  IF (output_fields(field_num)%local_output .AND. .NOT. output_fields(field_num)%need_compute) cycle
1870  ! Make sure that 1 file contains either time_average or instantaneous fields
1871  ! cannot have both time_average and instantaneous in 1 file
1872  IF ( .NOT.mix_snapshot_average_fields ) THEN
1873  IF ( (output_fields(field_num)%time_ops.NEQV.output_fields(field_num1)%time_ops) .AND.&
1874  & .NOT.output_fields(field_num1)%static .AND. .NOT.output_fields(field_num)%static) THEN
1875  IF ( mpp_pe() == mpp_root_pe() ) THEN
1876  ! <ERROR STATUS="FATAL">
1877  ! <files(file)%name> can NOT have BOTH time average AND instantaneous fields.
1878  ! Create a new file or set mix_snapshot_average_fields=.TRUE. in the namelist diag_manager_nml.
1879  ! </ERROR>
1880  CALL error_mesg('diag_util_mod::opening_file','file '//trim(files(file)%name)// &
1881  &' can NOT have BOTH time average AND instantaneous fields.'//&
1882  &' Create a new file or set mix_snapshot_average_fields=.TRUE. in the namelist diag_manager_nml.'&
1883  &, fatal)
1884  END IF
1885  END IF
1886  END IF
1887  ! check if any field has the same name as aux_name
1888  IF ( aux_present .AND. .NOT.match_aux_name ) THEN
1889  fieldname = output_fields(field_num)%output_name
1890  IF ( index(aux_name, trim(fieldname)) > 0 ) match_aux_name = .true.
1891  END IF
1892  ! check if any field has the same name as req_fields
1893  IF ( req_present .AND. .NOT.match_req_fields ) THEN
1894  fieldname = output_fields(field_num)%output_name
1895  is = 1; last = len_trim(req_fields)
1896  DO
1897  ind = index(req_fields(is:last),' ')
1898  IF (ind .eq. 0) ind = last-is+2
1899  ie = is+(ind-2)
1900  if (req_fields(is:ie) .EQ. trim(fieldname)) then
1901  match_req_fields = .true.
1902  !CALL error_mesg('diag_util_mod::opening_file','matched required field: '//TRIM(fieldname),NOTE)
1903  EXIT
1904  END IF
1905  is = is+ind
1906  if (is .GT. last) EXIT
1907  END DO
1908  END IF
1909 
1910  ! Put the time axis in the axis field
1911  num_axes = output_fields(field_num)%num_axes
1912  axes(1:num_axes) = output_fields(field_num)%axes(1:num_axes)
1913  IF ( .NOT.output_fields(field_num)%static ) THEN
1914  num_axes=num_axes+1
1915  axes(num_axes) = files(file)%time_axis_id
1916  END IF
1917  IF(output_fields(field_num)%time_average) THEN
1918  avg = avg_name
1919  ELSE IF(output_fields(field_num)%time_max) THEN
1920  avg = avg_name
1921  ELSE IF(output_fields(field_num)%time_min) THEN
1922  avg = avg_name
1923  ELSE
1924  avg = " "
1925  END IF
1926 !> Use the correct file object
1927  if (fnum_for_domain(file) == "2d") then
1928  fileob => fileobj(file)
1929  elseif (fnum_for_domain(file) == "nd") then
1930  fileob => fileobjnd(file)
1931  elseif (fnum_for_domain(file) == "ug") then
1932  fileob => fileobju(file)
1933  endif
1934  IF ( input_fields(input_field_num)%missing_value_present ) THEN
1935  IF ( len_trim(input_fields(input_field_num)%interp_method) > 0 ) THEN
1936  output_fields(field_num)%f_type = write_field_meta_data(files(file)%file_unit,&
1937  & output_fields(field_num)%output_name, axes(1:num_axes),&
1938  & input_fields(input_field_num)%units,&
1939  & input_fields(input_field_num)%long_name,&
1940  & input_fields(input_field_num)%range, output_fields(field_num)%pack,&
1941  & input_fields(input_field_num)%missing_value, avg_name = avg,&
1942  & time_method=output_fields(field_num)%time_method,&
1943  & standard_name = input_fields(input_field_num)%standard_name,&
1944  & interp_method = input_fields(input_field_num)%interp_method,&
1945  & attributes=output_fields(field_num)%attributes,&
1946  & num_attributes=output_fields(field_num)%num_attributes,&
1947  & use_ugdomain=files(file)%use_domainUG , &
1948  & fileob=fileob)
1949  ELSE
1950  output_fields(field_num)%f_type = write_field_meta_data(files(file)%file_unit,&
1951  & output_fields(field_num)%output_name, axes(1:num_axes),&
1952  & input_fields(input_field_num)%units,&
1953  & input_fields(input_field_num)%long_name,&
1954  & input_fields(input_field_num)%range, output_fields(field_num)%pack,&
1955  & input_fields(input_field_num)%missing_value, avg_name = avg,&
1956  & time_method=output_fields(field_num)%time_method,&
1957  & standard_name = input_fields(input_field_num)%standard_name,&
1958  & attributes=output_fields(field_num)%attributes,&
1959  & num_attributes=output_fields(field_num)%num_attributes,&
1960  & use_ugdomain=files(file)%use_domainUG , &
1961  & fileob=fileob)
1962 
1963  END IF
1964  ! NEED TO TAKE CARE OF TIME AVERAGING INFO TOO BOTH CASES
1965  ELSE
1966  IF ( len_trim(input_fields(input_field_num)%interp_method) > 0 ) THEN
1967  output_fields(field_num)%f_type = write_field_meta_data(files(file)%file_unit,&
1968  & output_fields(field_num)%output_name, axes(1:num_axes),&
1969  & input_fields(input_field_num)%units,&
1970  & input_fields(input_field_num)%long_name,&
1971  & input_fields(input_field_num)%range, output_fields(field_num)%pack,&
1972  & avg_name = avg,&
1973  & time_method=output_fields(field_num)%time_method,&
1974  & standard_name = input_fields(input_field_num)%standard_name,&
1975  & interp_method = input_fields(input_field_num)%interp_method,&
1976  & attributes=output_fields(field_num)%attributes,&
1977  & num_attributes=output_fields(field_num)%num_attributes,&
1978  & use_ugdomain=files(file)%use_domainUG , &
1979  & fileob=fileob)
1980 
1981  ELSE
1982  output_fields(field_num)%f_type = write_field_meta_data(files(file)%file_unit,&
1983  & output_fields(field_num)%output_name, axes(1:num_axes),&
1984  & input_fields(input_field_num)%units,&
1985  & input_fields(input_field_num)%long_name,&
1986  & input_fields(input_field_num)%range, output_fields(field_num)%pack,&
1987  & avg_name = avg,&
1988  & time_method=output_fields(field_num)%time_method,&
1989  & standard_name = input_fields(input_field_num)%standard_name,&
1990  & attributes=output_fields(field_num)%attributes,&
1991  & num_attributes=output_fields(field_num)%num_attributes,&
1992  & use_ugdomain=files(file)%use_domainUG , &
1993  & fileob=fileob)
1994 
1995  END IF
1996  END IF
1997  END DO nfields_loop
1998  ! If any of the fields in the file are time averaged, need to output the axes
1999  ! Use double precision since time axis is double precision
2000  IF ( time_ops ) THEN
2001  time_axis_id(1) = files(file)%time_axis_id
2002  files(file)%f_avg_start = write_field_meta_data(files(file)%file_unit,&
2003  & avg_name // '_T1', time_axis_id, time_units,&
2004  & "Start time for average period", pack=pack_size , &
2005  & fileob=fileob)
2006  files(file)%f_avg_end = write_field_meta_data(files(file)%file_unit,&
2007  & avg_name // '_T2', time_axis_id, time_units,&
2008  & "End time for average period", pack=pack_size , &
2009  & fileob=fileob)
2010  files(file)%f_avg_nitems = write_field_meta_data(files(file)%file_unit,&
2011  & avg_name // '_DT', time_axis_id,&
2012  & trim(time_unit_list(files(file)%time_units)),&
2013  & "Length of average period", pack=pack_size , &
2014  & fileob=fileob)
2015  END IF
2016 
2017  IF ( time_ops ) THEN
2018  time_axis_id(1) = files(file)%time_axis_id
2019  time_bounds_id(1) = files(file)%time_bounds_id
2020  CALL get_diag_axis( time_axis_id(1), time_name, time_units, time_longname,&
2021  & cart_name, dir, edges, domain, domainu, open_file_data)
2022  CALL get_diag_axis( time_bounds_id(1), timeb_name, timeb_units, timeb_longname,&
2023  & cart_name, dir, edges, domain, domainu, open_file_data)
2024  ! CF Compliance requires the unit on the _bnds axis is the same as 'time'
2025  files(file)%f_bounds = write_field_meta_data(files(file)%file_unit,&
2026  & trim(time_name)//'_bnds', (/time_bounds_id,time_axis_id/),&
2027  & time_units, trim(time_name)//' axis boundaries', pack=pack_size , &
2028  & fileob=fileob)
2029  END IF
2030  ! Let lower levels know that all meta data has been sent
2031  call done_meta_data(files(file)%file_unit)
2032 
2033  IF( aux_present .AND. .NOT.match_aux_name ) THEN
2034  ! <ERROR STATUS="WARNING">
2035  ! one axis has auxiliary but the corresponding field is NOT
2036  ! found in file <file_name>
2037  ! </ERROR>
2038  IF ( mpp_pe() == mpp_root_pe() ) CALL error_mesg('diag_util_mod::opening_file',&
2039  &'one axis has auxiliary but the corresponding field is NOT found in file '// &
2040  & trim(files(file)%name), warning)
2041  END IF
2042  IF( req_present .AND. .NOT.match_req_fields ) THEN
2043  ! <ERROR STATUS="FATAL">
2044  ! one axis has required fields but the corresponding field is NOT
2045  ! found in file <file_name>
2046  ! </ERROR>
2047  IF ( mpp_pe() == mpp_root_pe() ) CALL error_mesg('diag_util_mod::opening_file',&
2048  &'one axis has required fields ('//trim(req_fields)//') but the '// &
2049  &'corresponding fields are NOT found in file '//trim(files(file)%name), fatal)
2050  END IF
2051  ! Clean up pointer
2052  if (associated(fileob)) nullify(fileob)
2053  END SUBROUTINE opening_file
2054 
2055  !> @brief Write data out to file, and if necessary flush the buffers.
2056  SUBROUTINE diag_data_out(file, field, dat, time, final_call_in, static_write_in, filename_time)
2057  INTEGER, INTENT(in) :: file !< File ID.
2058  INTEGER, INTENT(in) :: field !< Field ID.
2059  REAL, DIMENSION(:,:,:,:), INTENT(inout) :: dat !< Data to write out.
2060  TYPE(time_type), INTENT(in) :: time !< Current model time.
2061  LOGICAL, OPTIONAL, INTENT(in):: final_call_in !< <TT>.TRUE.</TT> if this is the last write for file.
2062  LOGICAL, OPTIONAL, INTENT(in):: static_write_in !< <TT>.TRUE.</TT> if static fields are to be written to file.
2063  type(time_type), intent(in), optional :: filename_time !< Time used in setting the filename when
2064  !! writing periodic files
2065 
2066  LOGICAL :: final_call, do_write, static_write
2067  REAL :: dif, time_data(2, 1, 1, 1), dt_time(1, 1, 1, 1), start_dif, end_dif
2068  REAL :: time_in_file !< Time in file at the beginning of this call
2069 
2070  !< Save the current time in the file. If the time in the file is not the same as the
2071  !! current time, files(file)%rtime_current will be updated
2072  time_in_file = files(file)%rtime_current
2073 
2074  do_write = .true.
2075  final_call = .false.
2076  IF ( PRESENT(final_call_in) ) final_call = final_call_in
2077  static_write = .false.
2078  IF ( PRESENT(static_write_in) ) static_write = static_write_in
2079 !> dif is the time as a real that is evaluated
2080  dif = get_date_dif(time, get_base_time(), files(file)%time_units)
2081 
2082  ! get file_unit, open new file and close curent file if necessary
2083  IF ( .NOT.static_write .OR. files(file)%file_unit < 0 ) &
2084  CALL check_and_open(file, time, do_write, filename_time=filename_time)
2085  IF ( .NOT.do_write ) RETURN ! no need to write data
2086 !> Set up the time index and write the correct time value to the time array
2087  if (dif > files(file)%rtime_current) then
2088  files(file)%time_index = files(file)%time_index + 1
2089  files(file)%rtime_current = dif
2090  if (fnum_for_domain(file) == "2d") then
2091  call diag_write_time (fileobj(file), files(file)%rtime_current, files(file)%time_index, &
2092  time_name=fileobj(file)%time_name)
2093  elseif (fnum_for_domain(file) == "ug") then
2094  call diag_write_time (fileobju(file), files(file)%rtime_current, files(file)%time_index, &
2095  time_name=fileobju(file)%time_name)
2096  elseif (fnum_for_domain(file) == "nd") then
2097  call diag_write_time (fileobjnd(file), files(file)%rtime_current, files(file)%time_index, &
2098  time_name=fileobjnd(file)%time_name)
2099  else
2100  call error_mesg("diag_util_mod::diag_data_out","Error opening the file "//files(file)%name,fatal)
2101  endif
2102  elseif (dif < files(file)%rtime_current .and. .not.(static_write) ) then
2103  call error_mesg("diag_util_mod::diag_data_out","The time for the file "//trim(files(file)%name)//&
2104  " has gone backwards. There may be missing values for some of the variables",note)
2105  endif
2106 !> Write data
2107  call diag_field_write (output_fields(field)%output_name, dat, static_write, file, fileobju, &
2108  fileobj, fileobjnd, fnum_for_domain(file), time_in=files(file)%time_index)
2109  ! record number of bytes written to this file
2110  files(file)%bytes_written = files(file)%bytes_written +&
2111  & (SIZE(dat,1)*SIZE(dat,2)*SIZE(dat,3))*(8/output_fields(field)%pack)
2112  IF ( .NOT.output_fields(field)%written_once ) output_fields(field)%written_once = .true.
2113  ! *** inserted this line because start_dif < 0 for static fields ***
2114  IF ( .NOT.output_fields(field)%static ) THEN
2115  start_dif = get_date_dif(output_fields(field)%last_output, get_base_time(),files(file)%time_units)
2116  IF ( .NOT.mix_snapshot_average_fields ) THEN
2117  end_dif = get_date_dif(output_fields(field)%next_output, get_base_time(), files(file)%time_units)
2118  ELSE
2119  end_dif = dif
2120  END IF
2121  END IF
2122 
2123  if (files(file)%rtime_current > time_in_file) then !< If time was written in this call
2124  if (output_fields(field)%time_ops) then !< If this is a time_average field
2125  ! Output the axes if this is first time-averaged field
2126  time_data(1, 1, 1, 1) = start_dif
2127  call diag_field_write (files(file)%f_avg_start%fieldname, time_data(1:1,:,:,:), static_write, file, &
2128  fileobju, fileobj, fileobjnd, &
2129  fnum_for_domain(file), time_in=files(file)%time_index)
2130  time_data(2, 1, 1, 1) = end_dif
2131  call diag_field_write (files(file)%f_avg_end%fieldname, time_data(2:2,:,:,:), static_write, file, &
2132  fileobju, fileobj, fileobjnd, &
2133  fnum_for_domain(file), time_in=files(file)%time_index)
2134  ! Compute the length of the average
2135  dt_time(1, 1, 1, 1) = end_dif - start_dif
2136  call diag_field_write (files(file)%f_avg_nitems%fieldname, dt_time(1:1,:,:,:), static_write, file, &
2137  fileobju, fileobj, fileobjnd, &
2138  fnum_for_domain(file), time_in=files(file)%time_index)
2139  ! Include boundary variable for CF compliance
2140  call diag_field_write (files(file)%f_bounds%fieldname, time_data(1:2,:,:,:), static_write, file, &
2141  fileobju, fileobj, fileobjnd, &
2142  fnum_for_domain(file), time_in=files(file)%time_index)
2143  END IF
2144  END IF
2145 
2146  ! If write time is greater (equal for the last call) than last_flush for this file, flush it
2147  IF ( final_call ) THEN
2148  IF ( time >= files(file)%last_flush ) THEN
2149  files(file)%last_flush = time
2150  END IF
2151  ELSE
2152  IF ( time > files(file)%last_flush .AND. (flush_nc_files.OR.debug_diag_manager) ) THEN
2153  call diag_flush(file, fileobju, fileobj, fileobjnd, fnum_for_domain(file))
2154  files(file)%last_flush = time
2155  END IF
2156  END IF
2157  END SUBROUTINE diag_data_out
2158 
2159  !> @brief Checks if it is time to open a new file.
2160  !! @details Checks if it is time to open a new file. If yes, it first closes the
2161  !! current file, opens a new file and returns file_unit
2162  !! previous diag_manager_end is replaced by closing_file and output_setup by opening_file.
2163  SUBROUTINE check_and_open(file, time, do_write, filename_time)
2164  INTEGER, INTENT(in) :: file !<File ID.
2165  TYPE(time_type), INTENT(in) :: time !< Current model time.
2166  LOGICAL, INTENT(out) :: do_write !< <TT>.TRUE.</TT> if file is expecting more data to write,
2167  !! <TT>.FALSE.</TT> otherwise.
2168  TYPE(time_type), INTENT(in), optional :: filename_time !< Time used in setting the filename when
2169  !! writing periodic files
2170 
2171  IF ( time >= files(file)%start_time ) THEN
2172  IF ( files(file)%file_unit < 0 ) THEN ! need to open a new file
2173  CALL opening_file(file, time, filename_time=filename_time)
2174  do_write = .true.
2175  ELSE
2176  do_write = .true.
2177  IF ( time > files(file)%close_time .AND. time < files(file)%next_open ) THEN
2178  do_write = .false. ! file still open but receives NO MORE data
2179  ELSE IF ( time > files(file)%next_open ) THEN ! need to close current file and open a new one
2180  CALL write_static(file) ! write all static fields and close this file
2181  CALL opening_file(file, time, filename_time=filename_time)
2182  files(file)%time_index = 0 !< Reset the number of times in the files back to 0
2183  files(file)%start_time = files(file)%next_open
2184  files(file)%close_time =&
2185  & diag_time_inc(files(file)%start_time,files(file)%duration, files(file)%duration_units)
2186  files(file)%next_open =&
2187  & diag_time_inc(files(file)%next_open, files(file)%new_file_freq,&
2188  & files(file)%new_file_freq_units)
2189  IF ( files(file)%close_time > files(file)%next_open ) THEN
2190  ! <ERROR STATUS="FATAL">
2191  ! <file_name> has close time GREATER than next_open time,
2192  ! check file duration and frequency
2193  ! </ERROR>
2194  CALL error_mesg('diag_util_mod::check_and_open',&
2195  & files(file)%name// &
2196  & ' has close time GREATER than next_open time, check file duration and frequency',fatal)
2197  END IF
2198  END IF ! no need to open new file, simply return file_unit
2199  END IF
2200  ELSE
2201  do_write = .false.
2202  END IF
2203  END SUBROUTINE check_and_open
2204 
2205  !> @brief Output all static fields in this file
2206  SUBROUTINE write_static(file)
2207  INTEGER, INTENT(in) :: file !< File ID.
2208  INTEGER :: j, i, input_num
2209 
2210  DO j = 1, files(file)%num_fields
2211  i = files(file)%fields(j)
2212  input_num = output_fields(i)%input_field
2213  ! skip fields that were not registered
2214  IF ( .NOT.input_fields(input_num)%register ) cycle
2215  IF ( output_fields(i)%local_output .AND. .NOT. output_fields(i)%need_compute) cycle
2216  ! only output static fields here
2217  IF ( .NOT.output_fields(i)%static ) cycle
2218  CALL diag_data_out(file, i, output_fields(i)%buffer, files(file)%last_flush, .true., .true.)
2219  END DO
2220 !! New FMS_IO close
2221  ! File is stil open. This is to protect when the diag_table has no Fields
2222  ! going to this file, and it was never opened (b/c diag_data_out was not
2223  ! called)
2224  if (fnum_for_domain(file) == "2d" )then
2225  if (check_if_open(fileobj(file))) call close_file (fileobj(file) )
2226  elseif (fnum_for_domain(file) == "nd") then
2227  if (check_if_open(fileobjnd(file)) ) then
2228  call close_file (fileobjnd(file))
2229  endif
2230  elseif (fnum_for_domain(file) == "ug") then
2231  if (check_if_open(fileobju(file))) call close_file (fileobju(file))
2232  endif
2233  files(file)%file_unit = -1
2234  END SUBROUTINE write_static
2235 
2236  !> @brief Checks to see if <TT>output_name</TT> and <TT>output_file</TT> are unique in <TT>output_fields</TT>.
2237  SUBROUTINE check_duplicate_output_fields(err_msg)
2238  CHARACTER(len=*), INTENT(out), OPTIONAL :: err_msg !< Error message. If empty, then no duplicates found.
2239 
2240  INTEGER :: i, j, tmp_file
2241  CHARACTER(len=128) :: tmp_name
2242  CHARACTER(len=256) :: err_msg_local
2243 
2244  IF ( PRESENT(err_msg) ) err_msg=''
2245  ! Do the checking when more than 1 output_fileds present
2246  IF ( num_output_fields <= 1 ) RETURN
2247  err_msg_local = ''
2248 
2249  i_loop: DO i = 1, num_output_fields-1
2250  tmp_name = trim(output_fields(i)%output_name)
2251  tmp_file = output_fields(i)%output_file
2252  DO j = i+1, num_output_fields
2253  IF ( (tmp_name == trim(output_fields(j)%output_name)) .AND. &
2254  &(tmp_file == output_fields(j)%output_file)) THEN
2255  err_msg_local = ' output_field "'//trim(tmp_name)//&
2256  &'" duplicated in file "'//trim(files(tmp_file)%name)//'"'
2257  EXIT i_loop
2258  END IF
2259  END DO
2260  END DO i_loop
2261  IF ( err_msg_local /= '' ) THEN
2262  IF ( fms_error_handler(' ERROR in diag_table',err_msg_local,err_msg) ) RETURN
2263  END IF
2264  END SUBROUTINE check_duplicate_output_fields
2265 
2266  !> @brief Allocates the atttype in out_field
2267  SUBROUTINE attribute_init_field(out_field, err_msg)
2268  TYPE(output_field_type), INTENT(inout) :: out_field !< output field to allocate memory for attribute
2269  CHARACTER(LEN=*), INTENT(out), OPTIONAL :: err_msg !< Error message, passed back to calling function
2270 
2271  INTEGER :: istat
2272 
2273  ! Need to initialize err_msg if present
2274  IF ( PRESENT(err_msg) ) err_msg = ''
2275 
2276  ! Allocate memory for the attributes
2277  IF ( .NOT.allocated(out_field%attributes) ) THEN
2278  ALLOCATE(out_field%attributes(max_field_attributes), stat=istat)
2279  IF ( istat.NE.0 ) THEN
2280  ! <ERROR STATUS="FATAL">
2281  ! Unable to allocate memory for attribute <name> to module/input_field <module_name>/<field_name>
2282  ! </ERROR>
2283  IF ( fms_error_handler('diag_util_mod::attribute_init_field',&
2284  & 'Unable to allocate memory for attributes', err_msg) ) THEN
2285  RETURN
2286  END IF
2287  ELSE
2288  ! Set equal to 0. It will be increased when attributes added
2289  out_field%num_attributes = 0
2290  END IF
2291  END IF
2292  END SUBROUTINE attribute_init_field
2293 
2294  !> @brief Prepends the attribute value to an already existing attribute. If the
2295  !! attribute isn't yet defined, then creates a new attribute
2296  SUBROUTINE prepend_attribute_field(out_field, att_name, prepend_value, err_msg)
2297  TYPE(output_field_type), INTENT(inout) :: out_field !< output field that will get the attribute
2298  CHARACTER(len=*), INTENT(in) :: att_name !< Name of the attribute
2299  CHARACTER(len=*), INTENT(in) :: prepend_value !< Value to prepend
2300  CHARACTER(len=*), INTENT(out) , OPTIONAL :: err_msg !< Error message, passed back to calling routine
2301 
2302  INTEGER :: length, i, this_attribute
2303  CHARACTER(len=512) :: err_msg_local
2304 
2305  ! Initialize string characters
2306  err_msg_local=''
2307  IF ( PRESENT(err_msg) ) err_msg = ''
2308 
2309  ! Make sure the attributes for this out field have been initialized
2310  CALL attribute_init_field(out_field, err_msg_local)
2311  IF ( trim(err_msg_local) .NE. '' ) THEN
2312  IF ( fms_error_handler('diag_util_mod::prepend_attribute_field', trim(err_msg_local), err_msg) ) THEN
2313  RETURN
2314  END IF
2315  END IF
2316 
2317  ! Find if attribute exists
2318  this_attribute = 0
2319  DO i=1, out_field%num_attributes
2320  IF ( trim(out_field%attributes(i)%name) .EQ. trim(att_name) ) THEN
2321  this_attribute = i
2322  EXIT
2323  END IF
2324  END DO
2325 
2326  IF ( this_attribute > 0 ) THEN
2327  IF ( out_field%attributes(this_attribute)%type .NE. nf90_char ) THEN
2328  ! <ERROR STATUS="FATAL">
2329  ! Attribute <name> is not a character attribute.
2330  ! </ERROR>
2331  IF ( fms_error_handler('diag_util_mod::prepend_attribute_field', &
2332  & 'Attribute "'//trim(att_name)//'" is not a character attribute.',&
2333  & err_msg) ) THEN
2334  RETURN
2335  END IF
2336  END IF
2337  ELSE
2338  ! Defining a new attribute
2339  ! Increase the number of field attributes
2340  this_attribute = out_field%num_attributes + 1
2341  IF ( this_attribute .GT. max_field_attributes ) THEN
2342  ! <ERROR STATUS="FATAL">
2343  ! Number of attributes exceeds max_field_attributes for attribute <name>.
2344  ! Increase diag_manager_nml:max_field_attributes.
2345  ! </ERROR>
2346  IF ( fms_error_handler('diag_util_mod::prepend_attribute_field',&
2347  & 'Number of attributes exceeds max_field_attributes for attribute "'&
2348  & //trim(att_name)//'". Increase diag_manager_nml:max_field_attributes.',&
2349  & err_msg) ) THEN
2350  RETURN
2351  END IF
2352  ELSE
2353  out_field%num_attributes = this_attribute
2354  ! Set name and type
2355  out_field%attributes(this_attribute)%name = att_name
2356  out_field%attributes(this_attribute)%type = nf90_char
2357  ! Initialize catt to a blank string, as len_trim doesn't always work on an uninitialized string
2358  out_field%attributes(this_attribute)%catt = ''
2359  END IF
2360  END IF
2361 
2362  ! Check if string is already included, and return if found
2363  IF ( index(trim(out_field%attributes(this_attribute)%catt), trim(prepend_value)).EQ.0 ) THEN
2364  ! Check if new string length goes beyond the length of catt
2365  length = len_trim(trim(prepend_value)//" "//trim(out_field%attributes(this_attribute)%catt))
2366  IF ( length.GT.len(out_field%attributes(this_attribute)%catt) ) THEN
2367  ! <ERROR STATUS="FATAL">
2368  ! Prepend length for attribute <name> is longer than allowed.
2369  ! </ERROR>
2370  IF ( fms_error_handler('diag_util_mod::prepend_attribute_field',&
2371  & 'Prepend length for attribute "'//trim(att_name)//'" is longer than allowed.',&
2372  & err_msg) ) THEN
2373  RETURN
2374  END IF
2375  END IF
2376  ! Set fields
2377  out_field%attributes(this_attribute)%catt =&
2378  & trim(prepend_value)//' '//trim(out_field%attributes(this_attribute)%catt)
2379  out_field%attributes(this_attribute)%len = length
2380  END IF
2381  END SUBROUTINE prepend_attribute_field
2382  !> @brief Allocates the atttype in out_file
2383  SUBROUTINE attribute_init_file(out_file, err_msg)
2384  TYPE(file_type), INTENT(inout) :: out_file !< output file to allocate memory for attribute
2385  CHARACTER(LEN=*), INTENT(out), OPTIONAL :: err_msg !< Error message, passed back to calling function
2386 
2387  INTEGER :: istat
2388 
2389  ! Initialize err_msg
2390  IF ( PRESENT(err_msg) ) err_msg = ''
2391 
2392  ! Allocate memory for the attributes
2393  IF ( .NOT.allocated(out_file%attributes) ) THEN
2394  ALLOCATE(out_file%attributes(max_field_attributes), stat=istat)
2395  IF ( istat.NE.0 ) THEN
2396  ! <ERROR STATUS="FATAL">
2397  ! Unable to allocate memory for file attributes
2398  ! </ERROR>
2399  IF ( fms_error_handler('diag_util_mod::attribute_init_file', &
2400  & 'Unable to allocate memory for file attributes', err_msg) ) THEN
2401  RETURN
2402  END IF
2403  ELSE
2404  ! Set equal to 0. It will be increased when attributes added
2405  out_file%num_attributes = 0
2406  END IF
2407  END IF
2408  END SUBROUTINE attribute_init_file
2409 
2410  !> @brief Prepends the attribute value to an already existing attribute. If the
2411  !! attribute isn't yet defined, then creates a new attribute
2412  SUBROUTINE prepend_attribute_file(out_file, att_name, prepend_value, err_msg)
2413  TYPE(file_type), INTENT(inout) :: out_file !< output file that will get the attribute
2414  CHARACTER(len=*), INTENT(in) :: att_name !< Name of the attribute
2415  CHARACTER(len=*), INTENT(in) :: prepend_value !< Value to prepend
2416  CHARACTER(len=*), INTENT(out) , OPTIONAL :: err_msg !< Error message, passed back to calling routine
2417 
2418  INTEGER :: length, i, this_attribute
2419  CHARACTER(len=512) :: err_msg_local
2420 
2421  ! Initialize string variables
2422  err_msg_local = ''
2423  IF ( PRESENT(err_msg) ) err_msg = ''
2424 
2425  ! Make sure the attributes for this out file have been initialized
2426  CALL attribute_init_file(out_file, err_msg_local)
2427  IF ( trim(err_msg_local) .NE. '' ) THEN
2428  IF ( fms_error_handler('diag_util_mod::prepend_attribute_file', trim(err_msg_local), err_msg) ) THEN
2429  RETURN
2430  END IF
2431  END IF
2432 
2433  ! Find if attribute exists
2434  this_attribute = 0
2435  DO i=1, out_file%num_attributes
2436  IF ( trim(out_file%attributes(i)%name) .EQ. trim(att_name) ) THEN
2437  this_attribute = i
2438  EXIT
2439  END IF
2440  END DO
2441 
2442  IF ( this_attribute > 0 ) THEN
2443  IF ( out_file%attributes(this_attribute)%type .NE. nf90_char ) THEN
2444  ! <ERROR STATUS="FATAL">
2445  ! Attribute <name> is not a character attribute.
2446  ! </ERROR>
2447  IF ( fms_error_handler('diag_util_mod::prepend_attribute_file',&
2448  & 'Attribute "'//trim(att_name)//'" is not a character attribute.',&
2449  & err_msg) ) THEN
2450  RETURN
2451  END IF
2452  END IF
2453  ELSE
2454  ! Defining a new attribute
2455  ! Increase the number of file attributes
2456  this_attribute = out_file%num_attributes + 1
2457  IF ( this_attribute .GT. max_file_attributes ) THEN
2458  ! <ERROR STATUS="FATAL">
2459  ! Number of attributes exceeds max_file_attributes for attribute <name>.
2460  ! Increase diag_manager_nml:max_file_attributes.
2461  ! </ERROR>
2462  IF ( fms_error_handler('diag_util_mod::prepend_attribute_file',&
2463  & 'Number of attributes exceeds max_file_attributes for attribute "'&
2464  &//trim(att_name)//'". Increase diag_manager_nml:max_file_attributes.',&
2465  & err_msg) ) THEN
2466  RETURN
2467  END IF
2468  ELSE
2469  out_file%num_attributes = this_attribute
2470  ! Set name and type
2471  out_file%attributes(this_attribute)%name = att_name
2472  out_file%attributes(this_attribute)%type = nf90_char
2473  ! Initialize catt to a blank string, as len_trim doesn't always work on an uninitialized string
2474  out_file%attributes(this_attribute)%catt = ''
2475  END IF
2476  END IF
2477 
2478  ! Only add string only if not already defined
2479  IF ( index(trim(out_file%attributes(this_attribute)%catt), trim(prepend_value)).EQ.0 ) THEN
2480  ! Check if new string length goes beyond the length of catt
2481  length = len_trim(trim(prepend_value)//" "//trim(out_file%attributes(this_attribute)%catt))
2482  IF ( length.GT.len(out_file%attributes(this_attribute)%catt) ) THEN
2483  ! <ERROR STATUS="FATAL">
2484  ! Prepend length for attribute <name> is longer than allowed.
2485  ! </ERROR>
2486  IF ( fms_error_handler('diag_util_mod::prepend_attribute_file',&
2487  & 'Prepend length for attribute "'//trim(att_name)//'" is longer than allowed.',&
2488  & err_msg) ) THEN
2489  RETURN
2490  END IF
2491  END IF
2492  ! Set files
2493  out_file%attributes(this_attribute)%catt =&
2494  & trim(prepend_value)//' '//trim(out_file%attributes(this_attribute)%catt)
2495  out_file%attributes(this_attribute)%len = length
2496  END IF
2497  END SUBROUTINE prepend_attribute_file
2498 
2499  !> @brief Get the a diag_file's start_time as it is defined in the diag_table
2500  !! @return the start_time for the file
2501  function get_file_start_time(file_num) &
2502  result(start_time)
2503  integer, intent(in) :: file_num !< File number of the file to get the start_time from
2504 
2505  TYPE(time_type) :: start_time !< The start_time to return
2506 
2507  start_time = files(file_num)%start_time
2508  end function get_file_start_time
2509 END MODULE diag_util_mod
2510 !> @}
2511 ! close documentation grouping
subroutine, public get_compressed_axes_ids(id, r)
given an index of compressed-by-gathering axis, return an array of axes used in compression....
Definition: diag_axis.F90:1261
subroutine, public get_diag_axis(id, name, units, long_name, cart_name, direction, edges, Domain, DomainU, array_data, num_attributes, attributes, domain_position)
Return information about the axis with index ID.
Definition: diag_axis.F90:467
subroutine, public get_diag_axis_cart(id, cart_name)
Return the axis cartesian.
Definition: diag_axis.F90:565
integer function, public get_axis_global_length(id)
Return the global length of the axis.
Definition: diag_axis.F90:649
type(domain2d) function, public get_domain2d(ids)
Return the 2D domain for the axis IDs given.
Definition: diag_axis.F90:698
type(domainug) function, public get_domainug(id)
Retrun the 1D domain for the axis ID given.
Definition: diag_axis.F90:730
subroutine, public get_diag_axis_data(id, axis_data)
Return the axis data.
Definition: diag_axis.F90:574
character(len=128) function, public get_axis_reqfld(id)
Return the required field names for the axis.
Definition: diag_axis.F90:640
subroutine, public get_diag_axis_name(id, axis_name)
Return the short name of the axis.
Definition: diag_axis.F90:588
character(len=128) function, public get_axis_aux(id)
Return the auxiliary name for the axis.
Definition: diag_axis.F90:631
integer function, public diag_axis_init(name, array_data, units, cart_name, long_name, direction, set_name, edges, Domain, Domain2, DomainU, aux, req, tile_count, domain_position)
Initialize the axis, and return the axis ID.
Definition: diag_axis.F90:110
logical function, public axis_is_compressed(id)
given an axis, returns TRUE if the axis uses compression-by-gathering: that is, if this is an axis fo...
Definition: diag_axis.F90:1241
subroutine, public get_axes_shift(ids, ishift, jshift)
Return the value of the shift for the axis IDs given.
Definition: diag_axis.F90:822
subroutine, public get_diag_axis_domain_name(id, name)
Return the name of the axis' domain.
Definition: diag_axis.F90:601
type(domain1d) function, public get_domain1d(id)
Retrun the 1D domain for the axis ID given.
Definition: diag_axis.F90:685
integer function, public diag_subaxes_init(axis, subdata, start_indx, end_indx, domain_2d)
Create a subaxis on a parent axis.
Definition: diag_axis.F90:382
integer max_output_fields
Maximum number of output fields. Increase via diag_manager_nml.
Definition: diag_data.F90:357
integer pack_size
1 for double and 2 for float
Definition: diag_data.F90:408
integer function get_base_minute()
gets the module variable base_minute
Definition: diag_data.F90:551
integer function get_base_year()
gets the module variable base_year
Definition: diag_data.F90:519
integer num_input_fields
Number of input fields in use.
Definition: diag_data.F90:349
integer function get_base_hour()
gets the module variable base_hour
Definition: diag_data.F90:543
logical use_cmor
Indicates if we should overwrite the MISSING_VALUE to use the CMOR missing value.
Definition: diag_data.F90:372
logical flush_nc_files
Control if diag_manager will force a flush of the netCDF file on each write. Note: changing this to ....
Definition: diag_data.F90:365
integer, parameter glo_reg_val_alt
Alternate value used in the region specification of the diag_table to indicate to use the full axis i...
Definition: diag_data.F90:109
type(time_type) function get_base_time()
gets the module variable base_time
Definition: diag_data.F90:511
integer, parameter diag_field_not_found
Return value for a diag_field that isn't found in the diag_table.
Definition: diag_data.F90:112
integer max_out_per_in_field
Maximum number of output_fields per input_field. Increase via diag_manager_nml.
Definition: diag_data.F90:359
integer max_input_fields
Maximum number of input fields. Increase via diag_manager_nml.
Definition: diag_data.F90:358
logical region_out_use_alt_value
Will determine which value to use when checking a regional output if the region is the full axis or a...
Definition: diag_data.F90:377
integer num_output_fields
Number of output fields in use.
Definition: diag_data.F90:350
integer function get_base_day()
gets the module variable base_day
Definition: diag_data.F90:535
type(time_type) diag_init_time
Time diag_manager_init called. If init_time not included in diag_manager_init call,...
Definition: diag_data.F90:417
integer max_files
Maximum number of output files allowed. Increase via diag_manager_nml.
Definition: diag_data.F90:356
integer num_files
Number of output files currenly in use by the diag_manager.
Definition: diag_data.F90:348
integer max_file_attributes
Maximum number of user definable global attributes per file.
Definition: diag_data.F90:384
real(r8_kind), parameter cmor_missing_value
CMOR standard missing value.
Definition: diag_data.F90:111
logical prepend_date
Should the history file have the start date prepended to the file name. .TRUE. is only supported if t...
Definition: diag_data.F90:386
character(len=2), dimension(:), allocatable fnum_for_domain
If this file number in the array is for the "unstructured" or "2d" domain.
Definition: diag_data.F90:433
integer function get_base_month()
gets the module variable base_month
Definition: diag_data.F90:527
integer, parameter glo_reg_val
Value used in the region specification of the diag_table to indicate to use the full axis instead of ...
Definition: diag_data.F90:107
integer function get_base_second()
gets the module variable base_second
Definition: diag_data.F90:559
integer, parameter max_fields_per_file
Maximum number of fields per file.
Definition: diag_data.F90:89
integer max_field_attributes
Maximum number of user definable attributes per field. Liptak: Changed from 2 to 4 20170718.
Definition: diag_data.F90:382
Define the region for field output.
Definition: diag_data.F90:171
Type to define the diagnostic files that will be written as defined by the diagnostic table.
Definition: diag_data.F90:182
Type to hold the output field description.
Definition: diag_data.F90:249
subroutine, public get_local_indexes(latStart, latEnd, lonStart, lonEnd, istart, iend, jstart, jend)
Find the local start and local end indexes on the local PE for regional output.
Definition: diag_grid.F90:391
subroutine, public diag_field_write(varname, buffer, static, file_num, fileobjU, fileobj, fileobjND, fnum_for_domain, time_in)
Writes diagnostic data out using fms2_io routine.
subroutine, public diag_output_init(file_name, file_title, file_unit, domain, domainU, fileobj, fileobjU, fileobjND, fnum_domain, attributes)
Opens the output file.
Definition: diag_output.F90:89
subroutine, public diag_flush(file_num, fileobjU, fileobj, fileobjND, fnum_for_domain)
Flushes the file into disk.
subroutine, public done_meta_data(file_unit)
Writes axis data to file.
subroutine, public diag_write_time(fileob, rtime_value, time_index, time_name)
Writes the time data to the history file.
subroutine, public write_axis_meta_data(file_unit, axes, fileob, time_ops, time_axis_registered)
Write the axis meta data to file.
type(diag_fieldtype) function, public write_field_meta_data(file_unit, name, axes, units, long_name, range, pack, mval, avg_name, time_method, standard_name, interp_method, attributes, num_attributes, use_UGdomain, fileob)
Write the field meta data to file.
subroutine, public get_subfield_size(axes, outnum)
Get the size, start, and end indices for output fields.
Definition: diag_util.F90:134
subroutine opening_file(file, time, filename_time)
Open file for output, and write the meta data.
Definition: diag_util.F90:1635
subroutine, public init_input_field(module_name, field_name, tile_count)
Initialize the input field.
Definition: diag_util.F90:1310
subroutine, public init_file(name, output_freq, output_units, format, time_units, long_name, tile_count, new_file_freq, new_file_freq_units, start_time, file_duration, file_duration_units, filename_time_bounds)
Initialize the output file.
Definition: diag_util.F90:1064
subroutine, public check_bounds_are_exact_static(out_num, diag_field_id, err_msg)
Check if the array indices for output_fields(out_num) are equal to the output_fields(out_num)buffer u...
Definition: diag_util.F90:1017
integer function, public find_input_field(module_name, field_name, tile_count)
Return the field number for the given module name, field name, and tile number.
Definition: diag_util.F90:1291
subroutine, public write_static(file)
Output all static fields in this file.
Definition: diag_util.F90:2207
subroutine, public diag_data_out(file, field, dat, time, final_call_in, static_write_in, filename_time)
Write data out to file, and if necessary flush the buffers.
Definition: diag_util.F90:2057
subroutine, public init_output_field(module_name, field_name, output_name, output_file, time_method, pack, tile_count, local_coord)
Initialize the output field.
Definition: diag_util.F90:1343
type(time_type) function, public get_file_start_time(file_num)
Get the a diag_file's start_time as it is defined in the diag_table.
Definition: diag_util.F90:2503
logical function a_lessthan_b(a, b)
return true iff a<b.
Definition: diag_util.F90:819
integer function find_file(name, tile_count)
Return the file number for file name and tile.
Definition: diag_util.F90:1274
subroutine, public check_duplicate_output_fields(err_msg)
Checks to see if output_name and output_file are unique in output_fields.
Definition: diag_util.F90:2238
subroutine, public sync_file_times(file_id, init_time, err_msg)
Synchronize the file's start and close times with the model start and end times.
Definition: diag_util.F90:1244
subroutine fms_diag_check_out_of_bounds_r8(ofb, bounds, output_name, module_name, err_msg)
Checks if the array indices for output_field buffer (ofb) are outside the are outside the bounding bo...
Definition: diag_util.F90:906
character(len=1), public field_log_separator
separator used for csv-style log of registered fields set by nml in diag_manager init
Definition: diag_util.F90:116
subroutine, public check_out_of_bounds(out_num, diag_field_id, err_msg)
Checks if the array indices for output_fields(out_num) are outside the output_fields(out_num)buffer u...
Definition: diag_util.F90:842
subroutine, public fms_diag_check_bounds_are_exact_static(current_bounds, bounds, output_name, module_name, err_msg)
Check if the array indices specified in the bounding box "current_bounds" are equal to those specifie...
Definition: diag_util.F90:1039
logical function a_greaterthan_b(a, b)
return true iff a>b.
Definition: diag_util.F90:826
subroutine, public get_subfield_vert_size(axes, outnum)
Get size, start and end indices for output fields.
Definition: diag_util.F90:416
logical function compare_buffer_bounds_to_size(current_bounds, bounds, error_str, lowerb_comp, upperb_comp)
Compares the bounding indices of an array specified in "current_bounds" to the corresponding lower an...
Definition: diag_util.F90:767
subroutine, public log_diag_field_info(module_name, field_name, axes, long_name, units, missing_value, range, dynamic)
Writes brief diagnostic field info to the log file.
Definition: diag_util.F90:636
subroutine, public update_bounds(out_num, lower_i, upper_i, lower_j, upper_j, lower_k, upper_k)
Update the output_fields x, y, and z min and max boundaries (array indices) with the six specified bo...
Definition: diag_util.F90:747
subroutine prepend_attribute_file(out_file, att_name, prepend_value, err_msg)
Prepends the attribute value to an already existing attribute. If the attribute isn't yet defined,...
Definition: diag_util.F90:2413
subroutine prepend_attribute_field(out_field, att_name, prepend_value, err_msg)
Prepends the attribute value to an already existing attribute. If the attribute isn't yet defined,...
Definition: diag_util.F90:2297
subroutine, public diag_util_init()
Write the version number of this file to the log file.
Definition: diag_util.F90:124
subroutine attribute_init_file(out_file, err_msg)
Allocates the atttype in out_file.
Definition: diag_util.F90:2384
logical function a_noteq_b(a, b)
return true iff a /= b
Definition: diag_util.F90:833
subroutine attribute_init_field(out_field, err_msg)
Allocates the atttype in out_field.
Definition: diag_util.F90:2268
subroutine check_and_open(file, time, do_write, filename_time)
Checks if it is time to open a new file.
Definition: diag_util.F90:2164
subroutine, public check_bounds_are_exact_dynamic(out_num, diag_field_id, Time, err_msg)
This is an adaptor to the check_bounds_are_exact_dynamic_modern function to maintain an interface ser...
Definition: diag_util.F90:990
subroutine, public fms_diag_check_bounds_are_exact_dynamic(current_bounds, bounds, output_name, module_name, Time, field_prev_Time, err_msg)
Checks that array indices specified in the bounding box "current_bounds" are identical to those in th...
Definition: diag_util.F90:940
integer function get_index(number, array)
Find index i of array such that array(i) is closest to number.
Definition: diag_util.F90:547
subroutine fms_diag_check_out_of_bounds_r4(ofb, bounds, output_name, module_name, err_msg)
Checks if the array indices for output_fields(out_num) are outside the output_fields(out_num)buffer u...
Definition: diag_util.F90:874
Allocates the atttype in out_file.
Definition: diag_util.F90:97
Prepend a value to a string attribute in the output field or output file.
Definition: diag_util.F90:90
character(len=128) function, public get_time_string(filename, current_time)
This function determines a string based on current time. This string is used as suffix in output file...
real function, public get_date_dif(t2, t1, units)
Return the difference between two times in units.
subroutine, public write_version_number(version, tag, unit)
Prints to the log file (or a specified unit) the version id string and tag name.
Definition: fms.F90:701
logical function, public fms_error_handler(routine, message, err_msg)
Facilitates the control of fatal error conditions.
Definition: fms.F90:468
subroutine, public error_mesg(routine, message, level)
Print notes, warnings and error messages; terminates program for warning and error messages....
Definition: fms.F90:441
integer function mpp_get_tile_npes(domain)
Returns number of processors used on current tile.
subroutine mpp_get_domain_components(domain, x, y, tile_count)
Retrieve 1D components of 2D decomposition.
integer function, dimension(size(domain%tile_id(:))) mpp_get_tile_id(domain)
Returns the tile_id on current pe.
logical function mpp_mosaic_defined()
Accessor function for value of mosaic_defined.
integer function mpp_get_current_ntile(domain)
Returns number of tile on current pe.
integer function mpp_get_ntile_count(domain)
Returns number of tiles in mosaic.
These routines retrieve the axis specifications associated with the compute domains....
Modifies the extents (compute, data and global) of a given domain.
One dimensional domain used to manage shared data access between pes.
The domain2D type contains all the necessary information to define the global, compute and data domai...
Domain information for managing data on unstructured grids.
integer function mpp_npes()
Returns processor count for current pelist.
Definition: mpp_util.inc:420
integer function mpp_pe()
Returns processor ID.
Definition: mpp_util.inc:406
logical function, public leap_year(Time, err_msg)
Returns true if the year corresponding to the input time is a leap year (for default calendar)....
subroutine, public get_time(Time, seconds, days, ticks, err_msg)
Returns days and seconds ( < 86400 ) corresponding to a time. err_msg should be checked for any error...
subroutine, public get_date(time, year, month, day, hour, minute, second, tick, err_msg)
Gets the date for different calendar types. Given a time_interval, returns the corresponding date und...
type(time_type) function, public increment_time(Time, seconds, days, ticks, err_msg, allow_neg_inc)
Increments a time by seconds and days.
integer function, public get_calendar_type()
Returns default calendar type for mapping from time to date.
Type to represent amounts of time. Implemented as seconds and days to allow for larger intervals.
Data structure holding a 3D bounding box. It is commonlyused to represent the interval bounds or limi...