26 subroutine char_read_0d(fileobj, variable_name, buf, corner, append_error_msg, err, varid)
27 class(fmsnetcdffile_t),
intent(in) :: fileobj
28 character(len=*),
intent(in) :: variable_name
29 character(len=*),
intent(inout) :: buf
30 integer,
intent(in),
optional :: corner
32 character(len=200),
intent(in):: append_error_msg
33 integer,
intent(inout) :: err
34 integer,
intent(in) :: varid
36 integer,
dimension(2) :: start
38 character,
dimension(:),
allocatable :: charbuf
39 integer,
dimension(:),
allocatable :: dimsizes
43 ndims = get_variable_num_dimensions(fileobj, variable_name, broadcast=.false.)
44 allocate(dimsizes(ndims))
45 call get_variable_size(fileobj, variable_name, dimsizes, broadcast=.false.)
46 allocate(charbuf(dimsizes(1)))
48 if (ndims .eq. 2)
then
49 if (
present(corner))
then
53 elseif (ndims .gt. 2)
then
54 call error(
"Only scalar and 1d string values are currently supported: "//trim(append_error_msg))
56 err = nf90_get_var(fileobj%ncid, varid, charbuf, start=start, count=dimsizes)
57 if (len(buf) .lt. dimsizes(1))
then
58 call error(
"character buffer is too small; increase length: "//trim(append_error_msg))
62 if (charbuf(i) .eq. char(0))
then
72 subroutine char_read_1d(fileobj, variable_name, buf, c, append_error_msg, err, varid)
74 class(fmsnetcdffile_t),
intent(in) :: fileobj
75 character(len=*),
intent(in) :: variable_name
76 character(len=*),
dimension(:),
intent(inout) :: buf
78 integer,
dimension(2),
intent(in) :: c
79 character(len=200),
intent(in) :: append_error_msg
80 integer,
intent(inout) :: err
81 integer,
intent(in) :: varid
84 integer,
dimension(2) :: start
85 integer,
dimension(2) :: dimsizes
86 character,
dimension(:,:),
allocatable :: charbuf
87 character(len=1024) :: sbuf
91 ndims = get_variable_num_dimensions(fileobj, variable_name, broadcast=.false.)
92 if (ndims .ne. 2)
then
93 call error(trim(variable_name)//
"must be 2d dimensional in netcdf file.")
97 call get_variable_size(fileobj, variable_name, dimsizes, .false.)
98 dimsizes(2) = dimsizes(2) - start(2) + 1
99 call allocate_array(charbuf, dimsizes, initialize=.true.)
100 err = nf90_get_var(fileobj%ncid, varid, charbuf, start=start, count=dimsizes)
101 if (len(buf(1)) .lt. dimsizes(1))
then
102 call error(
"character buffer is too small; increase length: "//trim(append_error_msg))
104 if (
size(buf) .lt. dimsizes(2))
then
105 call error(
"incorrect buffer array size:: "//trim(append_error_msg))
107 do i = start(2), start(2)+dimsizes(2)-1
109 do j = 1, dimsizes(1)
110 if (charbuf(j, i-start(2)+1) .eq. char(0))
then
113 sbuf(j:j) = charbuf(j, i-start(2)+1)
115 call string_copy(buf(i-start(2)+1), sbuf)
125 class(fmsnetcdffile_t),
intent(in) :: fileobj
126 character(len=*),
intent(in) :: variable_name
127 class(*),
intent(inout) :: buf
128 integer,
intent(in),
optional :: unlim_dim_level
129 integer,
intent(in),
optional :: corner
131 logical,
intent(in),
optional :: broadcast
141 integer :: unlim_dim_index
142 integer,
dimension(1) :: c
143 character(len=1024),
dimension(1) :: buf1d
144 character(len=200) :: append_error_msg
146 append_error_msg =
"netcdf_read_data_0d: file:"//trim(fileobj%path)//
"- variable:"//trim(variable_name)
148 if (
present(broadcast))
then
154 if (
present(unlim_dim_level))
then
155 unlim_dim_index = get_variable_unlimited_dimension_index(fileobj, variable_name, &
157 if (unlim_dim_index .ne. 1)
then
158 call error(
"unlimited dimension must be the slowest varying dimension: "//trim(append_error_msg))
160 c(unlim_dim_index) = unlim_dim_level
162 if (fileobj%is_root)
then
163 varid = get_variable_id(fileobj%ncid, trim(variable_name), msg=append_error_msg)
165 type is (
integer(kind=i4_kind))
166 err = nf90_get_var(fileobj%ncid, varid, buf, start=c)
167 type is (
integer(kind=i8_kind))
168 err = nf90_get_var(fileobj%ncid, varid, buf, start=c)
169 type is (real(kind=r4_kind))
170 err = nf90_get_var(fileobj%ncid, varid, buf, start=c)
171 type is (real(kind=r8_kind))
172 err = nf90_get_var(fileobj%ncid, varid, buf, start=c)
173 type is (
character(len=*))
174 call char_read_0d(fileobj, variable_name, buf, corner, append_error_msg, err, varid)
176 call error(
"Unsupported variable type: "//trim(append_error_msg))
178 call check_netcdf_code(err, append_error_msg)
183 type is (
integer(kind=i4_kind))
184 call mpp_broadcast(buf, fileobj%io_root, pelist=fileobj%pelist)
185 type is (
integer(kind=i8_kind))
186 call mpp_broadcast(buf, fileobj%io_root, pelist=fileobj%pelist)
187 type is (real(kind=r4_kind))
188 call mpp_broadcast(buf, fileobj%io_root, pelist=fileobj%pelist)
189 type is (real(kind=r8_kind))
190 call mpp_broadcast(buf, fileobj%io_root, pelist=fileobj%pelist)
191 type is (
character(len=*))
192 call string_copy(buf1d(1), buf)
193 call mpp_broadcast(buf1d, len(buf1d(1)), fileobj%io_root, pelist=fileobj%pelist)
194 call string_copy(buf, buf1d(1))
196 call error(
"Unsupported variable type: "//trim(append_error_msg))
204 corner, edge_lengths, broadcast)
206 class(fmsnetcdffile_t),
intent(in) :: fileobj
207 character(len=*),
intent(in) :: variable_name
208 class(*),
dimension(:),
intent(inout) :: buf
210 integer,
intent(in),
optional :: unlim_dim_level
212 integer,
dimension(:),
intent(in),
optional :: corner
216 integer,
dimension(:),
intent(in),
optional :: edge_lengths
220 logical,
intent(in),
optional :: broadcast
230 integer :: unlim_dim_index
231 integer,
dimension(2) :: c
232 integer,
dimension(2) :: e
233 character(len=200) :: append_error_msg
235 append_error_msg =
"netcdf_read_data_1d: file:"//trim(fileobj%path)//
"- variable:"//trim(variable_name)
237 if (
present(broadcast))
then
243 if (
present(corner))
then
247 if (
present(edge_lengths))
then
248 e(1) = edge_lengths(1)
252 if (
present(unlim_dim_level))
then
253 unlim_dim_index = get_variable_unlimited_dimension_index(fileobj, variable_name, &
255 if (unlim_dim_index .ne. 2)
then
256 call error(
"unlimited dimension must be the slowest varying dimension: "//trim(append_error_msg))
258 c(unlim_dim_index) = unlim_dim_level
260 if (fileobj%is_root)
then
261 varid = get_variable_id(fileobj%ncid, trim(variable_name), msg=append_error_msg)
263 type is (
integer(kind=i4_kind))
264 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
265 type is (
integer(kind=i8_kind))
266 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
267 type is (real(kind=r4_kind))
268 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
269 type is (real(kind=r8_kind))
270 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
271 type is (
character(len=*))
272 call char_read_1d(fileobj, variable_name, buf, c, append_error_msg, err, varid)
274 call error(
"Unsupported variable type: "//trim(append_error_msg))
276 call check_netcdf_code(err, append_error_msg)
281 type is (
integer(kind=i4_kind))
282 call mpp_broadcast(buf,
size(buf), fileobj%io_root, pelist=fileobj%pelist)
283 type is (
integer(kind=i8_kind))
284 call mpp_broadcast(buf,
size(buf), fileobj%io_root, pelist=fileobj%pelist)
285 type is (real(kind=r4_kind))
286 call mpp_broadcast(buf,
size(buf), fileobj%io_root, pelist=fileobj%pelist)
287 type is (real(kind=r8_kind))
288 call mpp_broadcast(buf,
size(buf), fileobj%io_root, pelist=fileobj%pelist)
289 type is (
character(len=*))
290 call mpp_broadcast(buf, len(buf(1)), fileobj%io_root, pelist=fileobj%pelist)
292 call error(
"Unsupported variable type: "//trim(append_error_msg))
300 corner, edge_lengths, broadcast)
302 class(fmsnetcdffile_t),
intent(in) :: fileobj
303 character(len=*),
intent(in) :: variable_name
304 class(*),
dimension(:,:),
intent(inout) :: buf
306 integer,
intent(in),
optional :: unlim_dim_level
308 integer,
dimension(:),
intent(in),
optional :: corner
312 integer,
dimension(:),
intent(in),
optional :: edge_lengths
316 logical,
intent(in),
optional :: broadcast
326 integer :: unlim_dim_index
327 integer,
dimension(3) :: c
328 integer,
dimension(3) :: e
329 character(len=200) :: append_error_msg
331 append_error_msg =
"netcdf_read_data_2d: file:"//trim(fileobj%path)//
"- variable:"//trim(variable_name)
333 if (
present(broadcast))
then
339 if (
present(corner))
then
343 if (
present(edge_lengths))
then
344 e(1:2) = edge_lengths(1:2)
348 if (
present(unlim_dim_level))
then
349 unlim_dim_index = get_variable_unlimited_dimension_index(fileobj, variable_name, &
351 if (unlim_dim_index .ne. 3)
then
352 call error(
"unlimited dimension must be the slowest varying dimension: "//trim(append_error_msg))
354 c(unlim_dim_index) = unlim_dim_level
356 if(fileobj%use_collective)
then
357 varid = get_variable_id(fileobj%ncid, trim(variable_name), msg=append_error_msg)
360 err = nf90_var_par_access(fileobj%ncid, varid, nf90_collective)
361 call check_netcdf_code(err, append_error_msg)
363 type is (
integer(kind=i4_kind))
364 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
365 type is (
integer(kind=i8_kind))
366 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
367 type is (real(kind=r4_kind))
368 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
369 type is (real(kind=r8_kind))
370 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
372 call error(
"Unsupported variable type: "//trim(append_error_msg))
374 call check_netcdf_code(err, append_error_msg)
377 if (fileobj%is_root)
then
378 varid = get_variable_id(fileobj%ncid, trim(variable_name), msg=append_error_msg)
380 type is (
integer(kind=i4_kind))
381 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
382 type is (
integer(kind=i8_kind))
383 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
384 type is (real(kind=r4_kind))
385 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
386 type is (real(kind=r8_kind))
387 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
389 call error(
"Unsupported variable type: "//trim(append_error_msg))
391 call check_netcdf_code(err, append_error_msg)
396 type is (
integer(kind=i4_kind))
397 call mpp_broadcast(buf,
size(buf), fileobj%io_root, pelist=fileobj%pelist)
398 type is (
integer(kind=i8_kind))
399 call mpp_broadcast(buf,
size(buf), fileobj%io_root, pelist=fileobj%pelist)
400 type is (real(kind=r4_kind))
401 call mpp_broadcast(buf,
size(buf), fileobj%io_root, pelist=fileobj%pelist)
402 type is (real(kind=r8_kind))
403 call mpp_broadcast(buf,
size(buf), fileobj%io_root, pelist=fileobj%pelist)
405 call error(
"Unsupported variable type: "//trim(append_error_msg))
414 corner, edge_lengths, broadcast)
416 class(fmsnetcdffile_t),
intent(in) :: fileobj
417 character(len=*),
intent(in) :: variable_name
418 class(*),
dimension(:,:,:),
intent(inout) :: buf
420 integer,
intent(in),
optional :: unlim_dim_level
422 integer,
dimension(:),
intent(in),
optional :: corner
426 integer,
dimension(:),
intent(in),
optional :: edge_lengths
430 logical,
intent(in),
optional :: broadcast
440 integer :: unlim_dim_index
441 integer,
dimension(4) :: c
442 integer,
dimension(4) :: e
443 character(len=200) :: append_error_msg
445 append_error_msg =
"netcdf_read_data_3d: file:"//trim(fileobj%path)//
"- variable:"//trim(variable_name)
447 if (
present(broadcast))
then
453 if (
present(corner))
then
457 if (
present(edge_lengths))
then
458 e(1:3) = edge_lengths(1:3)
462 if (
present(unlim_dim_level))
then
463 unlim_dim_index = get_variable_unlimited_dimension_index(fileobj, variable_name, &
465 if (unlim_dim_index .ne. 4)
then
466 call error(
"unlimited dimension must be the slowest varying dimension: "//trim(append_error_msg))
468 c(unlim_dim_index) = unlim_dim_level
470 if(fileobj%use_collective)
then
471 varid = get_variable_id(fileobj%ncid, trim(variable_name), msg=append_error_msg)
474 err = nf90_var_par_access(fileobj%ncid, varid, nf90_collective)
475 call check_netcdf_code(err, append_error_msg)
477 type is (
integer(kind=i4_kind))
478 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
479 type is (
integer(kind=i8_kind))
480 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
481 type is (real(kind=r4_kind))
482 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
483 type is (real(kind=r8_kind))
484 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
486 call error(
"Unsupported variable type: "//trim(append_error_msg))
488 call check_netcdf_code(err, append_error_msg)
491 if (fileobj%is_root)
then
492 varid = get_variable_id(fileobj%ncid, trim(variable_name), msg=append_error_msg)
494 type is (
integer(kind=i4_kind))
495 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
496 type is (
integer(kind=i8_kind))
497 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
498 type is (real(kind=r4_kind))
499 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
500 type is (real(kind=r8_kind))
501 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
503 call error(
"Unsupported variable type: "//trim(append_error_msg))
505 call check_netcdf_code(err, append_error_msg)
510 type is (
integer(kind=i4_kind))
511 call mpp_broadcast(buf,
size(buf), fileobj%io_root, pelist=fileobj%pelist)
512 type is (
integer(kind=i8_kind))
513 call mpp_broadcast(buf,
size(buf), fileobj%io_root, pelist=fileobj%pelist)
514 type is (real(kind=r4_kind))
515 call mpp_broadcast(buf,
size(buf), fileobj%io_root, pelist=fileobj%pelist)
516 type is (real(kind=r8_kind))
517 call mpp_broadcast(buf,
size(buf), fileobj%io_root, pelist=fileobj%pelist)
519 call error(
"Unsupported variable type: "//trim(append_error_msg))
528 corner, edge_lengths, broadcast)
530 class(fmsnetcdffile_t),
intent(in) :: fileobj
531 character(len=*),
intent(in) :: variable_name
532 class(*),
dimension(:,:,:,:),
intent(inout) :: buf
534 integer,
intent(in),
optional :: unlim_dim_level
536 integer,
dimension(:),
intent(in),
optional :: corner
540 integer,
dimension(:),
intent(in),
optional :: edge_lengths
544 logical,
intent(in),
optional :: broadcast
554 integer :: unlim_dim_index
555 integer,
dimension(5) :: c
556 integer,
dimension(5) :: e
557 character(len=200) :: append_error_msg
559 append_error_msg =
"netcdf_read_data_4d: file:"//trim(fileobj%path)//
"- variable:"//trim(variable_name)
561 if (
present(broadcast))
then
567 if (
present(corner))
then
571 if (
present(edge_lengths))
then
572 e(1:4) = edge_lengths(1:4)
576 if (
present(unlim_dim_level))
then
577 unlim_dim_index = get_variable_unlimited_dimension_index(fileobj, variable_name, &
579 if (unlim_dim_index .ne. 5)
then
580 call error(
"unlimited dimension must be the slowest varying dimension: "//trim(append_error_msg))
582 c(unlim_dim_index) = unlim_dim_level
584 if (fileobj%is_root)
then
585 varid = get_variable_id(fileobj%ncid, trim(variable_name), msg=append_error_msg)
587 type is (
integer(kind=i4_kind))
588 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
589 type is (
integer(kind=i8_kind))
590 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
591 type is (real(kind=r4_kind))
592 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
593 type is (real(kind=r8_kind))
594 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
596 call error(
"Unsupported variable type: "//trim(append_error_msg))
598 call check_netcdf_code(err, append_error_msg)
603 type is (
integer(kind=i4_kind))
604 call mpp_broadcast(buf,
size(buf), fileobj%io_root, pelist=fileobj%pelist)
605 type is (
integer(kind=i8_kind))
606 call mpp_broadcast(buf,
size(buf), fileobj%io_root, pelist=fileobj%pelist)
607 type is (real(kind=r4_kind))
608 call mpp_broadcast(buf,
size(buf), fileobj%io_root, pelist=fileobj%pelist)
609 type is (real(kind=r8_kind))
610 call mpp_broadcast(buf,
size(buf), fileobj%io_root, pelist=fileobj%pelist)
612 call error(
"Unsupported variable type: "//trim(append_error_msg))
620 corner, edge_lengths, broadcast)
622 class(fmsnetcdffile_t),
intent(in) :: fileobj
623 character(len=*),
intent(in) :: variable_name
624 class(*),
dimension(:,:,:,:,:),
intent(inout) :: buf
626 integer,
intent(in),
optional :: unlim_dim_level
628 integer,
dimension(:),
intent(in),
optional :: corner
632 integer,
dimension(:),
intent(in),
optional :: edge_lengths
636 logical,
intent(in),
optional :: broadcast
646 integer :: unlim_dim_index
647 integer,
dimension(6) :: c
648 integer,
dimension(6) :: e
649 character(len=200) :: append_error_msg
651 append_error_msg =
"netcdf_read_data_5d: file:"//trim(fileobj%path)//
"- variable:"//trim(variable_name)
653 if (
present(broadcast))
then
659 if (
present(corner))
then
663 if (
present(edge_lengths))
then
664 e(1:5) = edge_lengths(1:5)
668 if (
present(unlim_dim_level))
then
669 unlim_dim_index = get_variable_unlimited_dimension_index(fileobj, variable_name, &
671 if (unlim_dim_index .ne. 6)
then
672 call error(
"unlimited dimension must be the slowest varying dimension: "//trim(append_error_msg))
674 c(unlim_dim_index) = unlim_dim_level
676 if (fileobj%is_root)
then
677 varid = get_variable_id(fileobj%ncid, trim(variable_name), msg=append_error_msg)
679 type is (
integer(kind=i4_kind))
680 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
681 type is (
integer(kind=i8_kind))
682 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
683 type is (real(kind=r4_kind))
684 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
685 type is (real(kind=r8_kind))
686 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
688 call error(
"Unsupported variable type: "//trim(append_error_msg))
690 call check_netcdf_code(err, append_error_msg)
695 type is (
integer(kind=i4_kind))
696 call mpp_broadcast(buf,
size(buf), fileobj%io_root, pelist=fileobj%pelist)
697 type is (
integer(kind=i8_kind))
698 call mpp_broadcast(buf,
size(buf), fileobj%io_root, pelist=fileobj%pelist)
699 type is (real(kind=r4_kind))
700 call mpp_broadcast(buf,
size(buf), fileobj%io_root, pelist=fileobj%pelist)
701 type is (real(kind=r8_kind))
702 call mpp_broadcast(buf,
size(buf), fileobj%io_root, pelist=fileobj%pelist)
704 call error(
"Unsupported variable type: "//trim(append_error_msg))
subroutine unpack_data_2d(fileobj, varid, varname, var_data)
subroutine unpack_data_1d(fileobj, varid, varname, var_data)
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 unpack_data_5d(fileobj, varid, varname, var_data)
subroutine unpack_data_0d(fileobj, varid, varname, var_data)
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_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 char_read_1d(fileobj, variable_name, buf, c, append_error_msg, err, varid)
Character read 1d function.
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 unpack_data_3d(fileobj, varid, varname, var_data)
subroutine char_read_0d(fileobj, variable_name, buf, corner, append_error_msg, err, varid)
Character read 0d function.
subroutine unpack_data_4d(fileobj, varid, varname, var_data)
subroutine netcdf_read_data_0d(fileobj, variable_name, buf, unlim_dim_level, corner, broadcast)
Read in data from a variable in a netcdf file.
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.