!### 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_fires !----------------------------------------------------------------------- ! purpose ! ------- ! perform c13 fire emissions ! ! interface ! --------- ! call emission_apply_fires ! ! method ! ------ ! subroutine is called from its general parent routine apply_emission !----------------------------------------------------------------------- use GO, only: gol, goErr, goPr, goLabel use Dims, only: nregions, nlon360,nlat180 use emission_data, only: gridbox_per_month_to_flux use chem_param , only: xmco2, xmc,xm13co2,xmc13 use global_types, only: emis_data use emission_common, only: enkf_is_forward, enkf_is_inverse, EmisInputDir use emission_common, only: c13_fires, c13_fires_m, c13_fires_n, flux_means_are_priors use emission_data, only: do_add_2d, flux_to_gridbox_per_month use chem_param, only: nmembersloc implicit none ! --- const ------------------------- character(len=*), parameter :: mname = 'emission_c13_fire' private ! public routines: public :: declare_emission_c13_fires, free_emission_c13_fires, calc_emission_c13_fires, emission_apply_c13_fires character(len=500) :: eps_prefix,fire_prefix, fire_postfix, fire_varname logical :: use_fires ! c13_fires_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_fires_regxmem contains subroutine declare_emission_c13_fires use GO, only : TrcFile, Init, Done, ReadRc use global_data, only : rcfile use dims, only: im,jm,nregions implicit none type(TrcFile) :: rcF integer :: status integer :: region, imr, jmr, n call Init( rcF, rcfile, status ) call ReadRc( rcF, 'emis.use.fires', use_fires ,status ,default=.false.) call ReadRc( rcF, 'fire.file.prefix', fire_prefix ,status,default=trim(EmisInputDir)//'/BBguido1x1' ) call ReadRc( rcF, 'fire.file.postfix', fire_postfix ,status ,default='.hdf') call ReadRc( rcF, 'fire.varname', fire_varname ,status ,default='bb') call ReadRc( rcF, 'eps.prefix', eps_prefix ,status) call Done( rcF, status ) if(enkf_is_forward()) then allocate(c13_fires_regxmem(nregions,nmembersloc)) do region=1,nregions imr=im(region) jmr=jm(region) do n=1,nmembersloc allocate(c13_fires_regxmem(region,n)%surf(imr,jmr)) end do enddo end if end subroutine declare_emission_c13_fires subroutine calc_emission_c13_fires use dims, only: nlon360, nlat180, im,jm, sec_month, newsrun, idate, newmonth, sec_year, itau, nread use emission_common, only: ReadFromFile3D, get_ntimesteps use emission_data, only: chardate use toolbox, only: coarsen_emission_2d, escape_tm use go_string, only: goNum2str use rs_read, only: da implicit none integer :: status,i,j,k integer :: region, imr, jmr, mem, ntimesteps integer, parameter :: annual_field = 0, rank2 = 2, level1 = 1 character(len=10) :: cdate real :: year2month real,dimension(nlon360,nlat180,248),save :: f0 real,dimension(4,nlon360,nlat180),save :: field2d_epsbio real,dimension(nlon360,nlat180) :: c13_fires_mem ! flux for a given ensemble member real,dimension(nlon360,nlat180),save :: alpha integer, parameter :: add_field=0 integer,save :: counter,counter2 character(len=*), parameter :: rname = mname//'/calc_emission_c13_fires' 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)')) call ReadFromFile3D(trim(fire_prefix)//trim(cdate)//trim(fire_postfix),trim(fire_varname),f0(:,:,:ntimesteps),status,netcdf4=.True.) call ReadFromFile3D(trim(eps_prefix)//trim(cdate)//'.nc','epsilon',field2d_epsbio,status,netcdf4=.True.) 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 alpha = 1+field2d_epsbio(counter2,:,:)/1000 c13_fires = f0(:,:,counter)*(0.011112*(1.+da/1000.))*alpha !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 fire emissions written by a forward code are non-zero even though we only read !WP! and coarsen these emissions each month ! accumulate weekly average if(allocated(c13_fires_m)) then c13_fires_m=c13_fires_m+c13_fires ! mol/m2/s c13_fires_n=c13_fires_n+1 endif if(newmonth.or.newsrun) then ! The rest only useful to the forward code if(enkf_is_inverse()) return ! Change units, and coarsen to the regions. do mem=1,nmembersloc c13_fires_mem=c13_fires*xmc13/1.e3 ! from mol/m2/s to kgC/m2/s call flux_to_gridbox_per_month(c13_fires_mem) ! to kgC/area/month call coarsen_emission_2d('c13_fires',360,180,c13_fires_mem,c13_fires_regxmem(:,mem),add_field) enddo endif ! if newmonth .or. newrun end subroutine calc_emission_c13_fires subroutine emission_apply_c13_fires(region) ! This routine should only be called by the forward enkf code. It ! applies the fires c13 flux by member and region to the ic13_fires ! 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_fires implicit none integer, intent(in) :: region integer :: mem if(enkf_is_forward()) then if(tracer_active(ic13_fires)) then call do_add_2d(region,ic13_fires,1,c13_fires_regxmem(region,1)%surf,xm13co2,xmc13) endif endif if(okdebug) write(*,*) 'end of emission_apply_c13_fires' end subroutine emission_apply_c13_fires subroutine free_emission_c13_fires use dims, only: nregions implicit none integer :: region,n if (enkf_is_forward()) then do region=1,nregions do n=1,nmembersloc deallocate(c13_fires_regxmem(region,n)%surf) end do enddo deallocate(c13_fires_regxmem) end if end subroutine free_emission_c13_fires end module emission_c13_fires