From 5a52f40416116a12d53c3ed776c846ec041c49ef Mon Sep 17 00:00:00 2001
From: s1924442 <b.arendarczyk@sms.ed.ac.uk>
Date: Fri, 12 Feb 2021 18:12:00 +0000
Subject: [PATCH] Further work on forestry and carbon components.

---
 GAMS/IntExtOpt.gms                            | 55 ++++++-------
 src/ac/ed/lurg/ModelMain.java                 |  2 -
 .../country/gams/GamsLocationOptimiser.java   | 77 ++++++++++++++++++-
 .../country/gams/GamsRasterOptimiser.java     | 16 ++--
 src/ac/ed/lurg/landuse/CarbonFluxItem.java    | 10 +--
 src/ac/ed/lurg/landuse/CarbonFluxReader.java  | 42 +++++-----
 src/ac/ed/lurg/landuse/LandUseItem.java       | 30 ++++++++
 src/ac/ed/lurg/landuse/WoodYieldItem.java     | 10 +--
 src/ac/ed/lurg/landuse/WoodYieldReader.java   | 42 +++++-----
 src/ac/ed/lurg/types/ForestType.java          | 16 ++++
 src/ac/ed/lurg/types/GamsLandCoverType.java   | 45 +++++++++++
 src/ac/ed/lurg/types/LandCoverType.java       |  3 +-
 12 files changed, 256 insertions(+), 92 deletions(-)
 create mode 100644 src/ac/ed/lurg/types/ForestType.java
 create mode 100644 src/ac/ed/lurg/types/GamsLandCoverType.java

diff --git a/GAMS/IntExtOpt.gms b/GAMS/IntExtOpt.gms
index 97955b75..1b94ca72 100644
--- a/GAMS/IntExtOpt.gms
+++ b/GAMS/IntExtOpt.gms
@@ -9,18 +9,13 @@
  SET feed_crop_less_pasture(feed_crop)                          / wheat, maize, rice, oilcrops, pulses, starchyRoots, fruitveg /;
  SET import_crop(all_types)   / monogastrics, ruminants,           wheat, maize, rice, oilcrops, pulses, starchyRoots, fruitveg, sugar, energycrops/;
 
- SET land_cover / cropland, pasture, timberForest, carbonForest, natural /;
- SET wood_producing / timberForest, carbonForest, natural /;
+ SET land_cover / cropland, pasture, timberForest, carbonForest, otherNatural /;
+ SET wood_producing / timberForest, carbonForest, otherNatural /;
  ALIAS (land_cover, land_cover_before);
- ALIAS (land_cover, land_cover_start);
- ALIAS (land_cover, land_cover_after);
-
 
  SET location;
  PARAMETER suitableLandArea(location)        areas of land in Mha;
  PARAMETER previousCropArea(crop, location)      areas for previous timestep in Mha;
- PARAMETER previousLandCoverArea(land_cover_before, land_cover_start, location) land cover area at timestep t-1 split by land cover at t-2;
- PARAMETER minimumLandCoverArea(land_cover, location) minimum land cover area to constrain conversion;
  PARAMETER previousFertIntensity(crop, location);
  PARAMETER previousIrrigIntensity(crop, location);
  PARAMETER previousOtherIntensity(crop, location);
@@ -63,8 +58,10 @@
  SCALAR domesticPriceMarkup                  factor price increased from cost of production;
  SCALAR maxLandExpansionRate                 max rate of country land expansion;
 
- PARAMETER carbonFluxRate(land_cover_before, land_cover_start, location)        carbon flux;
- PARAMETER woodYield(land_cover_before, land_cover_start, location)      wood yield;
+ PARAMETER previousLandCoverArea(land_cover, location) land cover area in Mha;
+* PARAMETER minimumLandCoverArea(land_cover, location) minimum land cover area to constrain conversion;
+ PARAMETER carbonFluxRate(land_cover_before, land_cover_after, location) carbon flux;
+ PARAMETER woodYield(land_cover_before, land_cover_after, location)      wood yield;
  SCALAR carbonPrice                                 price of carbon;
  SCALAR woodPrice                                   price of wood and timber;
 
@@ -75,7 +72,8 @@ $load previousCropArea, previousFertIntensity, previousIrrigIntensity, previousO
 $load yieldNone, yieldFertOnly, yieldIrrigOnly, yieldBoth, yieldShock
 $load fertParam, irrigParam, otherIParam, exportPrices, importPrices, maxNetImport, minNetImport, unhandledCropRate, setAsideRate, maxLandExpansionRate, subsidyRate
 $load meatEfficency, otherICost, irrigCost, irrigMaxRate, irrigConstraint, fertiliserUnitCost, domesticPriceMarkup, minDemandPerCereal, minDemandPerOilcrop, seedAndWasteRate
-$load carbonFluxRate, woodYield, carbonPrice, woodPrice
+$load previousLandCoverArea, carbonFluxRate, woodYield, carbonPrice, woodPrice
+*$load minimumLandCoverArea
 $gdxin
 
  SCALAR delta use to smooth power function see 7.5 www.gams.com dd docs solversconopt.pdf / 0.00000000001 /;
@@ -118,8 +116,7 @@ $gdxin
 *lcTransitionCost(lc_type_previous, lc_type, location) = carbonFluxRate(lc_type_previous, lc_type, location) * carbonPrice + agriExpansionCost(location);
 
  VARIABLES
-       cropArea(crop, location)               total area for crops - Mha
-       landCoverArea(land_cover_start, land_cover_after, location)  forest and natural area
+       cropArea(crop, location)             total area for crops - Mha
        fertI(crop, location)              fertilizer intensity for each crop - factor between 0 and 1
        irrigI(crop, location)             irrigation intensity for each crop - factor between 0 and 1
        otherIntensity(crop, location)
@@ -135,6 +132,8 @@ $gdxin
        cropDecrease(location)
        pastureIncrease(location)
        pastureDecrease(location)
+       landCoverArea(land_cover, location)  land cover area in Mha
+       landCoverChange(land_cover_before, land_cover_after, location) land cover change
        woodHarvest(location)              total wood harvested
        carbonFlux(location)               total carbon flux
        totalFeedDM
@@ -171,11 +170,12 @@ $gdxin
        PASTURE_INCREASE_CONV_CALC(location)
        PASTURE_DECREASE_CONV_CALC(location)
        PASTURE_TOTAL_CHANGE_CONSTRAINT(location)
-       PREVIOUS_LAND_COVER_CONSTRAINT(land_cover_start, location)  constraint so that previous land cover at this timestep equals land cover at previous timestep
-       CROPLAND_CONSTRAINT(location) cropland area equals sum of all crop areas
-       PASTURE_CONSTRAINT(location)  pasture area (land cover) equals pasture area (land use)
-       MINIMUM_LAND_COVER_CONSTRAINT(location) constraint on land cover conversion
-       TOTAL_LAND_COVER_CONSTRAINT(location)   total land cover equals suitable area
+       CROPLAND_LAND_COVER_CALC(location) cropland area equals sum of all crop areas
+       PASTURE_LAND_COVER_CALC(location)  pasture area (land cover) equals pasture area (land use)
+*       MINIMUM_LAND_COVER_CONSTRAINT(location) constraint on land cover conversion
+*       TOTAL_LAND_COVER_CONSTRAINT(location)   total land cover equals suitable area
+       LAND_COVER_CHANGE_CALC(land_cover, location)
+       LAND_COVER_CHANGE_CONSTRAINT(land_cover, location)
        WOOD_HARVEST_CALC(location)                          calc wood harvested
        CARBON_FLUX_CALC(location)                           calc carbon flux
        COST_EQ                                        total cost objective function;
@@ -217,8 +217,6 @@ $gdxin
 
  TOTAL_LAND_CHANGE_CONSTRAINT(location) .. suitableLandArea(location) =G= sum(crop_less_pasture, cropArea(crop_less_pasture, location)) / (1.0 - unhandledCropRate) + cropArea('pasture', location);
 
-
-
  MAX_NET_IMPORT_CONSTRAINT(import_crop) .. importAmount(import_crop) - exportAmount(import_crop) =L= maxNetImport(import_crop);
  MIN_NET_IMPORT_CONSTRAINT(import_crop) .. importAmount(import_crop) - exportAmount(import_crop) =G= minNetImport(import_crop);
  PASTURE_IMPORT_CONSTRAINT .. importAmount('pasture') =E= 0;
@@ -237,19 +235,22 @@ $gdxin
  PASTURE_DECREASE_CONV_CALC(location) .. pastureDecrease(location) =G= -(cropArea('pasture', location) - previousCropArea('pasture', location));
  PASTURE_TOTAL_CHANGE_CONSTRAINT(location) .. pastureIncrease(location) -pastureDecrease(location) =G= cropArea('pasture', location) - previousCropArea('pasture', location);
 
- PREVIOUS_LAND_COVER_CONSTRAINT(land_cover_start, location) .. sum(land_cover_after, landCoverArea(land_cover_start, land_cover_after, location)) =E= sum(land_cover_before, previousLandCoverArea(land_cover_before, land_cover_start, location));
+ CROPLAND_LAND_COVER_CAL(location) .. sum(land_cover, landCoverArea('cropland', location)) =E= sum(crop_less_pasture, cropArea(crop, location)) / (1.0 - unhandledCropRate);
+
+ PASTURE_LAND_COVER_CAL(location) .. sum(land_cover, landCoverArea('pasture', location)) =E= cropArea('pasture', location);
 
- CROPLAND_CONSTRAINT(location) .. sum(land_cover_start, landCoverArea(land_cover_start, 'cropland', location)) =E= sum(crop, cropArea(crop, location));
+* MINIMUM_LAND_COVER_CONSTRAINT(location, land_cover) .. landCoverArea(land_cover, location)) =G= minimumLandCover(land_cover_after, location);
 
- PASTURE_CONSTRAINT(location) .. sum(land_cover_start, landCoverArea(land_cover_start, 'pasture', location)) =E= cropArea('pasture', location);
+ LAND_COVER_CHANGE_CALC(land_cover, location) .. landCoverArea(land_cover, location) =E= previousLandCover(land_cover, location) +
+                     sum(land_cover_before, lcChange(land_cover_before, land_cover, location)) - sum(land_cover_after, lcChange(land_cover, land_cover_after, location));
 
- MINIMUM_LAND_COVER_CONSTRAINT(location, land_cover_after) .. sum(land_cover_start, landCoverArea(land_cover_start, land_cover_after, location)) =G= minimumLandCover(land_cover_after, location);
+ LAND_COVER_CHANGE_CONSTRAINT(land_cover, location) .. sum(land_cover_after, landCoverChange(land_cover, land_cover_after, location)) =L= previousLandCover(land_cover, location);
 
- TOTAL_LAND_COVER_CONSTRAINT(location) .. sum((land_cover_start, land_cover_after), landCoverArea(land_cover_start, land_cover_after, location)) =E= suitableArea(location);
+* TOTAL_LAND_COVER_CONSTRAINT(location) .. sum(land_cover, landCoverArea(land_cover, location)) =E= suitableArea(location);
 
- WOOD_HARVEST_CALC(location) .. woodHarvest(location) =E= sum((land_cover_start, land_cover_after), landCoverArea(land_cover_start, land_cover_after, location) * woodYield(land_cover_start, land_cover_after, location));
+ WOOD_HARVEST_CALC(location) .. woodHarvest(location) =E= sum((land_cover_before, land_cover_after, location), landCoverChange(land_cover_before, land_cover_after, location) * woodYield(land_cover_before, land_cover_after, location));
 
- CARBON_FLUX_CALC(location) .. carbonFlux(location) =E= sum((land_cover_start, land_cover_after), landCoverArea(land_cover_start, land_cover_after, location) * carbonFluxRate(land_cover_start, land_cover_after, location));
+ CARBON_FLUX_CALC(location) .. carbonFlux(location) =E= sum((land_cover_before, land_cover_after), landCoverChange(land_cover_before, land_cover_after, location) * carbonFluxRate(land_cover_before, land_cover_after, location));
 
 * CROPLAND_INCREASE_CONSTRAINT(location) .. cropIncrease(location)/(1.0 - unhandledCropRate) =L=  maxLandChangeRate * ( suitableLandArea(location) -sum(crop_less_pasture, previousCropArea(crop_less_pasture, location))/*(1.0 - unhandledCropRate) - previousCropArea('pasture', location) );
 * PASTURE_DECREASE_CONSTRAINT(location) .. pastureDecrease(location) =L=  maxLandChangeRate * previousCropArea('pasture', location);
@@ -277,7 +278,7 @@ $gdxin
  irrigI.L(crop, location) = previousIrrigIntensity(crop, location);
  otherIntensity.L(crop, location) = previousOtherIntensity(crop, location);
  cropArea.L(crop, location) = previousCropArea(crop, location);
- nonCropArea.L(wood_producing, location) = previousNonCropArea(wood_producing, location);
+ landCoverArea.L(land_cover, location) = previousLandCoverArea(land_cover, location);
  ruminantFeed.L(feed_crop) = previousRuminantFeed(feed_crop);
  monogastricFeed.L(feed_crop) = previousMonogastricFeed(feed_crop);
  importAmount.L(all_types) = previousImportAmount(all_types);
diff --git a/src/ac/ed/lurg/ModelMain.java b/src/ac/ed/lurg/ModelMain.java
index 5a678e9f..8525acd2 100644
--- a/src/ac/ed/lurg/ModelMain.java
+++ b/src/ac/ed/lurg/ModelMain.java
@@ -73,8 +73,6 @@ public class ModelMain {
 
 	private InternationalMarket internationalMarket;
 	private IrrigationRasterSet currentIrrigationData;
-	private CarbonFluxRasterSet currentCarbonFluxData;
-	private WoodYieldRasterSet currentWoodYieldData;
 	private RasterSet<LandUseItem> globalLandUseRaster;
 	private RasterSet<IntegerRasterItem> clusterIdRaster;
 
diff --git a/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java b/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java
index 256a88a5..302260cc 100644
--- a/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java
+++ b/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java
@@ -29,6 +29,8 @@ import ac.ed.lurg.landuse.IrrigationItem;
 import ac.ed.lurg.landuse.LandUseItem;
 import ac.ed.lurg.types.CommodityType;
 import ac.ed.lurg.types.CropType;
+import ac.ed.lurg.types.ForestType;
+import ac.ed.lurg.types.GamsLandCoverType;
 import ac.ed.lurg.types.LandCoverType;
 import ac.ed.lurg.types.Parameter;
 import ac.ed.lurg.types.YieldType;
@@ -88,8 +90,6 @@ public class GamsLocationOptimiser {
 
 		if (DEBUG) LogWriter.println("\nPrevious crop and land areas");	
 		GAMSParameter prevCropP = inDB.addParameter("previousCropArea", 2);
-		GAMSParameter prevLandCoverP = inDB.addParameter("previousLandCoverArea", 3);
-		GAMSParameter minLandCoverP = inDB.addParameter("minimumLandCoverArea", 3);
 		GAMSParameter prevFertIP = inDB.addParameter("previousFertIntensity", 2);
 		GAMSParameter prevIrrigIP = inDB.addParameter("previousIrrigIntensity", 2);
 		GAMSParameter prevOtherIP = inDB.addParameter("previousOtherIntensity", 2);
@@ -101,7 +101,10 @@ public class GamsLocationOptimiser {
 	
 		GAMSParameter landP = inDB.addParameter("suitableLandArea", 1);
 		GAMSParameter agriExpansionCostP = inDB.addParameter("agriExpansionCost", 1);		
-		GAMSParameter seedAndWasteRateP = inDB.addParameter("seedAndWasteRate", 1);		
+		GAMSParameter seedAndWasteRateP = inDB.addParameter("seedAndWasteRate", 1);	
+		
+		GAMSParameter prevLandCoverP = inDB.addParameter("previousLandCoverArea", 2);
+		//GAMSParameter minLandCoverP = inDB.addParameter("minimumLandCoverArea", 2);
 		
 		double totalAgriLand = 0;
 		double totalSuitable = 0;
@@ -162,6 +165,37 @@ public class GamsLocationOptimiser {
 				totalAgriLand += area;
 				
 			}
+			
+			/*
+			// Add forest land cover
+			for (ForestType ft : ForestType.values()) {
+				Vector<String> v = new Vector<String>();
+				v.add(ft.getName());
+				v.add(locString);
+				setGamsParamValue(prevLandCoverP.addRecord(v), landUseItem.getForestArea(ft), 2);
+			}
+			
+			// Add non-forest land cover. Need to do it like this because mapping is not 1:1 for natural and not all LandCoverTypes included.
+			Vector<String> v = new Vector<String>();
+			
+			v.add("cropland");
+			v.add(locString);
+			setGamsParamValue(prevLandCoverP.addRecord(v), landUseItem.getLandCoverArea(LandCoverType.CROPLAND), 2);
+
+			v.set(0, "pasture");
+			setGamsParamValue(prevLandCoverP.addRecord(v), landUseItem.getLandCoverArea(LandCoverType.PASTURE), 2);
+			
+			v.set(0, "otherNatural");
+			setGamsParamValue(prevLandCoverP.addRecord(v), landUseItem.getLandCoverArea(LandCoverType.OTHER_NATURAL) + landUseItem.getLandCoverArea(LandCoverType.UNMANAGED_FOREST), 2);
+			*/
+			
+			for (GamsLandCoverType lc : GamsLandCoverType.values()) {
+				Vector<String> v = new Vector<String>();
+				v.add(lc.getName());
+				v.add(locString);
+				setGamsParamValue(prevLandCoverP.addRecord(v), landUseItem.getGamsLandCoverArea(lc), 2);
+			}
+		
 		}
 		if (DEBUG) LogWriter.println(String.format("  Total agricultural %.1f,\t suitable %.1f", totalAgriLand, totalSuitable));
 
@@ -371,10 +405,12 @@ public class GamsLocationOptimiser {
 		GAMSParameter parmCroplandArea = outDB.getParameter("totalCropland");
 		GAMSParameter parmTotalArea = outDB.getParameter("totalArea");
 		GAMSParameter parmProdShock = outDB.getParameter("productionShock");
+		GAMSVariable varLandCoverArea = outDB.getVariable("landCoverArea");
 		
 		double totalCropArea = 0;
 		double totalPastureArea = 0;
 		double area, fertIntensity, irrigIntensity, otherIntensity = Double.NaN, ruminantFeed, monogastricFeed, netImport, netImportCost, yield, unitCost, prod, prodCost;
+		double landCoverArea;
 		
 		final LazyHashMap<Integer, LandUseItem> landUses = new LazyHashMap<Integer, LandUseItem>() { 
 			protected LandUseItem createValue() { return new LandUseItem(); }
@@ -420,7 +456,8 @@ public class GamsLocationOptimiser {
 				IrrigationItem irrigRefData = allIrrigationRefData.get(locId);
 				landUseItem.setIntensity(cropType, new Intensity(fertIntensity, irrigIntensity, otherIntensity, yield, unitCost, irrigRefData.getMaxIrrigAmount(cropType)));
 			}
-
+			
+			
 			double croplandArea = getParmValue(parmCroplandArea, locationName);
 			if (cropType.equals(CropType.PASTURE)) {
 				landUseItem.setLandCoverArea(LandCoverType.PASTURE, area);
@@ -432,7 +469,39 @@ public class GamsLocationOptimiser {
 			}
 			
 			landUseItem.setCropFraction(cropType, croplandArea > 0 ? area/croplandArea : 0);
+			
 		}
+		/*
+		for (GAMSVariableRecord rec : varLandCoverArea) {
+			String lc = rec.getKeys()[0];
+			String locationName = rec.getKeys()[1];
+			int locId = Integer.parseInt(locationName);
+			landCoverArea = rec.getLevel();
+			LandUseItem landUseItem = landUses.lazyGet(locId);
+			
+			switch(lc) {
+			case "forestT":
+				landUseItem.setForestArea(ForestType.TIMBER_FOREST, landCoverArea);
+				break;
+			case "forestC":
+				landUseItem.setForestArea(ForestType.CARBON_FOREST, landCoverArea);
+			}
+			
+		}
+		*/
+		for (GAMSVariableRecord rec : varLandCoverArea) {
+			String lc = rec.getKeys()[0];
+			String locationName = rec.getKeys()[1];
+			int locId = Integer.parseInt(locationName);
+			
+			landCoverArea = rec.getLevel();
+			
+			LandUseItem landUseItem = landUses.lazyGet(locId);
+			
+			landUseItem.setGamsLandCoverArea(GamsLandCoverType.getForName(lc), landCoverArea);
+			
+		}
+		
 
 		for (CropType meatTypes : CropType.getMeatTypes()) {
 			netImport = getParmValue(parmNetImports, meatTypes.getGamsName());
diff --git a/src/ac/ed/lurg/country/gams/GamsRasterOptimiser.java b/src/ac/ed/lurg/country/gams/GamsRasterOptimiser.java
index 292a38af..4511ed93 100644
--- a/src/ac/ed/lurg/country/gams/GamsRasterOptimiser.java
+++ b/src/ac/ed/lurg/country/gams/GamsRasterOptimiser.java
@@ -12,6 +12,7 @@ import ac.ed.lurg.landuse.IrrigationItem;
 import ac.ed.lurg.landuse.LandUseItem;
 import ac.ed.lurg.landuse.WoodYieldItem;
 import ac.ed.lurg.types.CropType;
+import ac.ed.lurg.types.GamsLandCoverType;
 import ac.ed.lurg.types.LandCoverType;
 import ac.ed.lurg.types.YieldType;
 import ac.ed.lurg.utils.LazyTreeMap;
@@ -189,14 +190,18 @@ public class GamsRasterOptimiser {
 			double pastureFromNatural = 0;
 			double cropFromNatural = 0;
 
+			// Gained pasture and lost cropland
 			if (totalPastureChange > 0 && totalCroplandChange < 0) {
 				pastureFromCrop = Math.min(totalPastureChange, -totalCroplandChange);
-
+				
+				// Gain in pasture greater than loss in cropland
 				if (totalPastureChange > -totalCroplandChange)
-					pastureFromNatural = totalPastureChange+totalCroplandChange;
+					pastureFromNatural = totalPastureChange+totalCroplandChange; // take pasture gain not met by cropland loss from natural
+				// Gain in pasture less than loss in cropland
 				else
-					cropFromNatural = totalPastureChange+totalCroplandChange;
+					cropFromNatural = totalPastureChange+totalCroplandChange; // (-ve) add cropland loss not met by pasture gain to natural 
 			}
+			// Gained cropland and lost pasture
 			else if (totalPastureChange < 0 && totalCroplandChange > 0) {
 				pastureFromCrop = -Math.min(-totalPastureChange, totalCroplandChange);
 
@@ -241,6 +246,7 @@ public class GamsRasterOptimiser {
 		return newLandUseRaster;
 	}
 
+
 	private void allocAllLandCropsTop(RasterSet<LandUseItem> newLandUseRaster, Set<RasterKey> keys, LandCoverType toLC, LandCoverType fromLC, double change, Integer locId) {
 		LandCoverType toType = toLC;
 		LandCoverType fromType = fromLC;
@@ -416,8 +422,8 @@ public class GamsRasterOptimiser {
 				}
 				
 				// Aggregate carbon fluxes and wood yields
-				for (LandCoverType lc_previous : LandCoverType.values()) {
-					for (LandCoverType lc_new : LandCoverType.values()) {
+				for (GamsLandCoverType lc_previous : GamsLandCoverType.values()) {
+					for (GamsLandCoverType lc_new : GamsLandCoverType.values()) {
 						aggCFlux.setCarbonFlux(lc_previous, lc_new, aggregateMean(aggCFlux.getCarbonFlux(lc_previous, lc_new), suitableAreaSoFar, 
 								carbonFluxItem.getCarbonFlux(lc_previous, lc_new), suitableAreaThisTime));
 						
diff --git a/src/ac/ed/lurg/landuse/CarbonFluxItem.java b/src/ac/ed/lurg/landuse/CarbonFluxItem.java
index e83ebf77..4aa1ae41 100644
--- a/src/ac/ed/lurg/landuse/CarbonFluxItem.java
+++ b/src/ac/ed/lurg/landuse/CarbonFluxItem.java
@@ -4,18 +4,18 @@ import ac.sac.raster.RasterItem;
 import java.util.HashMap;
 import java.util.Map;
 
-import ac.ed.lurg.types.LandCoverType;
+import ac.ed.lurg.types.GamsLandCoverType;
 
 public class CarbonFluxItem implements RasterItem {
-	Map<LandCoverType, Map<LandCoverType, Double>> carbonFluxes = new HashMap<LandCoverType, Map<LandCoverType, Double>>(); 
+	Map<GamsLandCoverType, Map<GamsLandCoverType, Double>> carbonFluxes = new HashMap<GamsLandCoverType, Map<GamsLandCoverType, Double>>(); 
 	
-	public void setCarbonFlux(LandCoverType previousLandCover, LandCoverType newLandCover, Double carbonFlux) {
-		Map<LandCoverType, Double> cfMap = new HashMap<LandCoverType, Double>();
+	public void setCarbonFlux(GamsLandCoverType previousLandCover, GamsLandCoverType newLandCover, double carbonFlux) {
+		Map<GamsLandCoverType, Double> cfMap = new HashMap<GamsLandCoverType, Double>();
 		cfMap.put(newLandCover, carbonFlux);		
 		carbonFluxes.put(previousLandCover, cfMap);
 	}
 	
-	public Double getCarbonFlux(LandCoverType previousLandCover, LandCoverType newLandCover) {
+	public Double getCarbonFlux(GamsLandCoverType previousLandCover, GamsLandCoverType newLandCover) {
 		return carbonFluxes.get(previousLandCover).get(newLandCover);
 	}
 
diff --git a/src/ac/ed/lurg/landuse/CarbonFluxReader.java b/src/ac/ed/lurg/landuse/CarbonFluxReader.java
index d23a771c..bac2e814 100644
--- a/src/ac/ed/lurg/landuse/CarbonFluxReader.java
+++ b/src/ac/ed/lurg/landuse/CarbonFluxReader.java
@@ -2,7 +2,7 @@ package ac.ed.lurg.landuse;
 
 import java.util.Map;
 
-import ac.ed.lurg.types.LandCoverType;
+import ac.ed.lurg.types.GamsLandCoverType;
 import ac.sac.raster.AbstractTabularRasterReader;
 import ac.sac.raster.RasterKey;
 import ac.sac.raster.RasterSet;
@@ -17,30 +17,30 @@ public class CarbonFluxReader extends AbstractTabularRasterReader<CarbonFluxItem
 	
 	@Override
 	protected void setData(RasterKey key, CarbonFluxItem item, Map<String, Double> rowValues) {
-		item.setCarbonFlux(LandCoverType.CROPLAND, LandCoverType.OTHER_NATURAL, getValueForCol(rowValues, "c2n"));
-		item.setCarbonFlux(LandCoverType.CROPLAND, LandCoverType.TIMBER_FOREST, getValueForCol(rowValues, "c2f"));
-		item.setCarbonFlux(LandCoverType.CROPLAND, LandCoverType.CARBON_FOREST, getValueForCol(rowValues, "c2x"));
-		item.setCarbonFlux(LandCoverType.CROPLAND, LandCoverType.PASTURE,       getValueForCol(rowValues, "c2p"));
+		item.setCarbonFlux(GamsLandCoverType.CROPLAND, GamsLandCoverType.OTHER_NATURAL, getValueForCol(rowValues, "c2n"));
+		item.setCarbonFlux(GamsLandCoverType.CROPLAND, GamsLandCoverType.TIMBER_FOREST, getValueForCol(rowValues, "c2f"));
+		item.setCarbonFlux(GamsLandCoverType.CROPLAND, GamsLandCoverType.CARBON_FOREST, getValueForCol(rowValues, "c2x"));
+		item.setCarbonFlux(GamsLandCoverType.CROPLAND, GamsLandCoverType.PASTURE,       getValueForCol(rowValues, "c2p"));
 		
-		item.setCarbonFlux(LandCoverType.PASTURE,  LandCoverType.OTHER_NATURAL,  getValueForCol(rowValues, "p2n"));
-		item.setCarbonFlux(LandCoverType.PASTURE,  LandCoverType.TIMBER_FOREST,  getValueForCol(rowValues, "p2f"));
-		item.setCarbonFlux(LandCoverType.PASTURE,  LandCoverType.CARBON_FOREST,  getValueForCol(rowValues, "p2x"));
-		item.setCarbonFlux(LandCoverType.PASTURE,  LandCoverType.CROPLAND,       getValueForCol(rowValues, "p2c"));
+		item.setCarbonFlux(GamsLandCoverType.PASTURE,  GamsLandCoverType.OTHER_NATURAL,  getValueForCol(rowValues, "p2n"));
+		item.setCarbonFlux(GamsLandCoverType.PASTURE,  GamsLandCoverType.TIMBER_FOREST,  getValueForCol(rowValues, "p2f"));
+		item.setCarbonFlux(GamsLandCoverType.PASTURE,  GamsLandCoverType.CARBON_FOREST,  getValueForCol(rowValues, "p2x"));
+		item.setCarbonFlux(GamsLandCoverType.PASTURE,  GamsLandCoverType.CROPLAND,       getValueForCol(rowValues, "p2c"));
 		
-		item.setCarbonFlux(LandCoverType.OTHER_NATURAL,  LandCoverType.PASTURE,  	   getValueForCol(rowValues, "n2p"));
-		item.setCarbonFlux(LandCoverType.OTHER_NATURAL,  LandCoverType.TIMBER_FOREST,  getValueForCol(rowValues, "n2f"));
-		item.setCarbonFlux(LandCoverType.OTHER_NATURAL,  LandCoverType.CARBON_FOREST,  getValueForCol(rowValues, "n2x"));
-		item.setCarbonFlux(LandCoverType.OTHER_NATURAL,  LandCoverType.CROPLAND,       getValueForCol(rowValues, "n2p"));
+		item.setCarbonFlux(GamsLandCoverType.OTHER_NATURAL,  GamsLandCoverType.PASTURE,  	   getValueForCol(rowValues, "n2p"));
+		item.setCarbonFlux(GamsLandCoverType.OTHER_NATURAL,  GamsLandCoverType.TIMBER_FOREST,  getValueForCol(rowValues, "n2f"));
+		item.setCarbonFlux(GamsLandCoverType.OTHER_NATURAL,  GamsLandCoverType.CARBON_FOREST,  getValueForCol(rowValues, "n2x"));
+		item.setCarbonFlux(GamsLandCoverType.OTHER_NATURAL,  GamsLandCoverType.CROPLAND,       getValueForCol(rowValues, "n2p"));
 		
-		item.setCarbonFlux(LandCoverType.TIMBER_FOREST,  LandCoverType.PASTURE,  	   getValueForCol(rowValues, "f2p"));
-		item.setCarbonFlux(LandCoverType.TIMBER_FOREST,  LandCoverType.OTHER_NATURAL,  getValueForCol(rowValues, "f2n"));
-		item.setCarbonFlux(LandCoverType.TIMBER_FOREST,  LandCoverType.CARBON_FOREST,  getValueForCol(rowValues, "f2x"));
-		item.setCarbonFlux(LandCoverType.TIMBER_FOREST,  LandCoverType.CROPLAND,       getValueForCol(rowValues, "f2c"));
+		item.setCarbonFlux(GamsLandCoverType.TIMBER_FOREST,  GamsLandCoverType.PASTURE,  	   getValueForCol(rowValues, "f2p"));
+		item.setCarbonFlux(GamsLandCoverType.TIMBER_FOREST,  GamsLandCoverType.OTHER_NATURAL,  getValueForCol(rowValues, "f2n"));
+		item.setCarbonFlux(GamsLandCoverType.TIMBER_FOREST,  GamsLandCoverType.CARBON_FOREST,  getValueForCol(rowValues, "f2x"));
+		item.setCarbonFlux(GamsLandCoverType.TIMBER_FOREST,  GamsLandCoverType.CROPLAND,       getValueForCol(rowValues, "f2c"));
 		
-		item.setCarbonFlux(LandCoverType.CARBON_FOREST,  LandCoverType.PASTURE,  	   getValueForCol(rowValues, "x2p"));
-		item.setCarbonFlux(LandCoverType.CARBON_FOREST,  LandCoverType.OTHER_NATURAL,  getValueForCol(rowValues, "x2n"));
-		item.setCarbonFlux(LandCoverType.CARBON_FOREST,  LandCoverType.TIMBER_FOREST,  getValueForCol(rowValues, "x2f"));
-		item.setCarbonFlux(LandCoverType.CARBON_FOREST,  LandCoverType.CROPLAND,       getValueForCol(rowValues, "x2c"));
+		item.setCarbonFlux(GamsLandCoverType.CARBON_FOREST,  GamsLandCoverType.PASTURE,  	   getValueForCol(rowValues, "x2p"));
+		item.setCarbonFlux(GamsLandCoverType.CARBON_FOREST,  GamsLandCoverType.OTHER_NATURAL,  getValueForCol(rowValues, "x2n"));
+		item.setCarbonFlux(GamsLandCoverType.CARBON_FOREST,  GamsLandCoverType.TIMBER_FOREST,  getValueForCol(rowValues, "x2f"));
+		item.setCarbonFlux(GamsLandCoverType.CARBON_FOREST,  GamsLandCoverType.CROPLAND,       getValueForCol(rowValues, "x2c"));
 		
                                                                                 
 	}
diff --git a/src/ac/ed/lurg/landuse/LandUseItem.java b/src/ac/ed/lurg/landuse/LandUseItem.java
index 8da6b60e..1729969e 100644
--- a/src/ac/ed/lurg/landuse/LandUseItem.java
+++ b/src/ac/ed/lurg/landuse/LandUseItem.java
@@ -8,6 +8,8 @@ import java.util.Map;
 import ac.ed.lurg.ModelConfig;
 import ac.ed.lurg.types.CropToDouble;
 import ac.ed.lurg.types.CropType;
+import ac.ed.lurg.types.ForestType;
+import ac.ed.lurg.types.GamsLandCoverType;
 import ac.ed.lurg.types.LandCoverType;
 import ac.ed.lurg.utils.Interpolator;
 import ac.sac.raster.InterpolatingRasterItem;
@@ -19,6 +21,8 @@ public class LandUseItem implements InterpolatingRasterItem<LandUseItem>, Serial
 	private Map<CropType, Double> cropFractions = new HashMap<CropType, Double>();
 	private Map<LandCoverType, Double> landCoverAreas = new HashMap<LandCoverType, Double>();
 	private Map<LandCoverType, Double> unprotectedAreas = new HashMap<LandCoverType, Double>();
+	private Map<ForestType, Double> forestAreas = new HashMap<ForestType, Double>();
+	private Map<GamsLandCoverType, Double> gamsLandCoverAreas = new HashMap<GamsLandCoverType, Double>();
 	private double protectedArea; //protected area in Mha
 	private double unavailableArea; //area unavailable due to altitude etc 
 	private double suitableArea;
@@ -425,6 +429,32 @@ public class LandUseItem implements InterpolatingRasterItem<LandUseItem>, Serial
 		
 		return total;
 	}
+	
+	public double getForestArea(ForestType forestType) {
+		return forestAreas.get(forestType);
+	}
+	
+	public void setForestArea(ForestType c, double d) {
+		if (Double.isNaN(d) || Double.isInfinite(d))
+			throw new RuntimeException("AreasItem for " + c + " is " + d);
+
+		double landCover = (d < 0.0) ? 0.0 : d;
+		
+		forestAreas.put(c, landCover);
+	}
+	
+	public double getGamsLandCoverArea(GamsLandCoverType landCoverType) {
+		return gamsLandCoverAreas.get(landCoverType);
+	}
+	
+	public void setGamsLandCoverArea(GamsLandCoverType c, double d) {
+		if (Double.isNaN(d) || Double.isInfinite(d))
+			throw new RuntimeException("AreasItem for " + c + " is " + d);
+
+		double landCover = (d < 0.0) ? 0.0 : d;
+		
+		gamsLandCoverAreas.put(c, landCover);
+	}
 
 	@Override
 	public String toString() {
diff --git a/src/ac/ed/lurg/landuse/WoodYieldItem.java b/src/ac/ed/lurg/landuse/WoodYieldItem.java
index ff7c23ae..04eef0a1 100644
--- a/src/ac/ed/lurg/landuse/WoodYieldItem.java
+++ b/src/ac/ed/lurg/landuse/WoodYieldItem.java
@@ -3,19 +3,19 @@ package ac.ed.lurg.landuse;
 import java.util.HashMap;
 import java.util.Map;
 
-import ac.ed.lurg.types.LandCoverType;
+import ac.ed.lurg.types.GamsLandCoverType;
 import ac.sac.raster.RasterItem;
 
 public class WoodYieldItem implements RasterItem {
-	Map<LandCoverType, Map<LandCoverType, Double>> woodYields = new HashMap<LandCoverType, Map<LandCoverType, Double>>(); 
+	Map<GamsLandCoverType, Map<GamsLandCoverType, Double>> woodYields = new HashMap<GamsLandCoverType, Map<GamsLandCoverType, Double>>(); 
 	
-	public void setWoodYield(LandCoverType previousLandCover, LandCoverType newLandCover, Double woodYield) {
-		Map<LandCoverType, Double> cfMap = new HashMap<LandCoverType, Double>();
+	public void setWoodYield(GamsLandCoverType previousLandCover, GamsLandCoverType newLandCover, double woodYield) {
+		Map<GamsLandCoverType, Double> cfMap = new HashMap<GamsLandCoverType, Double>();
 		cfMap.put(newLandCover, woodYield);		
 		woodYields.put(previousLandCover, cfMap);
 	}
 	
-	public Double getWoodYield(LandCoverType previousLandCover, LandCoverType newLandCover) {
+	public Double getWoodYield(GamsLandCoverType previousLandCover, GamsLandCoverType newLandCover) {
 		return woodYields.get(previousLandCover).get(newLandCover);
 	}
 }
diff --git a/src/ac/ed/lurg/landuse/WoodYieldReader.java b/src/ac/ed/lurg/landuse/WoodYieldReader.java
index 8a2522fa..d74880d1 100644
--- a/src/ac/ed/lurg/landuse/WoodYieldReader.java
+++ b/src/ac/ed/lurg/landuse/WoodYieldReader.java
@@ -2,7 +2,7 @@ package ac.ed.lurg.landuse;
 
 import java.util.Map;
 
-import ac.ed.lurg.types.LandCoverType;
+import ac.ed.lurg.types.GamsLandCoverType;
 import ac.sac.raster.AbstractTabularRasterReader;
 import ac.sac.raster.RasterKey;
 import ac.sac.raster.RasterSet;
@@ -16,30 +16,30 @@ public class WoodYieldReader extends AbstractTabularRasterReader<WoodYieldItem>
 	}
 
 	protected void setData(RasterKey key, WoodYieldItem item, Map<String, Double> rowValues) {
-		item.setWoodYield(LandCoverType.CROPLAND, LandCoverType.OTHER_NATURAL, getValueForCol(rowValues, "c2n"));
-		item.setWoodYield(LandCoverType.CROPLAND, LandCoverType.TIMBER_FOREST, getValueForCol(rowValues, "c2f"));
-		item.setWoodYield(LandCoverType.CROPLAND, LandCoverType.CARBON_FOREST, getValueForCol(rowValues, "c2x"));
-		item.setWoodYield(LandCoverType.CROPLAND, LandCoverType.PASTURE,       getValueForCol(rowValues, "c2p"));
+		item.setWoodYield(GamsLandCoverType.CROPLAND, GamsLandCoverType.OTHER_NATURAL, getValueForCol(rowValues, "c2n"));
+		item.setWoodYield(GamsLandCoverType.CROPLAND, GamsLandCoverType.TIMBER_FOREST, getValueForCol(rowValues, "c2f"));
+		item.setWoodYield(GamsLandCoverType.CROPLAND, GamsLandCoverType.CARBON_FOREST, getValueForCol(rowValues, "c2x"));
+		item.setWoodYield(GamsLandCoverType.CROPLAND, GamsLandCoverType.PASTURE,       getValueForCol(rowValues, "c2p"));
 		
-		item.setWoodYield(LandCoverType.PASTURE,  LandCoverType.OTHER_NATURAL,  getValueForCol(rowValues, "p2n"));
-		item.setWoodYield(LandCoverType.PASTURE,  LandCoverType.TIMBER_FOREST,  getValueForCol(rowValues, "p2f"));
-		item.setWoodYield(LandCoverType.PASTURE,  LandCoverType.CARBON_FOREST,  getValueForCol(rowValues, "p2x"));
-		item.setWoodYield(LandCoverType.PASTURE,  LandCoverType.CROPLAND,       getValueForCol(rowValues, "p2c"));
+		item.setWoodYield(GamsLandCoverType.PASTURE,  GamsLandCoverType.OTHER_NATURAL,  getValueForCol(rowValues, "p2n"));
+		item.setWoodYield(GamsLandCoverType.PASTURE,  GamsLandCoverType.TIMBER_FOREST,  getValueForCol(rowValues, "p2f"));
+		item.setWoodYield(GamsLandCoverType.PASTURE,  GamsLandCoverType.CARBON_FOREST,  getValueForCol(rowValues, "p2x"));
+		item.setWoodYield(GamsLandCoverType.PASTURE,  GamsLandCoverType.CROPLAND,       getValueForCol(rowValues, "p2c"));
 		
-		item.setWoodYield(LandCoverType.OTHER_NATURAL,  LandCoverType.PASTURE,  	   getValueForCol(rowValues, "n2p"));
-		item.setWoodYield(LandCoverType.OTHER_NATURAL,  LandCoverType.TIMBER_FOREST,  getValueForCol(rowValues, "n2f"));
-		item.setWoodYield(LandCoverType.OTHER_NATURAL,  LandCoverType.CARBON_FOREST,  getValueForCol(rowValues, "n2x"));
-		item.setWoodYield(LandCoverType.OTHER_NATURAL,  LandCoverType.CROPLAND,       getValueForCol(rowValues, "n2p"));
+		item.setWoodYield(GamsLandCoverType.OTHER_NATURAL,  GamsLandCoverType.PASTURE,  	   getValueForCol(rowValues, "n2p"));
+		item.setWoodYield(GamsLandCoverType.OTHER_NATURAL,  GamsLandCoverType.TIMBER_FOREST,  getValueForCol(rowValues, "n2f"));
+		item.setWoodYield(GamsLandCoverType.OTHER_NATURAL,  GamsLandCoverType.CARBON_FOREST,  getValueForCol(rowValues, "n2x"));
+		item.setWoodYield(GamsLandCoverType.OTHER_NATURAL,  GamsLandCoverType.CROPLAND,       getValueForCol(rowValues, "n2p"));
 		
-		item.setWoodYield(LandCoverType.TIMBER_FOREST,  LandCoverType.PASTURE,  	   getValueForCol(rowValues, "f2p"));
-		item.setWoodYield(LandCoverType.TIMBER_FOREST,  LandCoverType.OTHER_NATURAL,  getValueForCol(rowValues, "f2n"));
-		item.setWoodYield(LandCoverType.TIMBER_FOREST,  LandCoverType.CARBON_FOREST,  getValueForCol(rowValues, "f2x"));
-		item.setWoodYield(LandCoverType.TIMBER_FOREST,  LandCoverType.CROPLAND,       getValueForCol(rowValues, "f2c"));
+		item.setWoodYield(GamsLandCoverType.TIMBER_FOREST,  GamsLandCoverType.PASTURE,  	   getValueForCol(rowValues, "f2p"));
+		item.setWoodYield(GamsLandCoverType.TIMBER_FOREST,  GamsLandCoverType.OTHER_NATURAL,  getValueForCol(rowValues, "f2n"));
+		item.setWoodYield(GamsLandCoverType.TIMBER_FOREST,  GamsLandCoverType.CARBON_FOREST,  getValueForCol(rowValues, "f2x"));
+		item.setWoodYield(GamsLandCoverType.TIMBER_FOREST,  GamsLandCoverType.CROPLAND,       getValueForCol(rowValues, "f2c"));
 		
-		item.setWoodYield(LandCoverType.CARBON_FOREST,  LandCoverType.PASTURE,  	   getValueForCol(rowValues, "x2p"));
-		item.setWoodYield(LandCoverType.CARBON_FOREST,  LandCoverType.OTHER_NATURAL,  getValueForCol(rowValues, "x2n"));
-		item.setWoodYield(LandCoverType.CARBON_FOREST,  LandCoverType.TIMBER_FOREST,  getValueForCol(rowValues, "x2f"));
-		item.setWoodYield(LandCoverType.CARBON_FOREST,  LandCoverType.CROPLAND,       getValueForCol(rowValues, "x2c"));
+		item.setWoodYield(GamsLandCoverType.CARBON_FOREST,  GamsLandCoverType.PASTURE,  	   getValueForCol(rowValues, "x2p"));
+		item.setWoodYield(GamsLandCoverType.CARBON_FOREST,  GamsLandCoverType.OTHER_NATURAL,  getValueForCol(rowValues, "x2n"));
+		item.setWoodYield(GamsLandCoverType.CARBON_FOREST,  GamsLandCoverType.TIMBER_FOREST,  getValueForCol(rowValues, "x2f"));
+		item.setWoodYield(GamsLandCoverType.CARBON_FOREST,  GamsLandCoverType.CROPLAND,       getValueForCol(rowValues, "x2c"));
 		                                                                          
 	}
 }
diff --git a/src/ac/ed/lurg/types/ForestType.java b/src/ac/ed/lurg/types/ForestType.java
new file mode 100644
index 00000000..e591ca24
--- /dev/null
+++ b/src/ac/ed/lurg/types/ForestType.java
@@ -0,0 +1,16 @@
+package ac.ed.lurg.types;
+
+public enum ForestType {
+	TIMBER_FOREST("timberForest"),
+	CARBON_FOREST("carbonForest");
+	
+	private String name;
+	
+	ForestType(String name) {
+		this.name = name;
+	}
+	
+	public String getName() {
+		return name;
+	}
+}
diff --git a/src/ac/ed/lurg/types/GamsLandCoverType.java b/src/ac/ed/lurg/types/GamsLandCoverType.java
new file mode 100644
index 00000000..a986955b
--- /dev/null
+++ b/src/ac/ed/lurg/types/GamsLandCoverType.java
@@ -0,0 +1,45 @@
+package ac.ed.lurg.types;
+
+import java.util.HashMap;
+import java.util.Map;
+
+import ac.ed.lurg.utils.LogWriter;
+
+public enum GamsLandCoverType {
+	TIMBER_FOREST("timberForest"),
+	CARBON_FOREST("carbonForest"),
+	OTHER_NATURAL("otherNatural"),
+	PASTURE("pasture"),
+	CROPLAND("cropland");
+	
+	private String name;
+	
+	GamsLandCoverType(String name) {
+		this.name = name;
+	}
+	
+	public String getName() {
+		return name;
+	}
+	
+	private static final Map<String, GamsLandCoverType> nameCache = new HashMap<String, GamsLandCoverType>();
+    static {
+        for (GamsLandCoverType c : values()) {
+        	nameCache.put(c.getName(), c);
+       }
+    }
+    
+	private static GamsLandCoverType getFromCache(Map<String, GamsLandCoverType> cache, String name) {
+		GamsLandCoverType type = cache.get(name);
+		
+		if (type == null)
+			LogWriter.printlnError("Can't find Item for " + name);
+		
+		return type;
+	}
+    
+	public static GamsLandCoverType getForName(String name) {
+		return getFromCache(nameCache, name);
+	}
+}
+
diff --git a/src/ac/ed/lurg/types/LandCoverType.java b/src/ac/ed/lurg/types/LandCoverType.java
index 1bbc1c6f..6c0c8428 100644
--- a/src/ac/ed/lurg/types/LandCoverType.java
+++ b/src/ac/ed/lurg/types/LandCoverType.java
@@ -2,8 +2,7 @@ package ac.ed.lurg.types;
 
 public enum LandCoverType {
 	
-	TIMBER_FOREST ("timberForest", true),
-	CARBON_FOREST ("carbonForest", true),
+	MANAGED_FOREST ("managedForest", true),
 	UNMANAGED_FOREST ("unmanagedForest", true),
 	OTHER_NATURAL ("otherNatural", true),
 	CROPLAND ("cropland", false),
-- 
GitLab