#define IF_ERROR_RETURN(action) if (status> 0) then; TRACEBACK; action; return; end if #include "tm5.inc" !################################################################# module emission_c13_biodis !----------------------------------------------------------------------- ! purpose ! ------- ! perform c13 biosphere disequilibrium emissions ! ! EMISSIONS ! interface: call emission_apply_bio ! method: subroutine is called from its general parent routine apply_emission ! ! ! Table of contents: ! module emission_c13_biodis ! subroutine declare_emission_c13_biodis ! end subroutine declare_emission_c13_biodis ! subroutine calc_emission_c13_biodis ! end subroutine calc_emission_c13_biodis ! subroutine emission_apply_c13_biodis(region) ! end subroutine emission_apply_c13_biodis ! subroutine free_emission_c13_biodis ! end subroutine free_emission_c13_biodis ! end module emission_c13_biodis ! !----------------------------------------------------------------------- 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,gridbox_per_month_to_flux use emission_common, only: enkf_is_forward,enkf_is_inverse, EmisInputDir,C13InputDir use emission_common, only: c13_biodis, c13_biodis_m, c13_biodis_n, flux_means_are_priors, fluxmultiplier use chem_param, only: nmembersloc implicit none ! --- const ------------------------- character(len=*), parameter :: mname = 'emission_c13_biodis' private ! public routines: public :: declare_emission_c13_biodis, free_emission_c13_biodis, calc_emission_c13_biodis, emission_apply_c13_biodis logical :: use_fires character(len=500) :: terdis_prefix ! c13_biodis_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_biodis_regxmem contains subroutine declare_emission_c13_biodis 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_biodis' call Init( rcF, rcfile, status ) call ReadRc( rcF, 'emis.use.fires', use_fires ,status ,default=.false.) call Done( rcF, status ) call Set( t2m_dat(iglbsfc), status, used=.true. ) call Set( ssr_dat(iglbsfc), status, used=.true. ) call Init( rcF, rcfile, status ) call ReadRc( rcF, 'terdis.prefix', terdis_prefix ,status) call Done( rcF ,status) if(enkf_is_forward()) then allocate(c13_biodis_regxmem(nregions,nmembersloc)) do region=1,nregions imr=im(region) jmr=jm(region) do n=1,nmembersloc allocate(c13_biodis_regxmem(region,n)%surf(imr,jmr)) end do enddo end if end subroutine declare_emission_c13_biodis subroutine calc_emission_c13_biodis use dims, only : nlon360, nlat180, im,jm, sec_month, newsrun, idate, itau, staggered, nread, newmonth, iglbsfc use emission_common, only : ReadFromFile3D,ReadFromFile 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,xm13co2 use go_string, only : goNum2str implicit none integer :: status integer :: region, imr, jmr, mem integer, parameter :: annual_field = 0, rank2 = 2, level1 = 1 real,dimension(nlon360,nlat180) :: p0,r0 real,dimension(nlon360,nlat180,12),save :: field3d_dis real,dimension(nlon360,nlat180),save :: field2d_dis real,dimension(nlon360,nlat180) :: c13_biodis_mem ! flux for a given ensemble member integer, parameter :: add_field=0 real,parameter :: delta_atm=-8. real,pointer,dimension(:,:) :: t2m,ssr character(len=10) :: cdate character(len=*), parameter :: rname = mname//'/calc_emission_c13_biodis' !WP! CAUTION, this statement will prevent accumulation of flux1x1 means on timescale that are not integer multiples of !WP! nread. Perhaps in the future change this to hourly, or every time step for instance when interoplation of meteorology !WP! is used???! if(mod(itau,nread) /= 0) return ! only every nread hours if(newmonth) then cdate=trim(goNum2str(idate(1),'(i4.4)'))//trim(goNum2str(idate(2),'(i2.2)')) !Read biosphere disequilibrium write(gol,'(a)') 'BIODIS data to be read '//trim(cdate) ; call goPr call ReadFromFile(trim(terdis_prefix)//trim(cdate)//'.nc','ter_dis',field2d_dis,status,netcdf4=.True.) if (status /= 0) call escape_tm('stop in '//rname) endif write(gol,'(a)') 'CarbonTracker BIODIS data to be recalculated ' ; call goPr c13_biodis=field2d_dis(:,:) !kg 13CO2/cell/month call gridbox_per_month_to_flux(c13_biodis) ! to kg 13CO2 /m2/s c13_biodis=c13_biodis*1.e3/xm13co2 ! to mol/m2/s ! accumulate weekly average if(allocated(c13_biodis_m)) then if(flux_means_are_priors .eqv. .TRUE.) then c13_biodis_m=c13_biodis_m+c13_biodis ! mol/m2/s else c13_biodis_m=c13_biodis_m+c13_biodis endif c13_biodis_n=c13_biodis_n+1 endif 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_biodis_mem=c13_biodis !*fluxmultiplier(:,:,mem) ! parameter applied, mol/m2/s c13_biodis_mem=c13_biodis_mem*xmc13/1.e3 ! from mol/m2/s to kgC13/m2/s call flux_to_gridbox_per_month(c13_biodis_mem) ! to kgC13/area/month call coarsen_emission_2d('c13_biodis',360,180,c13_biodis_mem,c13_biodis_regxmem(:,mem),add_field) enddo end subroutine calc_emission_c13_biodis subroutine emission_apply_c13_biodis(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: myid,tracer_active use chem_param, only : ic13_biodis, xmc13, xm13co2 implicit none integer, intent(in) :: region integer :: mem character(len=*), parameter :: rname = mname//'/emission_apply_c13_biodis' if(enkf_is_forward()) then if(tracer_active(ic13_biodis)) then call do_add_2d(region,ic13_biodis,1,c13_biodis_regxmem(region,1)%surf,xm13co2,xmc13) endif endif if(okdebug) write(*,*) 'end of emission_apply_c13_biodis' end subroutine emission_apply_c13_biodis subroutine free_emission_c13_biodis use dims, only: nregions implicit none integer :: region,n character(len=*), parameter :: rname = mname//'/free_emission_c13_biodis' if(enkf_is_forward()) then do region=1,nregions do n=1,nmembersloc deallocate(c13_biodis_regxmem(region,n)%surf) end do enddo deallocate(c13_biodis_regxmem) endif end subroutine free_emission_c13_biodis end module emission_c13_biodis