initstatsdump Subroutine

public subroutine initstatsdump()

Uses

  • proc~~initstatsdump~~UsesGraph proc~initstatsdump initstatsdump decomp_2d decomp_2d proc~initstatsdump->decomp_2d module~modfields modfields proc~initstatsdump->module~modfields module~modglobal modglobal proc~initstatsdump->module~modglobal module~modmpi modmpi proc~initstatsdump->module~modmpi module~modstat_nc modstat_nc proc~initstatsdump->module~modstat_nc module~modfields->decomp_2d mpi mpi module~modmpi->mpi module~modstat_nc->module~modmpi netcdf netcdf module~modstat_nc->netcdf

Arguments

None

Calls

proc~~initstatsdump~~CallsGraph proc~initstatsdump initstatsdump mpi_bcast mpi_bcast proc~initstatsdump->mpi_bcast proc~define_nc define_nc proc~initstatsdump->proc~define_nc proc~ncinfo ncinfo proc~initstatsdump->proc~ncinfo proc~open_nc open_nc proc~initstatsdump->proc~open_nc proc~writestat_dims_nc writestat_dims_nc proc~initstatsdump->proc~writestat_dims_nc zend zend proc~initstatsdump->zend zstart zstart proc~initstatsdump->zstart nf90_def_var nf90_def_var proc~define_nc->nf90_def_var nf90_enddef nf90_enddef proc~define_nc->nf90_enddef nf90_inq_dimid nf90_inq_dimid proc~define_nc->nf90_inq_dimid nf90_inq_varid nf90_inq_varid proc~define_nc->nf90_inq_varid nf90_put_att nf90_put_att proc~define_nc->nf90_put_att nf90_redef nf90_redef proc~define_nc->nf90_redef proc~nchandle_error nchandle_error proc~define_nc->proc~nchandle_error nf90_create nf90_create proc~open_nc->nf90_create nf90_def_dim nf90_def_dim proc~open_nc->nf90_def_dim proc~open_nc->nf90_def_var proc~open_nc->nf90_enddef nf90_get_var nf90_get_var proc~open_nc->nf90_get_var proc~open_nc->nf90_inq_dimid proc~open_nc->nf90_inq_varid nf90_inquire nf90_inquire proc~open_nc->nf90_inquire nf90_inquire_dimension nf90_inquire_dimension proc~open_nc->nf90_inquire_dimension nf90_open nf90_open proc~open_nc->nf90_open proc~open_nc->nf90_put_att nf90_sync nf90_sync proc~open_nc->nf90_sync proc~writestat_dims_nc->nf90_inq_varid proc~writestat_dims_nc->nf90_inquire_dimension nf90_put_var nf90_put_var proc~writestat_dims_nc->nf90_put_var nf90_strerror nf90_strerror proc~nchandle_error->nf90_strerror

Called by

proc~~initstatsdump~~CalledByGraph proc~initstatsdump initstatsdump program~dalesurban DALESURBAN program~dalesurban->proc~initstatsdump

Namelists

Namelist NAMSTATSDUMP


Variables

Name Type Default Description
lydump logical .false.
tsample real None
klow integer None
khigh integer None
tstatsdump real None
lytdump logical .false.
ltkedump logical .false.
lxydump logical .false.
lxytdump logical .false.
ltdump logical .false.
ltreedump logical .false.
lmintdump logical .false.

Source Code

  subroutine initstatsdump
    use modmpi,   only : my_real,mpierr,comm3d,mpi_logical,mpi_integer,mpi_character,cmyid,cmyidx,cmyidy
    use modglobal,only : imax,jmax,kmax,cexpnr,ifnamopt,fname_options,ib,ie,jb,je,kb,ke,ladaptive,btime,&
                         nsv,lkslicedump,lislicedump,ljslicedump,ltreedump,ib,ie,islice,islicerank,isliceloc,jslice,jslicerank,jsliceloc
    use modstat_nc,only: open_nc, define_nc,ncinfo,writestat_dims_nc
    use modfields, only : ncstaty,ncstatyt,ncstattke,ncstatxy,ncstatkslice,ncstatislice,ncstatjslice,ncstatxyt,ncstatt,ncstattr,ncstatmint
    use decomp_2d, only : zstart, zend
    implicit none
    integer :: ierr

    namelist/NAMSTATSDUMP/ &
         lydump,tsample,klow,khigh,tstatsdump,lytdump,ltkedump,lxydump,lxytdump,ltdump,ltreedump,lmintdump    ! maybe removed; NAMSTATSDUMP is not in use anymore

    allocate(ncstaty(nstaty,4))
    allocate(ncstatyt(nstatyt,4))
    allocate(ncstattke(nstattke,4))
    allocate(ncstatxy(nstatxy,4))
    allocate(ncstatkslice(nstatkslice,4))
    allocate(ncstatislice(nstatislice,4))
    allocate(ncstatjslice(nstatjslice,4))
    allocate(ncstatxyt(nstatxyt,4))
    allocate(ncstatt(nstatt,4))
    allocate(ncstattr(nstattr,4))
    allocate(ncstatmint(nstatmint,4))

    klow=kb
    khigh=ke

    if(myid==0)then
       open(ifnamopt,file=fname_options,status='old',iostat=ierr)
       read (ifnamopt,NAMSTATSDUMP,iostat=ierr)
       if (ierr > 0) then
          write(0, *) 'ERROR: Problem in namoptions NAMSTATSDUMP'
          write(0, *) 'iostat error: ', ierr
          stop 1
       endif
       !write(6 ,NAMSTATSDUMP)
       close(ifnamopt)
    end if

    call MPI_BCAST(klow        ,1,MPI_INTEGER,0,comm3d,ierr) !have to do this? just want nc for first CPU
    call MPI_BCAST(khigh       ,1,MPI_INTEGER,0,comm3d,ierr)
    call MPI_BCAST(nstatt      ,1,MPI_INTEGER,0,comm3d,ierr)
    call MPI_BCAST(nstatmint      ,1,MPI_INTEGER,0,comm3d,ierr)
    ! call MPI_BCAST(nstaty      ,1,MPI_INTEGER,0,comm3d,ierr)
    call MPI_BCAST(ncstatyt    ,80,MPI_CHARACTER,0,comm3d,mpierr)
    call MPI_BCAST(ncstaty     ,80,MPI_CHARACTER,0,comm3d,mpierr)
    call MPI_BCAST(ncstattke   ,80,MPI_CHARACTER,0,comm3d,mpierr)
    call MPI_BCAST(ncstatxy    ,80,MPI_CHARACTER,0,comm3d,mpierr)
    call MPI_BCAST(ncstatxyt   ,80,MPI_CHARACTER,0,comm3d,mpierr)
    call MPI_BCAST(ncstatt     ,80,MPI_CHARACTER,0,comm3d,mpierr)
    call MPI_BCAST(ncstatmint     ,80,MPI_CHARACTER,0,comm3d,mpierr)
    !call MPI_BCAST(ltdump      ,1,MPI_LOGICAL,0,comm3d,ierr)      ! maybe removed; unnecessary broadcast; this variable already broadcasted in modstartup
    !call MPI_BCAST(ltreedump   ,1,MPI_LOGICAL,0,comm3d,ierr)      ! maybe removed; unnecessary broadcast; this variable already broadcasted in modstartup
    call MPI_BCAST(lmintdump      ,1,MPI_LOGICAL,0,comm3d,ierr)

    !> Generate y-averaged NetCDF: ydump.xxx.nc
    if(lydump) then

      yname(7:9) = cexpnr
      call ncinfo(tncstaty(1,:),'time'      ,'Time'                         ,'s'      ,'time')
      call ncinfo(ncstaty( 1,:),'uy'        ,'Streamwise velocity'          ,'m/s'    ,'m0tt')
      call ncinfo(ncstaty( 2,:),'vy'        ,'Spanwise velocity'            ,'m/s'    ,'t0tt')
      call ncinfo(ncstaty( 3,:),'wy'        ,'Vertical velocity'            ,'m/s'    ,'t0mt')
      call ncinfo(ncstaty( 4,:),'thly'      ,'Temperature'                  ,'K'      ,'t0tt')
      call ncinfo(ncstaty( 5,:),'qty'       ,'Moisture'                     ,'kg/kg'  ,'t0tt')
      call ncinfo(ncstaty( 6,:),'sca1y'     ,'Scalar field 1'               ,'kg/m^3' ,'t0tt')
      call ncinfo(ncstaty( 7,:),'sca2y'     ,'Scalar field 2'               ,'kg/m^3' ,'t0tt')
      call ncinfo(ncstaty( 8,:),'sca3y'     ,'Scalar field 3'               ,'kg/m^3' ,'t0tt')
      call ncinfo(ncstaty( 9,:),'upwpy'     ,'Turbulent mom. flux'          ,'m^2/s^2','m0mt')
      call ncinfo(ncstaty(10,:),'wpthlpy'   ,'Turbulent heat flux'          ,'K m/s'  ,'t0mt')
      call ncinfo(ncstaty(11,:),'usgsy'     ,'SGS mom. flux'                ,'m^2/s^2','m0mt')
      call ncinfo(ncstaty(12,:),'thlsgsy'   ,'SGS heat flux'                ,'K m/s'  ,'t0mt')
      call ncinfo(ncstaty(13,:),'uwyik'     ,'Advective mom. flux'          ,'m^2/s^2','m0mt')
      call ncinfo(ncstaty(14,:),'wthlyk'    ,'Advective heat flux'          ,'K m/s'  ,'t0mt')

      if (myid==0) then
        call open_nc(yname, ncidy, nrecy, n1=imax, n3=khigh-klow+1)
        if (nrecy==0) then
          call define_nc( ncidy, 1, tncstaty)
          call writestat_dims_nc(ncidy)
        end if
      call define_nc( ncidy, nstaty, ncstaty)
      endif !myid==0
    endif

    !> Generate time and y averaged NetCDF: ytdump.xxx.nc
    if (lytdump) then

      ytname(8:10) = cexpnr
      call ncinfo(tncstatyt(1,:),'time'       ,'Sampling time'             ,'s'       ,'time')
      call ncinfo(ncstatyt( 1,:),'uyt'        ,'Streamwise velocity'       ,'m/s'     ,'m0tt')
      call ncinfo(ncstatyt( 2,:),'vyt'        ,'Spanwise velocity'         ,'m/s'     ,'t0tt')
      call ncinfo(ncstatyt( 3,:),'wyt'        ,'Vertical velocity'         ,'m/s'     ,'t0mt')
      call ncinfo(ncstatyt( 4,:),'thlyt'      ,'Temperature'               ,'K'       ,'t0tt')
      call ncinfo(ncstatyt( 5,:),'qtyt'       ,'Moisture'                  ,'kg/kg'   ,'t0tt')
      call ncinfo(ncstatyt( 6,:),'sca1yt'     ,'Scalar field 1'            ,'kg/m^3'  ,'t0tt')
      call ncinfo(ncstatyt( 7,:),'sca2yt'     ,'Scalar field 2'            ,'kg/m^3'  ,'t0tt')
      call ncinfo(ncstatyt( 8,:),'sca3yt'     ,'Scalar field 3'            ,'kg/m^3'  ,'t0tt')

      call ncinfo(ncstatyt( 9,:),'upwpyt'     ,'Turbulent mom. flux'       ,'m^2/s^2' ,'m0mt')
      call ncinfo(ncstatyt( 10,:),'wpthlpyt'   ,'Turbulent heat flux'       ,'K m/s'   ,'t0mt')
      call ncinfo(ncstatyt( 11,:),'wpqtpyt'   ,'Turbulent moisture flux'   ,'kg/kg m/s','t0mt')
      call ncinfo(ncstatyt( 12,:),'wpsca1tpyt','Turbulent scalar flux'     ,'M m/s'   ,'t0mt')
      call ncinfo(ncstatyt( 13,:),'wpsca2tpyt','Turbulent scalar flux'     ,'M m/s'   ,'t0mt')
      call ncinfo(ncstatyt( 14,:),'wpsca3tpyt','Turbulent scalar flux'     ,'M m/s'   ,'t0mt')

      call ncinfo(ncstatyt( 15,:),'uwyt'      ,'Kinematic mom. flux'       ,'m^2/s^2' ,'m0mt')
      call ncinfo(ncstatyt( 16,:),'wthlyt'    ,'Kinematic heat flux'       ,'K m/s'   ,'t0mt')
      call ncinfo(ncstatyt( 17,:),'wqtyt'     ,'Kinematic moisture flux'   ,'K m/s'   ,'t0mt')
      call ncinfo(ncstatyt( 18,:),'wsca1yt'   ,'Kinematic scalar flux'     ,'K m/s'   ,'t0mt')
      call ncinfo(ncstatyt( 19,:),'wsca2yt'   ,'Kinematic scalar flux'     ,'K m/s'   ,'t0mt')
      call ncinfo(ncstatyt( 20,:),'wsca3yt'   ,'Kinematic scalar flux'     ,'K m/s'   ,'t0mt')

      call ncinfo(ncstatyt( 21,:),'upupyt'     ,'mom. variance'            ,'m^2/s^2' ,'m0tt')
      call ncinfo(ncstatyt( 22,:),'wpwpyt'     ,'mom. variance'            ,'m^2/s^2' ,'t0mt')
      call ncinfo(ncstatyt( 23,:),'thlpthlpyt' ,'temp. variance'           ,'K^2'     ,'t0tt')
      call ncinfo(ncstatyt( 24,:),'qtpqtpyt'   ,'moisture. variance'       ,'kg^2/kg^2','t0tt')
      call ncinfo(ncstatyt( 25,:),'sca1tpsca1pyt','scalar. variance'       ,'M^2'     ,'t0tt')
      call ncinfo(ncstatyt( 26,:),'sca2tpsca2pyt','scalar. variance'       ,'M^2'     ,'t0tt')
      call ncinfo(ncstatyt( 27,:),'sca3tpsca3pyt','scalar. variance'       ,'M^2'     ,'t0tt')

      call ncinfo(ncstatyt( 28,:),'usgsyt'    ,'SGS mom. flux'             ,'m^2/s^2' ,'m0mt')
      call ncinfo(ncstatyt( 29,:),'wsgsyt'    ,'SGS mom. flux'             ,'m^2/s^2' ,'t0mt')
      call ncinfo(ncstatyt( 30,:),'thlsgsyt'  ,'SGS heat flux'             ,'K m/s'   ,'t0mt')
      call ncinfo(ncstatyt( 31,:),'qtsgsyt'   ,'SGS moisture flux'         ,'kg/kg m/s','t0mt')
      call ncinfo(ncstatyt( 32,:),'sca1sgsyt' ,'SGS scalar flux'           ,'M m/s'   ,'t0mt')
      call ncinfo(ncstatyt( 33,:),'sca2sgsyt' ,'SGS scalar flux'           ,'M m/s'   ,'t0mt')
      call ncinfo(ncstatyt( 34,:),'sca3sgsyt' ,'SGS scalar flux'           ,'M m/s'   ,'t0mt')

      if (myid==0) then
        call open_nc(ytname, ncidyt, nrecyt, n1=imax, n3=khigh-klow+1)
        if (nrecyt==0) then
          call define_nc( ncidyt, 1, tncstatyt)
          call writestat_dims_nc(ncidyt)
        end if
        call define_nc( ncidyt, nstatyt, ncstatyt)
      endif !myid==0
    endif

    !> Generate y and x averaged NetCDF: xydump.xxx.nc
    if (lxydump) then

      xyname(8:10) = cexpnr
      call ncinfo(tncstatxy(1,:),'time'    ,'Time'                        ,'s'      ,'time')
      call ncinfo(ncstatxy( 1,:),'uxy'     ,'Streamwise velocity'         ,'m/s'    ,'tt'  )
      call ncinfo(ncstatxy( 2,:),'vxy'     ,'Spanwise velocity'           ,'m/s'    ,'tt'  )
      call ncinfo(ncstatxy( 3,:),'wxy'     ,'Vertical velocity'           ,'m/s'    ,'mt'  )
      call ncinfo(ncstatxy( 4,:),'thlxy'   ,'Temperature'                 ,'K'      ,'tt'  )
      call ncinfo(ncstatxy( 5,:),'qtxy'    ,'Moisture'                    ,'kg/kg'  ,'tt'  )
      call ncinfo(ncstatxy( 6,:),'pxy'     ,'Pressure'                    ,'kgm/s^2','tt'  )
      call ncinfo(ncstatxy( 7,:),'upwpxy'  ,'Mom. flux'                   ,'m^2/s^2','mt'  )
      call ncinfo(ncstatxy( 8,:),'wpthlpxy','Heat flux'                   ,'Km/s'   ,'mt'  )
      call ncinfo(ncstatxy( 9,:),'vpwpxy'  ,'Mom. flux'                   ,'Km/s'   ,'mt'  )
      call ncinfo(ncstatxy(10,:),'usgsxy'  ,'SGS mom. flux'               ,'m^2/s^2','mt'  )
      call ncinfo(ncstatxy(11,:),'thlsgsxy','SGS heat flux'               ,'Km/s'   ,'mt'  )
      call ncinfo(ncstatxy(12,:),'vsgsxy'  ,'SGS mom. flux'               ,'m^2/s^2','mt'  )
      call ncinfo(ncstatxy(13,:),'uwxyik'  ,'Advective mom. flux'         ,'m^2/s^2','mt'  )
      call ncinfo(ncstatxy(14,:),'wthlxy'  ,'Advective heat flux'         ,'K m/s'  ,'mt'  )
      call ncinfo(ncstatxy(15,:),'vwxy'    ,'Advective mom. flux'         ,'m^2/s^2','mt'  )
      if (myid==0) then
        call open_nc(xyname, ncidxy, nrecxy, n3=khigh-klow+1)
        if (nrecxy==0) then
          call define_nc( ncidxy, 1, tncstatxy)
          call writestat_dims_nc(ncidxy)
        end if
        call define_nc( ncidxy, nstatxy, ncstatxy)
      end if
    end if

    !> Generate time, y and x averaged NetCDF: xytdump.xxx.nc
    if (lxytdump) then

      xytname(9:11) = cexpnr
      call ncinfo(tncstatxyt(1,:),'time'      ,'Time'                        ,'s'      ,'time')
      call ncinfo(ncstatxyt( 1,:),'uxyt'       ,'Streamwise velocity'         ,'m/s'    ,'tt'  )
      call ncinfo(ncstatxyt( 2,:),'vxyt'       ,'Spanwise velocity'           ,'m/s'    ,'tt'  )
      call ncinfo(ncstatxyt( 3,:),'wxyt'       ,'Vertical velocity'           ,'m/s'    ,'mt'  )
      call ncinfo(ncstatxyt( 4,:),'thlxyt'     ,'Temperature'                 ,'K'      ,'tt'  )
      call ncinfo(ncstatxyt( 5,:),'qtxyt'      ,'Moisture'                    ,'kg/kg'  ,'tt'  )
      call ncinfo(ncstatxyt( 6,:),'pxyt'       ,'Pressure'                    ,'kgm/s^2','tt'  )
      call ncinfo(ncstatxyt( 7,:),'upwpxyt'    ,'Turbulent mom. flux'         ,'m^2/s^2','mt'  )
      call ncinfo(ncstatxyt( 8,:),'wpthlpxyt'  ,'Turbulent heat flux'         ,'K m/s'  ,'mt'  )
      call ncinfo(ncstatxyt( 9,:),'vpwpxyt'    ,'Turbulent mom. flux'         ,'m^2/s^2','mt'  )
      call ncinfo(ncstatxyt( 10,:),'upvpxyt'   ,'Turbulent mom. flux'         ,'m^2/s^2','mt'  )
      call ncinfo(ncstatxyt( 11,:),'uwxyt'     ,'Kinematic mom. flux'         ,'m^2/s^2','mt'  )
      call ncinfo(ncstatxyt( 12,:),'wthlxyt'   ,'Kinematic heat flux'         ,'K m/s'  ,'mt'  )
      call ncinfo(ncstatxyt( 13,:),'uvxyt'     ,'Kinematic mom. flux'         ,'m^2/s^2','mt'  )
      call ncinfo(ncstatxyt( 14,:),'vwxyt'     ,'Kinematic mom. flux'         ,'m^2/s^2','mt'  )
      call ncinfo(ncstatxyt( 15,:),'wwxyt'     ,'Kinematic mom. flux'         ,'m^2/s^2','mt'  )
      call ncinfo(ncstatxyt( 16,:),'usgsxyt'   ,'SGS mom. flux'               ,'m^2/s^2','mt'  )
      call ncinfo(ncstatxyt( 17,:),'thlsgsxyt' ,'SGS heat flux'               ,'K m/s'  ,'mt'  )
      call ncinfo(ncstatxyt( 18,:),'vsgsxyt'   ,'SGS mom. flux'               ,'K m/s'  ,'mt'  )
      call ncinfo(ncstatxyt( 19,:),'thlpthlptxy','Temp. variance'             ,'K^2'    ,'tt'  )
      call ncinfo(ncstatxyt( 20,:),'upuptxyc'  ,'u variance'                  ,'m^2/s^2','tt'  )
      call ncinfo(ncstatxyt( 21,:),'vpvptxyc'  ,'v variance'                  ,'m^2/s^2','tt'  )
      call ncinfo(ncstatxyt( 22,:),'wpwptxyc'  ,'w variance'                  ,'m^2/s^2','tt'  )
      call ncinfo(ncstatxyt( 23,:),'tketxyc'   ,'tke'                         ,'m^2/s^2','tt'  )

      if (myid==0) then
        call open_nc(xytname, ncidxyt, nrecxyt, n3=khigh-klow+1)
        if (nrecxyt==0) then
          call define_nc( ncidxyt, 1, tncstatxyt)
          call writestat_dims_nc(ncidxyt)
        end if
        call define_nc( ncidxyt, nstatxyt, ncstatxyt)
      end if
    end if

    !> Generate time averaged NetCDF: tdump.xxx.nc
    if (ltdump) then

      tname(7:9) = cmyidx
      tname(11:13) = cmyidy
      tname(15:17) = cexpnr
      call ncinfo(tncstatt(1,:),'time'      ,'Time'                        ,'s'      ,'time')
      call ncinfo(ncstatt( 1,:),'ut'        ,'Streamwise velocity'         ,'m/s'    ,'mttt'  )
      call ncinfo(ncstatt( 2,:),'vt'        ,'Spanwise velocity'           ,'m/s'    ,'tmtt'  )
      call ncinfo(ncstatt( 3,:),'wt'        ,'Vertical velocity'           ,'m/s'    ,'ttmt'  )
      call ncinfo(ncstatt( 4,:),'thlt'      ,'Temperature'                 ,'K'      ,'tttt'  )
      call ncinfo(ncstatt( 5,:),'qtt'       ,'Moisture'                    ,'kg/kg'  ,'tttt'  )
      call ncinfo(ncstatt( 6,:),'pt'        ,'Pressure'                    ,'kgm/s^2','tttt'  )
      call ncinfo(ncstatt( 7,:),'sca1t'     ,'Concentration field 1'       ,'g/m^3'  ,'tttt'  )
      call ncinfo(ncstatt( 8,:),'sca2t'     ,'Concentration field 2'       ,'g/m^3'  ,'tttt'  )
      call ncinfo(ncstatt( 9,:),'sca3t'     ,'Concentration field 3'       ,'g/m^3'  ,'tttt'  )
      call ncinfo(ncstatt(10,:),'sca4t'     ,'Concentration field 4'       ,'g/m^3'  ,'tttt'  )
      call ncinfo(ncstatt(11,:),'PSS'       ,'PSS defect'                  ,'gm/s'   ,'tttt'  )

      call ncinfo(ncstatt(12,:),'upwpt'     ,'Turbulent momentum flux'     ,'m^2/s^2','mtmt'  )
      call ncinfo(ncstatt(13,:),'vpwpt'     ,'Turbulent momentum flux'     ,'m^2/s^2','tmmt'  )
      call ncinfo(ncstatt(14,:),'upvpt'     ,'Turbulent momentum flux'     ,'m^2/s^2','mmtt'  )
      call ncinfo(ncstatt(15,:),'wpthlpt'   ,'Turbulent heat flux'         ,'K m/s'  ,'ttmt'  )
      call ncinfo(ncstatt(16,:),'wpsca1pt'  ,'Turbulent flux 1'            ,'gm/s'   ,'ttmt'  )
      call ncinfo(ncstatt(17,:),'wpsca2pt'  ,'Turbulent flux 2'            ,'gm/s'   ,'ttmt'  )
      call ncinfo(ncstatt(18,:),'wpsca3pt' ,'Turbulent flux 3'            ,'gm/s'   ,'ttmt'  )
      call ncinfo(ncstatt(19,:),'wpsca4pt'  ,'Turbulent flux 4'            ,'gm/s'   ,'ttmt'  )

      call ncinfo(ncstatt(20,:),'thlpthlpt','Temperature variance'        ,'K^2'    ,'tttt'  )
      call ncinfo(ncstatt(21,:),'upuptc'   ,'u variance'                  ,'m^2/s^2','tttt'  )
      call ncinfo(ncstatt(22,:),'vpvptc'   ,'v variance'                  ,'m^2/s^2','tttt'  )
      call ncinfo(ncstatt(23,:),'wpwptc'   ,'w variance'                  ,'m^2/s^2','tttt'  )
      call ncinfo(ncstatt(24,:),'tketc'    ,'TKE'                         ,'m^2/s^2','tttt'  )
      call ncinfo(ncstatt(25,:),'sca1psca1pt','Concentration variance 1'  ,'g^2/m^6','tttt'  )
      call ncinfo(ncstatt(26,:),'sca2psca2pt','Concentration variance 2'  ,'g^2/m^6','tttt'  )
      call ncinfo(ncstatt(27,:),'sca3psca3pt','Concentration variance 3'  ,'g^2/m^6','tttt'  )
      call ncinfo(ncstatt(28,:),'sca4psca4pt','Concentration variance 4'  ,'g^2/m^6','tttt'  )

      call ncinfo(ncstatt(29,:),'sv1sgs'   ,'SGS flux 1'                  ,'gm/s'   ,'ttmt'  )
      call ncinfo(ncstatt(30,:),'sv2sgs'   ,'SGS flux 2'                  ,'gm/s'   ,'ttmt'  )
      call ncinfo(ncstatt(31,:),'sv3sgs'   ,'SGS flux 3'                  ,'gm/s'   ,'ttmt'  )
      call ncinfo(ncstatt(32,:),'sv4sgs'   ,'SGS flux 4'                  ,'gm/s'   ,'ttmt'  )

      ! call ncinfo(ncstatt(33,:),'sca1t_max','Max concentration field 1'   ,'g/m^3'  ,'tttt'  )
      ! call ncinfo(ncstatt(34,:),'sca2t_max','Max concentration field 2'   ,'g/m^3'  ,'tttt'  )
      ! call ncinfo(ncstatt(35,:),'sca3t_max','Max concentration field 3'   ,'g/m^3'  ,'tttt'  )
      ! call ncinfo(ncstatt(36,:),'sca4t_max','Max concentration field 4'   ,'g/m^3'  ,'tttt'  )

!      if (myid==0) then
        call open_nc(tname, ncidt, nrect, n1=imax, n2=jmax, n3=khigh-klow+1)
        if (nrect==0) then
          call define_nc( ncidt, 1, tncstatt)
          call writestat_dims_nc(ncidt)
        end if
        call define_nc( ncidt, nstatt, ncstatt)
!      end if
    end if

    if (lmintdump) then

      mintname(10:12) = cmyidx
      mintname(14:16) = cmyidy
      mintname(18:20) = cexpnr
      call ncinfo(tncstatmint(1,:),'time'      ,'Time'                        ,'s'      ,'time')
      call ncinfo(ncstatmint( 1,:),'ut'        ,'Streamwise velocity'         ,'m/s'    ,'mttt'  )
      call ncinfo(ncstatmint( 2,:),'vt'        ,'Spanwise velocity'           ,'m/s'    ,'tmtt'  )
      call ncinfo(ncstatmint( 3,:),'wt'        ,'Vertical velocity'           ,'m/s'    ,'ttmt'  )
      call ncinfo(ncstatmint( 4,:),'thlt'      ,'Temperature'                 ,'K'      ,'tttt'  )
      call ncinfo(ncstatmint( 5,:),'qtt'       ,'Moisture'                    ,'kg/kg'  ,'tttt'  )
      call ncinfo(ncstatmint( 6,:),'pt'        ,'Pressure'                    ,'kgm/s^2','tttt'  )

    !      if (myid==0) then
        call open_nc(mintname, ncidmint, nrecmint, n1=imax, n2=jmax, n3=khigh-klow+1)
        if (nrecmint==0) then
          call define_nc( ncidmint, 1, tncstatmint)
          call writestat_dims_nc(ncidmint)
        end if
        call define_nc( ncidmint, nstatmint, ncstatmint)
    !      end if
    end if

    !> Generate time averaged NetCDF: treedump.xxx.nc
    if (ltreedump) then

      trname(10:12) = cmyidx
      trname(14:16) = cmyidy
      trname(18:20) = cexpnr
      call ncinfo(tncstattr(1,:),'time'      ,'Time'                        ,'s'      ,'time')
      call ncinfo(ncstattr( 1,:),'tr_u'      ,'Drag in x'                   ,'m/s^2'  ,'tttt'  )
      call ncinfo(ncstattr( 2,:),'tr_v'      ,'Drag in y'                   ,'m/s^2'  ,'tttt'  )
      call ncinfo(ncstattr( 3,:),'tr_w'      ,'Drag in z'                   ,'m/s^2'  ,'ttmt'  )
      call ncinfo(ncstattr( 4,:),'tr_thl'    ,'Temp source/ sink'           ,'K/s'    ,'tttt'  )
      call ncinfo(ncstattr( 5,:),'tr_qt'     ,'Moisture source sink'        ,'1/s'    ,'tttt'  )
      call ncinfo(ncstattr( 6,:),'tr_qtR'    ,'Moisture source sink'        ,'1/s'    ,'tttt'  )
      call ncinfo(ncstattr( 7,:),'tr_qtA'    ,'Moisture source sink'        ,'1/s'    ,'tttt'  )
      call ncinfo(ncstattr( 8,:),'tr_sv1'    ,'Scalar source sink'          ,'kg/m^3s','tttt'  )
      call ncinfo(ncstattr( 9,:),'tr_sv2'    ,'Scalar source sink'          ,'kg/m^3s','tttt'  )
      call ncinfo(ncstattr(10,:),'tr_omega'  ,'Decoupling factor'           ,'-'      ,'tttt'  )

!      if (myid==0) then
        call open_nc(trname, ncidtr, nrectr, n1=imax, n2=jmax, n3=khigh-klow+1)
        if (nrectr==0) then
          call define_nc( ncidtr, 1, tncstattr)
          call writestat_dims_nc(ncidtr)
        end if
        call define_nc( ncidtr, nstattr, ncstattr)
!      end if
    end if

    !> Generate time, y and x averaged NetCDF for tke budget: tkedump.xxx.nc
    if (ltkedump) then

      tkename(9:11) = cexpnr
      call ncinfo(tncstattke(1,:),'time' ,'Time'                                 ,'s'       ,'time')
      call ncinfo(ncstattke( 1,:),'p_b'  ,'p_bant production or consumption term', 'm^2/s^3','tt'  )
      call ncinfo(ncstattke( 2,:),'t_p'  ,'total viscous transport (?)'          , 'm^2/s^3','tt'  )
      call ncinfo(ncstattke( 3,:),'adv'  ,'Advection by mean wind'               , 'm^2/s^3','tt'  )
      call ncinfo(ncstattke( 4,:),'t_t'  ,'Total turb???'                        , 'm^2/s^3','tt'  )
      call ncinfo(ncstattke( 5,:),'t_sgs','total SGS  term'                      , 'm^2/s^3','tt'  )
      call ncinfo(ncstattke( 6,:),'p_t'  ,'Shear production term'                , 'm^2/s^3','tt'  )
      call ncinfo(ncstattke( 7,:),'t_v'  ,'Resolved viscous dissipation term'    , 'm^2/s^3','tt'  )
      call ncinfo(ncstattke( 8,:),'d_sgs','SGS dissipation term'                 , 'm^2/s^3','tt'  )

      if (myid==0) then
        call open_nc(tkename, ncidtke, nrectke, n3=khigh-klow+1)
        if (nrectke==0) then
          call define_nc( ncidtke, 1, tncstattke)
          call writestat_dims_nc(ncidtke)
        end if
        call define_nc( ncidtke, nstattke, ncstattke)
      endif !myid==0

    endif

    !> Generate sliced NetCDF: slicedump.xxx.xxx.nc
    if (lkslicedump) then

      kslicename(12:14) = cmyidx
      kslicename(16:18) = cmyidy
      kslicename(20:22) = cexpnr

      call ncinfo(tncstatkslice(1,:),'time'     ,'Time'   ,'s'   ,'time')
      call ncinfo(ncstatkslice( 1,:),'u_kslice'     ,'Streamwise velocity at kslice', '-', 'mt0t')
      call ncinfo(ncstatkslice( 2,:),'v_kslice'     ,'Spanwise velocity at kslice', '-', 'tm0t')
      call ncinfo(ncstatkslice( 3,:),'w_kslice'     ,'Vertical velocity at kslice', '-', 'tt0t')
      call ncinfo(ncstatkslice( 4,:),'thl_kslice'   ,'Potential temperature at kslice', '-', 'tt0t')
      call ncinfo(ncstatkslice( 5,:),'qt_kslice'    ,'Specific humidity at kslice', '-', 'tt0t')

      call open_nc(kslicename, ncidkslice, nreckslice, n1=imax, n2=jmax)

      if (nreckslice==0) then
        call define_nc( ncidkslice, 1, tncstatkslice)
        call writestat_dims_nc(ncidkslice)
      end if

      call define_nc( ncidkslice, nstatkslice, ncstatkslice)

    end if

    if (lislicedump) then

      islicename(12:14) = cmyidx
      islicename(16:18) = cmyidy
      islicename(20:22) = cexpnr

      call ncinfo(tncstatislice(1,:),'time'     ,'Time'   ,'s'   ,'time')
      call ncinfo(ncstatislice( 1,:),'u_islice'     ,'Streamwise velocity at islice', '-', '0ttt')
      call ncinfo(ncstatislice( 2,:),'v_islice'     ,'Spanwise velocity at islice', '-', '0mtt')
      call ncinfo(ncstatislice( 3,:),'w_islice'     ,'Vertical velocity at islice', '-', '0tmt')
      call ncinfo(ncstatislice( 4,:),'thl_islice'   ,'Potential temperature at islice', '-', '0ttt')
      call ncinfo(ncstatislice( 5,:),'qt_islice'    ,'Specific humidity at islice', '-', '0ttt')

      if ((islice >= zstart(1)) .and. (islice <= zend(1))) then
        islicerank = .true.
        isliceloc = islice - zstart(1) + 1
      else
        islicerank = .false.
      end if

      if (islicerank) then
        call open_nc(islicename, ncidislice, nrecislice, n2=jmax, n3=khigh-klow+1)
        if (nrecislice==0) then
          call define_nc( ncidislice, 1, tncstatislice)
          call writestat_dims_nc(ncidislice)
        end if
        call define_nc( ncidislice, nstatislice, ncstatislice)
      end if

    end if

    if (ljslicedump) then

      jslicename(12:14) = cmyidx
      jslicename(16:18) = cmyidy
      jslicename(20:22) = cexpnr

      call ncinfo(tncstatjslice(1,:),'time'     ,'Time'   ,'s'   ,'time')
      call ncinfo(ncstatjslice( 1,:),'u_jslice'     ,'Streamwise velocity at jslice', '-', 'm0tt')
      call ncinfo(ncstatjslice( 2,:),'v_jslice'     ,'Spanwise velocity at jslice', '-', 't0tt')
      call ncinfo(ncstatjslice( 3,:),'w_jslice'     ,'Vertical velocity at jslice', '-', 't0mt')
      call ncinfo(ncstatjslice( 4,:),'thl_jslice'   ,'Potential temperature at jslice', '-', 't0tt')
      call ncinfo(ncstatjslice( 5,:),'qt_jslice'    ,'Specific humidity at jslice', '-', 't0tt')

      if ((jslice >= zstart(2)) .and. (jslice <= zend(2))) then
        jslicerank = .true.
        jsliceloc = jslice - zstart(2) + 1
      else
        jslicerank = .false.
      end if

      if (jslicerank) then
         call open_nc(jslicename, ncidjslice, nrecjslice, n1=imax, n3=khigh-klow+1)
         if (nrecjslice==0) then
            call define_nc( ncidjslice, 1, tncstatjslice)
            call writestat_dims_nc(ncidjslice)
         end if
         call define_nc( ncidjslice, nstatjslice, ncstatjslice)
      end if

    end if

    !> Set times to zero so works for warm starts... could have issues with warmstarts here...
    tsamplep = 0.
    tstatsdumpp = 0.

  end subroutine initstatsdump