Use NetCDF Real, Parameter :: Miss=-99. Integer, Parameter :: Nmon=12, Nreg=11 Integer, Parameter :: NPeriod=3, NSce=2, Nmodels=7 ! 7 = 6 models + Ens !--------------------------------------------------------------------------------- Real, allocatable :: xlon(:), xlat(:), xtime(:) Integer :: Nlon, Nlat, Ntime, Nyear Logical :: Status Real, allocatable :: Mask(:,:), Regs(:,:) ! (Nlon,Nlat) ! Input Data (Reading from files): Real, allocatable :: Ob_T2m(:,:,:), Ob_Pr(:,:,:) ! (Nlon,Nlat,Ntime) - Aug 2020 Real, allocatable :: BL_T2m(:,:,:,:), BL_Pr(:,:,:,:) ! (Nlon,Nlat,Ntime,Nmodels) Real, allocatable :: MF_T2m(:,:,:,:,:), MF_Pr(:,:,:,:,:) ! (Nlon,Nlat,Ntime,Nmodels,NSce) Real, allocatable :: FF_T2m(:,:,:,:,:), FF_Pr(:,:,:,:,:) ! (Nlon,Nlat,Ntime,Nmodels,NSce) Real, allocatable :: T2m(:,:,:), Pr(:,:,:) !--------------------------------------------------------------------------------- ! Characteristics: ! Truong Trung binh thang Real, allocatable :: Ob_T2m_Mean(:,:,:), Ob_Pr_Mean(:,:,:) ! (Nlon,Nlat,Ntime) - Aug 2020 Real, allocatable :: BL_T2m_Mean(:,:,:,:), BL_Pr_Mean(:,:,:,:) !(Nlon,Nlat,Nmon,Nmodels) Real, allocatable :: MF_T2m_Mean(:,:,:,:,:), MF_Pr_Mean(:,:,:,:,:) !(Nlon,Nlat,Nmon,Nmodels,NSce) Real, allocatable :: FF_T2m_Mean(:,:,:,:,:), FF_Pr_Mean(:,:,:,:,:) !(Nlon,Nlat,Nmon,Nmodels,NSce) ! Truong Trung binh nam Real, allocatable :: Ob_T2m_Ann(:,:), Ob_Pr_Ann(:,:) !(Nlon,Nlat) - Aug 2020 Real, allocatable :: BL_T2m_Ann(:,:,:), BL_Pr_Ann(:,:,:) !(Nlon,Nlat,Nmodels) Real, allocatable :: MF_T2m_Ann(:,:,:,:), MF_Pr_Ann(:,:,:,:) !(Nlon,Nlat,Nmodels,NSce) Real, allocatable :: FF_T2m_Ann(:,:,:,:), FF_Pr_Ann(:,:,:,:) !(Nlon,Nlat,Nmodels,NSce) ! Trung binh thang cho tung vung KH Real, allocatable :: Ob_T2m_RMon(:,:), Ob_Pr_RMon(:,:) !(Nlon,Nlat) - Aug 2020 Real, allocatable :: BL_T2m_RMon(:,:,:), BL_Pr_RMon(:,:,:) !(Nreg,Nmon,Nmodels) Real, allocatable :: MF_T2m_Rmon(:,:,:,:), MF_Pr_RMon(:,:,:,:) !(Nreg,Nmon,Nmodels,NSce) Real, allocatable :: FF_T2m_Rmon(:,:,:,:), FF_Pr_RMon(:,:,:,:) !(Nreg,Nmon,Nmodels,NSce) ! Trung binh nam cho tung vung KH Real, allocatable :: Ob_T2m_RAnn(:), Ob_Pr_RAnn(:) !(Nreg) - Aug 2020 Real, allocatable :: BL_T2m_RAnn(:,:), BL_Pr_RAnn(:,:) !(Nreg,Nmodels) Real, allocatable :: MF_T2m_RAnn(:,:,:), MF_Pr_RAnn(:,:,:) !(Nreg,Nmodels,NSce) Real, allocatable :: FF_T2m_RAnn(:,:,:), FF_Pr_RAnn(:,:,:) !(Nreg,Nmodels,NSce) ! Changes of Annual mean: Real, allocatable :: D_BL_T2m_Ann(:,:,:), D_BL_Pr_Ann(:,:,:) !(Nlon,Nlat,Nmodels): BL-Obs - Aug 2020 Real, allocatable :: D_MF_T2m_Ann(:,:,:,:), D_MF_Pr_Ann(:,:,:,:) !(Nlon,Nlat,Nmodels,NSce) Real, allocatable :: D_FF_T2m_Ann(:,:,:,:), D_FF_Pr_Ann(:,:,:,:) !(Nlon,Nlat,Nmodels,NSce) ! Changes of Monthly mean: Real, allocatable :: D_BL_T2m_Mon(:,:,:,:), D_BL_Pr_Mon(:,:,:,:) !(Nlon,Nlat,Nmon,Nmodels): BL-Obs - Aug 2020 Real, allocatable :: D_MF_T2m_Mon(:,:,:,:,:), D_MF_Pr_Mon(:,:,:,:,:) !(Nlon,Nlat,Nmon,Nmodels,NSce) Real, allocatable :: D_FF_T2m_Mon(:,:,:,:,:), D_FF_Pr_Mon(:,:,:,:,:) !(Nlon,Nlat,Nmon,Nmodels,NSce) Real :: DelT(Nmon+1,Nreg,NPeriod-1,NSce,Nmodels) Real :: DelP(Nmon+1,Nreg,NPeriod-1,NSce,Nmodels) Real :: MeanT(Nmon+1,Nreg,NPeriod,NSce,Nmodels) Real :: MeanP(Nmon+1,Nreg,NPeriod,NSce,Nmodels) !--------------------------------------------------------------------------------- Character (Len=100) :: Fname, LandMaskFile, RegionFile Character (Len=100) :: Data_Dir, Model_nm(Nmodels), Model_Obs, Out_Dir Integer :: YSta(NPeriod), YEnd(NPeriod), itime1(NPeriod), itime2(NPeriod) Character (Len=5) :: Sce_nm(NSce), Period(NPeriod) !--------------------------------------------------------------------------------- Namelist /Infor/ Data_Dir, Period, YSta, YEnd, Nyear, Model_Obs, Model_nm, Sce_nm, LandMaskFile, RegionFile, Out_Dir Integer :: ilon,ilat,it,imon,iper,imod,isce Integer :: Nt0, it1, it2 ! Timeseries !--------------------------------------------------------------------------------- Open (1, file="Namelist_R_T_2020_Aug.txt", Status="Old") Read (1, Infor) Close (1) Call Read_Land_Mask Ntime = Nyear*Nmon itime1(1) = 1; itime2(1) = Ntime itime1(2) = (YSta(2) -YEnd(1)-1)*Nmon+1 ! itime1=(2036-2005-1)*12+1=361 itime2(2) = itime1(2)+Ntime-1 itime1(3) = (YSta(3) -YEnd(1)-1)*Nmon+1 ! itime2=(2070-2005-1)*12+1=769 itime2(3) = itime1(3)+Ntime-1 Allocate (Ob_T2m(Nlon,Nlat,Ntime), Ob_Pr(Nlon,Nlat,Ntime)) Allocate (BL_T2m(Nlon,Nlat,Ntime,Nmodels), BL_Pr(Nlon,Nlat,Ntime,Nmodels)) Allocate (MF_T2m(Nlon,Nlat,Ntime,Nmodels,NSce), MF_Pr(Nlon,Nlat,Ntime,Nmodels,NSce)) Allocate (FF_T2m(Nlon,Nlat,Ntime,Nmodels,NSce), FF_Pr(Nlon,Nlat,Ntime,Nmodels,NSce)) !Allocate (T2m(Nlon,Nlat,Ntime)); Allocate (Pr (Nlon,Nlat,Ntime)) Do iper=1, NPeriod if (iper==1) then ! BL: Obs +6 models + 1 Ens Fname=trim(Data_Dir)//trim(Period(iper))//"_"//trim(Model_Obs)//".nc" print*,trim(Fname); Call Read_T_R (Fname) Do ilon=1,Nlon Do ilat=1,Nlat Ob_T2m(ilon,ilat,1:Ntime)=T2m(ilon,ilat,1:Ntime) Ob_Pr (ilon,ilat,1:Ntime)=Pr (ilon,ilat,1:Ntime) Enddo Enddo Deallocate (T2m, Pr) Do imod=1, Nmodels Fname=trim(Data_Dir)//trim(Period(iper))//"_"//trim(Model_nm(imod))//".nc" Call Read_T_R (Fname) Do ilon=1,Nlon Do ilat=1,Nlat BL_T2m(ilon,ilat,1:Ntime,imod)=T2m(ilon,ilat,1:Ntime) BL_Pr (ilon,ilat,1:Ntime,imod)=Pr (ilon,ilat,1:Ntime) Enddo Enddo Deallocate (T2m, Pr) Enddo else Do isce=1,NSce Do imod=1, Nmodels Fname=trim(Data_Dir)//"SCE_"//trim(Sce_nm(isce))//"_"//trim(Model_nm(imod))//".nc" Call Read_T_R (Fname) it1=itime1(iper); it2=itime1(iper) if (iper==2) then ! MF Do ilon=1,Nlon Do ilat=1,Nlat MF_T2m(ilon,ilat,1:Ntime,imod,isce)=T2m(ilon,ilat,it1:it2) MF_Pr (ilon,ilat,1:Ntime,imod,isce)=Pr (ilon,ilat,it1:it2) Enddo Enddo !print*,MF_T2m(ilon,ilat,20:30,imod,isce) !MF_T2m(:,:,:,imod,isce)=T2m(:,:,:); MF_Pr(:,:,:,imod,isce)=Pr(:,:,:) else ! FF Do ilon=1,Nlon Do ilat=1,Nlat FF_T2m(ilon,ilat,1:Ntime,imod,isce)=T2m(ilon,ilat,it1:it2) FF_Pr (ilon,ilat,1:Ntime,imod,isce)=Pr (ilon,ilat,it1:it2) Enddo Enddo !FF_T2m(:,:,:,imod,isce)=T2m(:,:,:); FF_Pr(:,:,:,imod,isce)=Pr(:,:,:) endif Deallocate (T2m, Pr) Enddo Enddo endif Enddo print*,"Reading data ... Done !" !Do ilon=1,Nlon ! Do ilat=1, Nlat ! if (Mask(ilon,ilat)>0) then ! print*,ilon,ilat ! print*,Ob_T2m(ilon,ilat,1:Ntime) ! endif ! Enddo !Enddo !pause print*,"Calling Chars..." Call Chars print*,"Calling Chars...Done !" print*,"Calling Calc_Regional_Mean..." Call Calc_Regional_Mean !Call For_Can_Tho CONTAINS !--------------------------------------------------------------------------------- Subroutine Calc_Regional_Mean Integer :: ilon,ilat,it,imon,imon1,iper,imod,isce, ireg, icount Real :: X(Nlon,Nlat) Character (Len=100) :: Fname, Sce_nm, t_unit, var_ln, var_unit, lon_nm, lat_nm, z_nm, t_nm Character (Len=100) :: Fname1, Fname2 Real, Allocatable :: X2D1(:,:),X2D2(:,:), X3D1(:,:,:), X3D2(:,:,:), Xmod(:), Xmon(:), Xreg(:), XX(:) Real, Allocatable :: Y3D1(:,:,:), Y3D2(:,:,:) ! Dung de sap xep lai cac thang 11,12,1,2,..,10 Allocate (Ob_T2m_RMon(Nreg,Nmon), Ob_Pr_RMon(Nreg,Nmon)) ! Aug 2020 Allocate (BL_T2m_RMon(Nreg,Nmon,Nmodels), BL_Pr_RMon(Nreg,Nmon,Nmodels)) Allocate (MF_T2m_Rmon(Nreg,Nmon,Nmodels,NSce), MF_Pr_RMon(Nreg,Nmon,Nmodels,NSce)) Allocate (FF_T2m_Rmon(Nreg,Nmon,Nmodels,NSce), FF_Pr_RMon(Nreg,Nmon,Nmodels,NSce)) Allocate (Ob_T2m_RAnn(Nreg), Ob_Pr_RAnn(Nreg)) ! Aug 2020 Allocate (BL_T2m_RAnn(Nreg,Nmodels), BL_Pr_RAnn(Nreg,Nmodels)) Allocate (MF_T2m_RAnn(Nreg,Nmodels,NSce), MF_Pr_RAnn(Nreg,Nmodels,NSce)) Allocate (FF_T2m_RAnn(Nreg,Nmodels,NSce), FF_Pr_RAnn(Nreg,Nmodels,NSce)) Ob_T2m_RMon=Miss; Ob_Pr_RMon=Miss; Ob_T2m_RAnn=Miss; Ob_Pr_RAnn=Miss ! Aug 2020 BL_T2m_RMon=Miss; BL_Pr_RMon=Miss; MF_T2m_Rmon=Miss; MF_Pr_RMon=Miss; FF_T2m_Rmon=Miss; FF_Pr_RMon=Miss ! Obs Monthly Regional mean: ! Aug 2020 Do imon=1, Nmon Do ireg=1, Nreg X(1:Nlon,1:Nlat)=Ob_T2m_Mean(1:Nlon,1:Nlat,imon); Ob_T2m_RMon(ireg,imon)=Regional_Mean(X,ireg) X(1:Nlon,1:Nlat)=Ob_Pr_Mean(1:Nlon,1:Nlat,imon); Ob_Pr_RMon (ireg,imon)=Regional_Mean(X,ireg) Enddo Enddo ! BL: Monthly Regional mean Do imod=1,Nmodels Do imon=1, Nmon Do ireg=1, Nreg X(1:Nlon,1:Nlat)=BL_T2m_Mean(1:Nlon,1:Nlat,imon,imod); BL_T2m_RMon(ireg,imon,imod)=Regional_Mean(X,ireg) !+1) X(1:Nlon,1:Nlat)=BL_Pr_Mean (1:Nlon,1:Nlat,imon,imod); BL_Pr_RMon (ireg,imon,imod)=Regional_Mean(X,ireg) !+1) Enddo Enddo Enddo ! SCE: Monthly Regional mean Do isce=1, NSce Do imod=1,Nmodels Do imon=1, Nmon Do ireg=1, Nreg X(1:Nlon,1:Nlat)=MF_T2m_Mean(1:Nlon,1:Nlat,imon,imod,isce); MF_T2m_Rmon(ireg,imon,imod,isce)=Regional_Mean(X,ireg) !+1) X(1:Nlon,1:Nlat)=MF_Pr_Mean (1:Nlon,1:Nlat,imon,imod,isce); MF_Pr_RMon (ireg,imon,imod,isce)=Regional_Mean(X,ireg) !+1) X(1:Nlon,1:Nlat)=FF_T2m_Mean(1:Nlon,1:Nlat,imon,imod,isce); FF_T2m_Rmon(ireg,imon,imod,isce)=Regional_Mean(X,ireg) !+1) X(1:Nlon,1:Nlat)=FF_Pr_Mean (1:Nlon,1:Nlat,imon,imod,isce); FF_Pr_RMon (ireg,imon,imod,isce)=Regional_Mean(X,ireg) !+1) Enddo Enddo Enddo Enddo ! Obs Annual Regional mean ! Aug 2020 print*," Obs Annual Regional mean" Do ireg=1, Nreg X(1:Nlon,1:Nlat)=Ob_T2m_Ann(1:Nlon,1:Nlat); Ob_T2m_RAnn(ireg)=Regional_Mean(X,ireg) X(1:Nlon,1:Nlat)=Ob_Pr_Ann(1:Nlon,1:Nlat); Ob_Pr_RAnn(ireg)=Regional_Mean(X,ireg) Enddo ! BL: Annual Regional mean Do imod=1,Nmodels Do ireg=1, Nreg X(1:Nlon,1:Nlat)=BL_T2m_Ann(1:Nlon,1:Nlat,imod); BL_T2m_RAnn(ireg,imod)=Regional_Mean(X,ireg) !+1) X(1:Nlon,1:Nlat)=BL_Pr_Ann (1:Nlon,1:Nlat,imod); BL_Pr_RAnn (ireg,imod)=Regional_Mean(X,ireg) !+1) Enddo Enddo ! SCE: Annual Regional mean Do isce=1, NSce Do imod=1,Nmodels Do ireg=1, Nreg X(1:Nlon,1:Nlat)=MF_T2m_Ann(1:Nlon,1:Nlat,imod,isce); MF_T2m_RAnn(ireg,imod,isce)=Regional_Mean(X,ireg) !+1) X(1:Nlon,1:Nlat)=MF_Pr_Ann (1:Nlon,1:Nlat,imod,isce); MF_Pr_RAnn (ireg,imod,isce)=Regional_Mean(X,ireg) !+1) X(1:Nlon,1:Nlat)=FF_T2m_Ann(1:Nlon,1:Nlat,imod,isce); FF_T2m_RAnn(ireg,imod,isce)=Regional_Mean(X,ireg) !+1) X(1:Nlon,1:Nlat)=FF_Pr_Ann (1:Nlon,1:Nlat,imod,isce); FF_Pr_RAnn (ireg,imod,isce)=Regional_Mean(X,ireg) !+1) Enddo Enddo Enddo ! Write Obs regional data to file: Allocate (X3D1(Nreg,Nmodels,Nmon),X3D2(Nreg,Nmodels,Nmon), Xmod(Nmodels), Xmon(Nmon), Xreg(Nreg)) Allocate (Y3D1(Nreg,Nmodels,Nmon),Y3D2(Nreg,Nmodels,Nmon)) Fname2=trim(Out_Dir)//"Obs_BL_Reg_mean.nc" Do ireg=1,Nreg; Xreg(ireg)=ireg; Enddo Do imon=1,Nmon; Xmon(imon)=imon; Enddo Do imod=1,Nmodels; Xmod(imod)=imod; Enddo t_unit="months since 1986-01-01" var_ln="Monthly regional mean of Model Temp for period 1986-2005"; lon_nm="reg"; lat_nm="mon"; t_nm="model" Call wr_ncdf_3d (Fname2,Xreg,Xmon,Xmod,Nreg,Nmon,Nmodels,BL_T2m_RMon,"t_bl",t_unit, var_ln, "C", lon_nm, lat_nm, t_nm, Miss) var_ln="Monthly regional mean of Model Prec for period 1986-2005"; Call wr_ncdf_3d_apd (Fname2,BL_Pr_RMon,"p_bl",Nreg,Nmon,Nmodels, var_ln, "mm/mon", lon_nm, lat_nm, t_nm, Miss) var_ln="Monthly regional mean of Obs Temp for period 1986-2005"; lon_nm="reg"; lat_nm="mon" Call wr_ncdf_2d_apd (Fname2,Ob_T2m_RMon,"tobs",Nreg,Nmon, var_ln,"C", lon_nm, lat_nm, Miss) var_ln="Monthly regional mean of Obs Preci for period 1986-2005" Call wr_ncdf_2d_apd (Fname2,Ob_Pr_RMon,"pobs",Nreg,Nmon, var_ln, "mm/mon", lon_nm, lat_nm, Miss) var_ln="Annual regional mean of Temp for period 1986-2005" Call wr_ncdf_1d_apd (Fname2,Ob_T2m_RAnn,"t_ann", Nreg, var_ln, "C", lon_nm, Miss) var_ln="Annual regional mean of Preci for period 1986-2005" Call wr_ncdf_1d_apd (Fname2,Ob_Pr_RAnn,"p_ann", Nreg, var_ln, "mm/mon", lon_nm, Miss) ! Write to file: !Allocate (X3D1(Nreg,Nmodels,Nmon),X3D2(Nreg,Nmodels,Nmon), Xmod(Nmodels), Xmon(Nmon), Xreg(Nreg)) !Allocate (Y3D1(Nreg,Nmodels,Nmon),Y3D2(Nreg,Nmodels,Nmon)) Do imon=1,Nmon; Xmon(imon)=imon; Enddo Do imod=1,Nmodels; Xmod(imod)=imod; Enddo Do ireg=1,Nreg; Xreg(ireg)=ireg; Enddo t_unit="months since 1986-01-01"; var_unit="C"; lon_nm="reg"; lat_nm="model"; z_nm="scenario"; t_nm="mon" var_ln="Monthly regional mean of temperature for the period 1986-2005 projeted by different models" Fname=trim(Out_Dir)//"Monthly_Reg_Mean.nc" Fname1=trim(Out_Dir)//"Monthly_Reg_Mean_Reshape.nc" ! File co thu tu THANG la 11,12,1,2,...,10 Do ireg=1,Nreg Do imod=1,Nmodels X3D1(ireg,imod,1:Nmon)=BL_T2m_RMon(ireg,1:Nmon,imod) X3D2(ireg,imod,1:Nmon)=BL_Pr_RMon (ireg,1:Nmon,imod) Y3D1(ireg,imod,1:2) =BL_T2m_RMon(ireg,11:12,imod); Y3D1(ireg,imod,3:12) = BL_T2m_RMon(ireg,1:10,imod) Y3D2(ireg,imod,1:2) =BL_Pr_RMon (ireg,11:12,imod); Y3D2(ireg,imod,3:12) = BL_Pr_RMon (ireg,1:10,imod) Enddo Enddo Call wr_ncdf_3d (Fname,Xreg,Xmod,Xmon,Nreg,Nmodels,Nmon,X3D1,"t_bl",t_unit, var_ln, var_unit, lon_nm, lat_nm, t_nm, Miss) Call wr_ncdf_3d (Fname1,Xreg,Xmod,Xmon,Nreg,Nmodels,Nmon,Y3D1,"t_bl",t_unit, var_ln, var_unit, lon_nm, lat_nm, t_nm, Miss) ! Fname1 var_ln="Monthly regional mean of rainfall amount for the period 1986-2005" var_unit="mm/mon" Call wr_ncdf_3d_apd (Fname,X3D2,"p_bl",Nreg,Nmodels,Nmon, var_ln, var_unit, lon_nm, lat_nm, t_nm, Miss) Call wr_ncdf_3d_apd (Fname1,Y3D2,"p_bl",Nreg,Nmodels,Nmon, var_ln, var_unit, lon_nm, lat_nm, t_nm, Miss) ! Fname1 Do isce=1,NSce if (isce==1) then; Sce_nm="45"; else; Sce_nm="85"; endif Do ireg=1,Nreg Do imod=1,Nmodels X3D1(ireg,imod,1:Nmon)=MF_T2m_Rmon(ireg,1:Nmon,imod,isce) X3D2(ireg,imod,1:Nmon)=MF_Pr_RMon (ireg,1:Nmon,imod,isce) Y3D1(ireg,imod,1:2) =MF_T2m_Rmon(ireg,11:12,imod,isce); Y3D1(ireg,imod,3:12) = MF_T2m_Rmon(ireg,1:10,imod,isce) Y3D2(ireg,imod,1:2) =MF_Pr_RMon (ireg,11:12,imod,isce); Y3D2(ireg,imod,3:12) = MF_Pr_RMon (ireg,1:10,imod,isce) Enddo Enddo var_ln="Monthly regional mean of temperature for the period 2046-2065, RCP"//trim(Sce_nm); var_unit="C" Call wr_ncdf_3d_apd (Fname,X3D1,"t_mf"//trim(Sce_nm),Nreg,Nmodels,Nmon, var_ln, var_unit, lon_nm, lat_nm, t_nm, Miss) Call wr_ncdf_3d_apd (Fname1,Y3D1,"t_mf"//trim(Sce_nm),Nreg,Nmodels,Nmon, var_ln, var_unit, lon_nm, lat_nm, t_nm, Miss) ! Fname1 var_ln="Monthly regional mean of rainfall amount for the period 2046-2065, RCP"//trim(Sce_nm); var_unit="mm/mon" Call wr_ncdf_3d_apd (Fname,X3D2,"p_mf"//trim(Sce_nm),Nreg,Nmodels,Nmon, var_ln, var_unit, lon_nm, lat_nm, t_nm, Miss) Call wr_ncdf_3d_apd (Fname1,Y3D2,"p_mf"//trim(Sce_nm),Nreg,Nmodels,Nmon, var_ln, var_unit, lon_nm, lat_nm, t_nm, Miss) ! Fname1 Do ireg=1,Nreg Do imod=1,Nmodels X3D1(ireg,imod,1:Nmon)=FF_T2m_Rmon(ireg,1:Nmon,imod,isce) X3D2(ireg,imod,1:Nmon)=FF_Pr_RMon (ireg,1:Nmon,imod,isce) Y3D1(ireg,imod,1:2) =FF_T2m_Rmon(ireg,11:12,imod,isce); Y3D1(ireg,imod,3:12) = FF_T2m_Rmon(ireg,1:10,imod,isce) Y3D2(ireg,imod,1:2) =FF_Pr_RMon (ireg,11:12,imod,isce); Y3D2(ireg,imod,3:12) = FF_Pr_RMon (ireg,1:10,imod,isce) Enddo Enddo var_ln="Monthly regional mean of temperature for the period 2080-2099, RCP"//trim(Sce_nm); var_unit="C" Call wr_ncdf_3d_apd (Fname,X3D1,"t_ff"//trim(Sce_nm),Nreg,Nmodels,Nmon, var_ln, var_unit, lon_nm, lat_nm, t_nm, Miss) Call wr_ncdf_3d_apd (Fname1,Y3D1,"t_ff"//trim(Sce_nm),Nreg,Nmodels,Nmon, var_ln, var_unit, lon_nm, lat_nm, t_nm, Miss) ! Fname1 var_ln="Monthly regional mean of rainfall amount for the period 2080-2099, RCP"//trim(Sce_nm); var_unit="mm/mon" Call wr_ncdf_3d_apd (Fname,X3D2,"p_ff"//trim(Sce_nm),Nreg,Nmodels,Nmon, var_ln, var_unit, lon_nm, lat_nm, t_nm, Miss) Call wr_ncdf_3d_apd (Fname1,Y3D2,"p_ff"//trim(Sce_nm),Nreg,Nmodels,Nmon, var_ln, var_unit, lon_nm, lat_nm, t_nm, Miss) ! Fname1 Enddo Fname=trim(Out_Dir)//"Annual_Reg_Mean.nc" ! Baseline: var_ln="Annual regional mean of temperature for the period 1986-2005"; var_unit="C" Call wr_ncdf_2d (Fname,Xreg,Xmod,Nreg,Nmodels,BL_T2m_RAnn,"t_bl", var_ln, var_unit, lon_nm, lat_nm, Miss) var_ln="Annual regional mean of rainfall amount for the period 1986-2005"; var_unit="mm/mon" Call wr_ncdf_2d_apd (Fname,BL_Pr_RAnn,"p_bl",Nreg,Nmodels, var_ln, var_unit, lon_nm, lat_nm, Miss) ! Obs ann regional mean: Aug 2020 var_ln="Observed Annual regional mean of temperature for the period 1986-2005"; var_unit="C" Call wr_ncdf_1d_apd (Fname,Ob_T2m_RAnn,"tobs", Nreg, var_ln, var_unit, lon_nm, Miss) var_ln="Observed Annual regional mean of precipitation for the period 1986-2005"; var_unit="mm/yr" Call wr_ncdf_1d_apd (Fname,Ob_Pr_RAnn,"pobs", Nreg, var_ln, var_unit, lon_nm, Miss) Allocate (X2D1(Nreg,Nmodels), X2D2(Nreg,Nmodels)) Do isce=1,NSce if (isce==1) then; Sce_nm="45"; else; Sce_nm="85"; endif X2D1(:,:)=MF_T2m_RAnn(:,:,isce); X2D2(:,:)=MF_Pr_RAnn(:,:,isce) var_ln="Annual regional mean of temperature for the period 2046-2065, RCP"//trim(Sce_nm); var_unit="C" Call wr_ncdf_2d_apd (Fname,X2D1,"t_mf"//trim(Sce_nm),Nreg,Nmodels, var_ln, var_unit, lon_nm, lat_nm, Miss) var_ln="Monthly regional mean of rainfall amount for the period 2046-2065, RCP"//trim(Sce_nm); var_unit="mm/mon" Call wr_ncdf_2d_apd (Fname,X2D2,"p_mf"//trim(Sce_nm),Nreg,Nmodels, var_ln, var_unit, lon_nm, lat_nm, Miss) X2D1(:,:)=FF_T2m_RAnn(:,:,isce); X2D2(:,:)=FF_Pr_RAnn(:,:,isce) var_ln="Annual regional mean of temperature for the period 2080-2099, RCP"//trim(Sce_nm); var_unit="C" Call wr_ncdf_2d_apd (Fname,X2D1,"t_ff"//trim(Sce_nm),Nreg,Nmodels, var_ln, var_unit, lon_nm, lat_nm, Miss) var_ln="Monthly regional mean of rainfall amount for the period 2080-2099, RCP"//trim(Sce_nm); var_unit="mm/mon" Call wr_ncdf_2d_apd (Fname,X2D2,"p_ff"//trim(Sce_nm),Nreg,Nmodels, var_ln, var_unit, lon_nm, lat_nm, Miss) ENddo End Subroutine Calc_Regional_Mean !--------------------------------------------------------------------------------- Real Function Regional_Mean (X, ireg) Integer :: ilon,ilat,it,imon,imon1,iper,imod,isce, ireg, icount Real :: X(Nlon,Nlat), tmp Logical :: IsMask icount=0; tmp=0 Do ilon=1, Nlon Do ilat=1, Nlat !if (ireg>1) then; IsMask=(Mask(ilon,ilat)==ireg); else; IsMask=(Mask(ilon,ilat)>0); endif !if (IsMask) then !Mask(ilon,ilat)==ireg) then if (Regs(ilon,ilat)==ireg) then if (X(ilon,ilat)/=Miss) then icount=icount+1; tmp=tmp+X(ilon,ilat) endif endif Enddo Enddo if (icount>0) then Regional_Mean=tmp/real(icount) else Regional_Mean=Miss endif End Function Regional_Mean !--------------------------------------------------------------------------------- Real Function Get_Point_Value(X, KDo, VDo) Integer :: ilon, ilat, ilo, ila Real :: X(Nlon, Nlat), KDo, VDo Logical :: L1, L2 Do ilon=1, Nlon-1 Do ilat=1, Nlat-1 L1=KDo>=xlon(ilon).AND.KDo<=xlon(ilon+1) L2=VDo>=xlat(ilat).AND.VDo<=xlat(ilat+1) if (L1.AND.L2) then; ilo=ilon; ila=ilat; exit; endif Enddo Enddo Get_Point_Value=X(ilo,ila) End Function Get_Point_Value !--------------------------------------------------------------------------------- Subroutine For_Can_Tho Integer :: ilon,ilat,imon,imod,isce, it, iyr, imon1,it1, it2 Real :: X(Nlon,Nlat), Xmod(Nmodels), Xmon(Nmon), Xyrs(Nyear) Real, Parameter :: KDo=105.767, VDo=10.033 ! Can Tho Station Real, Allocatable :: BL_T(:,:), BL_P(:,:) ! (Nmodels,Ntime) Real, Allocatable :: MF_T(:,:,:), MF_P(:,:,:) ! (Nmodels,NSce,Ntime) Real, Allocatable :: FF_T(:,:,:), FF_P(:,:,:) ! (Nmodels,NSce,Ntime) Real, Allocatable :: BL_T_Mon(:,:), BL_P_Mon(:,:) ! (Nmodels,Nmon) Real, Allocatable :: MF_T_Mon(:,:,:), MF_P_Mon(:,:,:) ! (Nmodels,NSce,Nmon) Real, Allocatable :: FF_T_Mon(:,:,:), FF_P_Mon(:,:,:) ! (Nmodels,NSce,Nmon) Real, Allocatable :: BL_T_Yrs(:,:), BL_P_Yrs(:,:) ! (Nmodels,Nyear) Real, Allocatable :: MF_T_Yrs(:,:,:), MF_P_Yrs(:,:,:) ! (Nmodels,NSce,Nyear) Real, Allocatable :: FF_T_Yrs(:,:,:), FF_P_Yrs(:,:,:) ! (Nmodels,NSce,Nyear) Real, Allocatable :: X1(:),X2(:),Y1(:),Y2(:), X2Dm(:,:), X2Dy(:,:) Character (Len=100) :: Fname, Sce_nm, t_unit, var_ln, var_unit, lon_nm, lat_nm, z_nm, t_nm Allocate (BL_T(Nmodels,Ntime),BL_P(Nmodels,Ntime)) Allocate (MF_T(Nmodels,NSce,Ntime),MF_P(Nmodels,NSce,Ntime)) Allocate (FF_T(Nmodels,NSce,Ntime),FF_P(Nmodels,NSce,Ntime)) Allocate (BL_T_Mon(Nmodels,Nmon),BL_P_Mon(Nmodels,Nmon)) Allocate (MF_T_Mon(Nmodels,NSce,Nmon),MF_P_Mon(Nmodels,NSce,Nmon)) Allocate (FF_T_Mon(Nmodels,NSce,Nmon),FF_P_Mon(Nmodels,NSce,Nmon)) Allocate (BL_T_Yrs(Nmodels,Nyear),BL_P_Yrs(Nmodels,Nyear)) Allocate (MF_T_Yrs(Nmodels,NSce,Nyear),MF_P_Yrs(Nmodels,NSce,Nyear)) Allocate (FF_T_Yrs(Nmodels,NSce,Nyear),FF_P_Yrs(Nmodels,NSce,Nyear)) ! Input Data (Interpolating to station location from origin data): Do imod=1, Nmodels Do it=1, Ntime X(:,:)=BL_T2m(:,:,it,imod); BL_T(imod,it)=Get_Point_Value(X, KDo, VDo) X(:,:)=BL_Pr (:,:,it,imod); BL_P(imod,it)=Get_Point_Value(X, KDo, VDo) Enddo Enddo Do isce=1, NSce Do imod=1, Nmodels Do it=1, Ntime X(:,:)=MF_T2m(:,:,it,imod,isce); MF_T(imod,isce,it)=Get_Point_Value(X, KDo, VDo) X(:,:)=MF_Pr (:,:,it,imod,isce); MF_P(imod,isce,it)=Get_Point_Value(X, KDo, VDo) X(:,:)=FF_T2m(:,:,it,imod,isce); FF_T(imod,isce,it)=Get_Point_Value(X, KDo, VDo) X(:,:)=FF_Pr (:,:,it,imod,isce); FF_P(imod,isce,it)=Get_Point_Value(X, KDo, VDo) Enddo Enddo Enddo ! Monthly mean: Allocate (X2(Nyear),Y2(Nyear),X1(Nyear),Y1(Nyear)) ! Baseline: Do imod=1,Nmodels Do imon=1, Nmon iyr=0 Do it=1, Ntime imon1=mod(it,Nmon); if (imon1==0) imon1=Nmon if (imon1==imon) then iyr=iyr+1; X1(iyr)=BL_T(imod,it); X2(iyr)=BL_P(imod,it) endif Enddo BL_T_Mon(imod,imon)=Ave(X1,iyr,Miss); BL_P_Mon(imod,imon)=Ave(X2,iyr,Miss) Enddo Enddo ! Future: DO isce=1, NSce Do imod=1,Nmodels Do imon=1, Nmon iyr=0 Do it=1, Ntime imon1=mod(it,Nmon); if (imon1==0) imon1=Nmon if (imon1==imon) then iyr=iyr+1; X1(iyr)=MF_T(imod,isce,it); X2(iyr)=MF_P(imod,isce,it) Y1(iyr)=FF_T(imod,isce,it); Y2(iyr)=FF_P(imod,isce,it) endif Enddo MF_T_Mon(imod,isce,imon)=Ave(X1,iyr,Miss); MF_P_Mon(imod,isce,imon)=Ave(X2,iyr,Miss) FF_T_Mon(imod,isce,imon)=Ave(Y1,iyr,Miss); FF_P_Mon(imod,isce,imon)=Ave(Y2,iyr,Miss) Enddo Enddo Enddo Deallocate (X1,X2,Y1,Y2) ! Yearly: Allocate (X1(Nmon),X2(Nmon),Y1(Nmon),Y2(Nmon)) Do imod=1, Nmodels ! Baseline: Do iyr=1, Nyear it1=(iyr-1)*Nmon+1; it2=iyr*Nmon X1(1:Nmon)=BL_T(imod,it1:it2); X2(1:Nmon)=BL_P(imod,it1:it2) BL_T_Yrs(imod,iyr)=Ave(X1,Nmon,Miss); BL_P_Yrs(imod,iyr)=Ave(X2,Nmon,Miss)*Nmon Enddo ! Future: Do isce=1, NSce Do iyr=1, Nyear it1=(iyr-1)*Nmon+1; it2=iyr*Nmon X1(1:Nmon)=MF_T(imod,isce,it1:it2); X2(1:Nmon)=MF_P(imod,isce,it1:it2) Y1(1:Nmon)=FF_T(imod,isce,it1:it2); Y2(1:Nmon)=FF_P(imod,isce,it1:it2) MF_T_Yrs(imod,isce,iyr)=Ave(X1,Nmon,Miss); MF_P_Yrs(imod,isce,iyr)=Ave(X2,Nmon,Miss)*Nmon FF_T_Yrs(imod,isce,iyr)=Ave(Y1,Nmon,Miss); FF_P_Yrs(imod,isce,iyr)=Ave(Y2,Nmon,Miss)*Nmon Enddo Enddo Enddo Deallocate (X1,X2,Y1,Y2) ! Wite to files: Allocate (X2Dm(Nmodels,Nmon), X2Dy(Nmodels,Nyear)) Do imon=1,Nmon; Xmon(imon)=imon; Enddo Do imod=1,Nmodels; Xmod(imod)=imod; Enddo Do iyr=1, Nyear; Xyrs(iyr)=iyr; Enddo Fname=trim(Out_Dir)//"CT_Monthly_Mean.nc" var_ln="Monthly mean 2m-temperature for Can Tho"; var_unit="C"; lon_nm="Model"; lat_nm="Mon" Call wr_ncdf_2d (Fname,Xmod,Xmon,Nmodels,Nmon,BL_T_Mon,"t_bl", var_ln, var_unit, lon_nm, lat_nm, Miss) var_ln="Monthly mean rainfall for Can Tho"; var_unit="mm/mon" Call wr_ncdf_2d_apd (Fname,BL_P_Mon,"p_bl",Nmodels,Nmon, var_ln, var_unit, lon_nm, lat_nm, Miss) Do isce=1, Nsce if (isce==1) then; Sce_nm="45"; else; Sce_nm="85"; endif X2Dm(1:Nmodels,1:Nmon)=MF_T_Mon(1:Nmodels,isce,1:Nmon) var_ln="Monthly mean 2m-temperature for Can Tho"; var_unit="C" Call wr_ncdf_2d_apd (Fname,X2Dm,"t_mf"//trim(Sce_nm),Nmodels,Nmon, var_ln, var_unit, lon_nm, lat_nm, Miss) X2Dm(1:Nmodels,1:Nmon)=MF_P_Mon(1:Nmodels,isce,1:Nmon) var_ln="Monthly mean rainfall for Can Tho"; var_unit="mm/mon" Call wr_ncdf_2d_apd (Fname,X2Dm,"p_mf"//trim(Sce_nm),Nmodels,Nmon, var_ln, var_unit, lon_nm, lat_nm, Miss) X2Dm(1:Nmodels,1:Nmon)=FF_T_Mon(1:Nmodels,isce,1:Nmon) var_ln="Monthly mean 2m-temperature for Can Tho"; var_unit="C" Call wr_ncdf_2d_apd (Fname,X2Dm,"t_ff"//trim(Sce_nm),Nmodels,Nmon, var_ln, var_unit, lon_nm, lat_nm, Miss) X2Dm(1:Nmodels,1:Nmon)=FF_P_Mon(1:Nmodels,isce,1:Nmon) var_ln="Monthly mean rainfall for Can Tho"; var_unit="mm/mon" Call wr_ncdf_2d_apd (Fname,X2Dm,"p_ff"//trim(Sce_nm),Nmodels,Nmon, var_ln, var_unit, lon_nm, lat_nm, Miss) Enddo Fname=trim(Out_Dir)//"CT_Annual_Mean.nc" var_ln="Annual mean 2m-temperature for Can Tho"; var_unit="C"; lon_nm="Model"; lat_nm="year" Call wr_ncdf_2d (Fname,Xmod,Xyrs,Nmodels,Nyear,BL_T_Yrs,"t_bl", var_ln, var_unit, lon_nm, lat_nm, Miss) var_ln="Annual mean rainfall for Can Tho"; var_unit="mm/year" Call wr_ncdf_2d_apd (Fname,BL_P_Yrs,"p_bl",Nmodels,Nyear, var_ln, var_unit, lon_nm, lat_nm, Miss) Do isce=1, Nsce if (isce==1) then; Sce_nm="45"; else; Sce_nm="85"; endif X2Dy(1:Nmodels,1:Nyear)=MF_T_Yrs(1:Nmodels,isce,1:Nyear) var_ln="Annual mean 2m-temperature for Can Tho"; var_unit="C" Call wr_ncdf_2d_apd (Fname,X2Dy,"t_mf"//trim(Sce_nm),Nmodels,Nyear, var_ln, var_unit, lon_nm, lat_nm, Miss) X2Dy(1:Nmodels,1:Nyear)=MF_P_Yrs(1:Nmodels,isce,1:Nyear) var_ln="Annual mean rainfall for Can Tho"; var_unit="mm/year" Call wr_ncdf_2d_apd (Fname,X2Dy,"p_mf"//trim(Sce_nm),Nmodels,Nyear, var_ln, var_unit, lon_nm, lat_nm, Miss) X2Dy(1:Nmodels,1:Nyear)=FF_T_Yrs(1:Nmodels,isce,1:Nyear) var_ln="Annual mean 2m-temperature for Can Tho"; var_unit="C" Call wr_ncdf_2d_apd (Fname,X2Dy,"t_ff"//trim(Sce_nm),Nmodels,Nyear, var_ln, var_unit, lon_nm, lat_nm, Miss) X2Dy(1:Nmodels,1:Nyear)=FF_P_Yrs(1:Nmodels,isce,1:Nyear) var_ln="Annual mean rainfall for Can Tho"; var_unit="mm/year" Call wr_ncdf_2d_apd (Fname,X2Dy,"p_ff"//trim(Sce_nm),Nmodels,Nyear, var_ln, var_unit, lon_nm, lat_nm, Miss) Enddo Deallocate (X2Dm,X2Dy) End Subroutine For_Can_Tho !--------------------------------------------------------------------------------- Subroutine Chars Integer :: ilon,ilat,it,imon,imon1,iper,imod,isce Integer :: iyr Real, Allocatable :: X(:),Y(:),X1(:),Y1(:), X3D(:,:,:), X4D(:,:,:,:), Xmod(:), Xmon(:), XX(:) Real :: Tmp Character (Len=100) :: Fname, Sce_nm, t_unit, var_ln, var_unit, lon_nm, lat_nm, z_nm, t_nm !print*,Nlon,Nlat,Nmodels,Nmon,Ntime,NSce, Nyear; pause Allocate (Ob_T2m_Mean(Nlon,Nlat,Nmon), Ob_Pr_Mean(Nlon,Nlat,Nmon)) ! Aug 2020 Allocate (BL_T2m_Mean(Nlon,Nlat,Nmon,Nmodels),BL_Pr_Mean(Nlon,Nlat,Nmon,Nmodels)) Allocate (MF_T2m_Mean(Nlon,Nlat,Nmon,Nmodels,NSce),MF_Pr_Mean(Nlon,Nlat,Nmon,Nmodels,NSce)) Allocate (FF_T2m_Mean(Nlon,Nlat,Nmon,Nmodels,NSce),FF_Pr_Mean(Nlon,Nlat,Nmon,Nmodels,NSce)) Allocate (Ob_T2m_Ann(Nlon,Nlat), Ob_Pr_Ann(Nlon,Nlat)) ! Aug 2020 Allocate (BL_T2m_Ann(Nlon,Nlat,Nmodels),BL_Pr_Ann(Nlon,Nlat,Nmodels)) Allocate (MF_T2m_Ann(Nlon,Nlat,Nmodels,NSce),MF_Pr_Ann(Nlon,Nlat,Nmodels,NSce)) Allocate (FF_T2m_Ann(Nlon,Nlat,Nmodels,NSce),FF_Pr_Ann(Nlon,Nlat,Nmodels,NSce)) ! Changes of Annual mean: Allocate (D_BL_T2m_Ann(Nlon,Nlat,Nmodels), D_BL_Pr_Ann(Nlon,Nlat,Nmodels)) ! Aug 2020 Allocate (D_MF_T2m_Ann(Nlon,Nlat,Nmodels,NSce),D_MF_Pr_Ann(Nlon,Nlat,Nmodels,NSce)) Allocate (D_FF_T2m_Ann(Nlon,Nlat,Nmodels,NSce),D_FF_Pr_Ann(Nlon,Nlat,Nmodels,NSce)) ! Changes of Monthly mean: Allocate (D_BL_T2m_Mon(Nlon,Nlat,Nmon,Nmodels), D_BL_Pr_Mon(Nlon,Nlat,Nmon,Nmodels)) ! Aug 2020 Allocate (D_MF_T2m_Mon(Nlon,Nlat,Nmon,Nmodels,NSce),D_MF_Pr_Mon(Nlon,Nlat,Nmon,Nmodels,NSce)) Allocate (D_FF_T2m_Mon(Nlon,Nlat,Nmon,Nmodels,NSce),D_FF_Pr_Mon(Nlon,Nlat,Nmon,Nmodels,NSce)) Allocate (X(Nyear),Y(Nyear),X1(Nyear),Y1(Nyear), XX(Ntime)) ! Xmon(Nmon), XX(Ntime)) Ob_T2m_Mean=Miss; Ob_Pr_Mean=Miss ! Aug 2020 BL_T2m_Mean=Miss; BL_Pr_Mean=Miss; MF_T2m_Mean=Miss; MF_Pr_Mean=Miss; FF_T2m_Mean=Miss; FF_Pr_Mean=Miss Ob_T2m_Ann=Miss; Ob_Pr_Ann=Miss ! Aug 2020 BL_T2m_Ann=Miss; BL_Pr_Ann=Miss; MF_T2m_Ann=Miss; MF_Pr_Ann=Miss; FF_T2m_Ann=Miss; FF_Pr_Ann=Miss D_BL_T2m_Mon=Miss; D_BL_Pr_Mon=Miss ! Aug 2020 D_MF_T2m_Mon=Miss; D_MF_Pr_Mon=Miss; D_FF_T2m_Mon=Miss; D_FF_Pr_Mon=Miss D_BL_T2m_Ann=Miss; D_BL_Pr_Ann=Miss ! Aug 2020 D_MF_T2m_Ann=Miss; D_MF_Pr_Ann=Miss; D_FF_T2m_Ann=Miss; D_FF_Pr_Ann=Miss Do ilon=1, Nlon Do ilat=1, Nlat if (Mask(ilon,ilat)>0) then print*,ilon,ilat,Nmodels,Nmon,Ntime ! Obs: ! Aug 2020 Do imon=1, Nmon iyr=0 Do it=1, Ntime imon1=mod(it,Nmon); if (imon1==0) imon1=Nmon if (imon1==imon) then iyr=iyr+1; X(iyr)=Ob_T2m(ilon,ilat,it); Y(iyr)=Ob_Pr(ilon,ilat,it) endif Enddo Ob_T2m_Mean(ilon,ilat,imon)=Ave(X,iyr,Miss) Ob_Pr_Mean(ilon,ilat,imon) =Ave(Y,iyr,Miss) Enddo !print "(12F7.1)",Ob_T2m_Mean(ilon,ilat,:); pause !Obs Annual mean: XX(1:Ntime)=Ob_T2m(ilon,ilat,1:Ntime); Ob_T2m_Ann(ilon,ilat)=Ave(XX,Ntime,Miss) XX(1:Ntime)=Ob_Pr (ilon,ilat,1:Ntime); Ob_Pr_Ann (ilon,ilat)=Ave(XX,Ntime,Miss) if (Ob_Pr_Ann (ilon,ilat)/=Miss) Ob_Pr_Ann (ilon,ilat)=Ob_Pr_Ann (ilon,ilat)*Nmon ! Ann total rainfall ! Baseline: print*,"Baseline" Do imod=1,Nmodels Do imon=1, Nmon iyr=0 Do it=1, Ntime imon1=mod(it,Nmon); if (imon1==0) imon1=Nmon if (imon1==imon) then iyr=iyr+1; X(iyr)=BL_T2m(ilon,ilat,it,imod); Y(iyr)=BL_Pr(ilon,ilat,it,imod) !print*,"iyr=",iyr endif Enddo ! BL Monthly mean: BL_T2m_Mean(ilon,ilat,imon,imod)=Ave(X,iyr,Miss) BL_Pr_Mean(ilon,ilat,imon,imod)=Ave(Y,iyr,Miss) ! BL: Differences BL & Obs monthly: ! Aug 2020 D_BL_T2m_Mon(ilon,ilat,imon,imod)=BL_T2m_Mean(ilon,ilat,imon,imod)-Ob_T2m_Mean(ilon,ilat,imon) D_BL_Pr_Mon (ilon,ilat,imon,imod)=(BL_Pr_Mean(ilon,ilat,imon,imod)-Ob_Pr_Mean(ilon,ilat,imon))/ & Ob_Pr_Mean(ilon,ilat,imon)*100.0 Enddo ! BL Annual mean: XX(1:Ntime)=BL_T2m(ilon,ilat,1:Ntime,imod); BL_T2m_Ann(ilon,ilat,imod)=Ave(XX,Ntime,Miss) XX(1:Ntime)=BL_Pr (ilon,ilat,1:Ntime,imod); BL_Pr_Ann (ilon,ilat,imod)=Ave(XX,Ntime,Miss) if (BL_Pr_Ann (ilon,ilat,imod)/=Miss) BL_Pr_Ann (ilon,ilat,imod)=BL_Pr_Ann (ilon,ilat,imod)*Nmon ! Ann total rainfall ! BL: Differences BL & Obs Ann: D_BL_T2m_Ann(ilon,ilat,imod)=BL_T2m_Ann(ilon,ilat,imod)-Ob_T2m_Ann(ilon,ilat) D_BL_pr_Ann (ilon,ilat,imod)=(BL_Pr_Ann(ilon,ilat,imod)-Ob_Pr_Ann (ilon,ilat))/ & Ob_Pr_Ann (ilon,ilat)*100.0 Enddo ! Scenarios: print*,"Scenarios" Do isce=1,NSce Do imod=1,Nmodels Do imon=1, Nmon iyr=0 Do it=1, Ntime imon1=mod(it,Nmon); if (imon1==0) imon1=Nmon if (imon1==imon) then iyr=iyr+1; X(iyr) =MF_T2m(ilon,ilat,it,imod,isce); Y(iyr) =MF_Pr(ilon,ilat,it,imod,isce) X1(iyr)=FF_T2m(ilon,ilat,it,imod,isce); Y1(iyr)=FF_Pr(ilon,ilat,it,imod,isce) endif Enddo MF_T2m_Mean(ilon,ilat,imon,imod,isce) = Ave(X, iyr,Miss) MF_Pr_Mean (ilon,ilat,imon,imod,isce) = Ave(Y, iyr,Miss) FF_T2m_Mean(ilon,ilat,imon,imod,isce) = Ave(X1,iyr,Miss) FF_Pr_Mean (ilon,ilat,imon,imod,isce) = Ave(Y1,iyr,Miss) ! Changes of Monthly mean: D_MF_T2m_Mon(ilon,ilat,imon,imod,isce)= MF_T2m_Mean(ilon,ilat,imon,imod,isce)-BL_T2m_Mean(ilon,ilat,imon,imod) D_MF_Pr_Mon (ilon,ilat,imon,imod,isce)=(MF_Pr_Mean (ilon,ilat,imon,imod,isce)-BL_Pr_Mean (ilon,ilat,imon,imod))/ & BL_Pr_Mean (ilon,ilat,imon,imod)*100.0 D_FF_T2m_Mon(ilon,ilat,imon,imod,isce)= FF_T2m_Mean(ilon,ilat,imon,imod,isce)-BL_T2m_Mean(ilon,ilat,imon,imod) D_FF_Pr_Mon (ilon,ilat,imon,imod,isce)=(FF_Pr_Mean (ilon,ilat,imon,imod,isce)-BL_Pr_Mean (ilon,ilat,imon,imod))/ & BL_Pr_Mean (ilon,ilat,imon,imod)*100.0 Enddo XX(1:Ntime)=MF_T2m(ilon,ilat,1:Ntime,imod,isce); MF_T2m_Ann(ilon,ilat,imod,isce)=Ave(XX,Ntime,Miss) XX(1:Ntime)=MF_Pr (ilon,ilat,1:Ntime,imod,isce); MF_Pr_Ann (ilon,ilat,imod,isce)=Ave(XX,Ntime,Miss) if (MF_Pr_Ann (ilon,ilat,imod,isce)/=Miss) MF_Pr_Ann (ilon,ilat,imod,isce)=MF_Pr_Ann (ilon,ilat,imod,isce)*Nmon ! Ann total rainfall XX(1:Ntime)=FF_T2m(ilon,ilat,1:Ntime,imod,isce); FF_T2m_Ann(ilon,ilat,imod,isce)=Ave(XX,Ntime,Miss) XX(1:Ntime)=FF_Pr (ilon,ilat,1:Ntime,imod,isce); FF_Pr_Ann (ilon,ilat,imod,isce)=Ave(XX,Ntime,Miss) if (FF_Pr_Ann (ilon,ilat,imod,isce)/=Miss) FF_Pr_Ann (ilon,ilat,imod,isce)=FF_Pr_Ann (ilon,ilat,imod,isce)*Nmon ! Ann total rainfall ! Changes of Annual mean: D_MF_T2m_Ann(ilon,ilat,imod,isce)= MF_T2m_Ann(ilon,ilat,imod,isce)-BL_T2m_Ann(ilon,ilat,imod) D_MF_Pr_Ann (ilon,ilat,imod,isce)=(MF_Pr_Ann (ilon,ilat,imod,isce)-BL_Pr_Ann (ilon,ilat,imod))/ & BL_Pr_Ann (ilon,ilat,imod)*100.0 D_FF_T2m_Ann(ilon,ilat,imod,isce)= FF_T2m_Ann(ilon,ilat,imod,isce)-BL_T2m_Ann(ilon,ilat,imod) D_FF_Pr_Ann (ilon,ilat,imod,isce)=(FF_Pr_Ann (ilon,ilat,imod,isce)-BL_Pr_Ann (ilon,ilat,imod))/ & BL_Pr_Ann (ilon,ilat,imod)*100.0 Enddo Enddo print*,ilon,ilat endif Enddo Enddo print*,"Deallocate..."; Deallocate (X,Y,X1,Y1,XX) ! ,Xmon) print*,"Writing in Chars ..."; pause ! Write to files: Monthly Mean Allocate (X4D(Nlon,Nlat,Nmodels,Nmon), Xmod(Nmodels), Xmon(Nmon)) Do imon=1,Nmon; Xmon(imon)=imon; Enddo Do imod=1,Nmodels; Xmod(imod)=imod; Enddo ! BL Mean Temp Do ilon=1, Nlon Do ilat=1, Nlat Do imod=1,Nmodels X4D(ilon,ilat,imod,1:Nmon)=BL_T2m_Mean(ilon,ilat,1:Nmon,imod) Enddo Enddo Enddo t_unit="months since 1986-01-01"; var_unit="C"; lon_nm="lon"; lat_nm="lat"; z_nm="model"; t_nm="mon" var_ln="Monthly mean temperature for the period 1986-2005 projeted by different models" Fname=trim(Out_Dir)//"Monthly_Mean.nc" Call wr_ncdf_4d (Fname,xlon,xlat,Xmod,Xmon,nlon,nlat,Nmodels,Nmon, X4D,"t_bl", & t_unit, var_ln, var_unit, lon_nm, lat_nm, z_nm, t_nm, Miss) ! BL Mean Prec Do ilon=1, Nlon Do ilat=1, Nlat Do imod=1,Nmodels X4D(ilon,ilat,imod,1:Nmon)=BL_Pr_Mean(ilon,ilat,1:Nmon,imod) Enddo Enddo Enddo var_ln="Monthly mean rainfall amount for the period 1986-2005 projeted by different models"; var_unit="mm/mon" Call wr_ncdf_4d_apd (Fname, X4D,"p_bl", nlon,nlat,Nmodels,Nmon, & var_ln,var_unit, lon_nm, lat_nm, z_nm, t_nm,Miss) !--------------------------------------------- ! Obs Mean Temp - Aug 2020 var_ln="Observed monthly mean temperature for the period 1986-2005"; var_unit="C" Call wr_ncdf_3d_apd (Fname,Ob_T2m_Mean,"tobs",nlon,nlat,Nmon, var_ln, var_unit, lon_nm, lat_nm, t_nm, Miss) ! Obs monthly Prec - Aug 2020 var_ln="Observed monthly precipitation for the period 1986-2005"; var_unit="mm/mon" Call wr_ncdf_3d_apd (Fname,Ob_Pr_Mean,"pobs",nlon,nlat,Nmon, var_ln, var_unit, lon_nm, lat_nm, t_nm, Miss) !--------------------------------------------- print*,"Here ==>" Do isce=1,NSce if (isce==1) then; Sce_nm="45"; else; Sce_nm="85"; endif ! MF mean Temp Do ilon=1, Nlon Do ilat=1, Nlat Do imod=1,Nmodels X4D(ilon,ilat,imod,1:Nmon)=MF_T2m_Mean(ilon,ilat,1:Nmon,imod,isce) Enddo Enddo Enddo var_ln="Monthly mean temperature for the period 2046-2065 projeted by different models and scenario of RCP"//trim(Sce_nm) Call wr_ncdf_4d_apd (Fname, X4D, "t_mf"//trim(Sce_nm),nlon,nlat,Nmodels,Nmon, & var_ln,var_unit, lon_nm, lat_nm, z_nm, t_nm,Miss) ! MF mean Prec Do ilon=1, Nlon Do ilat=1, Nlat Do imod=1,Nmodels X4D(ilon,ilat,imod,1:Nmon)=MF_Pr_Mean(ilon,ilat,1:Nmon,imod,isce) Enddo Enddo Enddo var_ln="Monthly mean rainfall amount for the period 2046-2065 projeted by different models and scenario of RCP"//trim(Sce_nm) Call wr_ncdf_4d_apd (Fname, X4D, "p_mf"//trim(Sce_nm),nlon,nlat,Nmodels,Nmon, & var_ln,var_unit, lon_nm, lat_nm, z_nm, t_nm,Miss) ! FF mean Temp Do ilon=1, Nlon Do ilat=1, Nlat Do imod=1,Nmodels X4D(ilon,ilat,imod,1:Nmon)=FF_T2m_Mean(ilon,ilat,1:Nmon,imod,isce) Enddo Enddo Enddo var_ln="Monthly mean temperature for the period 2080-2099 projeted by different models and scenario of RCP"//trim(Sce_nm) Call wr_ncdf_4d_apd (Fname, X4D, "t_ff"//trim(Sce_nm),nlon,nlat,Nmodels,Nmon, & var_ln,var_unit, lon_nm, lat_nm, z_nm, t_nm,Miss) ! FF mean Prec Do ilon=1, Nlon Do ilat=1, Nlat Do imod=1,Nmodels X4D(ilon,ilat,imod,1:Nmon)=FF_Pr_Mean(ilon,ilat,1:Nmon,imod,isce) Enddo Enddo Enddo var_ln="Monthly mean rainfall amount for the period 2080-2099 projeted by different models and scenario of RCP"//trim(Sce_nm) Call wr_ncdf_4d_apd (Fname, X4D, "p_ff"//trim(Sce_nm),nlon,nlat,Nmodels,Nmon, & var_ln,var_unit, lon_nm, lat_nm, z_nm, t_nm,Miss) Enddo ! Write to files: Annual mean Fname=trim(Out_Dir)//"Annual_Mean.nc" Allocate (X3D(Nlon,Nlat,Nmodels)) ! Testing only: !======================================================== !X3D(:,:,:)=BL_T2m_Ann(:,:,:) !BL_T2m_Ann(:,:,1:3)=X3D(:,:,5:7); BL_T2m_Ann(:,:,4:7)=X3D(:,:,1:4) !X3D(:,:,:)=BL_Pr_Ann(:,:,:) !BL_Pr_Ann(:,:,1:3)=X3D(:,:,5:7); BL_Pr_Ann(:,:,4:7)=X3D(:,:,1:4) !X4D(:,:,:,:)=MF_T2m_Ann(:,:,:,:) !MF_T2m_Ann(:,:,1:3,:)=X4D(:,:,5:7,:); MF_T2m_Ann(:,:,4:7,:)=X4D(:,:,1:4,:) !X4D(:,:,:,:)=MF_Pr_Ann(:,:,:,:) !MF_Pr_Ann(:,:,1:3,:)=X4D(:,:,5:7,:); MF_Pr_Ann(:,:,4:7,:)=X4D(:,:,1:4,:) !X4D(:,:,:,:)=FF_T2m_Ann(:,:,:,:) !FF_T2m_Ann(:,:,1:3,:)=X4D(:,:,5:7,:); FF_T2m_Ann(:,:,4:7,:)=X4D(:,:,1:4,:) !X4D(:,:,:,:)=FF_Pr_Ann(:,:,:,:) !FF_Pr_Ann(:,:,1:3,:)=X4D(:,:,5:7,:); FF_Pr_Ann(:,:,4:7,:)=X4D(:,:,1:4,:) ! End of testing !======================================================== var_ln="Annual mean temperature for the period 1986-2005"; t_nm="model" X3D(:,:,:)=BL_T2m_Ann(:,:,:) Call wr_ncdf_3d (Fname,xlon,xlat,Xmod,nlon,nlat,Nmodels,X3D,"t_bl",t_unit, var_ln, var_unit, lon_nm, lat_nm, t_nm, Miss) !Call wr_ncdf_3d (Fname,xlon,xlat,Xmod,nlon,nlat,Nmodels,BL_T2m_Ann,"t_bl",t_unit, var_ln, var_unit, lon_nm, lat_nm, t_nm, Miss) var_ln="Annual total rainfall amount for the period 1986-2005" X3D(:,:,:)=BL_Pr_Ann(:,:,:) Call wr_ncdf_3d_apd (Fname,X3D,"p_bl",nlon,nlat,Nmodels, var_ln, var_unit, lon_nm, lat_nm, t_nm, Miss) ! Obs Ann: Aug 2020 var_ln="Observed annual mean temperature for the period 1986-2005" Call wr_ncdf_2d_apd (Fname, Ob_T2m_Ann,"tannob",nlon,nlat, var_ln, "C", lon_nm, lat_nm, Miss) var_ln="Observed annual prec for the period 1986-2005" Call wr_ncdf_2d_apd (Fname, Ob_Pr_Ann,"pannob",nlon,nlat, var_ln, "C", lon_nm, lat_nm, Miss) !Call wr_ncdf_3d_apd (Fname,BL_Pr_Ann,"p_bl",nlon,nlat,Nmodels, var_ln, var_unit, lon_nm, lat_nm, t_nm, Miss) !Allocate (X3D(Nlon,Nlat,Nmodels)) Do isce=1,NSce if (isce==1) then; Sce_nm="45"; else; Sce_nm="85"; endif X3D(:,:,:)=MF_T2m_Ann(:,:,:,isce) var_ln="Annual mean temperature for the period 2046-2065, RCP"//trim(Sce_nm); var_unit="C" Call wr_ncdf_3d_apd (Fname,X3D,"t_mf"//trim(Sce_nm),nlon,nlat,Nmodels, var_ln, var_unit, lon_nm, lat_nm, t_nm, Miss) X3D(:,:,:)=MF_Pr_Ann(:,:,:,isce) var_ln="Annual total rainfall amount for the period 2046-2065, RCP"//trim(Sce_nm); var_unit="mm/mon" Call wr_ncdf_3d_apd (Fname,X3D,"p_mf"//trim(Sce_nm),nlon,nlat,Nmodels, var_ln, var_unit, lon_nm, lat_nm, t_nm, Miss) X3D(:,:,:)=FF_T2m_Ann(:,:,:,isce) var_ln="Annual mean temperature for the period 2080-2099, RCP"//trim(Sce_nm); var_unit="C" Call wr_ncdf_3d_apd (Fname,X3D,"t_ff"//trim(Sce_nm),nlon,nlat,Nmodels, var_ln, var_unit, lon_nm, lat_nm, t_nm, Miss) X3D(:,:,:)=FF_Pr_Ann(:,:,:,isce) var_ln="Annual total rainfall amount for the period 2080-2099, RCP"//trim(Sce_nm); var_unit="mm/mon" Call wr_ncdf_3d_apd (Fname,X3D,"p_ff"//trim(Sce_nm),nlon,nlat,Nmodels, var_ln, var_unit, lon_nm, lat_nm, t_nm, Miss) Enddo ! Changes in Temp and Prec: ! Changes of Annual mean: Fname=trim(Out_Dir)//"Change_Annual_Mean.nc" Do isce=1,NSce if (isce==1) then; Sce_nm="45"; else; Sce_nm="85"; endif var_ln="Changes in Annual mean temperature for the period 2046-2065, RCP"//trim(Sce_nm); t_nm="model"; var_unit="C" X3D(:,:,:)=D_MF_T2m_Ann(:,:,:,isce) if (isce==1) then Call wr_ncdf_3d (Fname,xlon,xlat,Xmod,nlon,nlat,Nmodels,X3D,"t_mf"//trim(Sce_nm),t_unit, var_ln, var_unit, lon_nm, lat_nm, t_nm, Miss) else Call wr_ncdf_3d_apd (Fname,X3D,"t_mf"//trim(Sce_nm),nlon,nlat,Nmodels, var_ln, var_unit, lon_nm, lat_nm, t_nm, Miss) endif var_ln="Changes in Annual mean rainfall for the period 2046-2065, RCP"//trim(Sce_nm); t_nm="model"; var_unit="mm/mon" X3D(:,:,:)=D_MF_Pr_Ann(:,:,:,isce) Call wr_ncdf_3d_apd (Fname,X3D,"p_mf"//trim(Sce_nm),nlon,nlat,Nmodels, var_ln, var_unit, lon_nm, lat_nm, t_nm, Miss) var_ln="Changes in Annual mean temperature for the period 2080-2099, RCP"//trim(Sce_nm); t_nm="model"; var_unit="C" X3D(:,:,:)=D_FF_T2m_Ann(:,:,:,isce) Call wr_ncdf_3d_apd (Fname,X3D,"t_ff"//trim(Sce_nm),nlon,nlat,Nmodels, var_ln, var_unit, lon_nm, lat_nm, t_nm, Miss) var_ln="Changes in Annual mean rainfall for the period 2080-2099, RCP"//trim(Sce_nm); t_nm="model"; var_unit="mm/mon" X3D(:,:,:)=D_FF_Pr_Ann(:,:,:,isce) Call wr_ncdf_3d_apd (Fname,X3D,"p_ff"//trim(Sce_nm),nlon,nlat,Nmodels, var_ln, var_unit, lon_nm, lat_nm, t_nm, Miss) Enddo ! BL: Obs-BL - Aug 2020: D_BL_T2m_Ann var_ln="Differences in Annual mean temperature between Obs and BL for the period 1986-2005" Call wr_ncdf_3d_apd (Fname,D_BL_T2m_Ann,"delt",nlon,nlat,Nmodels, var_ln, "C", lon_nm, lat_nm, t_nm, Miss) var_ln="Differences in Annual preci between Obs and BL for the period 1986-2005" Call wr_ncdf_3d_apd (Fname,D_BL_PR_Ann,"delp",nlon,nlat,Nmodels, var_ln, "mm/yr",lon_nm, lat_nm, t_nm, Miss) ! Changes of Monthly mean: Fname=trim(Out_Dir)//"Change_Monthly_Mean.nc" z_nm="model"; t_nm="mon" Do isce=1,NSce if (isce==1) then; Sce_nm="45"; else; Sce_nm="85"; endif ! MF mean Temp Do ilon=1, Nlon Do ilat=1, Nlat Do imod=1,Nmodels X4D(ilon,ilat,imod,1:Nmon)=D_MF_T2m_Mon(ilon,ilat,1:Nmon,imod,isce) Enddo Enddo Enddo var_ln="Change in Monthly mean temperature for the period 2046-2065, RCP"//trim(Sce_nm); var_unit="C" if (isce==1) then Call wr_ncdf_4d (Fname,xlon,xlat,Xmod,Xmon,nlon,nlat,Nmodels,Nmon, X4D,"t_mf"//trim(Sce_nm), & t_unit, var_ln, var_unit, lon_nm, lat_nm, z_nm, t_nm, Miss) else Call wr_ncdf_4d_apd (Fname, X4D, "t_mf"//trim(Sce_nm),nlon,nlat,Nmodels,Nmon, & var_ln,var_unit, lon_nm, lat_nm, z_nm, t_nm,Miss) endif ! MF mean Prec Do ilon=1, Nlon Do ilat=1, Nlat Do imod=1,Nmodels X4D(ilon,ilat,imod,1:Nmon)=D_MF_Pr_Mon(ilon,ilat,1:Nmon,imod,isce) Enddo Enddo Enddo var_ln="Change in Monthly mean rainfall amount for the period 2046-2065,RCP"//trim(Sce_nm); var_unit="mm/mon" Call wr_ncdf_4d_apd (Fname, X4D, "p_mf"//trim(Sce_nm),nlon,nlat,Nmodels,Nmon, & var_ln,var_unit, lon_nm, lat_nm, z_nm, t_nm,Miss) ! FF mean Temp Do ilon=1, Nlon Do ilat=1, Nlat Do imod=1,Nmodels X4D(ilon,ilat,imod,1:Nmon)=D_FF_T2m_Mon(ilon,ilat,1:Nmon,imod,isce) Enddo Enddo Enddo var_ln="Change in Monthly mean temperature for the period 2080-2099, RCP"//trim(Sce_nm); var_unit="C" Call wr_ncdf_4d_apd (Fname, X4D, "t_ff"//trim(Sce_nm),nlon,nlat,Nmodels,Nmon, & var_ln,var_unit, lon_nm, lat_nm, z_nm, t_nm,Miss) ! FF mean Prec Do ilon=1, Nlon Do ilat=1, Nlat Do imod=1,Nmodels X4D(ilon,ilat,imod,1:Nmon)=D_FF_Pr_Mon(ilon,ilat,1:Nmon,imod,isce) Enddo Enddo Enddo var_ln="Change in Monthly mean rainfall amount for the period 2080-2099, RCP"//trim(Sce_nm); var_unit="mm/mon" Call wr_ncdf_4d_apd (Fname, X4D, "p_ff"//trim(Sce_nm),nlon,nlat,Nmodels,Nmon, & var_ln,var_unit, lon_nm, lat_nm, z_nm, t_nm,Miss) Enddo ! BL: Differences between Obs & BL - Aug 2020: print*,"Differences between Obs & BL - Aug 2020" var_ln="Differences of Monthly mean temperature between BL and Obs for the period 1986-2005" Call wr_ncdf_4d_apd (Fname, D_BL_T2m_Mon,"delt",nlon,nlat,Nmodels,Nmon, & var_ln,"C", lon_nm, lat_nm, z_nm, t_nm,Miss) var_ln="Differences of Monthly preci between BL and Obs for the period 1986-2005" Call wr_ncdf_4d_apd (Fname, D_BL_Pr_Mon,"delp",nlon,nlat,Nmodels,Nmon, & var_ln,"C", lon_nm, lat_nm, z_nm, t_nm,Miss) Deallocate (X3D, X4D, Xmod, Xmon) print*,"End of Differences between Obs & BL - Aug 2020" !Deallocate (X3D, X4D, Xmod, Xmon) End Subroutine Chars !--------------------------------------------------------------------------------- Subroutine Read_T_R (Fname) Character (Len=*) :: Fname Real, Allocatable :: xlon(:), xlat(:), xtime(:), T(:,:,:), R(:,:,:) Integer :: nlon, nlat, nt Call rd_ncdf_3d (trim(Fname), xlon, xlat, xtime, nlon, nlat, nt0, T2m, "t2m") Deallocate (xlon, xlat, xtime) Call rd_ncdf_3d (trim(Fname), xlon, xlat, xtime, nlon, nlat, nt0, Pr, "pr") Deallocate (xlon, xlat, xtime) !if (.not.Allocated(T2m)) Allocate (T2m(Nlon,Nlat,Nt)) !if (.not.Allocated(Pr )) Allocate (Pr (Nlon,Nlat,Nt)) !print*,"In file "//trim(Fname),": ",Nlon,Nlat,Nt; pause !T2m(:,:,:)=T(:,:,:); Pr(:,:,:)=R(:,:,:) !Deallocate (T,R) !Ntime = Nt !print*,"In file "//trim(Fname),": ",Nlon,Nlat,Nt; pause End Subroutine Read_T_R !--------------------------------------------------------------------------------- Subroutine Read_Land_Mask Integer :: nx, ny Real, Allocatable :: lonx(:), latx(:) Call rd_ncdf_2d (trim(LandMaskFile), xlon,xlat, nlon,nlat,Mask,"land") ! all=2..8=7 regs !Allocate (Regs(nlon,nlat)) !Regs=Mask Call rd_ncdf_2d (trim(RegionFile), lonx, latx, nx, ny, Regs, "reg") End Subroutine Read_Land_Mask !--------------------------------------------------------------------------------- Logical Function Empty_Ar1(X,N) Integer :: N,i Real :: X(N) Logical :: Empty Empty=.False. Do i=1,N if (X(i)==Miss.OR.IsNaN(X(i))) then Empty=.True.; Exit endif Enddo Empty_Ar1=Empty End Function Empty_Ar1 !--------------------------------------------------------------------------------- Real Function Ave(X, N, Miss) Integer :: N, i, icount Real :: X(N), Mean, Miss Mean=0.; icount=0 Do i=1, N if (X(i)/=Miss) then icount=icount+1; Mean=Mean+X(i) endif Enddo if (icount>0) then Mean=Mean/real(icount) else Mean=Miss endif Ave=Mean End Function Ave !--------------------------------------------------------------------------------- Include "/home/tanpv/Mylibs/sub_rd_ncdf_Tan.F90" Include "/home/tanpv/Mylibs/sub_wr_ncdf_Tan.F90" !--------------------------------------------------------------------------------- END