! #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_DATA ! ! !DESCRIPTION: Provides general methods and fields needed for the emission ! routines. ! ! emission_vdist_by_sector : distribute vertically emissions according to source sector ! msg_emis : print formated emission total ! do_add_2d : add 2D emission to tracer mass array ! do_add_3d : add 3D emission to tracer mass array ! do_add_3d_cycle : add 3D emission to tracer mass array with a daily cycle ! scale_cycle : evaluate a daily cycle ! distribute_l44 : scale input field with NH (SH) fires distributed over JAS (JFM) !\\ !\\ ! !INTERFACE: ! MODULE EMISSION_DATA ! ! !USES: ! use GO, only : gol, goPr, goErr use Dims, only : nregions, lm use global_types, only : emis_data, d2_data, d23_data, d3_data #ifdef with_m7 use mo_aero_m7, only : nmod use mo_aero, only : cmr_ff, cmr_bb, cmr_sk, cmr_sa, cmr_sc, facso2 use mo_aero, only : zbb_wsoc_perc, zbge_wsoc_perc !use mo_aero, only : zom2oc use global_types, only : modal_emissions #endif #ifdef with_budgets use budget_global,only : budg_dat, nzon_vg #endif IMPLICIT NONE ! ! !PUBLIC DATA MEMBERS & TYPES: ! public :: emis_data, d2_data, d23_data ! type(emis_data), target :: plandr(nregions) ! Land fraction type(emis_data), target :: emis2D(nregions) ! 2D emiss prior vertical redistribution type(emis_data), target :: emis_temp(nregions) ! temp array to store emission #ifdef with_m7 type(modal_emissions), dimension(nregions,nmod), target :: emis_mass type(modal_emissions), dimension(nregions,nmod), target :: emis_number #endif #ifdef with_m7 ! Count median / geometric mean radii of ! primary carbonaceous and sulfate emissions ! based on AeroCom-I recommendations (Dentener et al., ACP, 2006). ! The corresponding values for the M7 modes are given by Stier et al. (ACP, 2005). ! The same values have also been adopted in GLOMAP (Mann et al., 2010). ! ! Count median radii for carbonaceous aerosol emissions from Dentener et al., ! corresponding to sigma = 1.8: !real, parameter :: rad_emi_ff = 0.015e-6 !real, parameter :: rad_emi_bb = 0.04e-6 ! The corresponding values for sigma = 1.59 used in M7 ! are cmr_ff and cmr_bb, specified in mo_aero.F90. ! For comparison, Bond et al. (JGR, 2013) give number median radii ! between 25 and 40 nm for fresh BC in the urban areas ! of Tokyo, Nagoya, and Seoul, ! of 60 nm in plumes associated with wildfires, ! and about 15 nm from aircraft jet engines. ! These values are volume-equivalent radii (see their Fig. 4). ! ! According to the original paper by Schwarz et al. (GRL, 2008) ! the corresponding geometric standard deviation ! is sigma = 1.71 for the urban BC ! and 1.43 for the biomass burning aerosol. ! ! For BC in biomass burning plumes, ! Kondo et al (JGR, 2011) estimated ! number median radii in the range 68-70.5 nm (+- 6-8 nm) ! and geometric standard deviation between 1.32 and 1.36 (+- 0.01-0.04), ! for particles thickly coated by organics. ! ! Janhaell et al. (ACP, 2010) have compiled measurements of ! particle size in fresh biomass burning smoke from vegetation fires. ! They mention that particles from biomass burning are dominated ! by an accumulation mode. ! They also present a relation between the geometric mean diameter Dg ! and geometric standard deviation sigma for fresh smoke: ! Dg (um) = (584 +- 5) - (269 +-1) sigma ! This gives a geometric mean radius of 78 um for sigma = 1.59, ! in close agreement with the value used by Stier et al. ! ! Another way to account for differences in sigma ! is to modify the number median radius rg ! such that the number of particles (N) ! emitted for a certain mass (M) is the same. ! N is proportional to M/rv^3, ! where rv is the volume mean radius, ! which is related to rg by ! rv = rg * exp(1.5*(ln(sigma))^2). ! According to Janhall, fresh smoke has an average ! Dg of 117 +- 13 nm and sigma of 1.7 +- 0.1. ! At sigma=1.59, this would translate into Dg of about 129 nm, ! or rg of about 65 nm. ! For the estimates from Kondo et al. and Schwarz et al., ! this would give somewhat smaller rg values ! of about 58 and 53 nm, resp. ! ! Particles emitted by grass and savannah fires are generally ! somewhat smaller than those from wood burning. ! Janhaell et al. estimate that the mean emission radii ! for grass and savannah fires, resp., are 12.5 and 10 nm smaller. ! These differences is not accounted for in the model. ! ! In a later version of ECHAM-HAM particles the emission radius ! for biomass burning was reduced to the value for fossil fuel ! (Zhang et al., ACP, 2012). ! However, such as a small value seems inconsistent with measurements. ! ! In the CMIP6 emission data set, ! the contributions from solid biofuel combustion ! are included in the 2-D anthropogenic sectors, ! and provided separately in a supplementary data set. ! In CMIP6, biofuel is only non-zero ! for the energy, industry and residential and ! transportation sectors. ! Most of the residential emissions are due ! to solid biofuel combustion. ! The contribution to the industrial sector ! is only substantial in developing countries, ! while the contributions to the energy ! and transportation sectors are generally very small. ! ! For inventories other than CMIP6, ! no distinction between fossil and biofuel emissions ! is made in the model. ! ! Winijkul et al. (Atm. Env., 2015) have measured ! size distributions from energy-related combustion ! sources for the residential, industrial, power and ! transportation sectors. ! They give regional and global estimates of ! mass median diameters. ! In their supplementary material they give ! an overview of results from other studies. ! These results indicate that the size distribution ! for both biofuel and fossil fuel combustion sources ! are strongly dependent on the technique. ! Count median radii for the residential sector, ! presented in their Table S1, vary between about ! 15 nm for modern (improved) wood-fueled stoves, ! and 25-30 for fireplaces, ! to 160 nm for (regular) wood-fueled heating stoves, ! to 270 nm for traditional cookstoves. ! General fossil fuel (2-d sectors) ! emitted in the Aitken mode real, parameter :: rad_emi_ff_sol = cmr_ff real, parameter :: rad_emi_ff_insol = cmr_ff ! Energy sector ! emitted in the Aitken mode. real, parameter :: rad_emi_ene_sol = cmr_ff real, parameter :: rad_emi_ene_insol = cmr_ff ! Industry sector ! emitted in the Aitken mode. real, parameter :: rad_emi_ind_sol = cmr_ff real, parameter :: rad_emi_ind_insol = cmr_ff ! Transportation sector ! emitted in the Aitken mode. real, parameter :: rad_emi_tra_sol = cmr_ff real, parameter :: rad_emi_tra_insol = cmr_ff ! Shipping sector ! emitted in the Aitken mode. ! Currently set to cmr_ff, ! but could be changed. real, parameter :: rad_emi_shp_sol = cmr_ff real, parameter :: rad_emi_shp_insol = cmr_ff ! Aircraft sector ! emitted in the Aitken mode. ! Currently set to cmr_ff, ! but a smaller value could be used, ! e.g. 10 or 15 nm. real, parameter :: rad_emi_air_sol = cmr_ff real, parameter :: rad_emi_air_insol = cmr_ff ! Open biomass burning ! Soluble part emitted in the accumulation mode, ! insoluble part in the Aitken mode. ! Stier et al. apply cmr_ff to the insoluble part, ! but a slightly larger value seems more realistic, ! e.g. 40 nm. real, parameter :: rad_emi_bb_sol = cmr_bb real, parameter :: rad_emi_bb_insol = cmr_ff ! Solid biofuel combustion ! Currently treated as biomass burning, ! but could be changed (see comments above). ! The value cmr_bb = 75 nm corresponds ! with the value of 50 nm at sigma=2 ! assumed by Kondros et al. (ACP, 2015) ! for emissions from biofuel combusion ! in their baseline run. ! In one of their sensitivity runs, ! they increase it by a factor of 2. real, parameter :: rad_emi_bf_sol = cmr_bb real, parameter :: rad_emi_bf_insol = cmr_ff ! not used anymore: !real, parameter :: rad_soa = 0.01e-6 ! soa average radius ! assuming 3nm particle formation and growth to ! that size in half an hour ! Count median radii for sulfate aerosol emissions adapted to the M7 modes: real, parameter :: rad_so4_ait = cmr_sk ! aitken mode radius real, parameter :: rad_so4_acc = cmr_sa ! accumulation mode radius real, parameter :: rad_so4_coa = cmr_sc ! coarse mode radius ! Count median dry radii for sea salt emissions. ! These values have been updated ! following Vignati et al. (Atmos. Environ., 2010) ! For further explanations, see emission_ss.F90. !real, parameter :: radius_ssa = 0.0794e-6 !real, parameter :: radius_ssc = 0.63e-6 real, parameter :: radius_ssa = 0.09e-6 ! accumulation mode real, parameter :: radius_ssc = 0.794e-6 ! coarse mode ! Soluble fraction of POM mass ! According to Janhall et al. (ACP, 2010) ! 40 to 80% of the organic matter from ! vegetation fires is water soluble. ! Kondros et al. (ACP, 2015) use 80% ! for POM from biofuel combustion ! in their base run. ! For the moment, emisions from ! biofuel combustion and open biomass burning ! are treated in the same way. ! ! open biomass burning real, parameter :: frac_pom_sol_bb = zbb_wsoc_perc ! solid biofuel combustion real, parameter :: frac_pom_sol_bf = zbb_wsoc_perc !real, parameter :: frac_pom_sol_bf = 0.8 ! alternative value ! fossil fuel combustion !real, parameter :: frac_pom_sol_ff = 0.65 ! original value real, parameter :: frac_pom_sol_ff = 0.0 ! value since 2015 revision ! Soluble fraction of BC mass ! In the original code, ! fresh BC was assumed 100% insoluble, ! and therefore emitted into the Aitken mode. ! The new code allows to use a non-zero fraction, ! to account for non-resolved ageing close to the source. ! See e.g. Kondros et al., ! who use 50% for biofuel combustion ! in their base run. ! In the standard EMEP MSC-W model, ! 20% of the elemental carbon (EC) from ! anthropogenic sources ! (representative of fossil fuel combustion) ! and all EC from open biomass fires ! is assumed hygroscopic. ! ! open biomass burning !real, parameter :: frac_bc_sol_bb = 0.0 ! original value !real, parameter :: frac_bc_sol_bb = 0.5 real, parameter :: frac_bc_sol_bb = 0.95 ! To reduce the AOD over china and outflow region of ! Africa the water soluble fraction was increasde to 95% ! in preparation for CMIP6. ! solid biofuel combustion !real, parameter :: frac_bc_sol_bf = 0.0 ! original value !real, parameter :: frac_bc_sol_bf = 0.5 ! Aerocom real, parameter :: frac_bc_sol_bf = 0.95 ! TB: ! To reduce the AOD over china and outflow region of ! Africa the water soluble fraction was increasde to 95% ! in preparation for CMIP6. ! ! Some basiss for the choice can be found here: ! (e.g. Janhall et al., 2010; ! https://doi.org/10.5194/acp-10-1427-2010 ; Winijkul et al., 2015; ! https://doi.org/10.1016/j.atmosenv.2015.02.037; Li et al., 2009; ! https://pubs.acs.org/doi/abs/10.1021/es803330j). ! fossil fuel combustion real, parameter :: frac_bc_sol_ff = 0.0 ! Soluble fraction of surrogate SOA emissions ! POM from SOA is considered 65% soluble, ! as recommended by AeroCom. ! The paper by Kanakidou et al. (ACP, 2004) ! (mentioned in emission_pom.F90) ! does not seem to support 100% solubility, ! as was previously assumed. !real, parameter :: frac_soa_sol = 1.0 ! original value real, parameter :: frac_soa_sol = zbge_wsoc_perc ! value since 2015 revision ! Fraction of SOx mass emitted directly as sulfate real, parameter :: frac_so4=1.-facso2 #else real, parameter :: frac_so4=1.-0.975 #endif #ifdef with_m7 ! The value of 1.4 for the POM to OC mass ratio, ! set in mo_aero.F90, is an outdated estimate, ! see e.g. Turpin and Lim (Aerosol Sci. Technol., 2001) and ! Aiken et al. (Environ. Sci. Technol., 2008). ! Turpin and Lim estimate a ratio of 1.6 +- 0.2 for urban aerosol, ! and 2.1 +- 0.2 for aged (nonurban) aerosol; ! They also note that aerosols heavily impacted by woodsmoke can ! have an even higher ratio (2.2 to 2.6). ! According to Reid et al. (ACP, 2005), ! the POM to OC ratio in fresh biomass burning smoke is very uncertain, ! somewhere in between 1.4 and ~2. ! Aiken et al. measure ambient aerosol values ! between 1.25 for O/C = 0 to 2.44 for O/C = 1.0. ! For OA from biomass burning, they measure 1.56-1.70, ! lower than the estimates from Turpin and Lim (~2.0). ! They find the highest ratios for ! aged and freshly formed SOA (~2.4 and ~1.9, respectively) ! and lowest values for primary OA from urban combustion. ! Based on these studies, we apply different values ! for emissions from biomass burning versus other emissions, ! as is also done in some other models ! (see Table 1 in Tsigaridis et al., ACP, 2014), ! For SOA (seem emission_pom.F90) ! we apply a relatively high value valid for aged SOA. ! This compensates for the lack of SOA formation from isoprene, ! and improves the agreement with aerosol optical depth (AOD) ! derived from satellite observations (MODIS). ! ! We could go even further and apply different values for the ! water soluble and insoluble fractions (Turpin and Lim). ! ! It should be acknowledged that the representation of OA ! with a single tracer is very simplistic. ! In particular, increase in OA mass due to ageing ! is not properly accounted for. ! !real, parameter :: oc2pom = zom2oc !factor for conversion of OC mass to POM real, parameter :: oc2pom_ff = 1.6 ! fossil fuel real, parameter :: oc2pom_bf = 1.6 ! solid biofuel combustion real, parameter :: oc2pom_bb = 1.6 ! open biomass burning real, parameter :: oc2pom_soa = 2.4 ! SOA #endif ! CMIP6 - AR5 - EDGAR4 - RETRO - GFED3 - MACC - LPJ - HYMN - MEGAN logical :: LCMIP6, LCMIP6BMB, LCMIP6_CH4 logical :: LAR5, LAR5BMB, LEDGAR4, LRETROF, LGFED3, LMACC, LLPJ, LHYMN, LMACCITY, LMEGAN character(len=256) :: cmip6_ch4_dirname character(len=256) :: emis_input_dir character(len=256) :: emis_input_dir_cmip6 character(len=256) :: emis_input_dir_ar5 character(len=256) :: emis_input_dir_megan character(len=256) :: emis_input_dir_ed4 character(len=256) :: emis_input_dir_mac integer :: emis_input_year_nat integer :: emis_input_year_bc integer :: emis_input_year_oc integer :: emis_input_year_sox integer :: emis_input_year_nh3 integer :: emis_input_year_nox integer :: emis_input_year_co integer :: emis_input_year_nmvoc integer :: emis_input_year_ch4 character(len=256) :: emis_input_dir_gfed character(len=256) :: emis_input_dir_retro character(len=256) :: emis_input_dir_aerocom character(len=256) :: emis_input_dir_dms character(len=256) :: emis_input_dir_rn222 character(len=256) :: emis_input_dir_dust, emis_input_dust #ifdef with_ch4_emis character(len=256) :: emis_input_dir_natch4 #endif logical :: emis_ch4_single logical :: emis_ch4_fix3d real :: emis_ch4_fixed_ppb character(len=256) :: emis_zch4_fname logical :: ch4_fixyear integer :: ch4_year !------------------------------ ! SSP scenario name !------------------------------ character(len=14) :: ssp_name ! biomass burning levels: integer, parameter :: bb_nlev = 5 integer, parameter :: bb_lm = lm(1) ! to be reduced to strip stratosphere! integer, parameter :: bl_lm = lm(1) ! to be reduced to the BL (around 8) real, parameter :: bb_hlow (0:bb_nlev) = (/ 0., 100., 500., 1000., 2000., 3000./) real, parameter :: bb_hhigh(0:bb_nlev) = (/ 100., 500., 1000., 2000., 3000., 6000./) ! biomass burning levels in tropics (-20 - 20):: integer, parameter :: bb_nlev_trop = 3 real, parameter :: bb_hlow_trop (0:bb_nlev_trop) = (/ 0., 100., 500., 1000./) real, parameter :: bb_hhigh_trop(0:bb_nlev_trop) = (/ 100., 500., 1000., 2000./) ! ----------------------------------------------------------------------------------- ! Biomassburning time splitting factors (now globally, instead of by constituent) logical :: emis_bb_trop_cycle type bb_cycle_data real, pointer :: scalef(:) end type bb_cycle_data type(bb_cycle_data), dimension(nregions), target :: bb_cycle ! ----------------------------------------------------------------------------------- ! Budgets #ifdef with_budgets real, dimension(nregions) :: sum_emission type budemi_data real,dimension(:,:,:,:),pointer :: emi ! [i, j, nbud_vg, ntracet] end type budemi_data type(budemi_data),dimension(nregions),target :: budemi_dat #endif logical :: use_tiedkte ! read from rc file, used to scale LiNOx ! ---------------------------------------------------- ! Vertical splitting : used now for all sectors, for biomassburning and other categories ! (tables in Dentener, 2006, ACP) ! ! ATTENTION: changes here have impact on the routine emission_vdist_by_sector (contained) ! Note: vdist_* variable could be hold private ! integer, parameter :: vd_class_name_len = 12 ! biomassburning integer, parameter :: vdist_bb_nlev = 6 real,dimension(vdist_bb_nlev),parameter :: vdist_bb_hlow = (/ 0., 100., 500., 1000., 2000., 3000./) real,dimension(vdist_bb_nlev),parameter :: vdist_bb_hhigh = (/ 100., 500., 1000., 2000., 3000., 6000./) ! other sources (related to surface) integer, parameter :: vdist_nn_nlev = 6 real,dimension(vdist_nn_nlev),parameter :: vdist_nn_hlow = (/ 0., 30. , 100., 300., 600., 1000. /) real,dimension(vdist_nn_nlev),parameter :: vdist_nn_hhigh = (/ 30., 100., 300., 600., 1000., 2000. /) ! ! !PRIVATE TYPES: ! character(len=*), parameter, private :: mname = 'emission_data' ! ! !REVISION HISTORY: ! 1 Oct 2010 - Achim Strunk - adapted for AR5 ! 28 Jun 2012 - Ph. Le Sager - made everything public by default ! - adapted for lon-lat MPI domain decomposition ! !REMARKS: ! !EOP !------------------------------------------------------------------------ CONTAINS !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: EMISSION_VDIST_BY_SECTOR ! ! !DESCRIPTION: Return vertically distributed emissions depending on given ! sector and constituent. ! Distribution is done as 'compromise' between ! 1) Dentener et al., 2006, ACP, slightly modified for ! emissions supposed to enter the surface layer. ! 2) De Meij et al., ACP, 2006 ! 3) EMEP model (docu, model code, publications) ! 4) Bieser et al., GMDD, 2010 ! 5) Small updates to biomassburning following ! Huijnen et al., GMDD, 2010 !\\ !\\ ! !INTERFACE: ! SUBROUTINE EMISSION_VDIST_BY_SECTOR( splitname, constituent, region, emis2D, emis3D, status ) ! ! !USES: ! use toolbox, only : distribute_emis2D ! ! !INPUT PARAMETERS: ! character(len=vd_class_name_len), intent(in) :: splitname character(len=*), intent(in) :: constituent integer, intent(in) :: region type(emis_data), intent(in) :: emis2D ! ! !INPUT/OUTPUT PARAMETERS: ! type(d3_data), intent(inout) :: emis3D ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: status ! ! !REVISION HISTORY: ! 1 Oct 2010 - Achim Strunk - v0 ! ! !REMARKS: Splitting depending on constituent is still to code. ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/emission_vdist_by_sector' integer :: ilev character(len=2) :: splittype integer, parameter :: maxvdist_nlev = max(vdist_bb_nlev,vdist_nn_nlev) real, dimension(maxvdist_nlev) :: sfglobal, sfintrop, sfintemp ! --- begin -------------------------------------- ! control implicit settings here in this routine if( vdist_bb_nlev /= 6 .or. vdist_nn_nlev /= 6 ) then write(gol, '("ERROR: wrong number of layers: bb_nlev /= 6 or nn_lev /= 6")'); call goErr status=1; return end if ! initialise sfglobal = 0.0 sfintrop = 0.0 sfintemp = 0.0 ! ------------------------------------------------- ! This here is selection by emission source sectors ! ------------------------------------------------- select case( trim(splitname) ) case( 'combenergy' ) ! Combustion from energy production and transformation, power-plants ! no zonal differences, no dependency on species ! -- ! assumed to be related to EMEP SNAP category [1] ! AEROCOM equivalent: power-plants ! -- splittype = 'nn' sfglobal = (/ 0.0, 0.1, 0.7, 0.2, 0.0, 0.0 /) case( 'combrescom' ) ! Residential and commercial combustion ! no zonal differences, no dependency on species ! -- ! assumed to be related to EMEP SNAP category [2] ! AEROCOM equivalent: domestic/industry ! -- splittype = 'nn' sfglobal = (/ 0.4, 0.4, 0.2, 0.0, 0.0, 0.0 /) case( 'industry' ) ! Industrial processes and combustion ! no zonal differences, no dependency on species ! -- ! assumed to be an aggregate of EMEP SNAP categories [3,4,5] ! AEROCOM equivalent: industry ! -- splittype = 'nn' sfglobal = (/ 0.1, 0.2, 0.6, 0.1, 0.0, 0.0 /) case( 'waste' ) ! Waste treatment and disposal ! no zonal differences, no dependency on species ! -- ! assumed to be related to EMEP SNAP category [9] ! AEROCOM equivalent: - ! -- splittype = 'nn' sfglobal = (/ 0.1, 0.2, 0.4, 0.3, 0.0, 0.0 /) case( 'nearsurface' ) ! Near surface emissions ! no zonal differences, no dependency on species ! -- ! AR5: solvent production and use, agricultural waste burning, & ! ships, grassland fire, mineral dust (AEROCOM or Tegen-Vignati) ! EMEP equivalents: SNAP categories [6,8] ! AEROCOM equivalents: ships, agricultural waste ! -- splittype = 'nn' sfglobal = (/ 0.8, 0.2, 0.0, 0.0, 0.0, 0.0 /) case( 'surface' ) ! Surface emissions ! no zonal differences, no dependency on species ! -- ! AR5: agriculture, transport, soil, oceanic, biogenic ! EMEP equivalents: SNAP categories [7,8,10] ! AEROCOM equivalents: surface emissions from road/off-road, ships, agricultural waste ! -- splittype = 'nn' sfglobal = (/ 1.0, 0.0, 0.0, 0.0, 0.0, 0.0 /) !!$ case ( 'mindust' ) !!$ ! Mineral dust emissions !!$ ! no zonal differences !!$ ! -- !!$ ! AR5: surface dust emissions from either from AEROCOM or Tegen-Vignati !!$ ! assumed to be slightly higher than near-surface !!$ splittype = 'nn' !!$ sfglobal = (/ 0.8, 0.2, 0.0, 0.0, 0.0, 0.0 /) case( 'volcanic' ) ! Volcanic emissions ! -- ! AR5: natural SOx emissions (from MACC repository) ! EMEP equivalents: SNAP categories [11] ! AEROCOM equivalents: continuous: 2/3 to 1 of height of volcano top ! explosive: 0.5 to 1.5 km above volcano top ! Since no data is available about volcano top levels, we assume here, that ! a volcano top is usually 200-1000 m higher than the grid's surrounding terrain ! height, thus emissions are distributed from 100 - 2000 m almost homogeneously ! -- splittype = 'nn' sfglobal = (/ 0.0, 0.0, 0.1, 0.3, 0.4, 0.2 /) case( 'forestfire' ) ! forest fires: 6 layer model, different splitting for different regions ! dont distinguish between boreal Europe and boreal Canada, use Europe globally ! Dentener et al., modified in tropics (up to 2 km) with respect to Huijnen et al., GMD 2010, \ ! who are citing Labonne et al., 2007 ofr new information based on satellite data. splittype = 'bb' sfglobal = (/ 0.1, 0.1, 0.2, 0.2, 0.4, 0.0 /) sfintemp = (/ 0.2, 0.2, 0.2, 0.4, 0.0, 0.0 /) sfintrop = (/ 0.2, 0.2, 0.2, 0.4, 0.0, 0.0 /) case default write(gol, '("ERROR: wrong specification of emission sector ",a,"in vdist_by_sector called for constituent ",a)') & trim(splitname), trim(constituent); call goErr status=1; return end select ! -------------------- ! Do the splitting - ! no need to skip cases of fraction (sfglobal, sfintemp, ..) being zero: ! distribute_emis2D is optimized for those and will return right away ! -------------------- select case( splittype ) case( 'nn' ) do ilev = 1, vdist_nn_nlev call distribute_emis2D( region, emis2D, emis3D, vdist_nn_hlow(ilev), vdist_nn_hhigh(ilev), status, sfglobal(ilev) ) IF_NOTOK_RETURN(status=1) end do case( 'bb' ) do ilev = 1, vdist_bb_nlev ! boreal call distribute_emis2D( region, emis2D, emis3D, vdist_bb_hlow(ilev), vdist_bb_hhigh(ilev), status, & sfglobal(ilev), -90.,-61.) IF_NOTOK_RETURN(status=1) call distribute_emis2D( region, emis2D, emis3D, vdist_bb_hlow(ilev), vdist_bb_hhigh(ilev), status, & sfglobal(ilev), 60., 90.) IF_NOTOK_RETURN(status=1) ! temperate call distribute_emis2D( region, emis2D, emis3D, vdist_bb_hlow(ilev), vdist_bb_hhigh(ilev), status, & sfintemp(ilev), -60.,-31.) IF_NOTOK_RETURN(status=1) call distribute_emis2D( region, emis2D, emis3D, vdist_bb_hlow(ilev), vdist_bb_hhigh(ilev), status, & sfintemp(ilev), 30., 59.) IF_NOTOK_RETURN(status=1) ! tropics call distribute_emis2D( region, emis2D, emis3D, vdist_bb_hlow(ilev), vdist_bb_hhigh(ilev), status, & sfintrop(ilev), -30., 29.) IF_NOTOK_RETURN(status=1) end do case default write(gol, '("ERROR: wrong specification of splitting type in vdist_by_sector.")'); call goErr status=1; return end select ! OK status=0 END SUBROUTINE EMISSION_VDIST_BY_SECTOR !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: MSG_EMIS ! ! !DESCRIPTION: Provide output to verify the amount of emissions ! entering the calculations !\\ !\\ ! !INTERFACE: ! SUBROUTINE MSG_EMIS(year_month, provider, sector, msg_comp, msg_mass, summed_emissions) ! ! !USES: ! use dims, only: idate ! ! !INPUT PARAMETERS: ! integer, intent(in) :: year_month ! 1: year, 2: monthly amount character(len=*), intent(in) :: provider, sector, msg_comp real, intent(in) :: msg_mass, summed_emissions ! ! !REVISION HISTORY: ! 1 Oct 2010 - Achim Strunk - protex; increase space for output ! !EOP !------------------------------------------------------------------------ !BOC character(len=7),dimension(2), parameter:: ym=(/'Yearly ','Monthly'/) 1111 format(' EMISS-INFO - ',a7,' Emission from ',a8,' -> ',a18,1x,a6,' : mw=',f5.1,', month=',i2.2,', Sum [kg]=',1x,1pe11.4) 1112 format(' EMISS-INFO - ',a7,' Emission from ',a8,' -> ',a18,1x,a6,' : mw=',f5.1,', --------, Sum [kg]=',1x,1pe11.4) if (year_month==1) then write (gol,1112) ym(year_month), provider, sector, msg_comp, msg_mass, summed_emissions; call goPr else write (gol,1111) ym(year_month), provider, sector, msg_comp, msg_mass, idate(2), summed_emissions; call goPr endif END SUBROUTINE MSG_EMIS !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: SCALE_CYCLE ! ! !DESCRIPTION: Evaluates a daily cycle to be used for the biomass ! burning emissions. It is based on the isoprene daily ! cycle, except that the max peak is at 14:00h LT. ! Calculates factors to be multiplied with an ! emission field (e.g. co) in order to simulate a ! diurnal cycle in emissions (caused by solar dependency/ ! human activity), e.g. ! rm(i,j,1,isop) = rm(i,j,1,isop) + e_isop(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 only once, since ! the position relative to the sun is not relevant here. !\\ !\\ ! !INTERFACE: ! SUBROUTINE SCALE_CYCLE(ntim, scalef) ! ! !INPUT PARAMETERS: ! integer, intent(in) :: ntim ! ! !OUTPUT PARAMETERS: ! real, dimension(ntim), intent(out) :: scalef ! ! !REVISION HISTORY: ! 1 Oct 2010 - Achim Strunk - added without_diurnal_emission_cycle ! (to be removed again :-( ) ! 19 Oct 2010 - Narcisa Banda - changed without_diurnal_emission_cycle to with_.... ! 31 Jan 2013 - Ph. Le Sager - get rid of with_diurnal_emission_cycle! ! ! !REMARKS: ! - This routine is called only once (see emission.F90), and only if ! emis_bb_trop_cycle is .true.. The later is set with the ! input.emis.bb.dailycycle in the rc file (chem.input.rc) ! !EOP !------------------------------------------------------------------------ !BOC real :: pi, piby, obliq, deday, delta, dt, lon, hr, lat, houra, smm, w, sig integer :: i ! The calling routine (emission_declare) defines NTIM as follows: ! dtime=float(ndyn_max)/(2*tref(region)) ! timestep emissions ! ntim=86400/nint(dtime) ! number of timesteps in 24 hours. pi = acos(-1.0) piby = pi/180.0 ! obliq = 23.45 * piby w = 0.2 sig = 2.0 dt = 24./ntim !timestep in hours lon = -180.0*piby !calculate for dateline hr = - 0.5*dt !shift half an interval back if (hr .gt. 12) hr = hr - 24 ! hr = hr - 2 !shift two hours. This means that Local-day-time maximum will be at 14:00 h LT. smm = 0.0 do i=1,ntim hr = hr + dt houra = lon - pi + hr * (2.*pi/24.) scalef(i) = w + (1-w)*24.0/(sig*sqrt(2*pi))*exp(-0.5*(((hr-1.5)/sig)**2)) ! scalef(i) = max(1+(cos(houra)),0.0) smm = smm+scalef(i)/ntim end do if ( smm > 1e-5 ) then scalef(1:ntim) = scalef(1:ntim)/smm else scalef(1:ntim) = 1.0 end if END SUBROUTINE SCALE_CYCLE !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: DO_ADD_2D ! ! !DESCRIPTION: General routine to add a 2D emission field (given in kg ! /month) into tracer mass array. The mass increase is done ! at level L_LEVEL, and for tracer I_TRACER. ! ! Update accordingly the budget. !\\ !\\ ! !INTERFACE: ! SUBROUTINE DO_ADD_2D( region, i_tracer, l_level, is, js, emis_field, & mol_mass, mol_mass_e, status, xfrac ) ! ! !USES: ! use dims, only : CheckShape, tref, nsrce use dims, only : sec_month, im, jm, okdebug use global_data, only : mass_dat, region_dat use chem_param, only : ntracet, names use tm5_distgrid, only : get_distgrid, dgrid #ifdef with_budgets use budget_global, only : budg_dat, nzon_vg #endif ! ! !INPUT PARAMETERS: ! integer, intent(in) :: region integer, intent(in) :: i_tracer integer, intent(in) :: l_level integer, intent(in) :: is, js real, intent(in) :: emis_field(is:,js:) real, intent(in) :: mol_mass ! mol mass of tracer field real, intent(in) :: mol_mass_e ! mol mass of emission field (e.g. NOx emission ar in kgN/month) real, intent(in), optional :: xfrac ! partial addition ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: status ! ! !REVISION HISTORY: ! 27 Mar 2012 - P. Le Sager - adapted for lon-lat MPI domain decomposition ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/do_add_2d' integer :: i, j, nzone, nzone_v, i1, i2, j1, j2 real :: x, efrac, dtime, month2dt real :: mbef, maft, sume ! mass BEFore, AFTer, SUM Emiss (debug) real, dimension(:,:,:,:), pointer :: rm, rzm ! --- begin -------------------------------------- call GET_DISTGRID( dgrid(region), I_STRT=i1, I_STOP=i2, J_STRT=j1, J_STOP=j2 ) call CheckShape( (/i2-i1+1,j2-j1+1/), shape(emis_field), status ) IF_NOTOK_RETURN(status=1) dtime=float(nsrce)/(2*tref(region)) ! timestep emissions month2dt=dtime/sec_month ! conversion to emission per timestep rm => mass_dat(region)%rm rzm => mass_dat(region)%rzm if(present(xfrac)) then efrac = xfrac else efrac = 1.0 endif if (okdebug) then write (gol, '("--------------------------------------------------------------")') ; call goPr write (gol, '(" DO ADD 2D debug in region:",i3) ' ) region ; call goPr write (gol, '(" itracer :",i3,a8)' ) i_tracer, ' :'//names(i_tracer) ; call goPr write (gol, '(" xfrac :",f6.1)' ) efrac ; call goPr write (gol, '(" mol_mass :",f6.1 ) ' ) mol_mass ; call goPr write (gol, '(" mol_mass_e :",f6.1 )' ) mol_mass_e ; call goPr write (gol, '(" emisfield :",2e12.4 )') minval(emis_field),maxval(emis_field) ; call goPr sume = sum(emis_field)*month2dt*(mol_mass/mol_mass_e)*efrac write (gol, '(" sum emmis :",e12.4 )' ) sume ; call goPr mbef = sum(rm(i1:i2,j1:j2,:,i_tracer)) endif do j=j1,j2 do i=i1,i2 x = emis_field(i,j)*month2dt*(mol_mass/mol_mass_e)*efrac rm(i,j,l_level,i_tracer) = rm(i,j,l_level,i_tracer) + x #ifdef slopes ! injection of emissions at surface, thus vertical slope more negative ! (keep concentration at top of layer the same as before emission) if ( i_tracer <= ntracet ) then rzm(i,j,l_level,i_tracer) = rzm(i,j,l_level,i_tracer) - x end if #endif #ifdef with_budgets nzone = budg_dat(region)%nzong(i,j) !global budget nzone_v = nzon_vg(l_level) #ifndef without_emission budemi_dat(region)%emi(i,j,nzone_v,i_tracer) = & budemi_dat(region)%emi(i,j,nzone_v,i_tracer) + x/mol_mass*1e3 if(i_tracer == 1) sum_emission(region) = sum_emission(region) + x !in kg #endif #endif enddo enddo if (okdebug) then maft = sum(rm(i1:i2,j1:j2,:,i_tracer)) write (gol, '(" Added amount :",e12.4 )' ) maft-mbef ; call goPr write (gol, '(" Added amount can be different when inner zoom is present!" )' ) ; call goPr if (maft-mbef-sume > 1e-8*sume) then write (gol, '(" Warning: error in emission ")' ) ; call goErr endif endif nullify(rm) nullify(rzm) status=0 END SUBROUTINE DO_ADD_2D !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: DO_ADD_3D ! ! !DESCRIPTION: General routine to add a 3D emission field (given in kg ! /month) into tracer mass array. The mass increase is done ! for tracer I_TRACER. ! ! Update accordingly the budget. !\\ !\\ ! !INTERFACE: ! SUBROUTINE DO_ADD_3D( region, i_tracer, is, js, emis_field, mol_mass, mol_mass_e, status, xfrac) ! ! !USES: ! use dims, only : CheckShape, nsrce, tref use dims, only : sec_month, im, jm, lm, okdebug use global_data, only : mass_dat, region_dat use chem_param, only : ntracet, names use tm5_distgrid, only : get_distgrid, dgrid #ifdef with_budgets use budget_global, only : budg_dat, nzon_vg #endif ! ! !INPUT PARAMETERS: ! integer, intent(in) :: region ! region # integer, intent(in) :: i_tracer ! id of tracer to increase integer, intent(in) :: is, js real, intent(in) :: emis_field(is:,js:,:) real, intent(in) :: mol_mass, mol_mass_e real, intent(in), optional :: xfrac !partial addition ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: status ! ! !REVISION HISTORY: ! 27 Mar 2012 - P. Le Sager - adapted for lon-lat MPI domain decomposition ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/do_add_3d' integer :: i, j, l, nzone, nzone_v, i1,j1,i2,j2 integer, dimension(3) :: ubound_e real :: x, efrac, dtime, month2dt real, dimension(:,:,:,:), pointer :: rm real :: mbef, maft, sume ! debug ! --- begin -------------------------------------- call GET_DISTGRID( dgrid(region), I_STRT=i1, I_STOP=i2, J_STRT=j1, J_STOP=j2 ) status = 0 ! check arguments ubound_e = ubound(emis_field) if(ubound_e(1) /= i2) then write(gol,'(A)') 'do_add_3d: im emission not consistent' ; call goPr status = 1 endif if(ubound_e(2) /= j2) then write(gol,'(A)') 'do_add_3d: jm emission not consistent' ; call goPr status = 1 endif if(ubound_e(3) > lm(region)) then write(gol,'(A)') 'do_add_3d: lm emission not consistent' ; call goPr status = 1 endif IF_NOTOK_RETURN(status=1) rm => mass_dat(region)%rm dtime=float(nsrce)/(2*tref(region)) ! timestep emissions month2dt=dtime/sec_month ! conversion to emission per timestep if(present(xfrac)) then efrac = xfrac else efrac = 1.0 endif ! if (okdebug) then ! write (gol, '("--------------------------------------------------------------")') ; call goPr ! write (gol, '(" DO ADD 3D debug in region:",i3) ' ) region ; call goPr ! write (gol, '(" itracer :",i3,a8)' ) i_tracer, ' :'//names(i_tracer) ; call goPr ! write (gol, '(" xfrac :",f6.1)' ) efrac ; call goPr ! write (gol, '(" mol_mass :",f6.1 ) ' ) mol_mass ; call goPr ! write (gol, '(" mol_mass_e :",f6.1 )' ) mol_mass_e ; call goPr ! write (gol, '(" emisfield :",2e12.4 )') minval(emis_field), maxval(emis_field) ; call goPr ! write (gol, '(" ubound(3) :",i3 )' ) ubound_e(3) ; call goPr ! sume = sum(emis_field)*month2dt*(mol_mass/mol_mass_e)*efrac ! write (gol, '(" sume :", e12.4 )') sume ; call goPr ! mbef = sum(rm(i1:i2,j1:j2,:,i_tracer)) ! endif do l=1,ubound_e(3) do j=j1,j2 do i=i1,i2 x = emis_field(i,j,l)*month2dt*(mol_mass/mol_mass_e)*efrac rm(i,j,l,i_tracer) = rm(i,j,l,i_tracer) + x #ifdef with_budgets ! budget___ nzone=budg_dat(region)%nzong(i,j) !global budget nzone_v=nzon_vg(l) #ifndef without_emission budemi_dat(region)%emi(i,j,nzone_v,i_tracer) = & budemi_dat(region)%emi(i,j,nzone_v,i_tracer) + x/mol_mass*1e3 if(i_tracer ==1) sum_emission(region) = sum_emission(region) + x !in kg #endif #endif enddo enddo enddo !levels ! if (okdebug) then ! maft = sum(rm(i1:i2,j1:j2,:,i_tracer)) ! write (gol, '(" after-before :", e12.4 )' ) maft-mbef ; call goPr ! write (gol, '(" END debug output ")' ) ; call goPr ! endif nullify(rm) END SUBROUTINE DO_ADD_3D !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: DO_ADD_3D_CYCLE ! ! !DESCRIPTION: Routine to add a 3D emission field (given in kg/month), ! scaled with a daily cycle. To be used for biomass burning ! emissions. !\\ !\\ ! !INTERFACE: ! SUBROUTINE DO_ADD_3D_CYCLE( region, i_tracer, is, js, emis_field, curr_cycle, mol_mass, mol_mass_e, status, xfrac) ! ! !USES: ! use dims, only : CheckShape, tref, nsrce, itaur use dims, only : sec_month, im,jm,lm, okdebug use dims, only : ndyn_max,dx, xref, xbeg,dy,yref,ybeg use global_data, only : mass_dat, region_dat use chem_param, only : ntracet, names use tm5_distgrid, only : get_distgrid, dgrid #ifdef with_budgets use budget_global, only : budg_dat, nzon_vg #endif use datetime, only : tau2date use toolbox, only : dumpfield ! ! !INPUT PARAMETERS: ! integer, intent(in) :: region integer, intent(in) :: i_tracer integer, intent(in) :: is, js real, intent(in) :: emis_field(is:,js:,:) real, dimension(:), intent(in) :: curr_cycle real, intent(in) :: mol_mass ! first is the mass of tracer field real, intent(in) :: mol_mass_e ! mass of emission field (e.g. NOx emission ar in kgN/month) real, intent(in), optional :: xfrac ! partial addition ! ! !OUTPUT PARAMETERS: ! integer, intent(out) :: status ! ! !REVISION HISTORY: ! 1 Oct 2010 - Achim Strunk - v0 based on do_add_3d ! 20 Jun 2012 - Ph. Le Sager - adapted for lon-lat MPI domain decomposition ! ! !REMARKS: ! !EOP !------------------------------------------------------------------------ !BOC character(len=*), parameter :: rname = mname//'/do_add_3d_cycle' integer :: i, j, l, ipos, nzone, nzone_v integer,dimension(6) :: idater integer :: sec_in_day, ntim, itim, i1, j1, i2, j2 integer, dimension(3) :: ubound_e real :: x, xlon, xlat, efrac, dtime, month2dt, dtime2 real,dimension(:,:,:,:),pointer :: rm real :: mbef, maft, sume ! --- begin -------------------------------------- call GET_DISTGRID( dgrid(region), I_STRT=i1, I_STOP=i2, J_STRT=j1, J_STOP=j2 ) status = 0 ! check arguments ubound_e = ubound(emis_field) if(ubound_e(1) /= i2) then write(gol,'(A)') 'do_add_3d: im emission not consistent' ; call goPr status = 1 endif if(ubound_e(2) /= j2) then write(gol,'(A)') 'do_add_3d: jm emission not consistent' ; call goPr status = 1 endif if(ubound_e(3) > lm(region)) then write(gol,'(A)') 'do_add_3d: lm emission not consistent' ; call goPr status = 1 endif IF_NOTOK_RETURN(status=1) rm => mass_dat(region)%rm dtime = float(nsrce)/(2*tref(region)) !timestep emissions month2dt = dtime/sec_month ! conversion to emission per timestep dtime2 = float(ndyn_max)/(2*tref(region)) !timestep emissions based on non-CFL interference (CMK 05/2006) 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) !elapsed seconds this day itim = sec_in_day/nint(dtime2)+1 !time interval if(present(xfrac)) then efrac = xfrac else efrac = 1.0 endif ! if (okdebug) then ! write (gol, '("--------------------------------------------------------------")') ; call goPr ! write (gol, '(" DO ADD 3D CYCLE debug in region:",i3) ' ) region ; call goPr ! write (gol, '(" itracer :",i3,a8)' ) i_tracer, ' :'//names(i_tracer) ; call goPr ! write (gol, '(" xfrac :",f6.1)' ) efrac ; call goPr ! write (gol, '(" mol_mass :",f6.1 ) ' ) mol_mass ; call goPr ! write (gol, '(" mol_mass_e :", f6.1 )' ) mol_mass_e ; call goPr ! write (gol, '(" emisfield :", 2e12.4 )') minval(emis_field),maxval(emis_field) ; call goPr ! write (gol, '(" ubound(3) :", i3 )' ) ubound_e(3) ; call goPr ! sume = sum(emis_field)*month2dt*(mol_mass/mol_mass_e)*efrac ! write (gol, '(" sume :", e12.4 )' ) sume ; call goPr ! mbef = sum(rm(i1:i2,j1:j2,:,i_tracer)) ! write (gol,*) 'lbounds:', is,i1,js,j1 ; call goPr ! endif do l=1,ubound_e(3) do j=j1,j2 do i=i1,i2 xlat = ybeg(region) + (j-0.5)*dy/yref(region) if (xlat .gt. -20 .and. xlat .lt. 20) then ! apply emission cycle in tropics only ! itim = 1 and lon = -180 --->position 1 ! itim = ntim ant lon = -180 --->position ntim, etc. ! itim = 1 and lon = 0 ---->position ntim/2 xlon = xbeg(region) + (i-0.5)*dx/xref(region) ipos = 1 + mod(int((xlon+180.)*ntim/360.) + (itim-1),ntim) !position in array depending on time and lon. x = emis_field(i,j,l)*month2dt*(mol_mass/mol_mass_e)*efrac*curr_cycle(ipos) else x = emis_field(i,j,l)*month2dt*(mol_mass/mol_mass_e)*efrac endif rm(i,j,l,i_tracer) = rm(i,j,l,i_tracer) + x #ifdef with_budgets ! budget___ ! if(region==nregions) then !sub budget finest region ! nzone=nzon(i,j) ! nzone_v=nzon_v(l) ! budemi(nzone,nzone_v,i_tracer)=budemi(nzone,nzone_v,i_tracer)+x/mol_mass*1e3 !mole ! endif nzone=budg_dat(region)%nzong(i,j) !global budget nzone_v=nzon_vg(l) #ifndef without_emission budemi_dat(region)%emi(i,j,nzone_v,i_tracer) = & budemi_dat(region)%emi(i,j,nzone_v,i_tracer) + x/mol_mass*1e3 if(i_tracer ==1) sum_emission(region) = sum_emission(region) + x !in kg #endif #endif enddo enddo enddo !levels ! if (okdebug) then ! maft = sum(rm(i1:i2,j1:j2,:,i_tracer)) ! write (gol, '(" after-before :", e12.4 )' ) maft-mbef ; call goPr ! !if ( sume > tiny(sume) ) then ! ! This is wrong... ! ! if ( (maft-mbef)/sume < 0.5) then ! ! write (gol, '(" spurious do_add_3d_cycle.... ")' ) ; call goPr ! ! call dumpfield(1, 'spurious.hdf', im(region), jm(region), ubound_e(3), 1, emis_field, names(i_tracer)) ! ! endif ! !endif ! write (gol, '(" END debug output ")' ) ; call goPr ! endif nullify(rm) end subroutine do_add_3d_cycle !EOC !-------------------------------------------------------------------------- ! TM5 ! !-------------------------------------------------------------------------- !BOP ! ! !IROUTINE: DISTRIBUTE_L44 ! ! !DESCRIPTION: Scales input field (field2d) such that ! 1) NH temperate vegetation fires are distributed over the ! months: july-august-september ! 2) SH temperate vegetation fires are distributed over the ! months: January-February-March !\\ !\\ ! !INTERFACE: ! SUBROUTINE DISTRIBUTE_L44(month, imdim, jmdim, field2d) ! ! !USES: ! use dims, only : okdebug ! ! !INPUT PARAMETERS: ! integer, intent(in) :: month ! month integer, intent(in) :: jmdim, imdim ! dimension of grid ! ! !REVISION HISTORY: ! 1 Oct 2010 - Achim Strunk - protex ! !EOP !------------------------------------------------------------------------ !BOC real,dimension(imdim,jmdim) :: field2d real :: bb_nh,bb_sh ! multiplication factor for yearly temperate fires integer :: j,i bb_nh=0. bb_sh=0. if (month==7.or.month==8.or.month==9) bb_nh=3./12. if (month==1.or.month==2.or.month==3) bb_sh=3./12. do i=1, imdim do j=1, jmdim/2 field2d(i,j)=field2d(i,j)*bb_sh enddo do j=jmdim/2+1,jmdim field2d(i,j)=field2d(i,j)*bb_nh enddo enddo if(okdebug) then write(gol,*) 'l44 is distributed' ; call goPr endif END SUBROUTINE DISTRIBUTE_L44 !EOC END MODULE EMISSION_DATA