FMS  2025.04
Flexible Modeling System
netcdf_read_data.inc
1 !***********************************************************************
2 !* Apache License 2.0
3 !*
4 !* This file is part of the GFDL Flexible Modeling System (FMS).
5 !*
6 !* Licensed under the Apache License, Version 2.0 (the "License");
7 !* you may not use this file except in compliance with the License.
8 !* You may obtain a copy of the License at
9 !*
10 !* http://www.apache.org/licenses/LICENSE-2.0
11 !*
12 !* FMS is distributed in the hope that it will be useful, but WITHOUT
13 !* WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied;
14 !* without even the implied warranty of MERCHANTABILITY or FITNESS FOR A
15 !* PARTICULAR PURPOSE. See the License for the specific language
16 !* governing permissions and limitations under the License.
17 !***********************************************************************
18 !> @file
19 !> @brief Routines for retrieving variable values from netcdf files with different
20 !! dimension sizes for the @ref netcdf_read_data interface
21 
22 !> @addtogroup netcdf_io_mod
23 !> @{
24 
25 !> @brief Character read 0d function.
26 subroutine char_read_0d(fileobj, variable_name, buf, corner, append_error_msg, err, varid)
27  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
28  character(len=*), intent(in) :: variable_name !< Variable name.
29  character(len=*), intent(inout) :: buf !< Array that the data will be read into
30  integer, intent(in), optional :: corner !< Index of the string to read if the variable
31  !! contains a 1D array of strings.
32  character(len=200), intent(in):: append_error_msg !< Msg to be appended to FATAL error message
33  integer, intent(inout) :: err
34  integer, intent(in) :: varid
35 
36  integer, dimension(2) :: start
37  integer :: ndims
38  character, dimension(:), allocatable :: charbuf
39  integer, dimension(:), allocatable :: dimsizes
40  integer :: i
41 
42  start(:) = 1
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)))
47  charbuf(:) = ""
48  if (ndims .eq. 2) then
49  if (present(corner)) then
50  start(2) = corner
51  endif
52  dimsizes(2) = 1
53  elseif (ndims .gt. 2) then
54  call error("Only scalar and 1d string values are currently supported: "//trim(append_error_msg))
55  endif
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))
59  endif
60  buf = ""
61  do i = 1, dimsizes(1)
62  if (charbuf(i) .eq. char(0)) then
63  exit
64  endif
65  buf(i:i) = charbuf(i)
66  enddo
67  deallocate(charbuf)
68  deallocate(dimsizes)
69 end subroutine char_read_0d
70 
71 !> @brief Character read 1d function.
72 subroutine char_read_1d(fileobj, variable_name, buf, c, append_error_msg, err, varid)
73 
74  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
75  character(len=*), intent(in) :: variable_name !< Variable name.
76  character(len=*), dimension(:), intent(inout) :: buf !< Array that the data
77  !! will be read into.
78  integer, dimension(2), intent(in) :: c
79  character(len=200), intent(in) :: append_error_msg !< Msg to be appended to FATAL error message
80  integer, intent(inout) :: err
81  integer, intent(in) :: varid
82 
83  integer :: ndims
84  integer, dimension(2) :: start
85  integer, dimension(2) :: dimsizes
86  character, dimension(:,:), allocatable :: charbuf
87  character(len=1024) :: sbuf
88  integer :: i
89  integer :: j
90 
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.")
94  endif
95  start(1) = 1
96  start(2) = c(1)
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))
103  endif
104  if (size(buf) .lt. dimsizes(2)) then
105  call error("incorrect buffer array size:: "//trim(append_error_msg))
106  endif
107  do i = start(2), start(2)+dimsizes(2)-1
108  sbuf = ""
109  do j = 1, dimsizes(1)
110  if (charbuf(j, i-start(2)+1) .eq. char(0)) then
111  exit
112  endif
113  sbuf(j:j) = charbuf(j, i-start(2)+1)
114  enddo
115  call string_copy(buf(i-start(2)+1), sbuf)
116  enddo
117  deallocate(charbuf)
118 
119 end subroutine char_read_1d
120 
121 !> @brief Read in data from a variable in a netcdf file.
122 subroutine netcdf_read_data_0d(fileobj, variable_name, buf, unlim_dim_level, &
123  corner, broadcast)
124 
125  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
126  character(len=*), intent(in) :: variable_name !< Variable name.
127  class(*), intent(inout) :: buf !< Array that the data will be read into.
128  integer, intent(in), optional :: unlim_dim_level !< Unlimited dimension level.
129  integer, intent(in), optional :: corner !< Index of the string to read if the variable
130  !! contains a 1D array of strings.
131  logical, intent(in), optional :: broadcast !< Flag controlling whether or
132  !! not the data will be
133  !! broadcasted to non
134  !! "I/O root" ranks.
135  !! The broadcast will be done
136  !! by default.
137 
138  logical :: bcast
139  integer :: err
140  integer :: varid
141  integer :: unlim_dim_index
142  integer, dimension(1) :: c
143  character(len=1024), dimension(1) :: buf1d
144  character(len=200) :: append_error_msg !< Msg to be appended to FATAL error message
145 
146  append_error_msg = "netcdf_read_data_0d: file:"//trim(fileobj%path)//"- variable:"//trim(variable_name)
147 
148  if (present(broadcast)) then
149  bcast = broadcast
150  else
151  bcast = .true.
152  endif
153  c(1) = 1
154  if (present(unlim_dim_level)) then
155  unlim_dim_index = get_variable_unlimited_dimension_index(fileobj, variable_name, &
156  broadcast=bcast)
157  if (unlim_dim_index .ne. 1) then
158  call error("unlimited dimension must be the slowest varying dimension: "//trim(append_error_msg))
159  endif
160  c(unlim_dim_index) = unlim_dim_level
161  endif
162  if (fileobj%is_root) then
163  varid = get_variable_id(fileobj%ncid, trim(variable_name), msg=append_error_msg)
164  select type(buf)
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)
175  class default
176  call error("Unsupported variable type: "//trim(append_error_msg))
177  end select
178  call check_netcdf_code(err, append_error_msg)
179  call unpack_data_0d(fileobj, varid, variable_name, buf)
180  endif
181  if (bcast) then
182  select type(buf)
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))
195  class default
196  call error("Unsupported variable type: "//trim(append_error_msg))
197  end select
198  endif
199 end subroutine netcdf_read_data_0d
200 
201 
202 !> @brief Read in data from a variable in a netcdf file.
203 subroutine netcdf_read_data_1d(fileobj, variable_name, buf, unlim_dim_level, &
204  corner, edge_lengths, broadcast)
205 
206  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
207  character(len=*), intent(in) :: variable_name !< Variable name.
208  class(*), dimension(:), intent(inout) :: buf !< Array that the data
209  !! will be read into.
210  integer, intent(in), optional :: unlim_dim_level !< Unlimited dimension
211  !! level.
212  integer, dimension(:), intent(in), optional :: corner !< Array of starting
213  !! indices describing
214  !! where the data
215  !! will be read from.
216  integer, dimension(:), intent(in), optional :: edge_lengths !< The number of
217  !! elements that
218  !! will be read
219  !! in each dimension.
220  logical, intent(in), optional :: broadcast !< Flag controlling whether or
221  !! not the data will be
222  !! broadcasted to non
223  !! "I/O root" ranks.
224  !! The broadcast will be done
225  !! by default.
226 
227  logical :: bcast
228  integer :: err
229  integer :: varid
230  integer :: unlim_dim_index
231  integer, dimension(2) :: c
232  integer, dimension(2) :: e
233  character(len=200) :: append_error_msg !< Msg to be appended to FATAL error message
234 
235  append_error_msg = "netcdf_read_data_1d: file:"//trim(fileobj%path)//"- variable:"//trim(variable_name)
236 
237  if (present(broadcast)) then
238  bcast = broadcast
239  else
240  bcast = .true.
241  endif
242  c(:) = 1
243  if (present(corner)) then
244  c(1) = corner(1)
245  endif
246  e(:) = 1
247  if (present(edge_lengths)) then
248  e(1) = edge_lengths(1)
249  else
250  e(1) = size(buf)
251  endif
252  if (present(unlim_dim_level)) then
253  unlim_dim_index = get_variable_unlimited_dimension_index(fileobj, variable_name, &
254  broadcast=bcast)
255  if (unlim_dim_index .ne. 2) then
256  call error("unlimited dimension must be the slowest varying dimension: "//trim(append_error_msg))
257  endif
258  c(unlim_dim_index) = unlim_dim_level
259  endif
260  if (fileobj%is_root) then
261  varid = get_variable_id(fileobj%ncid, trim(variable_name), msg=append_error_msg)
262  select type(buf)
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)
273  class default
274  call error("Unsupported variable type: "//trim(append_error_msg))
275  end select
276  call check_netcdf_code(err, append_error_msg)
277  call unpack_data_1d(fileobj, varid, variable_name, buf)
278  endif
279  if (bcast) then
280  select type(buf)
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)
291  class default
292  call error("Unsupported variable type: "//trim(append_error_msg))
293  end select
294  endif
295 end subroutine netcdf_read_data_1d
296 
297 
298 !> @brief Read in data from a variable in a netcdf file.
299 subroutine netcdf_read_data_2d(fileobj, variable_name, buf, unlim_dim_level, &
300  corner, edge_lengths, broadcast)
301 
302  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
303  character(len=*), intent(in) :: variable_name !< Variable name.
304  class(*), dimension(:,:), intent(inout) :: buf !< Array that the data
305  !! will be read into.
306  integer, intent(in), optional :: unlim_dim_level !< Unlimited dimension
307  !! level.
308  integer, dimension(:), intent(in), optional :: corner !< Array of starting
309  !! indices describing
310  !! where the data
311  !! will be read from.
312  integer, dimension(:), intent(in), optional :: edge_lengths !< The number of
313  !! elements that
314  !! will be read
315  !! in each dimension.
316  logical, intent(in), optional :: broadcast !< Flag controlling whether or
317  !! not the data will be
318  !! broadcasted to non
319  !! "I/O root" ranks.
320  !! The broadcast will be done
321  !! by default.
322 
323  logical :: bcast
324  integer :: err
325  integer :: varid
326  integer :: unlim_dim_index
327  integer, dimension(3) :: c
328  integer, dimension(3) :: e
329  character(len=200) :: append_error_msg !< Msg to be appended to FATAL error message
330 
331  append_error_msg = "netcdf_read_data_2d: file:"//trim(fileobj%path)//"- variable:"//trim(variable_name)
332 
333  if (present(broadcast)) then
334  bcast = broadcast
335  else
336  bcast = .true.
337  endif
338  c(:) = 1
339  if (present(corner)) then
340  c(1:2) = corner(1:2)
341  endif
342  e(:) = 1
343  if (present(edge_lengths)) then
344  e(1:2) = edge_lengths(1:2)
345  else
346  e(1:2) = shape(buf)
347  endif
348  if (present(unlim_dim_level)) then
349  unlim_dim_index = get_variable_unlimited_dimension_index(fileobj, variable_name, &
350  broadcast=bcast)
351  if (unlim_dim_index .ne. 3) then
352  call error("unlimited dimension must be the slowest varying dimension: "//trim(append_error_msg))
353  endif
354  c(unlim_dim_index) = unlim_dim_level
355  endif
356  if(fileobj%use_collective) then
357  varid = get_variable_id(fileobj%ncid, trim(variable_name), msg=append_error_msg)
358  ! NetCDF does not have the ability to specify collective I/O at
359  ! the file basis so we must activate at the variable level
360  err = nf90_var_par_access(fileobj%ncid, varid, nf90_collective)
361  call check_netcdf_code(err, append_error_msg)
362  select type(buf)
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)
371  class default
372  call error("Unsupported variable type: "//trim(append_error_msg))
373  end select
374  call check_netcdf_code(err, append_error_msg)
375  call unpack_data_2d(fileobj, varid, variable_name, buf)
376  else
377  if (fileobj%is_root) then
378  varid = get_variable_id(fileobj%ncid, trim(variable_name), msg=append_error_msg)
379  select type(buf)
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)
388  class default
389  call error("Unsupported variable type: "//trim(append_error_msg))
390  end select
391  call check_netcdf_code(err, append_error_msg)
392  call unpack_data_2d(fileobj, varid, variable_name, buf)
393  endif
394  if (bcast) then
395  select type(buf)
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)
404  class default
405  call error("Unsupported variable type: "//trim(append_error_msg))
406  end select
407  endif
408  endif
409 end subroutine netcdf_read_data_2d
410 
411 
412 !> @brief Read in data from a variable in a netcdf file.
413 subroutine netcdf_read_data_3d(fileobj, variable_name, buf, unlim_dim_level, &
414  corner, edge_lengths, broadcast)
415 
416  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
417  character(len=*), intent(in) :: variable_name !< Variable name.
418  class(*), dimension(:,:,:), intent(inout) :: buf !< Array that the data
419  !! will be read into.
420  integer, intent(in), optional :: unlim_dim_level !< Unlimited dimension
421  !! level.
422  integer, dimension(:), intent(in), optional :: corner !< Array of starting
423  !! indices describing
424  !! where the data
425  !! will be read from.
426  integer, dimension(:), intent(in), optional :: edge_lengths !< The number of
427  !! elements that
428  !! will be read
429  !! in each dimension.
430  logical, intent(in), optional :: broadcast !< Flag controlling whether or
431  !! not the data will be
432  !! broadcasted to non
433  !! "I/O root" ranks.
434  !! The broadcast will be done
435  !! by default.
436 
437  logical :: bcast
438  integer :: err
439  integer :: varid
440  integer :: unlim_dim_index
441  integer, dimension(4) :: c
442  integer, dimension(4) :: e
443  character(len=200) :: append_error_msg !< Msg to be appended to FATAL error message
444 
445  append_error_msg = "netcdf_read_data_3d: file:"//trim(fileobj%path)//"- variable:"//trim(variable_name)
446 
447  if (present(broadcast)) then
448  bcast = broadcast
449  else
450  bcast = .true.
451  endif
452  c(:) = 1
453  if (present(corner)) then
454  c(1:3) = corner(1:3)
455  endif
456  e(:) = 1
457  if (present(edge_lengths)) then
458  e(1:3) = edge_lengths(1:3)
459  else
460  e(1:3) = shape(buf)
461  endif
462  if (present(unlim_dim_level)) then
463  unlim_dim_index = get_variable_unlimited_dimension_index(fileobj, variable_name, &
464  broadcast=bcast)
465  if (unlim_dim_index .ne. 4) then
466  call error("unlimited dimension must be the slowest varying dimension: "//trim(append_error_msg))
467  endif
468  c(unlim_dim_index) = unlim_dim_level
469  endif
470  if(fileobj%use_collective) then
471  varid = get_variable_id(fileobj%ncid, trim(variable_name), msg=append_error_msg)
472  ! NetCDF does not have the ability to specify collective I/O at
473  ! the file basis so we must activate at the variable level
474  err = nf90_var_par_access(fileobj%ncid, varid, nf90_collective)
475  call check_netcdf_code(err, append_error_msg)
476  select type(buf)
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)
485  class default
486  call error("Unsupported variable type: "//trim(append_error_msg))
487  end select
488  call check_netcdf_code(err, append_error_msg)
489  call unpack_data_3d(fileobj, varid, variable_name, buf)
490  else
491  if (fileobj%is_root) then
492  varid = get_variable_id(fileobj%ncid, trim(variable_name), msg=append_error_msg)
493  select type(buf)
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)
502  class default
503  call error("Unsupported variable type: "//trim(append_error_msg))
504  end select
505  call check_netcdf_code(err, append_error_msg)
506  call unpack_data_3d(fileobj, varid, variable_name, buf)
507  endif
508  if (bcast) then
509  select type(buf)
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)
518  class default
519  call error("Unsupported variable type: "//trim(append_error_msg))
520  end select
521  endif
522  endif
523 end subroutine netcdf_read_data_3d
524 
525 
526 !> @brief Read in data from a variable in a netcdf file.
527 subroutine netcdf_read_data_4d(fileobj, variable_name, buf, unlim_dim_level, &
528  corner, edge_lengths, broadcast)
529 
530  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
531  character(len=*), intent(in) :: variable_name !< Variable name.
532  class(*), dimension(:,:,:,:), intent(inout) :: buf !< Array that the data
533  !! will be read into.
534  integer, intent(in), optional :: unlim_dim_level !< Unlimited dimension
535  !! level.
536  integer, dimension(:), intent(in), optional :: corner !< Array of starting
537  !! indices describing
538  !! where the data
539  !! will be read from.
540  integer, dimension(:), intent(in), optional :: edge_lengths !< The number of
541  !! elements that
542  !! will be read
543  !! in each dimension.
544  logical, intent(in), optional :: broadcast !< Flag controlling whether or
545  !! not the data will be
546  !! broadcasted to non
547  !! "I/O root" ranks.
548  !! The broadcast will be done
549  !! by default.
550 
551  logical :: bcast
552  integer :: err
553  integer :: varid
554  integer :: unlim_dim_index
555  integer, dimension(5) :: c
556  integer, dimension(5) :: e
557  character(len=200) :: append_error_msg !< Msg to be appended to FATAL error message
558 
559  append_error_msg = "netcdf_read_data_4d: file:"//trim(fileobj%path)//"- variable:"//trim(variable_name)
560 
561  if (present(broadcast)) then
562  bcast = broadcast
563  else
564  bcast = .true.
565  endif
566  c(:) = 1
567  if (present(corner)) then
568  c(1:4) = corner(1:4)
569  endif
570  e(:) = 1
571  if (present(edge_lengths)) then
572  e(1:4) = edge_lengths(1:4)
573  else
574  e(1:4) = shape(buf)
575  endif
576  if (present(unlim_dim_level)) then
577  unlim_dim_index = get_variable_unlimited_dimension_index(fileobj, variable_name, &
578  broadcast=bcast)
579  if (unlim_dim_index .ne. 5) then
580  call error("unlimited dimension must be the slowest varying dimension: "//trim(append_error_msg))
581  endif
582  c(unlim_dim_index) = unlim_dim_level
583  endif
584  if (fileobj%is_root) then
585  varid = get_variable_id(fileobj%ncid, trim(variable_name), msg=append_error_msg)
586  select type(buf)
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)
595  class default
596  call error("Unsupported variable type: "//trim(append_error_msg))
597  end select
598  call check_netcdf_code(err, append_error_msg)
599  call unpack_data_4d(fileobj, varid, variable_name, buf)
600  endif
601  if (bcast) then
602  select type(buf)
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)
611  class default
612  call error("Unsupported variable type: "//trim(append_error_msg))
613  end select
614  endif
615 end subroutine netcdf_read_data_4d
616 
617 
618 !> @brief Read in data from a variable in a netcdf file.
619 subroutine netcdf_read_data_5d(fileobj, variable_name, buf, unlim_dim_level, &
620  corner, edge_lengths, broadcast)
621 
622  class(fmsnetcdffile_t), intent(in) :: fileobj !< File object.
623  character(len=*), intent(in) :: variable_name !< Variable name.
624  class(*), dimension(:,:,:,:,:), intent(inout) :: buf !< Array that the data
625  !! will be read into.
626  integer, intent(in), optional :: unlim_dim_level !< Unlimited dimension
627  !! level.
628  integer, dimension(:), intent(in), optional :: corner !< Array of starting
629  !! indices describing
630  !! where the data
631  !! will be read from.
632  integer, dimension(:), intent(in), optional :: edge_lengths !< The number of
633  !! elements that
634  !! will be read
635  !! in each dimension.
636  logical, intent(in), optional :: broadcast !< Flag controlling whether or
637  !! not the data will be
638  !! broadcasted to non
639  !! "I/O root" ranks.
640  !! The broadcast will be done
641  !! by default.
642 
643  logical :: bcast
644  integer :: err
645  integer :: varid
646  integer :: unlim_dim_index
647  integer, dimension(6) :: c
648  integer, dimension(6) :: e
649  character(len=200) :: append_error_msg !< Msg to be appended to FATAL error message
650 
651  append_error_msg = "netcdf_read_data_5d: file:"//trim(fileobj%path)//"- variable:"//trim(variable_name)
652 
653  if (present(broadcast)) then
654  bcast = broadcast
655  else
656  bcast = .true.
657  endif
658  c(:) = 1
659  if (present(corner)) then
660  c(1:5) = corner(1:5)
661  endif
662  e(:) = 1
663  if (present(edge_lengths)) then
664  e(1:5) = edge_lengths(1:5)
665  else
666  e(1:5) = shape(buf)
667  endif
668  if (present(unlim_dim_level)) then
669  unlim_dim_index = get_variable_unlimited_dimension_index(fileobj, variable_name, &
670  broadcast=bcast)
671  if (unlim_dim_index .ne. 6) then
672  call error("unlimited dimension must be the slowest varying dimension: "//trim(append_error_msg))
673  endif
674  c(unlim_dim_index) = unlim_dim_level
675  endif
676  if (fileobj%is_root) then
677  varid = get_variable_id(fileobj%ncid, trim(variable_name), msg=append_error_msg)
678  select type(buf)
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)
687  class default
688  call error("Unsupported variable type: "//trim(append_error_msg))
689  end select
690  call check_netcdf_code(err, append_error_msg)
691  call unpack_data_5d(fileobj, varid, variable_name, buf)
692  endif
693  if (bcast) then
694  select type(buf)
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)
703  class default
704  call error("Unsupported variable type: "//trim(append_error_msg))
705  end select
706  endif
707 end subroutine netcdf_read_data_5d
708 !> @}
subroutine unpack_data_2d(fileobj, varid, varname, var_data)
subroutine unpack_data_1d(fileobj, varid, varname, var_data)
Definition: unpack_data.inc:68
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)
Definition: unpack_data.inc:24
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.