#include <misc.h> 
#include <params.h>
#if ( defined SCAM )
#include <max.h>
#endif
module prescribed_aerosols 9,7
!-----------------------------------------------------------------------
!
! Purposes:
! read, store, interpolate, and return fields
! of aerosols to CAM. The initialization
! file (mass.nc) is assumed to be a monthly climatology
! of aerosols from MATCH (on a sigma pressure
! coordinate system).
! also provide a "background" aerosol field to correct
! for any deficiencies in the physical parameterizations
! This fields is a "tuning" parameter.
! Public methods:
! (1) - initialization
! read aerosol masses from external file
! also pressure coordinates
! convert from monthly average values to mid-month values
! (2) - interpolation (time and vertical)
! interpolate onto pressure levels of CAM
! interpolate to time step of CAM
! return mmr's of aerosols
! provide background aerosol based on "tauback"
! (3) - diagnostics
! write out various diagnostic fields
!
! Calling Heirarchy
! aerosol_initialize (public)
! add_aer_diag_fields
! get_aerosol (public)
! vert_interpolate
! background
! scale_aerosols
! aerosol_diagnostics (public)
! aerosol_indirect (public)
! get_rf_scales (public)
! get_int_scales (public)
!
!-----------------------------------------------------------------------
use shr_kind_mod, only: r8 => shr_kind_r8
use ppgrid, only: pcols, pver, pverp, begchunk, endchunk
use pmgrid, only: plon, plat, plev, plevp, masterproc
use phys_grid, only: get_ncols_p, scatter_field_to_chunk
use time_manager, only: get_curr_calday
use infnan, only: inf, bigint
use abortutils, only: endrun
#if ( defined SCAM )
use scamMod, only :initlonidx,initlatidx
use getnetcdfdata
#endif
implicit none
! naer_all is total number of species
! naer is number of species in climatology
! naer_all = naer + 1 (background "species") + 1 (volcanic)
integer, public, parameter :: naer_all = 12
integer, private, parameter :: naer = 10
! indices to aerosol array (species portion)
integer, public, parameter :: &
idxSUL = 1, &
idxSSLT = 2, &
idxOCPHO = 7, &
idxBCPHO = 8, &
idxOCPHI = 9, &
idxBCPHI = 10, &
idxBG = 11, &
idxVOLC = 12
! indices to sections of array that represent
! groups of aerosols
integer, public, parameter :: &
idxDUSTfirst = 3, &
numDUST = 4, &
idxCARBONfirst = 7, &
numCARBON = 4
! names of aerosols are they are represented in
! the climatology file and for the purposes
! of outputting mmr's. Appended '_V' indicates field has been vertically summed.
character(len=8), public, parameter :: aerosol_name(naer_all) = &
(/"MSUL_V "&
,"MSSLT_V "&
,"MDUST1_V"&
,"MDUST2_V"&
,"MDUST3_V"&
,"MDUST4_V"&
,"MOCPHO_V"&
,"MBCPHO_V"&
,"MOCPHI_V"&
,"MBCPHI_V"&
,"Backgrnd"&
,"MVOLC "/)
!
! default values for namelist variables
!
! compute carbons scaling from population
character(len=16), public :: scenario_carbon_scale = 'FIXED' ! 'FIXED' or 'RAMPED'
character(len=16), public :: scenario_prescribed_sulfur = 'FIXED' ! 'FIXED' or 'RAMPED'
character(len=16), public :: prescribed_sulfur = 'direct' ! 'off', 'passive' or 'direct'
integer, public :: rampyear_prescribed_sulfur = bigint ! this is the only acceptable value
! compute radiative forcing per aerosol
logical, public :: radforce = .false.
! compute stratospheric volcanic aerosol fields and forcing
logical, public :: strat_volcanic = .false.
! portion of each species group to use in computation
! of relative radiative forcing.
real(r8), public :: sulscl_rf = 0._r8 !
real(r8), public :: carscl_rf = 0._r8
real(r8), public :: ssltscl_rf = 0._r8
real(r8), public :: dustscl_rf = 0._r8
real(r8), public :: bgscl_rf = 0._r8
real(r8), public :: volcscl_rf = 0._r8
! "background" aerosol species mmr.
real(r8), public :: tauback = 0._r8
! portion of each species group to use in computation
! of aerosol forcing in driving the climate
real(r8), public :: sulscl = 1._r8
real(r8), public :: carscl = 1._r8
real(r8), public :: ssltscl = 1._r8
real(r8), public :: dustscl = 1._r8
real(r8), public :: volcscl = 1._r8
private
save
integer :: aernid = -1 ! netcdf id for aerosol file (init to invalid)
integer :: species_id(naer) = -1 ! netcdf_id of each aerosol species (init to invalid)
integer :: Mpsid ! netcdf id for MATCH PS
integer :: nm = 1 ! index to prv month in array. init to 1 and toggle between 1 and 2
integer :: np = 2 ! index to nxt month in array. init to 2 and toggle between 1 and 2
integer :: mo_nxt = bigint ! index to nxt month in file
real(r8) :: cdaym = inf ! calendar day of prv month
real(r8) :: cdayp = inf ! calendar day of next month
real(r8) :: Mid(12) ! Days into year for mid month date
data Mid/16.5, 46.0, 75.5, 106.0, 136.5, 167.0, 197.5, 228.5, 259.0, 289.5, 320.0, 350.5 /
public aerosol_initialize ! read from file, interpolate onto horiz grid
public get_aerosol ! interpolate onto pressure levels, and time
public aerosol_diagnostics ! write out various data about aerosols
public aerosol_indirect ! compute change in effective liqw radius from sulfate
public get_rf_scales ! compute scale factors for forcing computation
public get_int_scales ! compute scale factors for interactive comps
public aerint ! read next month info
!
! values read from file and temporary values used for interpolation
!
! AEROSOLc is:
! Cumulative Mass at midpoint of each month
! on CAM's horizontal grid (col)
! on MATCH's levels (lev)
! AEROSOLc(month_ind, col_ind, Match_lev_ind, chunk_ind, species_ind)
!
integer, parameter :: paerlev = 28 ! number of levels for aerosol fields (MUST = naerlev)
integer naerlev ! size of level dimension in MATCH data
real(r8) :: M_hybi(paerlev+1) ! MATCH hybi
real(r8) :: M_ps(plon,plat) ! surface pressure from MATCH file
real(r8), allocatable :: AEROSOLc(:,:,:,:,:) ! Aerosol cumulative mass from MATCH
real(r8), allocatable :: M_ps_cam_col(:,:,:) ! PS from MATCH on Cam Columns
real(r8), allocatable :: Match_ps_chunk(:,:) ! surface pressure array section at a particular time step
integer :: mxaerl ! Maximum level of background aerosol
contains
subroutine aerosol_initialize 1,49
!------------------------------------------------------------------
! Reads in:
! file from which to read aerosol Masses on CAM grid. Currently
! assumed to be MATCH ncep runs, averaged by month.
! NOTE (Data have been externally interpolated onto CAM grid
! and backsolved to provide Mid-month values)
!
! Populates:
! module variables:
! AEROSOLc(pcols,paerlev+1,begchunk:endchunk,naer,2))
! AEROSOLc( column_index
! , level_index (match levels)
! , chunk_index
! , species_index
! , month = 1:2 )
! M_hybi(level_index = Lev_MATCH) = pressure at mid-level.
! M_ps_cam_col(column,chunk,month) ! PS from MATCH on Cam Columns
!
! Method:
! read data from file
! allocate memory for storage of aerosol data on CAM horizontal grid
! distribute data to remote nodes
! populates the module variables
!
!------------------------------------------------------------------
use ioFileMod, only: getfil
use filenames, only: bndtvaer
use volcanicmass, only: volcanic_initialize
use carbonscales, only: init_scale
#if ( defined SPMD )
use mpishorthand
#endif
include 'netcdf.inc'
! local variables
integer :: naerlon
integer :: naerlat
integer :: naerlev
integer dateid ! netcdf id for date variable
integer secid ! netcdf id for seconds variable
integer londimid ! netcdf id for longitude dimension
integer latdimid ! netcdf id for latitude dimension
integer levdimid ! netcdf id for level dimension
integer lonid ! netcdf id for longitude variable
integer timesiz ! number of time samples (=12) in netcdf file
integer latid ! netcdf id for latitude variable
integer Mhybiid ! netcdf id for MATCH hybi
integer timeid ! netcdf id for time variable
integer dimids(nf_max_var_dims) ! variable shape
integer :: start(4) ! start vector for netcdf calls
integer :: kount(4) ! count vector for netcdf calls
integer mo ! month index
integer m ! constituent index
integer :: n ! loop index
integer :: i,j,k ! spatial indices
integer :: date_aer(12) ! Date on aerosol dataset (YYYYMMDD)
integer :: attnum ! attribute number
#if ( defined SCAM )
real(r8) :: coldata(paerlev) ! aerosol field read in from dataset
integer :: ret
#endif
integer mo_prv ! index to previous month
character(len=256) :: locfn ! netcdf local filename to open
!
! aerosol_data will be read in from the aerosol boundary dataset, then scattered to chunks
! after filling in the bottom level with zeros
!
real(r8) :: aerosol_data(plon,plat,paerlev) ! aerosol field read in from dataset
real(r8) :: aerosol_field(plon,paerlev+1,plat) ! aerosol field to be scattered
real(r8) :: caldayloc ! calendar day of current timestep
call t_startf ('aerosol_init')
call print_memusage ('Start aerosol_initialize')
if(scenario_carbon_scale.eq."RAMPED") then
write(6,*)"aerosol_init: carbon ramping overrides carscl"
else if(scenario_carbon_scale.eq."FIXED") then
write(6,*) "aerosol_init: using namelist value for carscl, or default"
else
call endrun ('scenario_carbon_scale must be either RAMPED or FIXED')
endif
if(scenario_prescribed_sulfur.eq."RAMPED") then
call endrun ('aerosol_init: using sulfur from external historical sequence is not implemented')
else if(scenario_prescribed_sulfur.eq."FIXED") then
write(6,*) "aerosol_init: Using sulfate from cyclic climatology"
if ( rampyear_prescribed_sulfur .ne. bigint ) then
call endrun ('aerosol_init: rampyear_prescribed_sulfur set but there is only one year in the climatology')
endif
else
call endrun ('scenario_prescribed_sulfur must be either RAMPED or FIXED')
endif
if(prescribed_sulfur == 'off') then
call endrun ('aerosol_init: prescribed sulfur = off is not implemented, try passive')
else if(.not.(prescribed_sulfur == 'direct' .or. prescribed_sulfur == 'passive')) then
write(6,*)
call endrun ('prescribed sulfur must be one of off, direct, or passive')
endif
!
! Allocate memory for dynamic arrays local to this module
!
allocate (AEROSOLc(pcols,paerlev+1,begchunk:endchunk,naer,2))
allocate (M_ps_cam_col(pcols,begchunk:endchunk,2))
allocate (Match_ps_chunk(pcols,begchunk:endchunk))
! TBH: HACK to avoid use of uninitialized values when ncols < pcols
AEROSOLc(:,:,:,:,:) = 0.
M_ps_cam_col(:,:,:) = 0.
Match_ps_chunk(:,:) = 0.
! TBH: END HACK
!
! Add fields to master list for output
!
call add_aer_diag_fields ()
!
! compute mxaerl
!
call setmxaerl ()
if (masterproc) then
!
! find and open file; abort if fail (getfil(,,0)).
!
call getfil (bndtvaer, locfn, 0)
call wrap_open (locfn, 0, aernid)
write(6,*)'aerosol_initialize: reading aerosol dataset.... If this seems to be taking too'
write(6,*)'long, perhaps the dataset is being read from an nfs-mounted filesystem.'
!
! First ensure dataset is CAM-ready
!
if (nf_inq_attid (aernid, nf_global, 'cam-ready', attnum) /= nf_noerr) then
write(6,*)'interpaerosols needs to be run to create a cam-ready aerosol dataset'
call endrun ('aerosol_initialize: global attribute cam-ready not found')
end if
!
! Get and check dimension info
!
call wrap_inq_dimid( aernid, 'lon', londimid )
call wrap_inq_dimid( aernid, 'lev', levdimid )
call wrap_inq_dimid( aernid, 'time', timeid )
call wrap_inq_dimid( aernid, 'lat', latdimid )
call wrap_inq_dimlen( aernid, londimid, naerlon )
call wrap_inq_dimlen( aernid, levdimid, naerlev )
call wrap_inq_dimlen( aernid, latdimid, naerlat )
call wrap_inq_dimlen( aernid, timeid, timesiz )
#if ( !defined SCAM )
if (naerlon /= plon .or. naerlat /= plat .or. naerlev /= paerlev .or. timesiz /= 12) then
write(6,*)'aerosol_initialize: grid mismatch'
write(6,*)'Model: plon, plat, naerlev, 12 =', plon, plat, naerlev, 12
write(6,*)'Dataset: naerlon, naerlat, paerlev, timesiz=', naerlon, naerlat, naerlev, timesiz
call endrun ()
end if
#endif
call wrap_inq_varid( aernid, 'date', dateid )
call wrap_inq_varid( aernid, 'datesec', secid )
do m = 1, naer
call wrap_inq_varid( aernid, TRIM(aerosol_name(m)), species_id(m))
end do
call wrap_inq_varid( aernid, 'lon', lonid )
call wrap_inq_varid( aernid, 'lat', latid )
!
! quick sanity check on one field
!
call wrap_inq_vardimid (aernid, species_id(1), dimids)
if ( (dimids(4) /= timeid) .or. &
(dimids(3) /= levdimid) .or. &
(dimids(2) /= latdimid) .or. &
(dimids(1) /= londimid) ) then
write(6,*)'AEROSOL READ: Data must be ordered time, lev, lat, lon'
write(6,*)'data are ordered as', dimids(4), dimids(3), dimids(2), dimids(1)
write(6,*)'data should be ordered as', timeid, levdimid, latdimid, londimid
call endrun ()
end if
!
! use hybi,PS from MATCH
!
call wrap_inq_varid( aernid, 'hybi', Mhybiid )
call wrap_inq_varid( aernid, 'PS', Mpsid )
!
! check dimension order for MATCH's surface pressure
!
call wrap_inq_vardimid (aernid, Mpsid, dimids)
if ( (dimids(3) /= timeid) .or. &
(dimids(2) /= latdimid) .or. &
(dimids(1) /= londimid) ) then
write(6,*)'AEROSOL READ: Pressure must be ordered time, lat, lon'
write(6,*)'data are ordered as', dimids(3), dimids(2), dimids(1)
write(6,*)'data should be ordered as', timeid, levdimid, latdimid, londimid
call endrun ()
end if
!
! read in hybi from MATCH
!
call wrap_get_var_realx (aernid, Mhybiid, M_hybi)
!
! Retrieve date and sec variables.
!
call wrap_get_var_int (aernid, dateid, date_aer)
if (timesiz < 12) then
write(6,*)'AEROSOL READ: When cycling aerosols, dataset must have 12 consecutive ', &
'months of data starting with Jan'
write(6,*)'Current dataset has only ',timesiz,' months'
call endrun ()
end if
do mo = 1,12
if (mod(date_aer(mo),10000)/100 /= mo) then
write(6,*)'AEROSOL READ: When cycling aerosols, dataset must have 12 consecutive ', &
'months of data starting with Jan'
write(6,*)'Month ',mo,' of dataset says date=',date_aer(mo)
call endrun ()
end if
end do
end if ! masterproc
!
! broadcast hybi to nodes
!
#if ( defined SPMD )
call mpibcast (M_hybi, paerlev+1, mpir8, 0, mpicom)
#endif
caldayloc = get_curr_calday ()
if (caldayloc < Mid(1)) then
mo_prv = 12
mo_nxt = 1
else if (caldayloc >= Mid(12)) then
mo_prv = 12
mo_nxt = 1
else
do i = 2 , 12
if (caldayloc < Mid(i)) then
mo_prv = i-1
mo_nxt = i
exit
end if
end do
end if
!
! Set initial calendar day values
!
cdaym = Mid(mo_prv)
cdayp = Mid(mo_nxt)
!
! Retrieve Aerosol Masses (kg/m^2 in each layer), transpose to model order (lon,lev,lat),
! then scatter to slaves.
!
if (nm /= 1 .or. np /= 2) call endrun ('AEROSOL_INITIALIZE: bad nm or np value')
do n=nm,np
if (n == 1) then
mo = mo_prv
else
mo = mo_nxt
end if
do m=1,naer
if (masterproc) then
#if ( defined SCAM )
call getncdata (aernid, initLatIdx, initLonIdx, mo, &
TRIM(aerosol_name(m)), coldata(:), RET)
aerosol_data(1,1,:)=coldata(:)
#else
start(:) = (/1,1,1,mo/)
kount(:) = (/plon,plat,paerlev,1/)
call wrap_get_vara_realx (aernid, species_id(m), start, kount, aerosol_data(1,1,1))
#endif
do j=1,plat
do k=1,paerlev
aerosol_field(:,k,j) = aerosol_data(:,j,k)
end do
aerosol_field(:,paerlev+1,j) = 0. ! value at bottom
end do
end if
call scatter_field_to_chunk (1, paerlev+1, 1, plon, aerosol_field, AEROSOLc(1,1,begchunk,m,n))
end do
!
! Retrieve PS from Match
!
if (masterproc) then
#if ( defined SCAM )
call getncdata (aernid, initLatIdx, initLonIdx, mo, &
'PS', M_ps(1,1), RET)
#else
start(:) = (/1,1,mo,-1/)
kount(:) = (/plon,plat,1,-1/)
call wrap_get_vara_realx (aernid, Mpsid, start, kount, M_ps(1,1))
#endif
end if
call scatter_field_to_chunk (1, 1, 1, plon, M_ps(1,1), M_ps_cam_col(1,begchunk,n))
end do ! n=nm,np (=1,2)
!
! read in volcanic grid data (structural only, no masses)
!
if (strat_volcanic) then
call volcanic_initialize()
endif
if ( scenario_carbon_scale == "RAMPED") then
call init_scale()
endif
call print_memusage('End aerosol_initialize')
call t_stopf ('aerosol_init')
return
end subroutine aerosol_initialize
subroutine get_aerosol(c, pint, AEROSOLt, scale) 2,10
!------------------------------------------------------------------
!
! Input:
! time at which aerosol mmrs are needed (get_curr_calday())
! chunk index
! CAM's vertical grid (pint)
!
! Output:
! values for Aerosol Mass Mixing Ratios at specified time
! on vertical grid specified by CAM (AEROSOLt)
!
! Method:
! first determine which indexs of aerosols are the bounding data sets
! interpolate both onto vertical grid aerm(),aerp().
! from those two, interpolate in time.
!
!------------------------------------------------------------------
use volcanicmass, only: get_volcanic_mass
use timeinterp, only: getfactors
!
! aerosol fields interpolated to current time step
! on pressure levels of this time step.
! these should be made read-only for other modules
! Is allocation done correctly here?
!
integer, intent(in) :: c ! Chunk Id.
real(r8), intent(in) :: pint(pcols,pverp) ! midpoint pres.
real(r8), intent(in) :: scale(naer_all) ! scale each aerosol by this amount
real(r8), intent(out) :: AEROSOLt(pcols, pver, naer_all) ! aerosols
!
! Local workspace
!
real(r8) caldayloc ! calendar day of current timestep
real(r8) fact1, fact2 ! time interpolation factors
integer i, k, j ! spatial indices
integer m ! constituent index
integer lats(pcols),lons(pcols) ! latitude and longitudes of column
integer ncol ! number of columns
real(r8) speciesmin(naer) ! minimal value for each species
!
! values before current time step "the minus month"
! aerosolm(pcols,pver) is value of preceeding month's aerosol mmr
! aerosolp(pcols,pver) is value of next month's aerosol mmr
! (think minus and plus or values to left and right of point to be interpolated)
!
real(r8) AEROSOLm(pcols,pver,naer) ! aerosol mmr from MATCH in column at previous (minus) month
!
! values beyond (or at) current time step "the plus month"
!
real(r8) AEROSOLp(pcols,pver,naer) ! aerosol mmr from MATCH in column at next (plus) month
logical error_found
!------------------------------------------------------------------
caldayloc = get_curr_calday ()
!
! Determine time interpolation factors. 1st arg says we are cycling 1 year of data
!
call getfactors (.true., mo_nxt, cdaym, cdayp, caldayloc, &
fact1, fact2, 'GET_AEROSOL:')
!
! interpolate (prv and nxt month) bounding datasets onto cam vertical grid.
! compute mass mixing ratios on CAMS's pressure coordinate
! for both the "minus" and "plus" months
!
ncol = get_ncols_p(c)
call vert_interpolate (M_ps_cam_col(1,c,nm), pint, nm, AEROSOLm, ncol, c)
call vert_interpolate (M_ps_cam_col(1,c,np), pint, np, AEROSOLp, ncol, c)
!
! Time interpolate.
!
do m=1,naer
do k=1,pver
do i=1,ncol
AEROSOLt(i,k,m) = AEROSOLm(i,k,m)*fact1 + AEROSOLp(i,k,m)*fact2
end do
end do
end do
do i=1,ncol
Match_ps_chunk(i,c) = M_ps_cam_col(i,c,nm)*fact1 + M_ps_cam_col(i,c,np)*fact2
end do
!
! get background aerosol (tuning) field
!
call background (c, ncol, pint, AEROSOLt(:, :, idxBG))
!
! find volcanic aerosol masses
!
if (strat_volcanic) then
call get_volcanic_mass(c, AEROSOLt(:,:,idxVOLC))
else
AEROSOLt(:,:,idxVOLC) = 0._r8
endif
!
! exit if mmr is negative (we have previously set
! cumulative mass to be a decreasing function.)
!
speciesmin(:) = 0. ! speciesmin(m) = 0 is minimum mmr for each species
error_found = .false.
do m=1,naer
do k=1,pver
do i=1,ncol
if (AEROSOLt(i, k, m) < speciesmin(m)) error_found = .true.
end do
end do
end do
if (error_found) then
do m=1,naer
do k=1,pver
do i=1,ncol
if (AEROSOLt(i, k, m) < speciesmin(m)) then
write(6,*) 'AEROSOL_INTERPOLATE: negative mass mixing ratio, exiting'
write(6,*) 'm, column, pver',m, i, k ,AEROSOLt(i, k, m)
call endrun ()
end if
end do
end do
end do
end if
!
! scale any AEROSOLS as required
!
call scale_aerosols (AEROSOLt, ncol, c, scale)
return
end subroutine get_aerosol
subroutine vert_interpolate (Match_ps, pint, n, AEROSOL_mmr, ncol, c) 5,6
!--------------------------------------------------------------------
! Input: match surface pressure, cam interface pressure,
! month index, number of columns, chunk index
!
! Output: Aerosol mass mixing ratio (AEROSOL_mmr)
!
! Method:
! interpolate column mass (cumulative) from match onto
! cam's vertical grid (pressure coordinate)
! convert back to mass mixing ratio
!
!--------------------------------------------------------------------
use physconst, only: gravit
real(r8), intent(out) :: AEROSOL_mmr(pcols,pver,naer) ! aerosol mmr from MATCH
real(r8), intent(in) :: Match_ps(pcols) ! surface pressure at a particular month
real(r8), intent(in) :: pint(pcols,pverp) ! interface pressure from CAM
integer, intent(in) :: ncol,c ! chunk index and number of columns
integer, intent(in) :: n ! prv or nxt month index
!
! Local workspace
!
integer m ! index to aerosol species
integer kupper(pcols) ! last upper bound for interpolation
integer i, k, kk, kkstart, kount ! loop vars for interpolation
integer isv, ksv, msv ! loop indices to save
logical bad ! indicates a bad point found
logical lev_interp_comp ! interpolation completed for a level
logical error_found
real(r8) AEROSOL(pcols,pverp,naer) ! cumulative mass of aerosol in column beneath upper
! interface of level in column at particular month
real(r8) dpl, dpu ! lower and upper intepolation factors
real(r8) v_coord ! vertical coordinate
real(r8) m_to_mmr ! mass to mass mixing ratio conversion factor
real(r8) AER_diff ! temp var for difference between aerosol masses
call t_startf ('vert_interpolate')
!
! Initialize index array
!
do i=1,ncol
kupper(i) = 1
end do
!
! assign total mass to topmost level
!
AEROSOL(:,1,:) = AEROSOLc(:,1,c,:,n)
!
! At every pressure level, interpolate onto that pressure level
!
do k=2,pver
!
! Top level we need to start looking is the top level for the previous k
! for all longitude points
!
kkstart = paerlev+1
do i=1,ncol
kkstart = min0(kkstart,kupper(i))
end do
kount = 0
!
! Store level indices for interpolation
!
! for the pressure interpolation should be comparing
! pint(column,lev) with M_hybi(lev)*M_ps_cam_col(month,column,chunk)
!
lev_interp_comp = .false.
do kk=kkstart,paerlev
if(.not.lev_interp_comp) then
do i=1,ncol
v_coord = pint(i,k)
if (M_hybi(kk)*Match_ps(i) .lt. v_coord .and. v_coord .le. M_hybi(kk+1)*Match_ps(i)) then
kupper(i) = kk
kount = kount + 1
end if
end do
!
! If all indices for this level have been found, do the interpolation and
! go to the next level
!
! Interpolate in pressure.
!
if (kount.eq.ncol) then
do m=1,naer
do i=1,ncol
dpu = pint(i,k) - M_hybi(kupper(i))*Match_ps(i)
dpl = M_hybi(kupper(i)+1)*Match_ps(i) - pint(i,k)
AEROSOL(i,k,m) = &
(AEROSOLc(i,kupper(i) ,c,m,n)*dpl + &
AEROSOLc(i,kupper(i)+1,c,m,n)*dpu)/(dpl + dpu)
enddo !i
end do
lev_interp_comp = .true.
end if
end if
end do
!
! If we've fallen through the kk=1,levsiz-1 loop, we cannot interpolate and
! must extrapolate from the bottom or top pressure level for at least some
! of the longitude points.
!
if(.not.lev_interp_comp) then
do m=1,naer
do i=1,ncol
if (pint(i,k) .lt. M_hybi(1)*Match_ps(i)) then
AEROSOL(i,k,m) = AEROSOLc(i,1,c,m,n)
else if (pint(i,k) .gt. M_hybi(paerlev+1)*Match_ps(i)) then
AEROSOL(i,k,m) = 0.0
else
dpu = pint(i,k) - M_hybi(kupper(i))*Match_ps(i)
dpl = M_hybi(kupper(i)+1)*Match_ps(i) - pint(i,k)
AEROSOL(i,k,m) = &
(AEROSOLc(i,kupper(i) ,c,m,n)*dpl + &
AEROSOLc(i,kupper(i)+1,c,m,n)*dpu)/(dpl + dpu)
end if
end do
end do
if (kount.gt.ncol) then
call endrun ('VERT_INTERPOLATE: Bad data: non-monotonicity suspected in dependent variable')
end if
end if
end do
! call t_startf ('vi_checks')
!
! aerosol mass beneath lowest interface (pverp) must be 0
!
AEROSOL(1:ncol,pverp,:) = 0.
!
! Set mass in layer to zero whenever it is less than
! 1.e-40 kg/m^2 in the layer
!
do m = 1, naer
do k = 1, pver
do i = 1, ncol
if (AEROSOL(i,k,m) < 1.e-40) AEROSOL(i,k,m) = 0.
end do
end do
end do
!
! Set mass in layer to zero whenever it is less than
! 10^-15 relative to column total mass
! convert back to mass mixing ratios.
! exit if mmr is negative
!
error_found = .false.
do m = 1, naer
do k = 1, pver
do i = 1, ncol
AER_diff = AEROSOL(i,k,m) - AEROSOL(i,k+1,m)
if( abs(AER_diff) < 1e-15*AEROSOL(i,1,m)) then
AER_diff = 0.
end if
m_to_mmr = gravit / (pint(i,k+1)-pint(i,k))
AEROSOL_mmr(i,k,m)= AER_diff * m_to_mmr
if (AEROSOL_mmr(i,k,m) < 0) error_found = .true.
end do
end do
end do
if (error_found) then
do m = 1, naer
do k = 1, pver
do i = 1, ncol
if (AEROSOL_mmr(i,k,m) < 0) then
write(6,*)'vert_interpolate: mmr < 0, m, col, lev, mmr',m, i, k, AEROSOL_mmr(i,k,m)
write(6,*)'vert_interpolate: aerosol(k),(k+1)',AEROSOL(i,k,m),AEROSOL(i,k+1,m)
write(6,*)'vert_interpolate: pint(k+1),(k)',pint(i,k+1),pint(i,k)
write(6,*)'n,c',n,c
call endrun()
end if
end do
end do
end do
end if
! call t_stopf ('vi_checks')
call t_stopf ('vert_interpolate')
return
end subroutine vert_interpolate
subroutine scale_aerosols(AEROSOLt, ncol, lchnk, scale) 1
!-----------------------------------------------------------------
! scale each species as determined by scale factors
!-----------------------------------------------------------------
integer, intent(in) :: ncol, lchnk ! number of columns and chunk index
real(r8), intent(in) :: scale(naer_all) ! scale each aerosol by this amount
real(r8), intent(inout) :: AEROSOLt(pcols, pver, naer_all) ! aerosols
integer m
do m = 1, naer_all
AEROSOLt(:ncol, :, m) = scale(m)*AEROSOLt(:ncol, :, m)
end do
return
end subroutine scale_aerosols
subroutine background(lchnk, ncol, pint, mmr) 1,3
!-----------------------------------------------------------------------
!
! Purpose:
! Set global mean tropospheric aerosol background (or tuning) field
!
! Method:
! Specify aerosol mixing ratio.
! Aerosol mass mixing ratio
! is specified so that the column visible aerosol optical depth is a
! specified global number (tauback). This means that the actual mixing
! ratio depends on pressure thickness of the lowest three atmospheric
! layers near the surface.
!
!-----------------------------------------------------------------------
use shr_kind_mod, only: r8 => shr_kind_r8
use aer_optics, only: kbg,idxVIS
use physconst, only: gravit
!-----------------------------------------------------------------------
implicit none
!-----------------------------------------------------------------------
#include <ptrrgrid.h>
!------------------------------Arguments--------------------------------
!
! Input arguments
!
integer, intent(in) :: lchnk ! chunk identifier
integer, intent(in) :: ncol ! number of atmospheric columns
real(r8), intent(in) :: pint(pcols,pverrp) ! Interface pressure (mks)
!
! Output arguments
!
real(r8), intent(out) :: mmr(pcols,pverr) ! "background" aerosol mass mixing ratio
!
!---------------------------Local variables-----------------------------
!
integer i ! Longitude index
integer k ! Level index
!
real(r8) mass2mmr ! Factor to convert mass to mass mixing ratio
real(r8) mass ! Mass of "background" aerosol as specified by tauback
!
!-----------------------------------------------------------------------
!
do i=1,ncol
mass2mmr = gravit / (pint(i,pverrp)-pint(i,pverrp-mxaerl))
do k=1,pverr
!
! Compute aerosol mass mixing ratio for specified levels (1.e3 factor is
! for units conversion of the extinction coefficiant from m2/g to m2/kg)
!
if ( k >= pverrp-mxaerl ) then
! kaervs is not consistent with the values in aer_optics
! this ?should? be changed.
! rhfac is also implemented differently
mass = tauback / (1.e3 * kbg(idxVIS))
mmr(i,k) = mass2mmr*mass
else
mmr(i,k) = 0._r8
endif
!
enddo
enddo
!
return
end subroutine background
subroutine aerosol_diagnostics(state, AEROSOL_mmr) 1,16
!-----------------------------------------------------------------
! write out the mass of aerosol in each layer
!-----------------------------------------------------------------
use physics_types, only: physics_state
use physconst, only: rga
use history, only: outfld
type(physics_state), intent(in) :: state ! state from physics code
real(r8), intent(in) :: AEROSOL_mmr(pcols, pver, naer_all) ! aerosols
integer lchnk,ncol ! chunk index and number of columns
real(r8) :: MAerosol(pcols, pver, naer_all) !mass of aerosols in each layer
! column totals of masses
real(r8) :: TMAllAerosols(pcols), TMCarbon(pcols), TMDust(pcols)
real(r8) :: TMAerosol(pcols,naer_all)
integer k,m ! level, constituent indices
TMAllAerosols = 0.0
TMCarbon = 0.0
TMDust = 0.0
TMAerosol = 0.0
ncol = state%ncol
lchnk = state%lchnk
!
! output surface pressure from MATCH
!
call outfld('PS_match',Match_ps_chunk(1,lchnk),pcols,lchnk)
!
! compute column mass of aerosols (and groups)
!
! first, mass in each layer
!
do m = 1, naer_all-1
MAerosol(:ncol,:,m) = AEROSOL_mmr(:ncol,:,m)*state%pdel(:ncol,:)*rga
enddo
MAerosol(:ncol,:,idxVOLC) = AEROSOL_mmr(:ncol,:,idxVOLC)
!
! now accumulate mass below upper interface bounding each layer
!
do k = pver-1,1,-1
MAerosol(:ncol,k,:) = MAerosol(:ncol,k+1,:) + MAerosol(:ncol,k,:)
enddo
!
! total mass in column is given by cumulative mass at level = 1
!
TMAerosol(:ncol,:) = MAerosol(:ncol,1,:)
do m = idxCARBONfirst, idxCARBONfirst+numCARBON-1
TMCarbon(:ncol) = TMCarbon(:ncol) + TMAerosol(:ncol,m)
enddo
do m = idxDUSTfirst, idxDUSTfirst+numDUST-1
TMDust(:ncol) = TMDust(:ncol) + TMAerosol(:ncol,m)
enddo
TMAllAerosols(:ncol)=TMCarbon(:ncol) &
+TMDust(:ncol) &
+TMAerosol(:ncol,idxSUL ) &
+TMAerosol(:ncol,idxSSLT) &
+TMAerosol(:ncol,idxVOLC)
!
! write out column totals of aerosol mmr's and groups
!
call outfld('MSUL_V ', MAerosol(:,:,idxSUL) ,pcols, lchnk)
call outfld('MSSLT_V ', MAerosol(:,:,idxSSLT) ,pcols, lchnk)
call outfld('MDUST1_V', MAerosol(:,:,idxDUSTfirst) ,pcols, lchnk)
call outfld('MDUST2_V', MAerosol(:,:,idxDUSTfirst+1) ,pcols, lchnk)
call outfld('MDUST3_V', MAerosol(:,:,idxDUSTfirst+2) ,pcols, lchnk)
call outfld('MDUST4_V', MAerosol(:,:,idxDUSTfirst+3) ,pcols, lchnk)
call outfld('MOCPHO_V', MAerosol(:,:,idxOCPHO) ,pcols, lchnk)
call outfld('MOCPHI_V', MAerosol(:,:,idxOCPHI) ,pcols, lchnk)
call outfld('MBCPHO_V', MAerosol(:,:,idxBCPHO) ,pcols, lchnk)
call outfld('MBCPHI_V', MAerosol(:,:,idxBCPHI) ,pcols, lchnk)
call outfld('MBG_V ', MAerosol(:,:,idxBG) ,pcols, lchnk)
call outfld('MVOLC ', MAerosol(:,:,idxVOLC) ,pcols, lchnk)
return
end subroutine aerosol_diagnostics
subroutine aerosol_indirect(ncol,lchnk,landfrac,pmid,t,qm1,cld,zm,rel) 2,9
!--------------------------------------------------------------
! Compute effect of sulfate on effective liquid water radius
! Method of Martin et. al.
!--------------------------------------------------------------
use constituents, only: ppcnst, cnst_get_ind
use history, only: outfld
#include <comctl.h>
integer, intent(in) :: ncol ! number of atmospheric columns
integer, intent(in) :: lchnk ! chunk identifier
real(r8), intent(in) :: landfrac(pcols) ! land fraction
real(r8), intent(in) :: pmid(pcols,pver) ! Model level pressures
real(r8), intent(in) :: t(pcols,pver) ! Model level temperatures
real(r8), intent(in) :: qm1(pcols,pver,ppcnst) ! Specific humidity and tracers
real(r8), intent(in) :: cld(pcols,pver) ! Fractional cloud cover
real(r8), intent(in) :: zm(pcols,pver) ! Height of midpoints (above surface)
real(r8), intent(in) :: rel(pcols,pver) ! liquid effective drop size (microns)
!
! local variables
!
real(r8) locrhoair(pcols,pver) ! dry air density [kg/m^3 ]
real(r8) lwcwat(pcols,pver) ! in-cloud liquid water path [kg/m^3 ]
real(r8) sulfmix(pcols,pver) ! sulfate mass mixing ratio [kg/kg ]
real(r8) so4mass(pcols,pver) ! sulfate mass concentration [g/cm^3 ]
real(r8) Aso4(pcols,pver) ! sulfate # concentration [#/cm^3 ]
real(r8) Ntot(pcols,pver) ! ccn # concentration [#/cm^3 ]
real(r8) relmod(pcols,pver) ! effective radius [microns]
real(r8) wrel(pcols,pver) ! weighted effective radius [microns]
real(r8) wlwc(pcols,pver) ! weighted liq. water content [kg/m^3 ]
real(r8) cldfrq(pcols,pver) ! frequency of occurance of...
! ! clouds (cld => 0.01) [fraction]
real(r8) locPi ! my piece of the pi
real(r8) Rdryair ! gas constant of dry air [J/deg/kg]
real(r8) rhowat ! density of water [kg/m^3 ]
real(r8) Acoef ! m->A conversion factor; assumes
! ! Dbar=0.10, sigma=2.0 [g^-1 ]
real(r8) rekappa ! kappa in evaluation of re(lmod)
real(r8) recoef ! temp. coeficient for calc of re(lmod)
real(r8) reexp ! 1.0/3.0
real(r8) Ntotb ! temp var to hold below cloud ccn
! -- Parameters for background CDNC (from `ambient' non-sulfate aerosols)...
real(r8) Cmarn ! Coef for CDNC_marine [cm^-3]
real(r8) Cland ! Coef for CDNC_land [cm^-3]
real(r8) Hmarn ! Scale height for CDNC_marine [m]
real(r8) Hland ! Scale height for CDNC_land [m]
parameter ( Cmarn = 50.0, Cland = 100.0 )
parameter ( Hmarn = 1000.0, Hland = 2000.0 )
real(r8) bgaer ! temp var to hold background CDNC
integer :: ixcldliq ! index of liquid cloud water
integer i,k ! loop indices
!
! Statement functions
!
logical land ! is this a column over land?
land(i) = nint(landfrac(i)).gt.0.5_r8
if (indirect) then
call endrun ('AEROSOL_INDIRECT: indirect effect is obsolete')
! ramping is not yet resolved so sulfmix is 0.
sulfmix(1:ncol,1:pver) = 0._r8
locPi = 3.141592654
Rdryair = 287.04
rhowat = 1000.0
Acoef = 1.2930E14
recoef = 3.0/(4.0*locPi*rhowat)
reexp = 1.0/3.0
call cnst_get_ind('CLDLIQ', ixcldliq)
do k=pver,1,-1
do i = 1,ncol
locrhoair(i,k) = pmid(i,k)/( Rdryair*t(i,k) )
lwcwat(i,k) = ( qm1(i,k,ixcldliq)/max(0.01_r8,cld(i,k)) )* &
locrhoair(i,k)
! NOTE: 0.001 converts kg/m3 -> g/cm3
so4mass(i,k) = sulfmix(i,k)*locrhoair(i,k)*0.001
Aso4(i,k) = so4mass(i,k)*Acoef
if (Aso4(i,k) <= 280.0) then
Aso4(i,k) = max(36.0_r8,Aso4(i,k))
Ntot(i,k) = -1.15E-3*Aso4(i,k)**2 + 0.963*Aso4(i,k)+5.30
rekappa = 0.80
else
Aso4(i,k) = min(1500.0_r8,Aso4(i,k))
Ntot(i,k) = -2.10E-4*Aso4(i,k)**2 + 0.568*Aso4(i,k)-27.9
rekappa = 0.67
end if
if (land(i)) then ! Account for local background aerosol;
bgaer = Cland*exp(-(zm(i,k)/Hland))
Ntot(i,k) = max(bgaer,Ntot(i,k))
else
bgaer = Cmarn*exp(-(zm(i,k)/Hmarn))
Ntot(i,k) = max(bgaer,Ntot(i,k))
end if
if (k == pver) then
Ntotb = Ntot(i,k)
else
Ntotb = Ntot(i,k+1)
end if
relmod(i,k) = (( (recoef*lwcwat(i,k))/(rekappa*Ntotb))**reexp)*10000.0
relmod(i,k) = max(4.0_r8,relmod(i,k))
relmod(i,k) = min(20.0_r8,relmod(i,k))
if (cld(i,k) >= 0.01) then
cldfrq(i,k) = 1.0
else
cldfrq(i,k) = 0.0
end if
wrel(i,k) = relmod(i,k)*cldfrq(i,k)
wlwc(i,k) = lwcwat(i,k)*cldfrq(i,k)
end do
end do
call outfld('MSO4 ',so4mass,pcols,lchnk)
call outfld('LWC ',lwcwat ,pcols,lchnk)
call outfld('CLDFRQ ',cldfrq ,pcols,lchnk)
call outfld('WREL ',wrel ,pcols,lchnk)
call outfld('WLWC ',wlwc ,pcols,lchnk)
write(6,*)'WARNING: indirect calculation has no effects'
else
do k = 1, pver
do i = 1, ncol
relmod(i,k) = rel(i,k)
end do
end do
endif
return
end subroutine aerosol_indirect
subroutine setmxaerl() 1
!------------------------------------
!
! set levels for background aerosol
!
!-----------------------------------
! need hypm from comhyb.h
#include <comhyb.h>
integer k ! index through vertical levels
!
! mxaerl = max number of levels (from bottom) for background aerosol
! Limit background aerosol height to regions below 900 mb
!
mxaerl = 0
do k=pver,1,-1
if (hypm(k) >= 9.e4) mxaerl = mxaerl + 1
end do
mxaerl = max(mxaerl,1)
if (masterproc) then
write(6,*)'AEROSOLS: Background aerosol will be limited to ', &
'bottom ',mxaerl,' model interfaces. Top interface is ', &
hypi(pverp-mxaerl),' pascals'
end if
return
end subroutine setmxaerl
subroutine add_aer_diag_fields() 1,30
!-------------------------------------------------------------
! make sure aerosols are in the master list so they can be written
!-------------------------------------------------------------
use history, only: phys_decomp, addfld
call addfld ('PS_match','N/m^2 ',1, 'I','Surface Pressure from aerosol climatology' ,phys_decomp)
call addfld ('frc_day ','None ',1, 'I','Portion of time column is lit' ,phys_decomp)
call addfld ('MSUL_V ','Kg/m^2 ',pver,'I','Mass of Sulfate in and below layer',phys_decomp)
call addfld ('MSSLT_V ','Kg/m^2 ',pver,'I','Mass of Sea Salt in and below layer',phys_decomp)
call addfld ('MDUST1_V','Kg/m^2 ',pver,'I','Mass of Dust bin 1 in and below layer',phys_decomp)
call addfld ('MDUST2_V','Kg/m^2 ',pver,'I','Mass of Dust bin 2 in and below layer',phys_decomp)
call addfld ('MDUST3_V','Kg/m^2 ',pver,'I','Mass of Dust bin 3 in and below layer',phys_decomp)
call addfld ('MDUST4_V','Kg/m^2 ',pver,'I','Mass of Dust bin 4 in and below layer',phys_decomp)
call addfld ('MOCPHO_V','Kg/m^2 ',pver,'I','Mass of OCPHO in and below layer',phys_decomp)
call addfld ('MOCPHI_V','Kg/m^2 ',pver,'I','Mass of OCPHI in and below layer',phys_decomp)
call addfld ('MBCPHO_V','Kg/m^2 ',pver,'I','Mass of BCPHO in and below layer',phys_decomp)
call addfld ('MBCPHI_V','Kg/m^2 ',pver,'I','Mass of BCPHI in and below layer',phys_decomp)
call addfld ('MBG_V ','Kg/m^2 ',pver,'I','Mass of Background Aerosol in and below layer',phys_decomp)
call addfld ('MVOLC ','Kg/m^2 ',pver,'I','Mass of Volcanic Aerosol in layer',phys_decomp)
call addfld ('SULOD_v ','None ',1, 'I','Sulfate Optical Depth in visible', phys_decomp)
call addfld ('SSLTOD_v','None ',1, 'I','Sea Salt Optical Depth in visible', phys_decomp)
call addfld ('CAROD_v ','None ',1, 'I','Carbon Optical Depth in visible', phys_decomp)
call addfld ('DUSTOD_v','None ',1, 'I','Dust Optical Depth in visible', phys_decomp)
call addfld ('BGOD_v ','None ',1, 'I','Background Aerosol Optical Depth in visible', phys_decomp)
call addfld ('VOLCOD_v','None ',1, 'I','Volcanic Aerosol Optical Depth in visible', phys_decomp)
call addfld ('AEROD_v ','None ',1, 'I','Total Aerosol Optical Depth in visible', phys_decomp)
call addfld ('AERSSA_v','None ',1, 'I','Total Aerosol Single Scattering Albedo in visible', phys_decomp)
call addfld ('AERASM_v','None ',1, 'I','Total Aerosol Asymmetry Parameter in visible', phys_decomp)
call addfld ('AERFWD_v','None ',1, 'I','Total Aerosol Forward Scattering in visible', phys_decomp)
!
! addfld calls for aerosol forcing-only calculations
!
if(radforce) then
call addfld ('FSNT_RF ','W/m^2 ',1, 'I','Total column absorbed solar flux (radforce)' ,phys_decomp)
call addfld ('FSNTC_RF','W/m^2 ',1, 'I','Clear sky total column absorbed solar flux (radforce)' ,phys_decomp)
call addfld ('FSNS_RF ','W/m^2 ',1, 'I','Surface absorbed solar flux (radforce)' ,phys_decomp)
call addfld ('FSNSC_RF','W/m^2 ',1, 'I','Clear sky surface absorbed solar flux (radforce)' ,phys_decomp)
call addfld ('QRS_RF ','K/s ',pver, 'I','Solar heating rate (radforce)' ,phys_decomp)
endif
return
end subroutine add_aer_diag_fields
subroutine get_rf_scales(scales) 1
real(r8), intent(out)::scales(naer_all) ! scale aerosols by this amount
integer i ! loop index
scales(idxBG) = bgscl_rf
scales(idxSUL) = sulscl_rf
scales(idxSSLT) = ssltscl_rf
do i = idxCARBONfirst, idxCARBONfirst+numCARBON-1
scales(i) = carscl_rf
enddo
do i = idxDUSTfirst, idxDUSTfirst+numDUST-1
scales(i) = dustscl_rf
enddo
scales(idxVOLC) = volcscl_rf
end subroutine get_rf_scales
subroutine get_int_scales(scales) 1
real(r8), intent(out)::scales(naer_all) ! scale each aerosol by this amount
integer i ! index through species
scales(idxBG) = 1._r8
scales(idxSUL) = sulscl
scales(idxSSLT) = ssltscl
do i = idxCARBONfirst, idxCARBONfirst+numCARBON-1
scales(i) = carscl
enddo
do i = idxDUSTfirst, idxDUSTfirst+numDUST-1
scales(i) = dustscl
enddo
scales(idxVOLC) = volcscl
return
end subroutine get_int_scales
subroutine aerint () 1,7
use carbonscales, only: get_carbonscale
implicit none
integer :: ntmp ! used in index swapping
integer :: start(4) ! start vector for netcdf calls
integer :: kount(4) ! count vector for netcdf calls
integer :: i,j,k ! spatial indices
integer :: m ! constituent index
real(r8) :: caldayloc ! calendar day of current timestep
real(r8) :: aerosol_data(plon,plat,paerlev) ! aerosol field read in from dataset
real(r8) :: aerosol_field(plon,paerlev+1,plat) ! aerosol field to be scattered
!
! determine if need to read in next month data
! also determine time interpolation factors
!
caldayloc = get_curr_calday ()
!
! If model time is past current forward ozone timeslice, then
! masterproc reads in the next timeslice for time interpolation. Messy logic is
! for interpolation between December and January (mo_nxt == 1). Just like
! oznint, sstint.
!
if (caldayloc > cdayp .and. .not. (mo_nxt == 1 .and. caldayloc > cdaym)) then
mo_nxt = mod(mo_nxt,12) + 1
cdaym = cdayp
cdayp = Mid(mo_nxt)
!
! Check for valid date info
!
if (.not. (mo_nxt == 1 .or. caldayloc <= cdayp)) then
call endrun ('AERINT: Non-monotonicity suspected in input aerosol data')
end if
ntmp = nm
nm = np
np = ntmp
do m=1,naer
if (masterproc) then
start(:) = (/1,1,1,mo_nxt/)
kount(:) = (/plon,plat,paerlev,1/)
call wrap_get_vara_realx (aernid, species_id(m), start, kount, aerosol_data(1,1,1))
do j=1,plat
do k=1,paerlev
aerosol_field(:,k,j) = aerosol_data(:,j,k)
end do
aerosol_field(:,paerlev+1,j) = 0. ! value at bottom
end do
end if
call scatter_field_to_chunk (1, paerlev+1, 1, plon, aerosol_field, AEROSOLc(1,1,begchunk,m,np))
end do
!
! Retrieve PS from Match
!
if (masterproc) then
start(:) = (/1,1,mo_nxt,-1/)
kount(:) = (/plon,plat,1,-1/)
call wrap_get_vara_realx (aernid, Mpsid, start, kount, M_ps(1,1))
write(6,*)'AERINT: Read aerosols data for julian day', Mid(mo_nxt)
end if
call scatter_field_to_chunk (1, 1, 1, plon, M_ps(1,1), M_ps_cam_col(1,begchunk,np))
end if
!
! Call get_carbonscale (if required) to set carbon scale
! at this particular time step
!
if(scenario_carbon_scale.eq."RAMPED") then
call get_carbonscale(carscl)
! if(masterproc) then
! write(6,*) "aerint: carscl" , carscl, 'caldayloc',caldayloc
! endif
endif
return
end subroutine aerint
end module prescribed_aerosols