FMS  2025.03
Flexible Modeling System
netcdf_io.F90
1 !***********************************************************************
2 !* GNU Lesser General Public License
3 !*
4 !* This file is part of the GFDL Flexible Modeling System (FMS).
5 !*
6 !* FMS is free software: you can redistribute it and/or modify it under
7 !* the terms of the GNU Lesser General Public License as published by
8 !* the Free Software Foundation, either version 3 of the License, or (at
9 !* your option) any later version.
10 !*
11 !* FMS is distributed in the hope that it will be useful, but WITHOUT
12 !* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
13 !* FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
14 !* for more details.
15 !*
16 !* You should have received a copy of the GNU Lesser General Public
17 !* License along with FMS. If not, see <http://www.gnu.org/licenses/>.
18 !***********************************************************************
19 !> @defgroup netcdf_io_mod netcdf_io_mod
20 !> @ingroup fms2_io
21 !> @brief Creates a basic netcdf type and routines to extend for other uses
22 !!
23 !> Handles type definitions and interfaces for netcdf I/O.
24 
25 !> @addtogroup netcdf_io_mod
26 !> @{
27 module netcdf_io_mod
28 #ifndef MAX_NUM_RESTART_VARS_
29 #define MAX_NUM_RESTART_VARS_ 250
30 #endif
31 use netcdf
32 use mpp_mod
33 use fms_io_utils_mod
34 use platform_mod
35 implicit none
36 private
37 
38 
39 !Module constants.
40 integer, parameter, public :: default_deflate_level = 0 !< The default (no compression) deflate level to use
41 integer, parameter :: variable_missing = -1
42 integer, parameter :: dimension_missing = -1
43 integer, parameter, public :: no_unlimited_dimension = -1 !> No unlimited dimension in file.
44 character(len=1), parameter :: missing_path = ""
45 integer, parameter :: missing_ncid = -1
46 integer, parameter :: missing_rank = -1
47 integer, parameter, public :: define_mode = 0
48 integer, parameter, public :: data_mode = 1
49 integer, parameter, public :: max_num_restart_vars = max_num_restart_vars_
50 integer, parameter, public :: unlimited = nf90_unlimited !> Wrapper to specify unlimited dimension.
51 integer, parameter :: dimension_not_found = 0
52 integer, parameter, public :: max_num_compressed_dims = 10 !> Maximum number of compressed
53  !! dimensions allowed.
54 integer, private :: fms2_ncchksz = -1 !< Chunksize (bytes) used in nc_open and nc_create
55 integer, private :: fms2_nc_format_param = -1 !< Netcdf format type param used in nc_create
56 character (len = 10), private :: fms2_nc_format !< Netcdf format type used in netcdf_file_open
57 integer, private :: fms2_header_buffer_val = -1 !< value used in NF__ENDDEF
58 integer, private :: fms2_deflate_level = default_deflate_level !< Netcdf deflate level to use in
59  !! nf90_def_var (integer between 1 to 9)
60 logical, private :: fms2_shuffle = .false. !< Flag indicating whether to use the netcdf shuffle filter
61 logical, private :: fms2_is_netcdf4 = .false. !< Flag indicating whether the default netcdf file format is netcdf4
62 
63 !> @}
64 
65 !> @brief information needed fr regional restart variables
66 !> @ingroup netcdf_io_mod
67 type, private :: bc_information
68  integer, dimension(:), allocatable :: indices !< Indices for the halo region for the variable
69  !! (starting x, ending x, starting y, ending y)
70  integer, dimension(:), allocatable :: global_size !< Size of the variable for each dimension
71  integer, dimension(:), allocatable :: pelist !< List of pelist that have the data for the variable
72  logical :: is_root_pe !< Flag indicating if this is the root_pe from the pelist
73  integer :: x_halo !< Number of halos in x
74  integer :: y_halo !< Number of halos in y
75  integer :: jshift !< Shift in the x axis (from center)
76  integer :: ishift !< Shift in the y axis (from center)
77  real(kind=r4_kind), dimension(:,:), allocatable :: globaldata2d_r4 !< 2d data pointer.
78  real(kind=r4_kind), dimension(:,:,:), allocatable :: globaldata3d_r4 !< 3d data pointer.
79  real(kind=r8_kind), dimension(:,:), allocatable :: globaldata2d_r8 !< 2d data pointer.
80  real(kind=r8_kind), dimension(:,:,:), allocatable :: globaldata3d_r8 !< 3d data pointer.
81  character(len=32) :: chksum !< The variable's checksum
82  logical :: data_on_file_root !< Flag indicating if the file root is part of the pelist that
83  !!contains data
84 endtype bc_information
85 
86 !> @brief Restart variable.
87 !> @ingroup netcdf_io_mod
88 type, private :: restartvariable_t
89  character(len=256) :: varname !< Variable name.
90  class(*), pointer :: data0d => null() !< Scalar data pointer.
91  class(*), dimension(:), pointer :: data1d => null() !< 1d data pointer.
92  class(*), dimension(:,:), pointer :: data2d => null() !< 2d data pointer.
93  class(*), dimension(:,:,:), pointer :: data3d => null() !< 3d data pointer.
94  class(*), dimension(:,:,:,:), pointer :: data4d => null() !< 4d data pointer.
95  class(*), dimension(:,:,:,:,:), pointer :: data5d => null() !< 5d data pointer.
96  logical :: was_read !< Flag to support legacy "query_initialized" feature, which
97  !! keeps track if a file was read.
98  logical :: is_bc_variable !< Flag indicating if variable is a bc_variable
99  type(bc_information) :: bc_info !< information about the boundary condition variable
100 endtype restartvariable_t
101 
102 !> @brief Compressed dimension.
103 !> @ingroup netcdf_io_mod
104 type, private :: compresseddimension_t
105  character(len=256) :: dimname !< Dimension name.
106  integer, dimension(:), allocatable :: npes_corner !< Array of starting
107  !! indices for each rank.
108  integer, dimension(:), allocatable :: npes_nelems !< Number of elements
109  !! associated with each
110  !! rank.
111  integer :: nelems !< Total size of the dimension.
112 endtype compresseddimension_t
113 
114 !> @brief information about the current dimensions for regional restart variables
115 !> @ingroup netcdf_io_mod
116 type, private :: dimension_information
117  integer, dimension(5) :: xlen !> The size of each unique x dimension
118  integer, dimension(5) :: ylen !> The size of each unique y dimension
119  integer, dimension(5) :: zlen !> The size of each unique z dimension
120  integer, dimension(3) :: cur_dim_len !> Number of unique:
121  !! cur_dim_len(1) : x dimensions
122  !! cur_dim_len(2) : y dimensions
123  !! cur_dim_len(3) : z dimensions
124 endtype dimension_information
125 
126 !> @brief Netcdf file type.
127 !> @ingroup netcdf_io_mod
128 type, public :: fmsnetcdffile_t
129  character(len=FMS_PATH_LEN) :: path !< File path.
130  logical :: is_readonly !< Flag telling if the file is readonly.
131  integer :: ncid !< Netcdf file id.
132  character(len=256) :: nc_format !< Netcdf file format.
133  logical :: is_netcdf4 !< Flag indicating if the netcdf file type is netcdf4
134  integer, dimension(:), allocatable :: pelist !< List of ranks who will
135  !! communicate.
136  integer :: io_root !< I/O root rank of the pelist.
137  logical :: is_root !< Flag telling if the current rank is the
138  !! I/O root.
139  logical :: is_restart !< Flag telling if the this file is a restart
140  !! file (that has internal pointers to data).
141  logical :: mode_is_append !! true if file is open in "append" mode
142  logical, allocatable :: is_open !< Allocated and set to true if opened.
143  type(restartvariable_t), dimension(:), allocatable :: restart_vars !< Array of registered
144  !! restart variables.
145  integer :: num_restart_vars !< Number of registered restart variables.
146  type(compresseddimension_t), dimension(:), allocatable :: compressed_dims !< "Compressed" dimension.
147  integer :: num_compressed_dims !< Number of compressed dimensions.
148  logical :: is_diskless !< Flag telling whether this is a diskless file.
149  character (len=20) :: time_name
150  type(dimension_information) :: bc_dimensions !<information about the current dimensions for regional
151  !! restart variables
152  logical :: use_collective = .false. !< Flag indicating if we should open the file for collective input
153  !! this should be set to .true. in the user application if they want
154  !! collective reads (put before open_file())
155  integer :: tile_comm=mpp_comm_null !< MPI communicator used for collective reads.
156  !! To be replaced with a real communicator at user request
157 
158 endtype fmsnetcdffile_t
159 
160 
161 !> @brief Range type for a netcdf variable.
162 !> @ingroup netcdf_io_mod
163 type, public :: valid_t
164  logical :: has_range !< Flag that's true if both min/max exist for a variable.
165  logical :: has_min !< Flag that's true if min exists for a variable.
166  logical :: has_max !< Flag that's true if max exists for a variable.
167  logical :: has_fill !< Flag that's true a user defined fill value.
168  logical :: has_missing !< Flag that's true a user defined missing value.
169  real(kind=r8_kind) :: fill_val !< Unpacked fill value for a variable.
170  real(kind=r8_kind) :: min_val !< Unpacked minimum value allowed for a variable.
171  real(kind=r8_kind) :: max_val !< Unpacked maximum value allowed for a variable.
172  real(kind=r8_kind) :: missing_val !< Unpacked missing value for a variable.
173 endtype valid_t
174 
175 
176 public :: netcdf_io_init
177 public :: netcdf_file_open
178 public :: netcdf_file_close
179 public :: netcdf_add_dimension
180 public :: netcdf_add_variable
182 public :: global_att_exists
183 public :: variable_att_exists
186 public :: get_global_attribute
187 public :: get_variable_attribute
188 public :: get_num_dimensions
189 public :: get_dimension_names
190 public :: dimension_exists
191 public :: is_dimension_unlimited
192 public :: get_dimension_size
193 public :: get_num_variables
194 public :: get_variable_names
195 public :: variable_exists
198 public :: get_variable_size
200 public :: netcdf_read_data
201 public :: netcdf_write_data
202 public :: compressed_write
203 public :: netcdf_save_restart
204 public :: netcdf_restore_state
205 public :: get_valid
206 public :: is_valid
208 public :: netcdf_file_open_wrap
209 public :: netcdf_file_close_wrap
210 public :: netcdf_add_variable_wrap
211 public :: netcdf_save_restart_wrap
212 public :: compressed_write_0d_wrap
213 public :: compressed_write_1d_wrap
214 public :: compressed_write_2d_wrap
215 public :: compressed_write_3d_wrap
216 public :: compressed_write_4d_wrap
217 public :: compressed_write_5d_wrap
218 public :: compressed_read_0d
219 public :: compressed_read_1d
220 public :: compressed_read_2d
221 public :: compressed_read_3d
222 public :: compressed_read_4d
223 public :: compressed_read_5d
235 public :: get_fill_value
236 public :: get_variable_sense
237 public :: get_variable_missing
238 public :: get_variable_units
239 public :: get_time_calendar
240 public :: is_registered_to_restart
241 public :: set_netcdf_mode
242 public :: check_netcdf_code
243 public :: check_if_open
244 public :: set_fileobj_time_name
245 public :: write_restart_bc
246 public :: read_restart_bc
247 public :: flush_file
248 
249 !> @ingroup netcdf_io_mod
251  module procedure netcdf_add_restart_variable_0d
252  module procedure netcdf_add_restart_variable_1d
253  module procedure netcdf_add_restart_variable_2d
254  module procedure netcdf_add_restart_variable_3d
255  module procedure netcdf_add_restart_variable_4d
256  module procedure netcdf_add_restart_variable_5d
257 end interface netcdf_add_restart_variable
258 
259 !> @ingroup netcdf_io_mod
261  module procedure netcdf_read_data_0d
262  module procedure netcdf_read_data_1d
263  module procedure netcdf_read_data_2d
264  module procedure netcdf_read_data_3d
265  module procedure netcdf_read_data_4d
266  module procedure netcdf_read_data_5d
267 end interface netcdf_read_data
268 
269 
270 !> @ingroup netcdf_io_mod
272  module procedure netcdf_write_data_0d
273  module procedure netcdf_write_data_1d
274  module procedure netcdf_write_data_2d
275  module procedure netcdf_write_data_3d
276  module procedure netcdf_write_data_4d
277  module procedure netcdf_write_data_5d
278 end interface netcdf_write_data
279 
280 
281 !> @ingroup netcdf_io_mod
283  module procedure compressed_write_0d
284  module procedure compressed_write_1d
285  module procedure compressed_write_2d
286  module procedure compressed_write_3d
287  module procedure compressed_write_4d
288  module procedure compressed_write_5d
289 end interface compressed_write
290 
291 
292 !> @ingroup netcdf_io_mod
294  module procedure register_global_attribute_0d
295  module procedure register_global_attribute_1d
296 end interface register_global_attribute
297 
298 
299 !> @ingroup netcdf_io_mod
301  module procedure register_variable_attribute_0d
302  module procedure register_variable_attribute_1d
303 end interface register_variable_attribute
304 
305 
306 !> @ingroup netcdf_io_mod
308  module procedure get_global_attribute_0d
309  module procedure get_global_attribute_1d
310 end interface get_global_attribute
311 
312 
313 !> @ingroup netcdf_io_mod
315  module procedure get_variable_attribute_0d
316  module procedure get_variable_attribute_1d
317 end interface get_variable_attribute
318 
319 
320 !> @ingroup netcdf_io_mod
322  module procedure scatter_data_bc_2d
323  module procedure scatter_data_bc_3d
324 end interface scatter_data_bc
325 
326 !> @ingroup netcdf_io_mod
327 interface gather_data_bc
328  module procedure gather_data_bc_2d
329  module procedure gather_data_bc_3d
330 end interface gather_data_bc
331 
332 !> The interface is needed to accomodate pgi because it can't handle class * and there was no other way around it
333 !> @ingroup netcdf_io_mod
334 interface is_valid
335  module procedure is_valid_r8
336  module procedure is_valid_r4
337 end interface is_valid
338 
339 !> @addtogroup netcdf_io_mod
340 !> @{
341 
342 contains
343 
344 !> @brief Accepts the namelist fms2_io_nml variables relevant to netcdf_io_mod
345 subroutine netcdf_io_init (chksz, header_buffer_val, netcdf_default_format, deflate_level, shuffle)
346 integer, intent(in) :: chksz !< Chunksize (bytes) used in nc_open and nc_create
347 character (len = 10), intent(in) :: netcdf_default_format !< Netcdf format type param used in nc_create
348 integer, intent(in) :: header_buffer_val !< Value used in NF__ENDDEF
349 integer, intent(in) :: deflate_level !< Netcdf deflate level to use in nf90_def_var
350  !! (integer between 1 to 9)
351 logical, intent(in) :: shuffle !< Flag indicating whether to use the netcdf shuffle filter
352 
353  fms2_ncchksz = chksz
354  fms2_deflate_level = deflate_level
355  fms2_shuffle = shuffle
356  fms2_is_netcdf4 = .false.
357  fms2_header_buffer_val = header_buffer_val
358  if (string_compare(netcdf_default_format, "64bit", .true.)) then
359  fms2_nc_format_param = nf90_64bit_offset
360  call string_copy(fms2_nc_format, "64bit")
361  elseif (string_compare(netcdf_default_format, "classic", .true.)) then
362  fms2_nc_format_param = nf90_classic_model
363  call string_copy(fms2_nc_format, "classic")
364  elseif (string_compare(netcdf_default_format, "netcdf4", .true.)) then
365  fms2_nc_format_param = nf90_netcdf4
366  fms2_is_netcdf4 = .true.
367  call string_copy(fms2_nc_format, "netcdf4")
368  else
369  call error("unrecognized netcdf file format "//trim(netcdf_default_format)// &
370  '. The acceptable values are "64bit", "classic", "netcdf4". Check fms2_io_nml: netcdf_default_format')
371  endif
372 
373 end subroutine netcdf_io_init
374 
375 !> @brief Check for errors returned by netcdf.
376 !! @internal
377 subroutine check_netcdf_code(err, msg)
378 
379  integer, intent(in) :: err !< Code returned by netcdf.
380  character(len=*), intent(in) :: msg !< Error message to be appended to the FATAL
381 
382  character(len=80) :: buf
383 
384  if (err .ne. nf90_noerr) then
385  buf = nf90_strerror(err)
386  call error(trim(buf)//": "//trim(msg))
387  endif
388 end subroutine check_netcdf_code
389 
390 
391 !> @brief Switch to the correct netcdf mode.
392 !! @internal
393 subroutine set_netcdf_mode(ncid, mode)
394 
395  integer, intent(in) :: ncid !< Netcdf file id.
396  integer, intent(in) :: mode !< Netcdf file mode.
397 
398  integer :: err
399 
400  if (mode .eq. define_mode) then
401  err = nf90_redef(ncid)
402  if (err .eq. nf90_eindefine .or. err .eq. nf90_eperm) then
403  return
404  endif
405  elseif (mode .eq. data_mode) then
406  if (fms2_header_buffer_val == -1) call error("set_netcdf_mode: fms2_header_buffer_val not set, call fms2_io_init")
407  err = nf90_enddef(ncid, h_minfree=fms2_header_buffer_val)
408  if (err .eq. nf90_enotindefine .or. err .eq. nf90_eperm) then
409  return
410  endif
411  else
412  call error("mode must be either define_mode or data_mode.")
413  endif
414  call check_netcdf_code(err, "set_netcdf_mode")
415 end subroutine set_netcdf_mode
416 
417 
418 !> @brief Get the id of a dimension from its name.
419 !! @return Dimension id, or dimension_missing if it doesn't exist.
420 !! @internal
421 function get_dimension_id(ncid, dimension_name, msg, allow_failure) &
422  result(dimid)
423 
424  integer, intent(in) :: ncid !< Netcdf file id.
425  character(len=*), intent(in) :: dimension_name !< Dimension name.
426  character(len=*), intent(in) :: msg !< Error message
427  logical, intent(in), optional :: allow_failure !< Flag that prevents
428  !! crash if dimension
429  !! does not exist.
430 
431  integer :: dimid
432 
433  integer :: err
434 
435  err = nf90_inq_dimid(ncid, trim(dimension_name), dimid)
436  if (present(allow_failure)) then
437  if (allow_failure .and. err .eq. nf90_ebaddim) then
438  dimid = dimension_missing
439  return
440  endif
441  endif
442  call check_netcdf_code(err, msg)
443 end function get_dimension_id
444 
445 
446 !> @brief Get the id of a variable from its name.
447 !! @return Variable id, or variable_missing if it doesn't exist.
448 !! @internal
449 function get_variable_id(ncid, variable_name, msg, allow_failure) &
450  result(varid)
451 
452  integer, intent(in) :: ncid !< Netcdf file object.
453  character(len=*), intent(in) :: variable_name !< Variable name.
454  character(len=*), intent(in) :: msg !< Error message
455  logical, intent(in), optional :: allow_failure !< Flag that prevents
456  !! crash if variable does
457  !! not exist.
458 
459  integer :: varid
460 
461  integer :: err
462 
463  err = nf90_inq_varid(ncid, trim(variable_name), varid)
464  if (present(allow_failure)) then
465  if (allow_failure .and. err .eq. nf90_enotvar) then
466  varid = variable_missing
467  return
468  endif
469  endif
470  call check_netcdf_code(err, msg)
471 end function get_variable_id
472 
473 
474 !> @brief Determine if an attribute exists.
475 !! @return Flag telling if the attribute exists.
476 !! @internal
477 function attribute_exists(ncid, varid, attribute_name, msg) &
478  result(att_exists)
479 
480  integer, intent(in) :: ncid !< Netcdf file id.
481  integer, intent(in) :: varid !< Variable id.
482  character(len=*), intent(in) :: attribute_name !< Attribute name.
483  character(len=*), intent(in), optional :: msg !< Error message
484 
485  logical :: att_exists
486 
487  integer :: err
488 
489  err = nf90_inquire_attribute(ncid, varid, trim(attribute_name))
490  if (err .eq. nf90_enotatt) then
491  att_exists = .false.
492  else
493  call check_netcdf_code(err, msg)
494  att_exists = .true.
495  endif
496 end function attribute_exists
497 
498 
499 !> @brief Get the type of a netcdf attribute.
500 !! @return The netcdf type of the attribute.
501 !! @internal
502 function get_attribute_type(ncid, varid, attname, msg) &
503  result(xtype)
504 
505  integer, intent(in) :: ncid !< Netcdf file id.
506  integer, intent(in) :: varid !< Variable id.
507  character(len=*), intent(in) :: attname !< Attribute name.
508  character(len=*), intent(in), optional :: msg !< Error message
509 
510  integer :: xtype
511 
512  integer :: err
513 
514  err = nf90_inquire_attribute(ncid, varid, attname, xtype=xtype)
515  call check_netcdf_code(err, msg)
516 end function get_attribute_type
517 
518 
519 !> @brief Get the type of a netcdf variable.
520 !! @return The netcdf type of the variable.
521 !! @internal
522 function get_variable_type(ncid, varid, msg) &
523  result(xtype)
524 
525  integer, intent(in) :: ncid !< Netcdf file id.
526  integer, intent(in) :: varid !< Variable id.
527  character(len=*), intent(in), optional :: msg !< Error message to append to netcdf error code
528 
529  integer :: xtype
530 
531  integer :: err
532 
533  err = nf90_inquire_variable(ncid, varid, xtype=xtype)
534  call check_netcdf_code(err, msg)
535 end function get_variable_type
536 
537 
538 !> @brief Open a netcdf file.
539 !! @return .true. if open succeeds, or else .false.
540 function netcdf_file_open(fileobj, path, mode, nc_format, pelist, is_restart, dont_add_res_to_filename) &
541  result(success)
542 
543  class(fmsnetcdffile_t), intent(inout) :: fileobj !< File object.
544  character(len=*), intent(in) :: path !< File path.
545  character(len=*), intent(in) :: mode !< File mode. Allowed values are:
546  !! "read", "append", "write", or
547  !! "overwrite".
548  character(len=*), intent(in), optional :: nc_format !< Netcdf format that
549  !! new files are written
550  !! as. Allowed values
551  !! are: "64bit", "classic",
552  !! or "netcdf4". Defaults to
553  !! "64bit". This overwrites
554  !! the value set in the fms2io
555  !! namelist
556  integer, dimension(:), intent(in), optional :: pelist !< List of ranks associated
557  !! with this file. If not
558  !! provided, only the current
559  !! rank will be able to
560  !! act on the file.
561  logical, intent(in), optional :: is_restart !< Flag telling if this file
562  !! is a restart file. Defaults
563  !! to false.
564  logical, intent(in), optional :: dont_add_res_to_filename !< Flag indicating not to add
565  !! ".res" to the filename
566  logical :: success
567 
568  integer :: nc_format_param
569  integer :: err
570  integer :: netcdf4 !< Query the file for the _IsNetcdf4 global attribute in the event
571  !! that the open for collective reads fails
572  character(len=FMS_PATH_LEN) :: buf !< File path with .res in the filename if it is a restart
573  character(len=FMS_PATH_LEN) :: buf2 !< File path with the filename appendix if there is one
574  logical :: is_res
575  logical :: dont_add_res !< flag indicated to not add ".res" to the filename
576 
577  if (allocated(fileobj%is_open)) then
578  if (fileobj%is_open) then
579  success = .true.
580  return
581  endif
582  endif
583  !< Only add ".res" to the file path if is_restart is set to true
584  !! and dont_add_res_to_filename is set to false.
585  is_res = .false.
586  if (present(is_restart)) then
587  is_res = is_restart
588  endif
589  dont_add_res = .false.
590  if (present(dont_add_res_to_filename)) then
591  dont_add_res = dont_add_res_to_filename
592  endif
593 
594  if (is_res .and. .not. dont_add_res) then
595  call restart_filepath_mangle(buf, trim(path))
596  else
597  call string_copy(buf, trim(path))
598  endif
599 
600  !< If it is a restart add the filename_appendix to the filename
601  if (is_res) then
602  call get_instance_filename(trim(buf), buf2)
603  else
604  call string_copy(buf2, trim(buf))
605  endif
606 
607  !Check if the file exists.
608  success = .true.
609  if (string_compare(mode, "read", .true.) .or. string_compare(mode, "append", .true.)) then
610  success = file_exists(buf2)
611  if (.not. success) then
612  return
613  endif
614  endif
615 
616  !Store properties in the derived type.
617  call string_copy(fileobj%path, trim(buf2))
618  if (present(pelist)) then
619  allocate(fileobj%pelist(size(pelist)))
620  fileobj%pelist(:) = pelist(:)
621  else
622  allocate(fileobj%pelist(1))
623  fileobj%pelist(1) = mpp_pe()
624  endif
625  fileobj%io_root = fileobj%pelist(1)
626  fileobj%is_root = mpp_pe() .eq. fileobj%io_root
627 
628  fileobj%is_netcdf4 = .false.
629  if (fms2_ncchksz == -1) call error("netcdf_file_open:: fms2_ncchksz not set, call fms2_io_init")
630  if (fms2_nc_format_param == -1) call error("netcdf_file_open:: fms2_nc_format_param not set, call fms2_io_init")
631 
632  if (present(nc_format)) then
633  if (string_compare(nc_format, "64bit", .true.)) then
634  nc_format_param = nf90_64bit_offset
635  elseif (string_compare(nc_format, "classic", .true.)) then
636  nc_format_param = nf90_classic_model
637  elseif (string_compare(nc_format, "netcdf4", .true.)) then
638  fileobj%is_netcdf4 = .true.
639  nc_format_param = nf90_netcdf4
640  else
641  call error("unrecognized netcdf file format: '"//trim(nc_format)//"' for file:"//trim(fileobj%path)//&
642  &"Check your open_file call, the acceptable values are 64bit, classic, netcdf4")
643  endif
644  call string_copy(fileobj%nc_format, nc_format)
645  else
646  call string_copy(fileobj%nc_format, trim(fms2_nc_format))
647  nc_format_param = fms2_nc_format_param
648  fileobj%is_netcdf4 = fms2_is_netcdf4
649  endif
650 
651  !Open the file with netcdf if this rank is the I/O root.
652  if (fileobj%is_root .and. .not.(fileobj%use_collective)) then
653  if (string_compare(mode, "read", .true.)) then
654  err = nf90_open(trim(fileobj%path), nf90_nowrite, fileobj%ncid, chunksize=fms2_ncchksz)
655  elseif (string_compare(mode, "append", .true.)) then
656  err = nf90_open(trim(fileobj%path), nf90_write, fileobj%ncid, chunksize=fms2_ncchksz)
657  elseif (string_compare(mode, "write", .true.)) then
658  err = nf90_create(trim(fileobj%path), ior(nf90_noclobber, nc_format_param), fileobj%ncid, chunksize=fms2_ncchksz)
659  elseif (string_compare(mode,"overwrite",.true.)) then
660  err = nf90_create(trim(fileobj%path), ior(nf90_clobber, nc_format_param), fileobj%ncid, chunksize=fms2_ncchksz)
661  else
662  call error("unrecognized file mode: '"//trim(mode)//"' for file:"//trim(fileobj%path)//&
663  &"Check your open_file call, the acceptable values are read, append, write, overwrite")
664  endif
665  call check_netcdf_code(err, "netcdf_file_open:"//trim(fileobj%path))
666  elseif(fileobj%use_collective .and. (fileobj%tile_comm /= mpp_comm_null)) then
667  if(string_compare(mode, "read", .true.)) then
668  ! Open the file for collective reads if the user requested that treatment in their application.
669  ! NetCDF does not have the ability to specify collective I/O at the file basis
670  ! so we must activate each variable in netcdf_read_data_2d() and netcdf_read_data_3d()
671  err = nf90_open(trim(fileobj%path), ior(nf90_nowrite, nf90_mpiio), fileobj%ncid, &
672  comm=fileobj%tile_comm, info=mpp_info_null)
673  if(err /= nf90_noerr) then
674  err = nf90_open(trim(fileobj%path), nf90_nowrite, fileobj%ncid)
675  err = nf90_get_att(fileobj%ncid, nf90_global, "_IsNetcdf4", netcdf4)
676  err = nf90_close(fileobj%ncid)
677  if(netcdf4 /= 1) then
678  call mpp_error(note,"netcdf_file_open: Open for collective read failed because the file is not &
679  netCDF-4 format."// &
680  " Falling back to parallel independent for file "// trim(fileobj%path))
681  fileobj%use_collective = .false.
682  fileobj%tile_comm = mpp_comm_null
683  else
684 #ifdef NO_NC_PARALLEL4
685  call mpp_error(fatal, "netCDF was not build with HDF5 parallel I/O features, "//&
686  "so collective netcdf io is not allowed. Please turn collective read off for file "//&
687  trim(fileobj%path))
688 #endif
689  endif
690  err = nf90_open(trim(fileobj%path), nf90_nowrite, fileobj%ncid, chunksize=fms2_ncchksz)
691  endif
692  elseif (string_compare(mode, "write", .true.)) then
693  call mpp_error(fatal,"netcdf_file_open: Attempt to create a file for collective write"// &
694  " This feature is not implemented"// trim(fileobj%path))
695  !err = nf90_create(trim(fileobj%path), ior(nf90_noclobber, nc_format_param), fileobj%ncid, &
696  ! comm=fileobj%tile_comm, info=MPP_INFO_NULL)
697  elseif (string_compare(mode,"overwrite",.true.)) then
698  call mpp_error(fatal,"netcdf_file_open: Attempt to create a file for collective overwrite"// &
699  " This feature is not implemented"// trim(fileobj%path))
700  !err = nf90_create(trim(fileobj%path), ior(nf90_clobber, nc_format_param), fileobj%ncid, &
701  ! comm=fileobj%tile_comm, info=MPP_INFO_NULL)
702  else
703  call error("unrecognized file mode: '"//trim(mode)//"' for file:"//trim(fileobj%path)//&
704  &"Check your open_file call, the acceptable values are read, append, write, overwrite")
705  endif
706  call check_netcdf_code(err, "netcdf_file_open:"//trim(fileobj%path))
707  else
708  fileobj%ncid = missing_ncid
709  endif
710 
711  fileobj%is_diskless = .false.
712 
713  !Allocate memory.
714  fileobj%is_restart = is_res
715  if (fileobj%is_restart) then
716  allocate(fileobj%restart_vars(max_num_restart_vars))
717  fileobj%num_restart_vars = 0
718  endif
719  fileobj%is_readonly = string_compare(mode, "read", .true.)
720  fileobj%mode_is_append = string_compare(mode, "append", .true.)
721  allocate(fileobj%compressed_dims(max_num_compressed_dims))
722  fileobj%num_compressed_dims = 0
723  ! Set the is_open flag to true for this file object.
724  if (.not.allocated(fileobj%is_open)) allocate(fileobj%is_open)
725  fileobj%is_open = .true.
726 
727  fileobj%bc_dimensions%xlen = 0
728  fileobj%bc_dimensions%ylen = 0
729  fileobj%bc_dimensions%zlen = 0
730  fileobj%bc_dimensions%cur_dim_len = 0
731 
732 end function netcdf_file_open
733 
734 
735 !> @brief Close a netcdf file.
736 subroutine netcdf_file_close(fileobj)
737 
738  class(fmsnetcdffile_t),intent(inout) :: fileobj !< File object.
739 
740  integer :: err
741  integer :: i
742 
743  if (fileobj%is_root) then
744  err = nf90_close(fileobj%ncid)
745  call check_netcdf_code(err, "netcdf_file_close:"//trim(fileobj%path))
746  endif
747  if (allocated(fileobj%is_open)) fileobj%is_open = .false.
748  fileobj%path = missing_path
749  fileobj%ncid = missing_ncid
750  if (allocated(fileobj%pelist)) then
751  deallocate(fileobj%pelist)
752  endif
753  fileobj%io_root = missing_rank
754  fileobj%is_root = .false.
755  if (allocated(fileobj%restart_vars)) then
756  deallocate(fileobj%restart_vars)
757  endif
758  fileobj%is_restart = .false.
759  fileobj%num_restart_vars = 0
760  do i = 1, fileobj%num_compressed_dims
761  if (allocated(fileobj%compressed_dims(i)%npes_corner)) then
762  deallocate(fileobj%compressed_dims(i)%npes_corner)
763  endif
764  if (allocated(fileobj%compressed_dims(i)%npes_nelems)) then
765  deallocate(fileobj%compressed_dims(i)%npes_nelems)
766  endif
767  enddo
768  if (allocated(fileobj%compressed_dims)) then
769  deallocate(fileobj%compressed_dims)
770  endif
771 end subroutine netcdf_file_close
772 
773 
774 !> @brief Get the index of a compressed dimension in a file object.
775 !! @return Index of the compressed dimension.
776 !! @internal
777 function get_compressed_dimension_index(fileobj, dim_name) &
778  result(dindex)
779 
780  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
781  character(len=*), intent(in) :: dim_name !< Dimension name.
782 
783  integer :: dindex
784  integer :: i
785 
786  dindex = dimension_not_found
787  do i = 1, fileobj%num_compressed_dims
788  if (string_compare(fileobj%compressed_dims(i)%dimname, dim_name)) then
789  dindex = i
790  return
791  endif
792  enddo
794 
795 
796 !> @brief Add a compressed dimension to a file object.
797 !! @internal
798 subroutine append_compressed_dimension(fileobj, dim_name, npes_corner, &
799  npes_nelems)
800 
801  class(fmsnetcdffile_t), intent(inout) :: fileobj !< File object.
802  character(len=*), intent(in) :: dim_name !< Dimension name.
803  integer, dimension(:), intent(in) :: npes_corner !< Array of starting
804  !! indices for each rank.
805  integer, dimension(:), intent(in) :: npes_nelems !< Number of elements
806  !! associated with each
807  !! rank.
808 
809  integer :: n
810 
811  if (get_compressed_dimension_index(fileobj, dim_name) .ne. dimension_not_found) then
812  call error("dimension "//trim(dim_name)//" already registered" &
813  //" to file "//trim(fileobj%path)//".")
814  endif
815  fileobj%num_compressed_dims = fileobj%num_compressed_dims + 1
816  n = fileobj%num_compressed_dims
817  if (n .gt. max_num_compressed_dims) then
818  call error("number of compressed dimensions exceeds limit.")
819  endif
820  call string_copy(fileobj%compressed_dims(n)%dimname, dim_name)
821  if (size(npes_corner) .ne. size(fileobj%pelist) .or. &
822  size(npes_nelems) .ne. size(fileobj%pelist)) then
823  call error("incorrect size for input npes_corner or npes_nelems arrays.")
824  endif
825  allocate(fileobj%compressed_dims(n)%npes_corner(size(fileobj%pelist)))
826  fileobj%compressed_dims(n)%npes_corner(:) = npes_corner(:)
827  allocate(fileobj%compressed_dims(n)%npes_nelems(size(fileobj%pelist)))
828  fileobj%compressed_dims(n)%npes_nelems(:) = npes_nelems(:)
829  fileobj%compressed_dims(n)%nelems = sum(fileobj%compressed_dims(n)%npes_nelems)
830 end subroutine append_compressed_dimension
831 
832 !> @brief Add a "compressed" unlimited dimension to a netcdf file.
833 !! @note Here compressed means that every rank has a different dimension_length
834 !! compressed. This was written specifically for the icebergs restarts.
835 subroutine register_unlimited_compressed_axis(fileobj, dimension_name, dimension_length)
836  class(fmsnetcdffile_t), intent(inout) :: fileobj !< File object.
837  character(len=*), intent(in) :: dimension_name !< Dimension name.
838  integer, intent(in) :: dimension_length !< Dimension length for the current rank
839 
840  integer, dimension(:), allocatable :: npes_start !< The starting index of the dimension for each of the PEs
841  integer, dimension(:), allocatable :: npes_count !< The size of the dimension for each of the PEs
842  integer :: i !< For do loops
843  integer :: err !< Netcdf error
844  integer :: dimid !< Netcdf id for the dimension
845 
846  !Gather all local dimension lengths on the I/O root pe.
847  allocate(npes_start(size(fileobj%pelist)))
848  allocate(npes_count(size(fileobj%pelist)))
849 
850  call mpp_gather((/dimension_length/),npes_count,pelist=fileobj%pelist)
851 
852  npes_start(1) = 1
853  do i = 1, size(fileobj%pelist)-1
854  npes_start(i+1) = npes_start(i) + npes_count(i)
855  enddo
856  call append_compressed_dimension(fileobj, dimension_name, npes_start, &
857  npes_count)
858 
859  if (fileobj%is_root .and. .not. fileobj%is_readonly) then
860  call set_netcdf_mode(fileobj%ncid, define_mode)
861  err = nf90_def_dim(fileobj%ncid, trim(dimension_name), unlimited, dimid)
862  call check_netcdf_code(err, "Netcdf_add_dimension: file:"//trim(fileobj%path)//" dimension name:"// &
863  & trim(dimension_name))
864  endif
866 
867 !> @brief Add a dimension to a file.
868 subroutine netcdf_add_dimension(fileobj, dimension_name, dimension_length, &
869  is_compressed)
870 
871  class(fmsnetcdffile_t), intent(inout) :: fileobj !< File object.
872  character(len=*), intent(in) :: dimension_name !< Dimension name.
873  integer, intent(in) :: dimension_length !< Dimension length.
874  logical, intent(in), optional :: is_compressed !< Changes the meaning of dim_len from
875  !! referring to the total size of the
876  !! dimension (when false) to the local
877  !! size for the current rank (when true).
878 
879  integer :: dim_len
880  integer, dimension(:), allocatable :: npes_start
881  integer, dimension(:), allocatable :: npes_count
882  integer :: i
883  integer :: err
884  integer :: dimid
885 
886  dim_len = dimension_length
887  if (present(is_compressed)) then
888  if (is_compressed) then
889  !Gather all local dimension lengths on the I/O root pe.
890  allocate(npes_start(size(fileobj%pelist)))
891  allocate(npes_count(size(fileobj%pelist)))
892  do i = 1, size(fileobj%pelist)
893  if (fileobj%pelist(i) .eq. mpp_pe()) then
894  npes_count(i) = dim_len
895  else
896  call mpp_recv(npes_count(i), fileobj%pelist(i), block=.false.)
897  call mpp_send(dim_len, fileobj%pelist(i))
898  endif
899  enddo
900  call mpp_sync_self(check=event_recv)
901  call mpp_sync_self(check=event_send)
902  npes_start(1) = 1
903  do i = 1, size(fileobj%pelist)-1
904  npes_start(i+1) = npes_start(i) + npes_count(i)
905  enddo
906  call append_compressed_dimension(fileobj, dimension_name, npes_start, &
907  npes_count)
908  dim_len = sum(npes_count)
909  endif
910  endif
911  if (fileobj%is_root .and. .not. fileobj%is_readonly) then
912  call set_netcdf_mode(fileobj%ncid, define_mode)
913  err = nf90_def_dim(fileobj%ncid, trim(dimension_name), dim_len, dimid)
914  call check_netcdf_code(err, "Netcdf_add_dimension: file:"//trim(fileobj%path)//" dimension name:"// &
915  & trim(dimension_name))
916  endif
917 end subroutine netcdf_add_dimension
918 
919 
920 !> @brief Add a compressed dimension.
921 subroutine register_compressed_dimension(fileobj, dimension_name, &
922  npes_corner, npes_nelems)
923 
924  class(fmsnetcdffile_t), intent(inout) :: fileobj !< File object.
925  character(len=*), intent(in) :: dimension_name !< Dimension name.
926  integer, dimension(:), intent(in) :: npes_corner !< Array of starting
927  !! indices for each rank.
928  integer, dimension(:), intent(in) :: npes_nelems !< Number of elements
929  !! associated with each
930  !! rank.
931 
932  integer :: dsize
933  integer :: fdim_size
934 
935  call append_compressed_dimension(fileobj, dimension_name, npes_corner, npes_nelems)
936  dsize = sum(npes_nelems)
937  if (fileobj%is_readonly) then
938  call get_dimension_size(fileobj, dimension_name, fdim_size, broadcast=.true.)
939  if (fdim_size .ne. dsize) then
940  call error("dimension "//trim(dimension_name)//" does not match" &
941  //" the size of the associated compressed axis.")
942  endif
943  else
944  call netcdf_add_dimension(fileobj, dimension_name, dsize)
945  endif
946 end subroutine register_compressed_dimension
947 
948 
949 !> @brief Add a variable to a file.
950 subroutine netcdf_add_variable(fileobj, variable_name, variable_type, dimensions, chunksizes)
951 
952  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
953  character(len=*), intent(in) :: variable_name !< Variable name.
954  character(len=*), intent(in) :: variable_type !< Variable type. Allowed
955  !! values are: "char", "int", "int64",
956  !! "float", or "double".
957  character(len=*), dimension(:), intent(in), optional :: dimensions !< Dimension names.
958  integer, optional, intent(in) :: chunksizes(:) !< netcdf chunksize to use for this variable
959  !! This feature is only
960  !! available for netcdf4 files
961  integer :: err
962  integer, dimension(:), allocatable :: dimids
963  integer :: vtype
964  integer :: varid
965  integer :: i
966  character(len=200) :: append_error_msg !< Msg to be appended to FATAL error message
967 
968  append_error_msg = "netcdf_add_variable: file:"//trim(fileobj%path)//" variable:"//trim(variable_name)
969 
970  if (fileobj%is_root) then
971  call set_netcdf_mode(fileobj%ncid, define_mode)
972  if (string_compare(variable_type, "int", .true.)) then
973  vtype = nf90_int
974  elseif (string_compare(variable_type, "int64", .true.)) then
975  if ( .not. fileobj%is_netcdf4) call error(trim(fileobj%path)//&
976  &": 64 bit integers are only supported with 'netcdf4' file format"//&
977  &". Set netcdf_default_format='netcdf4' in the fms2_io namelist OR "//&
978  &"add nc_format='netcdf4' to your open_file call")
979  vtype = nf90_int64
980  elseif (string_compare(variable_type, "float", .true.)) then
981  vtype = nf90_float
982  elseif (string_compare(variable_type, "double", .true.)) then
983  vtype = nf90_double
984  elseif (string_compare(variable_type, "char", .true.)) then
985  vtype = nf90_char
986  if (.not. present(dimensions)) then
987  call error("String variables require a string length dimension:"//trim(append_error_msg))
988  endif
989  else
990  call error("Unsupported variable type:"//trim(append_error_msg))
991  endif
992  if (present(dimensions)) then
993  allocate(dimids(size(dimensions)))
994  do i = 1, size(dimids)
995  dimids(i) = get_dimension_id(fileobj%ncid, trim(dimensions(i)),msg=append_error_msg)
996  enddo
997  if (fileobj%is_netcdf4) then
998  err = nf90_def_var(fileobj%ncid, trim(variable_name), vtype, dimids, varid, &
999  &deflate_level=fms2_deflate_level, shuffle=fms2_shuffle, chunksizes=chunksizes)
1000  else
1001  if (fms2_deflate_level .ne. default_deflate_level .or. fms2_shuffle .or. present(chunksizes)) &
1002  &call mpp_error(note,"Not able to use deflate_level or chunksizes if not using netcdf4"// &
1003  & " ignoring them")
1004  err = nf90_def_var(fileobj%ncid, trim(variable_name), vtype, dimids, varid)
1005  endif
1006  deallocate(dimids)
1007  else
1008  err = nf90_def_var(fileobj%ncid, trim(variable_name), vtype, varid)
1009  endif
1010  call check_netcdf_code(err, append_error_msg)
1011  endif
1012 end subroutine netcdf_add_variable
1013 
1014 
1015 !> @brief Given a compressed variable, get the index of the compressed
1016 !! dimension.
1017 !! @return Index of the compressed dimension or dimension_not_found.
1018 function get_variable_compressed_dimension_index(fileobj, variable_name, broadcast) &
1019  result(compressed_dimension_index)
1020 
1021  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
1022  character(len=*), intent(in) :: variable_name !< Variable name.
1023  logical, intent(in), optional :: broadcast !< Flag controlling whether or
1024  !! not the data will be
1025  !! broadcasted to non
1026  !! "I/O root" ranks.
1027  !! The broadcast will be done
1028  !! by default.
1029  integer, dimension(2) :: compressed_dimension_index
1030 
1031  integer :: ndims
1032  character(len=nf90_max_name), dimension(:), allocatable :: dim_names
1033  integer :: i
1034  integer :: j
1035 
1036  compressed_dimension_index = dimension_not_found
1037  if (fileobj%is_root) then
1038  ndims = get_variable_num_dimensions(fileobj, variable_name, broadcast=.false.)
1039  if (ndims .gt. 0) then
1040  allocate(dim_names(ndims))
1041  call get_variable_dimension_names(fileobj, variable_name, dim_names, broadcast=.false.)
1042  do i = 1, size(dim_names)
1043  j = get_compressed_dimension_index(fileobj,dim_names(i))
1044  if (j .ne. dimension_not_found) then
1045  compressed_dimension_index(1) = i
1046  compressed_dimension_index(2) = j
1047  exit
1048  endif
1049  enddo
1050  deallocate(dim_names)
1051  endif
1052  endif
1053  if (present(broadcast)) then
1054  if (.not. broadcast) then
1055  return
1056  endif
1057  endif
1058  call mpp_broadcast(compressed_dimension_index(1), fileobj%io_root, pelist=fileobj%pelist)
1059  call mpp_broadcast(compressed_dimension_index(2), fileobj%io_root, pelist=fileobj%pelist)
1061 
1062 
1063 !> @brief Add a restart variable to a FmsNetcdfFile_t type.
1064 !! @internal
1065 subroutine add_restart_var_to_array(fileobj, variable_name)
1066 
1067  class(fmsnetcdffile_t), intent(inout) :: fileobj !< Netcdf file object.
1068  character(len=*), intent(in) :: variable_name !< Variable name.
1069 
1070  integer :: i
1071 
1072  if (.not. fileobj%is_restart) then
1073  call error("file "//trim(fileobj%path)//" is not a restart file.")
1074  endif
1075  do i = 1, fileobj%num_restart_vars
1076  if (string_compare(fileobj%restart_vars(i)%varname, variable_name, .true.)) then
1077  call error("variable "//trim(variable_name)//" has already" &
1078  //" been added to restart file "//trim(fileobj%path)//".")
1079  endif
1080  enddo
1081  fileobj%num_restart_vars = fileobj%num_restart_vars + 1
1082  if (fileobj%num_restart_vars .gt. max_num_restart_vars) then
1083  call error("Number of restart variables exceeds limit.")
1084  endif
1085  call string_copy(fileobj%restart_vars(fileobj%num_restart_vars)%varname, &
1086  variable_name)
1087 end subroutine add_restart_var_to_array
1088 
1089 
1090 !> @brief Loop through registered restart variables and write them to
1091 !! a netcdf file.
1092 subroutine netcdf_save_restart(fileobj, unlim_dim_level)
1093 
1094  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
1095  integer, intent(in), optional :: unlim_dim_level !< Unlimited dimension
1096  !! level.
1097 
1098  integer :: i
1099 
1100  if (.not. fileobj%is_restart) then
1101  call error("write_restart:: file "//trim(fileobj%path)//" is not a restart file. &
1102  &Be sure the file was opened with is_restart=.true.")
1103  endif
1104  do i = 1, fileobj%num_restart_vars
1105  if (associated(fileobj%restart_vars(i)%data0d)) then
1106  call compressed_write(fileobj, fileobj%restart_vars(i)%varname, &
1107  fileobj%restart_vars(i)%data0d, &
1108  unlim_dim_level=unlim_dim_level)
1109  elseif (associated(fileobj%restart_vars(i)%data1d)) then
1110  call compressed_write(fileobj, fileobj%restart_vars(i)%varname, &
1111  fileobj%restart_vars(i)%data1d, &
1112  unlim_dim_level=unlim_dim_level)
1113  elseif (associated(fileobj%restart_vars(i)%data2d)) then
1114  call compressed_write(fileobj, fileobj%restart_vars(i)%varname, &
1115  fileobj%restart_vars(i)%data2d, &
1116  unlim_dim_level=unlim_dim_level)
1117  elseif (associated(fileobj%restart_vars(i)%data3d)) then
1118  call compressed_write(fileobj, fileobj%restart_vars(i)%varname, &
1119  fileobj%restart_vars(i)%data3d, &
1120  unlim_dim_level=unlim_dim_level)
1121  elseif (associated(fileobj%restart_vars(i)%data4d)) then
1122  call compressed_write(fileobj, fileobj%restart_vars(i)%varname, &
1123  fileobj%restart_vars(i)%data4d, &
1124  unlim_dim_level=unlim_dim_level)
1125  else
1126  call error("this branch should not be reached.")
1127  endif
1128  enddo
1129 end subroutine netcdf_save_restart
1130 
1131 
1132 !> @brief Loop through registered restart variables and read them from
1133 !! a netcdf file.
1134 subroutine netcdf_restore_state(fileobj, unlim_dim_level)
1135 
1136  type(fmsnetcdffile_t), intent(inout) :: fileobj !< File object.
1137  integer, intent(in), optional :: unlim_dim_level !< Unlimited dimension
1138  !! level.
1139 
1140  integer :: i
1141 
1142  if (.not. fileobj%is_restart) then
1143  call error("read_restart:: file "//trim(fileobj%path)//" is not a restart file. &
1144  &Be sure the file was opened with is_restart=.true.")
1145  endif
1146  do i = 1, fileobj%num_restart_vars
1147  if (associated(fileobj%restart_vars(i)%data0d)) then
1148  call netcdf_read_data(fileobj, fileobj%restart_vars(i)%varname, &
1149  fileobj%restart_vars(i)%data0d, &
1150  unlim_dim_level=unlim_dim_level, &
1151  broadcast=.true.)
1152  elseif (associated(fileobj%restart_vars(i)%data1d)) then
1153  call netcdf_read_data(fileobj, fileobj%restart_vars(i)%varname, &
1154  fileobj%restart_vars(i)%data1d, &
1155  unlim_dim_level=unlim_dim_level, &
1156  broadcast=.true.)
1157  elseif (associated(fileobj%restart_vars(i)%data2d)) then
1158  call netcdf_read_data(fileobj, fileobj%restart_vars(i)%varname, &
1159  fileobj%restart_vars(i)%data2d, &
1160  unlim_dim_level=unlim_dim_level, &
1161  broadcast=.true.)
1162  elseif (associated(fileobj%restart_vars(i)%data3d)) then
1163  call netcdf_read_data(fileobj, fileobj%restart_vars(i)%varname, &
1164  fileobj%restart_vars(i)%data3d, &
1165  unlim_dim_level=unlim_dim_level, &
1166  broadcast=.true.)
1167  elseif (associated(fileobj%restart_vars(i)%data4d)) then
1168  call netcdf_read_data(fileobj, fileobj%restart_vars(i)%varname, &
1169  fileobj%restart_vars(i)%data4d, &
1170  unlim_dim_level=unlim_dim_level, &
1171  broadcast=.true.)
1172  else
1173  call error("this branch should not be reached.")
1174  endif
1175  enddo
1176 end subroutine netcdf_restore_state
1177 
1178 
1179 !> @brief Determine if a global attribute exists.
1180 !! @return Flag telling if a global attribute exists.
1181 function global_att_exists(fileobj, attribute_name, broadcast) &
1182  result(att_exists)
1183 
1184  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
1185  character(len=*), intent(in) :: attribute_name !< Attribute name.
1186  logical, intent(in), optional :: broadcast !< Flag controlling whether or
1187  !! not the data will be
1188  !! broadcasted to non
1189  !! "I/O root" ranks.
1190  !! The broadcast will be done
1191  !! by default.
1192  logical :: att_exists
1193 
1194  if (fileobj%is_root) then
1195  att_exists = attribute_exists(fileobj%ncid, nf90_global, trim(attribute_name), &
1196  & msg="global_att_exists: file:"//trim(fileobj%path)//" attribute name:"//trim(attribute_name))
1197  endif
1198  if (present(broadcast)) then
1199  if (.not. broadcast) then
1200  return
1201  endif
1202  endif
1203  call mpp_broadcast(att_exists, fileobj%io_root, pelist=fileobj%pelist)
1204 end function global_att_exists
1205 
1206 
1207 !> @brief Determine if a variable's attribute exists.
1208 !! @return Flag telling if the variable's attribute exists.
1209 function variable_att_exists(fileobj, variable_name, attribute_name, &
1210  broadcast) &
1211  result(att_exists)
1212 
1213  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
1214  character(len=*), intent(in) :: variable_name !< Variable name.
1215  character(len=*), intent(in) :: attribute_name !< Attribute name.
1216  logical, intent(in), optional :: broadcast !< Flag controlling whether or
1217  !! not the data will be
1218  !! broadcasted to non
1219  !! "I/O root" ranks.
1220  !! The broadcast will be done
1221  !! by default.
1222  logical :: att_exists
1223 
1224  integer :: varid
1225 
1226  att_exists = .false.
1227  if (fileobj%is_root) then
1228  varid = get_variable_id(fileobj%ncid, trim(variable_name), &
1229  & msg="variable_att_exists: file:"//trim(fileobj%path)//"- variable:"//&
1230  &trim(variable_name))
1231  att_exists = attribute_exists(fileobj%ncid, varid, trim(attribute_name), &
1232  &msg="variable_att_exists: file:"//trim(fileobj%path)//" variable:"//trim(variable_name)//&
1233  &" attribute name:"//trim(attribute_name))
1234  endif
1235  if (present(broadcast)) then
1236  if (.not. broadcast) then
1237  return
1238  endif
1239  endif
1240  call mpp_broadcast(att_exists, fileobj%io_root, pelist=fileobj%pelist)
1241 end function variable_att_exists
1242 
1243 
1244 !> @brief Determine the number of dimensions in a file.
1245 !! @return The number of dimensions in the file.
1246 function get_num_dimensions(fileobj, broadcast) &
1247  result(ndims)
1248 
1249  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
1250  logical, intent(in), optional :: broadcast !< Flag controlling whether or
1251  !! not the data will be
1252  !! broadcasted to non
1253  !! "I/O root" ranks.
1254  !! The broadcast will be done
1255  !! by default.
1256  integer :: ndims
1257 
1258  integer :: err
1259 
1260  if (fileobj%is_root) then
1261  err = nf90_inquire(fileobj%ncid, ndimensions=ndims)
1262  call check_netcdf_code(err, "get_num_dimensions: file:"//trim(fileobj%path))
1263  endif
1264  if (present(broadcast)) then
1265  if (.not. broadcast) then
1266  return
1267  endif
1268  endif
1269  call mpp_broadcast(ndims, fileobj%io_root, pelist=fileobj%pelist)
1270 end function get_num_dimensions
1271 
1272 
1273 !> @brief Get the names of the dimensions in a file.
1274 subroutine get_dimension_names(fileobj, names, broadcast)
1275 
1276  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
1277  character(len=*), dimension(:), intent(inout) :: names !< Names of the
1278  !! dimensions.
1279  logical, intent(in), optional :: broadcast !< Flag controlling whether or
1280  !! not the data will be
1281  !! broadcasted to non
1282  !! "I/O root" ranks.
1283  !! The broadcast will be done
1284  !! by default.
1285 
1286  integer :: ndims
1287  integer :: i
1288  integer :: err
1289 
1290  if (fileobj%is_root) then
1291  ndims = get_num_dimensions(fileobj, broadcast=.false.)
1292  if (ndims .gt. 0) then
1293  if (size(names) .ne. ndims) then
1294  call error("'names' has to be the same size of the number of dimensions. &
1295  &Check your get_dimension_names call for file "//trim(fileobj%path))
1296  endif
1297  else
1298  call error("get_dimension_names: the file "//trim(fileobj%path)//" does not have any dimensions")
1299  endif
1300  names(:) = ""
1301  do i = 1, ndims
1302  err = nf90_inquire_dimension(fileobj%ncid, i, name=names(i))
1303  call check_netcdf_code(err, "get_dimension_names: file:"//trim(fileobj%path))
1304  enddo
1305  endif
1306  if (present(broadcast)) then
1307  if (.not. broadcast) then
1308  return
1309  endif
1310  endif
1311  call mpp_broadcast(ndims, fileobj%io_root, pelist=fileobj%pelist)
1312  if (.not. fileobj%is_root) then
1313  if (ndims .gt. 0) then
1314  if (size(names) .ne. ndims) then
1315  call error("'names' has to be the same size of the number of dimensions. &
1316  &Check your get_dimension_names call for file "//trim(fileobj%path))
1317  endif
1318  else
1319  call error("get_dimension_names: the file "//trim(fileobj%path)//" does not have any dimensions")
1320  endif
1321  names(:) = ""
1322  endif
1323  call mpp_broadcast(names, len(names(ndims)), fileobj%io_root, &
1324  pelist=fileobj%pelist)
1325 end subroutine get_dimension_names
1326 
1327 
1328 !> @brief Determine if a dimension exists.
1329 !! @return Flag telling if the dimension exists.
1330 function dimension_exists(fileobj, dimension_name, broadcast) &
1331  result(dim_exists)
1332 
1333  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
1334  character(len=*), intent(in) :: dimension_name !< Dimension name.
1335  logical, intent(in), optional :: broadcast !< Flag controlling whether or
1336  !! not the data will be
1337  !! broadcasted to non
1338  !! "I/O root" ranks.
1339  !! The broadcast will be done
1340  !! by default.
1341  logical :: dim_exists
1342 
1343  integer :: dimid
1344 
1345  if (fileobj%is_root) then
1346  dimid = get_dimension_id(fileobj%ncid, trim(dimension_name), &
1347  msg="dimension_exists: file:"//trim(fileobj%path)//" dimension:"//trim(dimension_name), &
1348  allow_failure=.true.)
1349  if (dimid .eq. dimension_missing) then
1350  dim_exists = .false.
1351  else
1352  dim_exists = .true.
1353  endif
1354  endif
1355  if (present(broadcast)) then
1356  if (.not. broadcast) then
1357  return
1358  endif
1359  endif
1360  call mpp_broadcast(dim_exists, fileobj%io_root, pelist=fileobj%pelist)
1361 end function dimension_exists
1362 
1363 
1364 !> @brief Determine where or not the dimension is unlimited.
1365 !! @return True if the dimension is unlimited, or else false.
1366 function is_dimension_unlimited(fileobj, dimension_name, broadcast) &
1367  result(is_unlimited)
1368 
1369  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
1370  character(len=*), intent(in) :: dimension_name !< Dimension name.
1371  logical, intent(in), optional :: broadcast !< Flag controlling whether or
1372  !! not the data will be
1373  !! broadcasted to non
1374  !! "I/O root" ranks.
1375  !! The broadcast will be done
1376  !! by default.
1377  logical :: is_unlimited
1378 
1379  character(len=200) :: append_error_msg !< Msg to be appended to FATAL error message
1380  integer :: dimid
1381  integer :: err
1382  integer :: ulim_dimid
1383 
1384  if (fileobj%is_root) then
1385  append_error_msg="is_dimension_unlimited: file:"//trim(fileobj%path)//&
1386  & " dimension_name:"//trim(dimension_name)
1387  dimid = get_dimension_id(fileobj%ncid, trim(dimension_name), msg=append_error_msg)
1388  err = nf90_inquire(fileobj%ncid, unlimiteddimid=ulim_dimid)
1389  call check_netcdf_code(err, append_error_msg)
1390  is_unlimited = dimid .eq. ulim_dimid
1391  endif
1392  if (present(broadcast)) then
1393  if (.not. broadcast) then
1394  return
1395  endif
1396  endif
1397  call mpp_broadcast(is_unlimited, fileobj%io_root, pelist=fileobj%pelist)
1398 end function is_dimension_unlimited
1399 
1400 
1401 !> @brief Get the name of the unlimited dimension.
1402 subroutine get_unlimited_dimension_name(fileobj, dimension_name, broadcast)
1403 
1404  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
1405  character(len=*), intent(out) :: dimension_name !< Dimension name.
1406  logical, intent(in), optional :: broadcast !< Flag controlling whether or
1407  !! not the data will be
1408  !! broadcasted to non
1409  !! "I/O root" ranks.
1410  !! The broadcast will be done
1411  !! by default.
1412 
1413  integer :: err
1414  integer :: dimid
1415  character(len=nf90_max_name), dimension(1) :: buffer
1416 
1417  dimension_name = ""
1418  if (fileobj%is_root) then
1419  err = nf90_inquire(fileobj%ncid, unlimiteddimid=dimid)
1420  call check_netcdf_code(err, "get_unlimited_dimension_name: file:"//trim(fileobj%path))
1421  err = nf90_inquire_dimension(fileobj%ncid, dimid, dimension_name)
1422  call check_netcdf_code(err, "get_unlimited_dimension_name: file:"//trim(fileobj%path))
1423  call string_copy(buffer(1), dimension_name)
1424  endif
1425  if (present(broadcast)) then
1426  if (.not. broadcast) then
1427  return
1428  endif
1429  endif
1430  call mpp_broadcast(buffer, nf90_max_name, fileobj%io_root, &
1431  pelist=fileobj%pelist)
1432  call string_copy(dimension_name, buffer(1))
1433 end subroutine get_unlimited_dimension_name
1434 
1435 
1436 !> @brief Get the length of a dimension.
1437 subroutine get_dimension_size(fileobj, dimension_name, dim_size, broadcast)
1438 
1439  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
1440  character(len=*), intent(in) :: dimension_name !< Dimension name.
1441  integer, intent(inout) :: dim_size !< Dimension size.
1442  logical, intent(in), optional :: broadcast !< Flag controlling whether or
1443  !! not the data will be
1444  !! broadcasted to non
1445  !! "I/O root" ranks.
1446  !! The broadcast will be done
1447  !! by default.
1448 
1449  integer :: dimid
1450  integer :: err
1451  character(len=200) :: append_error_msg !< Msg to be appended to FATAL error message
1452 
1453  if (fileobj%is_root) then
1454  append_error_msg = "get_dimension_size: file:"//trim(fileobj%path)//" dimension_name: "//trim(dimension_name)
1455  dimid = get_dimension_id(fileobj%ncid, trim(dimension_name), msg=append_error_msg)
1456  err = nf90_inquire_dimension(fileobj%ncid, dimid, len=dim_size)
1457  call check_netcdf_code(err, append_error_msg)
1458  endif
1459  if (present(broadcast)) then
1460  if (.not. broadcast) then
1461  return
1462  endif
1463  endif
1464  call mpp_broadcast(dim_size, fileobj%io_root, pelist=fileobj%pelist)
1465 end subroutine get_dimension_size
1466 
1467 
1468 !> @brief Determine the number of variables in a file.
1469 !! @return The number of variables in the file.
1470 function get_num_variables(fileobj, broadcast) &
1471  result(nvars)
1472 
1473  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
1474  logical, intent(in), optional :: broadcast !< Flag controlling whether or
1475  !! not the data will be
1476  !! broadcasted to non
1477  !! "I/O root" ranks.
1478  !! The broadcast will be done
1479  !! by default.
1480  integer :: nvars
1481 
1482  integer :: err
1483 
1484  if (fileobj%is_root) then
1485  err = nf90_inquire(fileobj%ncid, nvariables=nvars)
1486  call check_netcdf_code(err, "get_num_variables: file: "//trim(fileobj%path))
1487  endif
1488  if (present(broadcast)) then
1489  if (.not. broadcast) then
1490  return
1491  endif
1492  endif
1493  call mpp_broadcast(nvars, fileobj%io_root, pelist=fileobj%pelist)
1494 end function get_num_variables
1495 
1496 
1497 !> @brief Get the names of the variables in a file.
1498 subroutine get_variable_names(fileobj, names, broadcast)
1499 
1500  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
1501  character(len=*), dimension(:), intent(inout) :: names !< Names of the
1502  !! variables.
1503  logical, intent(in), optional :: broadcast !< Flag controlling whether or
1504  !! not the data will be
1505  !! broadcasted to non
1506  !! "I/O root" ranks.
1507  !! The broadcast will be done
1508  !! by default.
1509 
1510  integer :: nvars
1511  integer :: i
1512  integer :: err
1513 
1514  if (fileobj%is_root) then
1515  nvars = get_num_variables(fileobj, broadcast=.false.)
1516  if (nvars .gt. 0) then
1517  if (size(names) .ne. nvars) then
1518  call error("'names' has to be the same size of the number of variables. &
1519  &Check your get_variable_names call for file "//trim(fileobj%path))
1520  endif
1521  else
1522  call error("get_variable_names: the file "//trim(fileobj%path)//" does not have any variables")
1523  endif
1524  names(:) = ""
1525  do i = 1, nvars
1526  err = nf90_inquire_variable(fileobj%ncid, i, name=names(i))
1527  call check_netcdf_code(err, "get_variable_names: "//trim(fileobj%path))
1528  enddo
1529  endif
1530  if (present(broadcast)) then
1531  if (.not. broadcast) then
1532  return
1533  endif
1534  endif
1535  call mpp_broadcast(nvars, fileobj%io_root, pelist=fileobj%pelist)
1536  if (.not. fileobj%is_root) then
1537  if (nvars .gt. 0) then
1538  if (size(names) .ne. nvars) then
1539  call error("'names' has to be the same size of the number of variables. &
1540  &Check your get_variable_names call for file "//trim(fileobj%path))
1541  endif
1542  else
1543  call error("get_variable_names: the file "//trim(fileobj%path)//" does not have any variables")
1544  endif
1545  names(:) = ""
1546  endif
1547  call mpp_broadcast(names, len(names(nvars)), fileobj%io_root, &
1548  pelist=fileobj%pelist)
1549 end subroutine get_variable_names
1550 
1551 
1552 !> @brief Determine if a variable exists.
1553 !! @return Flag telling if the variable exists.
1554 function variable_exists(fileobj, variable_name, broadcast) &
1555  result(var_exists)
1556 
1557  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
1558  character(len=*), intent(in) :: variable_name !< Variable name.
1559  logical, intent(in), optional :: broadcast !< Flag controlling whether or
1560  !! not the data will be
1561  !! broadcasted to non
1562  !! "I/O root" ranks.
1563  !! The broadcast will be done
1564  !! by default.
1565  logical :: var_exists
1566 
1567  integer :: varid
1568 
1569  if (fileobj%is_root) then
1570  varid = get_variable_id(fileobj%ncid, trim(variable_name), &
1571  msg="variable_exists: file:"//trim(fileobj%path)//" variable:"//trim(variable_name), &
1572  allow_failure=.true.)
1573  var_exists = varid .ne. variable_missing
1574  endif
1575  if (present(broadcast)) then
1576  if (.not. broadcast) then
1577  return
1578  endif
1579  endif
1580  call mpp_broadcast(var_exists, fileobj%io_root, pelist=fileobj%pelist)
1581 end function variable_exists
1582 
1583 
1584 !> @brief Get the number of dimensions a variable depends on.
1585 !! @return Number of dimensions that the variable depends on.
1586 function get_variable_num_dimensions(fileobj, variable_name, broadcast) &
1587  result(ndims)
1588 
1589  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
1590  character(len=*), intent(in) :: variable_name !< Variable name.
1591  logical, intent(in), optional :: broadcast !< Flag controlling whether or
1592  !! not the data will be
1593  !! broadcasted to non
1594  !! "I/O root" ranks.
1595  !! The broadcast will be done
1596  !! by default.
1597  integer :: ndims
1598 
1599  integer :: varid
1600  integer :: err
1601  character(len=200) :: append_error_msg !< Msg to be appended to FATAL error message
1602 
1603 
1604  if (fileobj%is_root) then
1605  append_error_msg = "get_variable_num_dimension: file:"//trim(fileobj%path)//" variable: "//trim(variable_name)
1606  varid = get_variable_id(fileobj%ncid, trim(variable_name), msg=append_error_msg)
1607  err = nf90_inquire_variable(fileobj%ncid, varid, ndims=ndims)
1608  call check_netcdf_code(err, append_error_msg)
1609  endif
1610  if (present(broadcast)) then
1611  if (.not. broadcast) then
1612  return
1613  endif
1614  endif
1615  call mpp_broadcast(ndims, fileobj%io_root, pelist=fileobj%pelist)
1616 end function get_variable_num_dimensions
1617 
1618 
1619 !> @brief Get the name of a variable's dimensions.
1620 subroutine get_variable_dimension_names(fileobj, variable_name, dim_names, &
1621  broadcast)
1622 
1623  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
1624  character(len=*), intent(in) :: variable_name !< Variable name.
1625  character(len=*), dimension(:), intent(inout) :: dim_names !< Array of
1626  !! dimension
1627  !! names.
1628  logical, intent(in), optional :: broadcast !< Flag controlling whether or
1629  !! not the data will be
1630  !! broadcasted to non
1631  !! "I/O root" ranks.
1632  !! The broadcast will be done
1633  !! by default.
1634 
1635  integer :: varid
1636  integer :: err
1637  integer :: ndims
1638  integer,dimension(nf90_max_var_dims) :: dimids
1639  integer :: i
1640  character(len=200) :: append_error_msg !< Msg to be appended to FATAL error message
1641 
1642 
1643  if (fileobj%is_root) then
1644  append_error_msg = "get_variable_dimension_names: file:"//trim(fileobj%path)//" variable: "//trim(variable_name)
1645 
1646  varid = get_variable_id(fileobj%ncid, trim(variable_name), msg=append_error_msg)
1647  err = nf90_inquire_variable(fileobj%ncid, varid, ndims=ndims, &
1648  dimids=dimids)
1649  call check_netcdf_code(err, append_error_msg)
1650  if (ndims .gt. 0) then
1651  if (size(dim_names) .ne. ndims) then
1652  call error("'names' has to be the same size of the number of dimensions for the variable. &
1653  &Check your get_variable_dimension_names call for file "//trim(fileobj%path)// &
1654  " and variable:"//trim(variable_name))
1655  endif
1656  else
1657  call error("get_variable_dimension_names: the variable: "//trim(variable_name)//" in file: "//trim(fileobj%path)&
1658  & //" does not any dimensions. ")
1659  endif
1660  dim_names(:) = ""
1661  do i = 1, ndims
1662  err = nf90_inquire_dimension(fileobj%ncid, dimids(i), name=dim_names(i))
1663  call check_netcdf_code(err, append_error_msg)
1664  enddo
1665  endif
1666  if (present(broadcast)) then
1667  if (.not. broadcast) then
1668  return
1669  endif
1670  endif
1671  call mpp_broadcast(ndims, fileobj%io_root, pelist=fileobj%pelist)
1672  if (.not. fileobj%is_root) then
1673  if (ndims .gt. 0) then
1674  if (size(dim_names) .ne. ndims) then
1675  call error("'names' has to be the same size of the number of dimensions for the variable. &
1676  & Check your get_variable_dimension_names call for file "//trim(fileobj%path)// &
1677  " and variable:"//trim(variable_name))
1678  endif
1679  else
1680  call error("get_variable_dimension_names: the variable: "//trim(variable_name)//" in file: "//trim(fileobj%path)&
1681  & //" does not any dimensions. ")
1682  endif
1683  dim_names(:) = ""
1684  endif
1685  call mpp_broadcast(dim_names, len(dim_names(ndims)), fileobj%io_root, &
1686  pelist=fileobj%pelist)
1687 end subroutine get_variable_dimension_names
1688 
1689 
1690 !> @brief Get the size of a variable's dimensions.
1691 subroutine get_variable_size(fileobj, variable_name, dim_sizes, broadcast)
1692 
1693  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
1694  character(len=*), intent(in) :: variable_name !< Variable name.
1695  integer, dimension(:), intent(inout) :: dim_sizes !< Array of dimension
1696  !! sizes.
1697  logical, intent(in), optional :: broadcast !< Flag controlling whether or
1698  !! not the data will be
1699  !! broadcasted to non
1700  !! "I/O root" ranks.
1701  !! The broadcast will be done
1702  !! by default.
1703 
1704  integer :: varid
1705  integer :: err
1706  integer :: ndims
1707  integer,dimension(nf90_max_var_dims) :: dimids
1708  integer :: i
1709  character(len=200) :: append_error_msg !< Msg to be appended to FATAL error message
1710 
1711  if (fileobj%is_root) then
1712  append_error_msg = "get_variable_size: file:"//trim(fileobj%path)//" variable:"//trim(variable_name)
1713  varid = get_variable_id(fileobj%ncid, trim(variable_name), msg=append_error_msg)
1714  err = nf90_inquire_variable(fileobj%ncid, varid, ndims=ndims, dimids=dimids)
1715  call check_netcdf_code(err, append_error_msg)
1716  if (ndims .gt. 0) then
1717  if (size(dim_sizes) .ne. ndims) then
1718  call error("'dim_sizes' has to be the same size of the number of dimensions for the variable. &
1719  &Check your get_variable_size call for file "//trim(fileobj%path)// &
1720  " and variable:"//trim(variable_name))
1721  endif
1722  else
1723  call error("get_variable_size: the variable: "//trim(variable_name)//" in file: "//trim(fileobj%path)//&
1724  &" does not any dimensions. ")
1725  endif
1726  do i = 1, ndims
1727  err = nf90_inquire_dimension(fileobj%ncid, dimids(i), len=dim_sizes(i))
1728  call check_netcdf_code(err, append_error_msg)
1729  enddo
1730  endif
1731  if (present(broadcast)) then
1732  if (.not. broadcast) then
1733  return
1734  endif
1735  endif
1736  call mpp_broadcast(ndims, fileobj%io_root, pelist=fileobj%pelist)
1737  if (.not. fileobj%is_root) then
1738  if (ndims .gt. 0) then
1739  if (size(dim_sizes) .ne. ndims) then
1740  call error("'dim_sizes' has to be the same size of the number of dimensions for the variable. &
1741  &Check your get_variable_size call for file "//trim(fileobj%path)// &
1742  " and variable:"//trim(variable_name))
1743  endif
1744  else
1745  call error("get_variable_size: the variable: "//trim(variable_name)//" in file: "//trim(fileobj%path)//&
1746  &" does not any dimensions. ")
1747  endif
1748  endif
1749  call mpp_broadcast(dim_sizes, ndims, fileobj%io_root, pelist=fileobj%pelist)
1750 end subroutine get_variable_size
1751 
1752 
1753 !> @brief Get the index of a variable's unlimited dimensions.
1754 !! @return The index of the unlimited dimension, or else
1755 !! no_unlimited_dimension.
1756 function get_variable_unlimited_dimension_index(fileobj, variable_name, &
1757  broadcast) &
1758  result(unlim_dim_index)
1759 
1760  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
1761  character(len=*), intent(in) :: variable_name !< Variable name.
1762  logical, intent(in), optional :: broadcast !< Flag controlling whether or
1763  !! not the data will be
1764  !! broadcasted to non
1765  !! "I/O root" ranks.
1766  !! The broadcast will be done
1767  !! by default.
1768  integer :: unlim_dim_index
1769 
1770  integer :: ndims
1771  character(len=nf90_max_name), dimension(:), allocatable :: dim_names
1772  integer :: i
1773 
1774  unlim_dim_index = no_unlimited_dimension
1775  if (fileobj%is_root) then
1776  ndims = get_variable_num_dimensions(fileobj, variable_name, broadcast=.false.)
1777  allocate(dim_names(ndims))
1778  call get_variable_dimension_names(fileobj, variable_name, dim_names, &
1779  broadcast=.false.)
1780  do i = 1, size(dim_names)
1781  if (is_dimension_unlimited(fileobj, dim_names(i), .false.)) then
1782  unlim_dim_index = i
1783  exit
1784  endif
1785  enddo
1786  deallocate(dim_names)
1787  endif
1788  if (present(broadcast)) then
1789  if (.not. broadcast) then
1790  return
1791  endif
1792  endif
1793  call mpp_broadcast(unlim_dim_index, fileobj%io_root, pelist=fileobj%pelist)
1795 
1796 
1797 !> @brief Store the valid range for a variable.
1798 !! @return A ValidType_t object containing data about the valid
1799 !! range data for this variable can take.
1800 function get_valid(fileobj, variable_name) &
1801  result(valid)
1802 
1803  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
1804  character(len=*), intent(in) :: variable_name !< Variable name.
1805  type(valid_t) :: valid
1806 
1807  integer :: varid
1808  real(kind=r8_kind) :: scale_factor
1809  real(kind=r8_kind) :: add_offset
1810  real(kind=r8_kind), dimension(2) :: buffer
1811  integer :: xtype
1812  character(len=200) :: append_error_msg !< Msg to be appended to FATAL error message
1813 
1814  append_error_msg = "get_valid: file:"//trim(fileobj%path)
1815  if (fileobj%is_root) then
1816  varid = get_variable_id(fileobj%ncid, variable_name, msg=append_error_msg)
1817  valid%has_max = .false.
1818  valid%has_min = .false.
1819  valid%has_fill = .false.
1820  valid%has_missing = .false.
1821  valid%has_range = .false.
1822 
1823  !This routine makes use of netcdf's automatic type conversion to
1824  !store all range information in double precision.
1825  if (attribute_exists(fileobj%ncid, varid, "scale_factor", msg=append_error_msg)) then
1826  call get_variable_attribute(fileobj, variable_name, "scale_factor", scale_factor, &
1827  broadcast=.false.)
1828  else
1829  scale_factor = 1._r8_kind
1830  endif
1831  if (attribute_exists(fileobj%ncid, varid, "add_offset", msg=append_error_msg)) then
1832  call get_variable_attribute(fileobj, variable_name, "add_offset", add_offset, &
1833  broadcast=.false.)
1834  else
1835  add_offset = 0._r8_kind
1836  endif
1837 
1838  !valid%max_val and valid%min_val are defined by the "valid_range", "valid_min", and
1839  !"valid_max" variable attributes if they are present in the file. If either the maximum value
1840  !or minimum value is defined, valid%has_range is set to .true. (i.e. open ended ranges
1841  !are valid and should be tested within the is_valid function).
1842  if (attribute_exists(fileobj%ncid, varid, "valid_range", msg=append_error_msg)) then
1843  call get_variable_attribute(fileobj, variable_name, "valid_range", buffer, &
1844  broadcast=.false.)
1845  valid%max_val = buffer(2)*scale_factor + add_offset
1846  valid%has_max = .true.
1847  valid%min_val = buffer(1)*scale_factor + add_offset
1848  valid%has_min = .true.
1849  else
1850  if (attribute_exists(fileobj%ncid, varid, "valid_max", msg=append_error_msg)) then
1851  call get_variable_attribute(fileobj, variable_name, "valid_max", buffer(1), &
1852  broadcast=.false.)
1853  valid%max_val = buffer(1)*scale_factor + add_offset
1854  valid%has_max = .true.
1855  endif
1856  if (attribute_exists(fileobj%ncid, varid, "valid_min", msg=append_error_msg)) then
1857  call get_variable_attribute(fileobj, variable_name, "valid_min", buffer(1), &
1858  broadcast=.false.)
1859  valid%min_val = buffer(1)*scale_factor + add_offset
1860  valid%has_min = .true.
1861  endif
1862  endif
1863  valid%has_range = valid%has_min .or. valid%has_max
1864 
1865 
1866  !Get the missing value from the file if it exists.
1867  if (attribute_exists(fileobj%ncid, varid, "missing_value", msg=append_error_msg)) then
1868  call get_variable_attribute(fileobj, variable_name, "missing_value", buffer(1), &
1869  broadcast=.false.)
1870  valid%missing_val = buffer(1)*scale_factor + add_offset
1871  valid%has_missing = .true.
1872  endif
1873 
1874  !Get the fill value from the file if it exists.
1875  !If the _FillValue attribute is present and the maximum or minimum value is not defined,
1876  !then the maximum or minimum value will be determined by the _FillValue according to the NUG convention.
1877  !The NUG convention states that a positive fill value will be the exclusive upper
1878  !bound (i.e. valid values are less than the fill value), while a
1879  !non-positive fill value will be the exclusive lower bound (i.e. valis
1880  !values are greater than the fill value). As before, valid%has_range is true
1881  !if either a maximum or minimum value is set.
1882  if (attribute_exists(fileobj%ncid, varid, "_FillValue", msg=append_error_msg)) then
1883  call get_variable_attribute(fileobj, variable_name, "_FillValue", buffer(1), &
1884  broadcast=.false.)
1885  valid%fill_val = buffer(1)*scale_factor + add_offset
1886  valid%has_fill = .true.
1887  xtype = get_variable_type(fileobj%ncid, varid, msg=append_error_msg)
1888 
1889  if (.not. valid%has_range) then
1890  if (xtype .eq. nf90_short .or. xtype .eq. nf90_int) then
1891  if (buffer(1) .gt. 0) then
1892  valid%max_val = (buffer(1) - 1._r8_kind)*scale_factor + add_offset
1893  valid%has_max = .true.
1894  else
1895  valid%min_val = (buffer(1) + 1._r8_kind)*scale_factor + add_offset
1896  valid%has_min = .true.
1897  endif
1898  elseif (xtype .eq. nf90_float .or. xtype .eq. nf90_double) then
1899  if (buffer(1) .gt. 0) then
1900  valid%max_val = (nearest(nearest(buffer(1), -1._r8_kind), -1._r8_kind)) &
1901  *scale_factor + add_offset
1902  valid%has_max = .true.
1903  else
1904  valid%min_val = (nearest(nearest(buffer(1), 1._r8_kind), 1._r8_kind)) &
1905  *scale_factor + add_offset
1906  valid%has_min = .true.
1907  endif
1908  else
1909  call error("Unsupported variable type:"//trim(append_error_msg))
1910  endif
1911  valid%has_range = .true.
1912  endif
1913  endif
1914 
1915  endif
1916 
1917  call mpp_broadcast(valid%has_min, fileobj%io_root, pelist=fileobj%pelist)
1918  if (valid%has_min) then
1919  call mpp_broadcast(valid%min_val, fileobj%io_root, pelist=fileobj%pelist)
1920  endif
1921  call mpp_broadcast(valid%has_max, fileobj%io_root, pelist=fileobj%pelist)
1922  if (valid%has_max) then
1923  call mpp_broadcast(valid%max_val, fileobj%io_root, pelist=fileobj%pelist)
1924  endif
1925  call mpp_broadcast(valid%has_range, fileobj%io_root, pelist=fileobj%pelist)
1926 
1927  call mpp_broadcast(valid%has_fill, fileobj%io_root, pelist=fileobj%pelist)
1928  if (valid%has_fill) then
1929  call mpp_broadcast(valid%fill_val, fileobj%io_root, pelist=fileobj%pelist)
1930  endif
1931 
1932  call mpp_broadcast(valid%has_missing, fileobj%io_root, pelist=fileobj%pelist)
1933  if (valid%has_missing) then
1934  call mpp_broadcast(valid%missing_val, fileobj%io_root, pelist=fileobj%pelist)
1935  endif
1936 
1937 end function get_valid
1938 
1939 !> @brief Determine if a piece of (r4) data is "valid" (in the correct range.)
1940 !! @return A flag telling if the data element is "valid."
1941 elemental function is_valid_r8(datum, validobj) &
1942  result(valid_data)
1943 
1944  real(kind=r8_kind), intent(in) :: datum !< Unpacked data element.
1945  type(valid_t), intent(in) :: validobj !< Valid object.
1946  logical :: valid_data
1947 
1948  valid_data = check_if_valid(datum, validobj)
1949 
1950 end function is_valid_r8
1951 
1952 !> @brief Determine if a piece of (r8) data is "valid" (in the correct range.)
1953 !! @return A flag telling if the data element is "valid."
1954 elemental function is_valid_r4(datum, validobj) &
1955  result(valid_data)
1956 
1957  real(kind=r4_kind), intent(in) :: datum !< Unpacked data element.
1958  type(valid_t), intent(in) :: validobj !< Valid object.
1959  logical :: valid_data
1960 
1961  real(kind=r8_kind) :: rdatum
1962 
1963  !< Convert the data to r8 so it can be compared correctly since validobj%*_val are defined as r8
1964  rdatum = real(datum, kind=r8_kind)
1965  valid_data = check_if_valid(rdatum, validobj)
1966 
1967 end function is_valid_r4
1968 
1969 !> @brief Determine if a piece of data is "valid" (in the correct range.)
1970 !! @return A flag telling if the data element is "valid."
1971 elemental function check_if_valid(rdatum, validobj) &
1972  result(valid_data)
1973 
1974  real(kind=r8_kind), intent(in) :: rdatum !< packed data element.
1975  type(valid_t), intent(in) :: validobj !< Valid object.
1976  logical :: valid_data
1977 
1978  valid_data = .true.
1979  ! If the variable has a range (open or closed), valid values must be in that
1980  ! range.
1981  if (validobj%has_range) then
1982  if (validobj%has_min .and. .not. validobj%has_max) then
1983  valid_data = rdatum .ge. validobj%min_val
1984  elseif (validobj%has_max .and. .not. validobj%has_min) then
1985  valid_data = rdatum .le. validobj%max_val
1986  else
1987  valid_data = .not. (rdatum .lt. validobj%min_val .or. rdatum .gt. validobj%max_val)
1988  endif
1989  endif
1990  ! If the variable has a fill value or missing value, valid values must not be
1991  ! equal to either.
1992  if (validobj%has_fill .or. validobj%has_missing) then
1993  if (validobj%has_fill .and. .not. validobj%has_missing) then
1994  valid_data = rdatum .ne. validobj%fill_val
1995  elseif (validobj%has_missing .and. .not. validobj%has_fill) then
1996  valid_data = rdatum .ne. validobj%missing_val
1997  else
1998  valid_data = .not. (rdatum .eq. validobj%missing_val .or. rdatum .eq. validobj%fill_val)
1999  endif
2000  endif
2001 end function check_if_valid
2002 
2003 
2004 !> @brief Gathers a compressed arrays size and offset for each pe.
2005 subroutine compressed_start_and_count(fileobj, nelems, npes_start, npes_count)
2006 
2007  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
2008  integer, intent(in) :: nelems !< Number of elements on the current pe.
2009  integer, dimension(:), allocatable, intent(out) :: npes_start !< Offset for each pe.
2010  integer, dimension(:), allocatable, intent(out) :: npes_count !< Number of elements for
2011  !! each pe.
2012 
2013  integer :: i
2014 
2015  allocate(npes_start(size(fileobj%pelist)))
2016  allocate(npes_count(size(fileobj%pelist)))
2017  do i = 1, size(fileobj%pelist)
2018  if (fileobj%pelist(i) .eq. mpp_pe()) then
2019  npes_count(i) = nelems
2020  else
2021  call mpp_recv(npes_count(i), fileobj%pelist(i), block=.false.)
2022  call mpp_send(nelems, fileobj%pelist(i))
2023  endif
2024  enddo
2025  call mpp_sync_self(check=event_recv)
2026  call mpp_sync_self(check=event_send)
2027  npes_start(1) = 1
2028  do i = 1, size(fileobj%pelist)-1
2029  npes_start(i+1) = npes_start(i) + npes_count(i)
2030  enddo
2031 end subroutine compressed_start_and_count
2032 
2033 
2034 include "netcdf_add_restart_variable.inc"
2035 include "netcdf_read_data.inc"
2036 include "netcdf_write_data.inc"
2037 include "register_global_attribute.inc"
2038 include "register_variable_attribute.inc"
2039 include "get_global_attribute.inc"
2040 include "get_variable_attribute.inc"
2041 include "compressed_write.inc"
2042 include "compressed_read.inc"
2043 include "scatter_data_bc.inc"
2044 include "gather_data_bc.inc"
2045 include "unpack_data.inc"
2046 
2047 !> @brief Wrapper to distinguish interfaces.
2048 function netcdf_file_open_wrap(fileobj, path, mode, nc_format, pelist, is_restart, dont_add_res_to_filename) &
2049  result(success)
2050 
2051  type(fmsnetcdffile_t), intent(inout) :: fileobj !< File object.
2052  character(len=*), intent(in) :: path !< File path.
2053  character(len=*), intent(in) :: mode !< File mode. Allowed values are:
2054  !! "read", "append", "write", or
2055  !! "overwrite".
2056  character(len=*), intent(in), optional :: nc_format !< Netcdf format that
2057  !! new files are written
2058  !! as. Allowed values
2059  !! are: "64bit", "classic",
2060  !! or "netcdf4". Defaults to
2061  !! "64bit".
2062  integer, dimension(:), intent(in), optional :: pelist !< List of ranks associated
2063  !! with this file. If not
2064  !! provided, only the current
2065  !! rank will be able to
2066  !! act on the file.
2067  logical, intent(in), optional :: is_restart !< Flag telling if this file
2068  !! is a restart file. Defaults
2069  !! to false.
2070  logical, intent(in), optional :: dont_add_res_to_filename !< Flag indicating not to add
2071  !! ".res" to the filename
2072  logical :: success
2073 
2074  success = netcdf_file_open(fileobj, path, mode, nc_format, pelist, is_restart, dont_add_res_to_filename)
2075 end function netcdf_file_open_wrap
2076 
2077 
2078 !> @brief Wrapper to distinguish interfaces.
2079 subroutine netcdf_file_close_wrap(fileobj)
2080 
2081  type(fmsnetcdffile_t), intent(inout) :: fileobj !< File object.
2082 
2083  call netcdf_file_close(fileobj)
2084 end subroutine netcdf_file_close_wrap
2085 
2086 
2087 !> @brief Wrapper to distinguish interfaces.
2088 subroutine netcdf_add_variable_wrap(fileobj, variable_name, variable_type, dimensions, chunksizes)
2089 
2090  type(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
2091  character(len=*), intent(in) :: variable_name !< Variable name.
2092  character(len=*), intent(in) :: variable_type !< Variable type. Allowed
2093  !! values are: "int", "int64",
2094  !! "float", or "double".
2095  character(len=*), dimension(:), intent(in), optional :: dimensions !< Dimension names.
2096  integer, intent(in), optional :: chunksizes(:) !< netcdf chunksize to use for this variable (netcdf4 only)
2097 
2098  call netcdf_add_variable(fileobj, variable_name, variable_type, dimensions, chunksizes)
2099 end subroutine netcdf_add_variable_wrap
2100 
2101 !> @brief Wrapper to distinguish interfaces.
2102 subroutine netcdf_save_restart_wrap(fileobj, unlim_dim_level)
2103 
2104  type(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
2105  integer, intent(in), optional :: unlim_dim_level !< Unlimited dimension
2106  !! level.
2107 
2108  call netcdf_save_restart(fileobj, unlim_dim_level)
2109 end subroutine netcdf_save_restart_wrap
2110 
2111 
2112 !> @brief Returns a variable's fill value if it exists in the file.
2113 !! @return Flag telling if a fill value exists.
2114 function get_fill_value(fileobj, variable_name, fill_value, broadcast) &
2115  result(fill_exists)
2116 
2117  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
2118  character(len=*), intent(in) :: variable_name !< Variable name.
2119  class(*), intent(out) :: fill_value !< Fill value.
2120  logical, intent(in), optional :: broadcast !< Flag controlling whether or
2121  !! not the data will be
2122  !! broadcasted to non
2123  !! "I/O root" ranks.
2124  !! The broadcast will be done
2125  !! by default.
2126  logical :: fill_exists
2127 
2128  character(len=32), dimension(2) :: attribute_names
2129  logical :: bcast
2130  integer :: i
2131 
2132  fill_exists = .false.
2133  call string_copy(attribute_names(1), "_FillValue")
2134  call string_copy(attribute_names(2), "missing_value")
2135  if (present(broadcast)) then
2136  bcast = broadcast
2137  else
2138  bcast = .true.
2139  endif
2140  do i = 1, size(attribute_names)
2141  fill_exists = variable_att_exists(fileobj, variable_name, attribute_names(i), &
2142  broadcast=bcast)
2143  if (fill_exists) then
2144  call get_variable_attribute(fileobj, variable_name, attribute_names(i), &
2145  fill_value, broadcast=bcast)
2146  exit
2147  endif
2148  enddo
2149 end function get_fill_value
2150 
2151 
2152 function get_variable_sense(fileobj, variable_name) &
2153  result(variable_sense)
2154 
2155  class(fmsnetcdffile_t), intent(in) :: fileobj
2156  character(len=*), intent(in) :: variable_name
2157  integer :: variable_sense
2158 
2159  character(len=256) :: buf
2160 
2161  variable_sense = 0
2162  if (variable_att_exists(fileobj, variable_name, "positive")) then
2163  call get_variable_attribute(fileobj, variable_name, "positive", buf)
2164  if (string_compare(buf, "down")) then
2165  variable_sense = -1
2166  elseif (string_compare(buf, "up")) then
2167  variable_sense = 1
2168  endif
2169  endif
2170 end function get_variable_sense
2171 
2172 
2173 function get_variable_missing(fileobj, variable_name) &
2174  result(variable_missing)
2175 
2176  type(fmsnetcdffile_t), intent(in) :: fileobj
2177  character(len=*), intent(in) :: variable_name
2178  real(kind=r8_kind) :: variable_missing
2179 
2180  real(kind=r8_kind) :: variable_missing_1d(1) !< Workaround for pgi
2181 
2182  if (variable_att_exists(fileobj, variable_name, "_FillValue")) then
2183  call get_variable_attribute(fileobj, variable_name, "_FillValue", variable_missing_1d)
2184  elseif (variable_att_exists(fileobj, variable_name, "missing_value")) then
2185  call get_variable_attribute(fileobj, variable_name, "missing_value", variable_missing_1d)
2186  elseif (variable_att_exists(fileobj, variable_name, "missing")) then
2187  call get_variable_attribute(fileobj, variable_name, "missing", variable_missing_1d)
2188  else
2189  variable_missing_1d = mpp_fill_double
2190  endif
2191 
2192  variable_missing = variable_missing_1d(1)
2193 
2194 end function get_variable_missing
2195 
2196 
2197 subroutine get_variable_units(fileobj, variable_name, units)
2198 
2199  class(fmsnetcdffile_t), intent(in) :: fileobj
2200  character(len=*), intent(in) :: variable_name
2201  character(len=*), intent(out) :: units
2202 
2203  if (variable_att_exists(fileobj, variable_name, "units")) then
2204  call get_variable_attribute(fileobj, variable_name, "units", units)
2205  else
2206  units = "nounits"
2207  endif
2208 end subroutine get_variable_units
2209 
2210 
2211 subroutine get_time_calendar(fileobj, time_name, calendar_type)
2212 
2213  class(fmsnetcdffile_t), intent(in) :: fileobj
2214  character(len=*), intent(in) :: time_name
2215  character(len=*), intent(out) :: calendar_type
2216 
2217  if (variable_att_exists(fileobj, time_name, "calendar")) then
2218  call get_variable_attribute(fileobj, time_name, "calendar", calendar_type)
2219  elseif (variable_att_exists(fileobj, time_name, "calendar_type")) then
2220  call get_variable_attribute(fileobj, time_name, "calendar_type", calendar_type)
2221  else
2222  calendar_type = "unspecified"
2223  endif
2224 end subroutine get_time_calendar
2225 
2226 
2227 !> @brief Determine if a variable has been registered to a restart file..
2228 !! @return Flag telling if the variable has been registered to a restart file.
2229 function is_registered_to_restart(fileobj, variable_name) &
2230  result(is_registered)
2231 
2232  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
2233  character(len=*), intent(in) :: variable_name !< Variable name.
2234  logical :: is_registered
2235 
2236  integer :: i
2237 
2238  if (.not. fileobj%is_restart) then
2239  call error("file "//trim(fileobj%path)//" is not a restart file. &
2240  &Add is_restart=.true. to your open_file call")
2241  endif
2242  is_registered = .false.
2243  do i = 1, fileobj%num_restart_vars
2244  if (string_compare(fileobj%restart_vars(i)%varname, variable_name, .true.)) then
2245  is_registered = .true.
2246  exit
2247  endif
2248  enddo
2249 end function is_registered_to_restart
2250 
2251 
2252 function check_if_open(fileobj, fname) result(is_open)
2253  logical :: is_open !< True if the file in the file object is opened
2254  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
2255  character(len=*), intent(in), optional :: fname !< Optional filename for checking
2256 
2257  !Check if the is_open variable in the object has been allocated
2258  if (allocated(fileobj%is_open)) then
2259  is_open = fileobj%is_open !Return the value of the fileobj%is_open
2260  else
2261  is_open = .false. !If fileobj%is_open is not allocated, that the file has not been opened
2262  endif
2263 
2264  if (present(fname)) then
2265  !If the filename does not match the name in path,
2266  !then this is considered not open
2267  if (is_open .AND. trim(fname) .ne. trim(fileobj%path)) is_open = .false.
2268  endif
2269 end function check_if_open
2270 
2271 subroutine set_fileobj_time_name (fileobj,time_name)
2272  class(fmsnetcdffile_t), intent(inout) :: fileobj
2273  character(*),intent(in) :: time_name
2274  integer :: len_of_name
2275  len_of_name = len(trim(time_name))
2276  fileobj%time_name = ' '
2277  fileobj%time_name = time_name(1:len_of_name)
2278 ! if (.not. allocated(fileobj%time_name)) then
2279 ! allocate(character(len=len_of_name) :: fileobj%time_name)
2280 ! fileobj%time_name = time_name(1:len_of_name)
2281 ! else
2282 ! call error ("set_fileobj_time_name :: The time_name has already been set")
2283 ! endif
2284 end subroutine set_fileobj_time_name
2285 
2286 !> @brief Loop through the registered restart variables (including regional
2287 !! variables) and read them from the netcdf file
2288 subroutine read_restart_bc(fileobj, unlim_dim_level, ignore_checksum)
2289  class(fmsnetcdffile_t), intent(inout) :: fileobj !< File object
2290  integer, intent(in), optional :: unlim_dim_level !< Unlimited dimension
2291  !! level.
2292  logical, intent(in), optional :: ignore_checksum !< Checksum data integrity flag.
2293 
2294  integer :: i !< No description
2295 
2296  if (.not. fileobj%is_restart) then
2297  call error("file "//trim(fileobj%path)//" is not a restart file.")
2298  endif
2299 
2300  do i = 1, fileobj%num_restart_vars
2301  !> Go away if you are not in the pelist!
2302  if (.not.any(mpp_pe().eq.fileobj%restart_vars(i)%bc_info%pelist(:))) cycle
2303 
2304  !> The file's root pe reads the file and scatters it to the rest of the pes
2305  if (associated(fileobj%restart_vars(i)%data2d)) then
2306  call scatter_data_bc (fileobj, fileobj%restart_vars(i)%varname, &
2307  fileobj%restart_vars(i)%data2d, &
2308  fileobj%restart_vars(i)%bc_info, &
2309  unlim_dim_level = unlim_dim_level, &
2310  ignore_checksum=ignore_checksum)
2311  else if (associated(fileobj%restart_vars(i)%data3d)) then
2312  call scatter_data_bc (fileobj, fileobj%restart_vars(i)%varname, &
2313  fileobj%restart_vars(i)%data3d, &
2314  fileobj%restart_vars(i)%bc_info, &
2315  unlim_dim_level = unlim_dim_level, &
2316  ignore_checksum=ignore_checksum)
2317  endif
2318  end do
2319 
2320 
2321 end subroutine read_restart_bc
2322 
2323 !> @brief Loop through the registered restart variables (including regional
2324 !! variables) and write them to the netcdf file
2325 subroutine write_restart_bc(fileobj, unlim_dim_level)
2326  class(fmsnetcdffile_t), intent(inout) :: fileobj !< File object
2327  integer, intent(in), optional :: unlim_dim_level !< Unlimited dimension
2328  !! level.
2329  integer :: i !< No description
2330 
2331  if (.not. fileobj%is_restart) then
2332  call error("file "//trim(fileobj%path)//" is not a restart file. &
2333  &Add is_restart=.true. to your open_file call")
2334  endif
2335 
2336  !> Loop through the variables, root pe gathers the data from the other pes and writes out the checksum.
2337  !! Then loop through the variables again to write the data to the netcdf file.
2338  !! All the metadata should be written before the data to prevent netcdf from rewriting the file
2339  !! if the header is not big enough. That is why the two do loops are needed.
2340 
2341  do i = 1, fileobj%num_restart_vars
2342  !> Go away if you are not in the pelist!
2343  if (.not.any(mpp_pe().eq.fileobj%restart_vars(i)%bc_info%pelist(:))) cycle
2344 
2345  !> Go away if this is not a BC variable
2346  if (.not. fileobj%restart_vars(i)%is_bc_variable) cycle
2347 
2348  !> Root pe gathers the data from the other ranks, saves it in a buffer, and writes out the checksum.
2349  if (associated(fileobj%restart_vars(i)%data2d)) then
2350  call gather_data_bc(fileobj, fileobj%restart_vars(i)%data2d, fileobj%restart_vars(i)%bc_info)
2351  call register_variable_attribute(fileobj, fileobj%restart_vars(i)%varname, "checksum", &
2352  fileobj%restart_vars(i)%bc_info%chksum(1:len(fileobj%restart_vars(i)%bc_info%chksum)),&
2353  str_len=len(fileobj%restart_vars(i)%bc_info%chksum))
2354  else if (associated(fileobj%restart_vars(i)%data3d)) then
2355  call gather_data_bc(fileobj, fileobj%restart_vars(i)%data3d, fileobj%restart_vars(i)%bc_info)
2356  call register_variable_attribute(fileobj, fileobj%restart_vars(i)%varname, "checksum", &
2357  fileobj%restart_vars(i)%bc_info%chksum(1:len(fileobj%restart_vars(i)%bc_info%chksum)),&
2358  str_len=len(fileobj%restart_vars(i)%bc_info%chksum))
2359  endif
2360  enddo
2361 
2362  !> Write the data to the netcdf file
2363  do i = 1, fileobj%num_restart_vars
2364  if (allocated(fileobj%restart_vars(i)%bc_info%globaldata2d_r8 )) then
2365  call netcdf_write_data(fileobj, fileobj%restart_vars(i)%varname, &
2366  fileobj%restart_vars(i)%bc_info%globaldata2d_r8 , &
2367  unlim_dim_level=unlim_dim_level)
2368  deallocate(fileobj%restart_vars(i)%bc_info%globaldata2d_r8)
2369  else if (allocated(fileobj%restart_vars(i)%bc_info%globaldata2d_r4 )) then
2370  call netcdf_write_data(fileobj, fileobj%restart_vars(i)%varname, &
2371  fileobj%restart_vars(i)%bc_info%globaldata2d_r4 , &
2372  unlim_dim_level=unlim_dim_level)
2373  deallocate(fileobj%restart_vars(i)%bc_info%globaldata2d_r4)
2374  else if (allocated(fileobj%restart_vars(i)%bc_info%globaldata3d_r8 )) then
2375  call netcdf_write_data(fileobj, fileobj%restart_vars(i)%varname, &
2376  fileobj%restart_vars(i)%bc_info%globaldata3d_r8 , &
2377  unlim_dim_level=unlim_dim_level)
2378  deallocate(fileobj%restart_vars(i)%bc_info%globaldata3d_r8)
2379  else if (allocated(fileobj%restart_vars(i)%bc_info%globaldata3d_r4 )) then
2380  call netcdf_write_data(fileobj, fileobj%restart_vars(i)%varname, &
2381  fileobj%restart_vars(i)%bc_info%globaldata3d_r4 , &
2382  unlim_dim_level=unlim_dim_level)
2383  deallocate(fileobj%restart_vars(i)%bc_info%globaldata3d_r4 )
2384  endif
2385 
2386  enddo
2387 
2388 end subroutine write_restart_bc
2389 
2390 !> @brief flushes the netcdf file into disk
2391 subroutine flush_file(fileobj)
2392  class(fmsnetcdffile_t), intent(inout) :: fileobj !< FMS2_io fileobj
2393 
2394  integer :: err !< Netcdf error code
2395 
2396  if (fileobj%is_root) then
2397  err = nf90_sync(fileobj%ncid)
2398  call check_netcdf_code(err, "Flush_file: File:"//trim(fileobj%path))
2399  endif
2400 end subroutine flush_file
2401 
2402 end module netcdf_io_mod
2403 !> @}
2404 ! close documentation grouping
subroutine compressed_write_1d_wrap(fileobj, variable_name, cdata, unlim_dim_level, corner, edge_lengths)
Wrapper to distinguish interfaces.
subroutine register_variable_attribute_1d(fileobj, variable_name, attribute_name, attribute_value, str_len)
Add an attribute to a variable.
subroutine compressed_write_3d(fileobj, variable_name, cdata, unlim_dim_level, corner, edge_lengths)
For variables without a compressed dimension, this routine simply wraps netcdf_write data....
subroutine compressed_write_5d_wrap(fileobj, variable_name, cdata, unlim_dim_level, corner, edge_lengths)
Wrapper to distinguish interfaces.
subroutine compressed_write_0d(fileobj, variable_name, cdata, unlim_dim_level, corner)
For variables without a compressed dimension, this routine simply wraps netcdf_write data....
subroutine compressed_read_0d(fileobj, variable_name, cdata, unlim_dim_level, corner)
I/O domain reads in data from the netcdf file and broadcasts the data to the rest of the ranks....
subroutine compressed_write_5d(fileobj, variable_name, cdata, unlim_dim_level, corner, edge_lengths)
For variables without a compressed dimension, this routine simply wraps netcdf_write data....
subroutine compressed_read_5d(fileobj, variable_name, cdata, unlim_dim_level, corner, edge_lengths)
I/O domain reads in data from the netcdf file and broadcasts the data to the rest of the ranks....
subroutine compressed_read_3d(fileobj, variable_name, cdata, unlim_dim_level, corner, edge_lengths)
I/O domain reads in data from the netcdf file and broadcasts the data to the rest of the ranks....
subroutine compressed_write_0d_wrap(fileobj, variable_name, cdata, unlim_dim_level, corner)
Wrapper to distinguish interfaces.
subroutine compressed_read_1d(fileobj, variable_name, cdata, unlim_dim_level, corner, edge_lengths)
I/O domain reads in data from the netcdf file and broadcasts the data to the rest of the ranks....
subroutine compressed_read_2d(fileobj, variable_name, cdata, unlim_dim_level, corner, edge_lengths)
I/O domain reads in data from the netcdf file and broadcasts the data to the rest of the ranks....
subroutine compressed_write_4d(fileobj, variable_name, cdata, unlim_dim_level, corner, edge_lengths)
For variables without a compressed dimension, this routine simply wraps netcdf_write data....
subroutine compressed_write_2d_wrap(fileobj, variable_name, cdata, unlim_dim_level, corner, edge_lengths)
Wrapper to distinguish interfaces.
subroutine compressed_write_1d(fileobj, variable_name, cdata, unlim_dim_level, corner, edge_lengths)
For variables without a compressed dimension, this routine simply wraps netcdf_write data....
subroutine compressed_write_4d_wrap(fileobj, variable_name, cdata, unlim_dim_level, corner, edge_lengths)
Wrapper to distinguish interfaces.
subroutine compressed_read_4d(fileobj, variable_name, cdata, unlim_dim_level, corner, edge_lengths)
I/O domain reads in data from the netcdf file and broadcasts the data to the rest of the ranks....
subroutine compressed_write_2d(fileobj, variable_name, cdata, unlim_dim_level, corner, edge_lengths)
For variables without a compressed dimension, this routine simply wraps netcdf_write data....
subroutine register_variable_attribute_0d(fileobj, variable_name, attribute_name, attribute_value, str_len)
Add an attribute to a variable.
subroutine compressed_write_3d_wrap(fileobj, variable_name, cdata, unlim_dim_level, corner, edge_lengths)
Wrapper to distinguish interfaces.
subroutine, public get_instance_filename(name_in, name_out)
Adds the filename_appendix to name_in and sets it as name_out.
subroutine, public restart_filepath_mangle(dest, source)
Add ".res" to an input file path.
logical function, public file_exists(path)
Determine if a file exists.
logical function, public string_compare(string1, string2, ignore_case)
Compare strings.
subroutine mpp_sync_self(pelist, check, request, msg_size, msg_type)
This is to check if current PE's outstanding puts are complete but we can't use shmem_fence because w...
integer, parameter, public mpp_comm_null
MPP_COMM_NULL acts as an analagous mpp-macro for MPI_COMM_NULL to share with fms2_io NetCDF4 mpi-io....
Definition: mpp.F90:1359
integer, parameter, public mpp_info_null
MPP_INFO_NULL acts as an analagous mpp-macro for MPI_INFO_NULL to share with fms2_io NetCDF4 mpi-io....
Definition: mpp.F90:1350
integer function mpp_pe()
Returns processor ID.
Definition: mpp_util.inc:407
Perform parallel broadcasts.
Definition: mpp.F90:1105
Error handler.
Definition: mpp.F90:382
Gather data sent from pelist onto the root pe Wrapper for MPI_gather, can be used with and without in...
Definition: mpp.F90:698
Receive data from another PE.
Definition: mpp.F90:951
Send data to a receiving PE.
Definition: mpp.F90:1018
integer, private fms2_ncchksz
Chunksize (bytes) used in nc_open and nc_create.
Definition: netcdf_io.F90:54
subroutine, public netcdf_restore_state(fileobj, unlim_dim_level)
Loop through registered restart variables and read them from a netcdf file.
Definition: netcdf_io.F90:1135
subroutine append_compressed_dimension(fileobj, dim_name, npes_corner, npes_nelems)
Add a compressed dimension to a file object.
Definition: netcdf_io.F90:800
subroutine, public netcdf_file_close(fileobj)
Close a netcdf file.
Definition: netcdf_io.F90:737
subroutine netcdf_add_restart_variable_4d(fileobj, variable_name, vdata, dimensions, is_optional, chunksizes)
Add a restart variable to a netcdf file.
integer function, dimension(2) get_variable_compressed_dimension_index(fileobj, variable_name, broadcast)
Given a compressed variable, get the index of the compressed dimension.
Definition: netcdf_io.F90:1020
logical function, public netcdf_file_open(fileobj, path, mode, nc_format, pelist, is_restart, dont_add_res_to_filename)
Open a netcdf file.
Definition: netcdf_io.F90:542
elemental logical function is_valid_r8(datum, validobj)
Determine if a piece of (r4) data is "valid" (in the correct range.)
Definition: netcdf_io.F90:1943
logical function attribute_exists(ncid, varid, attribute_name, msg)
Determine if an attribute exists.
Definition: netcdf_io.F90:479
subroutine, public get_variable_size(fileobj, variable_name, dim_sizes, broadcast)
Get the size of a variable's dimensions.
Definition: netcdf_io.F90:1692
subroutine, public netcdf_add_dimension(fileobj, dimension_name, dimension_length, is_compressed)
Add a dimension to a file.
Definition: netcdf_io.F90:870
subroutine netcdf_add_restart_variable_0d(fileobj, variable_name, vdata, dimensions, is_optional, chunksizes)
Add a restart variable to a netcdf file.
type(valid_t) function, public get_valid(fileobj, variable_name)
Store the valid range for a variable.
Definition: netcdf_io.F90:1802
subroutine, public read_restart_bc(fileobj, unlim_dim_level, ignore_checksum)
Loop through the registered restart variables (including regional variables) and read them from the n...
Definition: netcdf_io.F90:2289
logical function, public is_dimension_unlimited(fileobj, dimension_name, broadcast)
Determine where or not the dimension is unlimited.
Definition: netcdf_io.F90:1368
character(len=10), private fms2_nc_format
Netcdf format type used in netcdf_file_open.
Definition: netcdf_io.F90:56
integer function, public get_variable_unlimited_dimension_index(fileobj, variable_name, broadcast)
Get the index of a variable's unlimited dimensions.
Definition: netcdf_io.F90:1759
logical function, public netcdf_file_open_wrap(fileobj, path, mode, nc_format, pelist, is_restart, dont_add_res_to_filename)
Wrapper to distinguish interfaces.
Definition: netcdf_io.F90:2050
integer function get_compressed_dimension_index(fileobj, dim_name)
Get the index of a compressed dimension in a file object.
Definition: netcdf_io.F90:779
subroutine, public compressed_start_and_count(fileobj, nelems, npes_start, npes_count)
Gathers a compressed arrays size and offset for each pe.
Definition: netcdf_io.F90:2006
subroutine scatter_data_bc_2d(fileobj, varname, vdata, bc_info, unlim_dim_level, ignore_checksum)
integer, private fms2_header_buffer_val
value used in NF__ENDDEF
Definition: netcdf_io.F90:57
integer, private fms2_nc_format_param
Netcdf format type param used in nc_create.
Definition: netcdf_io.F90:55
elemental logical function is_valid_r4(datum, validobj)
Determine if a piece of (r8) data is "valid" (in the correct range.)
Definition: netcdf_io.F90:1956
integer function, public get_variable_num_dimensions(fileobj, variable_name, broadcast)
Get the number of dimensions a variable depends on.
Definition: netcdf_io.F90:1588
subroutine, public get_dimension_size(fileobj, dimension_name, dim_size, broadcast)
Get the length of a dimension.
Definition: netcdf_io.F90:1438
subroutine netcdf_add_restart_variable_3d(fileobj, variable_name, vdata, dimensions, is_optional, chunksizes)
Add a restart variable to a netcdf file.
subroutine netcdf_write_data_4d(fileobj, variable_name, variable_data, unlim_dim_level, corner, edge_lengths)
Write data to a variable in a netcdf file.
subroutine netcdf_read_data_1d(fileobj, variable_name, buf, unlim_dim_level, corner, edge_lengths, broadcast)
Read in data from a variable in a netcdf file.
subroutine netcdf_add_restart_variable_5d_wrap(fileobj, variable_name, vdata, dimensions, is_optional, chunksizes)
Wrapper to distinguish interfaces.
elemental logical function check_if_valid(rdatum, validobj)
Determine if a piece of data is "valid" (in the correct range.)
Definition: netcdf_io.F90:1973
subroutine netcdf_add_restart_variable_1d_wrap(fileobj, variable_name, vdata, dimensions, is_optional, chunksizes)
Wrapper to distinguish interfaces.
logical function, public global_att_exists(fileobj, attribute_name, broadcast)
Determine if a global attribute exists.
Definition: netcdf_io.F90:1183
subroutine, public set_netcdf_mode(ncid, mode)
Switch to the correct netcdf mode.
Definition: netcdf_io.F90:394
subroutine netcdf_write_data_0d(fileobj, variable_name, variable_data, unlim_dim_level, corner)
Write data to a variable in a netcdf file.
subroutine register_global_attribute_0d(fileobj, attribute_name, attribute_value, str_len)
Add a global attribute.
subroutine get_variable_attribute_1d(fileobj, variable_name, attribute_name, attribute_value, broadcast)
Get the value of a variable's attribute.
subroutine netcdf_add_restart_variable_0d_wrap(fileobj, variable_name, vdata, dimensions, is_optional, chunksizes)
Wrapper to distinguish interfaces.
logical, private fms2_shuffle
Flag indicating whether to use the netcdf shuffle filter.
Definition: netcdf_io.F90:60
subroutine get_global_attribute_1d(fileobj, attribute_name, attribute_value, broadcast)
Get the value of a global attribute.
subroutine netcdf_add_restart_variable_3d_wrap(fileobj, variable_name, vdata, dimensions, is_optional, chunksizes)
Wrapper to distinguish interfaces.
subroutine netcdf_add_restart_variable_1d(fileobj, variable_name, vdata, dimensions, is_optional, chunksizes)
Add a restart variable to a netcdf file.
subroutine get_variable_attribute_0d(fileobj, variable_name, attribute_name, attribute_value, broadcast, reproduce_null_char_bug_flag)
Get the value of a variable's attribute.
integer function get_variable_id(ncid, variable_name, msg, allow_failure)
Get the id of a variable from its name.
Definition: netcdf_io.F90:451
subroutine register_restart_region_3d(fileobj, variable_name, vdata, indices, global_size, pelist, is_root_pe, x_halo, y_halo, jshift, ishift, is_optional)
Registers a regional 3D variable and stores the information needed.
subroutine register_global_attribute_1d(fileobj, attribute_name, attribute_value, str_len)
Add a global attribute.
logical function, public dimension_exists(fileobj, dimension_name, broadcast)
Determine if a dimension exists.
Definition: netcdf_io.F90:1332
logical function, public is_registered_to_restart(fileobj, variable_name)
Determine if a variable has been registered to a restart file..
Definition: netcdf_io.F90:2231
integer function get_dimension_id(ncid, dimension_name, msg, allow_failure)
Get the id of a dimension from its name.
Definition: netcdf_io.F90:423
subroutine, public register_compressed_dimension(fileobj, dimension_name, npes_corner, npes_nelems)
Add a compressed dimension.
Definition: netcdf_io.F90:923
logical function, public variable_att_exists(fileobj, variable_name, attribute_name, broadcast)
Determine if a variable's attribute exists.
Definition: netcdf_io.F90:1212
integer, private fms2_deflate_level
Netcdf deflate level to use in nf90_def_var (integer between 1 to 9)
Definition: netcdf_io.F90:58
subroutine netcdf_read_data_5d(fileobj, variable_name, buf, unlim_dim_level, corner, edge_lengths, broadcast)
Read in data from a variable in a netcdf file.
subroutine netcdf_add_restart_variable_2d(fileobj, variable_name, vdata, dimensions, is_optional, chunksizes)
Add a restart variable to a netcdf file.
subroutine netcdf_write_data_3d(fileobj, variable_name, variable_data, unlim_dim_level, corner, edge_lengths)
Write data to a variable in a netcdf file.
logical function, public get_fill_value(fileobj, variable_name, fill_value, broadcast)
Returns a variable's fill value if it exists in the file.
Definition: netcdf_io.F90:2116
subroutine gather_data_bc_3d(fileobj, vdata, bc_info)
gathers the 2d vdata from all of the relevant pes into the root_pe and saves it to a buffer.
subroutine, public get_variable_dimension_names(fileobj, variable_name, dim_names, broadcast)
Get the name of a variable's dimensions.
Definition: netcdf_io.F90:1622
subroutine, public check_netcdf_code(err, msg)
Check for errors returned by netcdf.
Definition: netcdf_io.F90:378
subroutine netcdf_write_data_1d(fileobj, variable_name, variable_data, unlim_dim_level, corner, edge_lengths)
Write data to a variable in a netcdf file.
subroutine, public flush_file(fileobj)
flushes the netcdf file into disk
Definition: netcdf_io.F90:2392
subroutine netcdf_add_restart_variable_4d_wrap(fileobj, variable_name, vdata, dimensions, is_optional, chunksizes)
Wrapper to distinguish interfaces.
subroutine get_global_attribute_0d(fileobj, attribute_name, attribute_value, broadcast)
Get the value of a global attribute.
subroutine netcdf_add_restart_variable_5d(fileobj, variable_name, vdata, dimensions, is_optional, chunksizes)
Add a restart variable to a netcdf file.
subroutine, public register_unlimited_compressed_axis(fileobj, dimension_name, dimension_length)
Add a "compressed" unlimited dimension to a netcdf file.
Definition: netcdf_io.F90:836
subroutine, public netcdf_save_restart(fileobj, unlim_dim_level)
Loop through registered restart variables and write them to a netcdf file.
Definition: netcdf_io.F90:1093
subroutine netcdf_read_data_3d(fileobj, variable_name, buf, unlim_dim_level, corner, edge_lengths, broadcast)
Read in data from a variable in a netcdf file.
subroutine, public netcdf_add_variable_wrap(fileobj, variable_name, variable_type, dimensions, chunksizes)
Wrapper to distinguish interfaces.
Definition: netcdf_io.F90:2089
subroutine netcdf_write_data_2d(fileobj, variable_name, variable_data, unlim_dim_level, corner, edge_lengths)
Write data to a variable in a netcdf file.
subroutine, public write_restart_bc(fileobj, unlim_dim_level)
Loop through the registered restart variables (including regional variables) and write them to the ne...
Definition: netcdf_io.F90:2326
integer function get_variable_type(ncid, varid, msg)
Get the type of a netcdf variable.
Definition: netcdf_io.F90:524
logical, private fms2_is_netcdf4
Flag indicating whether the default netcdf file format is netcdf4.
Definition: netcdf_io.F90:61
subroutine, public netcdf_save_restart_wrap(fileobj, unlim_dim_level)
Wrapper to distinguish interfaces.
Definition: netcdf_io.F90:2103
subroutine, public get_variable_names(fileobj, names, broadcast)
Get the names of the variables in a file.
Definition: netcdf_io.F90:1499
subroutine netcdf_read_data_2d(fileobj, variable_name, buf, unlim_dim_level, corner, edge_lengths, broadcast)
Read in data from a variable in a netcdf file.
subroutine netcdf_add_restart_variable_2d_wrap(fileobj, variable_name, vdata, dimensions, is_optional, chunksizes)
Wrapper to distinguish interfaces.
subroutine, public get_dimension_names(fileobj, names, broadcast)
Get the names of the dimensions in a file.
Definition: netcdf_io.F90:1275
subroutine gather_data_bc_2d(fileobj, vdata, bc_info)
gathers the 2d vdata from all of the relevant pes into the root_pe and saves it to a buffer.
logical function, public check_if_open(fileobj, fname)
Definition: netcdf_io.F90:2253
subroutine, public netcdf_io_init(chksz, header_buffer_val, netcdf_default_format, deflate_level, shuffle)
Accepts the namelist fms2_io_nml variables relevant to netcdf_io_mod.
Definition: netcdf_io.F90:346
integer function, public get_num_variables(fileobj, broadcast)
Determine the number of variables in a file.
Definition: netcdf_io.F90:1472
subroutine, public get_unlimited_dimension_name(fileobj, dimension_name, broadcast)
Get the name of the unlimited dimension.
Definition: netcdf_io.F90:1403
integer function get_attribute_type(ncid, varid, attname, msg)
Get the type of a netcdf attribute.
Definition: netcdf_io.F90:504
subroutine, public netcdf_add_variable(fileobj, variable_name, variable_type, dimensions, chunksizes)
Add a variable to a file.
Definition: netcdf_io.F90:951
subroutine scatter_data_bc_3d(fileobj, varname, vdata, bc_info, unlim_dim_level, ignore_checksum)
subroutine netcdf_write_data_5d(fileobj, variable_name, variable_data, unlim_dim_level, corner, edge_lengths)
Write data to a variable in a netcdf file.
subroutine register_restart_region_2d(fileobj, variable_name, vdata, indices, global_size, pelist, is_root_pe, x_halo, y_halo, jshift, ishift, is_optional)
Registers a regional 2D variable and stores the information needed.
subroutine, public netcdf_file_close_wrap(fileobj)
Wrapper to distinguish interfaces.
Definition: netcdf_io.F90:2080
subroutine netcdf_read_data_0d(fileobj, variable_name, buf, unlim_dim_level, corner, broadcast)
Read in data from a variable in a netcdf file.
logical function, public variable_exists(fileobj, variable_name, broadcast)
Determine if a variable exists.
Definition: netcdf_io.F90:1556
subroutine add_restart_var_to_array(fileobj, variable_name)
Add a restart variable to a FmsNetcdfFile_t type.
Definition: netcdf_io.F90:1066
integer function, public get_num_dimensions(fileobj, broadcast)
Determine the number of dimensions in a file.
Definition: netcdf_io.F90:1248
subroutine netcdf_read_data_4d(fileobj, variable_name, buf, unlim_dim_level, corner, edge_lengths, broadcast)
Read in data from a variable in a netcdf file.
The interface is needed to accomodate pgi because it can't handle class * and there was no other way ...
Definition: netcdf_io.F90:334
information needed fr regional restart variables
Definition: netcdf_io.F90:67
information about the current dimensions for regional restart variables
Definition: netcdf_io.F90:116
Range type for a netcdf variable.
Definition: netcdf_io.F90:163