25 integer,
intent(in) :: status
27 if (status /= nf90_noerr)
then 28 call die_msg(291021908, nf90_strerror(status))
40 integer,
intent(in) :: status
42 character(len=*),
intent(in) :: error_msg
44 if (status /= nf90_noerr)
then 45 call die_msg(701841139, trim(error_msg) &
46 //
" : " // trim(nf90_strerror(status)))
57 character(len=*),
intent(in) :: filename
59 integer,
intent(out) :: ncid
62 "opening " // trim(filename) //
" for reading")
72 character(len=*),
intent(in) :: filename
74 integer,
intent(out) :: ncid
77 "opening " // trim(filename) //
" for writing")
88 integer,
intent(in) :: ncid
100 integer,
intent(in) :: ncid
102 character(len=PMC_UUID_LEN),
intent(in) :: uuid
104 character(len=*),
intent(in) :: source
106 integer,
intent(in),
optional :: write_rank
108 integer,
intent(in),
optional :: write_n_proc
110 character(len=(len_trim(source) + 100)) :: history
111 integer :: use_rank, use_n_proc
114 call pmc_nc_check(nf90_put_att(ncid, nf90_global,
"title", &
115 trim(source) //
" output file"))
116 call pmc_nc_check(nf90_put_att(ncid, nf90_global,
"source", source))
117 call pmc_nc_check(nf90_put_att(ncid, nf90_global,
"UUID", uuid))
119 history((len_trim(history)+1):) = (
" created by " // trim(source))
120 call pmc_nc_check(nf90_put_att(ncid, nf90_global,
"history", history))
121 call pmc_nc_check(nf90_put_att(ncid, nf90_global,
"Conventions",
"CF-1.4"))
124 if (
present(write_rank))
then 125 use_rank = write_rank
129 if (
present(write_n_proc))
then 130 use_n_proc = write_n_proc
137 description=
"the process number (starting from 1) " &
138 //
"that output this data file")
140 description=
"total number of processes")
151 integer,
intent(in) :: ncid
153 real(kind=dp),
intent(out) :: var
155 character(len=*),
intent(in) :: name
158 logical,
optional,
intent(in) :: must_be_present
160 integer :: varid, status
161 logical :: use_must_be_present
163 if (
present(must_be_present))
then 164 use_must_be_present = must_be_present
166 use_must_be_present = .true.
168 status = nf90_inq_varid(ncid, name, varid)
169 if ((.not. use_must_be_present) .and. (status == nf90_enotvar))
then 176 "getting variable " // trim(name))
186 integer,
intent(in) :: ncid
188 integer,
intent(out) :: var
190 character(len=*),
intent(in) :: name
193 logical,
optional,
intent(in) :: must_be_present
195 integer :: varid, status
196 logical :: use_must_be_present
198 if (
present(must_be_present))
then 199 use_must_be_present = must_be_present
201 use_must_be_present = .true.
203 status = nf90_inq_varid(ncid, name, varid)
204 if ((.not. use_must_be_present) .and. (status == nf90_enotvar))
then 211 "getting variable " // trim(name))
221 integer,
intent(in) :: ncid
223 real(kind=dp),
intent(inout),
allocatable :: var(:)
225 character(len=*),
intent(in) :: name
228 logical,
optional,
intent(in) :: must_be_present
230 integer :: varid, status, dimids(1), size1
231 logical :: use_must_be_present
233 if (
present(must_be_present))
then 234 use_must_be_present = must_be_present
236 use_must_be_present = .true.
238 status = nf90_inq_varid(ncid, name, varid)
239 if ((.not. use_must_be_present) .and. (status == nf90_enotvar))
then 241 var = [
real(kind=dp)::]
246 "determining size of variable " // trim(name))
248 "determining size of dimension number " &
252 "getting variable " // trim(name))
262 integer,
intent(in) :: ncid
264 integer,
intent(inout),
allocatable :: var(:)
266 character(len=*),
intent(in) :: name
269 logical,
optional,
intent(in) :: must_be_present
271 integer :: varid, status, dimids(1), size1
272 logical :: use_must_be_present
274 if (
present(must_be_present))
then 275 use_must_be_present = must_be_present
277 use_must_be_present = .true.
279 status = nf90_inq_varid(ncid, name, varid)
280 if ((.not. use_must_be_present) .and. (status == nf90_enotvar))
then 287 "determining size of variable " // trim(name))
289 "determining size of dimension number " &
293 "getting variable " // trim(name))
303 integer,
intent(in) :: ncid
305 real(kind=dp),
intent(inout),
allocatable :: var(:,:)
307 character(len=*),
intent(in) :: name
310 logical,
optional,
intent(in) :: must_be_present
312 integer :: varid, status, dimids(2), size1, size2
313 logical :: use_must_be_present
315 if (
present(must_be_present))
then 316 use_must_be_present = must_be_present
318 use_must_be_present = .true.
320 status = nf90_inq_varid(ncid, name, varid)
321 if ((.not. use_must_be_present) .and. (status == nf90_enotvar))
then 323 var = reshape([
real(kind=dp)::], [0, 0])
328 "determining size of variable " // trim(name))
330 "determining size of dimension number " &
333 "determining size of dimension number " &
337 "getting variable " // trim(name))
347 integer,
intent(in) :: ncid
349 integer,
intent(inout),
allocatable :: var(:,:)
351 character(len=*),
intent(in) :: name
354 logical,
optional,
intent(in) :: must_be_present
356 integer :: varid, status, dimids(2), size1, size2
357 logical :: use_must_be_present
359 if (
present(must_be_present))
then 360 use_must_be_present = must_be_present
362 use_must_be_present = .true.
364 status = nf90_inq_varid(ncid, name, varid)
365 if ((.not. use_must_be_present) .and. (status == nf90_enotvar))
then 367 var = reshape([
integer::], [0, 0])
372 "determining size of variable " // trim(name))
374 "determining size of dimension number " &
377 "determining size of dimension number " &
381 "getting variable " // trim(name))
392 integer,
intent(in) :: ncid
394 integer,
intent(in) :: varid
396 character(len=*),
optional,
intent(in) :: unit
398 character(len=*),
optional,
intent(in) :: long_name
400 character(len=*),
optional,
intent(in) :: standard_name
402 character(len=*),
optional,
intent(in) :: description
404 if (
present(unit))
then 405 call pmc_nc_check(nf90_put_att(ncid, varid,
"unit", unit))
407 if (
present(long_name))
then 408 call pmc_nc_check(nf90_put_att(ncid, varid,
"long_name", long_name))
410 if (
present(standard_name))
then 411 call pmc_nc_check(nf90_put_att(ncid, varid,
"standard_name", &
414 if (
present(description))
then 415 call pmc_nc_check(nf90_put_att(ncid, varid,
"description", &
425 standard_name, description)
428 integer,
intent(in) :: ncid
430 real(kind=dp),
intent(in) :: var
432 character(len=*),
intent(in) :: name
434 character(len=*),
optional,
intent(in) :: unit
436 character(len=*),
optional,
intent(in) :: long_name
438 character(len=*),
optional,
intent(in) :: standard_name
440 character(len=*),
optional,
intent(in) :: description
442 integer :: varid, dimids(0)
445 call pmc_nc_check(nf90_def_var(ncid, name, nf90_double, dimids, varid))
458 standard_name, description)
461 integer,
intent(in) :: ncid
463 integer,
intent(in) :: var
465 character(len=*),
intent(in) :: name
467 character(len=*),
optional,
intent(in) :: unit
469 character(len=*),
optional,
intent(in) :: long_name
471 character(len=*),
optional,
intent(in) :: standard_name
473 character(len=*),
optional,
intent(in) :: description
475 integer :: varid, dimids(0)
478 call pmc_nc_check(nf90_def_var(ncid, name, nf90_int, dimids, varid))
494 integer,
intent(in) :: ncid
496 character(len=*),
intent(in) :: dim_name
498 integer,
intent(out) :: dimid
500 integer,
intent(in) :: dim_size
502 integer,
intent(in) :: array_dim
504 integer :: status, check_dim_size
505 character(len=NF90_MAX_NAME) :: check_name
507 status = nf90_inq_dimid(ncid, dim_name, dimid)
508 if (status == nf90_noerr)
then 509 call pmc_nc_check(nf90_inquire_dimension(ncid, dimid, check_name, &
511 call assert_msg(657263912, check_dim_size == dim_size, &
514 //
" of data array does not match size " &
516 // trim(dim_name) //
"' dim")
519 if (status /= nf90_ebaddim)
then 525 call pmc_nc_check(nf90_def_dim(ncid, dim_name, dim_size, dimid))
535 long_name, standard_name, description)
538 integer,
intent(in) :: ncid
540 real(kind=dp),
intent(in) :: var(:)
542 character(len=*),
intent(in) :: name
545 integer,
optional,
intent(in) :: dimids(1)
548 character(len=*),
optional,
intent(in) :: dim_name
550 character(len=*),
optional,
intent(in) :: unit
552 character(len=*),
optional,
intent(in) :: long_name
554 character(len=*),
optional,
intent(in) :: standard_name
556 character(len=*),
optional,
intent(in) :: description
558 integer :: varid, start(1), count(1), use_dimids(1)
560 if (
present(dimids))
then 562 elseif (
present(dim_name))
then 565 call die_msg(891890123,
"either dimids or dim_name must be present")
568 call pmc_nc_check(nf90_def_var(ncid, name, nf90_double, use_dimids, varid))
574 count = (/
size(var, 1) /)
576 start = start, count = count))
584 long_name, standard_name, description)
587 integer,
intent(in) :: ncid
589 integer,
intent(in) :: var(:)
591 character(len=*),
intent(in) :: name
594 integer,
optional,
intent(in) :: dimids(1)
597 character(len=*),
optional,
intent(in) :: dim_name
599 character(len=*),
optional,
intent(in) :: unit
601 character(len=*),
optional,
intent(in) :: long_name
603 character(len=*),
optional,
intent(in) :: standard_name
605 character(len=*),
optional,
intent(in) :: description
607 integer :: varid, start(1), count(1), use_dimids(1)
609 if (
present(dimids))
then 611 elseif (
present(dim_name))
then 614 call die_msg(464170526,
"either dimids or dim_name must be present")
617 call pmc_nc_check(nf90_def_var(ncid, name, nf90_int, use_dimids, varid))
623 count = (/
size(var, 1) /)
625 start = start, count = count))
633 dim_name_2, unit, long_name, standard_name, description)
636 integer,
intent(in) :: ncid
638 real(kind=dp),
intent(in) :: var(:,:)
640 character(len=*),
intent(in) :: name
643 integer,
optional,
intent(in) :: dimids(2)
646 character(len=*),
optional,
intent(in) :: dim_name_1
649 character(len=*),
optional,
intent(in) :: dim_name_2
651 character(len=*),
optional,
intent(in) :: unit
653 character(len=*),
optional,
intent(in) :: long_name
655 character(len=*),
optional,
intent(in) :: standard_name
657 character(len=*),
optional,
intent(in) :: description
659 integer :: varid, start(2), count(2), use_dimids(2)
661 if (
present(dimids))
then 663 elseif (
present(dim_name_1) .and.
present(dim_name_2))
then 668 "either dimids or both dim_name_1 and dim_name_2 must be present")
671 call pmc_nc_check(nf90_def_var(ncid, name, nf90_double, use_dimids, varid))
677 count = (/
size(var, 1),
size(var, 2) /)
679 start = start, count = count))
687 dim_name_2, unit, long_name, standard_name, description)
690 integer,
intent(in) :: ncid
692 integer,
intent(in) :: var(:,:)
694 character(len=*),
intent(in) :: name
697 integer,
optional,
intent(in) :: dimids(2)
700 character(len=*),
optional,
intent(in) :: dim_name_1
703 character(len=*),
optional,
intent(in) :: dim_name_2
705 character(len=*),
optional,
intent(in) :: unit
707 character(len=*),
optional,
intent(in) :: long_name
709 character(len=*),
optional,
intent(in) :: standard_name
711 character(len=*),
optional,
intent(in) :: description
713 integer :: varid, start(2), count(2), use_dimids(2)
715 if (
present(dimids))
then 717 elseif (
present(dim_name_1) .and.
present(dim_name_2))
then 722 "either dimids or both dim_name_1 and dim_name_2 must be present")
725 call pmc_nc_check(nf90_def_var(ncid, name, nf90_int, use_dimids, varid))
731 count = (/
size(var, 1),
size(var, 2) /)
733 start = start, count = count))
subroutine assert_msg(code, condition_ok, error_msg)
Errors unless condition_ok is true.
Random number generators.
subroutine pmc_nc_write_info(ncid, uuid, source, write_rank, write_n_proc)
Write basic information to a NetCDF file.
Wrapper functions for NetCDF. These all take a NetCDF ncid in data mode and return with it again in d...
subroutine pmc_nc_read_real_1d(ncid, var, name, must_be_present)
Read a simple real array from a NetCDF file.
integer function pmc_mpi_rank()
Returns the rank of the current process.
subroutine ensure_real_array_size(x, n, only_grow)
Allocate or reallocate the given array to ensure it is of the given size, preserving any data and/or ...
subroutine ensure_integer_array_size(x, n, only_grow)
Allocate or reallocate the given array to ensure it is of the given size, preserving any data and/or ...
subroutine pmc_nc_open_read(filename, ncid)
Open a NetCDF file for reading.
subroutine iso8601_date_and_time(date_time)
Current date and time in ISO 8601 format.
subroutine pmc_nc_close(ncid)
Close a NetCDF file.
subroutine pmc_nc_open_write(filename, ncid)
Open a NetCDF file for writing.
subroutine ensure_real_array_2d_size(x, n1, n2, only_grow)
Allocate or reallocate the given array to ensure it is of the given size, preserving any data and/or ...
subroutine pmc_nc_write_atts(ncid, varid, unit, long_name, standard_name, description)
Write attributes for a variable to a NetCDF file.
subroutine pmc_nc_read_integer(ncid, var, name, must_be_present)
Read a single integer from a NetCDF file.
subroutine pmc_nc_check(status)
Check the status of a NetCDF function call.
subroutine die_msg(code, error_msg)
Error immediately.
subroutine pmc_nc_write_integer(ncid, var, name, unit, long_name, standard_name, description)
Write a single integer to a NetCDF file.
subroutine pmc_nc_read_real(ncid, var, name, must_be_present)
Read a single real from a NetCDF file.
subroutine pmc_nc_write_real_2d(ncid, var, name, dimids, dim_name_1, dim_name_2, unit, long_name, standard_name, description)
Write a simple real 2D array to a NetCDF file.
subroutine pmc_nc_read_real_2d(ncid, var, name, must_be_present)
Read a simple real 2D array from a NetCDF file.
subroutine pmc_nc_read_integer_2d(ncid, var, name, must_be_present)
Read a simple integer 2D array from a NetCDF file.
subroutine ensure_integer_array_2d_size(x, n1, n2, only_grow)
Allocate or reallocate the given array to ensure it is of the given size, preserving any data and/or ...
subroutine pmc_nc_write_real(ncid, var, name, unit, long_name, standard_name, description)
Write a single real to a NetCDF file.
subroutine pmc_nc_write_real_1d(ncid, var, name, dimids, dim_name, unit, long_name, standard_name, description)
Write a simple real array to a NetCDF file.
subroutine pmc_nc_write_integer_2d(ncid, var, name, dimids, dim_name_1, dim_name_2, unit, long_name, standard_name, description)
Write a simple integer 2D array to a NetCDF file.
subroutine pmc_nc_ensure_dim(ncid, dim_name, dimid, dim_size, array_dim)
Create a dimension if necessary, or check its size if it already exists. In any case return the dimid...
subroutine pmc_nc_read_integer_1d(ncid, var, name, must_be_present)
Read a simple integer array from a NetCDF file.
character(len=pmc_util_convert_string_len) function integer_to_string(val)
Convert an integer to a string format.
integer function pmc_mpi_size()
Returns the total number of processes.
subroutine pmc_nc_check_msg(status, error_msg)
Check the status of a NetCDF function call and prints the given error message on failure.
Common utility subroutines.
subroutine pmc_nc_write_integer_1d(ncid, var, name, dimids, dim_name, unit, long_name, standard_name, description)
Write a simple integer array to a NetCDF file.