20#if defined PIO_LIB && defined DISTRIBUTE
26#if defined PIO_LIB && defined DISTRIBUTE
34#if defined PIO_LIB && defined DISTRIBUTE
80#if defined PIO_LIB && defined DISTRIBUTE
87#if defined PARALLEL_IO && defined DISTRIBUTE
91 & ncvarid, tindex, gtype, &
92 & LBi, UBi, LBj, UBj, LBk, UBk, Ascl, &
99 & MinValue, MaxValue)
RESULT (status)
104# if defined WRITE_WATER && defined MASKING
112 logical,
intent(in),
optional :: setfillval
114 integer,
intent(in) :: ng, model, ncid, ncvarid, tindex, gtype
115 integer,
intent(in) :: ifield
116 integer,
intent(in) :: lbi, ubi, lbj, ubj, lbk, ubk
118 integer,
intent(in),
optional :: extractfield
120 real(dp),
intent(in) :: ascl
124 real(r8),
intent(in) :: amask(lbi:,lbj:)
126 real(r8),
intent(in) :: adat(lbi:,lbj:,lbk:)
129 real(r8),
intent(in) :: amask(lbi:ubi,lbj:ubj)
131 real(r8),
intent(in) :: adat(lbi:ubi,lbj:ubj,lbk:ubk)
133 real(r8),
intent(out),
optional :: minvalue
134 real(r8),
intent(out),
optional :: maxvalue
140 integer :: extract_flag
141 integer :: i, ic, j, jc, k, kc, npts
142 integer :: imin, imax, jmin, jmax
143 integer :: ioff, joff, koff
144 integer :: istr, iend
145 integer :: ilen, isize, jlen, jsize, ijlen, klen
149 integer,
dimension(4) :: start, total
151 real(r8),
parameter :: inival = 0.0_r8
153 real(r8),
allocatable :: awrk(:)
166 SELECT CASE (abs(mytype))
211 IF (
PRESENT(setfillval))
THEN
223 IF (
PRESENT(extractfield))
THEN
224 extract_flag=extractfield
241 SELECT CASE (abs(mytype))
270 IF (.not.
allocated(awrk))
THEN
271 allocate ( awrk(npts) )
282 awrk(ic)=adat(i,j,k)*ascl
284 IF (abs(awrk(ic)).eq.0.0_r8)
THEN
289 IF ((amask(i,j).eq.0.0_r8).and.landfill)
THEN
308 status=nf90_put_var(ncid, ncvarid, awrk, start, total)
311# if defined WRITE_WATER && defined MASKING
321 SELECT CASE (abs(mytype))
350 IF (.not.
allocated(awrk))
THEN
351 allocate ( awrk(npts) )
365 awrk(ic)=adat(i,j,k)*ascl
367 IF (abs(awrk(ic)).eq.0.0_r8)
THEN
371 IF (amask(i,j).eq.0.0_r8)
THEN
380 CALL mp_collect (ng, model, npts, inival, awrk)
386 IF (awrk(i).lt.
spval)
THEN
403 status=nf90_put_var(ncid, ncvarid, awrk(istr:), start, total)
411 IF (
allocated(awrk))
THEN
423 & ncvarid, tindex, gtype, &
424 & LBi, UBi, LBj, UBj, LBk, UBk, Ascl, &
431 & MinValue, MaxValue)
RESULT (status)
451 logical,
intent(in),
optional :: setfillval
453 integer,
intent(in) :: ng, model, ncid, ncvarid, tindex, gtype
454 integer,
intent(in) :: ifield
455 integer,
intent(in) :: lbi, ubi, lbj, ubj, lbk, ubk
457 integer,
intent(in),
optional :: extractfield
459 real(dp),
intent(in) :: ascl
463 real(r8),
intent(in) :: amask(lbi:,lbj:)
465 real(r8),
intent(in) :: adat(lbi:,lbj:,lbk:)
468 real(r8),
intent(in) :: amask(lbi:ubi,lbj:ubj)
470 real(r8),
intent(in) :: adat(lbi:ubi,lbj:ubj,lbk:ubk)
472 real(r8),
intent(out),
optional :: minvalue
473 real(r8),
intent(out),
optional :: maxvalue
479 integer :: extract_flag
480 integer :: i, npts, tile
481 integer :: imin, imax, jmin, jmax, koff
482 integer :: ilen, jlen, klen, ijlen, mytype
485 integer,
dimension(4) :: start, total
487 real(r8),
dimension((Lm(ng)+2)*(Mm(ng)+2)*(UBk-LBk+1)) :: awrk
511 IF (
PRESENT(setfillval))
THEN
523 IF (
PRESENT(extractfield))
THEN
524 extract_flag=extractfield
533 stats % checksum=0_i8b
551 & gtype, ifield, tindex, &
552 & landfill, extract_flag, &
553 & lbi, ubi, lbj, ubj, lbk, ubk, &
558 & start, total, npts, awrk)
564 IF (
PRESENT(minvalue))
THEN
569 IF (abs(awrk(i)).lt.
spval)
THEN
570 minvalue=min(minvalue,awrk(i))
571 maxvalue=max(maxvalue,awrk(i))
582 status=nf90_put_var(ncid, ncvarid, awrk, start, total)
593 & lbi, ubi, lbj, ubj, lbk, ubk, &
600 WRITE (stdout,10) trim(
vname(1,ifield)), stats%min, stats%max, &
602 10
FORMAT (4x,
'- ',a,
':',t35,
'Min = ',1p,e15.8,
', Max = ', &
603 & 1p,e15.8,
', CheckSum = ',i0)
620#if defined PIO_LIB && defined DISTRIBUTE
624 & pioVar, tindex, pioDesc, &
625 & LBi, UBi, LBj, UBj, LBk, UBk, Ascl, &
632 & MinValue, MaxValue)
RESULT (status)
644 logical,
intent(in),
optional :: setfillval
646 integer,
intent(in) :: ng, model, tindex
647 integer,
intent(in) :: ifield
648 integer,
intent(in) :: lbi, ubi, lbj, ubj, lbk, ubk
650 integer,
intent(in),
optional :: extractfield
652 real(dp),
intent(in) :: ascl
656 real(r8),
intent(in) :: amask(lbi:,lbj:)
658 real(r8),
intent(in) :: adat(lbi:,lbj:,lbk:)
661 real(r8),
intent(in) :: amask(lbi:ubi,lbj:ubj)
663 real(r8),
intent(in) :: adat(lbi:ubi,lbj:ubj,lbk:ubk)
665 real(r8),
intent(out),
optional :: minvalue
666 real(r8),
intent(out),
optional :: maxvalue
668 TYPE (file_desc_t),
intent(inout) :: piofile
669 TYPE (io_desc_t),
intent(inout) :: piodesc
670 TYPE (my_vardesc),
intent(inout) :: piovar
674 logical :: landfill, lminmax
676 logical,
pointer :: lwater(:,:,:)
678 integer :: extract_flag
679 integer :: i, j, k, tile
680 integer :: imin, imax, jmin, jmax
681 integer :: cgrid, dkind, ghost, gtype
684 integer,
dimension(4) :: start, total
686 real(r8),
dimension(2) :: rbuffer
688 real(r4),
pointer :: awrk4(:,:,:)
689 real(r8),
pointer :: awrk8(:,:,:)
696 character (len= 3),
dimension(2) :: op_handle
711 SELECT CASE (abs(gtype))
731 IF (
PRESENT(minvalue))
THEN
733 IF (.not.
associated(lwater))
THEN
734 allocate ( lwater(lbi:ubi,lbj:ubj,lbk:ubk) )
744 IF (
PRESENT(setfillval))
THEN
756 IF (
PRESENT(extractfield))
THEN
757 extract_flag=extractfield
766 stats % checksum=0_i8b
781 IF (dkind.eq.pio_double)
THEN
782 IF (.not.
associated(awrk8))
THEN
783 allocate ( awrk8(lbi:ubi,lbj:ubj,lbk:ubk) )
790 awrk8(i,j,k)=adat(i,j,k)*ascl
792 IF ((amask(i,j).eq.0.0_r8).and.landfill)
THEN
794 IF (lminmax) lwater(i,j,k)=.false.
801 rbuffer(1)=minval(awrk8, mask=lwater)
802 rbuffer(2)=maxval(awrk8, mask=lwater)
805 IF (.not.
associated(awrk4))
THEN
806 allocate ( awrk4(lbi:ubi,lbj:ubj,lbk:ubk) )
813 awrk4(i,j,k)=real(adat(i,j,k)*ascl, r4)
815 IF ((amask(i,j).eq.0.0_r8).and.landfill)
THEN
816 awrk4(i,j,k)=real(
spval, r4)
823 rbuffer(1)=real(minval(awrk4, mask=lwater),r8)
824 rbuffer(2)=real(maxval(awrk4, mask=lwater),r8)
830 IF (tindex.gt.0)
THEN
831 CALL pio_setframe (piofile, &
833 & int(tindex, kind=pio_offset_kind))
838 IF (dkind.eq.pio_double)
THEN
839 CALL pio_write_darray (piofile, &
842 & awrk8(imin:imax,jmin:jmax,lbk:ubk), &
845 CALL pio_write_darray (piofile, &
848 & awrk4(imin:imax,jmin:jmax,lbk:ubk), &
865 & lbi, ubi, lbj, ubj, lbk, ubk, &
870 & extract_flag = extract_flag, &
873 WRITE (stdout,10) trim(
vname(1,ifield)), stats%min, stats%max, &
875 10
FORMAT (4x,
'- ',a,
':',t35,
'Min = ',1p,e15.8,
', Max = ', &
876 & 1p,e15.8,
', CheckSum = ',i0)
887 CALL mp_reduce (ng, model, 2, rbuffer, op_handle)
890 IF (
associated(lwater))
deallocate (lwater)
895 IF (dkind.eq.pio_double)
THEN
896 IF (
associated(awrk8))
deallocate (awrk8)
898 IF (
associated(awrk4))
deallocate (awrk4)
subroutine mp_gather3d(ng, model, lbi, ubi, lbj, ubj, lbk, ubk, tindex, gtype, ascl, amask, a, npts, awrk, setfillval)
subroutine mp_gather2d(ng, model, lbi, ubi, lbj, ubj, tindex, gtype, ascl, amask, a, npts, awrk, setfillval)
subroutine, public tile_bounds_1d(ng, tile, imax, istr, iend)
character(len=maxlen), dimension(6, 0:nv) vname
type(t_bounds), dimension(:), allocatable bounds
integer, parameter b3dvar
integer, parameter r3dvar
type(t_iobounds), dimension(:), allocatable iobounds
integer, parameter u3dvar
integer, parameter u2dvar
integer, parameter p2dvar
integer, parameter r2dvar
integer, parameter l3dvar
integer, parameter v2dvar
integer, parameter p3dvar
integer, parameter v3dvar
real(dp), parameter spval
integer function nf90_fwrite3d(ng, model, ncid, ifield, ncvarid, tindex, gtype, lbi, ubi, lbj, ubj, lbk, ubk, ascl, amask, adat, setfillval, extractfield, minvalue, maxvalue)
integer function pio_fwrite3d(ng, model, piofile, ifield, piovar, tindex, piodesc, lbi, ubi, lbj, ubj, lbk, ubk, ascl, amask, adat, setfillval, extractfield, minvalue, maxvalue)
subroutine, public pack_field3d(ng, model, tile, gtype, ifield, tindex, landfill, extract_flag, lbi, ubi, lbj, ubj, lbk, ubk, amask, ascl, adat, start, total, npts, awrk)
subroutine, public stats_3dfld(ng, tile, model, gtype, s, extract_flag, lbi, ubi, lbj, ubj, lbk, ubk, f, fmask, debug)