MODFLOW 6  version 6.7.0.dev3
USGS Modular Hydrologic Model
NCArrayReader.f90
Go to the documentation of this file.
1 !> @brief This module contains the NCArrayReaderModule
2 !!
3 !! This module defines the netcdf_array_load interface
4 !! which can read layered (UGRID) and non-layered (STRUCTURED)
5 !! netcdf arrays stored in modflow6 designated input variables.
6 !!
7 !<
9 
10  use kindmodule, only: dp, i4b, lgp
11  use constantsmodule, only: linelength
12  use simvariablesmodule, only: errmsg
18  use netcdfcommonmodule, only: nf_verify
19  use netcdf
20 
21  implicit none
22  private
23  public :: netcdf_array_load
24 
26  module procedure nc_array_load_int1d, nc_array_load_int2d, &
29  end interface netcdf_array_load
30 
31 contains
32 
33  !> @brief does the grid support per layer variables
34  !<
35  function is_layered(grid) result(layered)
36  character(len=*), intent(in) :: grid
37  logical(LGP) :: layered
38  select case (grid)
39  case ('LAYERED MESH')
40  layered = .true.
41  case ('STRUCTURED')
42  layered = .false.
43  case default
44  layered = .false.
45  end select
46  end function is_layered
47 
48  !> @brief Load NetCDF integer 1D array
49  !<
50  subroutine nc_array_load_int1d(int1d, mshape, idt, mf6_input, nc_vars, &
51  input_fname, iout, kper)
52  integer(I4B), dimension(:), pointer, contiguous, intent(in) :: int1d
53  integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape !< model shape
54  type(inputparamdefinitiontype), intent(in) :: idt !< input data type object describing this record
55  type(modflowinputtype), intent(in) :: mf6_input
56  type(ncpackagevarstype), pointer, intent(in) :: nc_vars
57  character(len=*), intent(in) :: input_fname
58  integer(I4B), intent(in) :: iout
59  integer(I4B), optional, intent(in) :: kper !< flag if set > 0 indicates ts
60  integer(I4B) :: varid, iper
61  logical(LGP) :: layered
62 
63  iper = 0
64  layered = (idt%layered .and. is_layered(nc_vars%grid))
65 
66  if (present(kper)) then
67  iper = kper
68  end if
69 
70  if (layered) then
71  if (iper > 0) then
72  call load_integer1d_layered_spd(int1d, mf6_input, mshape, idt, nc_vars, &
73  iper, input_fname)
74  else
75  call load_integer1d_layered(int1d, mf6_input, mshape, idt, nc_vars, &
76  input_fname)
77  end if
78  else
79  if (iper > 0) then
80  call load_integer1d_spd(int1d, mf6_input, mshape, idt, nc_vars, &
81  iper, input_fname)
82  else
83  varid = nc_vars%varid(idt%tagname)
84  call load_integer1d_type(int1d, mf6_input, mshape, idt, nc_vars, &
85  varid, input_fname)
86  end if
87  end if
88  end subroutine nc_array_load_int1d
89 
90  !> @brief Load NetCDF integer 2D array
91  !<
92  subroutine nc_array_load_int2d(int2d, mshape, idt, mf6_input, nc_vars, &
93  input_fname, iout)
94  integer(I4B), dimension(:, :), pointer, contiguous, intent(in) :: int2d
95  integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape !< model shape
96  type(inputparamdefinitiontype), intent(in) :: idt !< input data type object describing this record
97  type(modflowinputtype), intent(in) :: mf6_input
98  type(ncpackagevarstype), pointer, intent(in) :: nc_vars
99  character(len=*), intent(in) :: input_fname
100  integer(I4B), intent(in) :: iout
101  integer(I4B) :: varid
102  logical(LGP) :: layered
103 
104  layered = (idt%layered .and. is_layered(nc_vars%grid))
105 
106  if (layered) then
107  call load_integer2d_layered(int2d, mf6_input, mshape, idt, nc_vars, &
108  input_fname)
109  else
110  varid = nc_vars%varid(idt%tagname)
111  call load_integer2d_type(int2d, mf6_input, mshape, idt, nc_vars, &
112  varid, input_fname)
113  end if
114  end subroutine nc_array_load_int2d
115 
116  !> @brief Load NetCDF integer 3D array
117  !<
118  subroutine nc_array_load_int3d(int3d, mshape, idt, mf6_input, nc_vars, &
119  input_fname, iout)
120  integer(I4B), dimension(:, :, :), pointer, contiguous, intent(in) :: int3d
121  integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape !< model shape
122  type(inputparamdefinitiontype), intent(in) :: idt !< input data type object describing this record
123  type(modflowinputtype), intent(in) :: mf6_input
124  type(ncpackagevarstype), pointer, intent(in) :: nc_vars
125  character(len=*), intent(in) :: input_fname
126  integer(I4B), intent(in) :: iout
127  integer(I4B) :: varid
128  logical(LGP) :: layered
129 
130  layered = (idt%layered .and. is_layered(nc_vars%grid))
131 
132  if (layered) then
133  call load_integer3d_layered(int3d, mf6_input, mshape, idt, nc_vars, &
134  input_fname)
135  else
136  varid = nc_vars%varid(idt%tagname)
137  call load_integer3d_type(int3d, mf6_input, mshape, idt, nc_vars, &
138  varid, input_fname)
139  end if
140  end subroutine nc_array_load_int3d
141 
142  !> @brief Load NetCDF double 1D array
143  !<
144  subroutine nc_array_load_dbl1d(dbl1d, mshape, idt, mf6_input, nc_vars, &
145  input_fname, iout, kper, iaux)
146  real(DP), dimension(:), pointer, contiguous, intent(in) :: dbl1d
147  integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape !< model shape
148  type(inputparamdefinitiontype), intent(in) :: idt !< input data type object describing this record
149  type(modflowinputtype), intent(in) :: mf6_input
150  type(ncpackagevarstype), pointer, intent(in) :: nc_vars
151  character(len=*), intent(in) :: input_fname
152  integer(I4B), intent(in) :: iout
153  integer(I4B), optional, intent(in) :: kper !< flag if set > 0 indicates ts
154  integer(I4B), optional, intent(in) :: iaux
155  integer(I4B) :: varid, iper
156  logical(LGP) :: layered
157 
158  iper = 0
159  layered = (idt%layered .and. is_layered(nc_vars%grid))
160 
161  if (present(kper)) then
162  iper = kper
163  end if
164 
165  if (layered) then
166  if (iper > 0) then
167  call load_double1d_layered_spd(dbl1d, mf6_input, mshape, idt, nc_vars, &
168  iper, input_fname, iaux)
169  else
170  call load_double1d_layered(dbl1d, mf6_input, mshape, idt, nc_vars, &
171  input_fname)
172  end if
173  else
174  if (iper > 0) then
175  call load_double1d_spd(dbl1d, mf6_input, mshape, idt, nc_vars, &
176  iper, input_fname, iaux)
177  else
178  varid = nc_vars%varid(idt%tagname)
179  call load_double1d_type(dbl1d, mf6_input, mshape, idt, nc_vars, &
180  varid, input_fname)
181  end if
182  end if
183  end subroutine nc_array_load_dbl1d
184 
185  !> @brief Load NetCDF double 2D array
186  !<
187  subroutine nc_array_load_dbl2d(dbl2d, mshape, idt, mf6_input, nc_vars, &
188  input_fname, iout)
189  real(DP), dimension(:, :), pointer, contiguous, intent(in) :: dbl2d
190  integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape !< model shape
191  type(inputparamdefinitiontype), intent(in) :: idt !< input data type object describing this record
192  type(modflowinputtype), intent(in) :: mf6_input
193  type(ncpackagevarstype), pointer, intent(in) :: nc_vars
194  character(len=*), intent(in) :: input_fname
195  integer(I4B), intent(in) :: iout
196  integer(I4B) :: varid
197  logical(LGP) :: layered
198 
199  layered = (idt%layered .and. is_layered(nc_vars%grid))
200 
201  if (layered) then
202  call load_double2d_layered(dbl2d, mf6_input, mshape, idt, nc_vars, &
203  input_fname)
204  else
205  varid = nc_vars%varid(idt%tagname)
206  call load_double2d_type(dbl2d, mf6_input, mshape, idt, nc_vars, &
207  varid, input_fname)
208  end if
209  end subroutine nc_array_load_dbl2d
210 
211  !> @brief Load NetCDF double 3D array
212  !<
213  subroutine nc_array_load_dbl3d(dbl3d, mshape, idt, mf6_input, nc_vars, &
214  input_fname, iout)
215  real(DP), dimension(:, :, :), pointer, contiguous, intent(in) :: dbl3d
216  integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape !< model shape
217  type(inputparamdefinitiontype), intent(in) :: idt !< input data type object describing this record
218  type(modflowinputtype), intent(in) :: mf6_input
219  type(ncpackagevarstype), pointer, intent(in) :: nc_vars
220  character(len=*), intent(in) :: input_fname
221  integer(I4B), intent(in) :: iout
222  integer(I4B) :: varid
223  logical(LGP) :: layered
224 
225  layered = (idt%layered .and. is_layered(nc_vars%grid))
226 
227  if (layered) then
228  call load_double3d_layered(dbl3d, mf6_input, mshape, idt, nc_vars, &
229  input_fname)
230  else
231  varid = nc_vars%varid(idt%tagname)
232  call load_double3d_type(dbl3d, mf6_input, mshape, idt, nc_vars, &
233  varid, input_fname)
234  end if
235  end subroutine nc_array_load_dbl3d
236 
237  !> @brief load type 1d integer
238  !<
239  subroutine load_integer1d_type(int1d, mf6_input, mshape, idt, nc_vars, &
240  varid, input_fname)
241  integer(I4B), dimension(:), contiguous, pointer, intent(in) :: int1d
242  type(modflowinputtype), intent(in) :: mf6_input
243  integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape
244  type(inputparamdefinitiontype), intent(in) :: idt
245  type(ncpackagevarstype), pointer, intent(in) :: nc_vars
246  integer(I4B), intent(in) :: varid
247  character(len=*), intent(in) :: input_fname
248  integer(I4B), dimension(:), allocatable :: array_shape
249  integer(I4B), dimension(:, :, :), contiguous, pointer :: int3d_ptr
250  integer(I4B), dimension(:, :), contiguous, pointer :: int2d_ptr
251  integer(I4B) :: nvals
252 
253  ! initialize
254  nvals = 0
255 
256  if (idt%shape == 'NODES') then
257  ! set number of values
258  nvals = product(mshape)
259  if (size(mshape) == 3) then
260  int3d_ptr(1:mshape(3), 1:mshape(2), 1:mshape(1)) => int1d(1:nvals)
261  call nf_verify(nf90_get_var(nc_vars%ncid, varid, int3d_ptr), &
262  nc_vars%nc_fname)
263  else if (size(mshape) == 2) then
264  int2d_ptr(1:mshape(2), 1:mshape(1)) => int1d(1:nvals)
265  call nf_verify(nf90_get_var(nc_vars%ncid, varid, int2d_ptr), &
266  nc_vars%nc_fname)
267  else if (size(mshape) == 1) then
268  call nf_verify(nf90_get_var(nc_vars%ncid, varid, int1d), nc_vars%nc_fname)
269  end if
270  else
271  ! interpret shape
272  call get_shape_from_string(idt%shape, array_shape, mf6_input%mempath)
273  ! set nvals
274  nvals = array_shape(1)
275  ! read and set data
276  call nf_verify(nf90_get_var(nc_vars%ncid, varid, int1d), nc_vars%nc_fname)
277  end if
278  end subroutine load_integer1d_type
279 
280  !> @brief load type 1d double
281  !<
282  subroutine load_integer1d_spd(int1d, mf6_input, mshape, idt, nc_vars, &
283  iper, input_fname)
284  use constantsmodule, only: dnodata
285  use netcdfcommonmodule, only: ixstp
286  integer(I4B), dimension(:), contiguous, pointer, intent(in) :: int1d
287  type(modflowinputtype), intent(in) :: mf6_input
288  integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape
289  type(inputparamdefinitiontype), intent(in) :: idt
290  type(ncpackagevarstype), pointer, intent(in) :: nc_vars
291  integer(I4B), intent(in) :: iper
292  character(len=*), intent(in) :: input_fname
293  integer(I4B), dimension(:), allocatable :: layer_shape
294  integer(I4B) :: varid, nlay, ncpl, istp
295 
296  istp = ixstp()
297 
298  ! set varid
299  varid = nc_vars%varid(idt%tagname)
300 
301  call get_layered_shape(mshape, nlay, layer_shape)
302  ncpl = product(layer_shape)
303 
304  select case (idt%shape)
305  case ('NCPL', 'NAUX NCPL')
306  if (nc_vars%grid == 'STRUCTURED') then
307  call nf_verify(nf90_get_var(nc_vars%ncid, varid, int1d, &
308  start=(/1, 1, istp/), &
309  count=(/mshape(3), mshape(2), 1/)), &
310  nc_vars%nc_fname)
311  else if (nc_vars%grid == 'LAYERED MESH') then
312  call nf_verify(nf90_get_var(nc_vars%ncid, varid, int1d, &
313  start=(/1, istp/), count=(/ncpl, 1/)), &
314  nc_vars%nc_fname)
315  end if
316  case ('NODES', 'NAUX NODES')
317  write (errmsg, '(a,a,a)') &
318  'Timeseries netcdf input read not supported for full grid int1d &
319  &type ('//trim(idt%tagname)//').'
320  call store_error(errmsg)
321  call store_error_filename(input_fname)
322  case default
323  end select
324  end subroutine load_integer1d_spd
325 
326  !> @brief load type 1d integer layered
327  !<
328  subroutine load_integer1d_layered(int1d, mf6_input, mshape, idt, nc_vars, &
329  input_fname)
330  integer(I4B), dimension(:), contiguous, pointer, intent(in) :: int1d
331  type(modflowinputtype), intent(in) :: mf6_input
332  integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape
333  type(inputparamdefinitiontype), intent(in) :: idt
334  type(ncpackagevarstype), pointer, intent(in) :: nc_vars
335  character(len=*), intent(in) :: input_fname
336  integer(I4B), dimension(:), allocatable :: layer_shape
337  integer(I4B) :: nlay, varid
338  integer(I4B) :: k, ncpl
339  integer(I4B) :: index_start, index_stop
340  integer(I4B), dimension(:), contiguous, pointer :: int1d_ptr
341 
342  nullify (int1d_ptr)
343 
344  call get_layered_shape(mshape, nlay, layer_shape)
345 
346  ncpl = product(layer_shape)
347  index_start = 1
348  do k = 1, nlay
349  varid = nc_vars%varid(idt%tagname, layer=k)
350  index_stop = index_start + ncpl - 1
351  int1d_ptr(1:ncpl) => int1d(index_start:index_stop)
352  call nf_verify(nf90_get_var(nc_vars%ncid, varid, int1d_ptr), &
353  nc_vars%nc_fname)
354  index_start = index_stop + 1
355  end do
356  end subroutine load_integer1d_layered
357 
358  !> @brief load type 1d integer layered
359  !<
360  subroutine load_integer1d_layered_spd(int1d, mf6_input, mshape, idt, nc_vars, &
361  iper, input_fname)
362  use constantsmodule, only: dnodata
363  use netcdfcommonmodule, only: ixstp
364  integer(I4B), dimension(:), contiguous, pointer, intent(in) :: int1d
365  type(modflowinputtype), intent(in) :: mf6_input
366  integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape
367  type(inputparamdefinitiontype), intent(in) :: idt
368  type(ncpackagevarstype), pointer, intent(in) :: nc_vars
369  integer(I4B), intent(in) :: iper
370  character(len=*), intent(in) :: input_fname
371  integer(I4B), dimension(:), allocatable :: layer_shape
372  integer(I4B) :: nlay, varid
373  integer(I4B) :: ncpl, nvals, istp
374 
375  istp = ixstp()
376 
377  call get_layered_shape(mshape, nlay, layer_shape)
378  nvals = product(mshape)
379  ncpl = product(layer_shape)
380 
381  varid = nc_vars%varid(idt%tagname)
382  select case (idt%shape)
383  case ('NCPL', 'NAUX NCPL')
384  call nf_verify(nf90_get_var(nc_vars%ncid, varid, int1d, &
385  start=(/1, istp/), count=(/ncpl, 1/)), &
386  nc_vars%nc_fname)
387  case ('NODES', 'NAUX NODES')
388  call nf_verify(nf90_get_var(nc_vars%ncid, varid, int1d, &
389  start=(/1, istp/), count=(/nvals, 1/)), &
390  nc_vars%nc_fname)
391  case default
392  end select
393  end subroutine load_integer1d_layered_spd
394 
395  !> @brief load type 2d integer
396  !<
397  subroutine load_integer2d_type(int2d, mf6_input, mshape, idt, nc_vars, varid, &
398  input_fname)
399  integer(I4B), dimension(:, :), contiguous, pointer, intent(in) :: int2d
400  type(modflowinputtype), intent(in) :: mf6_input
401  integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape
402  type(inputparamdefinitiontype), intent(in) :: idt
403  type(ncpackagevarstype), pointer, intent(in) :: nc_vars
404  integer(I4B), intent(in) :: varid
405  character(len=*), intent(in) :: input_fname
406  integer(I4B), dimension(:), contiguous, pointer :: int1d_ptr
407  integer(I4B), dimension(:), allocatable :: array_shape
408  integer(I4B) :: ncpl, nlay
409 
410  nullify (int1d_ptr)
411 
412  if (nc_vars%grid == 'STRUCTURED') then
413  call nf_verify(nf90_get_var(nc_vars%ncid, varid, int2d), nc_vars%nc_fname)
414  else if (nc_vars%grid == 'LAYERED MESH') then
415  call get_layered_shape(mshape, nlay, array_shape)
416  ncpl = product(array_shape)
417  int1d_ptr(1:ncpl) => int2d(:, :)
418  call nf_verify(nf90_get_var(nc_vars%ncid, varid, int1d_ptr), &
419  nc_vars%nc_fname)
420  end if
421  end subroutine load_integer2d_type
422 
423  !> @brief load type 2d integer layered
424  !<
425  subroutine load_integer2d_layered(int2d, mf6_input, mshape, idt, nc_vars, &
426  input_fname)
427  integer(I4B), dimension(:, :), contiguous, pointer, intent(in) :: int2d
428  type(modflowinputtype), intent(in) :: mf6_input
429  integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape
430  type(inputparamdefinitiontype), intent(in) :: idt
431  type(ncpackagevarstype), pointer, intent(in) :: nc_vars
432  character(len=*), intent(in) :: input_fname
433  integer(I4B), dimension(:), allocatable :: layer_shape
434  integer(I4B) :: k
435  integer(I4B) :: ncpl, nlay, varid
436  integer(I4B), dimension(:), contiguous, pointer :: int1d_ptr
437 
438  nullify (int1d_ptr)
439 
440  if (size(mshape) == 3) then
441  write (errmsg, '(a,a,a)') &
442  'Layered netcdf read not supported for DIS int2d type ('// &
443  trim(idt%tagname)//').'
444  call store_error(errmsg)
445  call store_error_filename(input_fname)
446  else if (size(mshape) == 2) then
447  call get_layered_shape(mshape, nlay, layer_shape)
448  ncpl = layer_shape(1)
449  do k = 1, nlay
450  varid = nc_vars%varid(idt%tagname, layer=k)
451  int1d_ptr(1:ncpl) => int2d(1:ncpl, k)
452  call nf_verify(nf90_get_var(nc_vars%ncid, varid, int1d_ptr), &
453  nc_vars%nc_fname)
454  end do
455  end if
456  end subroutine load_integer2d_layered
457 
458  !> @brief load type 3d integer
459  !<
460  subroutine load_integer3d_type(int3d, mf6_input, mshape, idt, nc_vars, varid, &
461  input_fname)
462  integer(I4B), dimension(:, :, :), contiguous, pointer, intent(in) :: int3d
463  type(modflowinputtype), intent(in) :: mf6_input
464  integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape
465  type(inputparamdefinitiontype), intent(in) :: idt
466  type(ncpackagevarstype), pointer, intent(in) :: nc_vars
467  integer(I4B), intent(in) :: varid
468  character(len=*), intent(in) :: input_fname
469  call nf_verify(nf90_get_var(nc_vars%ncid, varid, int3d), nc_vars%nc_fname)
470  end subroutine load_integer3d_type
471 
472  !> @brief load type 3d integer layered
473  !<
474  subroutine load_integer3d_layered(int3d, mf6_input, mshape, idt, nc_vars, &
475  input_fname)
476  integer(I4B), dimension(:, :, :), contiguous, pointer, intent(in) :: int3d
477  type(modflowinputtype), intent(in) :: mf6_input
478  integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape
479  type(inputparamdefinitiontype), intent(in) :: idt
480  type(ncpackagevarstype), pointer, intent(in) :: nc_vars
481  character(len=*), intent(in) :: input_fname
482  integer(I4B), dimension(:), allocatable :: layer_shape
483  integer(I4B) :: k !, i, j
484  integer(I4B) :: ncpl, nlay, varid
485  integer(I4B) :: index_start, index_stop
486  integer(I4B), dimension(:), contiguous, pointer :: int1d_ptr
487 
488  nullify (int1d_ptr)
489  index_start = 1
490  call get_layered_shape(mshape, nlay, layer_shape)
491  ncpl = product(layer_shape)
492 
493  do k = 1, nlay
494  varid = nc_vars%varid(idt%tagname, layer=k)
495  index_stop = index_start + ncpl - 1
496  int1d_ptr(1:ncpl) => int3d(:, :, k:k)
497  call nf_verify(nf90_get_var(nc_vars%ncid, varid, int1d_ptr), &
498  nc_vars%nc_fname)
499  index_start = index_stop + 1
500  end do
501  end subroutine load_integer3d_layered
502 
503  !> @brief load type 1d double
504  !<
505  subroutine load_double1d_type(dbl1d, mf6_input, mshape, idt, nc_vars, &
506  varid, input_fname)
507  real(DP), dimension(:), contiguous, pointer, intent(in) :: dbl1d
508  type(modflowinputtype), intent(in) :: mf6_input
509  integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape
510  type(inputparamdefinitiontype), intent(in) :: idt
511  type(ncpackagevarstype), pointer, intent(in) :: nc_vars
512  integer(I4B), intent(in) :: varid
513  character(len=*), intent(in) :: input_fname
514  integer(I4B), dimension(:), allocatable :: array_shape
515  real(DP), dimension(:, :, :), contiguous, pointer :: dbl3d_ptr
516  real(DP), dimension(:, :), contiguous, pointer :: dbl2d_ptr
517  integer(I4B) :: nvals
518 
519  ! initialize
520  nvals = 0
521 
522  if (idt%shape == 'NODES') then
523  ! set number of values
524  nvals = product(mshape)
525  if (size(mshape) == 3) then
526  dbl3d_ptr(1:mshape(3), 1:mshape(2), 1:mshape(1)) => dbl1d(1:nvals)
527  call nf_verify(nf90_get_var(nc_vars%ncid, varid, dbl3d_ptr), &
528  nc_vars%nc_fname)
529  else if (size(mshape) == 2) then
530  dbl2d_ptr(1:mshape(2), 1:mshape(1)) => dbl1d(1:nvals)
531  call nf_verify(nf90_get_var(nc_vars%ncid, varid, dbl2d_ptr), &
532  nc_vars%nc_fname)
533  else if (size(mshape) == 1) then
534  call nf_verify(nf90_get_var(nc_vars%ncid, varid, dbl1d), nc_vars%nc_fname)
535  end if
536  else
537  ! interpret shape
538  call get_shape_from_string(idt%shape, array_shape, mf6_input%mempath)
539  ! set nvals
540  nvals = array_shape(1)
541  ! read and set data
542  call nf_verify(nf90_get_var(nc_vars%ncid, varid, dbl1d), nc_vars%nc_fname)
543  end if
544  end subroutine load_double1d_type
545 
546  !> @brief load type 1d double
547  !<
548  subroutine load_double1d_spd(dbl1d, mf6_input, mshape, idt, nc_vars, &
549  iper, input_fname, iaux)
550  use constantsmodule, only: dnodata
551  use netcdfcommonmodule, only: ixstp
552  real(DP), dimension(:), contiguous, pointer, intent(in) :: dbl1d
553  type(modflowinputtype), intent(in) :: mf6_input
554  integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape
555  type(inputparamdefinitiontype), intent(in) :: idt
556  type(ncpackagevarstype), pointer, intent(in) :: nc_vars
557  integer(I4B), intent(in) :: iper
558  character(len=*), intent(in) :: input_fname
559  integer(I4B), optional, intent(in) :: iaux
560  integer(I4B), dimension(:), allocatable :: layer_shape
561  real(DP), dimension(:, :, :), contiguous, pointer :: dbl3d
562  integer(I4B) :: varid, nlay, ncpl, nvals
563  integer(I4B) :: n, istp
564 
565  ! initialize
566  n = 0
567  istp = ixstp()
568 
569  ! set varid
570  if (present(iaux)) then
571  varid = nc_vars%varid(idt%tagname, iaux=iaux)
572  else
573  varid = nc_vars%varid(idt%tagname)
574  end if
575 
576  call get_layered_shape(mshape, nlay, layer_shape)
577  ncpl = product(layer_shape)
578  nvals = product(mshape)
579 
580  select case (idt%shape)
581  case ('NCPL', 'NAUX NCPL')
582  if (nc_vars%grid == 'STRUCTURED') then
583  call nf_verify(nf90_get_var(nc_vars%ncid, varid, dbl1d, &
584  start=(/1, 1, istp/), &
585  count=(/mshape(3), mshape(2), 1/)), &
586  nc_vars%nc_fname)
587  else if (nc_vars%grid == 'LAYERED MESH') then
588  call nf_verify(nf90_get_var(nc_vars%ncid, varid, dbl1d, &
589  start=(/1, istp/), count=(/ncpl, 1/)), &
590  nc_vars%nc_fname)
591  end if
592  case ('NODES', 'NAUX NODES')
593  if (nc_vars%grid == 'STRUCTURED') then
594  dbl3d(1:mshape(3), 1:mshape(2), 1:mshape(1)) => dbl1d(1:nvals)
595  call nf_verify(nf90_get_var(nc_vars%ncid, varid, dbl3d, &
596  start=(/1, 1, 1, istp/), &
597  count=(/mshape(3), mshape(2), mshape(1), &
598  1/)), nc_vars%nc_fname)
599  else if (nc_vars%grid == 'LAYERED MESH') then
600  call nf_verify(nf90_get_var(nc_vars%ncid, varid, dbl1d, &
601  start=(/1, istp/), count=(/nvals, 1/)), &
602  nc_vars%nc_fname)
603  end if
604  case default
605  end select
606  end subroutine load_double1d_spd
607 
608  !> @brief load type 1d double layered
609  !<
610  subroutine load_double1d_layered(dbl1d, mf6_input, mshape, idt, nc_vars, &
611  input_fname)
612  real(DP), dimension(:), contiguous, pointer, intent(in) :: dbl1d
613  type(modflowinputtype), intent(in) :: mf6_input
614  integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape
615  type(inputparamdefinitiontype), intent(in) :: idt
616  type(ncpackagevarstype), pointer, intent(in) :: nc_vars
617  character(len=*), intent(in) :: input_fname
618  integer(I4B), dimension(:), allocatable :: layer_shape
619  integer(I4B) :: nlay, varid
620  integer(I4B) :: k, ncpl
621  integer(I4B) :: index_start, index_stop
622  real(DP), dimension(:), contiguous, pointer :: dbl1d_ptr
623 
624  nullify (dbl1d_ptr)
625  index_start = 1
626  call get_layered_shape(mshape, nlay, layer_shape)
627  ncpl = product(layer_shape)
628 
629  do k = 1, nlay
630  varid = nc_vars%varid(idt%tagname, layer=k)
631  index_stop = index_start + ncpl - 1
632  dbl1d_ptr(1:ncpl) => dbl1d(index_start:index_stop)
633  call nf_verify(nf90_get_var(nc_vars%ncid, varid, dbl1d_ptr), &
634  nc_vars%nc_fname)
635  index_start = index_stop + 1
636  end do
637  end subroutine load_double1d_layered
638 
639  !> @brief load type 1d double layered
640  !<
641  subroutine load_double1d_layered_spd(dbl1d, mf6_input, mshape, idt, nc_vars, &
642  iper, input_fname, iaux)
643  use constantsmodule, only: dnodata
644  use netcdfcommonmodule, only: ixstp
645  real(DP), dimension(:), contiguous, pointer, intent(in) :: dbl1d
646  type(modflowinputtype), intent(in) :: mf6_input
647  integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape
648  type(inputparamdefinitiontype), intent(in) :: idt
649  type(ncpackagevarstype), pointer, intent(in) :: nc_vars
650  integer(I4B), intent(in) :: iper
651  character(len=*), intent(in) :: input_fname
652  integer(I4B), optional, intent(in) :: iaux
653  integer(I4B), dimension(:), allocatable :: layer_shape
654  integer(I4B) :: nlay, varid
655  integer(I4B) :: k, n, ncpl, idx, istp
656  real(DP), dimension(:), contiguous, pointer :: dbl1d_ptr
657 
658  istp = ixstp()
659 
660  call get_layered_shape(mshape, nlay, layer_shape)
661  ncpl = product(layer_shape)
662  allocate (dbl1d_ptr(ncpl))
663 
664  do k = 1, nlay
665  if (present(iaux)) then
666  varid = nc_vars%varid(idt%tagname, layer=k, iaux=iaux)
667  else
668  varid = nc_vars%varid(idt%tagname, layer=k)
669  end if
670  call nf_verify(nf90_get_var(nc_vars%ncid, varid, dbl1d_ptr, &
671  start=(/1, istp/), count=(/ncpl, 1/)), &
672  nc_vars%nc_fname)
673  if (idt%shape == 'NODES' .or. idt%shape == 'NAUX NODES') then
674  do n = 1, ncpl
675  idx = (k - 1) * ncpl + n
676  dbl1d(idx) = dbl1d_ptr(n)
677  end do
678  else if (idt%shape == 'NCPL' .or. idt%shape == 'NAUX NCPL') then
679  do n = 1, ncpl
680  dbl1d(n) = dbl1d_ptr(n)
681  end do
682  end if
683  end do
684 
685  ! cleanup
686  deallocate (dbl1d_ptr)
687  end subroutine load_double1d_layered_spd
688 
689  !> @brief load type 2d double
690  !<
691  subroutine load_double2d_type(dbl2d, mf6_input, mshape, idt, nc_vars, varid, &
692  input_fname)
693  real(DP), dimension(:, :), contiguous, pointer, intent(in) :: dbl2d
694  type(modflowinputtype), intent(in) :: mf6_input
695  integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape
696  type(inputparamdefinitiontype), intent(in) :: idt
697  type(ncpackagevarstype), pointer, intent(in) :: nc_vars
698  integer(I4B), intent(in) :: varid
699  character(len=*), intent(in) :: input_fname
700  real(DP), dimension(:), contiguous, pointer :: dbl1d_ptr
701  integer(I4B), dimension(:), allocatable :: array_shape
702  integer(I4B) :: ncpl, nlay
703 
704  nullify (dbl1d_ptr)
705 
706  if (nc_vars%grid == 'STRUCTURED') then
707  call nf_verify(nf90_get_var(nc_vars%ncid, varid, dbl2d), nc_vars%nc_fname)
708  else if (nc_vars%grid == 'LAYERED MESH') then
709  call get_layered_shape(mshape, nlay, array_shape)
710  ncpl = product(array_shape)
711  dbl1d_ptr(1:ncpl) => dbl2d(:, :)
712  call nf_verify(nf90_get_var(nc_vars%ncid, varid, dbl1d_ptr), &
713  nc_vars%nc_fname)
714  end if
715  end subroutine load_double2d_type
716 
717  !> @brief load type 2d double layered
718  !<
719  subroutine load_double2d_layered(dbl2d, mf6_input, mshape, idt, nc_vars, &
720  input_fname)
721  real(DP), dimension(:, :), contiguous, pointer, intent(in) :: dbl2d
722  type(modflowinputtype), intent(in) :: mf6_input
723  integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape
724  type(inputparamdefinitiontype), intent(in) :: idt
725  type(ncpackagevarstype), pointer, intent(in) :: nc_vars
726  character(len=*), intent(in) :: input_fname
727  integer(I4B), dimension(:), allocatable :: layer_shape
728  integer(I4B) :: k
729  integer(I4B) :: ncpl, nlay, varid
730  real(DP), dimension(:), contiguous, pointer :: dbl1d_ptr
731 
732  nullify (dbl1d_ptr)
733 
734  if (size(mshape) == 3) then
735  write (errmsg, '(a,a,a)') &
736  'Layered netcdf read not supported for DIS dbl2d type ('// &
737  trim(idt%tagname)//').'
738  call store_error(errmsg)
739  call store_error_filename(input_fname)
740  else if (size(mshape) == 2) then
741  call get_layered_shape(mshape, nlay, layer_shape)
742  ncpl = layer_shape(1)
743  do k = 1, nlay
744  varid = nc_vars%varid(idt%tagname, layer=k)
745  dbl1d_ptr(1:ncpl) => dbl2d(1:ncpl, k)
746  call nf_verify(nf90_get_var(nc_vars%ncid, varid, dbl1d_ptr), &
747  nc_vars%nc_fname)
748  end do
749  end if
750  end subroutine load_double2d_layered
751 
752  !> @brief load type 3d double
753  !<
754  subroutine load_double3d_type(dbl3d, mf6_input, mshape, idt, nc_vars, varid, &
755  input_fname)
756  real(DP), dimension(:, :, :), contiguous, pointer, intent(in) :: dbl3d
757  type(modflowinputtype), intent(in) :: mf6_input
758  integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape
759  type(inputparamdefinitiontype), intent(in) :: idt
760  type(ncpackagevarstype), pointer, intent(in) :: nc_vars
761  integer(I4B), intent(in) :: varid
762  character(len=*), intent(in) :: input_fname
763  !
764  call nf_verify(nf90_get_var(nc_vars%ncid, varid, dbl3d), nc_vars%nc_fname)
765  end subroutine load_double3d_type
766 
767  !> @brief load type 3d double layered
768  !<
769  subroutine load_double3d_layered(dbl3d, mf6_input, mshape, idt, nc_vars, &
770  input_fname)
771  real(DP), dimension(:, :, :), contiguous, pointer, intent(in) :: dbl3d
772  type(modflowinputtype), intent(in) :: mf6_input
773  integer(I4B), dimension(:), contiguous, pointer, intent(in) :: mshape
774  type(inputparamdefinitiontype), intent(in) :: idt
775  type(ncpackagevarstype), pointer, intent(in) :: nc_vars
776  character(len=*), intent(in) :: input_fname
777  integer(I4B), dimension(:), allocatable :: layer_shape
778  integer(I4B) :: k !, i, j
779  integer(I4B) :: ncpl, nlay, varid
780  integer(I4B) :: index_start, index_stop
781  real(DP), dimension(:), contiguous, pointer :: dbl1d_ptr
782 
783  nullify (dbl1d_ptr)
784 
785  call get_layered_shape(mshape, nlay, layer_shape)
786 
787  ncpl = product(layer_shape)
788  index_start = 1
789  do k = 1, nlay
790  varid = nc_vars%varid(idt%tagname, layer=k)
791  index_stop = index_start + ncpl - 1
792  dbl1d_ptr(1:ncpl) => dbl3d(:, :, k:k)
793  call nf_verify(nf90_get_var(nc_vars%ncid, varid, dbl1d_ptr), &
794  nc_vars%nc_fname)
795  index_start = index_stop + 1
796  end do
797  end subroutine load_double3d_layered
798 
799 end module ncarrayreadermodule
This module contains simulation constants.
Definition: Constants.f90:9
integer(i4b), parameter linelength
maximum length of a standard line
Definition: Constants.f90:45
real(dp), parameter dnodata
real no data constant
Definition: Constants.f90:95
This module contains the InputDefinitionModule.
This module defines variable data types.
Definition: kind.f90:8
This module contains the ModflowInputModule.
Definition: ModflowInput.f90:9
This module contains the NCArrayReaderModule.
subroutine load_double1d_layered(dbl1d, mf6_input, mshape, idt, nc_vars, input_fname)
load type 1d double layered
subroutine nc_array_load_int2d(int2d, mshape, idt, mf6_input, nc_vars, input_fname, iout)
Load NetCDF integer 2D array.
subroutine load_integer1d_layered_spd(int1d, mf6_input, mshape, idt, nc_vars, iper, input_fname)
load type 1d integer layered
subroutine load_double3d_layered(dbl3d, mf6_input, mshape, idt, nc_vars, input_fname)
load type 3d double layered
subroutine nc_array_load_dbl2d(dbl2d, mshape, idt, mf6_input, nc_vars, input_fname, iout)
Load NetCDF double 2D array.
subroutine load_integer3d_layered(int3d, mf6_input, mshape, idt, nc_vars, input_fname)
load type 3d integer layered
subroutine load_integer1d_type(int1d, mf6_input, mshape, idt, nc_vars, varid, input_fname)
load type 1d integer
subroutine load_double3d_type(dbl3d, mf6_input, mshape, idt, nc_vars, varid, input_fname)
load type 3d double
logical(lgp) function is_layered(grid)
does the grid support per layer variables
subroutine nc_array_load_int3d(int3d, mshape, idt, mf6_input, nc_vars, input_fname, iout)
Load NetCDF integer 3D array.
subroutine load_double2d_type(dbl2d, mf6_input, mshape, idt, nc_vars, varid, input_fname)
load type 2d double
subroutine load_integer1d_spd(int1d, mf6_input, mshape, idt, nc_vars, iper, input_fname)
load type 1d double
subroutine load_integer3d_type(int3d, mf6_input, mshape, idt, nc_vars, varid, input_fname)
load type 3d integer
subroutine load_integer2d_layered(int2d, mf6_input, mshape, idt, nc_vars, input_fname)
load type 2d integer layered
subroutine load_integer2d_type(int2d, mf6_input, mshape, idt, nc_vars, varid, input_fname)
load type 2d integer
subroutine nc_array_load_int1d(int1d, mshape, idt, mf6_input, nc_vars, input_fname, iout, kper)
Load NetCDF integer 1D array.
subroutine load_double1d_layered_spd(dbl1d, mf6_input, mshape, idt, nc_vars, iper, input_fname, iaux)
load type 1d double layered
subroutine load_integer1d_layered(int1d, mf6_input, mshape, idt, nc_vars, input_fname)
load type 1d integer layered
subroutine load_double1d_type(dbl1d, mf6_input, mshape, idt, nc_vars, varid, input_fname)
load type 1d double
subroutine load_double2d_layered(dbl2d, mf6_input, mshape, idt, nc_vars, input_fname)
load type 2d double layered
subroutine nc_array_load_dbl3d(dbl3d, mshape, idt, mf6_input, nc_vars, input_fname, iout)
Load NetCDF double 3D array.
subroutine nc_array_load_dbl1d(dbl1d, mshape, idt, mf6_input, nc_vars, input_fname, iout, kper, iaux)
Load NetCDF double 1D array.
subroutine load_double1d_spd(dbl1d, mf6_input, mshape, idt, nc_vars, iper, input_fname, iaux)
load type 1d double
This module contains the NCFileVarsModule.
Definition: NCFileVars.f90:7
This module contains the NetCDFCommonModule.
Definition: NetCDFCommon.f90:6
integer(i4b) function, public ixstp()
step index for timeseries data
subroutine, public nf_verify(res, nc_fname)
error check a netcdf-fortran interface call
This module contains simulation methods.
Definition: Sim.f90:10
subroutine, public store_error(msg, terminate)
Store an error message.
Definition: Sim.f90:92
subroutine, public store_error_filename(filename, terminate)
Store the erroring file name.
Definition: Sim.f90:203
This module contains simulation variables.
Definition: SimVariables.f90:9
character(len=maxcharlen) errmsg
error message string
This module contains the SourceCommonModule.
Definition: SourceCommon.f90:7
subroutine, public get_layered_shape(mshape, nlay, layer_shape)
subroutine, public get_shape_from_string(shape_string, array_shape, memoryPath)
derived type for storing input definition for a file
Type describing input variables for a package in NetCDF file.
Definition: NCFileVars.f90:22