!### 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_oceandis !----------------------------------------------------------------------- ! purpose ! ------- ! perform c13 ocean disequilibrium emissions ! ! EMISSIONS ! interface: call emission_apply_oceandis ! method: subroutine is called from its general parent routine apply_emission ! ! ! Table of contents: ! module emission_c13_oceandis ! subroutine declare_emission_c13_oceandis ! end subroutine declare_emission_c13_oceandis ! subroutine calc_emission_c13_oceandis ! end subroutine calc_emission_c13_oceandis ! subroutine emission_apply_c13_oceandis(region) ! end subroutine emission_apply_c13_oceandis ! subroutine free_emission_c13_oceandis ! end subroutine free_emission_c13_oceandis ! end module emission_c13_oceandis ! !----------------------------------------------------------------------- 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_oceandis, c13_oceandis_m, c13_oceandis_n, flux_means_are_priors, fluxmultiplier use chem_param, only: nmembersloc implicit none ! --- const ------------------------- character(len=*), parameter :: mname = 'emission_c13_oceandis' private ! public routines: public :: declare_emission_c13_oceandis, free_emission_c13_oceandis, calc_emission_c13_oceandis, emission_apply_c13_oceandis logical :: use_fires character(len=500) :: ocedis_prefix type(emis_data),dimension(:,:),allocatable,target :: c13_oceandis_regxmem contains subroutine declare_emission_c13_oceandis 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_oceandis' 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, 'ocedis.prefix', ocedis_prefix ,status) call Done( rcF ,status) if(enkf_is_forward()) then allocate(c13_oceandis_regxmem(nregions,nmembersloc)) do region=1,nregions imr=im(region) jmr=jm(region) do n=1,nmembersloc allocate(c13_oceandis_regxmem(region,n)%surf(imr,jmr)) end do enddo end if end subroutine declare_emission_c13_oceandis subroutine calc_emission_c13_oceandis 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_oceandis_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_oceandis' !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 ocean disequilibrium write(gol,'(a)') 'OCEANDIS data to be read '//trim(cdate) ; call goPr call ReadFromFile(trim(ocedis_prefix)//trim(cdate)//'.nc','oce_dis_clim',field2d_dis,status,netcdf4=.True.) if (status /= 0) call escape_tm('stop in '//rname) endif write(gol,'(a)') 'CarbonTracker OCEANDIS data to be recalculated ' ; call goPr c13_oceandis=field2d_dis(:,:) !kg 13CO2/cell/month call gridbox_per_month_to_flux(c13_oceandis) ! to kg 13CO2 /m2/s c13_oceandis=c13_oceandis*1.e3/xm13co2 ! to mol/m2/s ! accumulate weekly average if(allocated(c13_oceandis_m)) then if(flux_means_are_priors .eqv. .TRUE.) then c13_oceandis_m=c13_oceandis_m+c13_oceandis ! mol/m2/s else c13_oceandis_m=c13_oceandis_m+c13_oceandis !*fluxmultiplier(:,:,1) ! mol/m2/s, tracer 0 on PE0 is mean endif c13_oceandis_n=c13_oceandis_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_oceandis_mem=c13_oceandis ! *fluxmultiplier(:,:,mem) ! parameter applied, mol/m2/s c13_oceandis_mem=c13_oceandis_mem*xmc13/1.e3 ! from mol/m2/s to kgC13/m2/s call flux_to_gridbox_per_month(c13_oceandis_mem) ! to kgC13/area/month call coarsen_emission_2d('c13_oceandis',360,180,c13_oceandis_mem,c13_oceandis_regxmem(:,mem),add_field) enddo end subroutine calc_emission_c13_oceandis subroutine emission_apply_c13_oceandis(region) ! This routine should only be called by the forward enkf code. It ! applies the ocean c13 flux by member and region to the ic13_oceandis ! 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_oceandis, xmc13, xm13co2 implicit none integer, intent(in) :: region integer :: mem character(len=*), parameter :: rname = mname//'/emission_apply_c13_oceandis' if(enkf_is_forward()) then if(tracer_active(ic13_oceandis)) then call do_add_2d(region,ic13_oceandis,1,c13_oceandis_regxmem(region,1)%surf,xm13co2,xmc13) endif endif if(okdebug) write(*,*) 'end of emission_apply_c13_oceandis' end subroutine emission_apply_c13_oceandis subroutine free_emission_c13_oceandis use dims, only: nregions implicit none integer :: region,n character(len=*), parameter :: rname = mname//'/free_emission_c13_oceandis' if(enkf_is_forward()) then do region=1,nregions do n=1,nmembersloc deallocate(c13_oceandis_regxmem(region,n)%surf) end do enddo deallocate(c13_oceandis_regxmem) endif end subroutine free_emission_c13_oceandis end module emission_c13_oceandis