!### 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_bio !----------------------------------------------------------------------- ! purpose ! ------- ! perform c13 biosphere emissions ! ! EMISSIONS ! interface: call emission_apply_bio ! method: subroutine is called from its general parent routine apply_emission ! ! ! Table of contents: ! module emission_c13_bio ! subroutine declare_emission_c13_bio ! end subroutine declare_emission_c13_bio ! subroutine calc_emission_c13_bio ! end subroutine calc_emission_c13_bio ! subroutine emission_apply_c13_bio(region) ! end subroutine emission_apply_c13_bio ! subroutine free_emission_c13_bio ! end subroutine free_emission_c13_bio ! end module emission_c13_bio ! !----------------------------------------------------------------------- use GO, only: gol, goErr, goPr, goLabel use Dims, only: nregions, nlon360,nlat180 use Meteo, only: t2m_dat,ssr_dat use Meteo, only: Set use global_types, only: emis_data use emission_data, only: do_add_2d, flux_to_gridbox_per_month use emission_common, only: enkf_is_forward,enkf_is_inverse, EmisInputDir use emission_common, only: field2d_eps,field2d_eps_m,field2d_eps_n,co2_bio,c13_bio, c13_bio_m, c13_bio_n, flux_means_are_priors, fluxmultiplier,fluxmultiplier_epsbio,c13_sig use emission_common, only: c13_gpp, c13_gpp_m, c13_gpp_n use emission_common, only: c13_res, c13_res_m, c13_res_n use emission_common, only: field2d_eps,field2d_eps_m,field2d_eps_n,co2_gpp, co2_gpp_m, co2_gpp_n use chem_param, only: nmembersloc implicit none ! --- const ------------------------- character(len=*), parameter :: mname = 'emission_c13_bio' character(len=500) :: eps_prefix,bio_prefix, bio_prefix2, bio_postfix, bio_g_varname , bio_r_varname, bio_p_varname, bio_n_varname,bio_c13_n_varname,bio_w_varname private ! public routines: public :: declare_emission_c13_bio, free_emission_c13_bio, calc_emission_c13_bio, emission_apply_c13_bio logical :: use_fires ! c13_bio_regxmem is the result of the calculate routine. ! It holds the coarsened flux per ensemble member and region. type(emis_data),dimension(:,:),allocatable,target :: c13_bio_regxmem contains subroutine declare_emission_c13_bio use GO, only : TrcFile, Init, Done, ReadRc use global_data, only : rcfile use dims, only: im,jm,nregions, iglbsfc implicit none integer :: region, imr, jmr, n integer :: status type(TrcFile) :: rcF character(len=*), parameter :: rname = mname//'/declare_emission_c13_bio' call Init( rcF, rcfile, status ) call ReadRc( rcF, 'emis.use.fires', use_fires ,status ,default=.false.) bio_prefix = trim(EmisInputDir)//'/biofireparams_' call ReadRc( rcF, 'bio.file.prefix', bio_prefix ,status ,default=bio_prefix) call ReadRc( rcF, 'bio.file.prefix2', bio_prefix2 ,status ,default=bio_prefix2) call ReadRc( rcF, 'bio.file.postfix', bio_postfix ,status ,default='.hdf') call ReadRc( rcF, 'bio.r.varname', bio_r_varname ,status ,default='resp13') call ReadRc( rcF, 'bio.p.varname', bio_p_varname ,status ,default='gpp13') call ReadRc( rcF, 'bio.g.varname', bio_g_varname ,status ,default='gpp') call ReadRc( rcF, 'bio.c13.n.varname', bio_c13_n_varname ,status ,default='nep13') call ReadRc( rcF, 'bio.n.varname', bio_n_varname ,status ,default='nep') call ReadRc( rcF, 'bio.w.varname', bio_w_varname ,status ,default='wk') call ReadRc( rcF, 'eps.prefix', eps_prefix ,status) call Done( rcF, status ) call Set( t2m_dat(iglbsfc), status, used=.true. ) call Set( ssr_dat(iglbsfc), status, used=.true. ) if(enkf_is_forward()) then allocate(c13_bio_regxmem(nregions,nmembersloc)) do region=1,nregions imr=im(region) jmr=jm(region) do n=1,nmembersloc allocate(c13_bio_regxmem(region,n)%surf(imr,jmr)) end do enddo end if end subroutine declare_emission_c13_bio subroutine calc_emission_c13_bio use dims, only : nlon360, nlat180, im,jm, sec_month, newsrun, idate, itau, staggered, nread, newmonth, iglbsfc use emission_common, only : ReadFromFile3D, ReadFromFile, get_ntimesteps use toolbox, only : escape_tm use toolbox, only : coarsen_emission_2d use emission_data, only : chardate use Meteo, only : Set use chem_param, only : xmc,xmc13 use go_string, only : goNum2str use rs_read, only : da use partools, only : myid implicit none integer :: status,i,j,k integer :: region, imr, jmr, mem, ntimesteps integer, parameter :: annual_field = 0, rank2 = 2, level1 = 1 real,dimension(nlon360,nlat180,248),save :: p0,r0,g0,c13_signature real,dimension(nlon360,nlat180,248),save :: n0,n13c0,w0 real,dimension(nlon360,nlat180),save :: Delta,Delta2,alpha !,field2d_epsbio real,dimension(nlon360,nlat180) :: c13_bio_mem, c13_gpp_mem, c13_res_mem,co2_gpp_mem ! flux for a given ensemble member real,dimension(4,nlon360,nlat180),save :: field2d_epsbio integer, parameter :: add_field=0 real,pointer,dimension(:,:) :: t2m,ssr character(len=10) :: cdate integer,save :: counter,counter2 character(len=*), parameter :: rname = mname//'/calc_emission_c13_bio' if(mod(itau,nread) /= 0) return ! only every nread hours if(newmonth) then call get_ntimesteps(ntimesteps) cdate=trim(goNum2str(idate(1),'(i4.4)'))//trim(goNum2str(idate(2),'(i2.2)')) if(use_fires) then write (gol,'("[3-hourly biofireparams] filename : ",a,".")') (trim(bio_prefix)//trim(cdate)//trim(bio_postfix)) ; call goPr call ReadFromFile3D(trim(bio_prefix)//trim(cdate)//trim(bio_postfix),trim(bio_r_varname),r0(:,:,:ntimesteps),status,netcdf4=.True.) call ReadFromFile3D(trim(bio_prefix)//trim(cdate)//trim(bio_postfix),trim(bio_p_varname),p0(:,:,:ntimesteps),status,netcdf4=.True.) call ReadFromFile3D(trim(bio_prefix)//trim(cdate)//trim(bio_postfix),trim(bio_n_varname),n0(:,:,:ntimesteps),status,netcdf4=.True.) call ReadFromFile3D(trim(bio_prefix)//trim(cdate)//trim(bio_postfix),trim(bio_g_varname),g0(:,:,:ntimesteps),status,netcdf4=.True.) call ReadFromFile3D(trim(bio_prefix)//trim(cdate)//trim(bio_postfix),trim(bio_w_varname),w0(:,:,:ntimesteps),status,netcdf4=.True.) call ReadFromFile3D(trim(bio_prefix)//trim(cdate)//trim(bio_postfix),trim(bio_c13_n_varname),n13c0(:,:,:ntimesteps),status,netcdf4=.True.) call ReadFromFile3D(trim(eps_prefix)//trim(cdate)//'.nc','epsilon',field2d_epsbio,status,netcdf4=.True.) endif if (status /= 0) call escape_tm('stop in '//rname) endif counter = ((idate(3)-1)*8 + int(idate(4)/3) + 1) if (idate(2) .eq. 2 .and. idate(3) .eq. 29) then counter = counter - 8 endif counter2 = (int((idate(3)-1)/7) + 1) if (idate(3) .gt. 28) then counter2 = 4 endif co2_gpp=g0(:,:,counter) co2_bio=n0(:,:,counter) field2d_eps = field2d_epsbio(counter2,:,:) alpha = 1+field2d_eps/1000 c13_bio = co2_bio*(0.011112*(1.+da/1000.))*alpha ! accumulate weekly average if(allocated(c13_bio_m)) then if(flux_means_are_priors .eqv. .TRUE.) then c13_bio_m=c13_bio_m+c13_bio ! mol/m2/s c13_gpp_m=c13_gpp_m+co2_gpp*field2d_eps ! mol/m2/s co2_gpp_m=co2_gpp_m+co2_gpp ! mol/m2/s field2d_eps_m=field2d_eps_m+field2d_eps else c13_bio_m=c13_bio_m+c13_bio*fluxmultiplier(:,:,1) ! mol/m2/s, tracer 0 on PE0 is mean c13_gpp_m=c13_gpp_m+co2_gpp*fluxmultiplier_epsbio(:,:,1) ! mol/m2/s, tracer 0 on PE0 is mean co2_gpp_m=co2_gpp_m+co2_gpp !fluxmultiplier(:,:,1) ! mol/m2/s, tracer 0 on PE0 is mean field2d_eps_m=field2d_eps_m+field2d_eps*fluxmultiplier_epsbio(:,:,1) endif c13_bio_n=c13_bio_n+1 c13_gpp_n=c13_gpp_n+1 co2_gpp_n=co2_gpp_n+1 field2d_eps_n=field2d_eps_n+1 endif ! The rest only useful to the forward code if(enkf_is_inverse()) return ! Apply net flux scaling by ensemble member, change units, and ! coarsen to the regions. When parameters are used beyond simple ! scaling of the net flux, some piece of logic from above will ! have to be subsumed by this loop to get the parameter ! deviations. do mem=1,nmembersloc c13_bio_mem=c13_bio*fluxmultiplier(:,:,mem) ! parameter applied, mol/m2/s c13_bio_mem=c13_bio_mem*xmc13/1.e3 ! from mol/m2/s to kgC/m2/s call flux_to_gridbox_per_month(c13_bio_mem) ! to kgC13/area/month call coarsen_emission_2d('c13_bio',360,180,c13_bio_mem,c13_bio_regxmem(:,mem),add_field) enddo end subroutine calc_emission_c13_bio subroutine emission_apply_c13_bio(region) ! This routine should only be called by the forward enkf code. It ! applies the bio c13 flux by member and region to the ic13_bio ! tracer. In the inverse code, the calculated flux from each ! module is applied to the (single) c13 tracer. use dims, only: okdebug use ParTools, only: tracer_active use chem_param, only : ic13_bio, ic13_gpp, ic13_res, ico2_gpp, xmc13, xm13co2,xmc,xmco2 implicit none integer, intent(in) :: region integer :: mem character(len=*), parameter :: rname = mname//'/emission_apply_c13_bio' if(enkf_is_forward()) then if(tracer_active(ic13_bio)) then call do_add_2d(region,ic13_bio,1,c13_bio_regxmem(region,1)%surf,xm13co2,xmc13) endif endif if(okdebug) write(*,*) 'end of emission_apply_c13_bio' end subroutine emission_apply_c13_bio subroutine free_emission_c13_bio use dims, only: nregions implicit none integer :: region,n character(len=*), parameter :: rname = mname//'/free_emission_c13_bio' if(enkf_is_forward()) then do region=1,nregions do n=1,nmembersloc deallocate(c13_bio_regxmem(region,n)%surf) end do enddo deallocate(c13_bio_regxmem) endif end subroutine free_emission_c13_bio end module emission_c13_bio