Use NetCDF Real, Parameter :: Miss=-99.0 Integer, Parameter :: IMiss=-99, Nreg=7, Pentad=5 !Integer, Parameter :: Period=5, N5y=Period*73 ! 365 ! =5yrs x 73 Integer, Parameter :: Period=1, N5y=Period*73 ! =1yr x 73 Integer, Parameter :: Ncha=4 ! b(1), number of FD; b(2), mean duration of FD; b(3), mean severity of FD; b(4), mean speed of FD Integer, Allocatable :: time(:) Integer :: i, LenSt, sys Character (Len=255) :: OutDir, Fname, t_unit, DirIn, Str, time_unit, Fname1, Fname2 Character (Len=255), Allocatable :: Files(:) Character (Len=12) :: date0, timeres, idate0, idate1 Character (Len=1) :: Ch Integer :: yy1,mm1,dd1,hh1,yy2,mm2,dd2,hh2 Real*8 :: factor, offset, Missing, fill_value Real, Allocatable :: Prec(:,:,:), lon(:), lat(:), lev(:), xtime(:), Regs(:,:) Real, Allocatable, Dimension (:,:,:) :: sw1 ! (Nlon,Nlat,Ntime) Real, Allocatable, Dimension (:,:) :: sw1_Reg ! (Nreg,Ntime) Real, Allocatable, Dimension (:,:,:,:) :: psw5 ! (Nlon,Nlat,Nyrs,N5y) - Tinh cho tung thoi ky 5 nam Real, Allocatable, Dimension (:,:,:,:) :: FDCha5 ! (Nlon,Nlat,Ncha,Nyrs) Real, Allocatable, Dimension (:,:,:) :: psw5_Reg ! (Nreg,Nyrs,N5y) Real, Allocatable, Dimension (:,:,:) :: FDCha5_Reg ! (Nreg,Ncha,Nyrs) Real, Allocatable, Dimension (:,:,:) :: psw1 ! (Nlon,Nlat,Ntime) Real, Allocatable, Dimension (:,:,:) :: FDevt ! drt(Nlon,Nlat,Ntime) - 1-drought onset, 2-drought recovery, 0-no drought Real, Allocatable, Dimension (:,:,:) :: FDCha ! (Nlon,Nlat, Ncha) Real, Allocatable, Dimension (:,:) :: psw1_Reg ! (Nreg,Ntime) Real, Allocatable, Dimension (:,:) :: FDevt_Reg ! drt(Nreg,Ntime) Real, Allocatable, Dimension (:,:) :: FDCha_Reg ! (Nreg, Ncha) Real, Allocatable, Dimension (:,:,:,:) :: X4d Integer, Allocatable, Dimension (:) :: date_m1 ! (Nfcs) Integer :: Nlon, Nlat, Nlev, Ntime, Nyrs, Ndays, Npentad, Nlon0, Nlat0 ! Nfcs = So lan du bao, date_fcs = yyyymmdd cua ngay DB Integer :: imem, ifcs, iy, im, id !---------------------------------------------------------------------------------------------- Call Read_Regs ("/home/tanpv/Data/ERA5/Mask_VN_0.25.nc") !VN.grid_0.25.nc") Call Read_sw_percentile_data Call Calculating_FD_Chars OutDir="/home/tanpv/Progs/Flash_Dr/NC/" if (Period==5) then Call Write_data_to_file (trim(OutDir)//"FD_Char_1961-2020_5yrs.nc",trim(OutDir)//"FD_Char_1961-2020_5yrs_Reg.nc") else Call Write_data_to_file (trim(OutDir)//"FD_Char_1961-2020_1yrs.nc",trim(OutDir)//"FD_Char_1961-2020_1yrs_Reg.nc") endif CONTAINS !---------------------------------------------------------------------------------------------- Include "FL_Dr.f90" !---------------------------------------------------------------------------------------------- Subroutine Calculating_FD_Chars Integer :: ilon, ilat, it, ireg Real :: X(Ntime), drt(Ntime), b(Ncha) ! b(1), number of FD; b(2), mean duration of FD; b(3), mean severity of FD; b(4), mean speed of FD Allocate (FDCha5(Nlon,Nlat,Ncha,Nyrs)) Allocate (FDCha5_Reg(Nreg,Ncha,Nyrs)) FDCha5=Miss; FDCha5_Reg=Miss Do ilon=1,Nlon Do ilat=1,Nlat if (Regs(ilon,ilat)>0) then Do it=1,Nyrs X(1:Ntime)=psw5(ilon,ilat,it,1:Ntime) Call flashdrought(Ntime,X,b,drt) b(1) = b(1)/real(Pentad) ! number of FD per year FDCha5(ilon,ilat,1:Ncha,it) =b(1:Ncha) Enddo endif Enddo Enddo Do ireg=1,Nreg Do it=1,Nyrs X(1:Ntime)=psw5_Reg(ireg,it,1:Ntime) Call flashdrought(Ntime,X,b,drt) b(1) = b(1)/real(Pentad) ! number of FD per year FDCha5_Reg(ireg,1:Ncha,it) =b(1:Ncha) Enddo Enddo End Subroutine Calculating_FD_Chars !---------------------------------------------------------------------------------------------- Subroutine Read_Regs (Fname) Character (Len=*) :: Fname Integer :: ilon,ilat,Nlo,Nla, i1, i2,i3,i4,i5,i6,i7,i8,i9 Call rd_ncdf_2d (Fname,Lon,Lat,Nlon,Nlat,Regs,"reg") End Subroutine Read_Regs !---------------------------------------------------------------------------------------------- Subroutine Read_sw_percentile_data Integer :: Nt, it, j, iday, ilon,ilat, ifile, Nobs, i1,i2 Integer :: Nlo, Nla Real, Allocatable :: xlon(:),xlat(:),xtime(:),Y1d(:) Real :: tmp Fname="/home/tanpv/Progs/Flash_Dr/NC/Pentad_sw_1961-2020.nc" !Call rd_ncdf_3d (Fname,xlon,xlat,xtime,nlo,nla,nt,sw1,"sw") !Deallocate (xlon,xlat,xtime) Call rd_ncdf_3d (Fname,xlon,xlat,xtime,nlo,nla,nt,psw1,"swp") Deallocate (xlon,xlat,xtime) Fname="/home/tanpv/Progs/Flash_Dr/NC/Pentad_sw_1961-2020_Reg.nc" !Call rd_ncdf_2d (Fname,xlon,xtime,nlo,nt,sw1_Reg,"sw") !Deallocate (xlon,xtime) Call rd_ncdf_2d (Fname,xlon,xtime,nlo,nt,psw1_Reg,"swp") Deallocate (xlon,xtime) Nyrs=Nt/N5y Ntime=N5y Allocate (psw5(Nlon,Nlat,Nyrs,N5y),psw5_Reg(Nreg,Nyrs,N5y)) Do ilon=1,Nlon Do ilat=1,Nlat Do j=1,Nyrs i1=(j-1)*N5y+1; i2=j*N5y psw5(ilon,ilat,j,1:N5y)=psw1(ilon,ilat,i1:i2) Enddo Enddo Enddo Do ireg=1,Nreg Do j=1,Nyrs i1=(j-1)*N5y+1; i2=j*N5y psw5_Reg(ireg,j,1:N5y)=psw1_Reg(ireg,i1:i2) Enddo Enddo End Subroutine Read_sw_percentile_data !---------------------------------------------------------------------------------------------- Subroutine Write_data_to_file (Fname,Fnm1) ! ,X4d,X1d,Nlon,Nlat,Nt,Nfcs) Character (Len=*) :: Fname,Fnm1 ! Daily and Moving Pentad ! Dat(Nlon,Nlat,Nfiles,Ndays) Integer :: it,ifile,ireg,Ntime Real, Allocatable :: X3d(:,:,:),xday(:),X2d(:,:),xreg(:),xcha(:) Character (Len=50) :: t_unit, lon_nm, lat_nm, z_nm, t_nm Character (Len=10) :: var_unit Character (Len=150) :: var_ln(3)=(/"FD chars: b(1), number of FD per year; b(2), mean duration of FD; b(3), mean severity of FD; b(4), mean speed of FD", & "Daily volumetric soil water percemtile layer 2 (7 - 28cm)", & "Daily volumetric soil water percemtile layer 3 (28 - 100cm)"/) ! Allocate (FDCha5(Nlon,Nlat,Ncha,Nyrs)) ! Allocate (FDCha5_Reg(Nreg,Ncha,Nyrs)) Ntime=Nyrs Allocate (xday(Ntime),xreg(Nreg),xcha(Ncha)) Do ireg=1,Nreg; xreg(ireg)=ireg; Enddo Do it=1,Ntime; xday(it)=(it-1)*Period; Enddo Do it=1,Ncha; xcha(it)=it; Enddo var_unit="-"; t_unit="years since 1961-01-01 "; lon_nm="lon"; lat_nm="lat"; z_nm="cha"; t_nm="year" Call wr_ncdf_4d (Fname,Lon,Lat,xcha,xday,Nlon,Nlat,Ncha,Nyrs,FDCha5,"fd", & t_unit, trim(var_ln(1)), var_unit, lon_nm, lat_nm, z_nm, t_nm, Miss) ! Call wr_ncdf_3d (Fname,Lon,Lat,xcha,Nlon,Nlat,Ncha,FDCha,"fd",t_unit,trim(var_ln(1)),var_unit,lon_nm,lat_nm,t_nm,Miss) ! Call wr_ncdf_3d_apd (Fname,psw2,"sw2",Nlon,Nlat,Ntime,trim(var_ln(2)),var_unit,lon_nm,lat_nm,t_nm,Miss) ! Call wr_ncdf_3d_apd (Fname,psw3,"sw3",Nlon,Nlat,Ntime,trim(var_ln(3)),var_unit,lon_nm,lat_nm,t_nm,Miss) ! Areal mean of daily data: Do ireg=1,Nreg; xreg(ireg)=ireg; Enddo lon_nm="reg"; lat_nm="cha"; t_nm="year" Call wr_ncdf_3d (Fnm1,xreg,xcha,xday,Nreg,Ncha,Nyrs,FDCha5_Reg,"fd",t_unit,trim(var_ln(1)),var_unit,lon_nm,lat_nm,t_nm,Miss) ! Call wr_ncdf_2d (Fnm1,xreg,xcha,Nreg,Ncha,FDCha_Reg,"fd",trim(var_ln(1)),var_unit,"reg","cha",Miss) ! Call wr_ncdf_2d_apd (Fnm1,psw2_Reg,"sw2",Nreg,Ntime,"Areal mean of "//trim(var_ln(2)),var_unit,"reg","day",Miss) ! Call wr_ncdf_2d_apd (Fnm1,psw3_Reg,"sw3",Nreg,Ntime,"Areal mean of "//trim(var_ln(3)),var_unit,"reg","day",Miss) End Subroutine Write_data_to_file !---------------------------------------------------------------------------------------------- Include "/home/tanpv/Mylibs/sub_rd_ncdf_Tan.F90" Include "/home/tanpv/Mylibs/sub_wr_ncdf_Tan.F90" Include "/home/tanpv/Mylibs/Ave.f90" Include "/home/tanpv/Mylibs/Xmax.f90" Include "/home/tanpv/Mylibs/Xmin.f90" Include "/home/tanpv/Mylibs/Xsum.f90" Include "/home/tanpv/Mylibs/Sub_Percentile.f90" Include "/home/tanpv/Mylibs/Percentile_rank.f90" END