FMS 2025.01.02-dev
Flexible Modeling System
Loading...
Searching...
No Matches
netcdf_read_data.inc
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!> @file
20!> @brief Routines for retrieving variable values from netcdf files with different
21!! dimension sizes for the @ref netcdf_read_data interface
22
23!> @addtogroup netcdf_io_mod
24!> @{
25
26!> @brief Character read 0d function.
27subroutine char_read_0d(fileobj, variable_name, buf, corner, append_error_msg, err, varid)
28 class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
29 character(len=*), intent(in) :: variable_name !< Variable name.
30 character(len=*), intent(inout) :: buf !< Array that the data will be read into
31 integer, intent(in), optional :: corner !< Index of the string to read if the variable
32 !! contains a 1D array of strings.
33 character(len=200), intent(in):: append_error_msg !< Msg to be appended to FATAL error message
34 integer, intent(inout) :: err
35 integer, intent(in) :: varid
36
37 integer, dimension(2) :: start
38 integer :: ndims
39 character, dimension(:), allocatable :: charbuf
40 integer, dimension(:), allocatable :: dimsizes
41 integer :: i
42
43 start(:) = 1
44 ndims = get_variable_num_dimensions(fileobj, variable_name, broadcast=.false.)
45 allocate(dimsizes(ndims))
46 call get_variable_size(fileobj, variable_name, dimsizes, broadcast=.false.)
47 allocate(charbuf(dimsizes(1)))
48 charbuf(:) = ""
49 if (ndims .eq. 2) then
50 if (present(corner)) then
51 start(2) = corner
52 endif
53 dimsizes(2) = 1
54 elseif (ndims .gt. 2) then
55 call error("Only scalar and 1d string values are currently supported: "//trim(append_error_msg))
56 endif
57 err = nf90_get_var(fileobj%ncid, varid, charbuf, start=start, count=dimsizes)
58 if (len(buf) .lt. dimsizes(1)) then
59 call error("character buffer is too small; increase length: "//trim(append_error_msg))
60 endif
61 buf = ""
62 do i = 1, dimsizes(1)
63 if (charbuf(i) .eq. char(0)) then
64 exit
65 endif
66 buf(i:i) = charbuf(i)
67 enddo
68 deallocate(charbuf)
69 deallocate(dimsizes)
70end subroutine char_read_0d
71
72!> @brief Character read 1d function.
73subroutine char_read_1d(fileobj, variable_name, buf, c, append_error_msg, err, varid)
74
75 class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
76 character(len=*), intent(in) :: variable_name !< Variable name.
77 character(len=*), dimension(:), intent(inout) :: buf !< Array that the data
78 !! will be read into.
79 integer, dimension(2), intent(in) :: c
80 character(len=200), intent(in) :: append_error_msg !< Msg to be appended to FATAL error message
81 integer, intent(inout) :: err
82 integer, intent(in) :: varid
83
84 integer :: ndims
85 integer, dimension(2) :: start
86 integer, dimension(2) :: dimsizes
87 character, dimension(:,:), allocatable :: charbuf
88 character(len=1024) :: sbuf
89 integer :: i
90 integer :: j
91
92 ndims = get_variable_num_dimensions(fileobj, variable_name, broadcast=.false.)
93 if (ndims .ne. 2) then
94 call error(trim(variable_name)//"must be 2d dimensional in netcdf file.")
95 endif
96 start(1) = 1
97 start(2) = c(1)
98 call get_variable_size(fileobj, variable_name, dimsizes, .false.)
99 dimsizes(2) = dimsizes(2) - start(2) + 1
100 call allocate_array(charbuf, dimsizes, initialize=.true.)
101 err = nf90_get_var(fileobj%ncid, varid, charbuf, start=start, count=dimsizes)
102 if (len(buf(1)) .lt. dimsizes(1)) then
103 call error("character buffer is too small; increase length: "//trim(append_error_msg))
104 endif
105 if (size(buf) .lt. dimsizes(2)) then
106 call error("incorrect buffer array size:: "//trim(append_error_msg))
107 endif
108 do i = start(2), start(2)+dimsizes(2)-1
109 sbuf = ""
110 do j = 1, dimsizes(1)
111 if (charbuf(j, i-start(2)+1) .eq. char(0)) then
112 exit
113 endif
114 sbuf(j:j) = charbuf(j, i-start(2)+1)
115 enddo
116 call string_copy(buf(i-start(2)+1), sbuf)
117 enddo
118 deallocate(charbuf)
119
120end subroutine char_read_1d
121
122!> @brief Read in data from a variable in a netcdf file.
123subroutine netcdf_read_data_0d(fileobj, variable_name, buf, unlim_dim_level, &
124 corner, broadcast)
125
126 class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
127 character(len=*), intent(in) :: variable_name !< Variable name.
128 class(*), intent(inout) :: buf !< Array that the data will be read into.
129 integer, intent(in), optional :: unlim_dim_level !< Unlimited dimension level.
130 integer, intent(in), optional :: corner !< Index of the string to read if the variable
131 !! contains a 1D array of strings.
132 logical, intent(in), optional :: broadcast !< Flag controlling whether or
133 !! not the data will be
134 !! broadcasted to non
135 !! "I/O root" ranks.
136 !! The broadcast will be done
137 !! by default.
138
139 logical :: bcast
140 integer :: err
141 integer :: varid
142 integer :: unlim_dim_index
143 integer, dimension(1) :: c
144 character(len=1024), dimension(1) :: buf1d
145 character(len=200) :: append_error_msg !< Msg to be appended to FATAL error message
146
147 append_error_msg = "netcdf_read_data_0d: file:"//trim(fileobj%path)//"- variable:"//trim(variable_name)
148
149 if (present(broadcast)) then
150 bcast = broadcast
151 else
152 bcast = .true.
153 endif
154 c(1) = 1
155 if (present(unlim_dim_level)) then
156 unlim_dim_index = get_variable_unlimited_dimension_index(fileobj, variable_name, &
157 broadcast=bcast)
158 if (unlim_dim_index .ne. 1) then
159 call error("unlimited dimension must be the slowest varying dimension: "//trim(append_error_msg))
160 endif
161 c(unlim_dim_index) = unlim_dim_level
162 endif
163 if (fileobj%is_root) then
164 varid = get_variable_id(fileobj%ncid, trim(variable_name), msg=append_error_msg)
165 select type(buf)
166 type is (integer(kind=i4_kind))
167 err = nf90_get_var(fileobj%ncid, varid, buf, start=c)
168 type is (integer(kind=i8_kind))
169 err = nf90_get_var(fileobj%ncid, varid, buf, start=c)
170 type is (real(kind=r4_kind))
171 err = nf90_get_var(fileobj%ncid, varid, buf, start=c)
172 type is (real(kind=r8_kind))
173 err = nf90_get_var(fileobj%ncid, varid, buf, start=c)
174 type is (character(len=*))
175 call char_read_0d(fileobj, variable_name, buf, corner, append_error_msg, err, varid)
176 class default
177 call error("Unsupported variable type: "//trim(append_error_msg))
178 end select
179 call check_netcdf_code(err, append_error_msg)
180 call unpack_data_0d(fileobj, varid, variable_name, buf)
181 endif
182 if (bcast) then
183 select type(buf)
184 type is (integer(kind=i4_kind))
185 call mpp_broadcast(buf, fileobj%io_root, pelist=fileobj%pelist)
186 type is (integer(kind=i8_kind))
187 call mpp_broadcast(buf, fileobj%io_root, pelist=fileobj%pelist)
188 type is (real(kind=r4_kind))
189 call mpp_broadcast(buf, fileobj%io_root, pelist=fileobj%pelist)
190 type is (real(kind=r8_kind))
191 call mpp_broadcast(buf, fileobj%io_root, pelist=fileobj%pelist)
192 type is (character(len=*))
193 call string_copy(buf1d(1), buf)
194 call mpp_broadcast(buf1d, len(buf1d(1)), fileobj%io_root, pelist=fileobj%pelist)
195 call string_copy(buf, buf1d(1))
196 class default
197 call error("Unsupported variable type: "//trim(append_error_msg))
198 end select
199 endif
200end subroutine netcdf_read_data_0d
201
202
203!> @brief Read in data from a variable in a netcdf file.
204subroutine netcdf_read_data_1d(fileobj, variable_name, buf, unlim_dim_level, &
205 corner, edge_lengths, broadcast)
206
207 class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
208 character(len=*), intent(in) :: variable_name !< Variable name.
209 class(*), dimension(:), intent(inout) :: buf !< Array that the data
210 !! will be read into.
211 integer, intent(in), optional :: unlim_dim_level !< Unlimited dimension
212 !! level.
213 integer, dimension(:), intent(in), optional :: corner !< Array of starting
214 !! indices describing
215 !! where the data
216 !! will be read from.
217 integer, dimension(:), intent(in), optional :: edge_lengths !< The number of
218 !! elements that
219 !! will be read
220 !! in each dimension.
221 logical, intent(in), optional :: broadcast !< Flag controlling whether or
222 !! not the data will be
223 !! broadcasted to non
224 !! "I/O root" ranks.
225 !! The broadcast will be done
226 !! by default.
227
228 logical :: bcast
229 integer :: err
230 integer :: varid
231 integer :: unlim_dim_index
232 integer, dimension(2) :: c
233 integer, dimension(2) :: e
234 character(len=200) :: append_error_msg !< Msg to be appended to FATAL error message
235
236 append_error_msg = "netcdf_read_data_1d: file:"//trim(fileobj%path)//"- variable:"//trim(variable_name)
237
238 if (present(broadcast)) then
239 bcast = broadcast
240 else
241 bcast = .true.
242 endif
243 c(:) = 1
244 if (present(corner)) then
245 c(1) = corner(1)
246 endif
247 e(:) = 1
248 if (present(edge_lengths)) then
249 e(1) = edge_lengths(1)
250 else
251 e(1) = size(buf)
252 endif
253 if (present(unlim_dim_level)) then
254 unlim_dim_index = get_variable_unlimited_dimension_index(fileobj, variable_name, &
255 broadcast=bcast)
256 if (unlim_dim_index .ne. 2) then
257 call error("unlimited dimension must be the slowest varying dimension: "//trim(append_error_msg))
258 endif
259 c(unlim_dim_index) = unlim_dim_level
260 endif
261 if (fileobj%is_root) then
262 varid = get_variable_id(fileobj%ncid, trim(variable_name), msg=append_error_msg)
263 select type(buf)
264 type is (integer(kind=i4_kind))
265 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
266 type is (integer(kind=i8_kind))
267 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
268 type is (real(kind=r4_kind))
269 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
270 type is (real(kind=r8_kind))
271 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
272 type is (character(len=*))
273 call char_read_1d(fileobj, variable_name, buf, c, append_error_msg, err, varid)
274 class default
275 call error("Unsupported variable type: "//trim(append_error_msg))
276 end select
277 call check_netcdf_code(err, append_error_msg)
278 call unpack_data_1d(fileobj, varid, variable_name, buf)
279 endif
280 if (bcast) then
281 select type(buf)
282 type is (integer(kind=i4_kind))
283 call mpp_broadcast(buf, size(buf), fileobj%io_root, pelist=fileobj%pelist)
284 type is (integer(kind=i8_kind))
285 call mpp_broadcast(buf, size(buf), fileobj%io_root, pelist=fileobj%pelist)
286 type is (real(kind=r4_kind))
287 call mpp_broadcast(buf, size(buf), fileobj%io_root, pelist=fileobj%pelist)
288 type is (real(kind=r8_kind))
289 call mpp_broadcast(buf, size(buf), fileobj%io_root, pelist=fileobj%pelist)
290 type is (character(len=*))
291 call mpp_broadcast(buf, len(buf(1)), fileobj%io_root, pelist=fileobj%pelist)
292 class default
293 call error("Unsupported variable type: "//trim(append_error_msg))
294 end select
295 endif
296end subroutine netcdf_read_data_1d
297
298
299!> @brief Read in data from a variable in a netcdf file.
300subroutine netcdf_read_data_2d(fileobj, variable_name, buf, unlim_dim_level, &
301 corner, edge_lengths, broadcast)
302
303 class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
304 character(len=*), intent(in) :: variable_name !< Variable name.
305 class(*), dimension(:,:), intent(inout) :: buf !< Array that the data
306 !! will be read into.
307 integer, intent(in), optional :: unlim_dim_level !< Unlimited dimension
308 !! level.
309 integer, dimension(:), intent(in), optional :: corner !< Array of starting
310 !! indices describing
311 !! where the data
312 !! will be read from.
313 integer, dimension(:), intent(in), optional :: edge_lengths !< The number of
314 !! elements that
315 !! will be read
316 !! in each dimension.
317 logical, intent(in), optional :: broadcast !< Flag controlling whether or
318 !! not the data will be
319 !! broadcasted to non
320 !! "I/O root" ranks.
321 !! The broadcast will be done
322 !! by default.
323
324 logical :: bcast
325 integer :: err
326 integer :: varid
327 integer :: unlim_dim_index
328 integer, dimension(3) :: c
329 integer, dimension(3) :: e
330 character(len=200) :: append_error_msg !< Msg to be appended to FATAL error message
331
332 append_error_msg = "netcdf_read_data_2d: file:"//trim(fileobj%path)//"- variable:"//trim(variable_name)
333
334 if (present(broadcast)) then
335 bcast = broadcast
336 else
337 bcast = .true.
338 endif
339 c(:) = 1
340 if (present(corner)) then
341 c(1:2) = corner(1:2)
342 endif
343 e(:) = 1
344 if (present(edge_lengths)) then
345 e(1:2) = edge_lengths(1:2)
346 else
347 e(1:2) = shape(buf)
348 endif
349 if (present(unlim_dim_level)) then
350 unlim_dim_index = get_variable_unlimited_dimension_index(fileobj, variable_name, &
351 broadcast=bcast)
352 if (unlim_dim_index .ne. 3) then
353 call error("unlimited dimension must be the slowest varying dimension: "//trim(append_error_msg))
354 endif
355 c(unlim_dim_index) = unlim_dim_level
356 endif
357 if(fileobj%use_collective) then
358 varid = get_variable_id(fileobj%ncid, trim(variable_name), msg=append_error_msg)
359 ! NetCDF does not have the ability to specify collective I/O at
360 ! the file basis so we must activate at the variable level
361 err = nf90_var_par_access(fileobj%ncid, varid, nf90_collective)
362 call check_netcdf_code(err, append_error_msg)
363 select type(buf)
364 type is (integer(kind=i4_kind))
365 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
366 type is (integer(kind=i8_kind))
367 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
368 type is (real(kind=r4_kind))
369 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
370 type is (real(kind=r8_kind))
371 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
372 class default
373 call error("Unsupported variable type: "//trim(append_error_msg))
374 end select
375 call check_netcdf_code(err, append_error_msg)
376 call unpack_data_2d(fileobj, varid, variable_name, buf)
377 else
378 if (fileobj%is_root) then
379 varid = get_variable_id(fileobj%ncid, trim(variable_name), msg=append_error_msg)
380 select type(buf)
381 type is (integer(kind=i4_kind))
382 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
383 type is (integer(kind=i8_kind))
384 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
385 type is (real(kind=r4_kind))
386 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
387 type is (real(kind=r8_kind))
388 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
389 class default
390 call error("Unsupported variable type: "//trim(append_error_msg))
391 end select
392 call check_netcdf_code(err, append_error_msg)
393 call unpack_data_2d(fileobj, varid, variable_name, buf)
394 endif
395 if (bcast) then
396 select type(buf)
397 type is (integer(kind=i4_kind))
398 call mpp_broadcast(buf, size(buf), fileobj%io_root, pelist=fileobj%pelist)
399 type is (integer(kind=i8_kind))
400 call mpp_broadcast(buf, size(buf), fileobj%io_root, pelist=fileobj%pelist)
401 type is (real(kind=r4_kind))
402 call mpp_broadcast(buf, size(buf), fileobj%io_root, pelist=fileobj%pelist)
403 type is (real(kind=r8_kind))
404 call mpp_broadcast(buf, size(buf), fileobj%io_root, pelist=fileobj%pelist)
405 class default
406 call error("Unsupported variable type: "//trim(append_error_msg))
407 end select
408 endif
409 endif
410end subroutine netcdf_read_data_2d
411
412
413!> @brief Read in data from a variable in a netcdf file.
414subroutine netcdf_read_data_3d(fileobj, variable_name, buf, unlim_dim_level, &
415 corner, edge_lengths, broadcast)
416
417 class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
418 character(len=*), intent(in) :: variable_name !< Variable name.
419 class(*), dimension(:,:,:), intent(inout) :: buf !< Array that the data
420 !! will be read into.
421 integer, intent(in), optional :: unlim_dim_level !< Unlimited dimension
422 !! level.
423 integer, dimension(:), intent(in), optional :: corner !< Array of starting
424 !! indices describing
425 !! where the data
426 !! will be read from.
427 integer, dimension(:), intent(in), optional :: edge_lengths !< The number of
428 !! elements that
429 !! will be read
430 !! in each dimension.
431 logical, intent(in), optional :: broadcast !< Flag controlling whether or
432 !! not the data will be
433 !! broadcasted to non
434 !! "I/O root" ranks.
435 !! The broadcast will be done
436 !! by default.
437
438 logical :: bcast
439 integer :: err
440 integer :: varid
441 integer :: unlim_dim_index
442 integer, dimension(4) :: c
443 integer, dimension(4) :: e
444 character(len=200) :: append_error_msg !< Msg to be appended to FATAL error message
445
446 append_error_msg = "netcdf_read_data_3d: file:"//trim(fileobj%path)//"- variable:"//trim(variable_name)
447
448 if (present(broadcast)) then
449 bcast = broadcast
450 else
451 bcast = .true.
452 endif
453 c(:) = 1
454 if (present(corner)) then
455 c(1:3) = corner(1:3)
456 endif
457 e(:) = 1
458 if (present(edge_lengths)) then
459 e(1:3) = edge_lengths(1:3)
460 else
461 e(1:3) = shape(buf)
462 endif
463 if (present(unlim_dim_level)) then
464 unlim_dim_index = get_variable_unlimited_dimension_index(fileobj, variable_name, &
465 broadcast=bcast)
466 if (unlim_dim_index .ne. 4) then
467 call error("unlimited dimension must be the slowest varying dimension: "//trim(append_error_msg))
468 endif
469 c(unlim_dim_index) = unlim_dim_level
470 endif
471 if(fileobj%use_collective) then
472 varid = get_variable_id(fileobj%ncid, trim(variable_name), msg=append_error_msg)
473 ! NetCDF does not have the ability to specify collective I/O at
474 ! the file basis so we must activate at the variable level
475 err = nf90_var_par_access(fileobj%ncid, varid, nf90_collective)
476 call check_netcdf_code(err, append_error_msg)
477 select type(buf)
478 type is (integer(kind=i4_kind))
479 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
480 type is (integer(kind=i8_kind))
481 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
482 type is (real(kind=r4_kind))
483 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
484 type is (real(kind=r8_kind))
485 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
486 class default
487 call error("Unsupported variable type: "//trim(append_error_msg))
488 end select
489 call check_netcdf_code(err, append_error_msg)
490 call unpack_data_3d(fileobj, varid, variable_name, buf)
491 else
492 if (fileobj%is_root) then
493 varid = get_variable_id(fileobj%ncid, trim(variable_name), msg=append_error_msg)
494 select type(buf)
495 type is (integer(kind=i4_kind))
496 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
497 type is (integer(kind=i8_kind))
498 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
499 type is (real(kind=r4_kind))
500 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
501 type is (real(kind=r8_kind))
502 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
503 class default
504 call error("Unsupported variable type: "//trim(append_error_msg))
505 end select
506 call check_netcdf_code(err, append_error_msg)
507 call unpack_data_3d(fileobj, varid, variable_name, buf)
508 endif
509 if (bcast) then
510 select type(buf)
511 type is (integer(kind=i4_kind))
512 call mpp_broadcast(buf, size(buf), fileobj%io_root, pelist=fileobj%pelist)
513 type is (integer(kind=i8_kind))
514 call mpp_broadcast(buf, size(buf), fileobj%io_root, pelist=fileobj%pelist)
515 type is (real(kind=r4_kind))
516 call mpp_broadcast(buf, size(buf), fileobj%io_root, pelist=fileobj%pelist)
517 type is (real(kind=r8_kind))
518 call mpp_broadcast(buf, size(buf), fileobj%io_root, pelist=fileobj%pelist)
519 class default
520 call error("Unsupported variable type: "//trim(append_error_msg))
521 end select
522 endif
523 endif
524end subroutine netcdf_read_data_3d
525
526
527!> @brief Read in data from a variable in a netcdf file.
528subroutine netcdf_read_data_4d(fileobj, variable_name, buf, unlim_dim_level, &
529 corner, edge_lengths, broadcast)
530
531 class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
532 character(len=*), intent(in) :: variable_name !< Variable name.
533 class(*), dimension(:,:,:,:), intent(inout) :: buf !< Array that the data
534 !! will be read into.
535 integer, intent(in), optional :: unlim_dim_level !< Unlimited dimension
536 !! level.
537 integer, dimension(:), intent(in), optional :: corner !< Array of starting
538 !! indices describing
539 !! where the data
540 !! will be read from.
541 integer, dimension(:), intent(in), optional :: edge_lengths !< The number of
542 !! elements that
543 !! will be read
544 !! in each dimension.
545 logical, intent(in), optional :: broadcast !< Flag controlling whether or
546 !! not the data will be
547 !! broadcasted to non
548 !! "I/O root" ranks.
549 !! The broadcast will be done
550 !! by default.
551
552 logical :: bcast
553 integer :: err
554 integer :: varid
555 integer :: unlim_dim_index
556 integer, dimension(5) :: c
557 integer, dimension(5) :: e
558 character(len=200) :: append_error_msg !< Msg to be appended to FATAL error message
559
560 append_error_msg = "netcdf_read_data_4d: file:"//trim(fileobj%path)//"- variable:"//trim(variable_name)
561
562 if (present(broadcast)) then
563 bcast = broadcast
564 else
565 bcast = .true.
566 endif
567 c(:) = 1
568 if (present(corner)) then
569 c(1:4) = corner(1:4)
570 endif
571 e(:) = 1
572 if (present(edge_lengths)) then
573 e(1:4) = edge_lengths(1:4)
574 else
575 e(1:4) = shape(buf)
576 endif
577 if (present(unlim_dim_level)) then
578 unlim_dim_index = get_variable_unlimited_dimension_index(fileobj, variable_name, &
579 broadcast=bcast)
580 if (unlim_dim_index .ne. 5) then
581 call error("unlimited dimension must be the slowest varying dimension: "//trim(append_error_msg))
582 endif
583 c(unlim_dim_index) = unlim_dim_level
584 endif
585 if (fileobj%is_root) then
586 varid = get_variable_id(fileobj%ncid, trim(variable_name), msg=append_error_msg)
587 select type(buf)
588 type is (integer(kind=i4_kind))
589 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
590 type is (integer(kind=i8_kind))
591 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
592 type is (real(kind=r4_kind))
593 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
594 type is (real(kind=r8_kind))
595 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
596 class default
597 call error("Unsupported variable type: "//trim(append_error_msg))
598 end select
599 call check_netcdf_code(err, append_error_msg)
600 call unpack_data_4d(fileobj, varid, variable_name, buf)
601 endif
602 if (bcast) then
603 select type(buf)
604 type is (integer(kind=i4_kind))
605 call mpp_broadcast(buf, size(buf), fileobj%io_root, pelist=fileobj%pelist)
606 type is (integer(kind=i8_kind))
607 call mpp_broadcast(buf, size(buf), fileobj%io_root, pelist=fileobj%pelist)
608 type is (real(kind=r4_kind))
609 call mpp_broadcast(buf, size(buf), fileobj%io_root, pelist=fileobj%pelist)
610 type is (real(kind=r8_kind))
611 call mpp_broadcast(buf, size(buf), fileobj%io_root, pelist=fileobj%pelist)
612 class default
613 call error("Unsupported variable type: "//trim(append_error_msg))
614 end select
615 endif
616end subroutine netcdf_read_data_4d
617
618
619!> @brief Read in data from a variable in a netcdf file.
620subroutine netcdf_read_data_5d(fileobj, variable_name, buf, unlim_dim_level, &
621 corner, edge_lengths, broadcast)
622
623 class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
624 character(len=*), intent(in) :: variable_name !< Variable name.
625 class(*), dimension(:,:,:,:,:), intent(inout) :: buf !< Array that the data
626 !! will be read into.
627 integer, intent(in), optional :: unlim_dim_level !< Unlimited dimension
628 !! level.
629 integer, dimension(:), intent(in), optional :: corner !< Array of starting
630 !! indices describing
631 !! where the data
632 !! will be read from.
633 integer, dimension(:), intent(in), optional :: edge_lengths !< The number of
634 !! elements that
635 !! will be read
636 !! in each dimension.
637 logical, intent(in), optional :: broadcast !< Flag controlling whether or
638 !! not the data will be
639 !! broadcasted to non
640 !! "I/O root" ranks.
641 !! The broadcast will be done
642 !! by default.
643
644 logical :: bcast
645 integer :: err
646 integer :: varid
647 integer :: unlim_dim_index
648 integer, dimension(6) :: c
649 integer, dimension(6) :: e
650 character(len=200) :: append_error_msg !< Msg to be appended to FATAL error message
651
652 append_error_msg = "netcdf_read_data_5d: file:"//trim(fileobj%path)//"- variable:"//trim(variable_name)
653
654 if (present(broadcast)) then
655 bcast = broadcast
656 else
657 bcast = .true.
658 endif
659 c(:) = 1
660 if (present(corner)) then
661 c(1:5) = corner(1:5)
662 endif
663 e(:) = 1
664 if (present(edge_lengths)) then
665 e(1:5) = edge_lengths(1:5)
666 else
667 e(1:5) = shape(buf)
668 endif
669 if (present(unlim_dim_level)) then
670 unlim_dim_index = get_variable_unlimited_dimension_index(fileobj, variable_name, &
671 broadcast=bcast)
672 if (unlim_dim_index .ne. 6) then
673 call error("unlimited dimension must be the slowest varying dimension: "//trim(append_error_msg))
674 endif
675 c(unlim_dim_index) = unlim_dim_level
676 endif
677 if (fileobj%is_root) then
678 varid = get_variable_id(fileobj%ncid, trim(variable_name), msg=append_error_msg)
679 select type(buf)
680 type is (integer(kind=i4_kind))
681 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
682 type is (integer(kind=i8_kind))
683 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
684 type is (real(kind=r4_kind))
685 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
686 type is (real(kind=r8_kind))
687 err = nf90_get_var(fileobj%ncid, varid, buf, start=c, count=e)
688 class default
689 call error("Unsupported variable type: "//trim(append_error_msg))
690 end select
691 call check_netcdf_code(err, append_error_msg)
692 call unpack_data_5d(fileobj, varid, variable_name, buf)
693 endif
694 if (bcast) then
695 select type(buf)
696 type is (integer(kind=i4_kind))
697 call mpp_broadcast(buf, size(buf), fileobj%io_root, pelist=fileobj%pelist)
698 type is (integer(kind=i8_kind))
699 call mpp_broadcast(buf, size(buf), fileobj%io_root, pelist=fileobj%pelist)
700 type is (real(kind=r4_kind))
701 call mpp_broadcast(buf, size(buf), fileobj%io_root, pelist=fileobj%pelist)
702 type is (real(kind=r8_kind))
703 call mpp_broadcast(buf, size(buf), fileobj%io_root, pelist=fileobj%pelist)
704 class default
705 call error("Unsupported variable type: "//trim(append_error_msg))
706 end select
707 endif
708end subroutine netcdf_read_data_5d
709!> @}
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.