!### 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_ff !----------------------------------------------------------------------- ! purpose ! ------- ! perform c13 fossil fuel emissions ! ! EMISSIONS ! interface: call emission_apply_ff ! method: subroutine is called from its general parent routine apply_emission ! ! COVARIANCE ! interface: call covariance_c13_compute_ff ! method: subroutine is called from MakePriorUncertainty in enkf_tools.F90 ! ! Table of contents: ! module emission_c13_ff ! subroutine declare_emission_c13_ff ! end subroutine declare_emission_c13_ff ! subroutine calc_emission_c13_ff ! end subroutine calc_emission_c13_ff ! subroutine emission_apply_c13_ff(region) ! end subroutine emission_apply_c13_ff ! subroutine free_emission_c13_ff ! end subroutine free_emission_c13_ff ! end module emission_c13_ff ! !----------------------------------------------------------------------- use GO, only: gol, goErr, goPr, goLabel use Dims, only: nregions, nlon360,nlat180 use emission_data, only: gridbox_per_month_to_flux use global_types, only: emis_data use chem_param, only: nmembersloc use emission_common, only: enkf_is_forward, enkf_is_inverse, EmisInputDir use emission_common, only: co2_ff, c13_ff, c13_ff_m, c13_ff_n use emission_data, only: do_add_2d, flux_to_gridbox_per_month implicit none ! --- const ------------------------- character(len=*), parameter :: mname = 'emission_c13_ff' private character(len=2000) :: ff_prefix, ff_varname, ff_postfix ! public routines: public :: declare_emission_c13_ff, free_emission_c13_ff, calc_emission_c13_ff, emission_apply_c13_ff ! c13_ff_regxmem is the result of the calculate routine. ! It holds the coarsened flux per ensemble member and region. type(emis_data),allocatable,dimension(:,:),target :: c13_ff_regxmem contains subroutine declare_emission_c13_ff use dims, only: im,jm,nregions, idate use global_data, only: rcfile use GO, only: TrcFile, Init, Done, ReadRc implicit none integer :: status integer :: region, imr, jmr, n type(TrcFile) :: rcF character(len=*), parameter :: rname = mname//'/declare_emission_c13_ff' if(enkf_is_forward()) then allocate(c13_ff_regxmem(nregions, nmembersloc)) do region=1,nregions imr=im(region) jmr=jm(region) do n=1,nmembersloc allocate(c13_ff_regxmem(region,n)%surf(imr,jmr)) end do enddo end if call Init( rcF, rcfile, status ) IF_NOTOK_RETURN(status=1) call ReadRc( rcF, 'ff.file.prefix', ff_prefix ,status) IF_NOTOK_RETURN(status=1) call ReadRc( rcF, 'ff.file.postfix', ff_postfix ,status ) IF_NOTOK_RETURN(status=1) call ReadRc( rcF, 'ff.varname', ff_varname ,status) IF_NOTOK_RETURN(status=1) call Done( rcF, status ) IF_NOTOK_RETURN(status=1) end subroutine declare_emission_c13_ff subroutine calc_emission_c13_ff use dims, only: nlon360, nlat180, im,jm, sec_month, newsrun, idate, newmonth, sec_year, itau, nread use emission_common, only: ReadFromFile,ReadFromFile3D use emission_data, only: chardate use chem_param, only: xmc,xmc13 use toolbox, only: coarsen_emission_2d, escape_tm use go_string, only: goNum2str use global_data, only: rcfile use GO, only: TrcFile, Init, Done, ReadRc implicit none integer :: status integer :: region, imr, jmr integer, parameter :: annual_field = 0, rank2 = 2, level1 = 1 character(len=10) :: cdate real :: year2month real,dimension(nlon360,nlat180) :: c13_ff_mem ! flux for a given ensemble member real,dimension(:,:,:),allocatable :: emis3d integer, parameter :: add_field=0 integer :: mem type(TrcFile) :: rcF real,parameter :: deltaff=-30. character(len=*), parameter :: rname = mname//'/calc_emission_c13_ff' if(mod(itau,nread) /= 0) return ! only every nread hours if(newmonth) then allocate(emis3d(nlon360,nlat180,12)) cdate=trim(goNum2str(idate(1),'(i4.4)')) call ReadFromFile3D(trim(ff_prefix)//trim(cdate)//trim(ff_postfix),trim(ff_varname),emis3d,status,netcdf4=.True.) co2_ff = emis3d(:,:,idate(2))*1e3/xmc ! from kgC/m2/s to mol/m2/s deallocate(emis3d) if (status /= 0) call escape_tm('stop in '//rname) endif !WP! this block has been moved to be outside the (newmonth) flag to ensure that the averaging is done at each call. !WP! This ensures that the n-hourly average fossil fuel emissions written by a forward code are non-zero even though we only read !WP! and coarsen these emissions each month c13_ff=co2_ff*(0.011112*(1.+deltaff/1000.)) if(allocated(c13_ff_m)) then c13_ff_m=c13_ff_m+c13_ff ! mol/m2/s c13_ff_n=c13_ff_n+1 endif if(newmonth) then ! 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. do mem=1,nmembersloc c13_ff_mem=c13_ff ! mol/m2/s (this is where parameters would be applied for optimized fluxes) c13_ff_mem=c13_ff_mem*xmc13/1.e3 ! from mol/m2/s to kgC/m2/s call flux_to_gridbox_per_month(c13_ff_mem) ! to kgC/area/month call coarsen_emission_2d('c13_ff',360,180,c13_ff_mem,c13_ff_regxmem(:,mem),add_field) enddo endif end subroutine calc_emission_c13_ff subroutine emission_apply_c13_ff(region) ! This routine should only be called by the forward enkf code. It ! applies the ff c13 flux by member and region to the ic13ff ! 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_ff, xmc13, xm13co2 implicit none integer, intent(in) :: region integer :: mem character(len=*), parameter :: rname = mname//'/emission_apply_c13_ff' if(enkf_is_forward()) then if(tracer_active(ic13_ff)) then call do_add_2d(region,ic13_ff,1,c13_ff_regxmem(region,1)%surf,xm13co2,xmc13) endif endif if(okdebug) write(*,*) 'end of emission_apply_c13_ff' end subroutine emission_apply_c13_ff subroutine free_emission_c13_ff use dims, only: nregions implicit none integer :: region,n character(len=*), parameter :: rname = mname//'/free_emission_c13_ff' if(enkf_is_forward()) then do region=1,nregions do n=1,nmembersloc deallocate(c13_ff_regxmem(region,n)%surf) end do enddo deallocate(c13_ff_regxmem) end if end subroutine free_emission_c13_ff end module emission_c13_ff