program SPARRING ! reading from NetCDF files and passing the data to OASIS ! units for the climatic parameters: ! temperature (K) ! precipitation (m s-1) ! radiation (W m-2) ! Last update: 10 Dec 2015 - add 2 new fields for low veg. - see "ecev3" ! Now reads AMIP output. use mod_oasis use mpi implicit none ! integer,parameter :: MAXGRID = 35718 integer,parameter :: MAXGRID = 88838 integer,parameter :: num_years = 10 integer,parameter :: start_year = 1850 integer,parameter :: num_loops = 10 integer,parameter :: t_max_year = 365 integer,parameter :: t_max_leap = 366 character(len=6),parameter :: model_name = "SPARRI" double precision :: minhlai,maxhlai,minllai,maxllai double precision :: meanhlai,meanllai,meanvegl,meanvegh double precision :: maxcvh,maxcvl,maxtvl,maxtvh,cover ! PM, was SPATEMPE, SPAPRECI, SPARADIA, SPATOLAI character(len=128),parameter :: f_temp_name = "T2MVeg" integer :: f_temp_id = -77 double precision :: f_temp_data(MAXGRID,1) character(len=128),parameter :: f_prec_name = "TPVeg" integer :: f_prec_id = -77 double precision :: f_prec_data(MAXGRID,1) double precision :: f_precs_data(MAXGRID,1) double precision :: f_precc_data(MAXGRID,1) double precision :: f_snowfall_data(MAXGRID,1) character(len=128),parameter :: f_vegl_name = "CVLVeg" integer :: f_vegl_id = -77 double precision :: f_vegl_data(MAXGRID,1) character(len=128),parameter :: f_vegltype_name = "TVLVeg" integer :: f_vegltype_id = -77 double precision :: f_vegltype_data(MAXGRID,1) character(len=128),parameter :: f_vegh_name = "CVHVeg" integer :: f_vegh_id = -77 double precision :: f_vegh_data(MAXGRID,1) character(len=128),parameter :: f_snoc_name = "SDVeg" integer :: f_snoc_id = -77 double precision :: f_snoc_data(MAXGRID,1) character(len=128),parameter :: f_snod_name = "SDensVeg" integer :: f_snod_id = -77 double precision :: f_snod_data(MAXGRID,1) character(len=128),parameter :: f_st1l_name = "SoilTVeg.L001" integer :: f_st1l_id = -77 double precision :: f_st1l_data(MAXGRID,1) character(len=128),parameter :: f_st2l_name = "SoilTVeg.L002" integer :: f_st2l_id = -77 double precision :: f_st2l_data(MAXGRID,1) character(len=128),parameter :: f_st3l_name = "SoilTVeg.L003" integer :: f_st3l_id = -77 double precision :: f_st3l_data(MAXGRID,1) character(len=128),parameter :: f_st4l_name = "SoilTVeg.L004" integer :: f_st4l_id = -77 double precision :: f_st4l_data(MAXGRID,1) character(len=128),parameter :: f_sm1l_name = "SoilMVeg.L001" integer :: f_sm1l_id = -77 double precision :: f_sm1l_data(MAXGRID,1) character(len=128),parameter :: f_sm2l_name = "SoilMVeg.L002" integer :: f_sm2l_id = -77 double precision :: f_sm2l_data(MAXGRID,1) character(len=128),parameter :: f_sm3l_name = "SoilMVeg.L003" integer :: f_sm3l_id = -77 double precision :: f_sm3l_data(MAXGRID,1) character(len=128),parameter :: f_sm4l_name = "SoilMVeg.L004" integer :: f_sm4l_id = -77 double precision :: f_sm4l_data(MAXGRID,1) character(len=128),parameter :: f_radi_name = "SSRVeg" integer :: f_radi_id = -77 double precision :: f_radi_data(MAXGRID,1) character(len=128),parameter :: f_llai_name = "LAILVeg" integer :: f_llai_id = -77 double precision :: f_llai_data(MAXGRID,1) character(len=128),parameter :: f_hlai_name = "LAIHVeg" integer :: f_hlai_id = -77 double precision :: f_hlai_data(MAXGRID,1) character(len=128),parameter :: f_typh_name = "TypeHVeg" integer :: f_typh_id = -77 double precision :: f_typh_data(MAXGRID,1) character(len=128),parameter :: f_vggh_name = "FracHVeg" integer :: f_vggh_id = -77 double precision :: f_vggh_data(MAXGRID,1) ! ecev3 character(len=128),parameter :: f_typl_name = "TypeLVeg" integer :: f_typl_id = -77 double precision :: f_typl_data(MAXGRID,1) ! ecev3 character(len=128),parameter :: f_vggl_name = "FracLVeg" integer :: f_vggl_id = -77 double precision :: f_vggl_data(MAXGRID,1) integer :: comp_id = -77 integer :: part_id = -77 integer :: t,t_step,yr,lyr,t_max,t_step_full,cell,loopyr integer :: ierror integer::ix,iy ! character(len=*),parameter::nc_temp_file="/nobackup/rossby15/sm_paumi/preind/dpic_var167+1849.nc" ! character(len=*),parameter::nc_prec_file="/nobackup/rossby15/sm_paumi/preind/dpic_var142+1849.nc" ! character(len=*),parameter::nc_radi_file="/nobackup/rossby15/sm_paumi/preind/dpic_var169+1849.nc" ! VARIABLE FIELDS character(len=100) :: nc_temp_file character(len=100) :: nc_precs_file character(len=100) :: nc_precc_file character(len=100) :: nc_snowfall_file character(len=100) :: nc_st1l_file character(len=100) :: nc_st2l_file character(len=100) :: nc_st3l_file character(len=100) :: nc_st4l_file character(len=100) :: nc_sm1l_file character(len=100) :: nc_sm2l_file character(len=100) :: nc_sm3l_file character(len=100) :: nc_sm4l_file character(len=100) :: nc_radi_file ! FIXED FIELDS - T159!!!! character(len=*),parameter::nc_vegltype_file="/nobackup/rossby15/sm_paumi/preind/tvl.nc" character(len=*),parameter::nc_veghtype_file="/nobackup/rossby15/sm_paumi/preind/tvh.nc" ! STEM FOR VARIABLE FIELDS ! character(len=*),parameter::nc_file_path="/nobackup/rossby15/sm_paumi/ecev3/processed/" ! NEW, AMIP OUTPUT: character(len=*),parameter::nc_file_path="/nobackup/rossby18/sm_wyser/ecearth3-amip/processed/daily/" character(len=*),parameter::nc_file_tail="_dayavg.nc" logical :: isleapyear = .false. character yearstr*4 ! VARIABLE NAMES character(len=*),parameter::nc_temp_name="var167" character(len=*),parameter::nc_precs_name="var142" character(len=*),parameter::nc_precc_name="var143" character(len=*),parameter::nc_snowfall_name="var144" ! We do not have the SNOC and SNOD data, so we send fixed values to guess below. character(len=*),parameter::nc_st1l_name="var139" character(len=*),parameter::nc_st2l_name="var170" ! ST3L and ST4L - we do not have this data, so we'll send ST2L values instead - see below. character(len=*),parameter::nc_st3l_name="var183" character(len=*),parameter::nc_st4l_name="var236" character(len=*),parameter::nc_sm1l_name="var039" ! AMIP files have this character(len=*),parameter::nc_sm2l_name="var040" character(len=*),parameter::nc_sm3l_name="var041" character(len=*),parameter::nc_sm4l_name="var042" character(len=*),parameter::nc_sm1l_varname="var39" character(len=*),parameter::nc_sm2l_varname="var40" character(len=*),parameter::nc_sm3l_varname="var41" character(len=*),parameter::nc_sm4l_varname="var42" character(len=*),parameter::nc_radi_name="var176" character(len=*),parameter::nc_vegl_name="tvl" character(len=*),parameter::nc_vegh_name="tvh" ! FIELD IDS integer::nc_temp_fileid,nc_precs_fileid,nc_radi_fileid,nc_vegl_fileid,nc_vegh_fileid,nc_precc_fileid,nc_snowfall_fileid; integer::nc_temp_varid,nc_precs_varid,nc_radi_varid,nc_vegl_varid,nc_vegh_varid,nc_precc_varid,nc_snowfall_varid; integer::nc_vegltype_varid; ! We don't have a T255 file for vegltype yet integer::nc_st1l_fileid,nc_st2l_fileid,nc_st3l_fileid,nc_st4l_fileid; integer::nc_st1l_varid,nc_st2l_varid,nc_st3l_varid,nc_st4l_varid; integer::nc_sm1l_fileid,nc_sm2l_fileid,nc_sm3l_fileid,nc_sm4l_fileid; integer::nc_sm1l_varid,nc_sm2l_varid,nc_sm3l_varid,nc_sm4l_varid; integer::localcomm, cplcomm, icpl; ! *** START *** write (*,'(A)') "*II* sparring: Hello" write (*,'(A)') "*II* sparring: Now initialising Sparring using oasis_..." call oasis_init_comp(comp_id,model_name,ierror) write (*,'(A,I3)') "*II* sparring: oasis_init_comp returned ierror=",ierror ! call oasis_get_localcomm(localcomm, ierror) ! write (*,'(A,I3)') "*II* sparring: get_localcomm ierror=",ierror ! write (*,'(A,I12)') "*II* sparring: oasis_get_localcomm returned localcomm =",localcomm icpl = 1 call oasis_create_couplcomm(icpl,MPI_COMM_SELF, cplcomm, ierror) ! call oasis_create_couplcomm(icpl,localcomm, cplcomm, ierror) write (*,'(A,I3)') "*II* sparring: oasis_create_couplcomm ierror=",ierror write (*,'(A,I12)') "*II* sparring: oasis_create_couplcomm returned cplcomm =",cplcomm call oasis_def_partition(part_id,(/ 0,0,MAXGRID*1 /),ierror) write (*,'(A,I3)') "*II* sparring: oasis_def_partition returned part_id =",part_id write (*,'(A,I3)') "*II* sparring: oasis_def_partition returned ierror =",ierror call oasis_def_var( f_temp_id, & f_temp_name, & part_id, & (/ 2,1 /), & PRISM_Out, & (/ 1,MAXGRID,1,1 /), & PRISM_Real, & ierror ) write (*,'(A,I3)') "*II* sparring: oasis_def_var returned f_temp_id =",f_temp_id write (*,'(A,I3)') "*II* sparring: oasis_def_var returned ierror =",ierror call oasis_def_var( f_prec_id, & f_prec_name, & part_id, & (/ 2,1 /), & PRISM_Out, & (/ 1,MAXGRID,1,1 /), & PRISM_Real, & ierror ) write (*,'(A,I3)') "*II* sparring: oasis_def_var returned f_prec_id =",f_prec_id write (*,'(A,I3)') "*II* sparring: oasis_def_var returned ierror =",ierror call oasis_def_var( f_vegl_id, & f_vegl_name, & part_id, & (/ 2,1 /), & PRISM_Out, & (/ 1,MAXGRID,1,1 /), & PRISM_Real, & ierror ) write (*,'(A,I3)') "*II* sparring: oasis_def_var returned f_vegl_id =",f_vegl_id write (*,'(A,I3)') "*II* sparring: oasis_def_var returned ierror =",ierror call oasis_def_var( f_vegh_id, & f_vegh_name, & part_id, & (/ 2,1 /), & PRISM_Out, & (/ 1,MAXGRID,1,1 /), & PRISM_Real, & ierror ) write (*,'(A,I3)') "*II* sparring: oasis_def_var returned f_vegh_id =",f_vegh_id write (*,'(A,I3)') "*II* sparring: oasis_def_var returned ierror =",ierror call oasis_def_var( f_snoc_id, & f_snoc_name, & part_id, & (/ 2,1 /), & PRISM_Out, & (/ 1,MAXGRID,1,1 /), & PRISM_Real, & ierror ) write (*,'(A,I3)') "*II* sparring: oasis_def_var returned f_snoc_id =",f_snoc_id write (*,'(A,I3)') "*II* sparring: oasis_def_var returned ierror =",ierror call oasis_def_var( f_snod_id, & f_snod_name, & part_id, & (/ 2,1 /), & PRISM_Out, & (/ 1,MAXGRID,1,1 /), & PRISM_Real, & ierror ) write (*,'(A,I3)') "*II* sparring: oasis_def_var returned f_snod_id =",f_snod_id write (*,'(A,I3)') "*II* sparring: oasis_def_var returned ierror =",ierror call oasis_def_var( f_st1l_id, & f_st1l_name, & part_id, & (/ 2,1 /), & PRISM_Out, & (/ 1,MAXGRID,1,1 /), & PRISM_Real, & ierror ) write (*,'(A,I3)') "*II* sparring: oasis_def_var returned f_st1l_id =",f_st1l_id write (*,'(A,I3)') "*II* sparring: oasis_def_var returned ierror =",ierror call oasis_def_var( f_st2l_id, & f_st2l_name, & part_id, & (/ 2,1 /), & PRISM_Out, & (/ 1,MAXGRID,1,1 /), & PRISM_Real, & ierror ) write (*,'(A,I3)') "*II* sparring: oasis_def_var returned f_st2l_id =",f_st2l_id write (*,'(A,I3)') "*II* sparring: oasis_def_var returned ierror =",ierror call oasis_def_var( f_st3l_id, & f_st3l_name, & part_id, & (/ 2,1 /), & PRISM_Out, & (/ 1,MAXGRID,1,1 /), & PRISM_Real, & ierror ) write (*,'(A,I3)') "*II* sparring: oasis_def_var returned f_st3l_id =",f_st3l_id write (*,'(A,I3)') "*II* sparring: oasis_def_var returned ierror =",ierror call oasis_def_var( f_st4l_id, & f_st4l_name, & part_id, & (/ 2,1 /), & PRISM_Out, & (/ 1,MAXGRID,1,1 /), & PRISM_Real, & ierror ) write (*,'(A,I3)') "*II* sparring: oasis_def_var returned f_st4l_id =",f_st4l_id write (*,'(A,I3)') "*II* sparring: oasis_def_var returned ierror =",ierror call oasis_def_var( f_sm1l_id, & f_sm1l_name, & part_id, & (/ 2,1 /), & PRISM_Out, & (/ 1,MAXGRID,1,1 /), & PRISM_Real, & ierror ) write (*,'(A,I3)') "*II* sparring: oasis_def_var returned f_sm1l_id =",f_sm1l_id write (*,'(A,I3)') "*II* sparring: oasis_def_var returned ierror =",ierror call oasis_def_var( f_sm2l_id, & f_sm2l_name, & part_id, & (/ 2,1 /), & PRISM_Out, & (/ 1,MAXGRID,1,1 /), & PRISM_Real, & ierror ) write (*,'(A,I3)') "*II* sparring: oasis_def_var returned f_sm2l_id =",f_sm2l_id write (*,'(A,I3)') "*II* sparring: oasis_def_var returned ierror =",ierror call oasis_def_var( f_sm3l_id, & f_sm3l_name, & part_id, & (/ 2,1 /), & PRISM_Out, & (/ 1,MAXGRID,1,1 /), & PRISM_Real, & ierror ) write (*,'(A,I3)') "*II* sparring: oasis_def_var returned f_sm3l_id =",f_sm3l_id write (*,'(A,I3)') "*II* sparring: oasis_def_var returned ierror =",ierror call oasis_def_var( f_sm4l_id, & f_sm4l_name, & part_id, & (/ 2,1 /), & PRISM_Out, & (/ 1,MAXGRID,1,1 /), & PRISM_Real, & ierror ) write (*,'(A,I3)') "*II* sparring: oasis_def_var returned f_sm4l_id =",f_sm4l_id write (*,'(A,I3)') "*II* sparring: oasis_def_var returned ierror =",ierror call oasis_def_var( f_radi_id, & f_radi_name, & part_id, & (/ 2,1 /), & PRISM_Out, & (/ 1,MAXGRID,1,1 /), & PRISM_Real, & ierror ) write (*,'(A,I3)') "*II* sparring: oasis_def_var returned f_radi_id =",f_radi_id write (*,'(A,I3)') "*II* sparring: oasis_def_var returned ierror =",ierror ! VEGLTYPE call oasis_def_var( f_vegltype_id, & f_vegltype_name, & part_id, & (/ 2,1 /), & PRISM_Out, & (/ 1,MAXGRID,1,1 /), & PRISM_Real, & ierror ) write (*,'(A,I3)') "*II* sparring: f_vegltype_id =",f_vegltype_id write (*,'(A,I3)') "*II* sparring: returned ierror=",ierror ! LLAI call oasis_def_var( f_llai_id, & f_llai_name, & part_id, & (/ 2,1 /), & PRISM_In, & (/ 1,MAXGRID,1,1 /), & PRISM_Real, & ierror ) write (*,'(A,I3)') "*II* sparring: oasis_def_var returned f_llai_id =",f_llai_id write (*,'(A,I3)') "*II* sparring: oasis_def_var returned ierror =",ierror ! HLAI call oasis_def_var( f_hlai_id, & f_hlai_name, & part_id, & (/ 2,1 /), & PRISM_In, & (/ 1,MAXGRID,1,1 /), & PRISM_Real, & ierror ) write (*,'(A,I3)') "*II* sparring: oasis_def_var returned f_hlai_id =",f_hlai_id write (*,'(A,I3)') "*II* sparring: oasis_def_var returned ierror =",ierror ! TYPH call oasis_def_var( f_typh_id, & f_typh_name, & part_id, & (/ 2,1 /), & PRISM_In, & (/ 1,MAXGRID,1,1 /), & PRISM_Real, & ierror ) write (*,'(A,I3)') "*II* sparring: oasis_def_var returned f_typh_id =",f_typh_id write (*,'(A,I3)') "*II* sparring: oasis_def_var returned ierror =",ierror ! VEGH call oasis_def_var( f_vggh_id, & f_vggh_name, & part_id, & (/ 2,1 /), & PRISM_In, & (/ 1,MAXGRID,1,1 /), & PRISM_Real, & ierror ) write (*,'(A,I3)') "*II* sparring: oasis_def_var returned f_vggh_id =",f_vggh_id write (*,'(A,I3)') "*II* sparring: oasis_def_var returned ierror =",ierror ! ecev3 TYPL call oasis_def_var( f_typl_id, & f_typl_name, & part_id, & (/ 2,1 /), & PRISM_In, & (/ 1,MAXGRID,1,1 /), & PRISM_Real, & ierror ) write (*,'(A,I3)') "*II* sparring: oasis_def_var returned f_typl_id =",f_typl_id write (*,'(A,I3)') "*II* sparring: oasis_def_var returned ierror =",ierror ! ecev3 VEGL call oasis_def_var( f_vggl_id, & f_vggl_name, & part_id, & (/ 2,1 /), & PRISM_In, & (/ 1,MAXGRID,1,1 /), & PRISM_Real, & ierror ) write (*,'(A,I3)') "*II* sparring: oasis_def_var returned f_vggl_id =",f_vggl_id write (*,'(A,I3)') "*II* sparring: oasis_def_var returned ierror =",ierror write (*,'(A)') "*II* sparring: before call oasis_enddef(ierror)" call oasis_enddef(ierror) write (*,'(A,I3)') "*II* sparring: oasis_enddef returned ierror =",ierror ! FIXED FIELDS write (*,'(A)') "*II* sparring: Reading vegl and vegh" ! VEGL ! call read_fixedfiels(MAXGRID,1,nc_vegl_file,nc_vegl_name,& ! f_vegl_data,nc_vegl_fileid,nc_vegl_varid) ! VEGH ! call read_fixedfiels(MAXGRID,1,nc_vegh_file,nc_vegh_name,& ! f_vegh_data,nc_vegh_fileid,nc_vegh_varid) write (*,'(A)') "*II* sparring: Beginning time loop" t_step_full = 0 ! Repetitions of forcing do lyr = 1, num_loops ! Year loop do loopyr = start_year, start_year+num_years-1 yr = loopyr if (loopyr < 1870) yr = loopyr + 20 ! if (loopyr < 1870 .and. loopyr > 1859) yr = loopyr + 10 ! Leap year? isleapyear = .false. if (mod(yr,4) .eq. 0) isleapyear = .true. if (mod(yr,100) .eq. 0) isleapyear = .false. if (mod(yr,400) .eq. 0) isleapyear = .true. ! Create file names write(unit=yearstr, fmt='(I4)') yr nc_temp_file = nc_file_path//nc_temp_name//"_"//yearstr//nc_file_tail write (*,'(A)') nc_temp_file nc_precs_file = nc_file_path//nc_precs_name//"_"//yearstr//nc_file_tail write (*,'(A)') nc_precs_file nc_precc_file = nc_file_path//nc_precc_name//"_"//yearstr//nc_file_tail write (*,'(A)') nc_precc_file nc_snowfall_file=nc_file_path//nc_snowfall_name//"_"//yearstr//nc_file_tail write (*,'(A)') nc_snowfall_file nc_st1l_file=nc_file_path//nc_st1l_name//"_"//yearstr//nc_file_tail write (*,'(A)') nc_st1l_file nc_st2l_file=nc_file_path//nc_st2l_name//"_"//yearstr//nc_file_tail write (*,'(A)') nc_st2l_file !nc_st3l_file = nc_file_path // nc_st3l_name // "_" // yearstr !write (*,'(A)') nc_st3l_file !nc_st4l_file = nc_file_path // nc_st4l_name // "_" // yearstr // ".nc" !write (*,'(A)') nc_st4l_file nc_sm1l_file=nc_file_path//nc_sm1l_name//"_"//yearstr//nc_file_tail write (*,'(A)') nc_sm1l_file nc_sm2l_file = nc_file_path//nc_sm2l_name//"_"//yearstr//nc_file_tail write (*,'(A)') nc_sm2l_file nc_sm3l_file = nc_file_path//nc_sm3l_name//"_"//yearstr//nc_file_tail write (*,'(A)') nc_sm3l_file nc_sm4l_file = nc_file_path // nc_sm4l_name //"_"//yearstr//nc_file_tail write (*,'(A)') nc_sm4l_file nc_radi_file = nc_file_path // nc_radi_name//"_"//yearstr//nc_file_tail write (*,'(A)') nc_radi_file if (isleapyear .eq. .true.) then t_max = t_max_leap else t_max = t_max_year endif ! Day/6-hr loop for each year write (*,'(A,I5)') "*II* sparring: t_max =",t_max do t_step=0,t_max-1 t = t_step_full*86400 ! reading climate data from NetCDF file ! PM - changed to 6,12 from 2,6 write (*,'(A,I6,A,I12)') "*II* sparring before read_ifstest: Time step t = ",t_step," - time t=",t ! PM to _pi version call read_ifstest_pi(MAXGRID,1,t_step,t_max,nc_temp_file,nc_temp_name,& f_temp_data,nc_temp_fileid,nc_temp_varid) call read_ifstest_pi(MAXGRID,1,t_step,t_max,nc_precs_file,nc_precs_name,& f_precs_data,nc_precs_fileid,nc_precs_varid) call read_ifstest_pi(MAXGRID,1,t_step,t_max,nc_precc_file,nc_precc_name,& f_precc_data,nc_precc_fileid,nc_precc_varid) ! PM - changed to 6,12 from 2,6 write (*,'(A,I6,A,I12)') "*II* sparring after precip: Time step t = ",t_step," - time t=",t call read_ifstest_pi_depth(MAXGRID,1,t_step,t_max,nc_st1l_file,nc_st1l_name,& f_st1l_data,nc_st1l_fileid,nc_st1l_varid) call read_ifstest_pi_depth(MAXGRID,1,t_step,t_max,nc_st2l_file,nc_st2l_name,& f_st2l_data,nc_st2l_fileid,nc_st2l_varid) ! PM - changed to 6,12 from 2,6 !write (*,'(A,I6,A,I12)') "*II* sparring after st1l and st2l: Time step t = ",t_step," - time t=",t !call read_ifstest_pi(MAXGRID,1,t_step,t_max,nc_st3l_file,nc_st3l_name,& ! f_st3l_data,nc_st3l_fileid,nc_st3l_varid) !call read_ifstest_pi(MAXGRID,1,t_step,t_max,nc_st4l_file,nc_st4l_name,& ! f_st4l_data,nc_st4l_fileid,nc_st4l_varid) call read_ifstest_pi_depth(MAXGRID,1,t_step,t_max,nc_sm1l_file,nc_sm1l_varname,& f_sm1l_data,nc_sm1l_fileid,nc_sm1l_varid) call read_ifstest_pi_depth(MAXGRID,1,t_step,t_max,nc_sm2l_file,nc_sm2l_varname,& f_sm2l_data,nc_sm2l_fileid,nc_sm2l_varid) call read_ifstest_pi_depth(MAXGRID,1,t_step,t_max,nc_sm3l_file,nc_sm3l_varname,& f_sm3l_data,nc_sm3l_fileid,nc_sm3l_varid) call read_ifstest_pi_depth(MAXGRID,1,t_step,t_max,nc_sm4l_file,nc_sm4l_varname,& f_sm4l_data,nc_sm4l_fileid,nc_sm4l_varid) ! PM - changed to 6,12 from 2,6 write (*,'(A,I6,A,I12)') "*II* sparring after SM files: Time step t = ",t_step," - time t=",t call read_ifstest_pi(MAXGRID,1,t_step,t_max,nc_snowfall_file,nc_snowfall_name,& f_snowfall_data,nc_snowfall_fileid,nc_snowfall_varid) call read_ifstest_pi(MAXGRID,1,t_step,t_max,nc_radi_file,nc_radi_name,& f_radi_data,nc_radi_fileid,nc_radi_varid) ! Sum precipitation and snow components ! Replace ST3L and ST4L with values from ST2L do cell = 1, MAXGRID ! Total precip includes snow here f_prec_data(cell,1) = 4 * 1000.0 *(f_precs_data(cell,1) + f_precc_data(cell,1) + f_snowfall_data(cell,1)) ! AMIP data have units m/6h average, so *4 * ! 1000 to get mm/day ! Now * 1000 and / 1000 to go to kg m-2 day-1, ! and / 86400.0 to get to kg m-2 s-1 ! PM_Apr2012 - convert m to kg m-2 s-1 to mimic ECE - we'll convert back in LPJG f_prec_data(cell,1)=f_prec_data(cell,1)/86400.0 f_st3l_data(cell,1) = f_st2l_data(cell,1) f_st4l_data(cell,1) = f_st2l_data(cell,1) ! Snow f_snoc_data(cell,1) = 0.0 f_snod_data(cell,1) = 330.0 ! PM_Apr2012 - change J m-2 to W m-2 to be consistent with ECE. We'll remove the / in LPJG. f_radi_data(cell,1) = 4* f_radi_data(cell,1) / 86400.0 enddo do cell = 1, MAXGRID ! Enforce vegltype as grassland f_vegltype_data(cell,1) = 2; enddo ! PM - changed to 6,12 from 2,6 write (*,'(A,I6,A,I12)') "*II* sparring after ifstest: Time step t = ",t_step," - time t=",t ! ----------------------------------------------------------------------------------- ! *** PUT variables ! ----------------------------------------------------------------------------------- write (*,'(A,I6,A,I12)') "*II* sparringPUT : Time step t = ",t_step," - time t=",t write (*,'(A,I3)') "*II* sparring: calling oasis_put with f_temp_id =",f_temp_id call oasis_put(f_temp_id,t,f_temp_data,ierror) write (*,'(A,I3)') "*II* sparring: oasis_put returned ierror =",ierror write (*,'(A,I3)') "*II* sparring: calling oasis_put with f_radi_id =",f_radi_id call oasis_put(f_radi_id,t,f_radi_data,ierror) write (*,'(A,I3)') "*II* sparring: oasis_put returned ierror =",ierror write (*,'(A,I3)') "*II* sparring: calling oasis_put with f_sm1l_id =",f_sm1l_id call oasis_put(f_sm1l_id,t,f_sm1l_data,ierror) write (*,'(A,I3)') "*II* sparring: oasis_put returned ierror =",ierror write (*,'(A,I3)') "*II* sparring: calling oasis_put with f_sm2l_id =",f_sm2l_id call oasis_put(f_sm2l_id,t,f_sm2l_data,ierror) write (*,'(A,I3)') "*II* sparring: oasis_put returned ierror =",ierror write (*,'(A,I3)') "*II* sparring: calling oasis_put with f_sm3l_id =",f_sm3l_id call oasis_put(f_sm3l_id,t,f_sm3l_data,ierror) write (*,'(A,I3)') "*II* sparring: oasis_put returned ierror =",ierror write (*,'(A,I3)') "*II* sparring: calling oasis_put with f_sm4l_id =",f_sm4l_id call oasis_put(f_sm4l_id,t,f_sm4l_data,ierror) write (*,'(A,I3)') "*II* sparring: oasis_put returned ierror =",ierror write (*,'(A,I3)') "*II* sparring: calling oasis_put with f_st1l_id =",f_st1l_id call oasis_put(f_st1l_id,t,f_st1l_data,ierror) write (*,'(A,I3)') "*II* sparring: oasis_put returned ierror =",ierror write (*,'(A,I3)') "*II* sparring: calling oasis_put with f_st2l_id =",f_st2l_id call oasis_put(f_st2l_id,t,f_st2l_data,ierror) write (*,'(A,I3)') "*II* sparring: oasis_put returned ierror =",ierror write (*,'(A,I3)') "*II* sparring: calling oasis_put with f_st3l_id =",f_st3l_id call oasis_put(f_st3l_id,t,f_st3l_data,ierror) write (*,'(A,I3)') "*II* sparring: oasis_put returned ierror =",ierror write (*,'(A,I3)') "*II* sparring: calling oasis_put with f_st4l_id =",f_st4l_id call oasis_put(f_st4l_id,t,f_st4l_data,ierror) write (*,'(A,I3)') "*II* sparring: oasis_put returned ierror =",ierror ! Snow mass/unit surface (kg/m2) write (*,'(A,I3)') "*II* sparring: calling oasis_put with f_snoc_id =",f_snoc_id call oasis_put(f_snoc_id,t,f_snoc_data,ierror) write (*,'(A,I3)') "*II* sparring: oasis_put returned ierror =",ierror ! Snow density (kg/m3) write (*,'(A,I3)') "*II* sparring: calling oasis_put with f_snod_id =",f_snod_id call oasis_put(f_snod_id,t,f_snod_data,ierror) write (*,'(A,I3)') "*II* sparring: oasis_put returned ierror =",ierror ! Veg low write (*,'(A,I3)') "*II* sparring: calling oasis_put with f_vegl_id =",f_vegl_id call oasis_put(f_vegl_id,t,f_vegl_data,ierror) write (*,'(A,I3)') "*II* sparring: oasis_put returned ierror =",ierror ! Veg type low write (*,'(A,I3)') "*II* sparring: calling oasis_put with f_vegltype_id =",f_vegltype_id call oasis_put(f_vegltype_id,t,f_vegltype_data,ierror) write (*,'(A,I3)') "*II* sparring: oasis_put returned ierror =",ierror ! Veg high write (*,'(A,I3)') "*II* sparring: calling oasis_put with f_vegh_id =",f_vegh_id call oasis_put(f_vegh_id,t,f_vegh_data,ierror) write (*,'(A,I3)') "*II* sparring: oasis_put returned ierror =",ierror write (*,'(A,I3)') "*II* sparring: calling oasis_put with f_prec_id =",f_prec_id call oasis_put(f_prec_id,t,f_prec_data,ierror) write (*,'(A,I3)') "*II* sparring: oasis_put returned ierror =",ierror ! ----------------------------------------------------------------------------------- ! *** GET variables ! ----------------------------------------------------------------------------------- write (*,'(A,I6,A,I12)') "*II* sparringGET : Time step t = ",t_step," - time t=",t write (*,'(A,I3)') "*II* sparring: calling oasis_get with f_llai_id =",f_llai_id call oasis_get(f_llai_id,t,f_llai_data,ierror) write (*,'(A,I3)') "*II* sparring: oasis_get returned ierror =",ierror write (*,'(A,I3)') "*II* sparring: calling oasis_get with f_hlai_id =",f_hlai_id call oasis_get(f_hlai_id,t,f_hlai_data,ierror) write (*,'(A,I3)') "*II* sparring: oasis_get returned ierror =",ierror write (*,'(A,I3)') "*II* sparring: calling oasis_get with f_typh_id =",f_typh_id call oasis_get(f_typh_id,t,f_typh_data,ierror) write (*,'(A,I3)') "*II* sparring: oasis_get returned ierror =",ierror write (*,'(A,I3)') "*II* sparring: calling oasis_get with f_vggh_id =",f_vggh_id call oasis_get(f_vggh_id,t,f_vggh_data,ierror) write (*,'(A,I3)') "*II* sparring: oasis_get returned ierror =",ierror ! ecev3 write (*,'(A,I3)') "*II* sparring: calling oasis_get with f_typl_id =",f_typl_id call oasis_get(f_typl_id,t,f_typl_data,ierror) write (*,'(A,I3)') "*II* sparring: oasis_get returned ierror =",ierror ! ecev3 write (*,'(A,I3)') "*II* sparring: calling oasis_get with f_vggl_id =",f_vggl_id call oasis_get(f_vggl_id,t,f_vggl_data,ierror) write (*,'(A,I3)') "*II* sparring: oasis_get returned ierror =",ierror ! checking whether the incoming LAI data are ok: write(*,*)"*II* sparring: LLAI obtained from OASIS for test gridcell (g=8000) is ",f_llai_data(8000,1) write(*,*)"*II* sparring: HLAI obtained from OASIS for test gridcell (g=8000) is ",f_hlai_data(8000,1) write(*,*)"*II* sparring: TYPH obtained from OASIS for test gridcell (g=8000) is ",f_typh_data(8000,1) write(*,*)"*II* sparring: VEGH obtained from OASIS for test gridcell (g=8000) is ",f_vggh_data(8000,1) write(*,*)"*II* sparring: TYPL obtained from OASIS for test gridcell (g=8000) is ",f_typl_data(8000,1) write(*,*)"*II* sparring: VEGL obtained from OASIS for test gridcell (g=8000) is ",f_vggl_data(8000,1) write (*,'(A,I6,A,I12)') "*II* sparring DIAGS LAI check: Time step t = ",t_step," - time t=",t ! Check for min/max LAI minhlai=0 maxhlai=0 minllai=0 maxllai=0 meanhlai=0 meanllai=0 meanvegl=0 meanvegh=0 maxcvh=-1 maxcvl=-1 maxtvh=-1 maxtvl=-1 cover=0.0 do cell = 1, MAXGRID if (f_llai_data(cell,1)>maxllai)then maxllai=f_llai_data(cell,1) endif meanllai = meanllai + f_llai_data(cell,1)/MAXGRID if (f_llai_data(cell,1)maxhlai)then maxhlai=f_hlai_data(cell,1) endif if (f_hlai_data(cell,1)maxcvl)then maxcvl=f_vggl_data(cell,1) endif if (f_vggh_data(cell,1)>maxcvh)then maxcvh=f_vggh_data(cell,1) endif meanvegl = meanvegl + f_vegl_data(cell,1)/MAXGRID meanvegh = meanvegh + f_vegh_data(cell,1)/MAXGRID if (f_typl_data(cell,1)>maxtvl)then maxtvl=f_typl_data(cell,1) endif if (f_typh_data(cell,1)>maxtvh)then maxtvh=f_typh_data(cell,1) endif cover=f_vggl_data(cell,1)+f_vggh_data(cell,1) if (cover>1.001) then write(*,*)"*II* sparring: COVER ERROR!",cover endif enddo write(*,*)"*II* sparring DIAGS: MAX HLAI :",maxhlai write(*,*)"*II* sparring: MIN HLAI :",minhlai write(*,*)"*II* sparring: MAX LLAI :",maxllai write(*,*)"*II* sparring: MIN LLAI :",minllai write(*,*)"*II* sparring: MEAN LLAI :",meanllai write(*,*)"*II* sparring: MEAN HLAI :",meanhlai write(*,*)"*II* sparring: MEAN VEGL :",meanvegl write(*,*)"*II* sparring: MEAN VEGH :",meanvegh write(*,*)"*II* sparring: MAX VEGL :",maxcvl write(*,*)"*II* sparring: MAX VEGH :",maxcvh write(*,*)"*II* sparring: MAX TYPL :",maxtvl write(*,*)"*II* sparring: MAX TYPH :",maxtvh write (*,*) "*II* sparring: Finished dayloop!",t_step, t_max ! Increase the full simulation counter t_step_full = t_step_full + 1 ! End of dayloop enddo write (*,*) "*II* sparring: End of year!",loopyr ! End of year loop enddo write (*,'(A)') "*II* sparring: Finished time loop!" ! End of repetition loop enddo write (*,'(A)') "*II* sparring: oasis_terminate..." call oasis_terminate(ierror) write (*,'(A,I3)') "*II* sparring: oasis_terminate returned ierror =",ierror end program SPARRING ! PM - for pre-Industrial nc files. ! now reads new_ files ! --------------------------------------------------------------------------------------------------------------- subroutine read_ifstest_pi(NX,NY,t_step,t_max,filename,parname, & data_in,fileid,varid) ! NetCDF use netcdf implicit none integer::NX,NY,t_step,t_max character(len=*)::filename character(len=*)::parname double precision::data_in(NX,NY) integer::fileid integer::varid ! for older dpic_ files: !integer::start(2),count(2) !start=(/1,t_step+1/) !count=(/NX,NY/) ! for new_ files: integer::start(3),count(3) start=(/1,1,t_step+1/) count=(/NX,NY,1/) ! read file properties in first occasion if(t_step==0)then call check(nf90_open(filename,NF90_NOWRITE,fileid)) call check(nf90_inq_varid(fileid,parname,varid)) endif ! read field for given timestep call check(nf90_get_var(fileid,varid,data_in,start=start,count=count)) ! close NetCDF file if(t_step==t_max-1)then call check(nf90_close(fileid)) endif end ! PM - for pre-Industrial nc files with DEPTH ! --------------------------------------------------------------------------------------------------------------- subroutine read_ifstest_pi_depth(NX,NY,t_step,t_max,filename,parname, & data_in,fileid,varid) ! NetCDF use netcdf implicit none integer::NX,NY,t_step,t_max character(len=*)::filename character(len=*)::parname double precision::data_in(NX,NY) integer::fileid integer::varid ! double precision::data_in(NX,NY) ! integer::fileid ! integer::varid !integer::start(3),count(3) !start=(/1,1,t_step+1/) !count=(/NX,NY,1/) ! new_ files integer::start(4),count(4) start=(/1,1,1,t_step+1/) count=(/NX,NY,1,1/) ! read file properties in first occasion if(t_step==0)then call check(nf90_open(filename,NF90_NOWRITE,fileid)) call check(nf90_inq_varid(fileid,parname,varid)) endif ! read field for given timestep call check(nf90_get_var(fileid,varid,data_in,start=start,count=count)) ! close NetCDF file if(t_step==t_max-1)then call check(nf90_close(fileid)) endif end ! PM - for fixed nc files. ! --------------------------------------------------------------------------------------------------------------- subroutine read_fixedfiels(NX,NY,filename,parname,data_in,fileid,varid) ! NetCDF use netcdf implicit none integer::NX,NY character(len=*)::filename character(len=*)::parname double precision::data_in(NX,NY) integer::fileid integer::varid integer::start(2),count(2) start=(/1,1/) count=(/NX,NY/) ! read file properties in first occasion call check(nf90_open(filename,NF90_NOWRITE,fileid)) call check(nf90_inq_varid(fileid,parname,varid)) ! read field call check(nf90_get_var(fileid,varid,data_in,start=start,count=count)) ! close NetCDF file call check(nf90_close(fileid)) end ! --------------------------------------------------------------------------------------------------------------- subroutine read_ifstest(NX,NY,t_step,t_max,filename,parname, & data_in,fileid,varid) ! NetCDF use netcdf implicit none integer::NX,NY,t_step,t_max character(len=*)::filename character(len=*)::parname double precision::data_in(NX,NY) integer::fileid integer::varid integer::start(3),count(3) start=(/1,1,t_step+1/) count=(/NX,NY,1/) ! read file properties in first occasion if(t_step==0)then call check(nf90_open(filename,NF90_NOWRITE,fileid)) call check(nf90_inq_varid(fileid,parname,varid)) endif ! read field for given timestep call check(nf90_get_var(fileid,varid,data_in,start=start,count=count)) ! close NetCDF file if(t_step==t_max-1)then call check(nf90_close(fileid)) endif end ! --------------------------------------------------------------------------------------------------------------- subroutine check(status) use netcdf integer, intent ( in) :: status if(status /= nf90_noerr) then print *, trim(nf90_strerror(status)) stop "Stopped" end if end subroutine check ! ---------------------------------------------------------------------------------------------------------------