#include <misc.h> 
module acbnd 3,4
!-----------------------------------------------------------------------
! Purpose:
! This code does time interpolation for 3D boundary data in a netCDF
! file. Assumptions on the data in the netCDF file are:
! 1. Coordinates are ordered (lon,lev,lat,time)
! 2. The time coordinate is in days, and the data is assumed to be periodic
! annually.
!
! Author: B. Eaton
! Modified 11 Jan 2000: to use F90 interface and allow 4 byte reals
!-----------------------------------------------------------------------
use shr_kind_mod, only: r8 => shr_kind_r8
use ppgrid, only: pcols, pver, begchunk, endchunk
use pmgrid, only: masterproc, plon, plat
use abortutils, only: endrun
!nf90 use netcdf
implicit none
save
private
public :: &
acbndini, &! initialize acbnd module
acbndint, &! interpolate acbnd data to requested date/time
acbndget ! return latitude slice data at current date/time
#include <netcdf.inc>
real(r8), allocatable, dimension(:) :: &
time ! time coordinate (calander days + frac)
real(r8), allocatable, dimension(:,:,:,:,:) :: &
datin ! input data (pcols,pver,begchunk:endchunk,2,nvar)
real(r8), allocatable, dimension(:,:,:,:) :: &
dat ! interpolated data (pcols,pver,begchunk:endchunk,nvar)
character(len=128), allocatable, dimension(:) :: &
varnam ! names of variables to get from netCDF file
integer :: &
nvar, &! number of variables requested
ncid, &! ID for netCDF file
nrec, &! number of records (time samples)
lotim, &! time(lotim) .le. current time
hitim, &! current time .lt. time(hitim)
loin, &! index into input data array containing time(lotim) data
hiin, &! index into input data array containing time(hitim) data
start(4), &! start vector for netCDF hyperslabs
count(4) ! count vector for netCDF hyperslabs
!##############################################################################
contains
!##############################################################################
subroutine acbndini( ncfile, calday, xnvar, ncnam ) 1,22
!-----------------------------------------------------------------------
! Purpose:
! Open netCDF file containing annual cycle data. Initialize
! arrays with the data to be interpolated to the current time.
!
! It is assumed that the time coordinate is increasing and represents
! calendar days.
!
! Author: B. Eaton
!-----------------------------------------------------------------------
use error_messages, only: alloc_err, handle_ncerr, handle_err
use phys_grid, only: scatter_field_to_chunk
#if ( defined SPMD )
use mpishorthand
#endif
!-----------------------------------------------------------------------
#include <netcdf.inc>
!-----------------------------------------------------------------------
real(r8), intent(in) ::&
calday ! current time in calendar days + fraction.
integer, intent(in) ::&
xnvar ! number of requested variables in the netCDF file
character(len=*), intent(in) ::&
ncfile, &! netCDF data file
ncnam(xnvar) ! names of requested variables in the netCDF file
! Local variables:
integer ::&
i, &
istat, &
recid, &! record ID
vid ! variable ID
real(r8), allocatable :: datlatlon(:,:,:) ! used for netCDF input
!-----------------------------------------------------------------------
if ( xnvar .eq. 0 ) return
nvar = xnvar
! Allocate space for data.
allocate( varnam(nvar), stat=istat )
call alloc_err( istat, 'acbndini', 'varnam', nvar )
allocate( datin(pcols,pver,begchunk:endchunk,2,nvar), stat=istat )
call alloc_err( istat, 'acbndini', 'datin', &
pcols*pver*((endchunk-begchunk)+1)*2*nvar )
allocate( dat(pcols,pver,begchunk:endchunk,nvar), stat=istat )
call alloc_err( istat, 'acbndini', 'dat', &
pcols*pver*((endchunk-begchunk)+1)*nvar )
allocate( datlatlon(plon,pver,plat), stat=istat )
call alloc_err( istat, 'acbndini', 'datlatlon', &
plon*pver*plat )
datin(:,:,:,:,:) = 0._r8
dat(:,:,:,:) = 0._r8
do i = 1, nvar
varnam(i) = ncnam(i)
end do
if (masterproc) then
start(1) = 1
start(2) = 1
start(3) = 1
start(4) = 1
count(1) = plon
count(2) = pver
count(3) = plat
count(4) = 1
! Open file.
!nf90 call handle_ncerr( nf90_open( trim( ncfile ), NF_NOWRITE, ncid ), &
!nf90 'acbndini: error opening file '//trim( ncfile ) )
call handle_ncerr( nf_open( trim( ncfile ), NF_NOWRITE, ncid ), &
'acbndini: error opening file '//trim( ncfile ) )
! get the record id
!nf90 call handle_ncerr( nf90_inquire( ncid, unlimiteddimid=recid), &
!nf90 'acbndini: no record variables ' )
call handle_ncerr( nf_inq_unlimdim( ncid, recid), &
'acbndini: no record variables ' )
! Get size of unlimited dimension.
!nf90 call handle_ncerr( nf90_inquire_dimension( ncid, recid, len=nrec ), 'acbndini: ' )
call handle_ncerr( nf_inq_dimlen( ncid, recid, nrec ), 'acbndini: ' )
endif ! end of masterproc
#if ( defined SPMD )
! broadcast nrec to nodes
call mpibcast (nrec, 1, mpiint, 0, mpicom)
#endif
! Allocate space for time coordinate data.
allocate( time(nrec), stat=istat )
call alloc_err( istat, 'acbndini', 'time', nrec )
if (masterproc) then
! Get time coordinate.
!nf90 call handle_ncerr( nf90_inq_varid( ncid, 'time', vid ), &
!nf90 'acbndini: cannot find time coordinate variable' )
!nf90 call handle_ncerr( nf90_get_var( ncid, vid, time ), &
!nf90 'acbndini: error getting time coordinate data' )
call handle_ncerr( nf_inq_varid( ncid, 'time', vid ), &
'acbndini: cannot find time coordinate variable' )
call handle_ncerr( nf_get_var_double( ncid, vid, time ), &
'acbndini: error getting time coordinate data' )
endif ! end of masterproc
#if ( defined SPMD )
! broadcast time to nodes
call mpibcast (time, nrec, mpir8, 0, mpicom)
#endif
! Make sure the time coordinate looks like calander day, and is
! increasing.
call chktime( time, nrec )
! Find indices for time samples that bound the current time.
call findplb( time, nrec, calday, lotim )
hitim = mod( lotim, nrec ) + 1
! Read data.
loin = 1
hiin = 2
do i = 1, nvar
if (masterproc) then
!nf90 call handle_ncerr( nf90_inq_varid( ncid, varnam(i), vid ), &
!nf90 'acbndini: cannot find variable '//varnam(i) )
call handle_ncerr( nf_inq_varid( ncid, varnam(i), vid ), &
'acbndini: cannot find variable '//varnam(i) )
start(4) = lotim
!nf90 call handle_ncerr( nf90_get_var( ncid, vid, datlatlon, start, count ), &
!nf90 'acbndini: cannot read data for '//varnam(i) )
call handle_ncerr( nf_get_vara_double( ncid, vid, start, count, datlatlon ), &
'acbndini: cannot read data for '//varnam(i) )
endif ! end of masterproc
! scatter data to chunked data structures
call scatter_field_to_chunk (1, pver, 1, plon, datlatlon, &
datin(1,1,begchunk,loin,i))
if (masterproc) then
start(4) = hitim
!nf90 call handle_ncerr( nf90_get_var( ncid, vid, datlatlon, start, count ), &
!nf90 'acbndini: cannot read data for '//varnam(i) )
call handle_ncerr( nf_get_vara_double( ncid, vid, start, count, datlatlon ), &
'acbndini: cannot read data for '//varnam(i) )
endif ! end of masterproc
! scatter data to chunked data structures
call scatter_field_to_chunk (1, pver, 1, plon, datlatlon, &
datin(1,1,begchunk,hiin,i))
end do
if (masterproc) then
write(6,*)'acbndini: calendar day = ',calday,' : read data for days ', &
time(lotim), ' and ',time(hitim)
endif ! end of masterproc
deallocate( datlatlon, stat=istat )
call handle_err( istat, &
'ERROR deallocating memory for datlatlon in routine acbndini')
end subroutine acbndini
!#######################################################################
subroutine acbndint( calday ) 1,12
!-----------------------------------------------------------------------
! Purpose:
! Interpolate annual cycle surface flux data to the current time. Read
! in new time samples of the input data as necessary.
!
! Author: B. Eaton
!-----------------------------------------------------------------------
use error_messages, only: alloc_err, handle_ncerr, handle_err
use phys_grid, only: scatter_field_to_chunk
!-----------------------------------------------------------------------
#include <netcdf.inc>
!-----------------------------------------------------------------------
real(r8), intent(in) ::&
calday ! current time in calendar days + fraction.
! Local variables:
integer ::&
i, &
oldhitim, &
vid, ret
real(r8) ::&
dt, dt1, tint
! TBH: Get rid of all the replicated code for I/O and scattering!
real(r8), allocatable :: datlatlon(:,:,:) ! used for netCDF input
integer :: numpts
integer :: istat
!-----------------------------------------------------------------------
if ( nvar .eq. 0 ) return
! Check to see if model time is still bounded by dataset times.
oldhitim = hitim
call findplb( time, nrec, calday, lotim )
hitim = mod( lotim, nrec ) + 1
if ( hitim .ne. oldhitim ) then
! Read in new hitim data. Replace old lotim data.
loin = hiin
hiin = mod( loin, 2 ) + 1
endif
if ( hitim .ne. oldhitim ) then
allocate( datlatlon(plon,pver,plat), stat=istat )
call alloc_err( istat, 'acbndini', 'datlatlon', &
plon*pver*plat )
start(4) = hitim
do i = 1, nvar
if (masterproc) then
!nf90 call handle_ncerr( nf90_inq_varid( ncid, varnam(i), vid ), &
!nf90 'acbndint: cannot find variable '//varnam(i) )
!nf90 call handle_ncerr( nf90_get_var( ncid, vid, datlatlon, start, count ), &
!nf90 'acbndint: cannot read data for '//varnam(i) )
call handle_ncerr( nf_inq_varid( ncid, varnam(i), vid ), &
'acbndint: cannot find variable '//varnam(i) )
call handle_ncerr( nf_get_vara_double( ncid, vid, start, count, datlatlon ), &
'acbndint: cannot read data for '//varnam(i) )
endif ! end of masterproc
! scatter data to chunked data structures
call scatter_field_to_chunk (1, pver, 1, plon, datlatlon, &
datin(1,1,begchunk,hiin,i))
end do
if (masterproc) then
write(6,*)'acbndint: read data for day ',time(hitim)
endif ! end of masterproc
if ( lotim .ne. oldhitim ) then
! Read in new lotim data. Replace old hitim data.
start(4) = lotim
do i = 1, nvar
if (masterproc) then
!nf90 call handle_ncerr( nf90_inq_varid( ncid, varnam(i), vid ), &
!nf90 'acbndint: cannot find variable '//varnam(i) )
!nf90 call handle_ncerr( nf90_get_var( ncid, vid, datlatlon, start, count ), &
!nf90 'acbndint: cannot read data for '//varnam(i) )
call handle_ncerr( nf_inq_varid( ncid, varnam(i), vid ), &
'acbndint: cannot find variable '//varnam(i) )
call handle_ncerr( nf_get_vara_double( ncid, vid, start, count, datlatlon ), &
'acbndint: cannot read data for '//varnam(i) )
endif ! end of masterproc
! scatter data to chunked data structures
call scatter_field_to_chunk (1, pver, 1, plon, datlatlon, &
datin(1,1,begchunk,loin,i))
end do
if (masterproc) then
write(6,*)'acbndint: read data for day ',time(lotim)
endif ! end of masterproc
end if
deallocate( datlatlon, stat=istat )
call handle_err( istat, &
'ERROR deallocating memory for datlatlon in routine acbndini')
end if
! Linear interpolation... Start by computing the number of days between
! the lower and upper bounds, and days between
! the model time and lower bound.
if( time(hitim) .lt. time(lotim) )then
dt = 365. - time(lotim) + time(hitim)
if( calday .le. time(hitim) )then
dt1 = 365. - time(lotim) + calday
else
dt1 = calday - time(lotim)
end if
else
dt = time(hitim) - time(lotim)
dt1 = calday - time(lotim)
end if
tint = dt1/dt
numpts = size(dat) / nvar
do i = 1, nvar
call linintp( numpts, 0._r8, 1._r8, tint, &
datin(1,1,begchunk,loin,i), &
datin(1,1,begchunk,hiin,i), &
dat(1,1,begchunk,i) )
end do
end subroutine acbndint
!#######################################################################
subroutine acbndget( name, ncol, lchnk, slice ) 5,2
!-----------------------------------------------------------------------
! Purpose:
! Return annual cycle data for the requested latitudes and variable.
!
! Author: B. Eaton
!-----------------------------------------------------------------------
character(len=*), intent(in) :: name ! requested variable
integer, intent(in) :: ncol ! max index of 1st dimension
integer, intent(in) :: lchnk ! chunk index
real(r8), intent(out) :: slice(pcols,pver) ! requested data
! Local variables:
integer ::&
i, k, n
!-----------------------------------------------------------------------
if ( nvar .eq. 0 ) call endrun ('ACBNDGET')
do n = 1, nvar
if ( name .eq. trim(varnam(n)) ) then
do k = 1, pver
do i = 1, ncol
slice(i,k) = dat(i,k,lchnk,n)
end do
end do
return
end if
end do
call endrun ('ACBNDGET: '//name//' not found')
end subroutine acbndget
!#######################################################################
end module acbnd