setwd("/Users/new_user/Desktop/for_thay/APH_Char/")
library(ncdf4); library(ggplot2); library(reshape2)
library("gridExtra");library("cowplot"); library(ggpubr)

c.clas <- c("Ext","Mod","Sev")
cvar <- c("df","tdd","tds")

for (ivar in  c(cvar)){
     remove(draw.dat)
     remove(c)
     idem <- 1
     for (iclas in c(c.clas)) {
          remove(in.file); remove(idat); remove(sel.dat);remove(mel.dat)
          in.file <- paste("APH_Chars_Reg_",iclas,sep="")
          innc <- nc_open(paste(in.file,".nc",sep=""))
          idat <- ncvar_get(innc,ivar)
          reg <- round(ncvar_get(innc,"reg"),digits = 3); rownames(idat) <- paste("R",reg,sep=""); nreg <- dim(reg)
          time <- round(ncvar_get(innc,"decade"),digits = 3); colnames(idat) <- time; ntime <- dim(time)
          idx <- ncvar_get(innc,"idx"); ndx <- dim(idx)
          sel.dat <- idat[,,1] # selet idx=1 => SPI???
          mel.dat <- melt(sel.dat,id=colnames(sel.dat))
          mel.dat$clas <- iclas
          mel.dat$var <- ivar
          colnames(mel.dat) <- c("reg","mon","value","clas","var")
          if (idem == 1){
            draw.dat <- mel.dat
          } else {
            draw.dat <- rbind(draw.dat,mel.dat)
          }
          idem <- idem + 1
     }
     
     print(max(draw.dat$value))
     print(min(draw.dat$value))
     # for plot 
     
     creg <- rownames(idat)
     c.tile <- c("R1","R2","R3","R4","SEA")
     itmp <- 1
     
     remove(cmax,cmin,cint)
     if (ivar=="df"){
         cmax <- 8
         cmin <- 0
         cint <- 2
     }
     
     
     if (ivar=="tdd"){
       cmax <- 80
       cmin <- 0
       cint <- 20
     }
     
     if (ivar=="tds"){
       cmax <- 0
       cmin <- -120
       cint <- 40
     }
     
     for (ireg in c(creg)){
       tmp.dat <- subset(draw.dat,draw.dat$reg==ireg)
       tmp.dat$clas <- factor(tmp.dat$clas,levels=c("Mod","Sev","Ext"))
       
       p <-ggplot(data = tmp.dat,aes(x=mon,y=value,fill=clas,width=4.85,col=clas))+
         geom_bar(stat="identity", position=position_dodge())+
         #### for white-black figure    
         #scale_fill_manual(name="",values = c("#cccccc","#666666","#0d0d0d"))+
         #scale_color_manual(name="",values = c("#cccccc","#666666","#0d0d0d"))+
         #### for color firgure
         scale_fill_brewer(name="",palette = "Reds")+
         scale_color_brewer(name="",palette = "Reds")+
         scale_x_continuous(breaks=c(1962,1973,1984,1995,2005,2015),
                         labels=c("62-72","73-83","84-94",
                                   "95-05","06-15","62-15"))+
         scale_y_continuous(breaks=seq(cmin,cmax,cint),limits = c(cmin,cmax))+
         theme_classic()+ggtitle(c.tile[itmp])+
         theme(legend.position="bottom",
               axis.text.x=element_text(angle=0))+
         xlab("")+ylab("")
       
       p <- p+theme(legend.position="none")
       
       eval(parse(text=paste("p",itmp,"<-p",sep="")))
       itmp <- itmp + 1
       
     }
     
     
     c <- ggarrange(p1, p2, p3, p4, p5, nrow=5, 
                    common.legend = TRUE, legend="bottom", heights = c(rep(4,5),1.75))
     pdf(file=paste("APH_",ivar,".pdf",sep=""), width = 3.5, height = 7.0)
     print(c)
     dev.off()
     
}