!################################################################# ! ! user input ! !### macro's ##################################################### ! #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 user_input use GO, only : gol, goPr, goErr implicit none ! --- in/out -------------------------- private public :: user_input_start ! --- const -------------------------------- character(len=*), parameter :: mname = 'user_input' contains ! ============================================================= ! === ! === start file ! === ! ============================================================= subroutine user_input_start( status ) use GO , only : gol, goPr, goErr use dims , only : nregions use global_data, only : mass_dat use meteodata, only : m_dat ! --- in/out ------------------------------ integer, intent(out) :: status ! --- const ------------------------------ character(len=*), parameter :: rname = mname//'/user_input_start' ! --- local ------------------------------- integer :: n ! --- begin -------------------------------- ! info ... call start_co2(status) IF_NOTOK_RETURN(status=1) call start_c13(status) IF_NOTOK_RETURN(status=1) ! ok status=0 end subroutine user_input_start subroutine start_co2(status) use GO, only : TrcFile, Init, Done, ReadRc use file_hdf , only : THdfFile, TSds, Init, Done, ReadAttribute, ReadData, GetInfo use toolbox, only : coarsen_emission use global_data, only : mass_dat use global_data, only : rcfile use meteodata, only : m_dat use global_types,only : d3_data use ParTools, only : ntracetloc, ntracet_ar, myid use Dims, only : im,jm,lm,nregions, adv_scheme, nlon360,nlat180,at,bt use emission_common,only : EmisInputDir use chem_param, only : xmc, xmco2, xmair, names use GO, only : goLoCase use meteodata , only : lli, levi, sp_dat, m_dat use grid , only : TllGridInfo, TLevelInfo,Init,Done,Fill3D implicit none !______________const______________________________ character(len=*), parameter :: rname = mname//'/start_co2' !______________I/O______________________________ integer, intent(out) :: status !______________local______________________________ type(TrcFile) :: rcF real,dimension(:,:,:),pointer ::m real,dimension(:,:,:,:),pointer ::rm,rxm,rym,rzm real, dimension(:,:,:), pointer :: sp real,dimension(:,:,:,:),allocatable :: co2_file type(d3_data),dimension(nregions),target :: co2_init integer,parameter :: avg_field=1 integer :: dimheight, dimlat, dimlon, data_rank integer :: region,l,imr,jmr,n, offsetn integer,dimension(4) :: data_dims real :: co2_multiplier type(THdfFile) :: hdf type(TSds) :: sds_co2 character(len=200) :: filename, varname type(TllGridInfo) :: lli_in type(TLevelInfo) :: levi_in integer :: xbeg_in, ybeg_in real :: dx_in, dy_in integer :: im_in, jm_in, lm_in real, allocatable :: at_in(:), bt_in(:) ! 34 (from standard) layer defintition: real, parameter :: at_tm5_std_34(0:34) = (/ & 0., 7.367743, 210.393890, 855.361755, & 2063.779785, 3850.913330, 6144.314941, 8802.356445, & 11632.758789, 14411.124023, 16899.468750, 17961.357422, & 18864.750000, 19584.330078, 20097.402344, 20384.480469, & 20429.863281, 20222.205078, 19755.109375, 19027.695313, & 18045.183594, 16819.474609, 15379.805664, 13775.325195, & 12077.446289, 10376.126953, 8765.053711, 7306.631348, & 6018.019531, 3960.291504, 1680.640259, 713.218079, & 298.495789, 95.636963, 0./) real, parameter :: bt_tm5_std_34(0:34) = (/ & 1.00000000, 0.99401945, 0.97966272, 0.95182151, & 0.90788388, 0.84737492, 0.77159661, 0.68326861, & 0.58616841, 0.48477158, 0.38389215, 0.33515489, & 0.28832296, 0.24393314, 0.20247590, 0.16438432, & 0.13002251, 0.09967469, 0.07353383, 0.05169041, & 0.03412116, 0.02067788, 0.01114291, 0.00508112, & 0.00181516, 0.00046139, 0.00007582, 0.00000000, & 0.00000000, 0.00000000, 0.00000000, 0.00000000, & 0.00000000, 0.00000000, 0.00000000/) write (gol,*) 'Initializing global mixing ratios from file' ; call goPr if(ntracetloc == 0) return do region=1,nregions imr=im(region) jmr=jm(region) allocate(co2_init(region)%d3(imr,jmr,lm(1))) enddo ! Get filename and variable name from rc-file filename = trim(EmisInputDir)//'/co2_init_2007.hdf' varname = 'co2_init_2007' call Init( rcF, rcfile, status ) IF_NOTOK_RETURN(status=1) ! put out settings ? call ReadRc( rcF, 'init.co2.filename', filename, status, default=filename ) IF_ERROR_RETURN(status=1) call ReadRc( rcF, 'init.co2.varname', varname, status, default=varname ) IF_ERROR_RETURN(status=1) call ReadRc( rcF, 'init.co2.multiply', co2_multiplier, status, default=1.0 ) IF_ERROR_RETURN(status=1) ! close rcfile: call Done( rcF, status ) IF_NOTOK_RETURN(status=1) write (gol,*) '.... Opening file with initial CO2 mixing ratios for reading: '//trim(filename) ; call goPr call Init( hdf, trim(filename), 'read', status ) IF_NOTOK_RETURN(status=1) call Init( sds_co2, hdf, trim(varname), status ) IF_NOTOK_RETURN(status=1) call GetInfo( sds_co2, status,data_rank=data_rank, data_dims=data_dims ) IF_NOTOK_RETURN(status=1) if(data_rank /= 3.and.data_rank /= 4) then write (gol,*) 'Rank of requested dataset ('//trim(varname)//') is not 3 or 4 (time,lon,lat,levels) ...error' ; call goPr IF_NOTOK_RETURN(status=1) endif allocate(co2_file(data_dims(1),data_dims(2),data_dims(3),1) ) if(data_rank.eq.3) then write (gol,*) '.... Retrieving variable '//trim(varname) ; call goPr endif if(data_rank.eq.4) then write (gol,*) '.... Retrieving first time step of variable '//trim(varname) ; call goPr endif call ReadData( sds_co2, co2_file, status) IF_NOTOK_RETURN(status=1) call Done( sds_co2, status ) IF_NOTOK_RETURN(status=1) write (gol,*) '.... Success' ; call goPr call Done( hdf, status ) IF_NOTOK_RETURN(status=1) write (gol,*) '.... File closed ' ; call goPr !WP! Multiply with specified factor co2_file = co2_file * co2_multiplier !WP Now grid the input file to the dimensions of the current model run !IV coarsen_emission does not regrid in the vertical, therefore must be replaced by Fill3D !call coarsen_emission('co2_init',data_dims(1),data_dims(2),data_dims(3),co2_file(:,:,:,1),co2_init,avg_field) ! not needed im_in=data_dims(1) jm_in=data_dims(2) lm_in=data_dims(3) !WP! We are assuming a 1x1 degree input file for now, this can be changed in the future if(data_dims(1) /= 360) then write (gol,'("input from initial concentration file implemented for 1 degree longitude only ...")'); call goErr TRACEBACK; status=1; return endif if(data_dims(2) /= 180) then write (gol,'("input from initial concentration file implemented for 1 degree latitude only ...")'); call goErr TRACEBACK; status=1; return endif !WP! We are assuming a 1x1 degree input file for now, this can be changed in the future xbeg_in = -180.0 dx_in = 360.0/im_in ybeg_in = -90.0 dy_in = 180.0/jm_in call Init( lli_in, xbeg_in+0.5*dx_in, dx_in, im_in, & ybeg_in+0.5*dy_in, dy_in, jm_in, status ) IF_NOTOK_RETURN(status=1) allocate( at_in(0:lm_in) ) allocate( bt_in(0:lm_in) ) select case ( lm_in ) case ( lm(1) ) !WP! If vertical dimension equals current setup, use at, bt from current run at_in = at bt_in = bt levi_in = levi case ( 34 ) at_in = at_tm5_std_34 bt_in = bt_tm5_std_34 call Init( levi_in, lm_in, at_in, bt_in, status ) IF_NOTOK_RETURN(status=1) case default write (gol,'("input from initial concentration file implemented for 25 or 34 levels only ...")'); call goErr TRACEBACK; status=1; return end select deallocate( at_in ) deallocate( bt_in ) write (gol,*) '.... Global mean CO2 mixing ratio in input field : ',sum(co2_file)/(im_in*jm_in*lm_in) ; call goPr do region=1,nregions m => m_dat(region)%data rm => mass_dat(region)%rm_t rxm => mass_dat(region)%rxm_t rym => mass_dat(region)%rym_t rzm => mass_dat(region)%rzm_t sp => sp_dat(region)%data call Fill3D( lli(region), levi, 'n', sp(1:im(region),1:jm(region),1), co2_init(region)%d3, & lli_in, levi_in, co2_file(:,:,:,1), 'mass-aver', status ) write (gol,*) '.... Mixing ratios coarsened to TM5 grid ' ; call goPr write (gol,*) '.... Mean CO2 mixing ratio in TM5 field region nr:',region,sum(co2_init(region)%d3)/(im(region)*jm(region)*lm(region)) ; call goPr do n=1,ntracetloc ! There is a bug with TM5 using the reduced grid. Apparently ! negative concentrations are not permitted in this scenario. ! Thus for the forward simulation, we start each component ! with the background CO2 field, on which sinks can act. To ! get the component concentration, the bg field must be ! subtracted. !WP! Get the global tracer number, we only initialize ico2_bg (if we !WP! run fwd) or we initialize all co2 tracers (in case of inv) offsetn = sum(ntracet_ar(0:myid-1)) if (goLoCase(trim(names(offsetn+n)))=='co2'.or.goLoCase(trim(names(offsetn+n)))=='co2_bg') then if (region==1) then write (gol,'(a,i4,a,a)') '.... Found tracer ',offsetn+n,' with name: ',trim(names(offsetn+n)) ; call goPr write (gol,'(a)') '.... Initializing with CO2 from file' ; call goPr endif rm(1:im(region),1:jm(region),1:lm(region),n) = & co2_init(region)%d3*1.e-6*m(1:im(region),1:jm(region),1:lm(region))*xmco2/xmair ! else ! rm(:,:,:,:) = 0.0 endif enddo #ifdef slopes rxm(:,:,:,:) = 0.0 rym(:,:,:,:) = 0.0 rzm(:,:,:,:) = 0.0 #endif nullify(m) nullify(rm) nullify(rxm) nullify(rym) nullify(rzm) nullify(sp) enddo deallocate(co2_file) do region=1,nregions deallocate(co2_init(region)%d3) enddo write (gol,'(a,a)') rname ,' Done' ; call goPr end subroutine start_co2 subroutine start_c13(status) use GO, only : TrcFile, Init, Done, ReadRc use file_hdf , only : THdfFile, TSds, Init, Done, ReadAttribute, ReadData, GetInfo use toolbox, only : coarsen_emission use global_data, only : mass_dat use global_data, only : rcfile use meteodata, only : m_dat use global_types,only : d3_data use ParTools, only : ntracetloc, ntracet_ar, myid use Dims, only : im,jm,lm,nregions, adv_scheme, nlon360,nlat180,at,bt use emission_common,only : C13InputDir,EmisInputDir use chem_param, only : xmc13, xm13co2,xmco2, xmair, names use GO, only : goLoCase use meteodata , only : lli, levi, sp_dat, m_dat use grid , only : TllGridInfo, TLevelInfo,Init,Done,Fill3D implicit none !______________const______________________________ character(len=*), parameter :: rname = mname//'/start_c13' !______________I/O______________________________ integer, intent(out) :: status !______________local______________________________ type(TrcFile) :: rcF real,dimension(:,:,:),pointer ::m real,dimension(:,:,:,:),pointer ::rm,rxm,rym,rzm real, dimension(:,:,:), pointer :: sp real,dimension(:,:,:,:),allocatable :: c13_file,co2_file type(d3_data),dimension(nregions),target :: c13_init,co2_init integer,parameter :: avg_field=1 integer :: dimheight, dimlat, dimlon, data_rank integer :: region,l,imr,jmr,n, offsetn integer,dimension(4) :: data_dims real :: c13_multiplier,co2_multiplier type(THdfFile) :: hdf type(TSds) :: sds_c13,sds_co2 character(len=200) :: filename, varname type(TllGridInfo) :: lli_in type(TLevelInfo) :: levi_in integer :: xbeg_in, ybeg_in real :: dx_in, dy_in integer :: im_in, jm_in, lm_in real, allocatable :: at_in(:), bt_in(:) ! 34 (from standard) layer defintition: real, parameter :: at_tm5_std_34(0:34) = (/ & 0., 7.367743, 210.393890, 855.361755, & 2063.779785, 3850.913330, 6144.314941, 8802.356445, & 11632.758789, 14411.124023, 16899.468750, 17961.357422, & 18864.750000, 19584.330078, 20097.402344, 20384.480469, & 20429.863281, 20222.205078, 19755.109375, 19027.695313, & 18045.183594, 16819.474609, 15379.805664, 13775.325195, & 12077.446289, 10376.126953, 8765.053711, 7306.631348, & 6018.019531, 3960.291504, 1680.640259, 713.218079, & 298.495789, 95.636963, 0./) real, parameter :: bt_tm5_std_34(0:34) = (/ & 1.00000000, 0.99401945, 0.97966272, 0.95182151, & 0.90788388, 0.84737492, 0.77159661, 0.68326861, & 0.58616841, 0.48477158, 0.38389215, 0.33515489, & 0.28832296, 0.24393314, 0.20247590, 0.16438432, & 0.13002251, 0.09967469, 0.07353383, 0.05169041, & 0.03412116, 0.02067788, 0.01114291, 0.00508112, & 0.00181516, 0.00046139, 0.00007582, 0.00000000, & 0.00000000, 0.00000000, 0.00000000, 0.00000000, & 0.00000000, 0.00000000, 0.00000000/) write (gol,*) 'Initializing global mixing ratios from file' ; call goPr if(ntracetloc == 0) return do region=1,nregions imr=im(region) jmr=jm(region) allocate(co2_init(region)%d3(imr,jmr,lm(1))) enddo ! Get filename and variable name from rc-file filename = trim(EmisInputDir)//'/co2_init_2007.hdf' varname = 'co2_init_2007' call Init( rcF, rcfile, status ) IF_NOTOK_RETURN(status=1) ! put out settings ? call ReadRc( rcF, 'init.co2.filename', filename, status, default=filename ) IF_ERROR_RETURN(status=1) call ReadRc( rcF, 'init.co2.varname', varname, status, default=varname ) IF_ERROR_RETURN(status=1) call ReadRc( rcF, 'init.co2.multiply', co2_multiplier, status, default=1.0 ) IF_ERROR_RETURN(status=1) ! close rcfile: call Done( rcF, status ) IF_NOTOK_RETURN(status=1) write (gol,*) '.... Opening file with initial CO2 mixing ratios for reading: '//trim(filename) ; call goPr call Init( hdf, trim(filename), 'read', status ) IF_NOTOK_RETURN(status=1) call Init( sds_co2, hdf, trim(varname), status ) IF_NOTOK_RETURN(status=1) call GetInfo( sds_co2, status,data_rank=data_rank, data_dims=data_dims ) IF_NOTOK_RETURN(status=1) if(data_rank /= 3.and.data_rank /= 4) then write (gol,*) 'Rank of requested dataset ('//trim(varname)//') is not 3 or 4 (time,lon,lat,levels) ...error' ; call goPr IF_NOTOK_RETURN(status=1) endif allocate(co2_file(data_dims(1),data_dims(2),data_dims(3),1) ) if(data_rank.eq.3) then write (gol,*) '.... Retrieving variable '//trim(varname) ; call goPr endif if(data_rank.eq.4) then write (gol,*) '.... Retrieving first time step of variable '//trim(varname) ; call goPr endif call ReadData( sds_co2, co2_file, status) IF_NOTOK_RETURN(status=1) call Done( sds_co2, status ) IF_NOTOK_RETURN(status=1) write (gol,*) '.... Success' ; call goPr call Done( hdf, status ) IF_NOTOK_RETURN(status=1) write (gol,*) '.... File closed ' ; call goPr !WP! Multiply with specified factor co2_file = co2_file * co2_multiplier !WP Now grid the input file to the dimensions of the current model run !IV coarsen_emission does not regrid in the vertical, therefore must be replaced by Fill3D !call coarsen_emission('co2_init',data_dims(1),data_dims(2),data_dims(3),co2_file(:,:,:,1),co2_init,avg_field) ! not needed im_in=data_dims(1) jm_in=data_dims(2) lm_in=data_dims(3) !WP! We are assuming a 1x1 degree input file for now, this can be changed in the future if(data_dims(1) /= 360) then write (gol,'("input from initial concentration file implemented for 1 degree longitude only ...")'); call goErr TRACEBACK; status=1; return endif if(data_dims(2) /= 180) then write (gol,'("input from initial concentration file implemented for 1 degree latitude only ...")'); call goErr TRACEBACK; status=1; return endif !WP! We are assuming a 1x1 degree input file for now, this can be changed in the future xbeg_in = -180.0 dx_in = 360.0/im_in ybeg_in = -90.0 dy_in = 180.0/jm_in call Init( lli_in, xbeg_in+0.5*dx_in, dx_in, im_in, & ybeg_in+0.5*dy_in, dy_in, jm_in, status ) IF_NOTOK_RETURN(status=1) allocate( at_in(0:lm_in) ) allocate( bt_in(0:lm_in) ) select case ( lm_in ) case ( lm(1) ) !WP! If vertical dimension equals current setup, use at, bt from current run at_in = at bt_in = bt levi_in = levi case ( 34 ) at_in = at_tm5_std_34 bt_in = bt_tm5_std_34 call Init( levi_in, lm_in, at_in, bt_in, status ) IF_NOTOK_RETURN(status=1) case default write (gol,'("input from initial concentration file implemented for 25 or 34 levels only ...")'); call goErr TRACEBACK; status=1; return end select deallocate( at_in ) deallocate( bt_in ) ! print*,sum(co2_file),im_in,jm_in,lm_in write (gol,*) '.... Global mean CO2 mixing ratio in input field : ',sum(co2_file)/(im_in*jm_in*lm_in) ; call goPr do region=1,nregions m => m_dat(region)%data rm => mass_dat(region)%rm_t rxm => mass_dat(region)%rxm_t rym => mass_dat(region)%rym_t rzm => mass_dat(region)%rzm_t sp => sp_dat(region)%data call Fill3D( lli(region), levi, 'n', sp(1:im(region),1:jm(region),1), co2_init(region)%d3, & lli_in, levi_in, co2_file(:,:,:,1), 'mass-aver', status ) write (gol,*) '.... Mixing ratios coarsened to TM5 grid ' ; call goPr write (gol,*) '.... Mean CO2 mixing ratio in TM5 field region nr:',region,sum(co2_init(region)%d3)/(im(region)*jm(region)*lm(region)) ; call goPr #ifdef slopes rxm(:,:,:,:) = 0.0 rym(:,:,:,:) = 0.0 rzm(:,:,:,:) = 0.0 #endif nullify(m) nullify(rm) nullify(rxm) nullify(rym) nullify(rzm) nullify(sp) enddo write (gol,*) 'Initializing global 13C mixing ratios from file' ; call goPr if(ntracetloc == 0) return do region=1,nregions imr=im(region) jmr=jm(region) allocate(c13_init(region)%d3(imr,jmr,lm(1))) enddo ! Get filename and variable name from rc-file filename = trim(C13InputDir)//'/d13c_init_2007_offset2.hdf' varname = 'd13c_init_2007' call Init( rcF, rcfile, status ) IF_NOTOK_RETURN(status=1) ! put out settings ? call ReadRc( rcF, 'init.c13.filename', filename, status, default=filename ) IF_ERROR_RETURN(status=1) call ReadRc( rcF, 'init.c13.varname', varname, status, default=varname ) IF_ERROR_RETURN(status=1) call ReadRc( rcF, 'init.c13.multiply', c13_multiplier, status, default=1.0 ) IF_ERROR_RETURN(status=1) ! close rcfile: call Done( rcF, status ) IF_NOTOK_RETURN(status=1) write (gol,*) '.... Opening file with initial d13C ratios for reading: '//trim(filename) ; call goPr call Init( hdf, trim(filename), 'read', status ) IF_NOTOK_RETURN(status=1) call Init( sds_c13, hdf, trim(varname), status ) IF_NOTOK_RETURN(status=1) call GetInfo( sds_c13, status,data_rank=data_rank, data_dims=data_dims ) IF_NOTOK_RETURN(status=1) if(data_rank /= 3.and.data_rank /= 4) then write (gol,*) 'Rank of requested dataset ('//trim(varname)//') is not 3 or 4 (time,lon,lat,levels) ...error' ; call goPr IF_NOTOK_RETURN(status=1) endif allocate(c13_file(data_dims(1),data_dims(2),data_dims(3),1) ) if(data_rank.eq.3) then write (gol,*) '.... Retrieving variable '//trim(varname) ; call goPr endif if(data_rank.eq.4) then write (gol,*) '.... Retrieving first time step of variable '//trim(varname) ; call goPr endif call ReadData( sds_c13, c13_file, status) IF_NOTOK_RETURN(status=1) call Done( sds_c13, status ) IF_NOTOK_RETURN(status=1) write (gol,*) '.... Success' ; call goPr call Done( hdf, status ) IF_NOTOK_RETURN(status=1) write (gol,*) '.... File closed ' ; call goPr !WP! Multiply with specified factor c13_file = c13_file * c13_multiplier !WP Now grid the input file to the dimensions of the current model run !IV coarsen_emission does not regrid in the vertical, therefore must be replaced by Fill3D !call coarsen_emission('co2_init',data_dims(1),data_dims(2),data_dims(3),co2_file(:,:,:,1),co2_init,avg_field) ! not needed im_in=data_dims(1) jm_in=data_dims(2) lm_in=data_dims(3) !WP! We are assuming a 1x1 degree input file for now, this can be changed in the future if(data_dims(1) /= 360) then write (gol,'("input from initial concentration file implemented for 1 degree longitude only ...")'); call goErr TRACEBACK; status=1; return endif if(data_dims(2) /= 180) then write (gol,'("input from initial concentration file implemented for 1 degree latitude only ...")'); call goErr TRACEBACK; status=1; return endif !WP! We are assuming a 1x1 degree input file for now, this can be changed in the future xbeg_in = -180.0 dx_in = 360.0/im_in ybeg_in = -90.0 dy_in = 180.0/jm_in call Init( lli_in, xbeg_in+0.5*dx_in, dx_in, im_in, & ybeg_in+0.5*dy_in, dy_in, jm_in, status ) IF_NOTOK_RETURN(status=1) allocate( at_in(0:lm_in) ) allocate( bt_in(0:lm_in) ) select case ( lm_in ) case ( lm(1) ) !WP! If vertical dimension equals current setup, use at, bt from current run at_in = at bt_in = bt levi_in = levi case ( 34 ) at_in = at_tm5_std_34 bt_in = bt_tm5_std_34 call Init( levi_in, lm_in, at_in, bt_in, status ) IF_NOTOK_RETURN(status=1) case default write (gol,'("input from initial concentration file implemented for 25 or 34 levels only ...")'); call goErr TRACEBACK; status=1; return end select deallocate( at_in ) deallocate( bt_in ) !print*,sum(c13_file),im_in,jm_in,lm_in write (gol,*) '.... Global mean 13C mixing ratio in input field : ',sum(c13_file)/(im_in*jm_in*lm_in) ; call goPr do region=1,nregions m => m_dat(region)%data rm => mass_dat(region)%rm_t rxm => mass_dat(region)%rxm_t rym => mass_dat(region)%rym_t rzm => mass_dat(region)%rzm_t sp => sp_dat(region)%data call Fill3D( lli(region), levi, 'n', sp(1:im(region),1:jm(region),1), c13_init(region)%d3, & lli_in, levi_in, c13_file(:,:,:,1), 'mass-aver', status ) write (gol,*) '.... Mean 13C mixing ratio in TM5 field region nr:',region,sum(c13_init(region)%d3)/(im(region)*jm(region)*lm(region)) ; call goPr do n=1,ntracetloc ! There is a bug with TM5 using the reduced grid. Apparently ! negative concentrations are not permitted in this scenario. ! Thus for the forward simulation, we start each component ! with the background C13 field, on which sinks can act. To ! get the component concentration, the bg field must be ! subtracted. !WP! Get the global tracer number, we only initialize ic13_bg (if we !WP! run fwd) or we initialize all c13 tracers (in case of inv) offsetn = sum(ntracet_ar(0:myid-1)) !print*,myid,n,offsetn+n,trim(names(offsetn+n)) if (goLoCase(trim(names(offsetn+n)))=='co2c13'.or.goLoCase(trim(names(offsetn+n)))=='co2c13_bg') then if (region==1) then !print*,'.... Found tracer ',offsetn+n,' with name: ',trim(names(offsetn+n)) write (gol,'(a,i4,a,a)') '.... Found tracer ',offsetn+n,' with name: ',trim(names(offsetn+n)) ; call goPr write (gol,'(a)') '.... Initializing with C13 from file' ; call goPr endif rm(1:im(region),1:jm(region),1:lm(region),n) = & ((c13_init(region)%d3/1000.)+1.)*0.011112*co2_init(region)%d3*1.e-6*m(1:im(region),1:jm(region),1:lm(region))*xm13co2/xmair ! else ! rm(:,:,:,:) = 0.0 endif enddo #ifdef slopes rxm(:,:,:,:) = 0.0 rym(:,:,:,:) = 0.0 rzm(:,:,:,:) = 0.0 #endif nullify(m) nullify(rm) nullify(rxm) nullify(rym) nullify(rzm) nullify(sp) enddo deallocate(co2_file) deallocate(c13_file) if (myid==0) write (gol,'(a)') '.... All other tracers initialized to 0.0! ' ; call goPr do region=1,nregions deallocate(c13_init(region)%d3) enddo write (gol,'(a,a)') rname ,' Done' ; call goPr end subroutine start_c13 end module user_input