36 character(len=*),
intent(in) :: grid
37 logical(LGP) :: layered
51 input_fname, iout, kper)
52 integer(I4B),
dimension(:),
pointer,
contiguous,
intent(in) :: int1d
53 integer(I4B),
dimension(:),
contiguous,
pointer,
intent(in) :: mshape
57 character(len=*),
intent(in) :: input_fname
58 integer(I4B),
intent(in) :: iout
59 integer(I4B),
optional,
intent(in) :: kper
60 integer(I4B) :: varid, iper
61 logical(LGP) :: layered
64 layered = (idt%layered .and.
is_layered(nc_vars%grid))
66 if (
present(kper))
then
83 varid = nc_vars%varid(idt%tagname)
94 integer(I4B),
dimension(:, :),
pointer,
contiguous,
intent(in) :: int2d
95 integer(I4B),
dimension(:),
contiguous,
pointer,
intent(in) :: mshape
99 character(len=*),
intent(in) :: input_fname
100 integer(I4B),
intent(in) :: iout
101 integer(I4B) :: varid
102 logical(LGP) :: layered
104 layered = (idt%layered .and.
is_layered(nc_vars%grid))
110 varid = nc_vars%varid(idt%tagname)
120 integer(I4B),
dimension(:, :, :),
pointer,
contiguous,
intent(in) :: int3d
121 integer(I4B),
dimension(:),
contiguous,
pointer,
intent(in) :: mshape
125 character(len=*),
intent(in) :: input_fname
126 integer(I4B),
intent(in) :: iout
127 integer(I4B) :: varid
128 logical(LGP) :: layered
130 layered = (idt%layered .and.
is_layered(nc_vars%grid))
136 varid = nc_vars%varid(idt%tagname)
145 input_fname, iout, kper, iaux)
146 real(DP),
dimension(:),
pointer,
contiguous,
intent(in) :: dbl1d
147 integer(I4B),
dimension(:),
contiguous,
pointer,
intent(in) :: mshape
151 character(len=*),
intent(in) :: input_fname
152 integer(I4B),
intent(in) :: iout
153 integer(I4B),
optional,
intent(in) :: kper
154 integer(I4B),
optional,
intent(in) :: iaux
155 integer(I4B) :: varid, iper
156 logical(LGP) :: layered
159 layered = (idt%layered .and.
is_layered(nc_vars%grid))
161 if (
present(kper))
then
168 iper, input_fname, iaux)
176 iper, input_fname, iaux)
178 varid = nc_vars%varid(idt%tagname)
189 real(DP),
dimension(:, :),
pointer,
contiguous,
intent(in) :: dbl2d
190 integer(I4B),
dimension(:),
contiguous,
pointer,
intent(in) :: mshape
194 character(len=*),
intent(in) :: input_fname
195 integer(I4B),
intent(in) :: iout
196 integer(I4B) :: varid
197 logical(LGP) :: layered
199 layered = (idt%layered .and.
is_layered(nc_vars%grid))
205 varid = nc_vars%varid(idt%tagname)
215 real(DP),
dimension(:, :, :),
pointer,
contiguous,
intent(in) :: dbl3d
216 integer(I4B),
dimension(:),
contiguous,
pointer,
intent(in) :: mshape
220 character(len=*),
intent(in) :: input_fname
221 integer(I4B),
intent(in) :: iout
222 integer(I4B) :: varid
223 logical(LGP) :: layered
225 layered = (idt%layered .and.
is_layered(nc_vars%grid))
231 varid = nc_vars%varid(idt%tagname)
241 integer(I4B),
dimension(:),
contiguous,
pointer,
intent(in) :: int1d
243 integer(I4B),
dimension(:),
contiguous,
pointer,
intent(in) :: mshape
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
256 if (idt%shape ==
'NODES')
then
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), &
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), &
267 else if (
size(mshape) == 1)
then
268 call nf_verify(nf90_get_var(nc_vars%ncid, varid, int1d), nc_vars%nc_fname)
274 nvals = array_shape(1)
276 call nf_verify(nf90_get_var(nc_vars%ncid, varid, int1d), nc_vars%nc_fname)
286 integer(I4B),
dimension(:),
contiguous,
pointer,
intent(in) :: int1d
288 integer(I4B),
dimension(:),
contiguous,
pointer,
intent(in) :: mshape
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
299 varid = nc_vars%varid(idt%tagname)
302 ncpl = product(layer_shape)
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/)), &
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/)), &
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)//
').'
330 integer(I4B),
dimension(:),
contiguous,
pointer,
intent(in) :: int1d
332 integer(I4B),
dimension(:),
contiguous,
pointer,
intent(in) :: mshape
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
346 ncpl = product(layer_shape)
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), &
354 index_start = index_stop + 1
364 integer(I4B),
dimension(:),
contiguous,
pointer,
intent(in) :: int1d
366 integer(I4B),
dimension(:),
contiguous,
pointer,
intent(in) :: mshape
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
378 nvals = product(mshape)
379 ncpl = product(layer_shape)
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/)), &
387 case (
'NODES',
'NAUX NODES')
388 call nf_verify(nf90_get_var(nc_vars%ncid, varid, int1d, &
389 start=(/1, istp/), count=(/nvals, 1/)), &
399 integer(I4B),
dimension(:, :),
contiguous,
pointer,
intent(in) :: int2d
401 integer(I4B),
dimension(:),
contiguous,
pointer,
intent(in) :: mshape
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
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
416 ncpl = product(array_shape)
417 int1d_ptr(1:ncpl) => int2d(:, :)
418 call nf_verify(nf90_get_var(nc_vars%ncid, varid, int1d_ptr), &
427 integer(I4B),
dimension(:, :),
contiguous,
pointer,
intent(in) :: int2d
429 integer(I4B),
dimension(:),
contiguous,
pointer,
intent(in) :: mshape
432 character(len=*),
intent(in) :: input_fname
433 integer(I4B),
dimension(:),
allocatable :: layer_shape
435 integer(I4B) :: ncpl, nlay, varid
436 integer(I4B),
dimension(:),
contiguous,
pointer :: int1d_ptr
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)//
').'
446 else if (
size(mshape) == 2)
then
448 ncpl = layer_shape(1)
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), &
462 integer(I4B),
dimension(:, :, :),
contiguous,
pointer,
intent(in) :: int3d
464 integer(I4B),
dimension(:),
contiguous,
pointer,
intent(in) :: mshape
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)
476 integer(I4B),
dimension(:, :, :),
contiguous,
pointer,
intent(in) :: int3d
478 integer(I4B),
dimension(:),
contiguous,
pointer,
intent(in) :: mshape
481 character(len=*),
intent(in) :: input_fname
482 integer(I4B),
dimension(:),
allocatable :: layer_shape
484 integer(I4B) :: ncpl, nlay, varid
485 integer(I4B) :: index_start, index_stop
486 integer(I4B),
dimension(:),
contiguous,
pointer :: int1d_ptr
491 ncpl = product(layer_shape)
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), &
499 index_start = index_stop + 1
507 real(DP),
dimension(:),
contiguous,
pointer,
intent(in) :: dbl1d
509 integer(I4B),
dimension(:),
contiguous,
pointer,
intent(in) :: mshape
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
522 if (idt%shape ==
'NODES')
then
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), &
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), &
533 else if (
size(mshape) == 1)
then
534 call nf_verify(nf90_get_var(nc_vars%ncid, varid, dbl1d), nc_vars%nc_fname)
540 nvals = array_shape(1)
542 call nf_verify(nf90_get_var(nc_vars%ncid, varid, dbl1d), nc_vars%nc_fname)
549 iper, input_fname, iaux)
552 real(DP),
dimension(:),
contiguous,
pointer,
intent(in) :: dbl1d
554 integer(I4B),
dimension(:),
contiguous,
pointer,
intent(in) :: mshape
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
570 if (
present(iaux))
then
571 varid = nc_vars%varid(idt%tagname, iaux=iaux)
573 varid = nc_vars%varid(idt%tagname)
577 ncpl = product(layer_shape)
578 nvals = product(mshape)
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/)), &
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/)), &
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/)), &
612 real(DP),
dimension(:),
contiguous,
pointer,
intent(in) :: dbl1d
614 integer(I4B),
dimension(:),
contiguous,
pointer,
intent(in) :: mshape
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
627 ncpl = product(layer_shape)
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), &
635 index_start = index_stop + 1
642 iper, input_fname, iaux)
645 real(DP),
dimension(:),
contiguous,
pointer,
intent(in) :: dbl1d
647 integer(I4B),
dimension(:),
contiguous,
pointer,
intent(in) :: mshape
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
661 ncpl = product(layer_shape)
662 allocate (dbl1d_ptr(ncpl))
665 if (
present(iaux))
then
666 varid = nc_vars%varid(idt%tagname, layer=k, iaux=iaux)
668 varid = nc_vars%varid(idt%tagname, layer=k)
670 call nf_verify(nf90_get_var(nc_vars%ncid, varid, dbl1d_ptr, &
671 start=(/1, istp/), count=(/ncpl, 1/)), &
673 if (idt%shape ==
'NODES' .or. idt%shape ==
'NAUX NODES')
then
675 idx = (k - 1) * ncpl + n
676 dbl1d(idx) = dbl1d_ptr(n)
678 else if (idt%shape ==
'NCPL' .or. idt%shape ==
'NAUX NCPL')
then
680 dbl1d(n) = dbl1d_ptr(n)
686 deallocate (dbl1d_ptr)
693 real(DP),
dimension(:, :),
contiguous,
pointer,
intent(in) :: dbl2d
695 integer(I4B),
dimension(:),
contiguous,
pointer,
intent(in) :: mshape
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
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
710 ncpl = product(array_shape)
711 dbl1d_ptr(1:ncpl) => dbl2d(:, :)
712 call nf_verify(nf90_get_var(nc_vars%ncid, varid, dbl1d_ptr), &
721 real(DP),
dimension(:, :),
contiguous,
pointer,
intent(in) :: dbl2d
723 integer(I4B),
dimension(:),
contiguous,
pointer,
intent(in) :: mshape
726 character(len=*),
intent(in) :: input_fname
727 integer(I4B),
dimension(:),
allocatable :: layer_shape
729 integer(I4B) :: ncpl, nlay, varid
730 real(DP),
dimension(:),
contiguous,
pointer :: dbl1d_ptr
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)//
').'
740 else if (
size(mshape) == 2)
then
742 ncpl = layer_shape(1)
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), &
756 real(DP),
dimension(:, :, :),
contiguous,
pointer,
intent(in) :: dbl3d
758 integer(I4B),
dimension(:),
contiguous,
pointer,
intent(in) :: mshape
761 integer(I4B),
intent(in) :: varid
762 character(len=*),
intent(in) :: input_fname
764 call nf_verify(nf90_get_var(nc_vars%ncid, varid, dbl3d), nc_vars%nc_fname)
771 real(DP),
dimension(:, :, :),
contiguous,
pointer,
intent(in) :: dbl3d
773 integer(I4B),
dimension(:),
contiguous,
pointer,
intent(in) :: mshape
776 character(len=*),
intent(in) :: input_fname
777 integer(I4B),
dimension(:),
allocatable :: layer_shape
779 integer(I4B) :: ncpl, nlay, varid
780 integer(I4B) :: index_start, index_stop
781 real(DP),
dimension(:),
contiguous,
pointer :: dbl1d_ptr
787 ncpl = product(layer_shape)
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), &
795 index_start = index_stop + 1
This module contains simulation constants.
integer(i4b), parameter linelength
maximum length of a standard line
real(dp), parameter dnodata
real no data constant
This module defines variable data types.
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.
This module contains the NetCDFCommonModule.
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.
subroutine, public store_error(msg, terminate)
Store an error message.
subroutine, public store_error_filename(filename, terminate)
Store the erroring file name.
This module contains simulation variables.
character(len=maxcharlen) errmsg
error message string
This module contains the SourceCommonModule.
subroutine, public get_layered_shape(mshape, nlay, layer_shape)
subroutine, public get_shape_from_string(shape_string, array_shape, memoryPath)
Type describing input variables for a package in NetCDF file.