!### 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_co2_ff !----------------------------------------------------------------------- ! purpose ! ------- ! perform co2 fossil fuel emissions ! ! EMISSIONS ! interface: call emission_apply_ff ! method: subroutine is called from its general parent routine apply_emission ! ! COVARIANCE ! interface: call covariance_co2_compute_ff ! method: subroutine is called from MakePriorUncertainty in enkf_tools.F90 ! ! Table of contents: ! module emission_co2_ff ! subroutine declare_emission_co2_ff ! end subroutine declare_emission_co2_ff ! subroutine calc_emission_co2_ff ! end subroutine calc_emission_co2_ff ! subroutine emission_apply_co2_ff(region) ! end subroutine emission_apply_co2_ff ! subroutine free_emission_co2_ff ! end subroutine free_emission_co2_ff ! end module emission_co2_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, co2_ff_m, co2_ff_n use emission_data, only: do_add_2d, flux_to_gridbox_per_month implicit none ! --- const ------------------------- character(len=*), parameter :: mname = 'emission_co2_ff' private character(len=2000) :: ff_prefix, ff_varname, ff_postfix ! public routines: public :: declare_emission_co2_ff, free_emission_co2_ff, calc_emission_co2_ff, emission_apply_co2_ff type(emis_data),allocatable,dimension(:,:),target :: co2_ff_regxmem contains subroutine declare_emission_co2_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_co2_ff' if(enkf_is_forward()) then allocate(co2_ff_regxmem(nregions, nmembersloc)) do region=1,nregions imr=im(region) jmr=jm(region) do n=1,nmembersloc allocate(co2_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_co2_ff subroutine calc_emission_co2_ff use dims, only: nlon360, nlat180, im,jm, sec_month, newsrun, idate, newmonth, sec_year, itau, nread use emission_common, only: ReadFromFile3D use emission_data, only: chardate use chem_param, only: xmc 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) :: co2_ff_mem ! flux for a given ensemble member real,dimension(:,:,:),allocatable :: emis3d integer, parameter :: add_field=0 integer :: mem integer :: hid type(TrcFile) :: rcF character(len=*), parameter :: rname = mname//'/calc_emission_co2_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 ! accumulate weekly average ! ! Since FF fluxes are not optimized, the following ! if statement is extraneous; we want to use co2_ff ! regardless of the setting of flux_means_are_priors. ! -arj, 7 Feb 07. ! ! if(flux_means_are_priors .eqv. .TRUE.) then !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 if(allocated(co2_ff_m)) then co2_ff_m=co2_ff_m+co2_ff ! mol/m2/s co2_ff_n=co2_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 co2_ff_mem=co2_ff ! mol/m2/s (this is where parameters would be applied for optimized fluxes) co2_ff_mem=co2_ff_mem*xmc/1.e3 ! from mol/m2/s to kgC/m2/s call flux_to_gridbox_per_month(co2_ff_mem) ! to kgC/area/month call coarsen_emission_2d('co2_ff',360,180,co2_ff_mem,co2_ff_regxmem(:,mem),add_field) enddo endif end subroutine calc_emission_co2_ff subroutine emission_apply_co2_ff(region) ! This routine should only be called by the forward enkf code. It ! applies the ff co2 flux by member and region to the ico2ff ! tracer. In the inverse code, the calculated flux from each ! module is applied to the (single) co2 tracer. use dims, only: okdebug use ParTools, only: tracer_active use chem_param, only : ico2_ff, xmc, xmco2 implicit none integer, intent(in) :: region integer :: mem character(len=*), parameter :: rname = mname//'/emission_apply_co2_ff' if(enkf_is_forward()) then if(tracer_active(ico2_ff)) then call do_add_2d(region,ico2_ff,1,co2_ff_regxmem(region,1)%surf,xmco2,xmc) endif endif if(okdebug) write(*,*) 'end of emission_apply_co2_ff' end subroutine emission_apply_co2_ff subroutine free_emission_co2_ff use dims, only: nregions implicit none integer :: region,n character(len=*), parameter :: rname = mname//'/free_emission_co2_ff' if(enkf_is_forward()) then do region=1,nregions do n=1,nmembersloc deallocate(co2_ff_regxmem(region,n)%surf) end do enddo deallocate(co2_ff_regxmem) end if end subroutine free_emission_co2_ff end module emission_co2_ff