use.case2_map.R

############################################################################
############################################################################
###                                                                      ###
###                              SECTION 2:                              ###
###                              CAPRI MAP                               ###
###                                                                      ###
############################################################################
############################################################################


#################################################################
##                            Setup                            ##
#################################################################
# set gdx path
gdx.dir <- paste0(getwd(), "/inst/extdata/gdx")

# set output data path
outdata.dir <- paste0(getwd(), "/inst/outdata")

# set gdx file name which you want to load and visualization.
file1 <- "res_2_1230cap_after_2014_refdefaultA" # baseline
file2 <- "res_2_1230fta_import_bandefaultA" # scenario


#################################################################
##                             MAP                             ##
#################################################################

# set commodities
# CERE RAPE SOYA PULS
prod1 = c("CERE","RAPE","SOYA","PULS")
# MAIF OFAR GRAI GRAE
prod2 = c("MAIF", "OFAR","GRAI","GRAE")

prods = c(prod1,prod2)

map_data_use <- caprir::map_data %>%
  filter(!grepl("TR.*", CAPRI_NUTS_ID)) %>% # TR: Tuerkiye
  filter(!grepl("UK.*", CAPRI_NUTS_ID))  %>%# UK:
  filter(!grepl("PT20.*", CAPRI_NUTS_ID)) %>% # PT20/PT30
  filter(!grepl("PT30.*", CAPRI_NUTS_ID)) %>%# PT20/PT30
  filter(!grepl("ES70.*", CAPRI_NUTS_ID))  # ES70


# load baseline
benchmark <- capri_data(filename = paste0(gdx.dir,"/",file1),
                          selregion = "all",
                          seldim5 = "",
                          selcols = prods,
                          selrows = "LEVL" ,
                          simyear = "2030",
                          scenarioname = "benchmark")
# load scenario
scenario <- capri_data(filename = paste0(gdx.dir,"/",file2),
                         selregion = "all",
                         seldim5 = "",
                         selcols = prods,
                         selrows = "LEVL" ,
                         simyear = "2030",
                         scenarioname = "scenario")
# prod 1, percent change
p1 <- map_capri(baseline = benchmark, scenario = scenario, prods = prod1, comparison = TRUE, percent_change = TRUE, 11)
# prod 2, percent change
p2 <- map_capri(baseline = benchmark, scenario = scenario, prods = prod2, comparison = TRUE, percent_change = TRUE, 11)

p1
p2

# save as png
ggsave(paste0( outdata.dir,"/png/cereals.png"), plot = p1,width = 16, height = 9, dpi = 100)
ggsave(paste0( outdata.dir,"/png/MAIF.png"), plot = p2,width = 16, height = 9, dpi = 100)

## test case for package: sf
library(rgdal)
library(tidyverse)
library(sp)
library(sf)
# load kml file
patialpolydf <- readOGR("D:/public/yang/CAPRIR_Project/caprir/inst/extdata/map/capri_map.kml", "capri_map")

# patialpolydf@data <- patialpolydf@data %>%
  # left_join(mySel, by = c("Description" = "region")) %>% select(Name,Description,cols,min)

# Convert spatialpolydf to an sf object
spatialpolysf <- patialpolydf %>% st_as_sf()


new_sf <- spatialpolysf %>% filter(!grepl(paste(c("TR.*", "UK.*","PT20.*","PT30.*","ES70.*"), collapse = "|"), Description))

# observation_place <- patialpolysf@data


final <- left_join(new_sf, mySel, by = c("Description" = "region"))

# mySel$region <- sub("0+$", "",as.character(mySel$region))
#
# spatialpolysf <- observation_place %>%
#   left_join(mySel, by = c("Description" = "region"))


ofar <- final[final$cols =="OFAR",]
soya <- final[final$cols =="SOYA",]

plot(ofar["min"])
plot(soya["min"])

plot(nc["AREA"], key.pos = 4)