writestat_dims_nc Subroutine

public subroutine writestat_dims_nc(ncid)

Uses

  • proc~~writestat_dims_nc~~UsesGraph proc~writestat_dims_nc writestat_dims_nc module~modglobal modglobal proc~writestat_dims_nc->module~modglobal module~modmpi modmpi proc~writestat_dims_nc->module~modmpi mpi mpi module~modmpi->mpi

Arguments

Type IntentOptional Attributes Name
integer, intent(in) :: ncid

Calls

proc~~writestat_dims_nc~~CallsGraph proc~writestat_dims_nc writestat_dims_nc nf90_inq_varid nf90_inq_varid proc~writestat_dims_nc->nf90_inq_varid nf90_inquire_dimension nf90_inquire_dimension proc~writestat_dims_nc->nf90_inquire_dimension nf90_put_var nf90_put_var proc~writestat_dims_nc->nf90_put_var

Called by

proc~~writestat_dims_nc~~CalledByGraph proc~writestat_dims_nc writestat_dims_nc proc~initeb initEB proc~initeb->proc~writestat_dims_nc proc~initfielddump initfielddump proc~initfielddump->proc~writestat_dims_nc proc~initibm initibm proc~initibm->proc~writestat_dims_nc proc~initstatsdump initstatsdump proc~initstatsdump->proc~writestat_dims_nc program~dalesurban DALESURBAN program~dalesurban->proc~initeb program~dalesurban->proc~initfielddump program~dalesurban->proc~initibm program~dalesurban->proc~initstatsdump

Source Code

  subroutine writestat_dims_nc(ncid)
    use modglobal, only : xf,xh,yf,yh,dy,zf,zh,jmax,imax,kb,dxh
    use modmpi, only : myidx, myidy
    implicit none
    integer, intent(in) :: ncid
    integer             :: i=0,iret,length,varid
    integer :: dx

    dx = dxh(1) ! Assume equidistant grid
    !write(*,*) 'writestat_dims_nc'
    iret = nf90_inq_varid(ncid, 'xt', VarID)
    if (iret==0) iret=nf90_inquire_dimension(ncid, xtID, len=length)
    !if (iret==0) iret = nf90_put_var(ncid, varID, xf(1:length),(/1/))
    if (iret==0) iret = nf90_put_var(ncid, varID, (/(xf(i+imax*myidx),i=1,length)/),(/1/))
    !if (iret==0) iret = nf90_put_var(ncid, varID, (/(dx*(0.5+i)+myidx*imax*dx,i=0,length-1)/),(/1/))

    iret = nf90_inq_varid(ncid, 'xm', VarID)
    if (iret==0) iret=nf90_inquire_dimension(ncid, xmID, len=length)
    !if (iret==0) iret = nf90_put_var(ncid, varID, xh(1:length),(/1/))
    if (iret==0) iret = nf90_put_var(ncid, varID, (/(xh(i+imax*myidx),i=1,length)/),(/1/))
    !if (iret==0) iret = nf90_put_var(ncid, varID, (/(dx*i+myidx*imax*dx,i=0,length-1)/),(/1/))

    iret = nf90_inq_varid(ncid, 'yt', VarID)
    if (iret==0) iret=nf90_inquire_dimension(ncid, ytID, len=length)
    !if (iret==0) iret = nf90_put_var(ncid, varID, (/(dy*(0.5+i)+myidy*jmax*dy,i=0,length-1)/),(/1/))
    if (iret==0) iret = nf90_put_var(ncid, varID, (/(yf(i+jmax*myidy),i=1,length)/),(/1/))
    iret = nf90_inq_varid(ncid, 'ym', VarID)
    if (iret==0) iret=nf90_inquire_dimension(ncid, ymID, len=length)
    !if (iret==0) iret = nf90_put_var(ncid, varID, (/(dy*i+myidy*jmax*dy,i=0,length-1)/),(/1/))
    if (iret==0) iret = nf90_put_var(ncid, varID, (/(yh(i+jmax*myidy),i=1,length)/),(/1/))

    iret = nf90_inq_varid(ncid, 'zt', VarID)
    if (iret==0) iret=nf90_inquire_dimension(ncid,ztID, len=length)
    if (iret==0) iret = nf90_put_var(ncid, varID, zf(1:length),(/1/))  !ils13, 29.06.2017 zf starts at 0, not at 1
    iret = nf90_inq_varid(ncid, 'zm', VarID)
    if (iret==0) iret=nf90_inquire_dimension(ncid, zmID, len=length)
    if (iret==0) iret = nf90_put_var(ncid, varID, zh(1:length),(/1/))  !same for zh
    !if (isurf==1) then
      !iret = nf90_inq_varid(ncid, 'zts', VarID)
      !if (iret==0) iret = nf90_inquire_dimension(ncid, ztsID, len=length)
      !if (iret==0) iret = nf90_put_var(ncid, varID, zsoilc(1:length),(/1/))
    !end if

  end subroutine writestat_dims_nc