require(rgdal)
require(rmapshaper)
require(maptools)
require(data.table)
require(stringr)
require(dplyr)

#1 Prostorová data
setwd()
#Loading
#--------------------
povodi_0 <- readOGR('data/geo/E_HEIS$UPV_HLGP#P2$wm.shp', 'E_HEIS$UPV_HLGP#P2$wm')
reky_0 <- readOGR('data/geo/E_ISVS$UPOV_R.shp', 'E_ISVS$UPOV_R')
jezera_0 <- readOGR('data/geo/E_ISVS$UPOV_J.shp', 'E_ISVS$UPOV_J')
nadrze_0 <- readOGR('data/geo/nadrze.shp')

#Preparation of data
#--------------------
povodi_0 <- spTransform(povodi_0, CRS("+init=epsg:4326"))
povodi <- ms_simplify(povodi_0, keep_shapes = TRUE, keep = 0.05)
#--------------------
proj4string(reky_0) <- CRS("+title=Krovak JTSK +proj=krovak +lat_0=49.5 +lon_0=24.83333333333333 +alpha=30.28813975277778 +k=0.9999 +x_0=0 +y_0=0 +ellps=bessel +units=m +towgs84=570.8,85.7,462.8,4.998,1.587,5.261,3.56 +no_defs <>")
reky_0 <- spTransform(reky_0, CRS("+init=epsg:4326"))
reky <- ms_simplify(reky_0, keep_shapes = TRUE, keep = 0.10)
#--------------------
proj4string(jezera_0) <- CRS("+title=Krovak JTSK +proj=krovak +lat_0=49.5 +lon_0=24.83333333333333 +alpha=30.28813975277778 +k=0.9999 +x_0=0 +y_0=0 +ellps=bessel +units=m +towgs84=570.8,85.7,462.8,4.998,1.587,5.261,3.56 +no_defs <>")
jezera_0 <- spTransform(jezera_0, CRS("+init=epsg:4326"))
jezera <- ms_simplify(jezera_0, keep_shapes = TRUE, keep = 0.05)

proj4string(nadrze_0) <- CRS("+title=Krovak JTSK +proj=krovak +lat_0=49.5 +lon_0=24.83333333333333 +alpha=30.28813975277778 +k=0.9999 +x_0=0 +y_0=0 +ellps=bessel +units=m +towgs84=570.8,85.7,462.8,4.998,1.587,5.261,3.56 +no_defs <>")
nadrze <- spTransform(nadrze_0, CRS("+init=epsg:4326"))

kraje_0 <- readOGR('data/geo/admin/kraje.shp')
kraje <- spTransform(kraje_0, CRS("+init=epsg:4326"))
kraje <- ms_simplify(kraje_0, keep_shapes = TRUE, keep = 0.10)

a <- kraje_0$NAZEV
a[1] <- "Ústecký"
kraje_0$NAZEV <- kraje_0$NK
kraje_0$NAZEV <- a


okresy_0 <- readOGR('data/geo/admin/okresy.shp')
okresy_0 <- spTransform(okresy_0, CRS("+init=epsg:4326"))
okresy <- ms_simplify(okresy_0, keep_shapes = TRUE, keep = 0.10)

#Saving
#--------------------
#writeSpatialShape(povodi,"data/prep/povodi")
# writeSpatialShape(reky,"data/prep/reky")
# writeSpatialShape(jezera,"data/prep/jezera")

writeOGR(povodi, "data/prep", "povodi", driver="ESRI Shapefile", encoding  = "UTF-8")
writeOGR(reky, "data/prep", "reky", driver="ESRI Shapefile", encoding  = "UTF-8")
writeOGR(jezera, "data/prep", "jezera", driver="ESRI Shapefile", encoding  = "UTF-8")
writeOGR(nadrze, "data/prep", "nadrze", driver="ESRI Shapefile", encoding  = "UTF-8")
writeOGR(okresy, "data/prep", "okresy", driver="ESRI Shapefile", encoding  = "UTF-8")

kraje <- readOGR("data/prep/kraje.shp")
a <- as.character(kraje$NAZEV)
a[2] <- "Ústecký"
kraje <- ms_simplify(kraje, keep_shapes = TRUE, keep = 0.10)
kraje$NAZEV <- a 

writeSpatialShape(kraje,"data/prep/kraje")

povodi_III_0 <- readOGR('data/geo/A08_Povodi_III.shp')
povodi_III_0 <- spTransform(povodi_III_0, CRS("+init=epsg:4326"))
povodi_III <- ms_simplify(povodi_III_0, keep_shapes = TRUE, keep = 0.10)

writeSpatialShape(povodi_III,"data/prep/povodi_III")

#2 Měsiční bilance

#Loading
#--------------------
r = readRDS('BER_0010.rds')

#Preparation of data
#--------------------
i ='BER_0010.rds'
M = list()
for (i in dir()){
  r = readRDS(i)
  mr = melt(r, id.vars = 'DTM')
  m1 = mr[variable!='T', .(value = sum(value)), by = .(year(DTM), month(DTM), variable)]
  m2 = mr[variable=='T', .(value = mean(value)), by = .(year(DTM), month(DTM), variable)]
  m = rbind(m1,m2)
  m[, DTM:=as.Date(paste(year, month, 1, sep = '-'))]
  M[[length(M)+1]] = m
}

names(M) = gsub('\\.rds', '', dir())
BM = rbindlist(M, idcol = 'UPOV_ID')

#Saving
#--------------------

mesice <- c("Leden","Unor","Brezen","Duben","Kveten","Cerven",
            "Cervenec","Srpen","Zari","Rijen","Listopad","Prosinec")

seasons <- data.frame(month = c(1:12), seasons = c("Zima", "Zima", "Jaro", "Jaro", "Jaro", "Leto", "Leto",
                                                   "Leto","Podzim", "Podzim", "Podzim", "Zima"))

BM <- BM %>% left_join(seasons, by="month") 
BM$month2 <- as.factor(BM$month)
levels(BM$month2) <- mesice

BM <- BM %>% group_by(UPOV_ID, year, variable) %>% mutate(annual.avg = mean(value)) %>% ungroup

BM <- BM %>% group_by(UPOV_ID, variable) %>% mutate(mean_ep = mean(value)) %>%  ungroup

BM.long <- dcast(BM, month+year+UPOV_ID+DTM~variable)
BM.long$m <- as.factor(BM.long$month)
levels(BM.long$m) <- mesice

BM <- as.data.table(BM)

saveRDS(BM, file.path(.datadir, 'webapp_data/mbilan/bilan_month.rds'))
saveRDS(BM.long, file.path(.datadir, 'webapp_data/mbilan/bilan_month_long.rds'))
saveRDS(BMt, file.path(.datadir, 'webapp_data/mbilan/bilan_month_data_table.rds'))

#3 Denní průtoky
#--------------------
setwd("..")

#Loading
#--------------------

stanice_0 <- readOGR("data/chmu/156_stanic.shp")
seznam.st <- read.csv('data/chmu/156_stanic_seznam.csv',encoding = 'UTF-8', header = TRUE, sep = ";", 
                      colClasses = c("factor", "character", "character", "character", "character"))

QD <- read.table('data/chmu/QD_156_stanic.txt',encoding = 'UTF-8', header = TRUE, sep=',', 
                 colClasses = c("factor", "character", "numeric"), col.names = c("DBCN", "DTM", "value"))

#Preparation of data
#--------------------
colnames(seznam.st)[2] <- "NAZEV.STANICE"
colnames(seznam.st)[4] <- "OBDOBI.S.DATY.DENNICH.PRUTOKU"

stanice <- spTransform(stanice_0, CRS("+init=epsg:4326"))

QD$DTM <- as.Date(QD$DTM, format = "%d.%m.%Y")
QD <- merge(seznam.st,QD, by="DBCN")

QD <- select(QD, -OBDOBI.S.DATY.DENNICH.PRUTOKU, -NAZEV.STANICE, -TOK, -Plocha..km2.)

# QD <- readRDS(file.choose())

#Saving
#--------------------
saveRDS(QD, "QD.rds")
writeOGR(stanice, "data/prep", "stanice", driver="ESRI Shapefile", encoding  = "UTF-8")


#4 uzivani_ocistene

#Loading
#--------------------
u <- read.csv('data/vuv/uzivani_utvary_06_16.csv',encoding = 'UTF-8', header = TRUE, sep = ";")
u <- readRDS('data/uzivani/79_15/uzivani.rds')
#Preparation of data
#--------------------

u <- u[complete.cases(X, Y)]
u$DTM <- as.character(u$DTM)

colnames(u)[1] <- "ICOC"
u = melt(u, id.vars = c("ICOC", "JEV", "CZ_NACE", "POVODI", "NAZEV", "ROCNI.MNOZSTVI.tis.m3", "POVOLENE.MNOZSTVI.ROK.tis.m3", "SOUR_X",	"SOUR_Y",	"UPOV_ID",	"HEIS_POZN", "ROK"
))

u$variable <- gsub("MVM*", "", u$variable)
u$value <- gsub(",", ".", u$value)

u$DTM <- paste0("01",str_sub(paste0("00",u$variable),-2,-1), u$ROK, sep="-")
u$DTM <- as.Date(u$DTM, format = "%d%m%Y")
u$value <- as.numeric(u$value)

colnames(u)[13] <- "MESIC"
colnames(u)[c(8,9)] <- c("X", "Y")
colnames(u)[5] <- "NAZICO"

u <- u[!is.na(u$X)|!is.na(u$Y),] #vynechani bodu s neznamymi souradnicemi 

#Saving
#--------------------
setwd("C:/Users/Irina/Disk Google/1_ČZU/Sucho")
saveRDS(u, "data/uzivani/06_16/uzivani_NA.rds") #obsahuje NA
saveRDS(u, 'data/uzivani/06_16/uzivani_na_rm.rds') #vynechane NA

### doplneni zaznamu od bodu se stejnym ICOC
u <- readRDS('data/uzivani/06_16/uzivani_NA.rds')

u$na <- as.numeric(is.na(u$X))

u <- u %>% group_by(ICOC) %>% mutate(mean.X = mean(X, na.rm = T), mean.Y = mean(Y, na.rm = T)) %>% ungroup
u$X[is.na(u$X)] <- u$mean.X[is.na(u$X)]
u$Y[is.na(u$Y)] <- u$mean.Y[is.na(u$Y)]

u <- u[!is.na(u$X)|!is.na(u$Y),]

saveRDS(u, 'data/uzivani/06_16/uzivani_na_nahraz.rds') # NA nahrazene průměrem za ICOC

u <- select(u, -ROCNI.MNOZSTVI.tis.m3, -POVOLENE.MNOZSTVI.ROK.tis.m3 , -HEIS_POZN, -na, -CZ_NACE)

saveRDS(u, file.path(.datadir, 'webapp_data/uzivani/06_16/uzivani_na_nahraz.rds')) # NA nahrazene průměrem za ICOC zmenseny

#5 uzivani / UPOV_ID
#--------------------
require(rgeos)

povodi <- readOGR("data/prep/povodi.shp")
popis 
povodi <- sp::merge(povodi, popis, by='UPOV_ID')

u <- readRDS('data/uzivani/79_15/uzivani.rds')

u <- u[complete.cases(X, Y)]
xy <- SpatialPoints(u[, c("X", "Y")], proj4string = CRS("+init=epsg:2065"))
xy <- spTransform(xy, CRS("+init=epsg:4326") )

pocitadlo <- 0

for (i in povodi$UPOV_ID) {
  kde <- gIntersects(povodi[povodi$UPOV_ID==i,], xy, byid = TRUE)
  u$UPOV_ID[kde[,1]] <- i
  pocitadlo <- pocitadlo+1
  print(paste(i, (pocitadlo/1121)*100, '%'))
}

 saveRDS(u, 'data/uzivani/79_15/uzivani_upovid.rds')
 
 #6 indikatory
 #--------------------
 
 spi <- readRDS(file.path(.datadir, "indikatory/spi.Rds"))
 spi$month <- month(spi$DTM) 
 spi$year <- year(spi$DTM)

 spei <- readRDS(file.path(.datadir, "indikatory/spei.Rds"))
 
 spi <- select(spi, -IID)
 spei <- select(spei, -IID)
 
 saveRDS(spi, file.path(.datadir, "webapp_data/indikatory/spi.rds")) 
 saveRDS(spei, file.path(.datadir, "webapp_data/indikatory/spei.rds")) 
 
 spi <- as.data.table(spi)
 spei <- as.data.table(spei)
 
 .datadir <- "/home/irina/ownCloud/Shared/BILAN_UPOV/used_data"
 
 spi <- readRDS(file.path(.datadir, "webapp_data/spi.rds")) 
 spei <- readRDS(file.path(.datadir, "webapp_data/spei.rds")) 
 
 for(i in unique(spi$scale)){
  x <- spi[scale == i]
  saveRDS(x,
          paste0("/home/irysek/ownCloud/Shared/BILAN_UPOV/used_data/webapp_data/indikatory/", "spi_", i, ".rds"))
 }
 
 for(i in unique(spei$scale)){
   x <- spei[scale == i]
   saveRDS(x,
           paste0("/home/irysek/ownCloud/Shared/BILAN_UPOV/used_data/webapp_data/indikatory/", "spei_", i, ".rds"))
 }
 
 spei_filter <- spei[scale == input$krok]
 
 #7 pars
 #--------------------
 pars <- readRDS(file.path(.datadir, "pars/pars.rds"))
 
 data <- c()
 name <- names(pars)
 for(i in 1:length(pars)){
   tab <- pars[[i]]
   data <- rbind(data,
                 data.frame(UPOV_ID = rep(name[i],6), tab$pars))
   print(paste("Processing....",i))
   
 }
 
 saveRDS(data, file.path(.datadir, "webapp_data/pars/pars.rds"))
 
#8 popis txt to rds
#--------------------
 
popis <- read.table('used_data/webapp_data/E_ISVS$UTV_POV.txt',encoding = 'UTF-8', header = TRUE, sep=';')
popis <- select(popis, UPOV_ID, NAZ_UTVAR, NAZ_POVODI, NAZ_OBLAST, KTG_UPOV, U_PMU)
saveRDS(popis, file.path(.datadir, "webapp_data/popis.rds"))

#9 uzivani_leaflet
#--------------------

u <- readRDS(file.path(.datadir, "webapp_data/uzivani/06_16/uzivani_na_nahraz.rds"))
u <- as.data.table(u)
saveRDS(u, file.path(.datadir, "webapp_data/uzivani/06_16/uzivani_na_nahraz_dt.rds"))

  u_leaflet <- u %>% group_by(UPOV_ID, NAZICO, POVODI, JEV, X, Y, ICOC) %>% summarise(mean=mean(value))
  u_leaflet <- SpatialPointsDataFrame(u_leaflet[, c("X", "Y")], u_leaflet, proj4string = CRS("+init=epsg:2065"))
  u_leaflet <- spTransform(u_leaflet, CRS("+init=epsg:4326"))
  
saveRDS(u_leaflet, file.path(.datadir, "webapp_data/uzivani/u_leaflet.rds"))

#10 cara prekroceni rds
#--------------------

BM <- readRDS(file.path(.datadir, "webapp_data/mbilan/bilan_month.rds"))

cp <- BM %>% group_by(UPOV_ID, variable) %>% mutate(k = rank(-value), n = n()) %>% mutate(p_year = (k-0.3)/(n+0.4)) 

cp <- cp %>%  group_by(UPOV_ID, variable, seasons) %>% mutate(k = rank(-value), n = n()) %>% mutate(p_season = (k-0.3)/(n+0.4))

cp <- cp %>%  group_by(UPOV_ID, variable, month) %>% mutate(k = rank(-value), n = n()) %>% mutate(p_month = (k-0.3)/(n+0.4))

cp_final <- cp %>% ungroup() %>% select(UPOV_ID, variable, p_year, p_month, p_season, seasons, month2, value)

saveRDS(cp_final, file.path(.datadir, "webapp_data/to_from/cara_prekroceni.rds"))

cp <- as.data.table(cp)
saveRDS(cp, file.path(.datadir, "webapp_data/to_from/cara_prekroceni_dt.rds"))

#11 m-denní prutoky (chars)
#--------------------
chars<- readRDS(file.path(.datadir, "chmu/chars_mm.rds"))

chars <- chars %>% select(UPOV_ID, Qa, choices.mday)

saveRDS(chars, file.path(.datadir, "webapp_data/chmu/chars_mm.rds"))

#qmd

chars <- readRDS(file.path(.datadir, "webapp_data/chmu/chars_mm.rds"))

r = readRDS('BER_0010.rds')

i ='BER_0010.rds'
M = list()
for (i in dir()){
  r = readRDS(i)
  mr = melt(r, id.vars = 'DTM')
  m1 = mr[variable=='TRM_mm_d', .(value = mean(value)), by = .(year(DTM), month(DTM), variable)]
  m2 = mr[variable=='TRM_m3s', .(value = mean(value)), by = .(year(DTM), month(DTM), variable)]
  m = rbind(m1,m2)
  m[, DTM:=as.Date(paste(year, month, 1, sep = '-'))]
  M[[length(M)+1]] = m
}

names(M) = gsub('\\.rds', '', dir())
chars.sim = rbindlist(M, idcol = 'UPOV_ID')

chars.sim.mm0 <- chars.sim[variable=="TRM_mm_d",]

vektor <- gsub('\\.rds', '', dir())

new <- c()
for (i in vektor){
  new <- rbind(new,c(i,CatCa::mdq(chars.sim.mm0$value[chars.sim.mm0$UPOV_ID == i], return = 'vector')))
}

new <- as.data.frame(new)

for(j in c(2:14)){
  new[,j] <- as.numeric(new[,j])
  }

colnames(new)[1] <- "UPOV_ID"  

saveRDS(new, "chars_mm_sim.rds")


#12 validace
#--------------------

dbcn_to_upov <- readRDS(file.path(.datadir, "chmu/dbcn_uid_mon.rds"))

QD <- merge(dbcn_to_upov, QD, by = 'DBCN')
QD <- QD[order(QD$DTM),]
saveRDS(QD, file.choose())

mdta_val <- readRDS(file.path(.datadir,"chmu/mdta_2016.rds"))
mdta_val <- mdta_val[, c(1,3,4,6,7,8)]
mdta_val <- merge(mdta_val, dbcn_to_upov, by = 'DBCN')
saveRDS(mdta_val, file.path(.datadir, "webapp_data/chmu/mdta_2016.rds"))


#13 simulovana mesicni data validace
#--------------------

list_of_files <- dir(file.path(.datadir, "routed-0_bilan_bez_VN/"))

dta <- c()  

for(i in list_of_files){
  x <- readRDS(file.path(.datadir, "routed-0_bilan_bez_VN/",i))
  upov_name <- gsub(".rds","",i)
  x$UPOV_ID <- upov_name
  dta <- rbind(dta,x)
  print(paste("scream ",round(which(i== list_of_files)/length(list_of_files),4)*100,"%"))
}

saveRDS(dta, file.path(.datadir,"used_data/webapp_data/chmu/dta_routed0.rds"))

dta$m <- as.factor(substr(dta$DTM,6,7))
dta$y <- as.factor(substr(dta$DTM,1,4))

dta_m <- dta %>% group_by(m,y,UPOV_ID) %>% summarise(RM = mean(TRM_mm_d)) %>% mutate(DTM = as.Date(paste(y,m,"15",sep="-")))
new_dta_m <- merge(mdta_val, dta_m, by=c("DTM","UPOV_ID"))

saveRDS(new_dta_m, file.path(.datadir, "webapp_data/chmu/data_validace.rds"))