#include <misc.h>
#include <params.h>
module carbonscales 2,4
!-----------------------------------------------------------------------
!
! Purposes:
! Read in scaling dataset.
! Provide (linearly interpolated) scaling when requested
!
! Public routines:
! get_carbonscale
! init_scale
!
!-----------------------------------------------------------------------
use shr_kind_mod
, only: r8 => shr_kind_r8
use time_manager
, only: get_curr_calday, get_curr_date
use pmgrid
, only: masterproc
use abortutils
, only: endrun
implicit none
public init_scale ! read file data
public get_carbonscale ! return scale at current time step
private
save
real(r8),allocatable :: fscale(:)
real(r8),allocatable :: days(:) ! date converted to days
integer :: indextopreviousdate
integer :: numberoffiledata
contains
subroutine init_scale() 1,19
!-----------------------------------------------------------------------
! read scale file into module variables
!-----------------------------------------------------------------------
use ioFileMod
, only: getfil
use filenames
, only: bndtvcarbonscale
#if ( defined SPMD )
use mpishorthand
#endif
character(len=256) :: locfn
integer :: scalenid = -1 ! netcdf id for scale file
integer :: datedimid, scaleid, dateid
integer,allocatable :: date(:)
integer :: i
integer :: dimids(1)
integer :: yr, mon, day,sec
real(r8) :: calday ! calendar day of current timestep
real(r8) :: caldayloc ! calendar day of current timestep
if (masterproc) then
!
! find and open file; abort if fail (getfil(,,0)).
!
call getfil
(bndtvcarbonscale, locfn, 0)
call wrap_open
(locfn, 0, scalenid)
write(6,*)'init_scale: reading scale dataset'
!
! Get and check dimension info
!
call wrap_inq_dimid
( scalenid, 'date', datedimid )
call wrap_inq_dimlen
( scalenid, datedimid, numberoffiledata)
allocate( date(numberoffiledata) )
allocate( fscale(numberoffiledata) )
call wrap_inq_varid
( scalenid, 'carbonscale', scaleid )
call wrap_inq_varid
( scalenid, 'date', dateid )
call wrap_inq_vardimid
(scalenid, scaleid, dimids)
!
! read in Population and dates
!
call wrap_get_var_realx
(scalenid, scaleid, fscale)
call wrap_get_var_int
(scalenid, dateid, date)
#if ( defined SPMD )
call mpibcast
(numberoffiledata, 1, mpir8, 0, mpicom)
else
call mpibcast
(numberoffiledata, 1, mpir8, 0, mpicom)
allocate( date(numberoffiledata) )
allocate( fscale(numberoffiledata) )
#endif
endif
!
! broadcast scale and date to nodes
!
#if ( defined SPMD )
call mpibcast
(date, numberoffiledata, mpiint, 0, mpicom)
call mpibcast
(fscale, numberoffiledata, mpir8, 0, mpicom)
#endif
!
! convert date to days
!
allocate( days(numberoffiledata) )
do i = 1, numberoffiledata
call bnddyi
(date(i), 0, days(i))
days(i) = days(i) + (date(i) / 10000) * 365
enddo
deallocate( date )
calday = get_curr_calday ()
call get_curr_date
(yr, mon, day, sec)
caldayloc = calday + yr*365
if (caldayloc < days(1)) then
write(6,*) 'init_scale: calday, caldayloc, yr', calday, caldayloc, yr
write(6,*) 'init_scale: days(1)',days(1)
call endrun
('init_scale: calendar manager returns date previous to scale data')
endif
if (caldayloc >= days(numberoffiledata)) then
write(6,*) 'init_scale: calday, caldayloc, yr', calday, caldayloc, yr
write(6,*) 'init_scale: days(numberoffiledata)',days(numberoffiledata)
call endrun
('init_scale: calendar manager returns date after all scale data')
endif
!
! find data which bounds to the left
!
do indextopreviousdate = 1, numberoffiledata
if(caldayloc < days(indextopreviousdate)) exit
enddo
indextopreviousdate = indextopreviousdate - 1
write(6,*) 'init_scale: caldayloc',caldayloc,'dateindex,date',indextopreviousdate,days
end subroutine init_scale
subroutine get_carbonscale(scaleattime) 1,3
!-----------------------------------------------------------------------
! compute scale for time indicated by time_manager
!-----------------------------------------------------------------------
real(r8), intent(out) :: scaleattime
integer :: yr, mon, day,sec
real(r8) :: calday ! calendar day of current timestep
real(r8) :: caldayloc ! calendar day of current timestep
real(r8) :: deltat, fact2, fact1 ! time interpolation factors
real(r8) :: cdayp, cdaym
calday = get_curr_calday ()
call get_curr_date
(yr, mon, day, sec)
caldayloc = calday + yr*365
if (caldayloc >= days(indextopreviousdate+1)) then
indextopreviousdate = indextopreviousdate + 1
if (indextopreviousdate + 1 > numberoffiledata) then
write(6,*) 'get_carbonscale: date exceeds scale data', caldayloc, days(numberoffiledata)
call endrun
endif
endif
if ( (caldayloc < days(indextopreviousdate)) .or. caldayloc >= days(indextopreviousdate+1) ) then
write(6,*) 'simulation day, bounding days',caldayloc, days(indextopreviousdate), days(indextopreviousdate+1)
call endrun
('get_carbonscale: problem with data in scale file')
endif
deltat = days(indextopreviousdate+1) - days(indextopreviousdate)
fact1 = (days(indextopreviousdate+1) - caldayloc) / deltat
fact2 = (caldayloc - days(indextopreviousdate )) / deltat
scaleattime = fscale(indextopreviousdate ) * fact1 + &
fscale(indextopreviousdate+1) * fact2
end subroutine get_carbonscale
end module carbonscales