FMS  2025.04
Flexible Modeling System
fms_diag_axis_object.F90
1 !***********************************************************************
2 !* Apache License 2.0
3 !*
4 !* This file is part of the GFDL Flexible Modeling System (FMS).
5 !*
6 !* Licensed under the Apache License, Version 2.0 (the "License");
7 !* you may not use this file except in compliance with the License.
8 !* You may obtain a copy of the License at
9 !*
10 !* http://www.apache.org/licenses/LICENSE-2.0
11 !*
12 !* FMS is distributed in the hope that it will be useful, but WITHOUT
13 !* WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied;
14 !* without even the implied warranty of MERCHANTABILITY or FITNESS FOR A
15 !* PARTICULAR PURPOSE. See the License for the specific language
16 !* governing permissions and limitations under the License.
17 !***********************************************************************
18 
19 !> @defgroup fms_diag_axis_object_mod fms_diag_axis_object_mod
20 !> @ingroup diag_manager
21 !! @brief fms_diag_axis_object_mod stores the diag axis object, a diag domain
22 !! object, and a subaxis object.
23 
24 !> @file
25 !> @brief File for @ref diag_axis_object_mod
26 
27 !> @addtogroup fms_diag_axis_object_mod
28 !> @{
29 module fms_diag_axis_object_mod
30 #ifdef use_yaml
31  use mpp_domains_mod, only: domain1d, domain2d, domainug, mpp_get_compute_domain, center, &
32  & mpp_get_global_domain, north, east, mpp_get_tile_id, &
34  use platform_mod, only: r8_kind, r4_kind, i4_kind, i8_kind
35  use diag_data_mod, only: diag_atttype, max_axes, no_domain, two_d_domain, ug_domain, &
37  diag_null, index_gridtype, latlon_gridtype, pack_size_str, &
40  use mpp_mod, only: fatal, mpp_error, uppercase, mpp_pe, mpp_root_pe, stdout
41  use fms2_io_mod, only: fmsnetcdffile_t, fmsnetcdfdomainfile_t, fmsnetcdfunstructureddomainfile_t, &
42  & register_axis, register_field, register_variable_attribute, write_data
43  use fms_diag_yaml_mod, only: subregion_type, diag_yaml, max_subaxes
44  use diag_grid_mod, only: get_local_indices_cubesphere => get_local_indexes
45  use axis_utils2_mod, only: nearest_index
46  implicit none
47 
48  PRIVATE
49 
56 
57  !> @}
58 
59  !> @brief Type to hold the domain info for an axis
60  !! This type was created to avoid having to send in "Domain", "Domain2", "DomainUG" as arguments into subroutines
61  !! and instead only 1 class(diagDomain_t) argument can be send
62  !> @ingroup diag_axis_object_mod
64  contains
65  procedure :: set => set_axis_domain
66  procedure :: length => get_length
67  procedure :: get_ntiles
68  end type diagdomain_t
69 
70  !> @brief Type to hold the 1d domain
71  type, extends(diagdomain_t) :: diagdomain1d_t
72  type(domain1d) :: domain !< 1d Domain of the axis
73  end type
74 
75  !> @brief Type to hold the 2d domain
76  type, extends(diagdomain_t) :: diagdomain2d_t
77  type(domain2d) :: domain2 !< 2d Domain of an "X" or "Y" axis
78  end type
79 
80  !> @brief Type to hold the unstructured domain
81  type, extends(diagdomain_t) :: diagdomainug_t
82  type(domainug) :: domainug !< Domain of "U" axis
83  end type
84 
85  !> @brief Type to hold the diagnostic axis description.
86  !> @ingroup diag_axis_object_mod
88  INTEGER , private :: axis_id !< ID of the axis
89 
90  contains
91  procedure :: get_parent_axis_id
92  procedure :: get_subaxes_id
93  procedure :: get_axis_name
94  procedure :: is_z_axis
95  procedure :: write_axis_metadata
96  procedure :: write_axis_data
97  procedure :: add_structured_axis_ids
98  procedure :: get_structured_axis
99  procedure :: is_unstructured_grid
100  procedure :: get_edges_id
101  END TYPE fmsdiagaxis_type
102 
103  !> @brief Type to hold the diag_axis (either subaxis or a full axis)
104  !> @ingroup diag_axis_object_mod
106  class(fmsDiagAxis_type), allocatable :: axis
107  end type
108 
109  !> @brief Type to hold the subaxis
110  !> @ingroup diag_axis_object_mod
112  CHARACTER(len=:), ALLOCATABLE , private :: subaxis_name !< Name of the subaxis
113  INTEGER , private :: starting_index !< Starting index of the subaxis relative to the
114  !! parent axis
115  INTEGER , private :: ending_index !< Ending index of the subaxis relative to the
116  !! parent axis
117  INTEGER , private :: parent_axis_id !< Id of the parent_axis
118  INTEGER , private :: compute_idx(2) !< Starting and ending index of the compute domain
119  INTEGER, allocatable, private :: global_idx(:) !< Starting and ending index of the global domain
120  real(kind=r4_kind), allocatable, private :: zbounds(:) !< Bounds of the Z axis
121  contains
122  procedure :: fill_subaxis
123  procedure :: axis_length
124  procedure :: get_starting_index
125  procedure :: get_ending_index
126  procedure :: get_compute_indices
127  END TYPE fmsdiagsubaxis_type
128 
129  !> @brief Type to hold the diurnal axis
130  !> @ingroup diag_axis_object_mod
132  INTEGER , private :: ndiurnal_samples !< The number of diurnal samples
133  CHARACTER(len=:), ALLOCATABLE, private :: axis_name !< The diurnal axis name
134  CHARACTER(len=:), ALLOCATABLE, private :: long_name !< The longname of the diurnal axis
135  CHARACTER(len=:), ALLOCATABLE, private :: units !< The units
136  INTEGER , private :: edges_id !< The id of the diurnal edges
137  CHARACTER(len=:), ALLOCATABLE, private :: edges_name !< The name of the edges axis
138  CLASS(*), ALLOCATABLE, private :: diurnal_data(:) !< The diurnal data
139 
140  contains
141  procedure :: get_diurnal_axis_samples
142  procedure :: write_diurnal_metadata
143  END TYPE fmsdiagdiurnalaxis_type
144 
145  !> @brief Type to hold the diagnostic axis description.
146  !> @ingroup diag_axis_object_mod
148  CHARACTER(len=:), ALLOCATABLE, private :: axis_name !< Name of the axis
149  CHARACTER(len=:), ALLOCATABLE, private :: units !< Units of the axis
150  CHARACTER(len=:), ALLOCATABLE, private :: long_name !< Long_name attribute of the axis
151  CHARACTER(len=1) , private :: cart_name !< Cartesian name "X", "Y", "Z", "T", "U", "N"
152  CLASS(*), ALLOCATABLE, private :: axis_data(:) !< Data of the axis
153  CHARACTER(len=:), ALLOCATABLE, private :: type_of_data !< The type of the axis_data ("float" or "double")
154  !< TO DO this can be a dlinked to avoid having limits
155  integer, ALLOCATABLE, private :: subaxis(:) !< Array of subaxis
156  integer , private :: nsubaxis !< Number of subaxis
157  class(diagdomain_t),ALLOCATABLE, private :: axis_domain !< Domain
158  INTEGER , private :: type_of_domain !< The type of domain ("NO_DOMAIN", "TWO_D_DOMAIN",
159  !! or "UG_DOMAIN")
160  INTEGER , private :: length !< Global axis length
161  INTEGER , private :: direction !< Direction of the axis 0, 1, -1
162  INTEGER, ALLOCATABLE, private :: edges_id !< Axis ID for the edges axis
163  !! This axis will be written to the file
164  CHARACTER(len=:), ALLOCATABLE, private :: edges_name !< Name for the previously defined "edges axis"
165  !! This will be written as an attribute
166  CHARACTER(len=:), ALLOCATABLE, private :: aux !< Auxiliary name, can only be <TT>geolon_t</TT>
167  !! or <TT>geolat_t</TT>
168  CHARACTER(len=128) , private :: req !< Required field names.
169  INTEGER , private :: tile_count !< The number of tiles
170  TYPE(fmsdiagattribute_type),allocatable , private :: attributes(:) !< Array to hold user definable attributes
171  INTEGER , private :: num_attributes !< Number of defined attibutes
172  INTEGER , private :: domain_position !< The position in the doman (NORTH, EAST or CENTER)
173  integer, allocatable , private :: structured_ids(:) !< If the axis is in the unstructured grid,
174  !! this is the axis ids of the structured axis
175  CHARACTER(len=:), ALLOCATABLE, private :: set_name !< Name of the axis set. This is to distinguish
176  !! two axis with the same name
177 
178  contains
179 
180  PROCEDURE :: add_axis_attribute
181  PROCEDURE :: register => register_diag_axis_obj
182  PROCEDURE :: axis_length => get_axis_length
183  PROCEDURE :: set_edges
184  PROCEDURE :: set_axis_id
185  PROCEDURE :: get_compute_domain
186  PROCEDURE :: get_indices
187  PROCEDURE :: get_global_io_domain
188  PROCEDURE :: get_aux
189  PROCEDURE :: has_aux
190  PROCEDURE :: get_set_name
191  PROCEDURE :: has_set_name
192  PROCEDURE :: is_x_or_y_axis
193  PROCEDURE :: get_dim_size_layout
194  ! TO DO:
195  ! Get/has/is subroutines as needed
196  END TYPE fmsdiagfullaxis_type
197 
198  !> @addtogroup fms_diag_yaml_mod
199  !> @{
200  contains
201 
202  !!!!!!!!!!!!!!!!! DIAG AXIS PROCEDURES !!!!!!!!!!!!!!!!!
203  !> @brief Initialize the axis
204  subroutine register_diag_axis_obj(this, axis_name, axis_data, units, cart_name, long_name, direction,&
205  & set_name, Domain, Domain2, DomainU, aux, req, tile_count, domain_position, axis_length )
206  class(fmsdiagfullaxis_type),INTENT(inout):: this !< Diag_axis obj
207  CHARACTER(len=*), INTENT(in) :: axis_name !< Name of the axis
208  class(*), INTENT(in) :: axis_data(:) !< Array of coordinate values
209  CHARACTER(len=*), INTENT(in) :: units !< Units for the axis
210  CHARACTER(len=1), INTENT(in) :: cart_name !< Cartesian axis ("X", "Y", "Z", "T", "U", "N")
211  CHARACTER(len=*), INTENT(in), OPTIONAL :: long_name !< Long name for the axis.
212  CHARACTER(len=*), INTENT(in), OPTIONAL :: set_name !< Name of the parent axis, if it is a subaxis
213  INTEGER, INTENT(in), OPTIONAL :: direction !< Indicates the direction of the axis
214  TYPE(domain1d), INTENT(in), OPTIONAL :: Domain !< 1D domain
215  TYPE(domain2d), INTENT(in), OPTIONAL :: Domain2 !< 2D domain
216  TYPE(domainug), INTENT(in), OPTIONAL :: DomainU !< Unstructured domain
217  CHARACTER(len=*), INTENT(in), OPTIONAL :: aux !< Auxiliary name, can only be <TT>geolon_t</TT>
218  !! or <TT>geolat_t</TT>
219  CHARACTER(len=*), INTENT(in), OPTIONAL :: req !< Required field names.
220  INTEGER, INTENT(in), OPTIONAL :: tile_count !< Number of tiles
221  INTEGER, INTENT(in), OPTIONAL :: domain_position !< Domain position, "NORTH" or "EAST"
222  integer, intent(in), optional :: axis_length !< The length of the axis size(axis_data(:))
223 
224  this%axis_name = trim(axis_name)
225  this%units = trim(units)
226  this%cart_name = uppercase(cart_name)
227  call check_if_valid_cart_name(this%cart_name)
228 
229  if (present(long_name)) this%long_name = trim(long_name)
230 
231  select type (axis_data)
232  type is (real(kind=r8_kind))
233  allocate(real(kind=r8_kind) :: this%axis_data(axis_length))
234  this%axis_data = axis_data
235  this%length = axis_length
236  this%type_of_data = "double" !< This is what fms2_io expects in the register_field call
237  type is (real(kind=r4_kind))
238  allocate(real(kind=r4_kind) :: this%axis_data(axis_length))
239  this%axis_data = axis_data
240  this%length = axis_length
241  this%type_of_data = "float" !< This is what fms2_io expects in the register_field call
242  class default
243  call mpp_error(fatal, "The axis_data in your diag_axis_init call is not a supported type. &
244  & Currently only r4 and r8 data is supported.")
245  end select
246 
247  this%type_of_domain = no_domain
248  if (present(domain)) then
249  if (present(domain2) .or. present(domainu)) call mpp_error(fatal, &
250  "The presence of Domain with any other domain type is prohibited. "//&
251  "Check you diag_axis_init call for axis_name:"//trim(axis_name))
252  allocate(diagdomain1d_t :: this%axis_domain)
253  call this%axis_domain%set(domain=domain)
254  else if (present(domain2)) then
255  if (present(domainu)) call mpp_error(fatal, &
256  "The presence of Domain2 with any other domain type is prohibited. "//&
257  "Check you diag_axis_init call for axis_name:"//trim(axis_name))
258  allocate(diagdomain2d_t :: this%axis_domain)
259  call this%axis_domain%set(domain2=domain2)
260  this%type_of_domain = two_d_domain
261  else if (present(domainu)) then
262  allocate(diagdomainug_t :: this%axis_domain)
263  call this%axis_domain%set(domainu=domainu)
264  this%type_of_domain = ug_domain
265  endif
266 
267  this%tile_count = 1
268  if (present(tile_count)) this%tile_count = tile_count
269 
270  this%domain_position = center
271  if (present(domain_position)) this%domain_position = domain_position
272  call check_if_valid_domain_position(this%domain_position)
273 
274  this%direction = 0
275  if (present(direction)) this%direction = direction
276  call check_if_valid_direction(this%direction)
277 
278  if (present(aux)) this%aux = trim(aux)
279  if (present(req)) this%req = trim(req)
280  this%set_name = ""
281  if (present(set_name)) this%set_name = trim(set_name)
282 
283  if (max_subaxes .gt. 0) then
284  allocate(this%subaxis(max_subaxes))
285  this%subaxis = diag_null
286  endif
287 
288  this%nsubaxis = 0
289  this%num_attributes = 0
290  end subroutine register_diag_axis_obj
291 
292  !> @brief Add an attribute to an axis
293  subroutine add_axis_attribute(this, att_name, att_value)
294  class(fmsdiagfullaxis_type),INTENT(INOUT) :: this !< diag_axis obj
295  character(len=*), intent(in) :: att_name !< Name of the attribute
296  class(*), intent(in) :: att_value(:) !< The attribute value to add
297 
298  integer :: j !< obj%num_attributes (for less typing)
299 
300  if (.not. allocated(this%attributes)) &
301  allocate(this%attributes(max_axis_attributes))
302 
303  this%num_attributes = this%num_attributes + 1
304 
305  j = this%num_attributes
306  call this%attributes(j)%add(att_name, att_value)
307  end subroutine add_axis_attribute
308 
309  !> @brief Write the axis meta data to an open fileobj
310  subroutine write_axis_metadata(this, fms2io_fileobj, edges_in_file, parent_axis)
311  class(fmsdiagaxis_type), target, INTENT(IN) :: this !< diag_axis obj
312  class(fmsnetcdffile_t), INTENT(INOUT) :: fms2io_fileobj!< Fms2_io fileobj to write the data to
313  logical, INTENT(IN) :: edges_in_file !< .True. if the edges to this axis are
314  !! already in the file
315  class(fmsdiagaxis_type), OPTIONAL, target, INTENT(IN) :: parent_axis !< If the axis is a subaxis, axis object
316  !! for the parent axis (this will be used
317  !! to get some of the metadata info)
318 
319  character(len=:), ALLOCATABLE :: axis_edges_name !< Name of the edges, if it exist
320  character(len=:), pointer :: axis_name !< Name of the axis
321  integer :: axis_length !< Size of the axis
322  integer :: i !< For do loops
323  type(fmsdiagfullaxis_type), pointer :: diag_axis !< Local pointer to the diag_axis
324 
325  integer :: type_of_domain !< The type of domain the current axis is in
326  logical :: is_subaxis !< .true. if the axis is a subaxis
327  logical :: needs_domain_decomposition !< .True. if the axis needs the domain decomposition attribute
328  !! (i.e for "X" and "Y" subaxis)
329  integer :: domain_decomposition(4) !< indices of the global (1:2) and compute (3:4) domain for a "X" and "Y" subaxis
330 
331  is_subaxis = .false.
332  needs_domain_decomposition = .false.
333 
334  select type(this)
335  type is (fmsdiagfullaxis_type)
336  axis_name => this%axis_name
337  axis_length = this%length
338  diag_axis => this
339  type_of_domain = this%type_of_domain
340  type is (fmsdiagsubaxis_type)
341  is_subaxis = .true.
342  axis_name => this%subaxis_name
343  axis_length = this%ending_index - this%starting_index + 1
344  if (allocated(this%global_idx)) then
345  needs_domain_decomposition = .true.
346  domain_decomposition(1:2) = this%global_idx
347  domain_decomposition(3) = this%starting_index
348  domain_decomposition(4) = this%ending_index
349  endif
350  !< Get all the other information from the parent axis (i.e the cart_name, units, etc)
351  if (present(parent_axis)) then
352  select type(parent_axis)
353  type is (fmsdiagfullaxis_type)
354  diag_axis => parent_axis
355  end select
356  endif
357  type_of_domain = no_domain !< All subaxes are treated as non-domain decomposed (each rank writes it own file)
358  type is (fmsdiagdiurnalaxis_type)
359  call this%write_diurnal_metadata(fms2io_fileobj)
360  return
361  end select
362 
363  !< Add the axis as a dimension in the netcdf file based on the type of axis_domain and the fileobj type
364  select type (fms2io_fileobj)
365  !< The register_field calls need to be inside the select type block so that it can go inside the correct
366  !! register_field interface
367  type is (fmsnetcdffile_t)
368  !< Here the axis is not domain decomposed (i.e z_axis)
369  call register_axis(fms2io_fileobj, axis_name, axis_length)
370  call register_field(fms2io_fileobj, axis_name, diag_axis%type_of_data, (/axis_name/))
371  if (needs_domain_decomposition) then
372  call register_variable_attribute(fms2io_fileobj, axis_name, "domain_decomposition", &
373  domain_decomposition)
374  endif
375  type is (fmsnetcdfdomainfile_t)
376  select case (type_of_domain)
377  case (no_domain)
378  !< Here the fms2io_fileobj is domain decomposed, but the axis is not
379  !! Domain decomposed fileobjs can have axis that are not domain decomposed (i.e "Z" axis)
380  call register_axis(fms2io_fileobj, axis_name, axis_length)
381  call register_field(fms2io_fileobj, axis_name, diag_axis%type_of_data, (/axis_name/))
382  case (two_d_domain)
383  !< Here the axis is domain decomposed
384  call register_axis(fms2io_fileobj, axis_name, diag_axis%cart_name, domain_position=diag_axis%domain_position)
385  call register_field(fms2io_fileobj, axis_name, diag_axis%type_of_data, (/axis_name/))
386  end select
387  type is (fmsnetcdfunstructureddomainfile_t)
388  select case (type_of_domain)
389  case (ug_domain)
390  !< Here the axis is in a unstructured domain
391  call register_axis(fms2io_fileobj, axis_name)
392  call register_field(fms2io_fileobj, axis_name, diag_axis%type_of_data, (/axis_name/))
393  case default
394  !< Here the fms2io_fileobj is in the unstructured domain, but the axis is not
395  !< Unstructured domain fileobjs can have axis that are not domain decomposed (i.e "Z" axis)
396  call register_axis(fms2io_fileobj, axis_name, axis_length)
397  call register_field(fms2io_fileobj, axis_name, diag_axis%type_of_data, (/axis_name/))
398  end select
399  end select
400 
401  !< Write its metadata
402  if(allocated(diag_axis%long_name)) &
403  call register_variable_attribute(fms2io_fileobj, axis_name, "long_name", diag_axis%long_name, &
404  str_len=len_trim(diag_axis%long_name))
405 
406  if (diag_axis%cart_name .NE. "N") &
407  call register_variable_attribute(fms2io_fileobj, axis_name, "axis", diag_axis%cart_name, str_len=1)
408 
409  if (trim(diag_axis%units) .NE. "none") &
410  call register_variable_attribute(fms2io_fileobj, axis_name, "units", diag_axis%units, &
411  str_len=len_trim(diag_axis%units))
412 
413  select case (diag_axis%direction)
414  case (direction_up)
415  call register_variable_attribute(fms2io_fileobj, axis_name, "positive", "up", str_len=2)
416  case (direction_down)
417  call register_variable_attribute(fms2io_fileobj, axis_name, "positive", "down", str_len=4)
418  end select
419 
420  !< Ignore the edges attribute, if the edges are already in the file or if it is subaxis
421  if (.not. edges_in_file .and. allocated(diag_axis%edges_name) .and. .not. is_subaxis) then
422  call register_variable_attribute(fms2io_fileobj, axis_name, "edges", diag_axis%edges_name, &
423  str_len=len_trim(diag_axis%edges_name))
424  endif
425 
426  if(allocated(diag_axis%attributes)) then
427  do i = 1, diag_axis%num_attributes
428  select type (att_value => diag_axis%attributes(i)%att_value)
429  type is (character(len=*))
430  call register_variable_attribute(fms2io_fileobj, axis_name, diag_axis%attributes(i)%att_name, &
431  trim(att_value(1)), str_len=len_trim(att_value(1)))
432  class default
433  call register_variable_attribute(fms2io_fileobj, axis_name, diag_axis%attributes(i)%att_name, att_value)
434  end select
435  enddo
436  endif
437 
438  end subroutine write_axis_metadata
439 
440  !> @brief Write the axis data to an open fms2io_fileobj
441  subroutine write_axis_data(this, fms2io_fileobj, parent_axis)
442  class(fmsdiagaxis_type), target, INTENT(IN) :: this !< diag_axis obj
443  class(fmsnetcdffile_t), INTENT(INOUT) :: fms2io_fileobj!< Fms2_io fileobj to write the data to
444  class(fmsdiagaxis_type), OPTIONAL, target, INTENT(IN) :: parent_axis !< The parent axis if this is a subaxis
445 
446  integer :: i !< Starting index of a sub_axis
447  integer :: j !< Ending index of a sub_axis
448  integer :: global_io_index(2)!< Global io domain starting and ending index
449  select type(this)
450  type is (fmsdiagfullaxis_type)
451  call this%get_global_io_domain(global_io_index, fms2io_fileobj%is_file_using_netcdf_mpi())
452  call write_data(fms2io_fileobj, this%axis_name, this%axis_data(global_io_index(1):global_io_index(2)))
453  type is (fmsdiagsubaxis_type)
454  i = this%starting_index
455  j = this%ending_index
456 
457  if (present(parent_axis)) then
458  select type(parent_axis)
459  type is (fmsdiagfullaxis_type)
460  call write_data(fms2io_fileobj, this%subaxis_name, parent_axis%axis_data(i:j))
461  end select
462  endif
463  type is (fmsdiagdiurnalaxis_type)
464  call write_data(fms2io_fileobj, this%axis_name, this%diurnal_data)
465  end select
466  end subroutine write_axis_data
467 
468 
469  !> @brief Defined a new diurnal axis
470  subroutine define_diurnal_axis(diag_axis, naxis, n_diurnal_samples, is_edges)
471  class(fmsdiagaxiscontainer_type), target, intent(inout) :: diag_axis(:) !< Array of axis containers
472  integer, intent(inout) :: naxis !< Number of axis that have
473  !! been defined
474  integer, intent(in) :: n_diurnal_samples !< The number of diurnal samples
475  !! for the curent axis
476  logical, intent(in) :: is_edges !< Flag indicating if this is
477  !! an edge axis
478 
479  CHARACTER(32) :: axis_name !< name of the axis
480  CHARACTER(32) :: long_name !< long name of the axis
481  CHARACTER(32) :: edges_name !< name of the axis edge
482  CHARACTER(128) :: units !< units of the axis
483  real(kind=r8_kind), allocatable :: diurnal_data(:) !< Data for the axis
484  integer :: edges_id !< Id of the axis edge
485  integer :: i !< For do loops
486 
487  naxis = naxis + 1
488 
489  axis_name = ''
490  edges_name = ''
491  if (is_edges) then
492  WRITE (axis_name,'(a,i2.2)') 'time_of_day_edges_', n_diurnal_samples
493  long_name = "time of day edges"
494  allocate(diurnal_data(n_diurnal_samples + 1))
495  diurnal_data(1) = 0.0
496  edges_id = diag_null
497  do i = 1, n_diurnal_samples
498  diurnal_data(i+1) = 24.0* real(i)/n_diurnal_samples
499  enddo
500  else
501  WRITE (axis_name,'(a,i2.2)') 'time_of_day_', n_diurnal_samples
502  long_name = "time of day"
503  allocate(diurnal_data(n_diurnal_samples))
504  edges_id = naxis -1 !< The diurnal edges is the last defined axis
505  do i = 1, n_diurnal_samples
506  diurnal_data(i) = 24.0*(real(i)-0.5)/n_diurnal_samples
507  enddo
508  WRITE (edges_name,'(a,i2.2)') 'time_of_day_edges_', n_diurnal_samples
509  endif
510 
511  WRITE (units,11) 'hours', get_base_year(), get_base_month(), &
513 11 FORMAT(a,' since ',i4.4,'-',i2.2,'-',i2.2,' ',i2.2,':',i2.2,':',i2.2)
514 
515  allocate(fmsdiagdiurnalaxis_type :: diag_axis(naxis)%axis)
516  select type (diurnal_axis => diag_axis(naxis)%axis)
517  type is (fmsdiagdiurnalaxis_type)
518  diurnal_axis%axis_id = naxis
519  diurnal_axis%ndiurnal_samples = n_diurnal_samples
520  diurnal_axis%axis_name = trim(axis_name)
521  diurnal_axis%long_name = trim(long_name)
522  diurnal_axis%units = trim(units)
523  diurnal_axis%diurnal_data = diurnal_data
524  diurnal_axis%edges_id = edges_id
525  if (is_edges) &
526  WRITE (edges_name,'(a,i2.2)') 'time_of_day_edges_', n_diurnal_samples
527  diurnal_axis%edges_name = trim(edges_name)
528  end select
529  end subroutine define_diurnal_axis
530 
531  !< @brief Determine if the axis is in the unstructured grid
532  !! @return .True. if the axis is in unstructured grid
533  pure logical function is_unstructured_grid(this)
534  class(fmsdiagaxis_type), target, INTENT(in) :: this !< diag_axis obj
535 
536  is_unstructured_grid = .false.
537  select type (this)
538  type is (fmsdiagfullaxis_type)
539  is_unstructured_grid = trim(this%cart_name) .eq. "U"
540  end select
541  end function is_unstructured_grid
542 
543  !< @brief Adds the structured axis ids to the axis object
544  subroutine add_structured_axis_ids(this, axis_ids)
545  class(fmsdiagaxis_type), target, INTENT(inout) :: this !< diag_axis obj
546  integer, intent(in) :: axis_ids(2) !< axis ids to add to the axis object
547 
548  select type (this)
549  type is (fmsdiagfullaxis_type)
550  allocate(this%structured_ids(2))
551  this%structured_ids = axis_ids
552  end select
553  end subroutine add_structured_axis_ids
554 
555  !< @brief Get the structured axis ids from the axis object
556  !! @return the structured axis ids
557  pure function get_structured_axis(this) &
558  result(rslt)
559  class(fmsdiagaxis_type), target, INTENT(in) :: this !< diag_axis obj
560  integer :: rslt(2)
561 
562  rslt = diag_null
563  select type (this)
564  type is (fmsdiagfullaxis_type)
565  rslt = this%structured_ids
566  end select
567  end function get_structured_axis
568 
569 
570  !< @brief Get the edges_id of an axis_object
571  !! @return The edges_id of an axis object
572  pure integer function get_edges_id(this)
573  class(fmsdiagaxis_type), INTENT(in) :: this !< diag_axis obj
574 
575  get_edges_id = diag_null
576  select type (this)
577  type is (fmsdiagfullaxis_type)
578  if (allocated(this%edges_id)) get_edges_id = this%edges_id
579  end select
580  end function
581 
582  !> @brief Get the starting and ending indices of the global io domain of the axis
583  subroutine get_global_io_domain(this, global_io_index, use_collective_writes)
584  class(fmsdiagfullaxis_type), target, intent(in) :: this !< diag_axis obj
585  integer, intent(out) :: global_io_index(2) !< Global io domain starting and ending index
586  logical, intent(in) :: use_collective_writes !< .True. if using collective writes
587 
588  type(domain2d), pointer :: io_domain !< pointer to the io domain
589 
590  global_io_index(1) = 1
591  global_io_index(2) = this%length
592 
593  if (allocated(this%axis_domain)) then
594  select type(domain => this%axis_domain)
595  type is (diagdomain2d_t)
596  if (use_collective_writes) then
597  io_domain => domain%domain2
598  else
599  io_domain => mpp_get_io_domain(domain%domain2)
600  endif
601 
602  if (this%cart_name .eq. "X") then
603  call mpp_get_global_domain(io_domain, xbegin=global_io_index(1), xend=global_io_index(2), &
604  position=this%domain_position)
605  elseif (this%cart_name .eq. "Y") then
606  call mpp_get_global_domain(io_domain, ybegin=global_io_index(1), yend=global_io_index(2), &
607  position=this%domain_position)
608  endif
609  end select
610  endif
611  end subroutine get_global_io_domain
612 
613  !> @brief Get the length of the axis
614  !> @return axis length
615  function get_axis_length(this) &
616  result(axis_length)
617  class(fmsdiagfullaxis_type), intent(in) :: this !< diag_axis obj
618  integer :: axis_length
619 
620  !< If the axis is domain decomposed axis_length will be set to the length for the current PE:
621  if (allocated(this%axis_domain)) then
622  axis_length = this%axis_domain%length(this%cart_name, this%domain_position, this%length)
623  else
624  axis_length = this%length
625  endif
626 
627  end function
628 
629 
630  !> @brief Determine if an axis object has an auxiliary name
631  !! @return .true. if an axis object has an auxiliary name
632  pure function has_aux(this) &
633  result(rslt)
634  class(fmsdiagfullaxis_type), intent(in) :: this !< diag_axis obj
635  logical :: rslt
636 
637  rslt = .false.
638  if (allocated(this%aux)) rslt = trim(this%aux) .ne. ""
639  end function has_aux
640 
641  !> @brief Determine if an axis object has a set_name
642  !! @return .true. if an axis object has a set_name
643  pure function has_set_name(this) &
644  result(rslt)
645  class(fmsdiagfullaxis_type), intent(in) :: this !< diag_axis obj
646  logical :: rslt
647 
648  rslt = .false.
649  if (allocated(this%set_name)) rslt = trim(this%set_name) .ne. ""
650  end function has_set_name
651 
652  !> @brief Determine if an axis object is an x or y axis
653  !! @return .true. if an axis object is an x or y axis, optionally return a flag indicating which it is
654  function is_x_or_y_axis(this, x_or_y) &
655  result(rslt)
656  class(fmsdiagfullaxis_type), intent(in) :: this !< diag_axis obj
657  integer, optional, intent(inout) :: x_or_y !< returns is_x_axis if it is a x axis
658  !! is_y_axis if it is a y axis
659  logical :: rslt
660 
661  select case (trim(this%cart_name))
662  case ("X")
663  if (present(x_or_y)) x_or_y = is_x_axis
664  rslt = .true.
665  case ("Y")
666  if (present(x_or_y)) x_or_y = is_y_axis
667  rslt = .true.
668  case default
669  rslt = .false.
670  if (present(x_or_y)) x_or_y = diag_null
671  end select
672  end function is_x_or_y_axis
673 
674  !< @brief Get the global size of the axis, and the layout
675  !! It is assumed that this function is only called on "X" and "Y" axes
676  !! using the `is_x_or_y_axis` function from above
677  subroutine get_dim_size_layout(this, dim_size, layout)
678  class(fmsdiagfullaxis_type), intent(in) :: this !< diag_axis obj
679  integer, intent(out) :: dim_size !< Size of the dimension
680  integer, intent(out) :: layout !< Layout of the dimension
681 
682  integer :: nx, ny
683  integer :: layout_xy(2)
684 
685  select type (domain => this%axis_domain)
686  type is (diagdomain2d_t)
687  call mpp_get_global_domain(domain%Domain2, xsize=nx, ysize=ny)
688  call mpp_get_layout(domain%Domain2, layout_xy)
689 
690  if (this%cart_name .eq. "X") then
691  dim_size = nx
692  layout = layout_xy(1)
693  else if (this%cart_name .eq. "Y") then
694  dim_size = ny
695  layout = layout_xy(2)
696  endif
697  end select
698  end subroutine get_dim_size_layout
699 
700  !> @brief Get the set name of an axis object
701  !! @return the set name of an axis object
702  pure function get_set_name(this) &
703  result(rslt)
704  class(fmsdiagfullaxis_type), intent(in) :: this !< diag_axis obj
705  character(len=:), allocatable :: rslt
706 
707  rslt = this%set_name
708  end function get_set_name
709 
710  !> @brief Get the auxiliary name of an axis object
711  !! @return the auxiliary name of an axis object
712  pure function get_aux(this) &
713  result(rslt)
714  class(fmsdiagfullaxis_type), intent(in) :: this !< diag_axis obj
715  character(len=:), allocatable :: rslt
716 
717  rslt = this%aux
718  end function get_aux
719 
720  !> @brief Set the axis_id
721  subroutine set_axis_id(this, axis_id)
722  class(fmsdiagfullaxis_type), intent(inout) :: this !< diag_axis obj
723  integer, intent(in) :: axis_id !< Axis_id
724 
725  this%axis_id = axis_id
726 
727  end subroutine set_axis_id
728 
729  !> @brief Set the name and ids of the edges
730  subroutine set_edges(this, edges_name, edges_id)
731  class(fmsdiagfullaxis_type), intent(inout) :: this !< diag_axis obj
732  CHARACTER(len=*), intent(in) :: edges_name !< Name of the edges
733  integer, intent(in) :: edges_id !< Axis id of the edges
734 
735  !< Saving the name and the id of the edges axis because it will make it easier to use
736  !! downstream (i.e you need the edges name to write the attribute to the current axis,
737  !! and you need the edges id to add to the diag file object so that you can write the edges
738  !! to the file)
739  this%edges_name = edges_name
740  this%edges_id = edges_id
741  end subroutine set_edges
742 
743  !> @brief Determine if the subRegion is in the current PE.
744  !! If it is, determine the starting and ending indices of the current PE that belong to the subRegion
745  subroutine get_indices(this, compute_idx, corners_indices, starting_index, ending_index, need_to_define_axis)
746  class(fmsdiagfullaxis_type), intent(in) :: this !< diag_axis obj
747  integer, intent(in) :: compute_idx(:) !< Current PE's compute domain
748  class(*), intent(in) :: corners_indices(:) !< The indices of the corners of the subRegion
749  integer, intent(out) :: starting_index !< Starting index of the subRegion
750  !! for the current PE
751  integer, intent(out) :: ending_index !< Ending index of the subRegion
752  !! for the current PE
753  logical, intent(out) :: need_to_define_axis !< .true. if it is needed to define
754  !! an axis
755 
756  integer :: subregion_start !< Starting index of the subRegion
757  integer :: subregion_end !< Ending index of the subRegion
758 
759  !< Get the rectangular coordinates of the subRegion
760  !! If the subRegion is not rectangular, the points outside of the subRegion will be masked
761  !! out later
762  select type (corners_indices)
763  type is (integer(kind=i4_kind))
764  subregion_start = minval(corners_indices)
765  subregion_end = maxval(corners_indices)
766  end select
767 
768  !< Initiliaze the output
769  need_to_define_axis = .false.
770  starting_index = diag_null
771  ending_index = diag_null
772 
773  !< If the compute domain of the current PE is outisde of the range of sub_axis, return
774  if (compute_idx(1) < subregion_start .and. compute_idx(2) < subregion_start) return
775  if (compute_idx(1) > subregion_end .and. compute_idx(2) > subregion_end) return
776 
777  need_to_define_axis = .true.
778  if (compute_idx(1) >= subregion_start .and. compute_idx(2) >= subregion_end) then
779  !< In this case all the point of the current PE are inside the range of the sub_axis
780  starting_index = compute_idx(1)
781  ending_index = subregion_end
782  else if (compute_idx(1) >= subregion_start .and. compute_idx(2) <= subregion_end) then
783  !< In this case all the points of the current PE are valid up to the end point
784  starting_index = compute_idx(1)
785  ending_index = compute_idx(2)
786  else if (compute_idx(1) <= subregion_start .and. compute_idx(2) <= subregion_end) then
787  !< In this case all the points of the current PE are valid starting with t subregion_start
788  starting_index = subregion_start
789  ending_index = compute_idx(2)
790  else if (compute_idx(1) <= subregion_start .and. compute_idx(2) >= subregion_end) then
791  !< In this case only the points in the current PE ar valid
792  starting_index = subregion_start
793  ending_index = subregion_end
794  endif
795 
796  if (this%domain_position .ne. center) then
797  if (subregion_end - subregion_start + 1 .eq. 1) then
798  !< If your subregion consitsts of just 1 one, only include 1 PE
799  if (ending_index .eq. compute_idx(2)) need_to_define_axis = .false.
800  else
801  if (ending_index - starting_index + 1 .eq. 1) then
802  !< If the PEs section is only 1, only include 1 PE
803  if (starting_index .eq. compute_idx(2) .or. ending_index .eq. compute_idx(1)) &
804  need_to_define_axis = .false.
805  endif
806  endif
807  endif
808 
809  end subroutine get_indices
810 
811  !< Get the compute domain of the axis
812  subroutine get_compute_domain(this, compute_idx, need_to_define_axis, tile_number)
813  class(fmsdiagfullaxis_type), intent(in) :: this !< diag_axis obj
814  integer, intent(inout) :: compute_idx(:) !< Compute domain of the axis
815  logical, intent(out) :: need_to_define_axis !< .true. if it needed to define the axis
816  integer, optional, intent(in) :: tile_number !< The tile number of the axis
817 
818  !< Initialize the output
819  need_to_define_axis = .false.
820  compute_idx = diag_null
821 
822  if (.not. allocated(this%axis_domain)) then
823  !< If the axis is not domain decomposed, use the whole axis as the compute domain
824  if (this%cart_name .eq. "X" .or. this%cart_name .eq. "Y") then
825  compute_idx(1) = 1
826  compute_idx(2) = size(this%axis_data)
827  need_to_define_axis = .true.
828  endif
829  return
830  endif
831 
832  select type(domain => this%axis_domain)
833  type is (diagdomain2d_t)
834  if (present(tile_number)) then
835  !< If the tile number is present and the current PE is not on the tile, then there is no need
836  !! to define the axis
837  if (any(mpp_get_tile_id(domain%Domain2) .ne. tile_number)) then
838  need_to_define_axis = .false.
839  return
840  endif
841  endif
842 
843  !< Get the compute domain for the current PE if it is an "X" or "Y" axis
844  select case (this%cart_name)
845  case ("X")
846  call mpp_get_compute_domain(domain%Domain2, xbegin=compute_idx(1), xend=compute_idx(2), &
847  & position=this%domain_position)
848  need_to_define_axis = .true.
849  case ("Y")
850  call mpp_get_compute_domain(domain%Domain2, ybegin=compute_idx(1), yend=compute_idx(2), &
851  & position=this%domain_position)
852  need_to_define_axis = .true.
853  end select
854  end select
855 
856  end subroutine get_compute_domain
857 
858  !!!!!!!!!!!!!!!!!! SUB AXIS PROCEDURES !!!!!!!!!!!!!!!!!
859  !> @brief Fills in the information needed to define a subaxis
860  subroutine fill_subaxis(this, starting_index, ending_index, axis_id, parent_id, parent_axis_name, compute_idx, &
861  global_idx, zbounds, nz_subaxis)
862  class(fmsdiagsubaxis_type) , INTENT(INOUT) :: this !< diag_sub_axis obj
863  integer , intent(in) :: starting_index !< Starting index of the subRegion for the PE
864  integer , intent(in) :: ending_index !< Ending index of the subRegion for the PE
865  integer , intent(in) :: axis_id !< Axis id to assign to the subaxis
866  integer , intent(in) :: parent_id !< The id of the parent axis the subaxis belongs to
867  character(len=*) , intent(in) :: parent_axis_name !< Name of the parent_axis
868  integer , intent(in) :: compute_idx(2) !< Starting and ending index of
869  !! the axis's compute domain
870  integer, optional, intent(in) :: global_idx(2) !< Starting and ending index of
871  !! the axis's compute domain
872  real(kind=r4_kind), optional, intent(in) :: zbounds(2) !< Bounds of the z-axis
873  integer, optional, intent(in) :: nz_subaxis !< The number of z subaxis that have been defined
874  !! in the file
875 
876  integer :: nsubaxis !< The subaxis number in the axis name subXX
877  character(len=2) :: nsubaxis_char !< nsubaxis converted to a string
878 
879  nsubaxis = 1
880  if (present(nz_subaxis)) nsubaxis = nz_subaxis
881 
882  this%axis_id = axis_id
883  this%starting_index = starting_index
884  this%ending_index = ending_index
885  this%parent_axis_id = parent_id
886  write(nsubaxis_char, '(i2.2)') nsubaxis
887  this%subaxis_name = trim(parent_axis_name)//"_sub"//nsubaxis_char
888  this%compute_idx = compute_idx
889 
890  if (present(zbounds)) then
891  ! This is needed to avoid duplicating z sub axis!
892  allocate(this%zbounds(2))
893  this%zbounds = zbounds
894  endif
895 
896  if (present(global_idx)) then
897  ! This is needed for the "domain_decomposition" attribute which is needed for the combiner
898  allocate(this%global_idx(2))
899  this%global_idx = global_idx
900  endif
901  end subroutine fill_subaxis
902 
903  !> @brief Get the axis length of a subaxis
904  !> @return the axis length
905  function axis_length(this) &
906  result(res)
907  class(fmsdiagsubaxis_type) , INTENT(IN) :: this !< diag_sub_axis obj
908  integer :: res
909 
910  res = this%ending_index - this%starting_index + 1
911  end function
912 
913  !> @brief Accesses its member starting_index
914  !! @return a copy of the starting_index
915  function get_starting_index(this) result(indx)
916  class(fmsdiagsubaxis_type), intent(in) :: this !< diag_sub_axis object
917  integer :: indx !< Result to return
918  indx = this%starting_index
919  end function get_starting_index
920 
921  !> @brief Accesses its member ending_index
922  !! @return a copy of the ending_index
923  function get_ending_index(this) result(indx)
924  class(fmsdiagsubaxis_type), intent(in) :: this !< diag_sub_axis object
925  integer :: indx !< Result to return
926  indx = this%ending_index
927  end function get_ending_index
928 
929  !> @brief Accesses its member compute_indices
930  !! @return a copy of the ending_index
931  function get_compute_indices(this) result(indx)
932  class(fmsdiagsubaxis_type), intent(in) :: this !< diag_sub_axis object
933  integer :: indx(2) !< Result to return
934  indx = this%compute_idx
935  end function get_compute_indices
936 
937  !> @brief Get the ntiles in a domain
938  !> @return the number of tiles in a domain
939  function get_ntiles(this) &
940  result(ntiles)
941  class(diagdomain_t), INTENT(IN) :: this !< diag_axis obj
942 
943  integer :: ntiles
944 
945  select type (this)
946  type is (diagdomain2d_t)
947  ntiles = mpp_get_ntile_count(this%domain2)
948  end select
949  end function get_ntiles
950 
951  !> @brief Get the length of a 2D domain
952  !> @return Length of the 2D domain
953  function get_length(this, cart_axis, domain_position, global_length) &
954  result(length)
955  class(diagdomain_t), INTENT(IN) :: this !< diag_axis obj
956  character(len=*), INTENT(IN) :: cart_axis !< cart_axis of the axis
957  integer, INTENT(IN) :: domain_position !< Domain position (CENTER, NORTH, EAST)
958  integer, INTENT(IN) :: global_length !< global_length of the axis
959 
960  integer :: length
961 
962  select type (this)
963  type is(diagdomain2d_t)
964  if (trim(cart_axis) == "X") call mpp_get_compute_domain(this%Domain2, xsize=length, position=domain_position)
965  if (trim(cart_axis) == "Y") call mpp_get_compute_domain(this%Domain2, ysize=length, position=domain_position)
966  class default
967  !< If domain is 1D or UG, just set it to the global length
968  length = global_length
969  end select
970  end function get_length
971 
972  !!!!!!!!!!!!!!!!! FMS_DOMAIN PROCEDURES !!!!!!!!!!!!!!!!!
973 
974  !> @brief Set the axis domain
975  subroutine set_axis_domain(this, Domain, Domain2, DomainU)
976  class(diagdomain_t) :: this !< fms_domain obj
977  TYPE(domain1d), INTENT(in), OPTIONAL :: Domain !< 1d domain
978  TYPE(domain2d), INTENT(in), OPTIONAL :: Domain2 !< 2d domain
979  TYPE(domainug), INTENT(in), OPTIONAL :: DomainU !< Unstructured domain
980 
981  select type(this)
982  type is (diagdomain1d_t)
983  this%Domain = domain
984  type is (diagdomain2d_t)
985  this%Domain2 = domain2
986  type is (diagdomainug_t)
987  this%DomainUG = domainu
988  end select
989  end subroutine set_axis_domain
990 
991  !< @brief Allocates the array of axis/subaxis objects
992  !! @return true if there the aray of axis/subaxis objects is allocated
993  logical function fms_diag_axis_object_init(axis_array)
994  class(fmsdiagaxiscontainer_type) , allocatable, intent(inout) :: axis_array(:) !< Array of diag_axis
995 
996  if (allocated(axis_array)) call mpp_error(fatal, "The diag_axis containers is already allocated")
997  allocate(axis_array(max_axes))
998  !axis_array%axis_id = DIAG_NULL
999 
1000  fms_diag_axis_object_init = .true.
1001  end function fms_diag_axis_object_init
1002 
1003  !< @brief Deallocates the array of axis/subaxis objects
1004  !! @return false if the aray of axis/subaxis objects was allocated
1005  logical function fms_diag_axis_object_end(axis_array)
1006  class(fmsdiagaxiscontainer_type) , allocatable, intent(inout) :: axis_array(:) !< Array of diag_axis
1007 
1008  if (allocated(axis_array)) deallocate(axis_array)
1009  fms_diag_axis_object_end = .false.
1010 
1011  end function fms_diag_axis_object_end
1012 
1013  !< @brief Determine the axis name of an axis_object
1014  !! @return The name of the axis
1015  !! @note This function may be called from the field object (i.e. to determine the dimension names for io),
1016  !! The field object only contains the parent axis ids, because the subregion is defined in a per file basis,
1017  !! so the is_regional flag is needed so that the correct axis name can be used
1018  pure function get_axis_name(this, is_regional) &
1019  result(axis_name)
1020  class(fmsdiagaxis_type), intent(in) :: this !< Axis object
1021  logical, intent(in), optional :: is_regional !< Flag indicating if the axis is regional
1022 
1023  character(len=:), allocatable :: axis_name
1024 
1025  select type (this)
1026  type is (fmsdiagfullaxis_type)
1027  axis_name = this%axis_name
1028  if (present(is_regional)) then
1029  if (is_regional) then
1030  if (this%cart_name .eq. "X" .or. this%cart_name .eq. "Y") axis_name = axis_name//"_sub01"
1031  endif
1032  endif
1033  type is (fmsdiagsubaxis_type)
1034  axis_name = this%subaxis_name
1035  end select
1036  end function get_axis_name
1037 
1038  !< @brief Determine if the axis is a Z axis by looking at the cartesian name
1039  !! @return .True. if the axis is a Z axis
1040  pure logical function is_z_axis(this)
1041  class(fmsdiagaxis_type), intent(in) :: this !< Axis object
1042  is_z_axis = .false.
1043  select type (this)
1044  type is (fmsdiagfullaxis_type)
1045  if (this%cart_name .eq. "Z") is_z_axis = .true.
1046  end select
1047  end function
1048 
1049  !> @brief Check if a cart_name is valid and crashes if it isn't
1050  subroutine check_if_valid_cart_name(cart_name)
1051  character(len=*), intent(in) :: cart_name
1052 
1053  select case (cart_name)
1054  case ("X", "Y", "Z", "T", "U", "N")
1055  case default
1056  call mpp_error(fatal, "diag_axit_init: Invalid cart_name: "//cart_name//&
1057  "The acceptable values are X, Y, Z, T, U, N.")
1058  end select
1059  end subroutine check_if_valid_cart_name
1060 
1061  !> @brief Check if a domain_position is valid and crashes if it isn't
1062  subroutine check_if_valid_domain_position(domain_position)
1063  integer, INTENT(IN) :: domain_position
1064 
1065  select case (domain_position)
1066  case (center, north, east)
1067  case default
1068  call mpp_error(fatal, "diag_axit_init: Invalid domain_positon. &
1069  &The acceptable values are NORTH, EAST, CENTER")
1070  end select
1071  end subroutine check_if_valid_domain_position
1072 
1073  !> @brief Check if a direction is valid and crashes if it isn't
1074  subroutine check_if_valid_direction(direction)
1075  integer, INTENT(IN) :: direction
1076 
1077  select case(direction)
1078  case(-1, 0, 1)
1079  case default
1080  call mpp_error(fatal, "diag_axit_init: Invalid direction. &
1081  &The acceptable values are-1 0 1")
1082  end select
1083  end subroutine check_if_valid_direction
1084 
1085  !> @brief Loop through a variable's axis_id to determine and return the domain type and domain to use
1086  subroutine get_domain_and_domain_type(diag_axis, axis_id, domain_type, domain, var_name)
1087  class(fmsdiagaxiscontainer_type), target, intent(in) :: diag_axis(:) !< Array of diag_axis
1088  integer, INTENT(IN) :: axis_id(:) !< Array of axis ids
1089  integer, INTENT(OUT) :: domain_type !< fileobj_type to use
1090  CLASS(diagdomain_t), POINTER, INTENT(OUT) :: domain !< Domain
1091  character(len=*), INTENT(IN) :: var_name !< Name of the variable (for error messages)
1092 
1093  integer :: i !< For do loops
1094  integer :: j !< axis_id(i) (for less typing)
1095 
1096  domain_type = no_domain
1097  domain => null()
1098 
1099  do i = 1, size(axis_id)
1100  j = axis_id(i)
1101  select type (axis => diag_axis(j)%axis)
1102  type is (fmsdiagfullaxis_type)
1103  !< Check that all the axis are in the same domain
1104  if (domain_type .ne. axis%type_of_domain) then
1105  !< If they are different domains, one of them can be NO_DOMAIN
1106  !! i.e a variable can have axis that are domain decomposed (x,y) and an axis that isn't (z)
1107  if (domain_type .eq. no_domain .or. axis%type_of_domain .eq. no_domain ) then
1108  !< Update the domain_type and domain, if needed
1109  if ((axis%type_of_domain .eq. two_d_domain .and. size(axis_id) > 1) &
1110  & .or. axis%type_of_domain .eq. ug_domain) then
1111  domain_type = axis%type_of_domain
1112  domain => axis%axis_domain
1113  endif
1114  else
1115  call mpp_error(fatal, "The variable:"//trim(var_name)//" has axis that are not in the same domain")
1116  endif
1117  endif
1118  end select
1119  enddo
1120  end subroutine get_domain_and_domain_type
1121 
1122  !> @brief Fill in the subaxis object for a subRegion defined by index
1123  subroutine define_new_subaxis_index(parent_axis, subRegion, diag_axis, naxis, is_x_or_y, write_on_this_pe)
1124  class(fmsdiagaxiscontainer_type), target, intent(inout) :: diag_axis(:) !< Diag_axis object
1125  type(fmsdiagfullaxis_type), intent(inout) :: parent_axis !< axis object of the parent
1126  integer, intent(inout) :: naxis !< Number of axis registered
1127  type(subregion_type), intent(in) :: subregion !< SubRegion definition from the yaml
1128  integer, intent(in) :: is_x_or_y !< Flag indicating if it is
1129  !! a x or y axis
1130  logical, intent(out) :: write_on_this_pe !< .true. if the subregion
1131  !! is on this PE
1132  integer :: compute_idx(2) !< Indices of the compute domain
1133  integer :: global_idx(2) !< Indices of the "global" domain
1134  integer :: starting_index !< starting index of the subregion
1135  integer :: ending_index !< ending index of the subregion
1136 
1137  call parent_axis%get_compute_domain(compute_idx, write_on_this_pe, tile_number=subregion%tile)
1138  if (.not. write_on_this_pe) return
1139 
1140  !< Determine if the PE's compute domain is inside the subRegion
1141  !! If it is get the starting and ending indices for that PE
1142  call parent_axis%get_indices(compute_idx, subregion%corners(:,is_x_or_y), starting_index, ending_index, &
1143  write_on_this_pe)
1144 
1145  if (.not. write_on_this_pe) return
1146 
1147  select type(corners=> subregion%corners)
1148  type is (integer(kind=i4_kind))
1149  global_idx(1) = minval(corners(:,is_x_or_y))
1150  global_idx(2) = maxval(corners(:,is_x_or_y))
1151  end select
1152 
1153  !< If it made it to this point, the current PE is in the subRegion!
1154  call define_new_axis(diag_axis, parent_axis, naxis, parent_axis%axis_id, &
1155  starting_index, ending_index, compute_idx, global_idx)
1156 
1157  end subroutine define_new_subaxis_index
1158 
1159  !> @brief Fill in the subaxis object for a subRegion defined by lat lon
1160  subroutine define_new_subaxis_latlon(diag_axis, axis_ids, naxis, subRegion, is_cube_sphere, write_on_this_pe)
1161  class(fmsdiagaxiscontainer_type), target, intent(inout) :: diag_axis(:) !< Diag_axis object
1162  integer, INTENT(in) :: axis_ids(:) !< Array of axes_ids
1163  integer, intent(inout) :: naxis !< Number of axis registered
1164  type(subregion_type), intent(in) :: subregion !< SubRegion definition from the yaml
1165  logical, intent(in) :: is_cube_sphere !< .true. if this is a cubesphere
1166  logical, intent(out) :: write_on_this_pe !< .true. if the subregion
1167  !! is on this PE
1168 
1169  real :: lat(2) !< Starting and ending lattiude of the subRegion
1170  real :: lon(2) !< Starting and ending longitude or the subRegion
1171  integer :: lat_indices(2) !< Starting and ending latitude indices of the subRegion
1172  integer :: lon_indices(2) !< Starting and ending longitude indices of the subRegion
1173  integer :: compute_idx(2) !< Compute domain of the current axis
1174  integer :: starting_index(2) !< Starting index of the subRegion for the current PE for the "x" and "y"
1175  !! direction
1176  integer :: ending_index(2) !< Ending index of the subRegion for the current PE for the "x" and "y" direction
1177  logical :: need_to_define_axis(2) !< .true. if it is needed to define the subaxis for the "x" and "y" direction
1178  integer :: i !< For do loops
1179  integer :: parent_axis_ids(2) !< The axis id of the parent axis for the "x" and "y" direction
1180  logical :: is_x_y_axis !< .true. if the axis is x or y
1181  integer :: compute_idx_2(2, 2) !< Starting and ending indices of the compute domain for the "x" and "y" direction
1182  integer :: global_idx (2, 2) !< Starting and ending indices of the global domain for the "x" and "y" direction
1183 
1184  write_on_this_pe = .false.
1185  need_to_define_axis = .true.
1186  parent_axis_ids = diag_null
1187 
1188  !< Get the rectangular coordinates of the subRegion
1189  !! If the subRegion is not rectangular, the points outside of the subRegion will be masked
1190  !! out later
1191  select type (corners => subregion%corners)
1192  type is (real(kind=r4_kind))
1193  lon(1) = minval(corners(:,1))
1194  lon(2) = maxval(corners(:,1))
1195  lat(1) = minval(corners(:,2))
1196  lat(2) = maxval(corners(:,2))
1197  end select
1198 
1199  if_is_cube_sphere: if (is_cube_sphere) then
1200  !< Get the starting and ending indices of the subregion in the cubesphere relative to the global domain
1201  call get_local_indices_cubesphere(lat(1), lat(2), lon(1), lon(2),&
1202  & lon_indices(1), lon_indices(2), lat_indices(1), lat_indices(2))
1203  loop_over_axis_ids: do i = 1, size(axis_ids)
1204  select_axis_type: select type (parent_axis => diag_axis(axis_ids(i))%axis)
1205  type is (fmsdiagfullaxis_type)
1206  !< Get the PEs compute domain
1207  call parent_axis%get_compute_domain(compute_idx, is_x_y_axis)
1208 
1209  !< If this is not a "X" or "Y" axis go to the next axis
1210  if (.not. is_x_y_axis) cycle
1211 
1212  !< Determine if the PE's compute domain is inside the subRegion
1213  !! If it is get the starting and ending indices for that PE
1214  if (parent_axis%cart_name .eq. "X") then
1215  call parent_axis%get_indices(compute_idx, lon_indices, starting_index(1), ending_index(1), &
1216  need_to_define_axis(1))
1217  parent_axis_ids(1) = axis_ids(i)
1218  compute_idx_2(1,:) = compute_idx
1219  global_idx(1,:) = lon_indices
1220  else if (parent_axis%cart_name .eq. "Y") then
1221  call parent_axis%get_indices(compute_idx, lat_indices, starting_index(2), ending_index(2), &
1222  need_to_define_axis(2))
1223  parent_axis_ids(2) = axis_ids(i)
1224  compute_idx_2(2,:) = compute_idx
1225  global_idx(2,:) = lat_indices
1226  endif
1227  end select select_axis_type
1228  enddo loop_over_axis_ids
1229  else if_is_cube_sphere
1230  loop_over_axis_ids2: do i = 1, size(axis_ids)
1231  select type (parent_axis => diag_axis(axis_ids(i))%axis)
1232  type is (fmsdiagfullaxis_type)
1233  !< Get the PEs compute domain
1234  call parent_axis%get_compute_domain(compute_idx, is_x_y_axis)
1235 
1236  !< If this is not a "X" or "Y" axis go to the next axis
1237  if (.not. is_x_y_axis) cycle
1238 
1239  !< Get the starting and ending indices of the subregion relative to the global grid
1240  if (parent_axis%cart_name .eq. "X") then
1241  select type(adata=>parent_axis%axis_data)
1242  type is (real(kind=r8_kind))
1243  lon_indices(1) = nearest_index(real(lon(1), kind=r8_kind), adata)
1244  lon_indices(2) = nearest_index(real(lon(2), kind=r8_kind), adata)
1245  type is (real(kind=r4_kind))
1246  lon_indices(1) = nearest_index(real(lon(1), kind=r4_kind), adata)
1247  lon_indices(2) = nearest_index(real(lon(2), kind=r4_kind), adata)
1248  end select
1249  call parent_axis%get_indices(compute_idx, lon_indices, starting_index(1), ending_index(1), &
1250  need_to_define_axis(1))
1251  parent_axis_ids(1) = axis_ids(i)
1252  compute_idx_2(1,:) = compute_idx
1253  global_idx(1,:) = lon_indices
1254  else if (parent_axis%cart_name .eq. "Y") then
1255  select type(adata=>parent_axis%axis_data)
1256  type is (real(kind=r8_kind))
1257  lat_indices(1) = nearest_index(real(lat(1), kind=r8_kind), adata)
1258  lat_indices(2) = nearest_index(real(lat(2), kind=r8_kind), adata)
1259  type is (real(kind=r4_kind))
1260  lat_indices(1) = nearest_index(real(lat(1), kind=r4_kind), adata)
1261  lat_indices(2) = nearest_index(real(lat(2), kind=r4_kind), adata)
1262  end select
1263  call parent_axis%get_indices(compute_idx, lat_indices, starting_index(2), ending_index(2), &
1264  need_to_define_axis(2))
1265  parent_axis_ids(2) = axis_ids(i)
1266  compute_idx_2(2,:) = compute_idx
1267  global_idx(2,:) = lat_indices
1268  endif
1269  end select
1270  enddo loop_over_axis_ids2
1271  endif if_is_cube_sphere
1272 
1273  !< If the PE's compute is not inside the subRegion move to the next axis
1274  if (any(.not. need_to_define_axis )) return
1275 
1276  !< If it made it to this point, the current PE is in the subRegion!
1277  write_on_this_pe = .true.
1278 
1279  do i = 1, size(parent_axis_ids)
1280  if (parent_axis_ids(i) .eq. diag_null) cycle
1281  select type (parent_axis => diag_axis(parent_axis_ids(i))%axis)
1282  type is (fmsdiagfullaxis_type)
1283  call define_new_axis(diag_axis, parent_axis, naxis, parent_axis_ids(i), &
1284  starting_index(i), ending_index(i), compute_idx_2(i,:), global_idx(i,:))
1285  end select
1286  enddo
1287 
1288  end subroutine define_new_subaxis_latlon
1289 
1290  !> @brief Creates a new subaxis and fills it will all the information it needs
1291  subroutine define_new_axis(diag_axis, parent_axis, naxis, parent_id, &
1292  starting_index, ending_index, compute_idx, global_idx, new_axis_id, zbounds, &
1293  nz_subaxis)
1294 
1295  class(fmsdiagaxiscontainer_type), target, intent(inout) :: diag_axis(:) !< Diag_axis object
1296  class(fmsdiagfullaxis_type), intent(inout) :: parent_axis !< The parent axis
1297  integer, intent(inout) :: naxis !< The number of axis that
1298  !! have been defined
1299  integer, intent(in) :: parent_id !< Id of the parent axis
1300  integer, intent(in) :: starting_index !< PE's Starting index
1301  integer, intent(in) :: ending_index !< PE's Ending index
1302  integer, intent(in) :: compute_idx(2) !< Starting and ending index of
1303  !! the axis's compute domain
1304  integer, optional, intent(in) :: global_idx(2) !< Starting and ending index of
1305  !! the axis's global domain
1306  integer, optional, intent(out) :: new_axis_id !< Axis id of the axis this is creating
1307  real(kind=r4_kind), optional, intent(in) :: zbounds(2) !< Bounds of the Z axis
1308  integer, optional, intent(in) :: nz_subaxis !< The number of z subaxis that have
1309  !! been defined in the file
1310 
1311  naxis = naxis + 1 !< This is the axis id of the new axis!
1312 
1313  !< Add the axis_id of the new subaxis to the parent axis
1314  parent_axis%nsubaxis = parent_axis%nsubaxis + 1
1315  parent_axis%subaxis(parent_axis%nsubaxis) = naxis
1316 
1317  !< Allocate the new axis as a subaxis and fill it
1318  allocate(fmsdiagsubaxis_type :: diag_axis(naxis)%axis)
1319  diag_axis(naxis)%axis%axis_id = naxis
1320  if (present(new_axis_id)) new_axis_id = naxis
1321 
1322  select type (sub_axis => diag_axis(naxis)%axis)
1323  type is (fmsdiagsubaxis_type)
1324  call sub_axis%fill_subaxis(starting_index, ending_index, naxis, parent_id, &
1325  parent_axis%axis_name, compute_idx, global_idx=global_idx, zbounds=zbounds, nz_subaxis=nz_subaxis)
1326  end select
1327  end subroutine define_new_axis
1328 
1329  !< @brief Determine the parent_axis_id of a subaxis
1330  !! @return parent_axis_id if it is a subaxis and diag_null if is not a subaxis
1331  pure function get_parent_axis_id(this) &
1332  result(parent_axis_id)
1333 
1334  class(fmsdiagaxis_type), intent(in) :: this !< Axis Object
1335  integer :: parent_axis_id
1336 
1337  select type (this)
1338  type is (fmsdiagfullaxis_type)
1339  parent_axis_id = diag_null
1340  type is (fmsdiagsubaxis_type)
1341  parent_axis_id = this%parent_axis_id
1342  type is (fmsdiagdiurnalaxis_type)
1343  parent_axis_id = diag_null
1344  end select
1345 
1346  end function
1347 
1348  !< @brief Determine the most recent subaxis id in a diag_axis object
1349  !! @return the most recent subaxis id in a diag_axis object
1350  pure function get_subaxes_id(this) &
1351  result(sub_axis_id)
1352 
1353  class(fmsdiagaxis_type), intent(in) :: this !< Axis Object
1354  integer :: sub_axis_id
1355 
1356  sub_axis_id = this%axis_id
1357  select type (this)
1358  type is (fmsdiagfullaxis_type)
1359  if (this%cart_name .ne. "Z") sub_axis_id = this%subaxis(this%nsubaxis)
1360  end select
1361 
1362  end function
1363 
1364  !< @brief Parses the "compress" attribute to get the names of the two axis
1365  !! @return the names of the structured axis
1366  pure function parse_compress_att(compress_att) &
1367  result(axis_names)
1368  class(*), intent(in) :: compress_att(:) !< The compress attribute to parse
1369  character(len=120) :: axis_names(2)
1370 
1371  integer :: ios !< Errorcode after parsing the compress attribute
1372 
1373  select type (compress_att)
1374  type is (character(len=*))
1375  read(compress_att(1),*, iostat=ios) axis_names
1376  if (ios .ne. 0) axis_names = ""
1377  class default
1378  axis_names = ""
1379  end select
1380  end function parse_compress_att
1381 
1382  !< @brief Determine the axis id of a axis
1383  !! @return Axis id
1384  pure function get_axis_id_from_name(axis_name, diag_axis, naxis, set_name) &
1385  result(axis_id)
1386  class(fmsdiagaxiscontainer_type), intent(in) :: diag_axis(:) !< Array of axis object
1387  character(len=*), intent(in) :: axis_name !< Name of the axis
1388  integer, intent(in) :: naxis !< Number of axis that have been registered
1389  character(len=*), intent(in) :: set_name !< Name of the axis set
1390  integer :: axis_id
1391 
1392  integer :: i !< For do loops
1393 
1394  axis_id = diag_null
1395  do i = 1, naxis
1396  select type(axis => diag_axis(i)%axis)
1397  type is (fmsdiagfullaxis_type)
1398  if (trim(axis%axis_name) .eq. trim(axis_name)) then
1399  if (trim(axis%set_name) .eq. trim(set_name)) then
1400  axis_id = i
1401  return
1402  endif
1403  endif
1404  end select
1405  enddo
1406 
1407  end function get_axis_id_from_name
1408 
1409  !< @brief Get the number of diurnal samples for a diurnal axis
1410  !! @return The number of diurnal samples
1411  pure function get_diurnal_axis_samples(this) &
1412  result(n_diurnal_samples)
1413 
1414  class(fmsdiagdiurnalaxis_type), intent(in) :: this !< Axis Object
1415  integer :: n_diurnal_samples
1416 
1417  n_diurnal_samples = this%ndiurnal_samples
1418  end function get_diurnal_axis_samples
1419 
1420  !< @brief Writes out the metadata for a diurnal axis
1421  subroutine write_diurnal_metadata(this, fms2io_fileobj)
1422  class(fmsdiagdiurnalaxis_type), intent(in) :: this !< Diurnal axis Object
1423  class(fmsnetcdffile_t), intent(inout) :: fms2io_fileobj !< Fms2_io fileobj to write the data to
1424 
1425  call register_axis(fms2io_fileobj, this%axis_name, size(this%diurnal_data))
1426  call register_field(fms2io_fileobj, this%axis_name, pack_size_str, (/trim(this%axis_name)/))
1427  call register_variable_attribute(fms2io_fileobj, this%axis_name, "units", &
1428  &trim(this%units), str_len=len_trim(this%units))
1429  call register_variable_attribute(fms2io_fileobj, this%axis_name, "long_name", &
1430  &trim(this%long_name), str_len=len_trim(this%long_name))
1431  if (this%edges_id .ne. diag_null) &
1432  call register_variable_attribute(fms2io_fileobj, this%axis_name, "edges", &
1433  &trim(this%edges_name), str_len=len_trim(this%edges_name))
1434  end subroutine write_diurnal_metadata
1435 
1436  !> @brief Creates a new z subaxis to use
1437  subroutine create_new_z_subaxis(zbounds, var_axis_ids, diag_axis, naxis, file_axis_id, nfile_axis, nz_subaxis)
1438  real(kind=r4_kind), intent(in) :: zbounds(2) !< Bounds of the Z axis
1439  integer, intent(inout) :: var_axis_ids(:) !< The variable's axis_ids
1440  class(fmsdiagaxiscontainer_type), target, intent(inout) :: diag_axis(:) !< Array of diag_axis objects
1441  integer, intent(inout) :: naxis !< Number of axis that have been
1442  !! registered
1443  integer, intent(inout) :: file_axis_id(:) !< The file's axis_ids
1444  integer, intent(inout) :: nfile_axis !< Number of axis that have been
1445  !! defined in file
1446  integer, intent(inout) :: nz_subaxis !< The number of z subaxis currently
1447  !! defined in the file
1448 
1449  class(*), pointer :: zaxis_data(:) !< The data of the full zaxis
1450  integer :: subaxis_indices(2) !< The starting and ending indices of the subaxis relative to the full
1451  !! axis
1452  integer :: i !< For do loops
1453  integer :: subaxis_id !< The id of the new z subaxis
1454  logical :: axis_found !< Flag that indicated if the zsubaxis already exists
1455 
1456  !< Determine if the axis was already created
1457  axis_found = .false.
1458  do i = 1, nfile_axis
1459  select type (axis => diag_axis(file_axis_id(i))%axis)
1460  type is (fmsdiagsubaxis_type)
1461  if (axis%zbounds(1) .eq. zbounds(1) .and. axis%zbounds(2) .eq. zbounds(2)) then
1462  axis_found = .true.
1463  subaxis_id = file_axis_id(i)
1464  exit
1465  endif
1466  end select
1467  enddo
1468 
1469  !< Determine which of the variable's axis is the zaxis!
1470  do i = 1, size(var_axis_ids)
1471  select type (parent_axis => diag_axis(var_axis_ids(i))%axis)
1472  type is (fmsdiagfullaxis_type)
1473  if (parent_axis%cart_name .eq. "Z") then
1474  !< If the axis was previously defined set the var_axis_ids and leave
1475  if (axis_found) then
1476  var_axis_ids(i) = subaxis_id
1477  return
1478  endif
1479  zaxis_data => parent_axis%axis_data
1480 
1481  select type(zaxis_data)
1482  type is (real(kind=r4_kind))
1483  !TODO need to include the conversion to "real" because nearest_index doesn't take r4s and r8s
1484  subaxis_indices(1) = nearest_index(real(zbounds(1)), real(zaxis_data))
1485  subaxis_indices(2) = nearest_index(real(zbounds(2)), real(zaxis_data))
1486  type is (real(kind=r8_kind))
1487  subaxis_indices(1) = nearest_index(real(zbounds(1)), real(zaxis_data))
1488  subaxis_indices(2) = nearest_index(real(zbounds(2)), real(zaxis_data))
1489  end select
1490 
1491  nz_subaxis = nz_subaxis + 1
1492  call define_new_axis(diag_axis, parent_axis, naxis, parent_axis%axis_id, &
1493  &subaxis_indices(1), subaxis_indices(2), (/lbound(zaxis_data,1), ubound(zaxis_data,1)/), &
1494  &new_axis_id=subaxis_id, zbounds=zbounds, nz_subaxis=nz_subaxis)
1495  var_axis_ids(i) = subaxis_id
1496  return
1497  endif
1498  end select
1499  enddo
1500 
1501  end subroutine
1502 
1503  !> @brief Determine if the diag_axis(parent_axis_id) is the parent of diag_axis(axis_id)
1504  !! @return .True. if diag_axis(parent_axis_id) is the parent of diag_axis(axis_id)
1505  function is_parent_axis(axis_id, parent_axis_id, diag_axis) &
1506  result(rslt)
1507  integer, intent(in) :: axis_id !< Axis id to check
1508  integer, intent(in) :: parent_axis_id !< Axis id of the parent to check
1509  class(fmsdiagaxiscontainer_type), target, intent(in) :: diag_axis(:) !< Array of diag_axis objects
1510 
1511  logical :: rslt
1512 
1513  rslt = .false.
1514  select type(axis => diag_axis(axis_id)%axis)
1515  type is (fmsdiagsubaxis_type)
1516  if (axis%parent_axis_id .eq. parent_axis_id) rslt = .true.
1517  end select
1518  end function is_parent_axis
1519 
1520 #endif
1521 end module fms_diag_axis_object_mod
1522 !> @}
1523 ! close documentation grouping
integer, parameter direction_down
The axis points down if positive.
Definition: diag_data.F90:106
integer function get_base_minute()
gets the module variable base_minute
Definition: diag_data.F90:551
integer function get_base_year()
gets the module variable base_year
Definition: diag_data.F90:519
integer function get_base_hour()
gets the module variable base_hour
Definition: diag_data.F90:543
integer, parameter no_domain
Use the FmsNetcdfFile_t fileobj.
Definition: diag_data.F90:100
integer max_axis_attributes
Maximum number of user definable attributes per axis.
Definition: diag_data.F90:385
character(len=6) pack_size_str
Pack size as a string to be used in fms2_io register call set to "double" or "float".
Definition: diag_data.F90:409
integer max_axes
Maximum number of independent axes.
Definition: diag_data.F90:361
integer, parameter is_x_axis
integer indicating that it is a x axis
Definition: diag_data.F90:130
integer, parameter is_y_axis
integer indicating that it is a y axis
Definition: diag_data.F90:131
integer function get_base_day()
gets the module variable base_day
Definition: diag_data.F90:535
integer, parameter ug_domain
Use the FmsNetcdfUnstructuredDomainFile_t fileobj.
Definition: diag_data.F90:102
integer, parameter direction_up
The axis points up if positive.
Definition: diag_data.F90:105
integer function get_base_month()
gets the module variable base_month
Definition: diag_data.F90:527
integer function get_base_second()
gets the module variable base_second
Definition: diag_data.F90:559
integer, parameter two_d_domain
Use the FmsNetcdfDomainFile_t fileobj.
Definition: diag_data.F90:101
Attribute type for diagnostic fields.
Definition: diag_data.F90:157
Type to hold the attributes of the field/axis/file.
Definition: diag_data.F90:334
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
Add a dimension to a given file.
Definition: fms2_io.F90:189
Defines a new field within the given file Example usage:
Definition: fms2_io.F90:211
Write data to a defined field within a file Example usage:
Definition: fms2_io.F90:261
subroutine, public define_new_subaxis_latlon(diag_axis, axis_ids, naxis, subRegion, is_cube_sphere, write_on_this_pe)
Fill in the subaxis object for a subRegion defined by lat lon.
subroutine get_dim_size_layout(this, dim_size, layout)
pure character(len=:) function, allocatable get_axis_name(this, is_regional)
integer function get_length(this, cart_axis, domain_position, global_length)
Get the length of a 2D domain.
integer function get_ending_index(this)
Accesses its member ending_index.
logical function, public fms_diag_axis_object_end(axis_array)
subroutine check_if_valid_domain_position(domain_position)
Check if a domain_position is valid and crashes if it isn't.
logical function is_x_or_y_axis(this, x_or_y)
Determine if an axis object is an x or y axis.
pure logical function is_unstructured_grid(this)
subroutine set_axis_id(this, axis_id)
Set the axis_id.
pure logical function is_z_axis(this)
integer function axis_length(this)
Get the axis length of a subaxis.
pure integer function get_diurnal_axis_samples(this)
integer function get_ntiles(this)
Get the ntiles in a domain.
subroutine, public define_new_subaxis_index(parent_axis, subRegion, diag_axis, naxis, is_x_or_y, write_on_this_pe)
Fill in the subaxis object for a subRegion defined by index.
pure integer function, public get_axis_id_from_name(axis_name, diag_axis, naxis, set_name)
subroutine get_compute_domain(this, compute_idx, need_to_define_axis, tile_number)
subroutine, public define_diurnal_axis(diag_axis, naxis, n_diurnal_samples, is_edges)
Defined a new diurnal axis.
pure integer function, dimension(2) get_structured_axis(this)
pure logical function has_aux(this)
Determine if an axis object has an auxiliary name.
subroutine add_structured_axis_ids(this, axis_ids)
subroutine check_if_valid_cart_name(cart_name)
Check if a cart_name is valid and crashes if it isn't.
subroutine write_axis_metadata(this, fms2io_fileobj, edges_in_file, parent_axis)
Write the axis meta data to an open fileobj.
subroutine fill_subaxis(this, starting_index, ending_index, axis_id, parent_id, parent_axis_name, compute_idx, global_idx, zbounds, nz_subaxis)
Fills in the information needed to define a subaxis.
subroutine get_global_io_domain(this, global_io_index, use_collective_writes)
Get the starting and ending indices of the global io domain of the axis.
logical function, public fms_diag_axis_object_init(axis_array)
pure integer function get_subaxes_id(this)
integer function get_starting_index(this)
Accesses its member starting_index.
subroutine get_indices(this, compute_idx, corners_indices, starting_index, ending_index, need_to_define_axis)
Determine if the subRegion is in the current PE. If it is, determine the starting and ending indices ...
pure integer function get_parent_axis_id(this)
subroutine check_if_valid_direction(direction)
Check if a direction is valid and crashes if it isn't.
subroutine, public create_new_z_subaxis(zbounds, var_axis_ids, diag_axis, naxis, file_axis_id, nfile_axis, nz_subaxis)
Creates a new z subaxis to use.
subroutine register_diag_axis_obj(this, axis_name, axis_data, units, cart_name, long_name, direction, set_name, Domain, Domain2, DomainU, aux, req, tile_count, domain_position, axis_length)
Initialize the axis.
subroutine write_axis_data(this, fms2io_fileobj, parent_axis)
Write the axis data to an open fms2io_fileobj.
subroutine, public define_new_axis(diag_axis, parent_axis, naxis, parent_id, starting_index, ending_index, compute_idx, global_idx, new_axis_id, zbounds, nz_subaxis)
Creates a new subaxis and fills it will all the information it needs.
integer function get_axis_length(this)
Get the length of the axis.
pure character(len=:) function, allocatable get_set_name(this)
Get the set name of an axis object.
pure character(len=120) function, dimension(2), public parse_compress_att(compress_att)
subroutine set_edges(this, edges_name, edges_id)
Set the name and ids of the edges.
logical function, public is_parent_axis(axis_id, parent_axis_id, diag_axis)
Determine if the diag_axis(parent_axis_id) is the parent of diag_axis(axis_id)
subroutine add_axis_attribute(this, att_name, att_value)
Add an attribute to an axis.
subroutine write_diurnal_metadata(this, fms2io_fileobj)
pure logical function has_set_name(this)
Determine if an axis object has a set_name.
pure integer function get_edges_id(this)
integer function, dimension(2) get_compute_indices(this)
Accesses its member compute_indices.
subroutine, public get_domain_and_domain_type(diag_axis, axis_id, domain_type, domain, var_name)
Loop through a variable's axis_id to determine and return the domain type and domain to use.
subroutine set_axis_domain(this, Domain, Domain2, DomainU)
Set the axis domain.
pure character(len=:) function, allocatable get_aux(this)
Get the auxiliary name of an axis object.
integer function, dimension(size(domain%tile_id(:))) mpp_get_tile_id(domain)
Returns the tile_id on current pe.
integer function mpp_get_ntile_count(domain)
Returns number of tiles in mosaic.
type(domain2d) function, pointer mpp_get_io_domain(domain)
Set user stack size.
These routines retrieve the axis specifications associated with the compute domains....
These routines retrieve the axis specifications associated with the global domains....
Retrieve layout associated with a domain decomposition The 1D version of this call returns the number...
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 stdout()
This function returns the current standard fortran unit numbers for output.
Definition: mpp_util.inc:42
integer function mpp_pe()
Returns processor ID.
Definition: mpp_util.inc:406
Error handler.
Definition: mpp.F90:381
Type to hold the domain info for an axis This type was created to avoid having to send in "Domain",...
Type to hold the unstructured domain.
Type to hold the diagnostic axis description.
Type to hold the diag_axis (either subaxis or a full axis)
Type to hold the diagnostic axis description.
type to hold the sub region information about a file