! PROGRAM Verify_QM Use Netcdf Integer, Parameter :: Yr1=1982, Yr2=2011, Nmon=12, Nreg=7, NLead=5, Nvar=4, Nfile=9 Integer, Parameter :: Nyrs=(Yr2-Yr1+1) Real, Parameter :: Miss=-99. Character (Len=3) :: vname (Nvar)=(/"obs","mod","bjp","q_m"/) Character (Len=3) :: vname0(Nvar)=(/"obs","bjp","q_m","mod"/) Character (Len=255) :: Fname Real, Allocatable :: lon(:),lat(:),Regs(:,:), X3d(:,:,:) Real, Allocatable :: xlon(:),xlat(:),xlev(:),xtime(:) Integer :: Nlo,Nla,Nz,Nt Integer :: Nlon, Nlat, Ntime, Ngrd, ilon, ilat, it Integer :: Nmodel,Nele Real, Allocatable :: Dat_All(:,:,:,:,:) ! (Nlon,Nlat,Ntime,Nvar,NLead) Real, Allocatable :: X_inp(:,:), Y_inp(:), Yest(:,:), Yhq(:,:) Character (Len=5),Allocatable :: yvarn(:), Model_nm(:) !-------------------------------------------------------------------------------- ! Verification: Real, Allocatable, Dimension(:,:,:,:,:) :: Pr_ME, Pr_MAE, Pr_RME, Pr_RMAE, Pr_COR, Pr_Pval, Pr_AV ! (Nlon,Nlat,Nmon,NLead,Nvar) Integer :: Nyear, imon ! Tong so nam co SL va chi so thang !-------------------------------------------------------------------------------- ! ANN: character (len=250) :: input_file, output_file, save_file character (len=250) :: Dir, Dir_ANN_Out, Dir_Cli, Dir_ANN_Save, Dir_Veri, Dir_4Dr character (len=250) :: Dir_QM, Dir_BJP logical :: training_mode ! Use for calculating scaling factors Integer :: Nsample, Nepoc, Nsample1 Real,allocatable, Dimension(:,:) :: Inputs, Outputs, Inputs_X, Outputs_Y Real,allocatable, Dimension(:) :: min_in, max_in, min_out, max_out Real :: mean_error, mean_error_check !------------------------------------------------------------------------------- integer :: Nin, Nhidden, Nout, iin, ihid, iout integer :: i,j,k, iele, imod Real,allocatable, Dimension(:,:) :: weight_ih, weight_ho, moment_ih, moment_ho Real,allocatable, Dimension(:) :: fin, fhid, fout, ftarget, delta_o, delta_h ! Training Real :: learning_rate0, learning_rate, momentum, error, tmp, mintmp logical :: debug, Run_Cli !-------------------------------------------------------------------------------- Integer :: Op_Nhidden, Op_Nepoc ! Optimal values Real :: Op_learning_rate Real, Allocatable :: X4d(:,:,:,:), X5d(:,:,:,:,:), F1d(:),O1d(:), Score_Cor(:,:,:,:), Score_AV(:,:,:,:) !------------------- For ANN ------------------------------------------------------------- Namelist /Model_info_1/ Nmodel, Nele Namelist /Model_info_2/ Model_nm, yvarn !-------------------------------------------------------------------------------- Open (1, file="Namelist.txt",Status="old") Read (1, Model_info_1) Allocate (yvarn(Nele), Model_nm(Nmodel)) Read (1, Model_info_2) Fname="/home/tanpv/Data/AWC/Mask_VN_0.25_Corrected.nc" Call rd_ncdf_2d (Fname,lon,lat,Nlon,Nlat,Regs,"reg") ! 7 sub-regs mask Nyear=Nyrs Nvar1=Nvar-1 ! model, bjp, q_m Allocate (Pr_RME(Nlon,Nlat,NLead,Nmon,Nvar1),Pr_RMAE(Nlon,Nlat,NLead,Nmon,Nvar1),Pr_COR(Nlon,Nlat,NLead,Nmon,Nvar1),Pr_Pval(Nlon,Nlat,NLead,Nmon,Nvar1)) Allocate (Pr_ME (Nlon,Nlat,NLead,Nmon,Nvar1),Pr_MAE (Nlon,Nlat,NLead,Nmon,Nvar1),Pr_AV (Nlon,Nlat,NLead,Nmon,Nvar1)) Pr_ME=Miss; Pr_MAE=Miss; Pr_RME=Miss; Pr_RMAE=Miss; Pr_COR=Miss; Pr_Pval=Miss; Pr_AV=Miss Allocate (Score_Cor(Nreg,Nmon,NLead,Nvar1), Score_AV(Nreg,Nmon,NLead,Nvar1)) ! Number of gridpoints with Cor>0, AV>0 Dir_Veri ="/home/tanpv/Data/562_VNU/0.25_VN/Verify/" !Dir_QM ="/home/tanpv/Data/562_VNU/0.25_VN/QM_Ha/"//"/"//trim(yvarn(iele))//"/"//"Baseline/" Dir_BJP ="/home/tanpv/Data/562_VNU/0.25_VN/BJP/" !-------------------------------------------------------------------------------- Do iele=1,Nele Do imod=1,Nmodel Dir_QM="/home/tanpv/Data/562_VNU/0.25_VN/QM_Ha/"//trim(yvarn(iele))//"/"//"Baseline/" Call Read_all_data (Dir_BJP,Dir_QM,Model_nm(imod),yvarn(iele)) ! pause print*,"Deading data done for ",Model_nm(imod),yvarn(iele) !-------------------------------------------------------------------------------- Call Verification Fname=trim(Dir_Veri)//"Veri_Training_BMJ_QM_"//trim(yvarn(iele))//"_"//trim(Model_nm(imod))//".nc" Call Write_Verification (Fname,yvarn(iele),Model_nm(imod)) Open (3,file=trim(Dir_Veri)//"Score_Cor_AV_Training_"//trim(yvarn(iele))//"_"//trim(Model_nm(imod))//".txt",status="Unknown") Do ivar=1,Nvar1 write (3,"(A)") "Percentage of grid points with positive correlation and added value of "//trim(vname0(ivar+1)) Do ilead=1,NLead write (3,"(A,I1.1,A)") "Lead time = ",ilead," month(s)" write (3,"(6x,A6,12I6)") " Reg. ",(imon, imon=1,Nmon) Do ireg=1,Nreg write (3,"(A6,I6,12F6.1)") "Cor - ",ireg,(Score_Cor(ireg,imon,ilead,ivar),imon=1,Nmon) write (3,"(2A6,12F6.1)") "AV - "," ",(Score_AV (ireg,imon,ilead,ivar),imon=1,Nmon) Enddo Enddo Enddo Close (3) Enddo Enddo !-------------------------------------------------------------------------------- CONTAINS !-------------------------------------------------------------------------------- Subroutine Verification !(ilead,imon) Integer :: ilon,ilat,imon,imon1,ilead,ireg, ic, ic0, ic1, ic2, it Integer, Parameter :: N=Nyrs !-1 Real :: obs(N),bjp(N),model(N),q_m(N),Corr,Prob Do ilead=1,NLead Do ilon=1,Nlon Do ilat=1,Nlat if (Int(Regs(ilon,ilat))>0) then Do imon=1,Nmon ic=0 Do it=1,Ntime imon1=mod(it,Nmon); if (imon1==0) imon1=Nmon if (imon1==imon) then ic=ic+1 obs(ic) =Dat_All(ilon,ilat,it,1,ilead); model(ic)=Dat_All(ilon,ilat,it,2,ilead) bjp(ic) =Dat_All(ilon,ilat,it,3,ilead); q_m(ic)=Dat_All(ilon,ilat,it,4,ilead) endif Enddo if (ic>N) then; print*,"Something's wrong in Verification"; stop; endif Pr_ME (ilon,ilat,ilead,imon,1) = ME (bjp, obs,ic,Miss) Pr_ME (ilon,ilat,ilead,imon,3) = ME (model,obs,ic,Miss) Pr_ME (ilon,ilat,ilead,imon,2) = ME (q_m, obs,ic,Miss) Pr_MAE (ilon,ilat,ilead,imon,1) = MAE (bjp, obs,ic,Miss) Pr_MAE (ilon,ilat,ilead,imon,3) = MAE (model,obs,ic,Miss) Pr_MAE (ilon,ilat,ilead,imon,2) = MAE (q_m, obs,ic,Miss) Pr_RME (ilon,ilat,ilead,imon,1) = RME_1 (bjp, obs,ic,Miss) Pr_RME (ilon,ilat,ilead,imon,3) = RME_1 (model,obs,ic,Miss) Pr_RME (ilon,ilat,ilead,imon,2) = RME_1 (q_m, obs,ic,Miss) Pr_RMAE(ilon,ilat,ilead,imon,1) = RMAE_1(bjp, obs,ic,Miss) Pr_RMAE(ilon,ilat,ilead,imon,3) = RMAE_1(model,obs,ic,Miss) Pr_RMAE(ilon,ilat,ilead,imon,2) = RMAE_1(q_m, obs,ic,Miss) Call Correlation (bjp, obs,ic,Corr,Prob,Miss); Pr_COR(ilon,ilat,ilead,imon,1)=Corr; Pr_Pval(ilon,ilat,ilead,imon,1)=Prob Call Correlation (model,obs,ic,Corr,Prob,Miss); Pr_COR(ilon,ilat,ilead,imon,3)=Corr; Pr_Pval(ilon,ilat,ilead,imon,3)=Prob Call Correlation (q_m, obs,ic,Corr,Prob,Miss); Pr_COR(ilon,ilat,ilead,imon,2)=Corr; Pr_Pval(ilon,ilat,ilead,imon,2)=Prob Pr_AV (ilon,ilat,ilead,imon,1) = AV (bjp,model,obs,ic,Miss) ! Added value of bjp vs model Pr_AV (ilon,ilat,ilead,imon,2) = AV (q_m,model,obs,ic,Miss) ! Added value of q_m vs model Pr_AV (ilon,ilat,ilead,imon,3) = Miss ! Not use Enddo endif Enddo Enddo Enddo ! Do ivar=1,Nvar1 Do ilead=1,NLead Do imon=1,Nmon Do ireg=1,Nreg ic0=0; ic1=0; ic2=0 Do ilon=1,Nlon Do ilat=1,Nlat if (Int(Regs(ilon,ilat))==ireg) then ic0=ic0+1 if (Pr_COR (ilon,ilat,ilead,imon,ivar)>0.0) ic1=ic1+1 if (Pr_AV (ilon,ilat,ilead,imon,ivar)>0.0) ic2=ic2+1 endif Enddo Enddo Score_Cor(ireg,imon,ilead,ivar)=real(ic1)/real(ic0)*100.0 if (ivar