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(out) :: var(:)
225 character(len=*),
intent(in) :: name
228 logical,
optional,
intent(in) :: must_be_present
230 integer :: varid, status
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
246 "getting variable " // trim(name))
256 integer,
intent(in) :: ncid
258 integer,
intent(out) :: var(:)
260 character(len=*),
intent(in) :: name
263 logical,
optional,
intent(in) :: must_be_present
265 integer :: varid, status
266 logical :: use_must_be_present
268 if (present(must_be_present))
then
269 use_must_be_present = must_be_present
271 use_must_be_present = .true.
273 status = nf90_inq_varid(ncid, name, varid)
274 if ((.not. use_must_be_present) .and. (status == nf90_enotvar))
then
281 "getting variable " // trim(name))
291 integer,
intent(in) :: ncid
293 real(kind=dp),
intent(out) :: var(:,:)
295 character(len=*),
intent(in) :: name
298 logical,
optional,
intent(in) :: must_be_present
300 integer :: varid, status
301 logical :: use_must_be_present
303 if (present(must_be_present))
then
304 use_must_be_present = must_be_present
306 use_must_be_present = .true.
308 status = nf90_inq_varid(ncid, name, varid)
309 if ((.not. use_must_be_present) .and. (status == nf90_enotvar))
then
316 "getting variable " // trim(name))
326 integer,
intent(in) :: ncid
328 integer,
intent(out) :: var(:,:)
330 character(len=*),
intent(in) :: name
333 logical,
optional,
intent(in) :: must_be_present
335 integer :: varid, status
336 logical :: use_must_be_present
338 if (present(must_be_present))
then
339 use_must_be_present = must_be_present
341 use_must_be_present = .true.
343 status = nf90_inq_varid(ncid, name, varid)
344 if ((.not. use_must_be_present) .and. (status == nf90_enotvar))
then
351 "getting variable " // trim(name))
362 integer,
intent(in) :: ncid
364 integer,
intent(in) :: varid
366 character(len=*),
optional,
intent(in) :: unit
368 character(len=*),
optional,
intent(in) :: long_name
370 character(len=*),
optional,
intent(in) :: standard_name
372 character(len=*),
optional,
intent(in) :: description
374 if (present(unit))
then
375 call
pmc_nc_check(nf90_put_att(ncid, varid,
"unit", unit))
377 if (present(long_name))
then
378 call
pmc_nc_check(nf90_put_att(ncid, varid,
"long_name", long_name))
380 if (present(standard_name))
then
381 call
pmc_nc_check(nf90_put_att(ncid, varid,
"standard_name", &
384 if (present(description))
then
385 call
pmc_nc_check(nf90_put_att(ncid, varid,
"description", &
395 standard_name, description)
398 integer,
intent(in) :: ncid
400 real(kind=dp),
intent(in) :: var
402 character(len=*),
intent(in) :: name
404 character(len=*),
optional,
intent(in) :: unit
406 character(len=*),
optional,
intent(in) :: long_name
408 character(len=*),
optional,
intent(in) :: standard_name
410 character(len=*),
optional,
intent(in) :: description
412 integer :: varid, dimids(0)
415 call
pmc_nc_check(nf90_def_var(ncid, name, nf90_double, dimids, varid))
428 standard_name, description)
431 integer,
intent(in) :: ncid
433 integer,
intent(in) :: var
435 character(len=*),
intent(in) :: name
437 character(len=*),
optional,
intent(in) :: unit
439 character(len=*),
optional,
intent(in) :: long_name
441 character(len=*),
optional,
intent(in) :: standard_name
443 character(len=*),
optional,
intent(in) :: description
445 integer :: varid, dimids(0)
448 call
pmc_nc_check(nf90_def_var(ncid, name, nf90_int, dimids, varid))
464 integer,
intent(in) :: ncid
466 character(len=*),
intent(in) :: dim_name
468 integer,
intent(out) :: dimid
470 integer,
intent(in) :: dim_size
472 integer,
intent(in) :: array_dim
474 integer :: status, check_dim_size
475 character(len=1000) :: check_name
477 status = nf90_inq_dimid(ncid, dim_name, dimid)
478 if (status == nf90_noerr)
then
479 call
pmc_nc_check(nf90_inquire_dimension(ncid, dimid, check_name, &
481 call
assert_msg(657263912, check_dim_size == dim_size, &
484 //
" of data array does not match size " &
486 // trim(dim_name) //
"' dim")
489 if (status /= nf90_ebaddim)
then
495 call
pmc_nc_check(nf90_def_dim(ncid, dim_name, dim_size, dimid))
505 long_name, standard_name, description)
508 integer,
intent(in) :: ncid
510 real(kind=dp),
intent(in) :: var(:)
512 character(len=*),
intent(in) :: name
515 integer,
optional,
intent(in) :: dimids(1)
518 character(len=*),
optional,
intent(in) :: dim_name
520 character(len=*),
optional,
intent(in) :: unit
522 character(len=*),
optional,
intent(in) :: long_name
524 character(len=*),
optional,
intent(in) :: standard_name
526 character(len=*),
optional,
intent(in) :: description
528 integer :: varid, start(1), count(1), use_dimids(1)
530 if (present(dimids))
then
532 elseif (present(dim_name))
then
535 call
die_msg(891890123,
"either dimids or dim_name must be present")
538 call
pmc_nc_check(nf90_def_var(ncid, name, nf90_double, use_dimids, varid))
544 count = (/
size(var, 1) /)
546 start = start, count = count))
554 long_name, standard_name, description)
557 integer,
intent(in) :: ncid
559 integer,
intent(in) :: var(:)
561 character(len=*),
intent(in) :: name
564 integer,
optional,
intent(in) :: dimids(1)
567 character(len=*),
optional,
intent(in) :: dim_name
569 character(len=*),
optional,
intent(in) :: unit
571 character(len=*),
optional,
intent(in) :: long_name
573 character(len=*),
optional,
intent(in) :: standard_name
575 character(len=*),
optional,
intent(in) :: description
577 integer :: varid, start(1), count(1), use_dimids(1)
579 if (present(dimids))
then
581 elseif (present(dim_name))
then
584 call
die_msg(464170526,
"either dimids or dim_name must be present")
587 call
pmc_nc_check(nf90_def_var(ncid, name, nf90_int, use_dimids, varid))
593 count = (/
size(var, 1) /)
595 start = start, count = count))
603 dim_name_2, unit, long_name, standard_name, description)
606 integer,
intent(in) :: ncid
608 real(kind=dp),
intent(in) :: var(:,:)
610 character(len=*),
intent(in) :: name
613 integer,
optional,
intent(in) :: dimids(2)
616 character(len=*),
optional,
intent(in) :: dim_name_1
619 character(len=*),
optional,
intent(in) :: dim_name_2
621 character(len=*),
optional,
intent(in) :: unit
623 character(len=*),
optional,
intent(in) :: long_name
625 character(len=*),
optional,
intent(in) :: standard_name
627 character(len=*),
optional,
intent(in) :: description
629 integer :: varid, start(2), count(2), use_dimids(2)
631 if (present(dimids))
then
633 elseif (present(dim_name_1) .and. present(dim_name_2))
then
638 "either dimids or both dim_name_1 and dim_name_2 must be present")
641 call
pmc_nc_check(nf90_def_var(ncid, name, nf90_double, use_dimids, varid))
647 count = (/
size(var, 1),
size(var, 2) /)
649 start = start, count = count))
657 dim_name_2, unit, long_name, standard_name, description)
660 integer,
intent(in) :: ncid
662 integer,
intent(in) :: var(:,:)
664 character(len=*),
intent(in) :: name
667 integer,
optional,
intent(in) :: dimids(2)
670 character(len=*),
optional,
intent(in) :: dim_name_1
673 character(len=*),
optional,
intent(in) :: dim_name_2
675 character(len=*),
optional,
intent(in) :: unit
677 character(len=*),
optional,
intent(in) :: long_name
679 character(len=*),
optional,
intent(in) :: standard_name
681 character(len=*),
optional,
intent(in) :: description
683 integer :: varid, start(2), count(2), use_dimids(2)
685 if (present(dimids))
then
687 elseif (present(dim_name_1) .and. present(dim_name_2))
then
692 "either dimids or both dim_name_1 and dim_name_2 must be present")
695 call
pmc_nc_check(nf90_def_var(ncid, name, nf90_int, use_dimids, varid))
701 count = (/
size(var, 1),
size(var, 2) /)
703 start = start, count = count))
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 die_msg(code, error_msg)
Error immediately.
subroutine iso8601_date_and_time(date_time)
Current date and time in ISO 8601 format.
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_check(status)
Check the status of a NetCDF function call.
subroutine pmc_nc_read_integer(ncid, var, name, must_be_present)
Read a single integer from a NetCDF file.
subroutine assert_msg(code, condition_ok, error_msg)
Errors unless condition_ok is true.
subroutine pmc_nc_read_integer_1d(ncid, var, name, must_be_present)
Read a simple integer array from a NetCDF file.
Random number generators.
integer function pmc_mpi_size()
Returns the total number of processes.
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.
Common utility subroutines.
subroutine pmc_nc_read_real_1d(ncid, var, name, must_be_present)
Read a simple real array from 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.
integer function pmc_mpi_rank()
Returns the rank of the current process.
character(len=pmc_util_convert_string_len) function integer_to_string(val)
Convert an integer to a string format.
subroutine pmc_nc_write_info(ncid, uuid, source, write_rank, write_n_proc)
Write basic information to 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_write_real(ncid, var, name, unit, long_name, standard_name, description)
Write a single real to a NetCDF file.
subroutine pmc_nc_close(ncid)
Close 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...
Wrapper functions for NetCDF. These all take a NetCDF ncid in data mode and return with it again in d...
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_open_read(filename, ncid)
Open a NetCDF file for reading.
subroutine pmc_nc_open_write(filename, ncid)
Open a NetCDF file for writing.
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_read_integer_2d(ncid, var, name, must_be_present)
Read a simple integer 2D array from a NetCDF file.
subroutine pmc_nc_check_msg(status, error_msg)
Check the status of a NetCDF function call and prints the given error message on failure.