#include <misc.h>
#include <preproc.h>
module inicFileMod 1,5
!-----------------------------------------------------------------------
!BOP
!
! !MODULE: inicFileMod
!
! !DESCRIPTION:
! Read CLM initial data netCDF files (used only in CAM perpetual mode)
!
! !USES:
use shr_kind_mod
, only : r8 => shr_kind_r8
use spmdMod
, only : masterproc
use abortutils
, only : endrun
use clm_varctl
, only : iulog
use ncdio
! !PUBLIC TYPES:
implicit none
save
!
! !PUBLIC MEMBER FUNCTIONS:
public :: inicfile_perp
!
! !REVISION HISTORY:
! Jan/2004 Mariana Vertenstein: Creation
! Jan/2010 Li Xu: Modified to correct ncd_iolocal and snow_fraction
!
!EOP
!-----------------------------------------------------------------------
contains
!------------------------------------------------------------------------
!BOP
!
! !IROUTINE: inicfile_perp
!
! !INTERFACE:
subroutine inicfile_perp 1,38
!
! !DESCRIPTION:
! Read perpetual initial data fields
!
! !USES:
use clmtype
use clm_varctl
, only : finidat
use clm_varpar
, only : nlevsno, nlevgrnd, nlevsoi, nlevlak
use clm_varcon
, only : denice, denh2o, zlnd, isturb
use fileutils
, only : getfil
use decompMod
, only : get_proc_bounds, get_proc_global
use fileutils
, only : getfil
use shr_sys_mod
, only : shr_sys_flush
!
! !ARGUMENTS:
implicit none
!
! !REVISION HISTORY:
!
!
! !LOCAL VARIABLES:
!EOP
integer :: ncid ! netCDF dataset id
integer :: j,c,l ! indices
integer :: begp, endp ! per-proc beginning and ending pft indices
integer :: begc, endc ! per-proc beginning and ending column indices
integer :: begl, endl ! per-proc beginning and ending landunit indices
integer :: begg, endg ! per-proc gridcell ending gridcell indices
integer :: numg ! total number of gridcells across all processors
integer :: numl ! total number of landunits across all processors
integer :: numc ! total number of column across all processors
integer :: nump ! total number of pfts across all processors
logical :: readvar ! determine if variable is on initial file
type(gridcell_type), pointer :: gptr ! pointer to gridcell derived subtype
type(landunit_type), pointer :: lptr ! pointer to landunit derived subtype
type(column_type) , pointer :: cptr ! pointer to column derived subtype
type(pft_type) , pointer :: pptr ! pointer to pft derived subtype
integer , pointer :: clandunit(:) ! landunit index associated with each column
integer , pointer :: itypelun(:) ! landunit type
character(len=256), save :: loc_fni ! local finidat name
logical, save :: opened_finidat = .false. ! true => finidat was opened for read
character(len= 32) :: subname='inicfile_perp' ! subroutine name
!------------------------------------------------------------------------
! Assign local pointers to derived subtypes components (landunit-level)
itypelun => clm3%g%l%itype
! Assign local pointers to derived subtypes components (column-level)
clandunit => clm3%g%l%c%landunit
! Set pointers into derived type
gptr => clm3%g
lptr => clm3%g%l
cptr => clm3%g%l%c
pptr => clm3%g%l%c%p
! Determine necessary processor subgrid bounds
call get_proc_bounds
(begg, endg, begl, endl, begc, endc, begp, endp)
call get_proc_global
(numg, numl, numc, nump)
! Read initial dataset
if (masterproc) then
if (.not. opened_finidat) then
call getfil
(finidat, loc_fni, 0)
call check_ret
(nf_open(loc_fni, nf_nowrite, ncid), subname)
write(iulog,*)trim(subname),': opened netcdf file ',loc_fni
call shr_sys_flush
(iulog)
call check_dim
(ncid, 'gridcell', numg)
call check_dim
(ncid, 'landunit', numl)
call check_dim
(ncid, 'column' , numc)
call check_dim
(ncid, 'pft' , nump)
call check_dim
(ncid, 'levsno' , nlevsno)
call check_dim
(ncid, 'levgrnd' , nlevgrnd)
call check_dim
(ncid, 'levlak' , nlevlak)
opened_finidat = .true.
else
call check_ret
(nf_open(loc_fni, nf_nowrite, ncid), subname)
end if
end if
call ncd_iolocal
(varname='ZSNO', data=cptr%cps%z, dim1name=namec, dim2name='levsno', &
lowerb2=-nlevsno+1, upperb2=0, ncid=ncid, flag='read', readvar=readvar)
if (.not. readvar) call endrun
()
call ncd_iolocal
(varname='DZSNO', data=cptr%cps%dz, dim1name=namec, dim2name='levsno', &
lowerb2=-nlevsno+1, upperb2=0, ncid=ncid, flag='read', readvar=readvar)
if (.not. readvar) call endrun
()
call ncd_iolocal
(varname='ZISNO', data=cptr%cps%zi, dim1name=namec, dim2name='levsno', &
lowerb2=-nlevsno, upperb2=-1, ncid=ncid, flag='read', readvar=readvar)
if (.not. readvar) call endrun
()
call ncd_iolocal
(varname='H2OSNO', data=cptr%cws%h2osno, dim1name=namec, &
ncid=ncid, flag='read', readvar=readvar)
if (.not. readvar) call endrun
()
call ncd_iolocal
(varname='SNOWDP', data=cptr%cps%snowdp, dim1name=namec, &
ncid=ncid, flag='read', readvar=readvar)
if (.not. readvar) call endrun
()
call ncd_iolocal
(varname='SNLSNO', data=cptr%cps%snl, dim1name=namec, &
ncid=ncid, flag='read', readvar=readvar)
if (.not. readvar) call endrun
()
call ncd_iolocal
(varname='H2OSOI_LIQ', data=cptr%cws%h2osoi_liq, dim1name=namec, dim2name='levtot', &
ncid=ncid, flag='read', readvar=readvar)
if (.not. readvar) call endrun
()
call ncd_iolocal
(varname='H2OSOI_ICE', data=cptr%cws%h2osoi_ice, dim1name=namec, dim2name='levtot', &
ncid=ncid, flag='read', readvar=readvar)
if (.not. readvar) call endrun
()
if (masterproc) then
call check_ret
(nf_close(ncid), subname)
end if
! Determine volumetric soil water
do j = 1,nlevsoi
do c = begc,endc
l = cptr%landunit(c)
if (.not. lptr%lakpoi(l)) then
cptr%cws%h2osoi_vol(c,j) = cptr%cws%h2osoi_liq(c,j)/(cptr%cps%dz(c,j)*denh2o) &
+cptr%cws%h2osoi_ice(c,j)/(cptr%cps%dz(c,j)*denice)
end if
end do
end do
! ============================================================================
! ! Fraction of soil covered by snow (Z.-L. Yang U. Texas)
! ============================================================================
do c = begc,endc
l = clandunit(c)
if (itypelun(l) == isturb) then
! Urban landunit use Bonan 1996 (LSM Technical Note)
cptr%cps%frac_sno(c) = min( cptr%cps%snowdp(c)/0.05_r8, 1._r8)
else
! snow cover fraction in Niu et al. 2007
cptr%cps%frac_sno(c) = 0.0_r8
if ( cptr%cps%snowdp(c) > 0.0_r8 ) then
cptr%cps%frac_sno(c) = tanh(cptr%cps%snowdp(c)/(2.5_r8*zlnd* &
(min(800._r8,cptr%cws%h2osno(c)/cptr%cps%snowdp(c))/100._r8)**1._r8) )
endif
end if
end do
end subroutine inicfile_perp
end module inicFileMod