#include "cppdefs.h" MODULE nf_fwrite4d_mod ! !git $Id$ !svn $Id: nf_fwrite4d.F 1190 2023-08-18 19:51:09Z arango $ !================================================== Hernan G. Arango === ! Copyright (c) 2002-2023 The ROMS/TOMS Group ! ! Licensed under a MIT/X style license ! ! See License_ROMS.md ! !======================================================================= ! ! ! This function writes out a generic floating point 4D array into an ! ! output NetCDF file. ! ! ! ! On Input: ! ! ! ! ng Nested grid number ! ! model Calling model identifier ! ! ncid NetCDF file ID ! ! ifield Field metadata index (integer) ! #if defined PIO_LIB && defined DISTRIBUTE !or pioFile PIO file descriptor structure, TYPE(file_desc_t) ! ! pioFile%fh file handler ! ! pioFile%iosystem IO system descriptor (struct) ! #endif ! ncvarid NetCDF variable ID ! #if defined PIO_LIB && defined DISTRIBUTE !or pioVar PIO variable descriptor structure, TYPE(My_VarDesc) ! ! pioVar%vd variable descriptor TYPE(Var_Desc_t)! ! pioVar%dkind variable data kind ! ! pioVar%gtype variable C-gridtype ! #endif ! tindex NetCDF time record index to write ! ! gtype Grid type. If negative, only write water points ! #if defined PIO_LIB && defined DISTRIBUTE !or pioDesc IO data decomposition descriptor, TYPE(IO_desc_t) ! #endif ! LBi I-dimension Lower bound ! ! UBi I-dimension Upper bound ! ! LBj J-dimension Lower bound ! ! UBj J-dimension Upper bound ! ! LBk K-dimension Lower bound ! ! UBk K-dimension Upper bound ! ! LBt Time-dimension Lower bound ! ! UBt Time-dimension Upper bound ! ! Amask land/Sea mask, if any (real) ! ! Ascl Factor to scale field before writing (real) ! ! Adat Field to write out (real) ! ! SetFillVal Logical switch to set fill value in land areas ! ! (OPTIONAL) ! ! ! ! On Output: ! ! ! ! status Error flag (integer) ! ! MinValue Minimum value (real, OPTIONAL) ! ! MaxValue Maximum value (real, OPTIONAL) ! ! ! #ifdef POSITIVE_ZERO ! Starting F95 zero values can be signed (-0 or +0) following the ! ! IEEE 754 floating point standard. This may produce different ! ! output data in serial and parallel applications. Since comparing ! ! serial and parallel output is essential for tracking parallel ! ! partition bugs, "positive zero" is enforced. ! ! ! #endif !======================================================================= ! USE mod_param USE mod_parallel USE mod_ncparam USE mod_scalars ! implicit none ! INTERFACE nf_fwrite4d MODULE PROCEDURE nf90_fwrite4d #if defined PIO_LIB && defined DISTRIBUTE MODULE PROCEDURE pio_fwrite4d #endif END INTERFACE nf_fwrite4d ! CONTAINS #if defined PARALLEL_IO && defined DISTRIBUTE ! !*********************************************************************** FUNCTION nf90_fwrite4d (ng, model, ncid, ifield, & & ncvarid, tindex, gtype, & & LBi, UBi, LBj, UBj, LBk, UBk, LBt, UBt, & & Ascl, & # ifdef MASKING & Amask, & # endif & Adat, SetFillVal, & & MinValue, MaxValue) RESULT (status) !*********************************************************************** ! USE mod_netcdf # ifdef MASKING ! USE distribute_mod, ONLY : mp_collect # endif ! ! Imported variable declarations. ! logical, intent(in), optional :: SetFillVal ! integer, intent(in) :: ng, model, ncid, ncvarid, tindex, gtype integer, intent(in) :: ifield integer, intent(in) :: LBi, UBi, LBj, UBj, LBk, UBk, LBt, UBt ! real(dp), intent(in) :: Ascl ! # ifdef ASSUMED_SHAPE # ifdef MASKING real(r8), intent(in) :: Amask(LBi:,LBj:) # endif real(r8), intent(in) :: Adat(LBi:,LBj:,LBk:,LBt:) # else # ifdef MASKING real(r8), intent(in) :: Amask(LBi:UBi,LBj:UBj) # endif real(r8), intent(in) :: Adat(LBi:UBi,LBj:UBj,LBk:UBk,LBt:UBt) # endif real(r8), intent(out), optional :: MinValue real(r8), intent(out), optional :: MaxValue ! ! Local variable declarations. ! # ifdef MASKING logical :: LandFill ! # endif integer :: i, ic, j, jc, k, kc, l, lc, Npts integer :: Imin, Imax, Jmin, Jmax, Kmin, Kmax integer :: Ioff, Joff, Koff, Loff integer :: Istr, Iend integer :: Ilen, Isize, Jlen, Jsize, IJlen, IJKlen, Klen, Llen integer :: MyType integer :: status integer, dimension(5) :: start, total ! real(r8), parameter :: IniVal = 0.0_r8 real(r8), allocatable :: Awrk(:) ! !----------------------------------------------------------------------- ! Set starting and ending indices to process. !----------------------------------------------------------------------- ! status=nf90_noerr ! ! Set first and last grid point according to staggered C-grid ! classification. ! MyType=gtype ! SELECT CASE (ABS(MyType)) CASE (p2dvar, p3dvar) Imin=BOUNDS(ng)%Istr (MyRank) Imax=BOUNDS(ng)%IendR(MyRank) Jmin=BOUNDS(ng)%Jstr (MyRank) Jmax=BOUNDS(ng)%JendR(MyRank) Isize=IOBOUNDS(ng)%xi_psi Jsize=IOBOUNDS(ng)%eta_psi CASE (r2dvar, r3dvar) Imin=BOUNDS(ng)%IstrR(MyRank) Imax=BOUNDS(ng)%IendR(MyRank) Jmin=BOUNDS(ng)%JstrR(MyRank) Jmax=BOUNDS(ng)%JendR(MyRank) Isize=IOBOUNDS(ng)%xi_rho Jsize=IOBOUNDS(ng)%eta_rho CASE (u2dvar, u3dvar) Imin=BOUNDS(ng)%Istr (MyRank) Imax=BOUNDS(ng)%IendR(MyRank) Jmin=BOUNDS(ng)%JstrR(MyRank) Jmax=BOUNDS(ng)%JendR(MyRank) Isize=IOBOUNDS(ng)%xi_u Jsize=IOBOUNDS(ng)%eta_u CASE (v2dvar, v3dvar) Imin=BOUNDS(ng)%IstrR(MyRank) Imax=BOUNDS(ng)%IendR(MyRank) Jmin=BOUNDS(ng)%Jstr (MyRank) Jmax=BOUNDS(ng)%JendR(MyRank) Isize=IOBOUNDS(ng)%xi_v Jsize=IOBOUNDS(ng)%eta_v CASE DEFAULT Imin=BOUNDS(ng)%IstrR(MyRank) Imax=BOUNDS(ng)%IendR(MyRank) Jmin=BOUNDS(ng)%JstrR(MyRank) Jmax=BOUNDS(ng)%JendR(MyRank) Isize=IOBOUNDS(ng)%xi_rho Jsize=IOBOUNDS(ng)%eta_rho END SELECT ! Ilen=Imax-Imin+1 Jlen=Jmax-Jmin+1 Klen=UBk-LBk+1 Llen=UBt-LBt+1 # ifdef MASKING ! ! Set switch to replace land areas with fill value, spval. ! IF (PRESENT(SetFillVal)) THEN LandFill=SetFillVal ELSE LandFill=tindex.gt.0 END IF # endif ! !----------------------------------------------------------------------- ! Parallel I/O: Pack tile data into 1D array in column-major order. # ifdef MASKING ! Overwrite masked points with special value. # endif !----------------------------------------------------------------------- ! IF (gtype.gt.0) THEN ! ! Set offsets due the NetCDF dimensions. Recall that some output ! variables not always start at one. ! SELECT CASE (ABS(MyType)) CASE (p2dvar, p3dvar) Ioff=0 Joff=0 CASE (r2dvar, r3dvar) Ioff=1 Joff=1 CASE (u2dvar, u3dvar) Ioff=0 Joff=1 CASE (v2dvar, v3dvar) Ioff=1 Joff=0 CASE DEFAULT Ioff=1 Joff=1 END SELECT IF (LBk.eq.0) THEN Koff=1 ELSE Koff=0 END IF IF (LBt.eq.0) THEN Loff=1 ELSE Loff=0 END IF ! ! Allocate and initialize scratch work array. ! IJlen=Ilen*Jlen IJKlen=IJlen*Klen Npts=IJKlen*Llen IF (.not.allocated(Awrk)) THEN allocate ( Awrk(Npts) ) Awrk=IniVal END IF ! ! Pack and scale tile data. ! ic=0 DO l=LBt,UBt DO k=LBk,UBk DO j=Jmin,Jmax DO i=Imin,Imax ic=ic+1 Awrk(ic)=Adat(i,j,k,l)*Ascl # ifdef POSITIVE_ZERO IF (ABS(Awrk(ic)).eq.0.0_r8) THEN Awrk(ic)=0.0_r8 ! impose positive zero END IF # endif # ifdef MASKING IF ((Amask(i,j).eq.0.0_r8).and.LandFill) THEN Awrk(ic)=spval END IF # endif END DO END DO END DO END DO ! ! Write out data: all parallel nodes write their own packed tile data. ! start(1)=Imin+Ioff total(1)=Ilen start(2)=Jmin+Joff total(2)=Jlen start(3)=LBk+Koff total(3)=Klen start(4)=LBt+Loff total(4)=Llen start(5)=tindex total(5)=1 status=nf90_put_var(ncid, ncvarid, Awrk, start, total) END IF # if defined WRITE_WATER && defined MASKING ! !----------------------------------------------------------------------- ! Parallel I/O: Remove land points and pack tile data into 1D array. !----------------------------------------------------------------------- ! IF (gtype.lt.0) THEN ! ! Set offsets due array packing into 1D array in column-major order. ! SELECT CASE (ABS(MyType)) CASE (p2dvar, p3dvar) Ioff=0 Joff=1 CASE (r2dvar, r3dvar) Ioff=1 Joff=0 CASE (u2dvar, u3dvar) Ioff=0 Joff=0 CASE (v2dvar, v3dvar) Ioff=1 Joff=1 CASE DEFAULT Ioff=1 Joff=0 END SELECT IF (LBk.eq.0) THEN Koff=0 ELSE Koff=1 END IF IF (LBt.eq.0) THEN Loff=0 ELSE Loff=1 END IF ! ! Allocate and initialize scratch work array. ! IJlen=Isize*Jsize IJKlen=IJlen*Klen Npts=IJKlen*Llen IF (.not.allocated(Awrk)) THEN allocate ( Awrk(Npts) ) Awrk=IniVal END IF ! ! Scale and gather data from all spawned nodes. Store data into a 1D ! global array, packed in column-major order. Flag land point with ! special value. ! DO l=LBt,UBt lc=(l-Loff)*IJKlen DO k=LBk,UBk kc=(k-Koff)*IJlen+lc DO j=Jmin,Jmax jc=(j-Joff)*Isize+kc DO i=Imin,Imax ic=i+Ioff+jc Awrk(ic)=Adat(i,j,k,l)*Ascl # ifdef POSITIVE_ZERO IF (ABS(Awrk(ic)).eq.0.0_r8) THEN Awrk(ic)=0.0_r8 ! impose positive zero END IF # endif IF (Amask(i,j).eq.0.0_r8) THEN Awrk(ic)=spval END IF END DO END DO END DO END DO ! ! Global reduction of work array. ! CALL mp_collect (ng, model, Npts, IniVal, Awrk) ! ! Remove land points. ! ic=0 DO i=1,Npts IF (Awrk(i).lt.spval) THEN ic=ic+1 Awrk(ic)=Awrk(i) END IF END DO Npts=ic ! ! Write out data: all parallel nodes write a section of the packed ! data. ! CALL tile_bounds_1d (ng, MyRank, Npts, Istr, Iend) ! start(1)=Istr total(1)=Iend-Istr+1 start(2)=tindex total(2)=1 status=nf90_put_var(ncid, ncvarid, Awrk(Istr:), start, total) END IF # endif ! !----------------------------------------------------------------------- ! Deallocate scratch work array. !----------------------------------------------------------------------- ! IF (allocated(Awrk)) THEN deallocate (Awrk) END IF ! RETURN END FUNCTION nf90_fwrite4d #else ! !*********************************************************************** FUNCTION nf90_fwrite4d (ng, model, ncid, ifield, & & ncvarid, tindex, gtype, & & LBi, UBi, LBj, UBj, LBk, UBk, LBt, UBt, & & Ascl, & # ifdef MASKING & Amask, & # endif & Adat, SetFillVal, & & MinValue, MaxValue) RESULT (status) !*********************************************************************** ! USE mod_netcdf # ifdef DISTRIBUTE ! USE distribute_mod, ONLY : mp_bcasti # ifdef INLINE_2DIO USE distribute_mod, ONLY : mp_gather2d # else USE distribute_mod, ONLY : mp_gather3d # endif # endif # ifdef OUTPUT_STATS USE stats_mod, ONLY : stats_4dfld # endif ! ! Imported variable declarations. ! logical, intent(in), optional :: SetFillVal ! integer, intent(in) :: ng, model, ncid, ncvarid, tindex, gtype integer, intent(in) :: ifield integer, intent(in) :: LBi, UBi, LBj, UBj, LBk, UBk, LBt, UBt ! real(dp), intent(in) :: Ascl ! # ifdef ASSUMED_SHAPE # ifdef MASKING real(r8), intent(in) :: Amask(LBi:,LBj:) # endif real(r8), intent(in) :: Adat(LBi:,LBj:,LBk:,LBt:) # else # ifdef MASKING real(r8), intent(in) :: Amask(LBi:UBi,LBj:UBj) # endif real(r8), intent(in) :: Adat(LBi:UBi,LBj:UBj,LBk:UBk,LBt:UBt) # endif real(r8), intent(out), optional :: MinValue real(r8), intent(out), optional :: MaxValue ! ! Local variable declarations. ! # ifdef MASKING logical :: LandFill ! # endif integer :: i, j, k, ic, fourth, Npts, tile integer :: Imin, Imax, Jmin, Jmax, Kmin, Kmax, Koff, Loff integer :: Ilen, Jlen, Klen, IJlen, MyType integer :: status integer, dimension(5) :: start, total ! # if defined INLINE_2DIO && defined DISTRIBUTE real(r8), dimension((Lm(ng)+2)*(Mm(ng)+2)) :: Awrk # else real(r8), dimension((Lm(ng)+2)*(Mm(ng)+2)*(UBk-LBk+1)) :: Awrk # endif # ifdef OUTPUT_STATS ! TYPE (T_STATS) :: Stats # endif ! !----------------------------------------------------------------------- ! Set starting and ending indices to process. !----------------------------------------------------------------------- ! status=nf90_noerr ! ! Set first and last grid point according to staggered C-grid ! classification. Set loops offsets. ! MyType=gtype ! SELECT CASE (ABS(MyType)) CASE (p2dvar, p3dvar) Imin=IOBOUNDS(ng)%ILB_psi Imax=IOBOUNDS(ng)%IUB_psi Jmin=IOBOUNDS(ng)%JLB_psi Jmax=IOBOUNDS(ng)%JUB_psi CASE (r2dvar, r3dvar) Imin=IOBOUNDS(ng)%ILB_rho Imax=IOBOUNDS(ng)%IUB_rho Jmin=IOBOUNDS(ng)%JLB_rho Jmax=IOBOUNDS(ng)%JUB_rho CASE (u2dvar, u3dvar) Imin=IOBOUNDS(ng)%ILB_u Imax=IOBOUNDS(ng)%IUB_u Jmin=IOBOUNDS(ng)%JLB_u Jmax=IOBOUNDS(ng)%JUB_u CASE (v2dvar, v3dvar) Imin=IOBOUNDS(ng)%ILB_v Imax=IOBOUNDS(ng)%IUB_v Jmin=IOBOUNDS(ng)%JLB_v Jmax=IOBOUNDS(ng)%JUB_v CASE DEFAULT Imin=IOBOUNDS(ng)%ILB_rho Imax=IOBOUNDS(ng)%IUB_rho Jmin=IOBOUNDS(ng)%JLB_rho Jmax=IOBOUNDS(ng)%JUB_rho END SELECT ! Ilen=Imax-Imin+1 Jlen=Jmax-Jmin+1 Klen=UBk-LBk+1 IJlen=Ilen*Jlen ! IF (LBk.eq.0) THEN Koff=0 ELSE Koff=1 END IF ! IF (LBt.eq.0) THEN Loff=1 ELSE Loff=0 END IF # ifdef MASKING ! ! Set switch to replace land areas with fill value, spval. ! IF (PRESENT(SetFillVal)) THEN LandFill=SetFillVal ELSE LandFill=tindex.gt.0 END IF # endif ! ! If appropriate, initialize minimum and maximum processed values. ! IF (PRESENT(MinValue)) THEN MinValue=spval MaxValue=-spval END IF ! ! Initialize local array to avoid denormalized numbers. This ! facilitates processing and debugging. ! Awrk=0.0_r8 # ifdef OUTPUT_STATS ! ! Initialize output variable statistics structure. ! Stats % checksum=0_i8b Stats % count=0 Stats % min=spval Stats % max=-spval Stats % avg=0.0_r8 Stats % rms=0.0_r8 # endif # ifdef DISTRIBUTE ! !----------------------------------------------------------------------- ! If distributed-memory set-up, collect tile data from all spawned ! nodes and store it into a global scratch 1D array, packed in column- ! major order. # ifdef MASKING # ifdef WRITE_WATER ! Remove land points and pack water points into 1D-array. # else ! Overwrite masked points with special value. # endif # endif !----------------------------------------------------------------------- ! DO fourth=LBt,UBt # ifdef INLINE_2DIO ! ! Process the data as 2D slices (level by level) to reduce memory ! requirements. ! DO k=LBk,UBk CALL mp_gather2d (ng, model, LBi, UBi, LBj, UBj, & & tindex, gtype, Ascl, & # ifdef MASKING & Amask, & # endif & Adat(:,:,k,fourth), Npts, Awrk, SetFillVal) # else ! ! Process the data as 3D slices. ! CALL mp_gather3d (ng, model, LBi, UBi, LBj, UBj, LBk, UBk, & & tindex, gtype, Ascl, & # ifdef MASKING & Amask, & # endif & Adat(:,:,:,fourth), Npts, Awrk, SetFillVal) # endif ! !----------------------------------------------------------------------- ! If applicable, compute output field minimum and maximum values. !----------------------------------------------------------------------- ! IF (PRESENT(MinValue)) THEN IF (OutThread) THEN DO i=1,Npts IF (ABS(Awrk(i)).lt.spval) THEN MinValue=MIN(MinValue,Awrk(i)) MaxValue=MAX(MaxValue,Awrk(i)) END IF END DO END IF END IF ! !----------------------------------------------------------------------- ! Write output buffer into NetCDF file. !----------------------------------------------------------------------- ! IF (OutThread) THEN IF (gtype.gt.0) THEN start(1)=1 total(1)=Ilen start(2)=1 total(2)=Jlen # ifdef INLINE_2DIO start(3)=k-Koff+1 total(3)=1 # else start(3)=1 total(3)=Klen # endif start(4)=fourth+Loff total(4)=1 start(5)=tindex total(5)=1 # ifdef MASKING ELSE start(1)=1+(fourth+Loff-1)*Npts total(1)=Npts start(2)=tindex total(2)=1 # endif END IF # ifdef POSITIVE_ZERO DO ic=1,Npts IF (ABS(Awrk(ic)).eq.0.0_r8) THEN Awrk(ic)=0.0_r8 ! impose positive zero END IF END DO # endif status=nf90_put_var(ncid, ncvarid, Awrk, start, total) END IF # ifdef INLINE_2DIO END DO # endif END DO # else ! !----------------------------------------------------------------------- ! If serial or shared-memory applications and serial output, pack data ! into a global 1D array in column-major order. # ifdef MASKING # ifdef WRITE_WATER ! Remove land points and pack water points into 1D-array. # else ! Overwrite masked points with special value. # endif # endif !----------------------------------------------------------------------- ! ! Process data as 3D slices. ! DO fourth=LBt,UBt IF (gtype.gt.0) THEN ic=0 Npts=IJlen*Klen DO k=LBk,UBk DO j=Jmin,Jmax DO i=Imin,Imax ic=ic+1 Awrk(ic)=Adat(i,j,k,fourth)*Ascl # ifdef MASKING IF ((Amask(i,j).eq.0.0_r8).and.LandFill) THEN Awrk(ic)=spval END IF # endif END DO END DO END DO # ifdef MASKING ELSE Npts=0 DO k=LBk,UBk DO j=Jmin,Jmax DO i=Imin,Imax IF (Amask(i,j).gt.0.0_r8) THEN Npts=Npts+1 Awrk(Npts)=Adat(i,j,k,fourth)*Ascl END IF END DO END DO END DO # endif END IF ! !----------------------------------------------------------------------- ! If applicable, compute output field minimum and maximum values. !----------------------------------------------------------------------- ! IF (PRESENT(MinValue)) THEN IF (OutThread) THEN DO i=1,Npts IF (ABS(Awrk(i)).lt.spval) THEN MinValue=MIN(MinValue,Awrk(i)) MaxValue=MAX(MaxValue,Awrk(i)) END IF END DO END IF END IF ! !----------------------------------------------------------------------- ! Write output buffer into NetCDF file. !----------------------------------------------------------------------- ! IF (OutThread) THEN IF (gtype.gt.0) THEN start(1)=1 total(1)=Ilen start(2)=1 total(2)=Jlen start(3)=1 total(3)=Klen start(4)=fourth+Loff total(4)=1 start(5)=tindex total(5)=1 # ifdef MASKING ELSE start(1)=1+(fourth+Loff-1)*Npts total(1)=Npts start(2)=tindex total(2)=1 # endif END IF # ifdef POSITIVE_ZERO DO ic=1,Npts IF (ABS(Awrk(ic)).eq.0.0_r8) THEN Awrk(ic)=0.0_r8 ! impose positive zero END IF END DO # endif status=nf90_put_var(ncid, ncvarid, Awrk, start, total) END IF END DO # endif # ifdef OUTPUT_STATS ! !----------------------------------------------------------------------- ! Compute and report output field statistics. !----------------------------------------------------------------------- ! # ifdef DISTRIBUTE tile=MyRank # else tile=-1 # endif CALL stats_4dfld (ng, tile, model, gtype, Stats, & & LBi, UBi, LBj, UBj, LBk, UBk, LBt, UBt, & & Adat, & # ifdef MASKING & Fmask = Amask, & # endif & debug = .FALSE.) IF (OutThread) THEN WRITE (stdout,10) TRIM(Vname(1,ifield)), Stats%min, Stats%max, & & Stats%checksum 10 FORMAT (4x,'- ',a,':',t35,'Min = ',1p,e15.8,', Max = ', & & 1p,e15.8,', CheckSum = ',i0) END IF # endif # ifdef DISTRIBUTE ! !----------------------------------------------------------------------- ! Broadcast IO error flag to all nodes. !----------------------------------------------------------------------- ! CALL mp_bcasti (ng, model, status) # endif ! RETURN END FUNCTION nf90_fwrite4d #endif #if defined PIO_LIB && defined DISTRIBUTE ! !*********************************************************************** FUNCTION pio_fwrite4d (ng, model, pioFile, ifield, & & pioVar, tindex, pioDesc, & & LBi, UBi, LBj, UBj, LBk, UBk, LBt, UBt, & & Ascl, & # ifdef MASKING & Amask, & # endif & Adat, SetFillVal, & & MinValue, MaxValue) RESULT (status) !*********************************************************************** ! USE mod_pio_netcdf ! USE distribute_mod, ONLY : mp_reduce # ifdef OUTPUT_STATS USE stats_mod, ONLY : stats_4dfld # endif ! ! Imported variable declarations. ! logical, intent(in), optional :: SetFillVal ! integer, intent(in) :: ng, model, tindex integer, intent(in) :: ifield integer, intent(in) :: LBi, UBi, LBj, UBj, LBk, UBk, LBt, UBt ! real(dp), intent(in) :: Ascl ! # ifdef ASSUMED_SHAPE # ifdef MASKING real(r8), intent(in) :: Amask(LBi:,LBj:) # endif real(r8), intent(in) :: Adat(LBi:,LBj:,LBk:,LBt:) # else # ifdef MASKING real(r8), intent(in) :: Amask(LBi:UBi,LBj:UBj) # endif real(r8), intent(in) :: Adat(LBi:UBi,LBj:UBj,LBk:UBk,LBt:UBt) # endif real(r8), intent(out), optional :: MinValue real(r8), intent(out), optional :: MaxValue ! TYPE (File_desc_t), intent(inout) :: pioFile TYPE (IO_Desc_t), intent(inout) :: pioDesc TYPE (My_VarDesc), intent(inout) :: pioVar ! ! Local variable declarations. ! logical :: LandFill, Lminmax logical, pointer :: Lwater(:,:,:,:) ! integer :: i, j, k, l, tile integer :: Imin, Imax, Jmin, Jmax integer :: Cgrid, dkind, ghost, gtype integer :: status integer, dimension(5) :: start, total ! real(r8), dimension(2) :: rbuffer real(r4), pointer :: Awrk4(:,:,:,:) real(r8), pointer :: Awrk8(:,:,:,:) # ifdef OUTPUT_STATS ! TYPE (T_STATS) :: Stats # endif ! character (len= 3), dimension(2) :: op_handle ! !----------------------------------------------------------------------- ! Set starting and ending indices to process. !----------------------------------------------------------------------- ! status=PIO_noerr ! ! Set first and last grid point according to staggered C-grid ! classification. Set loops offsets. ! ghost=0 dkind=pioVar%dkind gtype=pioVar%gtype ! SELECT CASE (ABS(gtype)) CASE (p3dvar) Cgrid=1 CASE (l3dvar, l4dvar, r3dvar) Cgrid=2 CASE (u3dvar) Cgrid=3 CASE (v3dvar) Cgrid=4 CASE DEFAULT Cgrid=2 END SELECT ! Imin=BOUNDS(ng)%Imin(Cgrid,ghost,MyRank) Imax=BOUNDS(ng)%Imax(Cgrid,ghost,MyRank) Jmin=BOUNDS(ng)%Jmin(Cgrid,ghost,MyRank) Jmax=BOUNDS(ng)%Jmax(Cgrid,ghost,MyRank) ! ! Set switch to compute minimum and maximum values. ! IF (PRESENT(MinValue)) THEN Lminmax=.TRUE. IF (.not.associated(Lwater)) THEN allocate ( Lwater(LBi:UBi,LBj:UBj,LBk:UBk,LBt:UBt) ) Lwater=.TRUE. END IF ELSE Lminmax=.FALSE. END IF # ifdef MASKING ! ! Set switch to replace land areas with fill value, spval. ! IF (PRESENT(SetFillVal)) THEN LandFill=SetFillVal ELSE LandFill=tindex.gt.0 END IF # endif # ifdef OUTPUT_STATS ! ! Initialize output variable statistics structure. ! Stats % checksum=0_i8b Stats % count=0 Stats % min=spval Stats % max=-spval Stats % avg=0.0_r8 Stats % rms=0.0_r8 # endif ! !----------------------------------------------------------------------- ! Write out data into NetCDF file. !----------------------------------------------------------------------- ! ! Allocate, initialize and load data into local array used for ! writing. Overwrite masked points with special value. ! IF (dkind.eq.PIO_double) THEN ! double precision IF (.not.associated(Awrk8)) THEN allocate ( Awrk8(LBi:UBi,LBj:UBj,LBk:UBk,LBt:UBt) ) Awrk8=0.0_r8 END IF ! DO l=LBt,UBt DO k=LBk,UBk DO j=Jmin,Jmax DO i=Imin,Imax Awrk8(i,j,k,l)=Adat(i,j,k,l)*Ascl # ifdef MASKING IF ((Amask(i,j).eq.0.0_r8).and.LandFill) THEN Awrk8(i,j,k,l)=spval IF (Lminmax) Lwater(i,j,k,l)=.FALSE. END IF # endif END DO END DO END DO END DO IF (Lminmax) THEN rbuffer(1)=MINVAL(Awrk8, MASK=Lwater) rbuffer(2)=MAXVAL(Awrk8, MASK=Lwater) END IF ELSE ! single precision IF (.not.associated(Awrk4)) THEN allocate ( Awrk4(LBi:UBi,LBj:UBj,LBk:UBk,LBt:UBt) ) Awrk4=0.0_r8 END IF ! DO l=LBt,UBt DO k=LBk,UBk DO j=Jmin,Jmax DO i=Imin,Imax Awrk4(i,j,k,l)=REAL(Adat(i,j,k,l)*Ascl, r4) # ifdef MASKING IF ((Amask(i,j).eq.0.0_r8).and.LandFill) THEN Awrk4(i,j,k,l)=REAL(spval, r4) END IF # endif END DO END DO END DO END DO IF (Lminmax) THEN rbuffer(1)=REAL(MINVAL(Awrk4, MASK=Lwater),r8) rbuffer(2)=REAL(MAXVAL(Awrk4, MASK=Lwater),r8) END IF END IF ! ! Set unlimited time dimension record to write, if any. ! IF (tindex.gt.0) THEN CALL PIO_setframe (pioFile, & & pioVar%vd, & & INT(tindex, kind=PIO_OFFSET_KIND)) END IF ! ! Write out data into NetCDF. ! IF (dkind.eq.PIO_double) THEN ! double precision CALL PIO_write_darray (pioFile, & & pioVar%vd, & & pioDesc, & & Awrk8(Imin:Imax,Jmin:Jmax, & & LBk:UBk,LBt:UBt), & & status) ELSE ! single precision CALL PIO_write_darray (pioFile, & & pioVar%vd, & & piodesc, & & Awrk4(Imin:Imax,Jmin:Jmax, & & LBk:UBk,LBt:UBt), & & status) END IF # ifdef OUTPUT_STATS ! !----------------------------------------------------------------------- ! Compute and report output field statistics. !----------------------------------------------------------------------- ! # ifdef DISTRIBUTE tile=MyRank # else tile=-1 # endif CALL stats_4dfld (ng, tile, model, gtype, Stats, & & LBi, UBi, LBj, UBj, LBk, UBk, LBt, UBt, & & Adat, & # ifdef MASKING & Fmask = Amask, & # endif & debug = .FALSE.) IF (OutThread) THEN WRITE (stdout,10) TRIM(Vname(1,ifield)), Stats%min, Stats%max, & & Stats%checksum 10 FORMAT (4x,'- ',a,':',t35,'Min = ',1p,e15.8,', Max = ', & & 1p,e15.8,', CheckSum = ',i0) END IF # endif ! !----------------------------------------------------------------------- ! If applicable, compute global minimum and maximum values. !----------------------------------------------------------------------- ! IF (Lminmax) THEN op_handle(1)='MIN' op_handle(2)='MAX' CALL mp_reduce (ng, model, 2, rbuffer, op_handle) MinValue=rbuffer(1) MaxValue=rbuffer(2) IF (associated(Lwater)) deallocate (Lwater) END IF ! ! Deallocate local array. ! IF (dkind.eq.PIO_double) THEN IF (associated(Awrk8)) deallocate (Awrk8) ELSE IF (associated(Awrk4)) deallocate (Awrk4) END IF ! RETURN END FUNCTION pio_fwrite4d #endif ! END MODULE nf_fwrite4d_mod