all crops done. treenuts, veggies, EU subsetting. Function R script created....

all crops done. treenuts, veggies, EU subsetting. Function R script created. Suitability file created and transformed, linked to zones AEZ ADM0 Soilagri.

Functions_cifos/function_collection.R

+easypackages::packages(c("tidyverse", 'sf', 'raster','sp','dplyr','spatialEco',
+                         'exactextractr','tibble', 'here', "readr",                         
+                         "readxl", "foreign", "data.table", "spData", "purrr", "rgdal", "lwgeom",
+                         "plyr"))
+
+
+# Crop variable conversion ------------------------------------------------
+# Calculate production from yields and area - here: suitable yields
+FunProdCalcSuit = function(Area,Yield){
+  bind_rows(Area,Yield) %>% 
+    dplyr::select(-soil_code, -id_suit) %>% 
+    pivot_longer(cols = "wheat":"rest", 
+                 names_to = "crops") %>% 
+    dplyr::select(-unit) %>% 
+    pivot_wider(names_from = rec_type, values_from = value) %>% 
+    mutate(SPR_ton = SAR*YIA/1000) %>% 
+    dplyr::rename(SAR_ha = SAR, 
+                  YIA_kgha = YIA)
+}
+
+
+# Aggregation functions ---------------------------------------------------
+# Aggregating suitability files per zone, crop and tech_type
+Fun_aggregate = function(suitability_dat){
+  suitability_dat %>% 
+    dplyr::rename(tech_type = level) %>% 
+    dplyr::select(ADM0_CODE, AEZ, Soil_agri, tech_type, crops, SAR_ha, SPR_ton) %>% 
+    dplyr::group_by(ADM0_CODE, AEZ, Soil_agri, tech_type, crops) %>% 
+    dplyr::summarise(SAR_ha = sum(SAR_ha), 
+                     SPR_ton = sum(SPR_ton))
+} 

Suitability_SPAM/Suitability_SPAM.R

+easypackages::packages(c("tidyverse", 'sf', 'raster','sp','dplyr','spatialEco',
+                         'exactextractr','tibble', 'here', "readr",                         
+                         "readxl", "foreign", "data.table", "spData", "purrr", "rgdal", "lwgeom",
+                         "plyr"))
+
+# Sourcing functions from external R script 
+source("Functions_cifos/function_collection.R")
+
+# Loading the data(was handed over from Ulrike July 2021)  -------------------------------------------------------
+# Area
+suit_p_all_H_SAR = read_csv(here::here("Input_data", "cropsuit_p_all_SAR_H_cell5m.csv"))
+suit_p_all_L_SAR = read_csv(here::here("Input_data", "cropsuit_p_all_SAR_L_cell5m.csv"))
+suit_p_all_I_SAR = read_csv(here::here("Input_data", "cropsuit_p_all_SAR_I_cell5m.csv"))
+
+# Yields
+suit_p_all_H_YIA = read_csv(here::here("Input_data", "cropsuit_p_all_YIA_H_cell5m.csv"))
+suit_p_all_L_YIA = read_csv(here::here("Input_data", "cropsuit_p_all_YIA_L_cell5m.csv"))
+suit_p_all_I_YIA = read_csv(here::here("Input_data", "cropsuit_p_all_YIA_I_cell5m.csv"))
+
+# Spam files for spatial vars
+dat_SPAM <- read_csv(here::here("Input_data", "spam2010V2r0_global_P_TA.csv")) 
+
+# Combined datasets for different intensity levels with yields, areas and production 
+Suit_high = FunProdCalcSuit(suit_p_all_H_SAR,suit_p_all_H_YIA)
+Suit_low = FunProdCalcSuit(suit_p_all_L_SAR,suit_p_all_L_YIA)
+Suit_irr = FunProdCalcSuit(suit_p_all_I_SAR,suit_p_all_I_YIA)
+
+
+# Cell5m - Making it spatial ----------------------------------------------
+# Loading the different cell5m files 
+
+# This file covers the whole globe and can be used to spatially join the suitability tables by cell5m (https://dataverse.harvard.edu/dataset.xhtml?persistentId=doi:10.7910/DVN/MZLXVQ)
+cell5m_global = st_read("Input_data/Mappings/cell5m/Cell5mPoints_Global.shp", 
+                        agr= c(grid_code="identity")) #%>% rename(cell5m = grid_code)
+
+# Cell5m from suitability table (Ulrike) + pulling cell5m into vector
+cell5m_Suit = Suit_high %>% distinct_at(vars(cell5m), .keep_all = FALSE) %>% pull()
+
+## Making the suitability cell5m spatial by filtering  
+Suit_cell5m_dist = cell5m_global %>%  filter(grid_code %in% cell5m_Suit)
+
+# cell5m = st_read("Input_data/cell5m.shp", agr= c(cell5m="identity"))
+zones = st_read("Input_data/Zone_ADM0_SOILagri_AEZ.shp", agr= c(AEZ="constant", Soil_agri="constant", ADM0_CODE="constant"))
+
+join_point = st_join(Suit_cell5m_dist, zones, join = st_intersects) #result in points
+st_write(join_point, "Input_data/join_point.shp")
+
+join_polygon = st_join(zones, Suit_cell5m_dist, join = st_intersects) #left_join; result in multipolygon -> downside: not covering all the points. 
+st_write(join_polygon, "Input_data/join_polygon.shp")
+
+# Create final zonation for the suitability maps
+suit_zone_cell5m = join_polygon %>% 
+  st_set_geometry(NULL) %>% as_tibble() %>% 
+  dplyr::rename("cell5m" = "grid_code")
+
+write_csv(suit_zone_cell5m, here::here("Input_data", "suit_zone_cell5m.csv"))
+
+# bind suitability and zonation 
+suit_H = Suit_high %>% right_join(., suit_zone_cell5m, by="cell5m")
+suit_L = Suit_low %>% right_join(., suit_zone_cell5m, by="cell5m")  
+suit_I = Suit_irr %>% right_join(., suit_zone_cell5m, by="cell5m")
+
+# write_csv(suit_H, here::here("Input_data", "suit_H_zone.csv"))
+# write_csv(suit_L, here::here("Input_data", "suit_L_zone.csv"))
+# write_csv(suit_I, here::here("Input_data", "suit_I_zone.csv"))
+
+# aggregating per zone, crop, tech-type and create one single df
+suit_agg = list(suit_H, suit_L, suit_I) %>% 
+  lapply(Fun_aggregate) %>% 
+  do.call(rbind.data.frame, .) %>% 
+  as_tibble()
+
+# Calculate the yields 
+suit_allvars = suit_agg %>% 
+  mutate(SYI_kgha = round((SPR_ton*1000)/SAR_ha),1) %>% 
+  mutate_at(vars(SYI_kgha), ~replace(., is.nan(.), 0))
+# write_csv(suit_allvars, here::here("Input_data","suit_allvars.csv"))
+
+
+# Mapping -----------------------------------------------------------------
+# CIFOS crop 
+
+
+# country names 
+
+
+
+# EU27 - subsetting -------------------------------------------------------
+
+
+
+