! #define TRACEBACK write (gol,'("in ",a," (",a,", line",i5,")")') rname, __FILE__, __LINE__; call goErr #define IF_NOTOK_RETURN(action) if (status/=0) then; TRACEBACK; action; return; end if #define IF_ERROR_RETURN(action) if (status> 0) then; TRACEBACK; action; return; end if ! #include "tm5.inc" ! !----------------------------------------------------------------------------- ! TM5 ! !----------------------------------------------------------------------------- !BOP ! ! !MODULE: EMISSION_TERP ! ! !DESCRIPTION: data and methods for terpene emissions. !\\ !\\ ! !INTERFACE: ! MODULE EMISSION_TERP ! ! !USES: ! use GO, only : gol, goPr, goErr, goBug use tm5_distgrid, only : dgrid, get_distgrid, scatter, gather use partools, only : isRoot, par_broadcast use dims, only : nregions, idate, okdebug use global_types, only : emis_data, d3_data, isop_data use emission_read, only : used_providers_terp, has_terp_emis IMPLICIT NONE PRIVATE ! ! !PUBLIC MEMBER FUNCTIONS: ! public :: Emission_terp_Init public :: Emission_terp_Done public :: emission_terp_Declare public :: Emission_terp_Apply ! ! !PRIVATE DATA MEMBERS: ! character(len=*), parameter :: mname = 'emission_terp' type(isop_data),dimension(nregions),target :: terp_dat !flags the daily calculation of terpene distribution func. integer,dimension(nregions) :: day_terp type( emis_data ), dimension(:,:), allocatable :: terp_emis_2d type( d3_data ), dimension(:,:), allocatable :: terp_emis_3d integer :: terp_2dsec, terp_3dsec logical, allocatable :: has_data_2d(:), has_data_3d(:) ! ! !REVISION HISTORY: ! 1 Oct 2010 - Achim Strunk - overhaul for AR5 ! 28 Jun 2012 - Ph. Le Sager - adapted for lon-lat MPI domain decomposition ! ! !REMARKS: ! !EOP !------------------------------------------------------------------------ CONTAINS !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: EMISSION_TERP_INIT ! ! !DESCRIPTION: Allocate memory !\\ !\\ ! !INTERFACE: ! SUBROUTINE EMISSION_TERP_INIT( status ) ! ! !USES: ! use dims, only : ndyn_max, nsrce, tref, im, jm, lm, idate, nlon360, nlat180 use dims, only : sec_month use global_data, only : rcfile use emission_read, only : providers_def, numb_providers, sectors_def ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: status ! ! !REVISION HISTORY: ! 02 Apr 2014 - v0 ! ! !REMARKS: ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/emission_terp_init' integer :: region, i1, i2, j1, j2 integer :: lmr, lsec, ntim, lprov real :: dtime ! --- begin ----------------------------------------- status = 0 if(.not. has_terp_emis) return terp_2dsec = 0 terp_3dsec = 0 do lprov = 1, numb_providers if (count(used_providers_terp == providers_def(lprov)%name)/=0) then ! limit AR5 to the 2D biomass burning sector if ( trim(providers_def(lprov)%name) == 'AR5' ) then terp_2dsec = terp_2dsec + count( sectors_def%prov == 'AR5 ' & .and. sectors_def%catname == 'biomassburning') else ! use all sectors terp_2dsec = terp_2dsec + providers_def(lprov)%nsect2d terp_3dsec = terp_3dsec + providers_def(lprov)%nsect3d end if endif enddo allocate( terp_emis_2d( nregions, terp_2dsec ) ) allocate( terp_emis_3d( nregions, terp_3dsec ) ) allocate( has_data_2d(terp_2dsec)) ; has_data_2d=.false. allocate( has_data_3d(terp_3dsec)) ; has_data_3d=.false. ! allocate information arrays (2d and 3d) do region=1,nregions CALL GET_DISTGRID( dgrid(region), I_STRT=i1, I_STOP=i2, J_STRT=j1, J_STOP=j2 ) lmr = lm(region) do lsec=1,terp_2dsec allocate( terp_emis_2d(region,lsec)%surf(i1:i2, j1:j2) ) end do do lsec=1,terp_3dsec allocate( terp_emis_3d(region,lsec)%d3(i1:i2, j1:j2,lmr) ) end do dtime=float(ndyn_max)/(2*tref(region)) ! timestep emissions ntim=86400/nint(dtime) ! number of timesteps in 24 hours for this region allocate(terp_dat(region)%scalef_isop(ntim,j1:j2)) ! allocate enough memory! enddo ! ok status = 0 END SUBROUTINE EMISSION_TERP_INIT !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: EMISSION_terp_DONE ! ! !DESCRIPTION: Free memory !\\ !\\ ! !INTERFACE: ! SUBROUTINE EMISSION_TERP_DONE( STATUS ) ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: status ! ! !REVISION HISTORY: ! 1 Oct 2010 - Achim Strunk - adapted to new structures ! 28 Jun 2012 - Ph. Le Sager - adapted for lon-lat MPI domain decomposition ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/Emission_terp_Done' integer :: region, lsec ! --- begin -------------------------------------- status = 0 if(.not. has_terp_emis) return do region = 1, nregions do lsec=1,terp_2dsec deallocate( terp_emis_2d(region,lsec)%surf ) end do do lsec=1,terp_3dsec deallocate( terp_emis_3d(region,lsec)%d3 ) end do deallocate(terp_dat(region)%scalef_isop ) end do deallocate( terp_emis_2d, terp_emis_3d ) deallocate( has_data_2d, has_data_3d ) status = 0 END SUBROUTINE EMISSION_TERP_DONE !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: EMISSION_TERP_DECLARE ! ! !DESCRIPTION: Opens, reads and evaluates input files (per month). ! Provides emissions on 2d/3d-arrays which are then added ! to mixing ratios in routine *apply. !\\ !\\ ! !INTERFACE: ! SUBROUTINE EMISSION_TERP_DECLARE( status ) ! ! !USES: ! use toolbox, only : coarsen_emission use dims, only : im, jm, lm, idate, sec_month, nlon360, nlat180, iglbsfc use chem_param, only : xmterp use emission_data, only : msg_emis, LCMIP6BMB, LAR5BMB, LMEGAN ! ---------------- CMIP6 - AR5 - MACC -------------------- use emission_data, only : emis_input_year_nmvoc, emis_input_year_nat use emission_data, only : emis_input_dir_mac use emission_data, only : emis_input_dir_megan use emission_data, only : emis_input_dir_retro use emission_data, only : emis_input_dir_gfed use emission_read, only : emission_ar5_regrid_aircraft use emission_read, only : emission_cmip6_ReadSector use emission_read, only : emission_cmip6bmb_ReadSector use emission_read, only : emission_ar5_ReadSector use emission_read, only : emission_macc_ReadSector use emission_read, only : emission_megan_ReadSector use emission_read, only : emission_gfed_ReadSector use emission_read, only : emission_retro_ReadSector use emission_read, only : sectors_def, numb_sectors use emission_read, only : ar5_dim_3ddata ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: status ! ! !REVISION HISTORY: ! 1 Oct 2010 - Achim Strunk - adapted for AR5 ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/emission_terp_declare' integer :: region, hasData integer, parameter :: add_field=0 integer, parameter :: amonth=2 integer :: imr, jmr, lmr, lsec ! AR5 real,dimension(:,:,:),allocatable :: field3d, field3d2 type(d3_data) :: emis3d, work(nregions) type(emis_data) :: wrk2D(nregions) integer :: seccount2d, seccount3d ! --- begin ----------------------------------------- status = 0 if(.not. has_terp_emis) return write(gol,'(" EMISS-INFO ------------- read TERPENE emissions -------------")'); call goPr ! flag calculation of daily terpene emission distribution day_terp(1:nregions) = -1 do region = 1, nregions do lsec=1,terp_2dsec terp_emis_2d(region,lsec)%surf = 0.0 end do do lsec=1,terp_3dsec terp_emis_3d(region,lsec)%d3 = 0.0 end do end do ! global arrays for coarsening do region = 1, nregions if (isRoot)then allocate(work(region)%d3(im(region),jm(region),ar5_dim_3ddata)) else allocate(work(region)%d3(1,1,1)) end if enddo do region = 1, nregions wrk2D(region)%surf => work(region)%d3(:,:,1) end do ! -------------------------------- ! do a loop over available sectors ! -------------------------------- ! count 2d and 3d sectors seccount2d = 0 seccount3d = 0 ! always allocate here 3d data set (for 2d sectors it will be filled in first layer only) if (isRoot) then allocate( field3d( nlon360, nlat180, ar5_dim_3ddata ) ) ; field3d = 0.0 else allocate( field3d( 1, 1, 1 ) ) end if sec : do lsec = 1, numb_sectors ! screen out unused providers if (count(used_providers_terp.eq.sectors_def(lsec)%prov).eq.0) cycle ! Explicitly skip AR5 non-fires and EDGAR and HYMN inventories ! Same w/ ED41 and ED42 (but already screened out when building providers list anyway) if ( sectors_def(lsec)%prov == 'AR5 ' .and. sectors_def(lsec)%catname /= 'biomassburning') cycle if ( sectors_def(lsec)%prov == 'ED41 ' .or. & sectors_def(lsec)%prov == 'ED42 ' .or. & sectors_def(lsec)%prov == 'HYMN ') cycle field3d = 0.0 if( sectors_def(lsec)%f3d ) then seccount3d = seccount3d + 1 else seccount2d = seccount2d + 1 end if if (isRoot) then ! READ select case( trim(sectors_def(lsec)%prov) ) case( 'CMIP6BMB' ) call emission_cmip6bmb_ReadSector( 'NMVOC-'//'C10H16', emis_input_year_nmvoc, idate(2), lsec, field3d, status ) IF_NOTOK_RETURN(status=1;deallocate(field3d)) case( 'AR5' ) call emission_ar5_ReadSector( 'terpenes', emis_input_year_nmvoc, idate(2), lsec, field3d, status ) IF_NOTOK_RETURN(status=1) case( 'MACC' ) ! skip if MEGAN is used, screen out sectors w/o terpenes else if( (trim(sectors_def(lsec)%name) .eq. 'emiss_bio') .and. (.not. LMEGAN) ) then if (trim(sectors_def(lsec)%catname) .eq. 'natural') then call emission_macc_ReadSector( emis_input_dir_mac, 'TERPENES', emis_input_year_nat, idate(2), & '0.5x0.5_kg.nc', sectors_def(lsec)%name, 'kg / s', field3d, status ) IF_NOTOK_RETURN(status=1) else call emission_macc_ReadSector( emis_input_dir_mac, 'TERPENES', emis_input_year_nmvoc, idate(2), & '0.5x0.5_kg.nc', sectors_def(lsec)%name, 'kg / s', field3d, status ) IF_NOTOK_RETURN(status=1) endif endif case('GFEDv3') call emission_gfed_ReadSector( emis_input_dir_gfed, 'terpenes', emis_input_year_nmvoc, idate(2), & sectors_def(lsec)%name, 'kg / s', field3d(:,:,1), status ) IF_NOTOK_RETURN(status=1) case('RETRO') call emission_retro_ReadSector( emis_input_dir_retro, 'MONOTERPENES', emis_input_year_nmvoc, idate(2), & sectors_def(lsec)%name, 'kg / s', field3d(:,:,1), status ) IF_NOTOK_RETURN(status=1) case('MEGAN') call emission_megan_ReadSector( emis_input_dir_megan, 'monoterpenes', emis_input_year_nat, idate(2), & sectors_def(lsec)%name, 'kg / s', field3d(:,:,1), status ) IF_NOTOK_RETURN(status=1) case('DUMMY') case default write(gol,*) "Error in buidling list of providers USED_PROVIDERS_TERP" ; call goErr status=1; TRACEBACK; return end select ! nothing found??? if( sum(field3d) < 100.*TINY(1.0) ) then if (okdebug) then write(gol,'("EMISS-INFO - no terp emissions found for ",a," ",a," for month ",i2 )') & trim(sectors_def(lsec)%prov), trim(sectors_def(lsec)%name), idate(2) ; call goPr endif hasData=0 else if (okdebug) then write(gol,'("EMISS-INFO - found terp emissions for ",a," ",a," for month ",i2 )') & trim(sectors_def(lsec)%prov), trim(sectors_def(lsec)%name), idate(2) ; call goPr endif ! scale from kg/s to kg/month field3d = field3d * sec_month ! kg / month hasData=1 end if end if call Par_broadcast(hasData, status) IF_NOTOK_RETURN(status=1) if (hasData == 0) then cycle sec else if ( sectors_def(lsec)%f3d ) then has_data_3d(seccount3d)=.true. else has_data_2d(seccount2d)=.true. end if end if ! Distinguish 2d/3d sectors if( sectors_def(lsec)%f3d ) then write(gol,'("EMISS-ERROR - Unexpected 3D data: code needed")'); call goErr status=1; TRACEBACK; return else ! ar5_sector is 2d ! --------------------------- ! 2d data (Anthropogenic, Ships, Biomassburning) ! --------------------------- if (isRoot) then ! print total & regrid call msg_emis( amonth, trim(sectors_def(lsec)%prov), sectors_def(lsec)%name, 'TERP', xmterp, sum(field3d(:,:,1)) ) call coarsen_emission( 'terp '//sectors_def(lsec)%name, nlon360, nlat180, field3d(:,:,1), wrk2D, add_field, status ) IF_NOTOK_RETURN(status=1) end if do region = 1, nregions call scatter(dgrid(region), terp_emis_2d(region,seccount2d)%surf, work(region)%d3(:,:,1), 0, status) IF_NOTOK_RETURN(status=1) end do end if end do sec ! sectors deallocate( field3d ) do region = 1, nregions if (associated(wrk2D(region)%surf)) nullify(wrk2D(region)%surf) deallocate( work(region)%d3 ) end do ! ok status = 0 END SUBROUTINE EMISSION_TERP_DECLARE !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: EMISSION_TERP_APPLY ! ! !DESCRIPTION: Take monthly emissions, and ! - split them vertically ! - apply time splitting factors ! - add them up (add_terp) !\\ !\\ ! !INTERFACE: ! SUBROUTINE EMISSION_TERP_APPLY( region, status ) ! ! !USES: ! use dims, only : itaur, nsrce, tref use dims, only : im, jm, lm, mlen, ndyn_max use dims, only : yref, dy, ybeg use datetime, only : tau2date, get_day use emission_data, only : emission_vdist_by_sector, bb_cycle use emission_data, only : emis_bb_trop_cycle use chem_param, only : iterp, xmterp use emission_data, only : do_add_3d, do_add_3d_cycle use emission_read, only : sectors_def, numb_sectors use tracer_data, only : tracer_print USE partools, ONLY : isRoot, par_reduce, par_reduce_element ! ! !INPUT PARAMETERS: ! integer, intent(in) :: region ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: status ! ! !REVISION HISTORY: ! 1 Oct 2010 - Achim Strunk - adapted for AR5 ! 27 Apr 2012 - Ph. Le Sager - adapted for lon-lat MPI domain decomposition ! 2 Apr 2014 - J. E. Williams - terpene emissions introduced ! ! !REMARKS: ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/emission_terp_apply' integer, dimension(6) :: idater real :: dtime, fraction integer :: imr, jmr, lmr, lsec integer :: seccount2d, seccount3d type(d3_data) :: emis3d real :: dlatr integer :: day, ntim, i1, i2, j1, j2 ! --- begin ----------------------------------------- status = 0 if(.not. has_terp_emis) return if( okdebug ) then write(gol,*) 'start of emission_apply_terp'; call goPr end if ! --- WHERE, WHEN CALL TAU2DATE( itaur(region), idater) CALL GET_DISTGRID( dgrid(region), I_STRT=i1, I_STOP=i2, J_STRT=j1, J_STOP=j2 ) day = GET_DAY(idater(2),idater(3),mlen) dtime = float(nsrce)/(2*tref(region)) !emissions are added in two steps...XYZECCEZYX. if(okdebug) then write(gol,*)'emission_terp_apply in region ',region,' at date: ',idater, ' with time step:', dtime call goPr end if ! --- DIURNAL CYCLE if ( day /= day_terp(region) ) then dlatr = dy/yref(region) ! should be the same as lli%dlat_deg dtime=float(ndyn_max)/(2*tref(region)) ! timestep emissions ntim=86400/nint(dtime) ! number of timesteps in 24 hours. if ( okdebug ) then write (gol,*) 'day',day,'ntim',ntim,'jm',jm(region),'ybeg',& ybeg(region),'dlatr',dlatr,'sum(mlen)',sum(mlen),nsrce,dtime call goPr end if ! NOTE : latitude_start must be an integer. This is ok as long as ! the regions ybeg are integer. call scale_terp(day, ntim, j1, terp_dat(region)%scalef_isop, ybeg(region), dlatr, sum(mlen)) day_terp(region) = day !avoid recalculation during this day! endif ! --- ADD EMISSIONS ! get a working structure for 3d emissions allocate( emis3d%d3(i1:i2,j1:j2,lm(region)) ) ; emis3d%d3 = 0.0 ! default: no additional splitting fraction = 1.0 ! count 2d and 3d sectors seccount2d = 0 seccount3d = 0 ! cycle over sectors do lsec = 1, numb_sectors ! screen out unused providers if (count(used_providers_terp.eq.sectors_def(lsec)%prov).eq.0) cycle ! add up contribution from all proc ! screen out sectors which do not contain terpene emissions ! if ( sectors_def(lsec)%prov == 'AR5 ' .and. sectors_def(lsec)%catname /= 'biomassburning') cycle if ( sectors_def(lsec)%prov == 'ED41 ' .or. & sectors_def(lsec)%prov == 'ED42 ' .or. & sectors_def(lsec)%prov == 'HYMN ') cycle ! distinguish between 2d/3d sectors if( sectors_def(lsec)%f3d ) then seccount3d = seccount3d + 1 if (.not.has_data_3d(seccount3d)) cycle emis3d%d3 = terp_emis_3d(region,seccount3d)%d3 else seccount2d = seccount2d + 1 if (.not.has_data_2d(seccount2d)) cycle emis3d%d3 = 0.0 ! vertically distribute according to sector call emission_vdist_by_sector( sectors_def(lsec)%vdisttype, 'TERP', region, terp_emis_2d(region,seccount2d), & emis3d, status ) IF_NOTOK_RETURN(status=1) endif ! add dataset according to sector and category if( trim(sectors_def(lsec)%catname) == "natural" ) then call do_add_terp( region, i1, j1, emis3d%d3, status) ! special routine IF_NOTOK_RETURN(status=1) elseif ( trim(sectors_def(lsec)%catname) == "biomassburning" ) then if( emis_bb_trop_cycle) then call do_add_3d_cycle( region, iterp, i1, j1, emis3d%d3, bb_cycle(region)%scalef, xmterp, xmterp, status, fraction ) IF_NOTOK_RETURN(status=1) else call do_add_3d( region, iterp, i1, j1, emis3d%d3, xmterp, xmterp, status) IF_NOTOK_RETURN(status=1) end if endif end do deallocate( emis3d%d3 ) if(okdebug) then write(gol,*) 'end of emission_terp_apply' ; call goPr end if ! OK status = 0 END SUBROUTINE EMISSION_TERP_APPLY !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: SCALE_TERP ! ! !DESCRIPTION: calculates factors to be multiplied with an emission field ! (e.g. terp) in order to simulate a diurnal cycle in ! emissions (caused by solar dependency), e.g. : ! ! rm(i,j,1,terp) = rm(i,j,1,terp) + e_terp(i,j)*scalef(ipos,j)*dt ! ! with ipos depending on time and longitude. ! ! The scalefactor is calculated for -180 longitude and the ! mean value for ntim timesteps during a day is scaled to 1. ! The routine should be called every day, since the position ! relative to the sun changes. The scaling is based on the ! cos(zenith) which is 1 for overhead sun and 0 at sunset. ! For the polar night, the values are set to 1.0 (even distribution) ! One thing is sure afterwards: we know that e_terp kg/s is added. !\\ !\\ ! !INTERFACE: ! SUBROUTINE SCALE_TERP(day, ntim, j1, scalef, lat_start, dlat, idayy) ! ! !INPUT PARAMETERS: ! integer, intent(in) :: day ! day number based on 365 days a year integer, intent(in) :: ntim ! number of timesteps during one day integer, intent(in) :: j1 ! start index for latitudes integer, intent(in) :: lat_start ! southern edge of the domain (integer!) real, intent(in) :: dlat ! increment (real!) integer, intent(in) :: idayy ! # days this year ! ! !OUTPUT PARAMETERS: ! real, dimension(:,j1:), intent(out) :: scalef ! ! !REVISION HISTORY: ! 10 May 2012 - P. Le Sager - adapted for lon-lat MPI domain decomposition ! ! !REMARKS: ! !EOP !------------------------------------------------------------------------ !BOC real :: pi, piby, obliq, deday, delta, dt, lon, hr, lat, houra, smm,scale_test integer :: i, j pi = acos(-1.0) piby = pi/180.0 obliq = 23.45 * piby deday = 4.88 + 2*pi/idayy * day delta = asin(sin(obliq)*sin(deday)) dt = 24./ntim ! timestep in hours lon = -180.0*piby ! calculate for dateline do j=j1,ubound(scalef,2) ! over latitudes hr = - 0.5*dt ! shift half a interval back lat = (lat_start + (j-0.5)*dlat)*piby smm = 0.0 do i=1,ntim hr = hr + dt houra = lon - pi + hr * (2.*pi/24.) !checking if the sun is above horizon scale_test = max(sin(delta)*sin(lat) + cos(delta)*cos(lat)*cos(houra),0.0) scalef(i,j) = cos(houra) !smm = smm+scalef(i,j)/ntim smm = smm+scale_test/ntim end do if ( smm > 1e-5 ) then scalef(1:ntim,j) = scalef(1:ntim,j)+1 else scalef(1:ntim,j) = 1.0 end if end do END SUBROUTINE SCALE_TERP !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: DO_ADD_terp ! ! !DESCRIPTION: Add local time splitted emissions to chemical array. The ! difference with "do_add_3d" is the number of levels this is ! applied to: 4 instead of all levels. !\\ !\\ ! !INTERFACE: ! SUBROUTINE DO_ADD_TERP(region, i1, j1, em_terp, status) ! ! !USES: ! use dims, only : tref, nsrce, sec_month,dx,xref,xbeg, ndyn_max,itaur use global_data, only : mass_dat use chem_param, only : iterp, xmterp #ifdef with_budgets use budget_global, only : budg_dat, nzon_vg use emission_data, only : sum_emission, budemi_dat #endif use datetime, only : tau2date ! ! !INPUT PARAMETERS: ! integer, intent(in) :: region, i1, j1 real, intent(in) :: em_terp(i1:,j1:,:) ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: status ! ! !REVISION HISTORY: ! 03 Apr 2014 - J. E. Williams - copied from do_add_isop and adapted for terpenes ! ! !REMARKS: ! - may be be wise to add a level_limit keyword to the original "do_add_3d" ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/do_add_terp' real :: dtime, month2dt, x, xlon, dtime2 integer :: i, j, l, nzone, nzone_v, i2, j2 integer :: sec_in_day,ntim,itim,ipos integer,dimension(6)::idater real, dimension(:,:,:,:), pointer :: rm ! --- begin ----------------------------------------- rm => mass_dat(region)%rm CALL GET_DISTGRID( dgrid(region), I_STOP=i2, J_STOP=j2 ) dtime = float(nsrce) / (2*tref(region)) ! timestep emissions dtime2 = float(ndyn_max) / (2*tref(region)) ! timestep emissions based on non-CFL interference (CMK 05/2006) month2dt = dtime/sec_month ! conversion to emission per timestep (CFl induced timestep) ntim = 86400/nint(dtime2) ! number of timesteps in 24 hours based on non-CFL interference (CMK 05/2006) call tau2date(itaur(region),idater) sec_in_day = idater(4)*3600 + idater(5)*60 + idater(6) ! elapsd seconds this day itim = sec_in_day/nint(dtime2)+1 ! # time interval CMK 05/2006 ! ! Terpene emissions only at the surface ! do l= 1,4 do j=j1,j2 do i=i1,i2 xlon = xbeg(region) + (i-0.5)*dx/xref(region) ! itim = 1 and lon = -180 --->position 1 ! itim = ntim ant lon = -180 --->position ntim, etc. ! itim = 1 and lon = 0 ---->position ntim/2 ipos = 1 + mod(int((xlon+180.)*ntim/360.) + (itim-1),ntim) !position in array depending on time and lon. x=em_terp(i,j,l)*terp_dat(region)%scalef_isop(ipos,j)*month2dt !x=em_terp(i,j,l)*month2dt !write(1111,*) terp_dat(region)%scalef_isop(ipos,j),x ! x = em_terp(i,j,l)*month2dt rm(i,j,l,iterp)=rm(i,j,l,iterp)+x #ifdef with_budgets nzone=budg_dat(region)%nzong(i,j) !global budget nzone_v=nzon_vg(l) budemi_dat(region)%emi(i,j,nzone_v,iterp) = & budemi_dat(region)%emi(i,j,nzone_v,iterp) + x/xmterp*1e3 if(iterp ==1) sum_emission(region) = sum_emission(region) + x !in kg #endif enddo enddo enddo nullify(rm) status=0 END SUBROUTINE DO_ADD_TERP !EOC END MODULE EMISSION_TERP