!### macros ##################################################### #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" !################################################################# module emission_c13_all !----------------------------------------------------------------------- ! purpose ! ------- ! perform c13 emissions all ! ! interface_ ! --------- ! call emission_apply_all ! ! method ! ------ ! subroutine is called from its general parent routine apply_emission !----------------------------------------------------------------------- use chem_param , only: xm13co2, xmc13 use global_types, only: emis_data use Dims, only: nregions, nlon360,nlat180 use emission_data, only: do_add_2d, flux_to_gridbox_per_month use emission_common, only: field2d_eps,c13_bio, c13_gpp,c13_gpp_m,c13_gpp_n,co2_gpp,co2_gpp_m,co2_gpp_n, c13_res, co2_bio,co2_fires,c13_sig,c13_ocean, c13_ff, c13_fires, c13_biodis, c13_oceandis, c13_all, c13_all_m, c13_all_n,flux_means_are_priors,fluxmultiplier,fluxmultiplier_biodis,fluxmultiplier_epsbio implicit none private ! public routines: public :: emission_apply_c13_all, declare_emission_c13_all, free_emission_c13_all, calc_emission_c13_all public :: c13_flux, emis_data, pseudo type(emis_data),dimension(:,:),allocatable,target :: c13_flux logical :: flux_unallocated=.true. logical :: pseudo=.false. contains subroutine declare_emission_c13_all use dims, only: im,jm,nregions use ParTools, only: ntracetloc use chem_param, only: nmembersloc implicit none integer :: region, imr, jmr, n allocate(c13_flux(nregions,ntracetloc)) do region=1,nregions imr=im(region) jmr=jm(region) ! do n=1,nmembersloc !,ntracetloc do n=1,ntracetloc allocate(c13_flux(region,n)%surf(imr,jmr)) end do enddo end subroutine declare_emission_c13_all subroutine calc_emission_c13_all use GO, only: gol, goErr, goPr, goLabel use dims, only: okdebug,nlon360, nlat180, im,jm, sec_month, newsrun, idate, newday,nregions, itau, staggered, nread use toolbox, only: coarsen_emission_2d use ParTools, only: tracer_active, ntracetloc use chem_param, only: xmc13,nmembersloc use rs_read, only: da use partools, only : myid,root_t implicit none integer :: st,region, imr, jmr, n, mem,tr integer, parameter :: add_field=0 real,dimension(nlon360,nlat180),save :: alpha if(mod(itau,nread) /= 0) return ! only every nread hours !WP! put full flux into emission array for each member do mem=1,ntracetloc alpha = 1+(field2d_eps*fluxmultiplier_epsbio(:,:,mem))/1000. c13_bio = co2_bio*(0.011112*(1. + da/1000.))*alpha c13_all=((c13_bio+c13_ocean)*fluxmultiplier(:,:,mem)+c13_ff+c13_fires+c13_biodis+c13_oceandis) c13_all=c13_all*xmc13/1.e3 ! full emissions in kg 13C/m2/s call flux_to_gridbox_per_month(c13_all) ! to kgC/area/month call coarsen_emission_2d('c13',360,180,c13_all,c13_flux(:,mem),add_field) end do ! accumulate weekly average if(flux_means_are_priors .eqv. .TRUE.) then if(allocated(c13_all_m)) then c13_all_m=c13_all_m+(c13_bio+c13_ocean+c13_ff+c13_fires+c13_biodis+c13_oceandis) c13_all_n=c13_all_n+1 end if else if(allocated(c13_all_m)) then alpha = 1+(field2d_eps*fluxmultiplier_epsbio(:,:,1))/1000 c13_bio = co2_bio*(0.011112*(1. + da/1000.))*alpha c13_all_m=c13_all_m+((c13_bio+c13_ocean)*fluxmultiplier(:,:,1)+c13_ff+c13_fires+c13_biodis+c13_oceandis) c13_all_n=c13_all_n+1 end if end if end subroutine calc_emission_c13_all subroutine emission_apply_c13_all(region) ! use GO, only: goLoCase use dims, only: okdebug use ParTools, only: tracer_active, ntracetloc, myid, ntracet_ar use chem_param, only: nmembersloc,ic13_bio,names implicit none ! io integer, intent(in) :: region ! local integer :: n, offsetn=0 ! start if(myid > 0) offsetn=sum(ntracet_ar(0:myid-1)) if (goLoCase(trim(names(offsetn+1)))=='co2c13'.or.goLoCase(trim(names(offsetn+1)))=='co2c13_bg') then do n=1,ntracetloc call do_add_2d(region,offsetn+n,1,c13_flux(region,n)%surf,xm13co2,xmc13) enddo endif if(okdebug) write(*,*) 'end of emission_apply_c13_casm' end subroutine emission_apply_c13_all !================================================================== subroutine free_emission_c13_all use dims, only: nregions_max use ParTools, only: ntracetloc implicit none integer :: region, n do region=1,nregions do n=1,ntracetloc deallocate(c13_flux(region,n)%surf) end do enddo end subroutine free_emission_c13_all end module emission_c13_all