diff --git a/CommonFunctions.R b/CommonFunctions.R
index 29e2a336dd0adecbf419b8fb6594f2a30f2db62f..d9e204fc973c05ec3f26740132f797f47c5b0fdb 100644
--- a/CommonFunctions.R
+++ b/CommonFunctions.R
@@ -36,19 +36,26 @@ itemGroupMapping = data.table(
bioenergyUse = c(rep(FALSE, 8), rep(TRUE, 10))
)
-rasterFromFile = function(lcType, fileName, offset=0.0) {
- plumDt=data.table(read.table(fileName, header=TRUE))
- plumDt[, c('Lon', 'Lat') := list(Lon-offset, Lat-offset)]
- rasterFromXYZ(plumDt[, list(Lon, Lat, get(lcType))], crs = "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
+getItemDetails = function() {
+ itemDetails = fread(paste0(data_dir_root, "/HALF/itemsInputs.csv"))
+ itemDetails[, allowedForVeggie:=TRUE]
+ itemDetails[animal == TRUE & !Item %in% c("Eggs", "Milk - Excluding Butter"), allowedForVeggie:=FALSE]
+ itemDetails
}
-
getCountryCodes = function(includeWorld=TRUE) {
dt = fread(paste0(data_dir_root, "/SSP/country_codes2.csv"), header=TRUE)[, list(Country, CountryCode, Region, IncomeGroup)]
if (includeWorld == TRUE) dt = rbind(dt, data.table(Country = "World", CountryCode = "WLD", Region = "World", IncomeGroup="variable"))
dt
}
+##### Plot helper functions
+rasterFromFile = function(lcType, fileName, offset=0.0) {
+ plumDt=data.table(read.table(fileName, header=TRUE))
+ plumDt[, c('Lon', 'Lat') := list(Lon-offset, Lat-offset)]
+ rasterFromXYZ(plumDt[, list(Lon, Lat, get(lcType))], crs = "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
+}
+
########## Country data
readCountryData = function() {
########## Land
diff --git a/PLUM/Cluster/ProcessClusterRes.R b/PLUM/Cluster/ProcessClusterRes.R
index 8c179da1c9b43995bfa0e1eb64366735dbb8cef0..3c5ae25544a1c154ef26f97ec6f9a4f4b7fce9ba 100644
--- a/PLUM/Cluster/ProcessClusterRes.R
+++ b/PLUM/Cluster/ProcessClusterRes.R
@@ -198,7 +198,7 @@ plotConeLcHind(filename="~/Documents/LURG/PLUM/offline_paper/hind1970/hc/lc_conc
-plotConeRcp = function (ensemble, outBase=baseOutputDir, baseYr=2010, endYr=2100, outfilename, saveToFile=FALSE) {
+plotConeRcpOrSsp = function (ensemble, outBase=baseOutputDir, baseYr=2010, endYr=2100, outfilename, saveToFile=FALSE) {
lc=fread(file.path(outBase, "lc_concat.txt"))
if (ensemble %like% 'rcp') {
@@ -311,6 +311,6 @@ plotConeRcp = function (ensemble, outBase=baseOutputDir, baseYr=2010, endYr=2100
layout_matrix = rbind(c(1,2, 3,4), c(5,6,7,8), c(9, 10,11,11)))
if (saveToFile) dev.off()
}
-#plotConeRcp('rcp', outBase="~/Documents/LURG/PLUM/offline_paper/rcp", outfilename="rcp_comp_over_time.pdf", saveToFile=TRUE)
-plotConeRcp('ssp11', outBase="~/Documents/LURG/PLUM/SamsCouplingPaper", outfilename="ssp.pdf", saveToFile=TRUE)
+#plotConeRcpOrSsp('rcp', outBase="~/Documents/LURG/PLUM/offline_paper/rcp", outfilename="rcp_comp_over_time.pdf", saveToFile=TRUE)
+plotConeRcpOrSsp('ssp11', outBase="~/Documents/LURG/PLUM/SamsCouplingPaper", outfilename="ssp.pdf", saveToFile=TRUE)
diff --git a/PLUM/SetupData/PreparePLUMData.R b/PLUM/SetupData/PreparePLUMData.R
index 9d67045e1e1ccf1a2049de1b8a51994779e43ce2..0137a0aa49413b9ec5f2260a933954287b2b6454 100644
--- a/PLUM/SetupData/PreparePLUMData.R
+++ b/PLUM/SetupData/PreparePLUMData.R
@@ -37,7 +37,7 @@ writePlumFiles = function(dt, year, output_dir) {
popForG = processedFull$histDt[Year == 2010 & Item == "Cereals" & !is.na(population), list(Country, popForGrouping=population)]
newBaseConsumption = merge(newBaseConsumption, popForG, by='Country', all.x=TRUE)
newBaseConsumption[is.na(popForGrouping), popForGrouping:=basePop]
- write.table(newBaseConsumption[Country != "Sudan (former)"], file.path(output_dir,"base_consump.csv"), sep=",", row.names=FALSE, quote=FALSE)
+ write.table(newBaseConsumption[!Country %in% c("World", "Sudan (former)")], file.path(output_dir,"base_consump.csv"), sep=",", row.names=FALSE, quote=FALSE)
# Net imports and production (of crops)
consumpDt = dt$com_bal[Item %in% itemGroupMapping$combalItem & Year == year, list(Country, Year, Item, net_import=import-export, prod)]
diff --git a/PLUM/SetupData/ReadFaoForPlum-notHindcast.R b/PLUM/SetupData/ReadFaoForPlum-notHindcast.R
index decdaad4c0273f3ee5c4990becf6a3f5d1a3213c..120e572af0aee52c64b80a52fbc6bd4fee77f852 100644
--- a/PLUM/SetupData/ReadFaoForPlum-notHindcast.R
+++ b/PLUM/SetupData/ReadFaoForPlum-notHindcast.R
@@ -5,6 +5,9 @@ require(zoo)
data_dir_root = "~/Documents/LURG/Data"
fao_data_dir = file.path(data_dir_root, "FAOStat-Dec2015")
+getItemGroupMapping = function() {
+ itemGroupMapping
+}
readLivestockNumbers = function() {
animal_prod_data = fread(file.path(fao_data_dir,"Production_LivestockPrimary_E_All_Data_Norm.csv"))
@@ -84,7 +87,7 @@ disaggregateCountries = function(Dt,wbCountryDt, colsToChange){
getHistoricDataset = function(commodityDt, countryDt) {
- h_dt = commodityDt[Item %in% itemGroupMapping$combalItem,
+ h_dt = commodityDt[Item %in% getItemGroupMapping()$combalItem,
list(Country, Year, Item,
supply=supply,
net_import=(import-export),
@@ -92,7 +95,7 @@ getHistoricDataset = function(commodityDt, countryDt) {
feed_rate=feed/supply,
other_rate=other/supply)]
- h_dt = merge(h_dt, itemGroupMapping[, list(plumDemandItem, Item=combalItem, fcr)], by=c('Item'))
+ h_dt = merge(h_dt, getItemGroupMapping()[, list(plumDemandItem, Item=combalItem, fcr)], by=c('Item'))
itemDetails = fread(paste0(data_dir_root, "/HALF/itemsInputs-PLUM.csv"))
h_dt = merge(h_dt, itemDetails, by='Item', all.x=TRUE)
@@ -201,7 +204,7 @@ processData = function(minYear=1960, year=2010){
"Soyabean Cake", "Sunflowerseed Oil", "Sunflowerseed Cake", "Palmkernel Oil"))
# losses and seed by PLUM crop
- print (merge(com_bal_adj, itemGroupMapping[, list(plumCropItem, Item=combalItem, fcr)], by=c('Item'))[Country == "World" & Year == 2010,
+ print (merge(com_bal_adj, getItemGroupMapping()[, list(plumCropItem, Item=combalItem, fcr)], by=c('Item'))[Country == "World" & Year == year,
list(loss_seed_rate=sum(seed+waste)/sum(supply)*100), by= plumCropItem])
country_data = readCountryData()
@@ -212,7 +215,8 @@ processData = function(minYear=1960, year=2010){
projectedConsumption = getProjectAllConsumption(baseConsumption,ssp[Year >= year],com_curves)
livestockNumbers = readLivestockNumbers()
- return(list(ssp=ssp,histDt=histDt,baseConsumption=baseConsumption,com_curves=com_curves,projectedConsumption=projectedConsumption,com_bal=com_bal_adj, livestockNumbers=livestockNumbers))
+ return(list(ssp=ssp,histDt=histDt, baseConsumption=baseConsumption, com_curves=com_curves, country_data=country_data,
+ projectedConsumption=projectedConsumption, com_bal=com_bal_adj, livestockNumbers=livestockNumbers))
}
processedFull = processData()
diff --git a/PLUM/SetupData/UclProjections.R b/PLUM/SetupData/UclProjections.R
new file mode 100644
index 0000000000000000000000000000000000000000..d49ffaa249fc12510501cceb6f7b2e801bb98f5b
--- /dev/null
+++ b/PLUM/SetupData/UclProjections.R
@@ -0,0 +1,65 @@
+fao_data_dir = file.path(data_dir_root, "FAOStat-Sep2018")
+
+getItemGroupMapping = function() {
+ rbind(itemGroupMapping,
+ data.table(plumDemandItem=c("Friut & Veg"), plumCropItem=c("Friut & Veg"), cropProdItem=c("Vegetables Primary", "Fruit Primary"),
+ combalItem=c("Vegetables", "Fruits - Excluding Wine")), fill=TRUE)
+}
+
+processedFull = processData(year=2013)
+itemDetails = getItemDetails()
+
+consumpChange = merge(processedFull$projectedConsumption[, list(Country, Year, plumCropItem=Item, Model, Scenario, population, cpcAdj=projectedCpc/baseCpc)],
+ itemDetails,
+ by=c("plumCropItem"), allow.cartesian=TRUE)
+
+baseConsump = merge(processedFull$com_bal[Item %in% itemDetails$Item & Year == 2013],
+ processedFull$country_data[, list(Country, Year, population)],
+ by=c('Country', 'Year'))
+baseConsump = baseConsump[, list(Country, Item, baseCpc=food/population, baseFood=food)]
+
+dt = merge(consumpChange, baseConsump, by=c("Country", "Item"))[, list(Country, Year, Item, SspModel=Model, Scenario, population, cpcAdj, baseCpc, baseFood, cpc=baseCpc*cpcAdj, food=baseCpc*cpcAdj*population)]
+
+write.table(dt, "~/Download/foodItemProjections.csv", row.names=FALSE, sep=",")
+
+
+# FIGURE PLOT
+baseFigDir="~/Downloads"
+
+compDt=processedFull$histDt[Year <= 2013 & !is.na(cpc) & !is.na(population) & Country %in% getCountryCodes()$Country, list(Country, Year, Item, FAO=cpc*population*avg_fcr_dm_adj)]
+
+proj=processedFull$projectedConsumption
+proj=proj[!is.na(projectedCpc) & !is.na(population) & Country %in% getCountryCodes()$Country, list(Country, Year, Item, Projected=projectedCpc*population*avg_fcr_dm_adj, Model, Scenario)]
+projDt=merge(compDt[Year>=2010], proj, by=c('Country', 'Year', 'Item'), all.y=TRUE)
+projDt[!is.na(FAO) & is.na(Projected), Projected:=FAO]
+
+projDt=rbind(projDt[, list(Country, Year, Item, Model, Scenario, Projected)],
+ projDt[, list(Country="World", Projected=sum(Projected)), by=c('Year', 'Item', 'Model', 'Scenario')])
+compMelt=melt(compDt, id.vars=c('Country', 'Year', 'Item'), variable.name='Type')
+allDtToPlot = rbind(compMelt, projDt[, list(Country, Year, Item,
+ Type=paste(gsub('(SSP[0-9])_.*', '\\1', Scenario), gsub('([A-Z]+) .*', '\\1', Model)), value=Projected)])
+
+check=compDt[Country == "World" & Year == 2010]
+check
+
+plotDemandProjections = function(dt, fileName, countries=c("World", "United States of America", "India", "Brazil", "China")) {
+ dt2 = dt[Country %in% countries]
+ faoLabel = "Historic (FAO)"
+ dt2[Type %like% 'OECD', Type:=gsub('(.*) OECD', '\\1 ', Type)]
+ dt2[Type=="FAO", Type:=faoLabel]
+ dt2[, Type:=factor(Type, levels=c(faoLabel, sort(unique(as.character(dt2[!Type %in% c(faoLabel)]$Type)))))]
+ ggplot(dt2, aes(x=Year, y=value, colour=Type, size=Type, linetype=Type, order = Type)) +
+ geom_line() + facet_grid(Country ~ Item, scale="free") + scale_y_log10() +
+ # scale_colour_manual(values = c("black", "red", "#FFB72BCC","#FFB72BCC", "#96DA00CC","#96DA00CC", "#00EABECC", "#00EABECC", "#00D9FFCC", "#00D9FFCC", "#F9A3FFCC", "#F9A3FFCC"), guide = guide_legend(nrow = 2)) + # hcl(h=seq(50,290,60), l=c(80), c=100, alpha=0.8)
+ scale_colour_manual(values = c("black", "red", "#FFB72BCC","#96DA00CC","#00EABECC", "#00D9FFCC", "#F9A3FFCC")) + # hcl(h=seq(50,290,60), l=c(80), c=100, alpha=0.8)
+ scale_size_manual(values = c(0.55, 0.55, rep(.45,10))) +
+ # scale_linetype_manual(values = c(rep("solid", 2), rep(c("solid", "dashed"),5))) + with both OECD and IIASA projections
+ scale_linetype_manual(values = rep("solid", 7)) +
+ theme_minimal() +
+ theme(panel.border = element_rect(colour = "grey10", fill=NA, size=0.3), #axis.line = element_line(size = 0.3, colour ="grey10"),
+ legend.position="bottom", legend.title=element_blank()) +
+ labs(x="Year", y="Total consumption (Mt dry matter/year)") +
+ ggsave(file=fileName,width=9,height=9)
+}
+
+plotDemandProjections(allDtToPlot[!Type %like% 'IIASA'], file.path(baseFigDir,"demand_countries_1961_2100.pdf"))
diff --git a/altLU/AltLURead.R b/altLU/AltLURead.R
index 7ae8be8a5ba73e60d8e90dc76a6164864d447f55..194a59b947eb096e9fd80a60dc34a149e11c5eb0 100644
--- a/altLU/AltLURead.R
+++ b/altLU/AltLURead.R
@@ -17,13 +17,6 @@ getFeedConvRatio = function(generateRandom=FALSE, lowerB=-0.2, upperB=0.2) {
fcrDt
}
-getItemDetails = function() {
- itemDetails = fread(paste0(data_dir_root, "/HALF/itemsInputs.csv"))
- itemDetails[, allowedForVeggie:=TRUE]
- itemDetails[animal == TRUE & !Item %in% c("Eggs", "Milk - Excluding Butter"), allowedForVeggie:=FALSE]
- itemDetails
-}
-
########## Food supply
readFoodSupply = function() {
processFoodSupply = function(fao_file) {
@@ -306,7 +299,3 @@ con_prod = allocateMeatExports(con_prod)
write.table(merge(itemDetails, con_prod[Country == "World" & Year == 2011, list(Item, mj_per_kg, protein_g_per_kg)], by='Item'),
"/Users/peteralexander/Documents/LURG/Alt\ Scenarios/figures/items.csv", row.names=FALSE, sep=",")
-
-
-write.table(merge(itemDetails, con_prod[Country == "United Kingdom" & Year == 2011], by='Item'),
- "/Users/peteralexander/Documents/LURG/Alt\ Scenarios/figures/UK_con_prod.csv", row.names=FALSE, sep=",")